Regression

• Depth data from monocular images
• Mapping terrain
• Mapping an environmental quantity (such as temperature)
• Classification

1. \(y = w_0 + w_1 x_1 + \cdots + w_n x_n\). The goal is to find the weights \(w\).
2. By setting \(x_0 = 1\), we can write the model as \(y = w^T x\).

1. Stack all the \(x\) vectors: \(X = \begin{bmatrix}x_1 & \dots & x_n\end{bmatrix}\)
2. Stack all the measurements \(y = \begin{bmatrix}y_1 & \dots y_n \end{bmatrix}\)

3. The overall equation is

   \begin{equation} y = X^T w \end{equation}
4. Note that in general, \(X^T \in \mathbb{R}^{n \times D}\)

The regression equation \(y = X^T w\) cannot be solved precisely because there are more measurements than independent variables. Instead we try to minimize the square error:

\begin{equation} w_\ast = \mathrm{arg}\min_w ||X^Tw - y||_2^2 \end{equation}

So \(w_\ast = (X X^T)^{-1} X y\) when \(n > D\) (there are more measurements then the dimension of the data \(x\)).

1. The model: \(y_i = f(x_i) + \epsilon_i\), where \(f(x) = x^T w\) and \(\epsilon_i \sim \mathcal{N}(0, \sigma^2)\).

   • Noise is independent of the function value and other noise values: that is

     \begin{equation} \epsilon = \begin{bmatrix}\epsilon_1 \\ \vdots \\ \epsilon_n\end{bmatrix} \sim \mathcal{N}(0, \sigma^2 I_{n \times n}) \end{equation}
   • The goal of the regression is to find \(f(x)\), which in this case consists of finding the weights \(w \in \mathbb{R}^D\)

   • Without seeing any data, the weights are randomly distributed:

     \begin{equation} w \sim \mathcal{N}(0, \Sigma) \end{equation}

1. Write the posterior distribution of the weights using Bayes rule:

   \begin{equation} p(w | y, X) = \frac{p(y | X, w) p(w)}{p(y | X)} \end{equation}
   • We can apply Bayes rule directly, because we know \(p(w)\) and can compute \(p(y | X, w)\) from the formula for the model.
2. The likelihood is \(p(y | X, w) = \mathcal{N}(X^T w, \sigma^2 I_{n \times n})\).

   • This formula follows from the independence of each \(\epsilon_i\) and because \(y_i\) is a linear transform of a Gaussian random variable (\(\epsilon\)):

     \begin{align} \begin{bmatrix} y_1 \\ \vdots \\ y_n \end{bmatrix} &= \begin{bmatrix} x_1^T \\ \vdots \\ x_n^T \end{bmatrix} w + \begin{bmatrix} \epsilon_1 \\ \vdots \\ \epsilon_n \end{bmatrix}\\ y &= X^T w + \epsilon \end{align}

3. The normalization constant (marginal likelihood) can be computed as

   \begin{align} p(y | X) &= \int_{-\infty}^{\infty}p(y, w | x)dw \\ & = \int_{-\infty}^{\infty}p(y | x, w)p(w)dw \end{align}

4. After substituting the distributions into Bayes rule, we find

   \begin{equation} p(w | y, X) \sim \mathcal{N}(\bar{w}, A^{-1}), \end{equation}

   with

   \begin{equation} A = \frac{1}{\sigma^2}X X^T + \Sigma^{-1} \end{equation}

   and

   \begin{equation} \bar{w}= \frac{1}{\sigma^2}A^{-1} X y. \end{equation}

5. To make predictions from the model—that is, what is \(y_\ast = f(x_\ast)\) given \(x_\ast\)— use the probability distribution

   \begin{align} p(y_\ast | x_\ast, X, y) &= \int_{-\infty}^{\infty}p(y_\ast, w | x_\ast, X, y) dw \\ &=\int_{-\infty}^{\infty}p(y_\ast | x_\ast, w)p(w | X, y) dw\\ &= \mathcal{N}(\frac{1}{\sigma^2}x_\ast^T A^{-1} X y, x_\ast^T A^{-1}x_\ast)\\ &= \mathcal{N}(x_\ast^T \bar{w}, x_\ast^T A^{-1}x_\ast) \end{align}
   • The mean of the prediction is just the weighted sum of elements in the input vector, where the weights are the most probable weights
   • Variance grows with magnitude of input

• In its original form, \(p(f_\ast | x_\ast, X, y)\) requires inverting an \(N \times N\) matrix
  • With some manipulation we can make it so that we only need to invert an \(n \times n\) matrix.
• \(A \in \mathbb{R}^{N \times N}\), that is its dimensions scale with the number of features

• Let \(K = \Phi^T \Sigma \Phi\). Then note that

  \begin{equation} A^{-1} \Phi = \Sigma \Phi (K + \sigma^2 I)^{-1} \end{equation}

  and (from the Matrix Inversion Lemma a.k.a Woodbury Matrix Identity)

  \begin{equation} \phi_\ast^T A^{-1} \phi_\ast = \phi_\ast^T \Sigma \phi - \phi_\ast \Sigma \Phi(K + \sigma^2 I)^{-1}\Phi^T \Sigma \phi_\ast, \end{equation}

• After applying these transformations, the estimate becomes

  \begin{equation} \mathcal{N}(\phi_\ast^T\Sigma\Phi(K + \sigma^2I)^{-1}y, \phi_\ast^T\Sigma\phi_\ast - \phi_\ast^T\Sigma\Phi(K + \sigma^2)^{-1}\Phi^T\Sigma\phi_\ast \end{equation}
• Hand-wave thought about dimensionality:
  • We have an infinite number of basis functions \(N \to \infty\) with a finite number of inputs \(n\).
• The covariance kernel \(k(x, x') = \phi(x)^T \Sigma \phi(x')\).
• Note that \(k(x, x')\) is an inner product since \(\Sigma > 0\).
• The kernel trick:
  • A computational trick
  • If we know \(k(x, x')\), we can replace all these feature vector evaluations directly by using the kernel function on the original inputs.

Given the mean and covariance of a Gaussian process (and we'll assume mean 0 here), we can write a distribution for the process, partitioning into points we have measured and points at which we want a prediction

\begin{equation} \begin{bmatrix} f_\ast y \\ \end{bmatrix} \sim \mathcal{N}(0, \begin{bmatrix} K(X_\ast, X_\ast) & K(X_\ast, X) \\ K(X, X_\ast) & K(X, X) \end{bmatrix}), \end{equation}

where \(K(X, X_\ast)\) is the \(n \times n_\ast\) covariance matrix formed by \(k(x, x_\ast)\) at all training and test point combinations.

The above is just a jointly Gaussian distribution, with the prediction points and input points given. Thus, we can simply apply the formula for a conditional Gaussian to predict everything at \(f_\ast\).

\begin{equation} \mathcal{N}(K(X_\ast, X)(K(X, X) + \sigma^2 I)^{-1}y, K(X_\ast, X_\ast) - K(X_\ast, X)(K(X, X) + \sigma^2 I)^{-1}K(X, X_\ast)) \end{equation}