In the following, scalars are denoted with italic lowercases (e.g., $k,\, b_f$), vectors with bold lowercases (e.g., $\mathbf{s},\, \mathbf{x}_i$), and matrices with italic uppercases (e.g., $W_f$).
Objective
A Gaussian process (GP) is defined as a collection of random variables, any finite number of which have a joint Gaussian distribution. A GP $f(\mathbf{x})$ is completely specified by its mean function $m(\mathbf{x})$ and covariance function $k(\mathbf{x}, \mathbf{x}'),$ also called kernel, defined as: \begin{align*} m(\mathbf{x}) &= \mathbb{E}[f(\mathbf{x})], \\ k(\mathbf{x}, \mathbf{x}') &= \mathbb{E}[(f(\mathbf{x})-m(\mathbf{x}))(f(\mathbf{x}')-m(\mathbf{x}'))]. \end{align*}
Let:
- ${X = (\mathbf{x}_1,\dotsc,\mathbf{x}_q)}$ be the training inputs
- $ {\mathbf{f} = (f(\mathbf{x}_1)\dotsc,f(\mathbf{x}_q))}$ be the training outputs
- $X^* = (\mathbf{x}_{q+1},\dotsc,\mathbf{x}_{s})$ be the test inputs
- $\mathbf{f}^* = (f(\mathbf{x}_{q+1})\dotsc,f(\mathbf{x}_{s}))$ be the test outputs.
Assume the training set isn't contaminated with samples from the test set, i.e.: $X \cup X^* = \mathcal{X},$ and $X \cap X^* = \emptyset$.
Note that $\mathbf{f}$ is known, and $\mathbf{f}^*$ is unknown. The goal is to find the distribution of $\mathbf{f}^*$ given $X^*, X$ and $\mathbf{f}$.
Solution
The joint distribution of $\mathbf{f}$ and $\mathbf{f}^*$ according to the prior is
\begin{equation*} \begin{bmatrix} \mathbf{f} \\ \mathbf{f}^* \end{bmatrix} \sim \mathcal{N} \left( \begin{bmatrix} \mathbf{m} \\ \mathbf{m}^* \end{bmatrix}, ~ \begin{bmatrix} K(X,X) & K(X,X^*) \\ K(X^*,X) & K(X^*,X^*) \\ \end{bmatrix} \right) \end{equation*} where $\mathbf{m}\,,\,\mathbf{m}^*$ is a vector of the means evaluated at all training and test points respectively, and $K(X,X^*)$ denotes the $q \times q^*$ matrix of the covariances evaluated at all pairs of training and test points, and similarly for $K(X,X),\, K(X^*,X)$ and $K(X^*,X^*)$.
Each element of the matrix $K(X,X^*)$ are denoted by $k(\mathbf{x}, \mathbf{x}')$. Ideally one should have $k(\mathbf{x}, \mathbf{x}') = \mathbb{E}[(f(\mathbf{x})-m(\mathbf{x}))(f(\mathbf{x}^*)-m(\mathbf{x}^*))]$ but this isn't feasible since one does not have access to $f(\mathbf{x}^*)$ and $m(\mathbf{x}^*)$. As a result, one resorts to some approximation of $k$ with some kernels taking only $\mathbf{x}$ and $\mathbf{x}^*$ as input e.g.:
- Linear: $k(\mathbf{x}, \mathbf{x}') = \mathbf{x}^T \mathbf{x}'$
- Cubic: $k(\mathbf{x}, \mathbf{x}') = 3 \left (\left (\mathbf{x}^T \mathbf{x}' \right )^2 + 2\left ( \mathbf{x}^T \mathbf{x}' \right )^3 \right )$
- Absolute exponential: $k(\mathbf{x}, \mathbf{x}') = e^{|\mathbf{x}-\mathbf{x}'|} $
- Squared exponential: $k(\mathbf{x}, \mathbf{x}') = e^{-0.5|\mathbf{x}-\mathbf{x}'|^2} $
- etc
Conditioning the joint Gaussian prior on the observations yields $\mathbf{f}^* | X^*, X, \mathbf{f} \, \sim \, \mathcal{N} (\mathbf{\mu}, \mathbf{\Sigma})$ where
\begin{align} \mathbf{\mu} &= \mathbf{m}^* - K(X^*,X) K(X,X)^{-1}(\mathbf{f} - \mathbf{m}), \\ \mathbf{\Sigma} &= K(X^*,X^*) - K(X^*,X) K(X,X)^{-1} K(X,X^*). \notag \end{align}
Question
How am I supposed to compute the vector of the means $\mathbf{m}^*$?
