Gaussian-process regression and prediction

Let $\mathcal{X}$ be a set, let $m:\mathcal{X}\to \mathbb{R}$, and let $k:\mathcal{X}\times \mathcal{X}\to \mathbb{R}$ be symmetric and positive semidefinite, meaning that for every finite ordered list $t_1,\dots ,t_r$ in $\mathcal{X}$ and every $a\in {\mathbb{R}}^r$, ${\sum }_{i,j=1}^r{a}_i{a}_{j}k(t_i,t_j)\ge 0.$ For finite ordered lists $A=(a_i{)}_{i=1}^r$ and $B=(b_j{)}_{j=1}^s$ in $\mathcal{X}$, define $K_{AB}=(k(a_i,b_j){)}_{i,j};$ repeated indices are allowed.

Given a real Gaussian process $f=(f_t{)}_{t\in \mathcal{X}}$ with mean $m(t)=\mathbb{E}[f_t]$ and covariance $k(s,t)=\mathbb{E}[(f_s-m(s))(f_t-m(t))]$, and a nonempty finite training list $X=(x_i{)}_{i=1}^n$ in $\mathcal{X}$, define $f_X=(f_{x_i}{)}_{i=1}^n$ and $m_X=(m(x_i){)}_{i=1}^n$. For a centered Gaussian noise vector $\epsilon \sim N_n(0,R)$ with symmetric positive-semidefinite covariance $R$ that is independent of $f$, define the observation vector $Y=f_X+\epsilon$ and the matrix $A_X=K_{XX}+R$.

For a nonempty finite test list $X_{\ast }=(x_j^{\ast }{)}_{j=1}^p$ in $\mathcal{X}$, define the latent test vector $f_{\ast }=f_{X_{\ast }}$. If a centered Gaussian vector ${\eta }_{\ast }\sim N_p(0,R_{\ast })$ with symmetric positive-semidefinite covariance $R_{\ast }$ is independent of the pair $(f,\epsilon )$, define the noisy test vector $Y_{\ast }=f_{\ast }+{\eta }_{\ast }$.

A Gaussian process has multivariate-normal finite evaluation vectors, including singular laws and repeated indices.

Finite Gaussian laws allow positive-semidefinite, possibly singular covariance matrices.

Independence means independence of the sigma-algebras generated by the random elements.

Rasmussen and Williams, GPML, Chapter 2 §2.2, printed pp. 13–18

Bishop, PRML, §6.4.2, printed pp. 306–310

diagram failed to render:
\begin{tikzpicture}[>=stealth]
  \node[draw,rounded corners] (prior) at (0,0) {$f\sim\mathrm{GP}(m,k)$};
  \node[draw] (train) at (-2,-1.5) {$f_X$};
  \node[draw] (test) at (3,-1.5) {$f_*$};
  \node[draw] (noise) at (-2,-3) {$\varepsilon$};
  \node[draw,rounded corners] (obs) at (0,-3) {$Y=f_X+\varepsilon$};
  \node[draw] (testnoise) at (3,-3) {$\eta_*$};
  \node[draw,rounded corners] (testobs) at (5,-3) {$Y_*=f_*+\eta_*$};
  \draw[->] (prior) -- (train);
  \draw[->] (prior) -- (test);
  \draw[->] (train) -- (obs);
  \draw[->] (noise) -- (obs);
  \draw[->] (test) -- (testobs);
  \draw[->] (testnoise) -- (testobs);
\end{tikzpicture}

The single prior node represents one process, so $f_X$ and $f_{\ast }$ may be correlated. Arrows into $Y$ and $Y_{\ast }$ denote vector addition. The training noise $\epsilon$ is independent of $f$, and the test noise ${\eta }_{\ast }$ is independent of the pair $(f,\epsilon )$; the diagram records no posterior independence claim.

Let $H$ be a real inner-product space and $V\subseteq H$ a finite-dimensional subspace. For each $f\in H$ there is a unique $p\in V$ such that $f-p$ is orthogonal to $V$, meaning $⟨f-p,v⟩=0$ for every $v\in V$. For every $v\in V$, $\mathrm{\parallel }f-v{\mathrm{\parallel }}^2=\mathrm{\parallel }f-p{\mathrm{\parallel }}^2+\mathrm{\parallel }p-v{\mathrm{\parallel }}^2\mathrm{.}$ In particular, $\mathrm{\parallel }f{\mathrm{\parallel }}^2=\mathrm{\parallel }p{\mathrm{\parallel }}^2+\mathrm{\parallel }f-p{\mathrm{\parallel }}^2$.

Every finite-dimensional inner-product space has an orthonormal basis, with the empty basis in dimension zero.

Treil, Linear Algebra Done Wrong, Chapter 5 §5.3.2

In the finite observation model, let $K=K_{XX}$ and let $R$ be symmetric and positive semidefinite. Then $A=K+R$ is symmetric and positive semidefinite. Its ordinary inverse exists if and only if $A$ is positive definite, equivalently if and only if all singular values of $A$ are positive. Thus the standard inverse formulas are defined exactly under $A\succ 0$; no ordinary inverse formula applies when $A$ is singular.

The finite Gram and observed matrices are K_XX and A_X=K_XX+R, with positive-semidefinite kernel and noise covariances.

Rasmussen and Williams, GPML, Chapter 4 §4.1, printed pp. 79–80

Bishop, PRML, §6.4.2, printed pp. 306–309

Assume AC. Let $\mathcal{X}$ be a set, let $f=(f_t{)}_{t\in \mathcal{X}}$ be a real Gaussian process with mean $m$ and symmetric positive-semidefinite covariance kernel $k$, and let $X=(x_i{)}_{i=1}^n\in {\mathcal{X}}^n$ with $n\ge 1$. For independent observation noise $\epsilon \sim N_n(0,R)$ with symmetric positive-semidefinite covariance $R$, set $Y=f_X+\epsilon$, $m_X=(m(x_i))$, $K_{XX}=(k(x_i,x_j))$, and $A_X=K_{XX}+R$. Suppose $A_X$ is positive definite. Define the log marginal likelihood at $y\in {\mathbb{R}}^n$ as the logarithm of the Lebesgue density of $Y$ at $y$ under this prior-predictive model. Then $Y\sim N_n(m_X,A_X)$ and $\log p(y\mid X,m,k,R)=-\frac{1}{2}(y-m_X{)}^T{A}_X^{-1}(y-m_X)-\frac{1}{2}\log \det (A_X)-\frac{n}{2}\log (2\pi )\mathrm{.}$ This formula is asserted only for the nonsingular observed-covariance case $A_X\succ 0$.

Rasmussen and Williams, GPML, Chapter 2 §§2.2–2.3, equations (2.28)–(2.30), printed pp. 18–20

Bishop, PRML, §6.4.2, equations (6.61)–(6.62), printed pp. 306–307

Assume AC. Let $\mathcal{X}$ be a set, $m:\mathcal{X}\to \mathbb{R}$, and $k:\mathcal{X}\times \mathcal{X}\to \mathbb{R}$ a symmetric kernel whose every finite Gram quadratic form on points of $\mathcal{X}$ is nonnegative. Let $f=(f_t{)}_{t\in \mathcal{X}}$ be the canonical Gaussian process with $m(t)=\mathbb{E}[f_t]$ and $k(s,t)=\mathbb{E}[(f_s-m(s))(f_t-m(t))]$. For finite ordered lists $A=(a_i{)}_{i=1}^r$ and $B=(b_j{)}_{j=1}^s$ in $\mathcal{X}$, set $K_{AB}=(k(a_i,b_j))$; repetitions are allowed. Take $X=(x_i{)}_{i=1}^n\in {\mathcal{X}}^n$ and $X_{\ast }=(x_j^{\ast }{)}_{j=1}^p\in {\mathcal{X}}^p$ with $n,p\ge 1$, and write $m_X=(m(x_i))$, $m_{\ast }=(m(x_j^{\ast }))$, $f_X=(f_{x_i})$, and $f_{\ast }=(f_{x_j^{\ast }})$. Let $\epsilon \sim N_n(0,R)$ be independent of $f$, where $R$ is symmetric positive semidefinite, and set $Y=f_X+\epsilon$ and $A_X=K_{XX}+R$. Suppose $A_X$ is positive definite. Then for every $y\in {\mathbb{R}}^n$ the conditional law of $f_{\ast }$ given $Y=y$ has mean $m_{\ast }+K_{\ast X}{A}_X^{-1}(y-m_X)$ and covariance $K_{\ast \ast }-K_{\ast X}{A}_X^{-1}{K}_{X\ast }\mathrm{.}$ This covariance is positive semidefinite and may be singular. If ${\eta }_{\ast }\sim N_p(0,R_{\ast })$ is independent of the pair $(f,\epsilon )$, where $R_{\ast }$ is symmetric positive semidefinite, then the conditional law of $f_{\ast }+{\eta }_{\ast }$ has the same mean and covariance with $R_{\ast }$ added. The displayed kernel is defined for every $y$, although regular conditional versions are determined only $P_Y$-almost surely.

Rasmussen and Williams, GPML, Chapter 2 §2.2, equations (2.18)–(2.24), printed pp. 15–18

Bishop, PRML, §6.4.2, equations (6.61)–(6.68), printed pp. 306–309

A real Gaussian process is a family whose finite evaluation vectors are multivariate normal. Let $T$ be a set and let $f=(f_t{)}_{t\in T}$ be a real Gaussian process with mean $m(t)=\mathbb{E}{f}_t$ and symmetric positive-semidefinite covariance kernel $k(s,t)=\mathbb{E}[(f_s-m(s))(f_t-m(t))]$; positive semidefinite means every finite Gram quadratic form is nonnegative. For $n\ge 1$ and a finite training list $X=(x_i{)}_{i=1}^n$ in $T$, set $f_X=(f_{x_i})$, $K_{XX}=(k(x_i,x_j))$, and let $\epsilon$ be centered Gaussian observation noise independent of $f$, with symmetric positive-semidefinite covariance $R$. Set $Y=f_X+\epsilon$ and $A_X=K_{XX}+R\succ 0$. Write the dense Cholesky factorization as $A_X=L{L}^T$, where $L$ is lower triangular with positive diagonal. In the scalar-operation model, where each addition, subtraction, multiplication, division, or square root counts as one operation, factorization costs $n^3\mathrm{/}3+O(n^2)$ operations. For a test input $x_{\ast }\in T$, define $k_{x_{\ast }X}=(k(x_{\ast },x_i){)}_{i=1}^n$ and $k_{X{x}_{\ast }}=k_{x_{\ast }X}^T$. After factorization, the latent predictive variance $k(x_{\ast },x_{\ast })-k_{x_{\ast }X}{A}_X^{-1}{k}_{X{x}_{\ast }}$ requires $O(n^2)$ operations, where $O(n^2)$ denotes an asymptotic operation count as $n$ grows. These are arithmetic counts, not accuracy guarantees. For an invertible real $n\times n$ matrix $A$, define its spectral condition number by ${\kappa }_2(A)=\parallel A{\parallel }_2\parallel A^{-1}{\parallel }_2$, where $\parallel A{\parallel }_2={\sup }_{\parallel x{\parallel }_2=1}\parallel Ax{\parallel }_2$ for $x\in {\mathbb{R}}^n$. For $0<ϵ<1$, $A_{ϵ}=\mathrm{diag}(1,ϵ)$ has condition number $1\mathrm{/}ϵ$, while $E=\mathrm{diag}(0,-ϵ)$ has relative spectral-norm size $\parallel E{\parallel }_2\mathrm{/}\parallel A_{ϵ}{\parallel }_2=ϵ$ even though $A_{ϵ}+E$ is singular.

The negative log marginal likelihood need not be convex in a covariance hyperparameter. On the one-point set $T=\{t\}$, take the centered Gaussian prior $f_t\sim N(0,\theta )$ with $\theta \ge 0$, covariance kernel $k_{\theta }(t,t)=\theta$, and independent observation noise $\epsilon \sim N(0,1)$. For observation $Y=f_t+\epsilon$ at $y=0$, the observed covariance is $A_{\theta }=1+\theta >0$ and the negative log marginal likelihood is $L(\theta )=\frac{1}{2}\log (1+\theta )+\frac{1}{2}\log (2\pi )\mathrm{.}$ For every $\theta >0$, its second derivative is $L^{\mathrm{\prime }\mathrm{\prime }}(\theta )=-1\mathrm{/}(2(1+\theta {)}^2)<0$, so $L$ is not convex on its parameter domain.

The finite Gaussian-process regression posterior formula uses Gaussian observation noise. Write $\varphi (z)=(2\pi {)}^{-1\mathrm{/}2}{e}^{-z^2\mathrm{/}2}$ for the standard normal density and $\mathrm{\Phi }(z)={\int }_{-\mathrm{\infty }}^z\varphi (t)dt$ for its cumulative distribution function. With one latent value $z\sim N(0,1)$ and one positive label, the logistic likelihood $P(y=1\mid z)=\sigma (z)=(1+e^{-z}{)}^{-1}$ or probit likelihood $P(y=1\mid z)=\mathrm{\Phi }(z)$ gives a non-Gaussian posterior: the posterior log density has nonzero third derivative at $z=\log 3$ for logistic and $z=0$ for probit. Both likelihoods lie in $(0,1)$, so multiplying by the standard normal density gives a finite, positive normalizing integral. Finally, a Gaussian process specifies finite-dimensional laws but need not have continuous paths. A continuous modification means a process on the same probability space that agrees almost surely at each fixed input and has almost surely continuous sample paths. Under AC, on $T=[0,1]$ the equality kernel $k(s,t)=1$ when $s=t$ and $0$ otherwise has no continuous modification.

Rasmussen and Williams, GPML, Chapter 2 §§2.2–2.3, printed pp. 13–21

Rasmussen and Williams, GPML, Chapter 3 §§3.3–3.4.1, printed pp. 39–43

Rasmussen and Williams, GPML, Chapter 4 §§4.1–4.1.1, printed pp. 79–81

Bishop, PRML, §6.4.2, printed pp. 306–310

Driscoll and Braun, Fundamentals of Numerical Computation, §§2.5, 2.9

Driscoll and Braun, Fundamentals of Numerical Computation, §2.9.3

Assume AC. Let $T$ be a set, let $k:T\times T\to \mathbb{R}$ be a symmetric kernel whose every finite Gram quadratic form is nonnegative, and let $m:T\to \mathbb{R}$. Let $H$ be a real Hilbert space of functions on $T$ such that, for every $x\in T$, the section $k_x:T\to \mathbb{R}$, $k_x(t)=k(t,x)$, belongs to $H$ and satisfies $h(x)=⟨h,k_x{⟩}_H$ for every $h\in H$. Assume $m\in H$. For a nonempty finite ordered list $X=(x_i{)}_{i=1}^n$ (repetitions allowed), let $K=K_{XX}=(k(x_i,x_j))$, $m_X=(m(x_i))$, and take $y\in {\mathbb{R}}^n$ and $\sigma >0$. On the affine space $m+H$, define $J(g)={\sum }_{i=1}^n(y_i-g(x_i){)}^2+{\sigma }^2\parallel g-m{\parallel }_H^2\mathrm{.}$ Then $J$ has a unique minimizer $g_{\ast }=m+{\sum }_{i=1}^n{\alpha }_i{k}_{x_i},\alpha =(K+{\sigma }^2{I}_n{)}^{-1}(y-m_X)\mathrm{.}$ Its evaluation vector is $(g_{\ast }(x_i){)}_{i=1}^n=m_X+K(K+{\sigma }^2{I}_n{)}^{-1}(y-m_X)\mathrm{.}$ Under AC, the mean and kernel determine a canonical cylinder-space Gaussian-process prior. If it is observed at $X$ with independent $N_n(0,{\sigma }^2{I}_n)$ noise, this vector is the conditional mean of the latent evaluation vector given the observation $Y=y$. The Hilbert-space existence and reproducing-section properties are hypotheses here; the theorem makes no general kernel-to-Hilbert-space existence claim.

Rasmussen and Williams, GPML, Chapter 6 §§6.1–6.4, printed pp. 129–150

Hastie, Tibshirani, and Friedman, ESL, 2nd ed., §5.8, printed pp. 167–173

Treil, Linear Algebra Done Wrong, Chapter 5 §§5.3.2, 5.4.1