Machine Learning 1, companion to the lecture on Gaussian processes
We know how to put a Gaussian distribution on a vector. A Gaussian process puts one on a function. This page gets there in small steps: a vector is a function of its index, a function is a vector with one entry for every point, and observing random smooth functions tells us what covariance matrix such a distribution needs.
A vector f ∈ ℝN is a list of N numbers, and we tell them apart by their index i = 1, …, N. So the vector is really a rule that takes an index and returns a number, i ↦ fi. That is a function, just one whose domain is the finite set {1, …, N}. Drawn with the index on a horizontal axis it looks like any other plot of a function, except that it is only defined at a few points.
Now let the index be a coordinate. Take a function f(x) on [0, 1] and read it off on a grid x1, …, xN. The values form a vector
f = ( f(x1), f(x2), …, f(xN) )T ∈ ℝN
whose entries are indexed by grid points rather than by integers. Refine the grid and the vector gets longer. In the limit there is one entry for every x in [0, 1], infinitely many of them, and the vector is the function itself. A finite vector is a function on a finite set of indices, and a continuous function is a vector on an infinite set of coordinates.
Suppose we observe smooth functions that come out differently every time, say the temperature over a day, one curve per day. Read each one off on the same grid and we get a new vector every day. The vectors are random because the curves are.
Start with a single grid point xi. Every observed curve gives one value f(xi), so the value at a fixed location is a random variable, and collecting many observations shows its distribution. For the curves on this page that distribution is Gaussian. That is a property of how these curves were made, not of random functions in general: other ways of generating random curves give other shapes. We assume it from here on and come back to it in the definition at the end.
Now take two grid points. Both values are random, and together they form a pair of random variables with a joint distribution. The two are not independent. When a curve is high at xi it is also high at the neighboring grid point, because the curve is smooth and the two points are close.
We can measure this. With M observed vectors f(1), …, f(M), the sample covariance between entries i and j is
Σ̂ij = 1/(M − 1) ∑m ( fi(m) − f̄i ) ( fj(m) − f̄j )
For neighboring grid points it is large, and it falls off as the two points move apart. Doing this for every pair at once gives a covariance matrix with a band along the diagonal and values close to zero away from it.
Compare this with the random vectors we already know how to draw. With f ∼ 𝒩(0, I) every entry is independent of every other entry. Plotted on the grid and connected, the vector is noise, and refining the grid only adds more independent numbers. Its covariance matrix is the identity, which no collection of smooth curves would ever produce.
A Gaussian is fixed by its mean and its covariance. So instead of measuring the covariance from curves, we can use it to make curves. Drawing vectors from 𝒩(0, Σ̂) with the measured Σ̂ gives vectors that, plotted on the grid, look like the curves we observed. With few observations Σ̂ is noisy and the samples are rougher. This suggests describing a distribution over functions by a covariance matrix.
The measured covariance belongs to one grid, and a function is not a vector of that many numbers. To get a distribution over functions we need a covariance matrix on every grid, so we have to say how the measured pattern carries over to a finer one.
One way is by grid steps. Entries one step apart get the correlation that neighbors had on the grid we measured on, two steps apart what two steps had, and so on. But one grid step is a step of 1/(N − 1) in x, so on a finer grid the same rule correlates points that are closer together, and the curves get wigglier. Put differently, the covariance is then one fixed matrix whose entries depend on i − j, and a grid of N points uses its first N rows and columns. More grid points means more entries of the same random sequence, squeezed into [0, 1]. Every grid gives a different kind of curve, and nothing settles down as N grows.
The other way is by distance in x. Points a distance d apart get the correlation that points a distance d apart had before. That makes the covariance one function of the two positions, defined on the whole square [0, 1]2, and every grid reads it off over that whole square, only at a different resolution. That function is the kernel,
Σij = k(xi, xj) = exp( −(xi − xj)2 / 2ℓ2 )
with length scale ℓ = 0.1 in the figure. The rows and columns of the covariance matrix are now indexed by coordinates too, the way the entries of the vector are. On any grid the matrix is the kernel evaluated on that grid, the Gram matrix K from the lecture on kernel methods, and refining the grid only adds rows and columns. A function is a vector with one entry for every x, and in the same way the kernel is a covariance matrix with one entry for every pair (x, x′).
Refining the grid only adds rows and columns because of the marginalization property of Gaussians. If a vector is Gaussian, any sub-vector is Gaussian too, with the mean and covariance obtained by keeping the matching entries and deleting the rest. Here that means something useful. The distribution of the values at a few chosen points does not depend on how many other grid points there were. It is 𝒩(0, K) with K built on just those points.
So we never have to say what the grid is. For any finite set of points we know the joint distribution of the function values there, and the answers for different sets agree with each other. That is all we need to call it a distribution over functions.
Definition (Gaussian process)
A Gaussian process is a collection of random variables f(x), one for every input x, any finite number of which is jointly Gaussian distributed. It is fixed by a mean function and a covariance function, the kernel,
m(x) = 𝔼[ f(x) ], k(x, x′) = cov( f(x), f(x′) )
and we write f(·) ∼ GP( m(·), k(·, ·) ). For points x1, …, xN the vector of function values is 𝒩(m, K) with mi = m(xi) and Kij = k(xi, xj).
On this page the mean function is zero everywhere and the kernel is the squared exponential above. The mean says around which curve the draws scatter. The kernel says what the draws look like: how smooth they are, how quickly they wiggle, how far they stray from the mean.
How the curves on this page are made. The observed functions are themselves draws from the Gaussian process with the kernel above, so what the figure measures can be compared with the kernel that produced it. Every curve, observed or sampled, is a single draw on a grid of 193 points. We draw a vector z ∼ 𝒩(0, I), compute the Cholesky factor K = L LT, and set f = L z, which has covariance L LT = K. A coarser grid shows a subset of the same 193 values, which by the marginalization property is a correct sample on the coarse grid.
On a grid of N points the independent sample is w = LN−1 fN, where LN is the Cholesky factor of K on that grid. This w is exactly 𝒩(0, I), and the samples from the measured covariance and from the kernel are the same N numbers multiplied by the Cholesky factor of Σ̂ or of K. A small number, 10−6, is added to the diagonal before factoring. Without it the matrix is too close to singular for the factorization to finish in floating point.