TL;DR
- ベイズ的な考えでニューラルネットワークを解きたいので、重み$w$の事前分布や条件付き確率$p(t|x)$に分布を導入する
- 分布を導入するとニューラルネットワークの非線形性の影響で数学的に解けなくなる。ラプラス近似などの近似を導入することで予測分布やハイパーパラメータの推定式を導出できる
5.7 ベイズニューラルネットワーク
先の節ではニューラルネットのパラメータを尤度が最大（誤差が最小）となるように、すなわち最尤推定で求めた。
ここからはベイズの枠組みで推定を行おう。すなわち、パラメータの事前分布を定義し、事後分布を求める。
3.3節で見たように、事後分布から予測分布（新しい$x$に対する$t$の予測分布）を求めるには重み$w$での積分演算が必要になる（⇒part10）が、多層ネットワークでは関数がパラメータに極度に非線形となるため、もはや数学的に解を求めることができない。
そこで、事後分布をそのモードを中心とするガウス分布で近似するラプラス近似（⇒part17）などの近似を用いながら推定を行っていく。
5.7.1 パラメータの事後分布
1次元の連続な目標変数$t$を入力ベクトル$x$から予測する問題を考えよう。
条件付き確率$ p(t|x) $はガウス分布と仮定し、その平均はニューラルネットの出力$ y(x,w) $で精度（分散の逆数）は$\beta$とする。すなわち、
$$
p(t| x,w,\beta ) = N(t | y(x,w), \beta^{-1} ) \tag{5.161}
$$
また、重み$w$に関する事前分布もガウス分布を仮定し：
$$
p(w | \alpha) = N(w | 0, \alpha^{-1} I ) = \left( \frac{ \alpha }{ 2 \pi} \right)^{M/2} \exp
\left\{ - \frac{ \alpha }{2} w^T w \right\}
\tag{5.162}
$$
とする。
独立同分布による$N$回の観測値$ x_1, \cdots, x_N $をデータ集合とし、それぞれの目標値集合を$ D=\{t_1, \cdots, t_N \} $とすれば、尤度関数は
$$
p(D|w, \beta) = \prod_{n=1}^N N(t_n | y(x_n, w), \beta^{-1} )
= \prod_{n=1}^N \left( \frac{\beta}{2 \pi} \right)^{1/2} \exp
\left\{
- \frac{\beta}{2} (t_n - y(x_n, w) )
\right\}
\tag{5.163}
$$
で与えられ、尤度関数の事後分布は
\begin{align}
p(w | D, \alpha, \beta) &\propto p(w|\alpha) p(D|w, \beta) \tag{5.164} \\
&\propto \exp \left\{
- \prod_{n=1}^N \left[
\frac{\beta}{2} (t_n - y(x_n, w) )
\right] - \frac{\alpha}{2} w^T w
\right\}
\end{align}
となる。expの中身が$w$に関する二次式であればガウス分布である（⇒part4）が、$y(x_n,w)$が$w$に非線形に依存し二次式にならないため、式(5.164)の事後分布はガウス分布にならない。
そこで、ラプラス近似（⇒part17）を用いて事後分布をガウス分布で近似しよう。
そのためには、まず事後分布の（局所的）最大値$w_{MAP}$を見つける必要があり、これは反復的数値最適化法が必要らしい。
下式の事後分布の対数の$\alpha$と$\beta$を固定とみなして、共役勾配法などのアルゴリズムを使えば、事後分布の最大値$w_{MAP}$を見つけることができる（らしい。共役勾配法の解説はテキストにないので、そういうものだと思おう）。
$$
\ln p(w|D) = - \frac{\alpha}{2} w^T w - \frac{\beta}{2} \sum_{n=1}^N \{ y(x_n, w) - t_n \}^2 + const. \tag{5.165}
$$
近似ガウス分布の精度行列$A$は過去の議論より事後分布の負の対数尤度の2回微分の行列により与えられ、
$$
A = - \nabla \nabla \ln p(w|D, \alpha, \beta) = \alpha I + \beta H \tag{5.166}
$$
である。ただし、$H$は二乗和誤差関数$ \{ y(x_n, w) - t_n \}^2 $の$w$の各要素による2階微分からなるヘッセ行列である。
これを利用して事後分布を近似したガウス分布は
$$
q(w|D) = N(w| w_{MAP}, A^{-1} ) \tag{5.167}
$$
で与えられる。この事後分布を$w$に関して周辺化することで予測分布が得られ、
$$
p(t|x, D) = \int p(t|x, w) q(w|D) dw \tag{5.168}
$$
となる。これでめでたしといきたいが、事後分布をガウス分布で近似しても$ p(t|x, w) $に$w$の非線形関数である$y(w,x)$が含まれている（ガウス分布$ p(t|x, w) $の平均$y(w,x)$が非線形に変化する）ことから、この積分は解析的に扱いづらい。
そこで、$y(w,x)$を$ w_{MAP} $のまわりでテイラー展開して線形項だけで表す1と
\begin{align}
y(x,w) \simeq y(x, w_{MAP}) + \mathbf{g}^T (w - w_{MAP} ) \tag{5.169} \\
\mathbf{g} = \nabla_w y(x,w) |_{w=w_{MAP} } \tag{5.170}
\end{align}
となる。よって：
\begin{align}
p(t| x,w,\beta ) &= N(t | y(x,w), \beta^{-1} ) \tag{5.161} \\
& \simeq N(t | y(x, w_{MAP}) + \mathbf{g}^T (w - w_{MAP} ), \beta^{-1} ) \tag{5.171}
\end{align}
式(5.171)、式(5.167)の結果を使って予測分布を得るための積分（式(5.168)）は2.3.3節で導出した結果を用いると簡単に計算できる（詳細は後ほど：演習5.38）：
\begin{align}
p(t | x, D, \alpha, \beta) &= N(t | y(x, w_{MAP}), \sigma^2(x) ) \tag{5.172} \\
\sigma^2(x) &= \beta^{-1} + \mathbf{g}^T A^{-1} \mathbf{g} \tag{5.173}
\end{align}
この式より、予測分布は平均がネットワーク関数$ y(x, w_{MAP}) $で与えられるガウス分布であることが分かる。
分散は2つの項を持ち、1つめは目標変数における内在的なノイズに起因する。2つめは$x$に依存する項で、$w$の不確定性を表している2。
演習5.38：式(5.172-5.173)の導出
2.3.3節（part5）で導出したガウス分布の周辺分布に関する公式を用いればよい。
\begin{align}
p(x) &= N(x| \mu, \Lambda^{-1} ) \tag{2.113} \\
p(y|x) &= N(y | Ax + b, L^{-1} ) \tag{2.114} \\
p(y) &= N(y | A \mu + b, L^{-1} + A \Lambda^{-1} A^T )\tag{2.115}
\end{align}
- 求めたい予測分布$ p(t | x, D, \alpha, \beta) $を$p(y)$
- $q(w|D)$を$ p(x) $
- $p(t|x, w)$を$ p(y|x) $
と考えて、式(2.113-2.114)の$\mu, \Lambda^{-1}, A, x, b, L^{-1}$を読み替えれば計算できる。
- $x$は$w$, $y$は$t$
- $\mu$は$w_{MAP}$, $\Lambda^{-1}$は$A^{-1}$
- $L^{-1}$は$\beta^{-1}$
- $x$が$w$なことを踏まえると$A$は$\mathbf{g}^T $
- $b$は$ y(x, w_{MAP} ) - \mathbf{g}^T w_{MAP} $となるので、予測分布の$ A \mu + b $は$ \mathbf{g}^T w_{MAP} + y(x, w_{MAP} ) - \mathbf{g}^T w_{MAP} = y(x, w_{MAP} ) $となる
5.7.2 超パラメータ最適化
ここまでハイパーパラメータである$ \alpha, \beta $は固定で既知のものと扱ってきた。ここからは3.5節のエビデンス理論（⇒part12）を使ってハイパーパラメータの値を選ぶ実用的な手続きを得る、というのが本節の目的。
数学的な指針は3.5節と同じ、かつ式変形が大変なので簡単にだけ見ておこう。
ハイパーパラメータの周辺尤度（エビデンスともいう）はネットワークの重み$w$に関する積分で得られ、
$$
p(D | \alpha, \beta) = \int p(D | w, \beta) p(w | \alpha) dw \tag{5.174}
$$
である。これの対数は
\begin{align}
\ln p(D | \alpha, \beta) &\simeq -E(w_{MAP}) - \frac{1}{2} \ln |A| + \frac{W}{2} \ln \alpha + \frac{N}{2} \ln \beta - \frac{N}{2} \ln (2 \pi) \tag{5.175} \\
E(w_{MAP}) &= \frac{\beta}{2} \sum_{n=1}^N \{ y(x_n, w_{MAP}) - t_n \}^2 + \frac{\alpha}{2} w_{MAP}^T w_{MAP} \tag{5.176}
\end{align}
ここで$W$は$w$に含まれるパラメータの総数である。
エビデンス理論においては$ \ln p(D | \alpha, \beta) $を最大化して$\alpha, \beta$の点推定を行う。
推定式は
\begin{align}
\alpha &= \frac{ \gamma }{w_{MAP}^T w_{MAP}} \tag{5.178} \\
\gamma &= \sum_{i=1}^W \frac{\lambda_i}{ \alpha + \lambda_i } \tag{5.179} \\
\frac{1}{ \beta } &= \frac{1}{N - \gamma} \sum_{n=1}^N \{ y(x_n, w_{MAP} ) - t_n \}^2 \tag{5.180}
\end{align}
ここで$ \lambda_i $は誤差関数の2階微分からなるヘッセ行列の固有値である。
$ \alpha, \beta $, 事後分布は互いに依存しているため、$ \alpha, \beta $の再推定と事後分布の更新は交互に行わなければならない。
5.7.3 クラス分類のためのベイズニューラルネットワーク
これまでニューラルネットワーク回帰モデルをベイズの枠組みで解いてきた。
ここではクラス分類に適用できるようにこの枠組みを修正する。
1つのロジスティックシグモイド出力を持つネットワークによる2クラス分類問題を考える。
この節も数学的な指針は4.5節と同じ、かつ式変形が大変なので簡単にだけ見ておこう。
このモデルの対数尤度関数は
$$
\ln p(D|w) = \sum_{n=1}^N \{ t_n \ln y_n + (1 - t_n) \ln (1 - y_n) \} \tag{5.181}
$$
で与えられる。ここで、$ t_n \in \{0,1\} $は目標値であり、$ y_n = y(x_n, w) $である。
まずハイパーパラメータ$\alpha$を適当な値に初期化し、対数事後分布を最大化して$w$を決める。
これは正則化誤差関数
$$
E(w) = - \ln p(D|w) + \frac{\alpha}{2} w^T w \tag{5.182}
$$
を最小化することと等価であり、誤差逆伝播法で求めることができる。求まった重みを$w_{MAP}$とする。
ハイパーパラメータ$ \alpha $を最適化するには周辺尤度を最大化する。周辺尤度は式変形により
\begin{align}
\ln p(D|\alpha) &\simeq -E(w_{MAP}) - \frac{1}{2} \ln |A| + \frac{W}{2} \ln \alpha \tag{5.183} \\
E(w_{MAP}) &= - \sum_{n=1}^N \{ t_n \ln y_n + (1 - t_n) \ln (1 - y_n) \} + \frac{\alpha}{2} w_{MAP}^T w_{MAP} \tag{5.184}
\end{align}
次に予測分布について。
予測分布を得るための積分がネットワーク関数の非線形性により計算困難であること、ロジスティック回帰モデルでは出力が(0,1)の範囲に制限されるため、回帰の場合のようにネットワーク出力を線形近似するのは不適切であること、ロジスティックシグモイドを持つガウス分布の畳み込み計算は困難であることから、予測分布の導出にはいくつかの近似が必要となる。