Kernels for Gaussian Processes

Kernels for Gaussian Processes

A kernel is a symmetric positive semi-definite function $K:\Omega \times \Omega \to \mathbb{R}$ that specifies the covariance structure of a Gaussian process: $K(x,x^{\prime })=\text{Cov}(f(x),f(x^{\prime }))$. Positive semi-definite means that for any finite collection of points $u_1,\dots ,u_N\in \Omega$, the Gram matrix $\mathbf{K}=(K(u_i,u_j){)}_{i,j}$ is positive semi-definite. This is the only restriction: every valid kernel defines a valid GP.

The kernel encodes every prior belief about the function — smoothness, stationarity, length-scale, periodicity, boundary conditions. Canonical kernels on $\Omega =[0,\infty )$:

• Brownian motion: $K(s,t)=\min (s,t)$. Samples are continuous but nowhere differentiable; $f(0)=0$.
• Scaled BM plus constant: $K(s,t)=C+{\tau }^2\min (s,t)$. Adds an $\mathcal{N}(0,C)$ intercept so $f(0)\ne 0$; with $C\to \infty$ the intercept is flat.
• Integrated Brownian Motion (IBM): $K(s,t)=\frac{1}{2}\max (s,t)\min (s,t{)}^2-\frac{1}{6}\min (s,t{)}^3$. Samples are once-differentiable; $f(0)=f^{\prime }(0)=0$.
• Scaled IBM plus linear term: $K(s,t)=C(1+st)+{\tau }^2\left[\frac{1}{2}\max (s,t)\min (s,t{)}^2-\frac{1}{6}\min (s,t{)}^3\right]$. Adds $\mathcal{N}(0,C)$ priors on intercept and slope; used in 3.10.2 GP Regression because its posterior mean is a natural cubic spline.

Canonical kernels on $\Omega \subseteq {\mathbb{R}}^d$:

• Linear: $K(\mathbf{x},{\mathbf{x}}^{\prime })={\mathbf{x}}^{\top }{\mathbf{x}}^{\prime }$. GP is the Bayesian linear-regression function.
• Squared-exponential (RBF): $K(\mathbf{x},{\mathbf{x}}^{\prime })={\sigma }_f^2\exp \left(-\Vert \mathbf{x}-{\mathbf{x}}^{\prime }{\Vert }^2/(2{\ell }^2)\right)$. Samples are infinitely differentiable; $\ell$ is the length-scale.
• Matérn: smoothness-tunable generalization of RBF. At half-integer smoothness parameters, has closed-form simple expressions.

Valid kernels form a cone under three operations: sums $K_1+K_2$, products $K_1\cdot K_2$, and scalar multiples $\alpha K$ (for $\alpha \ge 0$). This kernel algebra allows building complex priors from simple parts.

Prerequisites

• 3.10 Gaussian Processes — parent index; the definition of a GP.
• 6.2 Multivariate Gaussian Distribution — what “positive semi-definite covariance” means.
• 4.6 Brownian Motion — the underlying stochastic process for most non-stationary kernels in this note.

Intuition

A kernel $K(x,x^{\prime })$ answers one question: if $f(x)$ takes a value above its mean, how confident am I that $f(x^{\prime })$ does too? Large $K(x,x^{\prime })$ means the function values are tightly correlated; small $K(x,x^{\prime })$ means they are nearly independent. The kernel shape translates prior knowledge about the function into correlations between function values at pairs of inputs.

Three axes govern kernel design. Smoothness determines how many times sample paths can be differentiated: Brownian samples are continuous but not differentiable, IBM samples are once differentiable, RBF samples are $C^{\infty }$. Stationarity determines whether $K(x,x^{\prime })$ depends only on the separation $x-x^{\prime }$ or on the absolute positions; stationary kernels (RBF, Matérn) treat the domain homogeneously, non-stationary kernels (Brownian, IBM) let the variance grow away from the origin. Boundary conditions fall out of the kernel: Brownian motion pins $f(0)=0$, IBM pins both $f(0)=0$ and $f^{\prime }(0)=0$, and adding a linear term with flat priors on the coefficients releases these constraints.

The positive semi-definiteness requirement is not a technicality — it is the statement that every finite-dimensional marginal of the induced GP is a proper Gaussian. Violate it and the implied “distribution” assigns negative variance to some linear combination of function values, which is nonsense.

Properties

(Kernel Algebra) Sums, products, and nonneg. scalar multiples of valid kernels are valid kernels.

If $K_1$ and $K_2$ are positive semi-definite, so are $K_1+K_2$, $K_1\cdot K_2$, and $\alpha K_1$ for $\alpha \ge 0$. This is the operational fact that lets one build elaborate kernels from simple parts. Summing kernels yields additive GPs — $f=f_1+f_2$ with $f_j\sim \text{GP}(0,K_j)$ independent. Multiplying kernels combines effects, for example trend times periodic component. The scaled-IBM-plus-linear-term kernel $K(s,t)=C(1+st)+{\tau }^2{K}_I(s,t)$ is literally the sum of the “constant $+$ linear” kernel $C(1+st)$ and the IBM kernel, reflecting the decomposition $f(x)={\beta }_0+{\beta }_{1}x+\tau I(x)$ with independent components.

(Brownian Motion Kernel) Continuous non-differentiable samples; variance grows linearly.

The kernel $K(s,t)=\min (s,t)$ is the covariance of standard Brownian motion $W_t$ on $[0,\infty )$. It enforces $f(0)=0$ almost surely (since $K(0,0)=0$) and gives $\text{Var}(f(t))=t$, so the variance grows linearly with $t$. Positive semi-definiteness follows directly from the independent-increments representation: any linear combination $\sum a_{i}f(t_i)$ rewrites as a linear combination of independent Gaussian increments. Sample paths are continuous but nowhere differentiable — a consequence of the square-root scaling of increments, $W_{t+h}-W_t\sim \mathcal{N}(0,h)$, which has magnitude $\sqrt{h}$ rather than $h$. This makes Brownian motion a very rough prior, appropriate when the function is expected to fluctuate on every scale. See 3.10.3 Integrated Brownian Motion for a smoother alternative.

(Integrated Brownian Motion Kernel) Once-differentiable samples; cubic-spline posterior mean.

The IBM kernel is

$$K_I(s,t)=\frac{1}{2}\max (s,t)\min (s,t{)}^2-\frac{1}{6}\min (s,t{)}^3=\frac{1}{6}\min (s,t{)}^2(3\max (s,t)-\min (s,t)).$$

It is obtained by integrating Brownian motion once: $I(t)={\int }_0^t{W}_{s}ds$. Since the integral of a continuous path is differentiable, IBM samples are $C^1$. The kernel pins $f(0)=0$ and $f^{\prime }(0)=0$, so in practice one uses the scaled-IBM-plus-linear-term kernel to release these constraints. The headline property: with the scaled-IBM-plus-linear-term kernel, the posterior mean for the GP regression problem is the natural cubic spline interpolant — the same smoother produced by penalized-splines methods in frequentist nonparametric regression. Full details in 3.10.3 Integrated Brownian Motion.

(Squared-Exponential / RBF Kernel) Infinitely smooth stationary default.

The squared-exponential kernel

$$K(\mathbf{x},{\mathbf{x}}^{\prime })={\sigma }_f^2\exp \left(-\frac{\Vert \mathbf{x}-{\mathbf{x}}^{\prime }{\Vert }^2}{2{\ell }^2}\right)$$

depends only on the Euclidean distance $\Vert \mathbf{x}-{\mathbf{x}}^{\prime }\Vert$, so it is stationary and isotropic. Samples are infinitely differentiable — in fact analytic — which makes RBF the default kernel for smooth regression problems. The length-scale $\ell$ controls how rapidly correlations decay: small $\ell$ gives wiggly samples, large $\ell$ gives nearly-constant samples. The signal scale ${\sigma }_f^2$ controls the overall variance. In high dimensions the RBF is often equipped with automatic relevance determination (ARD) weights, giving a separate length-scale per dimension: $K(\mathbf{x},{\mathbf{x}}^{\prime })={\sigma }_f^2\exp (-{\sum }_d(x_d-x_d^{\prime }{)}^2/(2{\ell }_d^2))$, which shrinks irrelevant dimensions to zero influence.

(Linear Kernel) GP with a linear kernel is Bayesian linear regression.

For $K(\mathbf{x},{\mathbf{x}}^{\prime })={\mathbf{x}}^{\top }{\mathbf{x}}^{\prime }$, the induced GP is $f(\mathbf{x})={\mathbf{x}}^{\top }\mathbf{\beta }$ with $\mathbf{\beta }\sim \mathcal{N}(0,\mathbf{I})$. The covariance reproduces: $\text{Cov}({\mathbf{x}}^{\top }\mathbf{\beta },{{\mathbf{x}}^{\prime }}^{\top }\mathbf{\beta })={\mathbf{x}}^{\top }{\mathbf{x}}^{\prime }$. This collapses GP regression to 3.8.1 Bayesian Linear Regression with a standard-normal prior. The linear kernel is degenerate — only $d$-dimensional, not infinite-dimensional — so its Gram matrices have rank at most $d$. This observation generalizes: a kernel $K(\mathbf{x},{\mathbf{x}}^{\prime })=\mathbf{\varphi }(\mathbf{x}{)}^{\top }\mathbf{\varphi }({\mathbf{x}}^{\prime })$ for any feature map $\mathbf{\varphi }$ is a linear kernel in $\mathbf{\varphi }$-space, and the GP is $f(\mathbf{x})=\mathbf{\varphi }(\mathbf{x}{)}^{\top }\mathbf{\beta }$.

(Non-Stationarity of BM/IBM Kernels) Variance depends on position, not just separation.

Brownian and Integrated Brownian kernels have $K(s,s)=s$ and $K(s,s)=s^3/3$ respectively — the variance grows with position, so the prior is non-stationary. This is usually undesirable when the domain is a bounded interval with no special “origin,” but it is exactly what is wanted when $x=0$ is a meaningful initial time and the function accumulates uncertainty as $x$ grows. For Bayesian nonparametric regression on a generic interval, shift-invariance is usually preferred and the RBF or Matérn kernel is the default choice. Non-stationary kernels shine in time-series settings where the origin is a real anchor.

(Mercer Decomposition) Every kernel has a spectral eigenfunction decomposition.

By Mercer’s theorem, a continuous positive semi-definite kernel on a compact domain admits

$$K(x,x^{\prime })=\sum _{k=1}^{\infty }{\lambda }_k{\varphi }_k(x){\varphi }_k(x^{\prime }),$$

where $\{{\varphi }_k\}$ is an orthonormal basis of eigenfunctions and ${\lambda }_k\ge 0$. A GP with this kernel can be written as $f(x)={\sum }_k\sqrt{{\lambda }_k}{Z}_k{\varphi }_k(x)$ with $Z_k\stackrel{\text{iid}}{\sim }\mathcal{N}(0,1)$. This connects GPs to Bayesian linear regression in an infinite-dimensional feature space and underlies the theoretical analysis of GP regression (effective degrees of freedom, rates of convergence). For the Brownian motion kernel on $[0,1]$, the eigenfunctions are sines and the eigenvalues decay as $1/k^2$. For the IBM kernel, eigenvalues decay as $1/k^4$ — faster decay reflects smoother samples.

(Hyperparameter Estimation) Kernel parameters are learned via marginal likelihood.

Kernels typically depend on free hyperparameters: length-scale $\ell$, signal variance ${\sigma }_f^2$, ratio $\gamma =\tau /\sigma$. These are estimated by maximizing (or placing priors on) the marginal likelihood $p(\mathbf{y}\mid \mathbf{\theta })=\mathcal{N}(\mathbf{y};0,{\mathbf{K}}_{\mathbf{\theta }}+{\sigma }^2\mathbf{I})$. On a log scale, the marginal likelihood is

$$\log p(\mathbf{y}\mid \mathbf{\theta })=-\frac{1}{2}{\mathbf{y}}^{\top }({\mathbf{K}}_{\mathbf{\theta }}+{\sigma }^2\mathbf{I}{)}^{-1}\mathbf{y}-\frac{1}{2}\log |{\mathbf{K}}_{\mathbf{\theta }}+{\sigma }^2\mathbf{I}|-\frac{n}{2}\log 2\pi ,$$

which has a built-in Occam penalty: the log-determinant term shrinks with overly flexible kernels. See 3.7 Marginal Likelihood for why this counts as automatic model selection.

Examples

(Easy) Verify positive semi-definiteness of the Brownian motion kernel.

Show that $K(s,t)=\min (s,t)$ is positive semi-definite on $[0,\infty )$.

Solution. Take any $t_1<\cdots <t_N$ and coefficients $a_1,\dots ,a_N$. Set ${\alpha }_i=W_{t_i}-W_{t_{i-1}}$ (with $t_0=0$), so $W_{t_i}={\alpha }_1+\cdots +{\alpha }_i$ and the ${\alpha }_j$ are independent with ${\alpha }_j\sim \mathcal{N}(0,t_j-t_{j-1})$. Then

$$\sum _{i,j}{a}_i{a}_{j}K(t_i,t_j)=\text{Var}\left(\sum _i{a}_i{W}_{t_i}\right)=\text{Var}\left(\sum _i{a}_i\sum _{j\le i}{\alpha }_j\right)=\sum _j(t_j-t_{j-1})(\sum _{i\ge j}{a}_i{)}^2\ge 0.$$

The equality on the right uses independence of the ${\alpha }_j$. This confirms $K$ is positive semi-definite.

(Medium) Build a kernel with a periodic and trend component.

Suppose $f:[0,10]\to \mathbb{R}$ is believed to consist of a smooth long-term trend plus a periodic oscillation of period roughly 1. Construct a kernel.

Solution. Use the sum $K=K_{\text{trend}}+K_{\text{periodic}}$ with

$$K_{\text{trend}}(s,t)={\sigma }_1^2\exp (-(s-t{)}^2/(2{\ell }_1^2))\text{(RBF with large length-scale ℓ1≈5)}$$ $$K_{\text{periodic}}(s,t)={\sigma }_2^2\exp \left(-\frac{2{\sin }^2(\pi (s-t))}{{\ell }_2^2}\right)\text{(periodic kernel with period 1)}.$$

The trend component captures slow variation; the periodic component captures the one-unit cycle. Because valid kernels are closed under sums, $K$ is positive semi-definite. Fitting to data requires learning ${\sigma }_1^2,{\ell }_1,{\sigma }_2^2,{\ell }_2$ via marginal likelihood — see 3.10.2 GP Regression.

(Hard) GP with RBF kernel on one point recovers the Normal-Normal model.

Show that if $n=1$, observing $y_1=f(x_1)+ϵ$ under a zero-mean GP with RBF kernel recovers the univariate 3.5.2 Normal-Normal Model for $f(x_1)$.

Solution. With $n=1$, the Gram matrix is the scalar $\mathbf{K}=K(x_1,x_1)={\sigma }_f^2$. The prior on $f(x_1)$ is $\mathcal{N}(0,{\sigma }_f^2)$. Observing $y_1\mid f(x_1)\sim \mathcal{N}(f(x_1),{\sigma }^2)$, the posterior is

$$f(x_1)\mid y_1\sim \mathcal{N}\left(\frac{{\sigma }_f^2}{{\sigma }_f^2+{\sigma }^2}{y}_1,\frac{{\sigma }_f^2{\sigma }^2}{{\sigma }_f^2+{\sigma }^2}\right).$$

This is exactly the Normal-Normal posterior for unknown mean with ${\mu }_0=0$, ${\tau }_0^2={\sigma }_f^2$, noise ${\sigma }^2$. The posterior mean shrinks $y_1$ toward the prior mean (zero) with weight ${\sigma }^2/({\sigma }_f^2+{\sigma }^2)$. At a new test point $x\ne x_1$, the prediction is ${\mathbf{k}}^{\top }(\mathbf{K}+{\sigma }^2{)}^{-1}{y}_1=(K(x,x_1)/({\sigma }_f^2+{\sigma }^2))y_1$, which interpolates between zero (far from $x_1$) and the shrunken observation (near $x_1$).

See Also

• 3.10 Gaussian Processes — parent index.
• 3.10.2 GP Regression — posterior formulas and hyperparameter estimation.
• 3.10.3 Integrated Brownian Motion — detailed treatment of the IBM kernel.
• 6.2 Multivariate Gaussian Distribution — conditioning formula that drives all posterior calculations.
• 4.6 Brownian Motion — underlying stochastic process for the canonical non-stationary kernels.