Gaussian processes
A practical guide to function priors, kernel assumptions, predictive uncertainty, and the cubic-cost trade-off behind Gaussian processes.
What you'll learn
- Explain a Gaussian process as a distribution over functions rather than a single fitted curve.
- Choose RBF, Matérn, and periodic kernels deliberately, including the effect of length-scale.
- Compute and interpret predictive mean, latent variance, and noisy-observation variance.
- Recognize cubic scaling and choose inducing points or a different model when it matters.
- Connect GP uncertainty to Bayesian optimisation and acquisition functions.
Before you start
Gaussian processes
At 3:07 a.m., an engineer has twelve expensive battery tests on a screen.
The tests cover temperatures from 20°C to 50°C. A model predicts that battery life at 42°C is 8.1 hours, with a narrow interval. That is useful. Then somebody asks about 85°C, where there are no measurements. The model still returns 8.1 hours, decorated with several decimal places.
A line gives you a prediction. A tree gives you a prediction. Neither automatically tells you, “This part is supported by data; that part is a guess.”
That distinction matters whenever a wrong prediction is expensive: in laboratory experiments, robot control, geospatial mapping, and hyperparameter tuning, where one new evaluation may cost 20 minutes or $500.
A Gaussian process, usually shortened to GP, treats the unknown regression function itself as random. It begins with many plausible curves, then uses observations to rule most of them out.
The result is a predictive distribution: a mean, a variance, and, underneath, correlations between predictions.
A distribution over functions
Ordinary regression asks for one best function. A GP asks:
Before seeing the data, which functions seem plausible, and how should the data change those beliefs?
A GP is written as f ~ GP(m, k). m is the mean function and k is the covariance function, or kernel. Formally, every finite collection of function values has a joint Gaussian distribution.
Imagine input locations x = 0, x = 1, and x = 2. The GP gives a three-dimensional Gaussian distribution for f(0), f(1), and f(2). Draw one triple and you get three points. Repeat this over a dense grid, connecting neighboring points, and you get one possible curve. Draw again and you get another.
The kernel controls whether those curves look like smooth hills, rough paths, or repeating waves. The input need not be time: it can be temperature, latitude, a robot’s joint angle, or a vector of hyperparameters. “Process” means a collection of random variables indexed by inputs.
With no observations, the GP is a prior over curves. After observing data, Bayes’ rule produces a posterior. Curves near reliable observations become more plausible; curves that disagree lose probability.
The kernel is the encoded assumption
The kernel repeatedly answers:
If I know the function here, how much should that tell me about the function there?
A large covariance means two function values should move together under the prior. A covariance near zero means observing one tells you little about the other. The kernel is therefore the model’s assumption about geometry.
The RBF kernel
The radial basis function (RBF) kernel is the common smooth default:
k(x, x') = sigma_f^2 exp(-(x - x')^2 / (2 l^2))
x and x' are inputs, sigma_f^2 sets the vertical scale of variation, and l is the length-scale. A shorter length-scale makes correlation decay faster, so the function can change sharply between neighborhoods. A longer one keeps distant points related and forces slower variation.
Length-scale is not a hard radius. At two length-scales, the model has not forgotten the other point completely; it considers the relationship weak.
RBF sample functions are extremely smooth. That suits a temperature response that changes gradually, but not a signal with corners, jumps, or rough physical variation. A model can draw a beautiful curve through data containing a regime change. Beauty is not evidence.
Matérn and periodic kernels
A Matérn kernel adds a smoothness parameter, usually nu. Smaller nu allows rougher sample functions; larger nu makes them smoother, approaching RBF behavior. Values such as nu = 1.5 and nu = 2.5 are useful intermediate choices.
A periodic kernel says that inputs separated by a period behave alike:
k(x, x') = sigma_f^2 exp(-2 sin^2(pi |x - x'| / p) / l^2)
Here, p is the period. Use it for time of day, day of week, or a genuinely repeating physical cycle. A pure periodic kernel repeats forever; multiplying it by an RBF creates a quasi-periodic model whose cycles gradually lose their relationship.
Kernels can also be added and multiplied. A sum can represent a smooth trend plus a periodic effect. A product can make one component modulate another. For example:
ConstantKernel * RBF + WhiteKernel
represents a smooth latent signal plus independent observation noise.
Hyperparameters are commonly learned by maximizing the log marginal likelihood. Length-scale, signal variance, and noise variance can trade off, so multiple optimizer restarts are often worthwhile.
Predictive mean and variance
Assume noisy observations:
y_i = f(x_i) + epsilon_i
The noise terms are independent Gaussian variables with variance sigma_n^2. Let K contain k(x_i, x_j) and define:
A = K + sigma_n^2 I
The noise term says that two measurements at the same input need not be identical.
For a new input x_*, let k_* be the vector of covariances between it and every training input. Under a Gaussian likelihood, the latent-function posterior is Gaussian:
mu_* = m_* + k_*^T A^-1 (y - m)
sigma_latent^2 = k(x_*, x_*) - k_*^T A^-1 k_*
The data influence the new point through k_*. Strong covariance transfers information; weak covariance does not. The subtracted term is uncertainty removed by conditioning. Near observations it is large, while far away it approaches zero for decaying kernels, so variance returns toward the prior variance.
For a noisy future observation, add measurement noise:
sigma_observation^2 = sigma_latent^2 + sigma_n^2
Latent uncertainty concerns the underlying response. Observation uncertainty concerns the response plus instrument or process noise. At an observed input, a noiseless GP can drive latent variance to zero; a noisy GP should retain some uncertainty.
A numerical example
Take one observation: x = 0, y = 2.
Use a zero prior mean, unit signal variance, RBF length-scale l = 1, and noise variance sigma_n^2 = 0.25. With one training point:
K = [1]
and
A = [1.25]
At x_* = 0.5, the covariance is:
k_* = exp(-0.5^2 / 2) = exp(-0.125), about 0.882.
The posterior mean is:
0 + 0.882 / 1.25 * 2 = 1.412
The latent variance is:
1 - 0.882^2 / 1.25 = 0.377
The latent standard deviation is about 0.614.
At x_* = 2, the covariance is only:
exp(-2^2 / 2) = exp(-2), about 0.135.
The posterior mean becomes about 0.217, and the latent variance becomes about 0.985.
| Test input | Covariance to data | Posterior mean | Latent variance | Latent standard deviation |
|---|---|---|---|---|
0.5 | 0.882 | 1.412 | 0.377 | 0.614 |
2.0 | 0.135 | 0.217 | 0.985 | 0.993 |
Near the observation, the model follows the evidence. Two length-scales away, it mostly returns to the zero prior and admits that it knows little.
A 95 percent latent interval at x = 0.5 is approximately 1.412 +/- 1.203, or [0.209, 2.615]. At x = 2, it is approximately 0.217 +/- 1.946, or [-1.729, 2.163].
For a future noisy observation, the standard deviations would be about 0.792 and 1.111, because noise variance 0.25 is added. The widening follows directly from conditioning a joint Gaussian distribution; it is not an added uncertainty heuristic.
Making a GP behave in production
Start with input geometry. Scale continuous features using statistics learned from the training data, and fit that preprocessing inside each training fold. Choose a kernel that describes a plausible function, add a noise term unless measurements are exact, and inspect learned length-scales and noise levels. A parameter at the edge of its allowed range deserves investigation.
Evaluate distributions as distributions. RMSE measures the center; log predictive density rewards a distribution that puts probability in the right places. For observed held-out targets, use the noisy predictive distribution, with sigma_observation^2 = sigma_latent^2 + sigma_n^2. Evaluate latent intervals only against known latent targets. Coverage estimates also have sampling uncertainty, especially when nearby test cases are correlated.
A GP provides a joint posterior over test points, not just separate intervals. Nearby predictions can be highly correlated; that covariance matters when selecting a batch of experiments or constructing simultaneous bands.
The cubic wall
An exact GP stores an n by n covariance matrix. Training solves a dense linear system, normally with Cholesky factorization.
- Time:
O(n^3) - Memory:
O(n^2)
At n = 10,000, the matrix has 100 million entries, about 800 MB in 64-bit floating point before factorization and work arrays. Hyperparameter optimization may repeat the factorization. This is why a GP that is instant on 500 rows can become painful on 5,000.
Inducing points and other escapes
An inducing-point GP introduces m pseudo-inputs, with m much smaller than n. Many variational methods have training cost roughly:
O(n m^2 + m^3)
instead of O(n^3). The inducing points can be learned or selected to cover the input space. Too few or poorly placed points blur local structure, so compare predictive accuracy and coverage against an exact GP on a smaller representative subset.
Random Fourier features approximate stationary kernels with finite features, turning the problem into a Bayesian or regularized linear model. Nyström methods use representative points to approximate the covariance. Both scale further, but they change the model.
Bayesian optimisation: uncertainty becomes an action
Suppose a model-training run takes 30 minutes. After eight learning-rate trials, the next choice should consider both promising regions and regions where the model is uncertain. That is Bayesian optimisation.
A GP acts as a surrogate for the expensive objective. After each evaluation, update the GP, compute its posterior mean and latent standard deviation, score candidates with an acquisition function, evaluate the best candidate, and repeat.
For maximisation, an upper confidence bound is:
UCB(x) = mu(x) + sqrt(beta) sigma_latent(x)
The mean rewards exploitation; latent standard deviation rewards exploration. Expected improvement asks:
EI(x) = E[max(f(x) - f_best, 0)]
These acquisitions concern the latent objective, not measurement noise. With a WhiteKernel, a library’s predictive standard deviation may include observation noise. Feeding it directly into EI or UCB can make the optimiser chase noisy measurements. Use latent variance separately or a noise-aware acquisition. With noisy evaluations, do not automatically treat the largest observed value as the best latent value; it may be a lucky measurement.
The connection is strongest when the search space is small or moderate, continuous, and meaningfully smooth. A GP over a learning rate on its raw scale is usually inferior to one over log10(learning_rate). One-hot categories may create meaningless distances.
When a GP is the wrong tool
Exact inference is unsuitable for large n. In high-dimensional inputs, distance-based kernels may find that nearly every point is far from every other point. ARD can identify weak dimensions, but adds hyperparameters and can be unstable with limited data.
Stationary RBF and Matérn kernels assume the same relationship throughout the input space. Thresholds or regime changes may require non-stationary kernels, local GPs, input warps, or a different model. Extrapolation is another trap: a zero-mean decaying-kernel GP usually reverts toward zero outside the observed range while its variance grows. Encode a defensible trend, collect data in the new regime, or refuse to extrapolate.
A GP is most attractive when data are limited, inputs have meaningful geometry, and uncertainty changes what you do next. Otherwise, boosted trees, linear models, splines, state-space models, quantile regression, or conformal prediction may be simpler and more reliable.
What to remember
- A GP is a probability distribution over functions; observations condition it into a posterior.
- The kernel is the model’s geometry: RBF is very smooth, Matérn permits roughness, and periodic means repetition.
- Length-scale controls how quickly information travels through input space.
- Variance grows away from data because covariance with observations fades.
- Exact inference costs
O(n^3)time andO(n^2)memory.
Quick check
Practice this in an interview
All questionsRidge regression is maximum a posteriori estimation for a linear model with Gaussian observation noise and a zero-mean Gaussian prior on the coefficients. The regularization strength is the noise variance divided by the prior variance, subject to the scaling convention used in the Ridge objective.
A GMM fits a weighted set of Gaussian distributions and gives each observation posterior membership probabilities, while k-means assigns it to one nearest centroid using squared Euclidean distance. Prefer a GMM for overlapping, elliptical, or unequal-variance clusters and density-based decisions; prefer k-means when speed, simplicity, and hard assignments matter.
EM fits a GMM by alternating two steps: the E-step computes each point's responsibility (posterior probability) under each Gaussian using current parameters, and the M-step updates the means, covariances, and mixing weights to maximize the expected log-likelihood given those responsibilities. It iterates until the likelihood converges. Because the objective is non-convex, EM only reaches a local optimum, so initialization and multiple restarts matter.
Work top-down: start at the model layer with quantization, distillation, or routing cheaper models for easy requests, since model choices drive every downstream cost. Then optimize the runtime with batching, caching, and techniques like prompt caching for LLMs, and finally match infrastructure to the load using autoscaling on queue depth and spot or batch capacity. Track cost per token or per prediction alongside latency percentiles and accuracy so optimizations never silently degrade quality.