ガウス過程回帰モデル

カーネル法を用いた事後予測分布を、さらに別の角度から解釈する。

関数 $f: \mathcal X \to \R$ について、その出力 $f_i = f(\bm x_i)$ を任意の個数だけ纏めたベクトル $\bm f \in \R^d$ が多変量正規分布

\bm f \sim \mathcal N_d(\bm f | \bm m, \bm K)

に従い、その平均ベクトル $\bm m$ および共分散行列 $\bm K$ がそれぞれ

\begin{aligned} \bm m ={}& \left[\begin{darray}{c} m(\bm x_i) \end{darray}\right]_{i} && \in \R^d, \\ \bm K ={}& \left[\begin{darray}{c} k(\bm x_i, \bm x_j) \end{darray}\right]_{ij} && \in \R^{d \times d} \end{aligned}

で表されるとする。このようなとき、関数 $f$ はガウス過程 (Gaussian process) に従うといい、

f \sim \mathcal{GP}(m, k)

と表す。 $m$ を平均関数 (mean function)、 $k$ を分散関数 (variance function) などという。

この特別な場合として、平均関数が $0$ である場合を考えよう。

\begin{aligned} &&& f \sim \mathcal{GP}(0, k) \\ &\text{i.e.}&& \bm f \sim \mathcal N_d(\bm y | \bm 0, \bm V_y) \end{aligned}

関数 $f$ の任意個数の出力がこの分布に従うということから、 $\bm f$ に更に $d_\ast$ 個の出力値 $\bm f_\ast$ を追加しても、分布は同じ形になる。

p(\bm f_\ast, \bm f | \bm X_\ast, \bm X) = \mathcal N\left( \left[\begin{darray}{c} \bm f \\ \bm f_\ast \end{darray}\right] \middle| \left[\begin{darray}{c} \bm 0 \\ \bm 0 \end{darray}\right], \left[\begin{darray}{c} \bm K & \bm K_\ast \\ \bm K_\ast^\mathsf{T} & \bm K_{\ast\ast} \end{darray}\right] \right)

このとき、条件付き正規分布の性質により、 $\bm f$ が得られたという条件のもとで $\bm f_\ast$ が従う分布は次のように書ける。

\begin{aligned} p(\bm f_\ast | \bm X_\ast, \bm X, \bm f) = \mathcal N(\bm y_\ast | \bm m'_{f_\ast | f}, \bm V'_{f_\ast | f}) \end{aligned}

\begin{aligned} \bm m_{f_\ast | f} ={}& \bm K_\ast^\mathsf{T} \bm K^{-1} \bm y \\ \bm V_{f_\ast | f} ={}& \bm K_{\ast\ast} - \bm K_\ast^\mathsf{T} \bm K^{-1} \bm K_\ast \\ \end{aligned}

ただし行列 $\bm K, \bm K_\ast, \bm K_{\ast \ast}$ は次のとおりである。

\begin{aligned} \bm K ={}& \left[\begin{darray}{c} k_{ij} \end{darray}\right]_{ij} && \in \R^{d \times d} \\ \bm K_{\ast} ={}& \left[\begin{darray}{c} k_{ij^\ast} \end{darray}\right]_{ij} && \in \R^{d \times d^\ast} \\ \bm K_{\ast\ast} ={}& \left[\begin{darray}{c} k_{i^\ast j^\ast} \end{darray}\right]_{ij} && \in \R^{d^{\ast} \times d^{\ast}} \end{aligned}

k_{ij} = k(\bm x_i, \bm x_j)

観測値 $\bm y$ は、関数の出力値 $\bm f$ に観測ノイズ $\bm \varepsilon \sim \mathcal N_d(\bm \varepsilon | \bm 0, \bm \sigma^2 \bm I_d)$ を追加したものと考える。すると次のように書ける。

p(\bm y_\ast, \bm y | \bm X_\ast, \bm X) = \mathcal N\left( \left[\begin{darray}{c} \bm y \\ \bm y_\ast \end{darray}\right] \middle| \left[\begin{darray}{c} \bm 0 \\ \bm 0 \end{darray}\right], \left[\begin{darray}{c} \bm K & \bm K_\ast \\ \bm K_\ast^\mathsf{T} & \bm K_{\ast\ast} \end{darray}\right] + \sigma^2 \bm I_{d + d^\ast} \right)

このとき、 $\bm y$ が得られたのち $\bm y_\ast$ が従う分布は

\begin{aligned} p(\bm y_\ast | \bm X_\ast, \bm X, \bm y) = \mathcal N(\bm y_\ast | \bm m'_{y_\ast | y}, \bm V'_{y_\ast | y}) \end{aligned}

\begin{aligned} \bm m_{y_\ast | y} ={}& \bm K_\ast^\mathsf{T} (\bm K + \sigma^2 \bm I_d)^{-1} \bm y \\ \bm V_{y_\ast | y} ={}& \bm K_{\ast\ast} + \sigma^2 \bm I_{d^\ast} - \bm K_\ast^\mathsf{T} (\bm K + \sigma^2 \bm I_d)^{-1} \bm K_\ast \\ \end{aligned}

と書ける。

面白いことに、このような設定で計算された $\bm m_{y_\ast|y}$ と $\bm V_{y_\ast|y}$ はそれぞれ、 📄 事後予測分布とWoodburyの公式で計算された $\bm m_{y_\ast | y}$ と $\bm V_{y_\ast|y}$ とまったく同じ形となる。

p(\bm y_\ast | \bm X_\ast, \bm X, \bm y) = \mathcal N_{d^\ast} (\bm y_\ast | \bm m_y, \bm V_y)

\begin{aligned} \bm m_{y_\ast | y} ={}& \bm K_\ast^\mathsf{T} (\bm K + \sigma^2 \bm I_d)^{-1} \bm y \\ \bm V_{y_\ast | y} ={}& \bm K_{\ast\ast} + \sigma^2 \bm I_{d^\ast} - \bm K_\ast^\mathsf{T} (\bm K + \sigma^2 \bm I_d)^{-1} \bm K_\ast \\ \end{aligned}

\left\{\begin{aligned} \bm K ={}& \bm X \bm V_0 \bm X^\mathsf{T} \\ \bm K_{\ast} ={}& \bm X \bm V_0 \bm X_\ast^\mathsf{T} \\ \bm K_{\ast\ast} ={}& \bm X_\ast \bm V_0 \bm X_\ast^\mathsf{T} \\ \end{aligned}\right.

この理屈から、ガウス過程を用いた回帰モデルはきわめて表現力の高いモデルとして注目されており、ガウス過程回帰モデル (Gaussian process regressor) と呼ばれている。

ガウス過程回帰モデルは、欠損データの内挿のための回帰モデルや、ベイズ最適化の代理モデルに用いられることが多い。