Gaussian Processes: Why Kernels Are Priors Over Functions, Not Just Similarity Scores

Most introductions to kernels describe them as similarity measures. Two points that are close get a high kernel value, two points that are far apart get a low one, and a learning algorithm uses those similarities to interpolate. That description is not wrong, but it hides the more powerful idea behind Gaussian processes. A kernel is not just a way to compare inputs. It is a prior over functions.

Once you see a kernel that way, many design choices that look arbitrary become principled. The lengthscale is a belief about how quickly a function can change. The choice between a squared exponential kernel and a Matern kernel is a belief about smoothness. The noise term is a belief about how much of the observed variation is measurement error rather than signal. Gaussian process regression is what happens when you take those beliefs seriously, combine them with data, and let Bayes' rule do the rest.

From a prior over weights to a prior over functions

Start with Bayesian linear regression. Suppose we model an unknown function as a linear combination of basis functions:

$$f(x) = \phi(x)^{\top} w$$

Put a Gaussian prior on the weights, $w \sim \mathcal{N}(0, \Sigma_p)$ , and assume Gaussian observation noise. The posterior over weights is Gaussian, and so is the predictive distribution at a new point. That is useful, but it forces us to commit to a fixed set of basis functions.

A Gaussian process removes that commitment. Instead of placing a prior on a finite weight vector, we place a prior directly on the function itself. A Gaussian process is a collection of random variables, one for every possible input x, such that any finite subset has a joint Gaussian distribution. It is fully specified by a mean function and a covariance function:

$$f \sim \mathcal{GP}(m(x), k(x, x'))$$

Here $m(x)$ is usually taken to be zero after centring the data, and $k(x, x')$ is the kernel. For any finite set of inputs, the function values are jointly Gaussian with mean given by m and covariance given by k evaluated at every pair of points.

This is the key shift. We no longer ask what weights are plausible. We ask what functions are plausible, before seeing any data. The kernel answers that question.

What a kernel actually encodes

If $k(x, x')$ is large, then $f(x)$ and $f(x')$ are strongly correlated under the prior. Knowing one tells you a lot about the other. If it is small, they are almost independent. So the kernel controls correlation structure, and correlation structure controls the shape of typical functions drawn from the prior.

Consider the squared exponential, or radial basis function, kernel:

$$k(x, x') = \sigma_f^{2} \exp\left(-\frac{|x - x'|^{2}}{2\ell^{2}}\right)$$

The signal variance $\sigma_f^{2}$ controls the typical amplitude of functions: how far they stray from the mean. The lengthscale $\ell$ controls how quickly correlation decays with distance. A short lengthscale allows rapid variation, so sampled functions look wiggly. A long lengthscale forces slow variation, so sampled functions look smooth and almost flat over short distances.

Draw samples from a Gaussian process with a short lengthscale and you get jagged, nervous functions. Draw samples with a long lengthscale and you get gentle, slowly varying curves. Neither set of samples is "right" in the abstract. Each corresponds to a different prior belief, and the data will favour one over the other through the marginal likelihood, discussed below.

Smoothness is a separate belief from lengthscale, and this is where the Matern family matters. The squared exponential kernel implies functions that are infinitely differentiable, which is often smoother than reality. The Matern kernel with parameter $\nu$ controls differentiability directly. With $\nu = 1/2$ , samples are continuous but not differentiable, like a random walk. With $\nu = 3/2$ or $\nu = 5/2$ , samples are once or twice differentiable. Many practitioners default to Matern 5/2 for real data precisely because infinite smoothness is a strong and often false assumption.

Other kernels encode other structure. A periodic kernel encodes the belief that the function repeats. A linear kernel encodes the belief that the underlying trend is linear. Kernels can be added and multiplied, and the resulting kernel corresponds to a prior that combines those structures: for example, a smooth trend plus a periodic component plus noise. Kernel design is prior design.

Gaussian process regression in closed form

Suppose we observe noisy targets $y = f(x) + \epsilon$ , with $\epsilon \sim \mathcal{N}(0, \sigma_n^{2})$ . Let X be the training inputs, y the training targets, and X_{*} the test inputs. Under the Gaussian process prior, training and test function values are jointly Gaussian:

$$\begin{bmatrix} y \ f_{} \end{bmatrix} \sim \mathcal{N}\left(0, \begin{bmatrix} K + \sigma_n^{2} I & K_{} \ K_{}^{\top} & K_{*} \end{bmatrix}\right)$$

where K is the kernel matrix on training points, $K_{}$ is the cross-covariance between training and test points, and $K_{*}$ is the kernel matrix on test points.

Conditioning a Gaussian on the observed y gives the posterior predictive distribution, also Gaussian, with mean and covariance:

$$\bar{f}{*} = K{*}^{\top} (K + \sigma_n^{2} I)^{-1} y$$

$$\mathrm{cov}(f_{}) = K_{} - K_{}^{\top} (K + \sigma_n^{2} I)^{-1} K_{*}$$

Read these equations carefully, because they explain the behaviour of Gaussian processes intuitively. The predictive mean is a weighted combination of the observed targets, with weights determined by kernel similarity to the test point and by the inverse training covariance. Points far from all training data get predictions close to the prior mean, because their kernel values are small. The predictive variance starts at the prior variance $K_{**}$ and is reduced by an amount that depends on how well the training data inform the test point. Far from the data, uncertainty returns to the prior level. Near the data, it shrinks, but never below the noise floor.

That calibrated uncertainty is the main reason to use a Gaussian process instead of, say, kernel ridge regression with the same kernel. The mean predictions are closely related, but only the Gaussian process gives a principled variance that grows sensibly outside the data range. That matters when the prediction will drive a decision, such as choosing the next experiment in Bayesian optimisation or deciding where a model is safe to trust.

Learning hyperparameters from the marginal likelihood

Hyperparameters such as lengthscale, signal variance, and noise variance are usually learned by maximising the log marginal likelihood:

$$\log p(y \mid X) = -\frac{1}{2} y^{\top} (K + \sigma_n^{2} I)^{-1} y - \frac{1}{2} \log |K + \sigma_n^{2} I| - \frac{n}{2} \log 2\pi$$

The three terms tell a clean story. The first term rewards fit to the data. The second term penalises complexity: a kernel that can explain many different datasets gets a larger determinant and a larger penalty. The third term is a normalising constant. Maximising this objective therefore trades fit against complexity automatically, without a separate validation set. This is sometimes called the Bayesian Occam's razor.

In practice, the marginal likelihood surface can have multiple local optima, especially with several hyperparameters. A short lengthscale with low noise can explain the data by fitting every wiggle, while a long lengthscale with higher noise explains it as a smooth trend plus noise. Both may be plausible, and optimisation from different starting points can land in different optima. It is worth trying several restarts and inspecting the fitted kernel, rather than trusting a single optimisation run blindly.

A practical example

In scikit-learn, Gaussian process regression is a few lines. The code below fits a process with a Matern kernel plus a white noise kernel, optimises the hyperparameters, and plots predictions with uncertainty bands:

import numpy as np
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import Matern, WhiteKernel, ConstantKernel

X = np.linspace(0, 10, 20).reshape(-1, 1)
rng = np.random.default_rng(0)
y = np.sin(X).ravel() + rng.normal(0, 0.1, size=X.shape[0])

kernel = ConstantKernel(1.0) * Matern(length_scale=1.0, nu=2.5) + WhiteKernel(0.1)
gp = GaussianProcessRegressor(kernel=kernel, n_restarts_optimizer=5, random_state=0)
gp.fit(X, y)

X_test = np.linspace(-2, 12, 200).reshape(-1, 1)
mean, std = gp.predict(X_test, return_std=True)
print(gp.kernel_)

Two things are worth noting in this example. First, extrapolating to -2 or 12, outside the training range 0 to 10, the mean reverts toward zero and the standard deviation grows. That is the prior reasserting itself, and it is honest behaviour. Second, the fitted kernel printed at the end is the learned prior: it tells you what lengthscale and noise level the data supported. If the fitted lengthscale is implausibly small, the model may be fitting noise. If it is implausibly large, the kernel family may be too rigid.

Limitations and common pitfalls

Gaussian processes are not a universal solution, and their failures are instructive.

First, computation scales cubically. Training requires inverting an n by n matrix, which costs $O(n^{3})$ time and $O(n^{2})$ memory. Beyond a few thousand points, exact inference becomes impractical, and sparse approximations with inducing points, or alternative methods entirely, are needed. This is the main reason Gaussian processes are popular for small-data problems like Bayesian optimisation and less common for large-scale deep learning datasets.

Second, standard kernels are stationary: they depend only on the difference $x - x'$ , not on absolute position. That means the same lengthscale and variance apply everywhere. If the true function is smooth in one region and wiggly in another, a stationary kernel will compromise, and neither region will be modelled well. Non-stationary kernels and input warping exist, but they add complexity.

Third, kernels struggle in high dimensions. Distance concentrates as dimension grows, so kernel values become less discriminative, and learning a separate lengthscale per dimension (automatic relevance determination) helps only up to a point. In very high dimensions, a Gaussian process with a generic kernel is rarely competitive with a model that learns its own representation, such as a neural network.

Fourth, the prior can mislead. If you choose a squared exponential kernel for a function that is actually rough, the model will be overconfident between data points, because it believes the function must be smooth. Uncertainty quantification is only as good as the prior. Checking the model by holding out data, or by comparing its predictive intervals against actual errors, is essential rather than optional.

Finally, do not confuse the kernel trick with the Gaussian process prior. Support vector machines and kernel ridge regression use kernels to define a feature space and a regulariser. Gaussian processes use the same kernels to define a distribution over functions. The maths overlaps, but the interpretation, and the uncertainty estimates that follow from it, are different.

Takeaway

A kernel is a statement about what functions you expect before seeing data: how large they are, how fast they vary, how smooth they are, and what structure they contain. Gaussian process regression is the machinery that turns that statement, combined with observations, into predictions with honest uncertainty.

When you choose a kernel, you are not just choosing a similarity function. You are choosing a prior. Choose it deliberately, check what it implies by sampling from it, fit its hyperparameters with care, and validate its uncertainty on held-out data. Done that way, a Gaussian process is one of the clearest examples in machine learning of a model whose assumptions are visible, testable, and useful.

---

Connect with me on LinkedIn: linkedin.com/in/umeshtharukaofficial