prml5.tex

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 第5章
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\setcounter{chapter}{4}
\chapter{「ニューラルネットワーク」の補足}
PRML5章では, 本来別の記号を割り当てるべきところに同じ記号を用いることがあり, 初読時には混乱しやすい.
慣れてしまえば読むのは難しくないが, ここでは出来るだけ区別してみる.

\newcommand{\ww}[2]{w_{#1}^{(#2)}}
\newcommand{\calRR}[1]{{\cal R}\left\{#1 \right\}}
\newcommand{\deltaj}{\delta_j^{(1)}}
\newcommand{\deltak}{\delta_k^{(2)}}
\newcommand{\wmap}{w_{\rm MAP}}
\newcommand{\amap}{a_{\rm MAP}}

\section{フィードフォワードネットワーク関数}
3章, 4章でやったモデルは基底関数$\phi_j$とパラメータ$w_j$の線形和を非線形活性化関数に入れたものだった.
ここではそれを拡張する.
$x_1, \ldots, x_D$を入力変数とし$x_0=1$をバイアス項（定数項）に対応する変数, $w_{ji}^{(1)}$をパラメータとして
$$
\hat{a}_j = \sum_{i=0}^D \ww{ji}{1} x_i
$$
とする.
PRMLでは上記の$\hat{a}_j$を$a_j$と書いているがすぐあとに出てくる$a_k$とは無関係である.
ここでは異なることを強調するために$\hat{a}_j$とする.

$\hat{a}_j$を活性化関数$h$で変換する.
$$
z_j=h(\hat{a}_j).
$$
$h$としてはロジスティックシグモイドなどのシグモイド関数が用いられる.
これらの線形和をとって出力ユニット活性を求める. $z_0=1$をバイアス項に対応する変数として
$$
a_k=\sum_{j=0}^M \ww{kj}{2} z_j.
$$

この出力ユニット活性を活性化関数を通してネットワークの出力$y_k$とする.
2クラス分類問題ならロジスティックシグモイド関数を使う.
$$
y_k=y_k(x,w)=\sigma(a_k).
$$
ここで$w$は$\{\ww{ji}{1}, \ww{kj}{2}\}$をまとめたベクトルである.
これらの式を組み合わせると
$$
y_k=\sigma\left(\sum_{j=0}^M w_{kj}^{(2)} h \left(\sum_{i=0}^D w_{ji}^{(1)} x_i \right) \right).
$$
\section{ネットワーク訓練}\label{ch5_net}
回帰問題を考える. 入力ベクトル$x$と$K$次元の目標変数$t$があり,
$x$, $w$における$t$の条件付き確率が精度$\beta I$のガウス分布とする.
$N$個の同時独立分布$x=\{x_1,\ldots,x_N\}$と対応する目標値$t=\{t_1,\ldots,t_N\}$を用意し, 出力ユニットの活性化関数を恒等写像として
$y_n=y(x_n,w)$とする.
$$
p(t|x,w)=\prod_{n=1}^N \calN(t_n|y(x_n,w),\beta^{-1}I).
$$
対数をとると
\begin{eqnarray}\label{cp5_1}
\log p(t|x,w) &=& -\sum_n \left(\frac{1}{2}\trans{(t_n-y_n)}(\beta I)(t_n-y_n)\right) - \sum_n \frac{K}{2}\log(2\pi) - \sum_n \frac{1}{2}\log|\beta^{-1}I| \nonumber \\
&=&-\frac{\beta}{2}\sum_n||t_n-y(x_n,w)||^2-\frac{DN}{2}\log(2\pi)+\frac{NK}{2}\log \beta.
\end{eqnarray}
そうするとこの関数の$w$についての最大化は最初の項の最小化, つまり
$$
E(w)=\frac{1}{2}\sum_n ||y(x_n,w)-t_n||^2
$$
の最小化と同等である. 最小値を与える$w(=w_{\text ML})$をなんらかの方法で求める.
その値を式(\ref{cp5_1})に代入して$\beta$で微分して0とおくと
$$
-\frac{1}{2}\sum_n||t_n-y(x_n,w_{\text ML})||^2+\frac{NK}{2}\frac{1}{\beta}=0.
$$
よって
$$
\frac{1}{\beta_{\text ML}}=\frac{1}{NK}\sum_n||t_n-y(x_n,w_{\text ML})||^2.
$$

\subsection{問題に応じた関数の選択}
回帰問題を考える. $a_k$の活性化関数を恒等写像にとる.
すると二乗和誤差関数の微分は
$$
\diff{a_k}{E}=y_k-t_k.
$$
クラス分類問題でも同様の関係式が成り立つことを確認しよう.

目標変数$t$が$t=1$でクラス$C_1$, $t=0$でクラス$C_2$を表す2クラス分類問題を考える.
活性化関数をロジスティックシグモイド関数に選ぶ.
$$
y=\sigma(a)=\frac{1}{1+\exp(-a)}.
$$
この微分は$dy/da=\sigma(a)(1-\sigma(a))=y(1-y)$であった.
$p(t=1|x)=y(x,w)$, $p(t=0|x)=1-y(x,w)$なので
$$
p(t|x,w)=y(x,w)^t \left(1-y(x,w)\right)^{1-t}.
$$
よって4章と同様にして交差エントロピー誤差関数は
$$
E(w)=-\sum_{n=1}^N \left(t_n \log y_n + (1-t_n) \log (1-y_n)\right).
$$
これを$a_k$で微分すると
\begin{eqnarray*}
\diff{a_k}{E} &=& -\left(t_k\frac{y_k(1-y_k)}{y_k} + (1-t_k)\frac{-y_k(1-y_k)}{1-y_k}\right)\\
              &=& -(t_k - t_k y_k - y_k + y_k t_k)\\
              &=& y_k - t_k.
\end{eqnarray*}
$K$個の2クラス分類問題を考える. それぞれの活性化関数がロジスティックシグモイド関数とする.
$$
p(t|x,w)=\prod_{k=1}^K y_k(x,w)^{t_k} (1-y_k(x,w))^{1-t_k}.
$$
$y_{nk}=y_k(x_n,w)$とし$n$番目の入力$x_n$に対する目標変数を$t_{nk}$で表す.
$t_{nk}\in\{0,1\}$であり, $\sum_k t_{nk}=1$である.
$$
E(w)=-\prod_{n,k} (t_{nk} \log y_{nk}+(1-t_{nk})\log (1-y_{nk})).
$$
$y_{nj}$に対応する$a$を$a_{nj}$とすると
$$
\diff{a_{nj}}{E(w)}=-(t_{nj}(1-y_{nj})+(1-t_{nj})(-y_{nj}))=y_{nj}-t_{nj}.
$$

最後に$K$クラス分類問題を考える. 同様に$n$番目の入力$x_n$に対する目標変数$t_{nk}$で表す.
$y_k(x_n,w)$を$t_{nk}$が1となる確率$p(t_{nk}=1|x_n)$とみなす.
$$
E(w) = -\log p(t|x,w)=-\sum_{n,k} t_{nk} \log y_k(x_n,w).
$$
活性化関数はソフトマックス関数で
$$
y_k(x,w)=\frac{\exp a_k(x,w)}{\sum_j \exp(a_j(x,w))}
$$
のとき
$$
\diff{a_j}{y_k}=y_k(\delta_{kj}-y_j).
$$
よって$a_{nj}=a_j(x_n,w)$とすると
$$
\dif{a_{nj}}\log y_k(x_n,w)=\frac{1}{y_{nk}}\diff{a_{nj}}{y_{nk}}=\delta_{kj}-y_{nj}.
$$
よって
$$
\diff{a_{nj}}{E}=-\sum_k t_{nk}(\delta_{kj}-y_{nj})=-t_{nj}+\left(\sum_k t_{nk}\right)y_{nj}=y_{nj}-t_{nj}.
$$
\section{局所二次近似}\label{ch5_loc}
スカラー$x$についての関数$E(x)$の$x=a$におけるテイラー展開を2次の項で打ち切った近似式は
$$
E(x) \approx E(a) + E'(a)(x-a) + \half E''(a) (x-a)^2
$$
であった. これを$n$変数関数$E(x_1, \ldots, x_n)$に拡張する.
$x=\trans{(x_1, \ldots, x_n)}$, $a=\trans{(a_1, \ldots, a_n)}$とおいて
\begin{eqnarray*}
E(x)
 &\approx& E(a) + \trans{(x-a)}\left(\dif{x_i}E\right) + \half \trans{(x-a)}\left(\ddiff{x_i}{x_j}{E}\right)(x-a)\\
 &=&       E(a) + \trans{(x-a)} (\nabla E) + \half \quads{H(E)}{(x-a)}
\end{eqnarray*}
となる.
$E(x)$が$x=a$の付近で極小ならば, そこでの勾配$\nabla E$は0なので
$$
E(x) \approx E(a) + \half \quads{H(E)}{(x-a)}.
$$
$H(E)$は対称行列なので\ref{pos_sym_matrix}節の議論より対角化することで
$$
E(x) \approx E(a) + \half \sum_i \lambda_i y_i^2
$$
の形にできる. そして$E(x)$が$x=a$の付近で極小となるのは$H(E)>0$（正定値）であるときとわかる.

なお, $H(f)=\nabla^2 f = \nabla(\nabla f)$という表記をすることがある. 微分作用素$\nabla$を2回するので2乗の形をしている.
ただ$\nabla f$が縦ベクトルならもう一度$\nabla$をするときは結果が行列になるように,
入力ベクトルの転置を取って作用するとみなす.
$n^2$次元の長いベクトルになるわけではない.
\section{誤差関数微分の評価}
与えられたネットワークに対して誤差関数の変化の割合を調べる.
この節ではどの変数がどの変数に依存しているか気をつけて微分する必要がある.
誤差関数が訓練集合の各データに対する誤差の和で表せると仮定する：
$$
E(w)=\sum_{n=1}^N E_n(w).
$$
一般のフィードフォワードネットワークで
\begin{eqnarray}\label{ch5_2}
a_j=\sum_i w_{ji} z_i, \quad z_j=h(a_j)
\end{eqnarray}
とする. 入力$z_i$が出力ユニット$a_j$に影響を与え, その$a_j$が非線形活性化関数$h()$を通して$z_j$に影響を与える.
ある特定のパターン$E_n$の重み$w_{ji}$に関する微分を考える.
以下, 特定のパターンを固定することで$E_n$以外の添え字の$n$を省略する.
式(\ref{ch5_2})のように$E_n$は非線形活性化関数$h$の変数$a_j$を通して$w_{ji}$に依存している.
$$
\diff{w_{ji}}{E_n}=\diff{a_j}{E_n}\diff{w_{ji}}{a_j}.
$$
$a_j$は$w_{ji}$に関しては線形なので
$$
\diff{w_{ji}}{a_j}=z_i.
$$
誤差と呼ばれる記号$\delta_j = \partial E_n / \partial a_j$を導入すると
$$
\diff{w_{ji}}{E_n}=\delta_j z_i
$$
とかける. $\delta_j$は$h()$が正準連結関数の場合は\ref{ch5_net}節での考察により
$$
\delta_j=\diff{a_j}{E_n}=y_j-t_j
$$
で計算できる.
ユニット$j$につながっているユニット$k$を通して$E_n$への$a_j$の影響があると考えると
\begin{eqnarray}\label{cp5_3}
\diff{a_j}{E_n}=\sum_k \diff{a_k}{E_n}\diff{a_j}{a_k}.
\end{eqnarray}
ここで$a_j$と$a_k$は$z$を経由して関係していると考えているので$\partial a_k / \partial a_j=\delta_{jk}$（クロネッカーのデルタ）にはならないことに注意する.
実際,
$$
a_k=\sum_i w_{ki} h(a_i)
$$
より
$$
\diff{a_j}{a_k}=w_{kj} h'(a_j).
$$
これを式(\ref{cp5_3})に代入して
$$
\delta_j=\diff{a_j}{E_n}=\sum_k \delta_k w_{kj} h'(a_j)=h'(a_j)\sum_k w_{kj} \delta_k.
$$
\section{外積による近似}
$$
E(w)=\half \sum_{n=1}^N \left(y_n-t_n\right)^2
$$
のときのヘッセ行列は
$$
H=H(E) = \sum_n \outp{(\nabla y_n)} + \sum_n (y_n-t_n)H(y_n).
$$
一般には成立しないがもしよく訓練された状態で$y_n$が目標値$t_n$に十分近ければ第2項を無視できる.
その場合$b_n=\nabla y_n=\nabla a_n$（活性化関数が恒等写像なので）とおくと
$$
H \approx \sum_n \outp{b_n}.
$$
これをLevenberg-Marquardt近似という.
$$
E=\half \iint (y(x,w)-t)^2p(x,t)\,dxdt
$$
のときのヘッセ行列を考えてみると
$$
\diff{w_i}{E}=\iint (y-t)\diff{w_i}{y} p(x,t)\,dxdt.
$$
\begin{eqnarray*}
\ddiff{w_i}{w_j}{E}&=&\iint \left( \diff{w_j}{y}\diff{w_i}{y} + (y-t)\ddiff{w_i}{w_j}{y} \right)p(x,t)\,dxdt.\\
 && （p(x,t)=p(t|x)p(x)より）\\
 &=&\int \diff{w_j}{y}\diff{w_i}{y}\left( \int p(t|x)\,dt\right) p(x)\,dx\\
 &+& \int \ddiff{w_i}{w_j}{y}\left(\int (y-t)p(t|x)\,dt\right)p(x)\,dx\\
 && （第1項のカッコ内は1. \quad 第2項はy(x)=\int t p(t|x)\,dtを使うと0）\\
 &=&\int \diff{w_j}{y}\diff{w_i}{y}p(x)\,dx.
\end{eqnarray*}
ロジスティックシグモイドのときは
\begin{eqnarray*}
\nabla E(w)&=&\sum_n \diff{a_n}{E} \nabla a_n=-\sum_n\left(\frac{t_n y_n(1-y_n)}{y_n}-\frac{(1-t_n)y_n (1-y_n)}{1-y_n}\right)\nabla a_n\\
&=& \sum_n (y_n-t_n)\nabla a_n.
\end{eqnarray*}
よって$y_n \approx t_n$なら
\begin{eqnarray*}
\nabla^2 E(w) &=& \sum_n \diff{a_n}{y_n} \outp{\nabla a_n}+\sum_n (y_n-t_n)\nabla^2 a_n\\
 &\approx& \sum_n y_n(1-y_n) \outp{\nabla a_n}.
\end{eqnarray*}

\section{ヘッセ行列の厳密な評価}
この節は計算は難しくはないが, 記号がややこしいので書いてみる.
変数の関係式は
$\hat{a}_j=\sum_i \ww{ji}{1}x_i$, $z_j=h(\hat{a}_j)$, $a_k=\sum_j \ww{kj}{2} z_j$, $y_k=a_k$である.
PRMLの$a_j$と$a_k$は違う対象であることに注意する. ここでは$a_j$の代わりに$\hat{a}_j$を使う.

添え字の$i$, $i'$は入力, $j$, $j'$は隠れユニット, $k$, $k'$は出力である.
また
$$
\delta_k=\diff{a_k}{E_n}, \quad M_{kk'}=\ddiff{a_k}{a_{k'}}{E_n}
$$
という記号を導入する.
$E_n$以外の添え字$n$を省略する.
\subsection{両方の重みが第2層にある}
$$
\diff{\ww{kj}{2}}{a_k}=z_j, \quad \diff{\ww{kj}{2}}{E_n}=\diff{\ww{kj}{2}}{a_k}\diff{a_k}{E_n}=z_j \delta_k.
$$
よって
$$
\ddiff{\ww{kj}{2}}{\ww{k'j'}{2}}{E_n}=\diff{\ww{k'j'}{2}}{a_k'}\dif{a_{k'}}\left(\diff{\ww{kj}{2}}{E_n}\right)=z_{j'}z_j \diff{a_{k'}}{\delta_k}=z_j z_{j'} M_{kk'}.
$$
\subsection{両方の重みが第1層にある}
$$
\diff{\hat{a}_j}{a_k}=\ww{kj}{2}h'(\hat{a}_j), \quad \diff{\ww{ji}{1}}{\hat{a}_j}=x_i
$$
より
$$
\diff{\ww{ji}{1}}{E_n}=\diff{\ww{ji}{1}}{\hat{a}_j}\diff{\hat{a}_j}{E_n}=x_i\sum_k \diff{\hat{a}_j}{a_k}\diff{a_k}{E_n}=x_i\sum_k \ww{kj}{2}h'(\hat{a}_j)\delta_k.
$$
$$
\ddiff{\ww{ji}{1}}{\ww{j'i'}{1}}{E_n}=\diff{\ww{j'i'}{1}}{\hat{a}_{j'}}\dif{\hat{a}_{j'}}\left(\diff{\ww{ji}{1}}{E_n}\right)=x_i x_{i'}\underbrace{\dif{\hat{a}_{j'}}\left(h'(\hat{a}_j)\sum_k \ww{kj}{2}\delta_k \right)}_{=:A}.
$$
$j=j'$のとき
$$
A=h''(\hat{a}_{j'})\sum_k \ww{kj}{2}\delta_k + B, \quad B:=h'(\hat{a}_j)\dif{\hat{a}_{j'}}\left(\sum_k \ww{kj}{2}\delta_k\right).
$$
$j\ne j'$のとき
$$
B= h'(\hat{a}_j) \sum_{k'} \diff{\hat{a}_j}{a_{k'}}\dif{a_{k'}}\left(\sum_k \ww{kj}{2}\delta_k\right)
 = \sum_{k,k'} h'(\hat{a}_j)h'(\hat{a}_{j'})\ww{k'j'}{2}\ww{kj}{2}M_{kk'}.
$$
二つをまとめて
$$
\ddiff{\ww{ji}{1}}{\ww{j'i'}{1}}{E_n}=x_i x_{i'}\left\{h''(\hat{a}_{j'})\delta_{jj'}\sum_k \ww{kj}{2}\delta_k+h'(\hat{a}_j)h'(\hat{a}_{j'})\sum_{k,k'}\ww{kj}{2}\ww{k'j'}{2}M_{kk'}\right\}.
$$
$\delta_{jj'}$はクロネッカーのデルタ.
\subsection{重みが別々の層に一つずつある}
$$
\diff{\ww{kj'}{2}}{E_n}=z_{j'}\delta_k, \quad \ddiff{\ww{ji}{1}}{\ww{kj'}{2}}{E_n}=\diff{\ww{ji}{1}}{\hat{a}_j} \underbrace{\dif{\hat{a}_j}(z_{j'}\delta_k)}_{=:A}, \quad \diff{\ww{ji}{1}}{\hat{a}_j}=x_i.
$$
$j=j'$のとき
$$
A=h'(\hat{a}_{j'})\delta_k+B, \quad B:=z_{j'}\diff{\hat{a}_j}{\delta_k}.
$$
$j\ne j'$のとき
$$
B=z_{j'}\sum_{k'} \diff{\hat{a}_j}{a_{k'}}\diff{a_{k'}}{\delta_k}=z_{j'}\sum_{k'} \ww{k'j}{2}h'(\hat{a}_j)M_{kk'}.
$$
よって
$$
\ddiff{\ww{ji}{1}}{\ww{kj'}{2}}{E_n}=x_i h'(\hat{a}_j)\left\{\delta_{jj'}\delta_k+z_{j'}\sum_{k'} \ww{k'j}{2}M_{kk'}\right\}.
$$

\section{ヘッセ行列の積の高速な計算}
応用面を考えると最終的に必要なものはヘッセ行列$H$そのものではなくあるベクトル$v$と$H$の積であることが多い.
$H$を計算せず直接$\trans{v}H=\trans{v}\nabla \nabla$を計算するために, 左半分だけを取り出して$\calRR{\cdot}=\trans{v}\nabla$という記法を導入する.
\ref{ch5_loc}節の終わりに書いたようにこの$\nabla$は入力が縦ベクトルなら転置を取ってから作用するとみなす.
なお, $v$に依存するものをあたかも依存しないかのように$\calRR{\cdot}$と書いてしまうのは筋がよいとは思わない.

簡単な例を見てみよう. 2変数関数$y=f(x_1, x_2)$について
$$
\calRR{\cdot}=(v_1,v_2)\nabla = (v_1,v_2)\vvec{\dif{x_1}}{\dif{x_2}}.
$$
よって
\begin{eqnarray*}
\calRR{x_1}&=&(v_1,v_2)\vvec{1}{0}=v_1,\\
\calRR{x_2}&=&(v_1,v_2)\vvec{0}{1}=v_2,
\end{eqnarray*}
これを, $\calRR{}$は入力値の$x_i$をその添え字に対応する$v_i$に置き換える作用と考えることにする.
$\calRR{}$は$x_i$について明らかに線形, つまり
$$
\calRR{ax_1+bx_2}=a v_1 + b v_2 = a\calRR{x_1}+b\calRR{x_2}.
$$

前節と同じ2層ネットワークで考えてみる.
$a_j=\sum_i \ww{ji}{1}x_i$, $z_j=h(a_j)$, $y_k=\sum_j \ww{kj}{2} z_j$である.
$a_k$の代わりに$y_k$を使うので$\hat{a}_j$ではなくPRMLと同じ$a_j$にする.
PRMLでは$w_{ji}$の肩の添え字を省略しているが念のためここではつけておく.

$\ww{ji}{1}$に対応する値を$v_{ji}$とすると$\ww{ji}{1}$の線形和である$a_j$について
$$
\calRR{a_j} =\sum_i x_i \calRR{\ww{ji}{1}} 	= \sum_i v_{ji} x_i.
$$
$$
\calRR{z_j}=\trans{v}\nabla \left(\sum_j \ww{ji}{1} h(a_j)\right)=\trans{v}\left(\diff{a_j}{h(a_j)} \nabla a_j\right)=h'(a_j)\calRR{a_j}.
$$
\begin{eqnarray*}
\calRR{y_k}
 &=&\trans{v^{(2)}}\left(\nabla \left(\sum_j \ww{kj}{2}z_j\right)\right)=\trans{v^{(2)}}\left(\sum_j \left(\nabla \ww{kj}{2}\right)z_j+\sum_j \ww{kj}{2}\nabla z_j\right)\\
 &=&\sum_j v_{kj}^{(2)} z_j+\sum_j\ww{kj}{2}\calRR{z_j}.
\end{eqnarray*}
なんとなくルールが見えてきたであろう.
$\calRR{\cdot}$は$\calRR{w}=v$という記号の置き換え以外は積や合成関数の微分のルールの形に従っている（もともと微分作用素を用いて定義しているので当然ではあるが）.

逆伝播の式：
\begin{eqnarray*}
\deltak&=&y_k-t_k,\\
\deltaj&=&h'(a_j)\sum_k\ww{kj}{2}\deltak
\end{eqnarray*}
で考えてみると
$$
\calRR{\deltak}=\calRR{y_k}.
$$
$$
\calRR{\deltaj}
 =h''(a_j)\calRR{a_j}\left(\sum_k \ww{kj}{2}\deltak \right)+h'(a_j)\left(\sum_k v_{kj}^{(2)} \deltak + \sum_k \ww{kj}{2} \calRR{\deltak}\right).
$$
誤差の微分の式:
\begin{eqnarray*}
\diff{\ww{kj}{2}}{E}&=&\deltak z_j.\\
\diff{\ww{jk}{1}}{E}&=&\deltaj x_i.
\end{eqnarray*}
より
\begin{eqnarray*}
\calRR{\diff{\ww{kj}{2}}{E}}&=&\calRR{\deltak}z_j+\deltak \calRR{z_j}.\\
\calRR{\diff{\ww{jk}{1}}{E}}&=&x_i\calRR{\deltaj}.
\end{eqnarray*}

\section{ソフト重み共有}
ネットワークの, あるグループに属する重みを等しくすることで複雑さを減らす手法がある.
しかし重みが等しいという制約は厳しい.
ソフト重み共有はその制約を外し, 代わりに正則化項を追加することで,
あるグループに属する重みが似た値をとれるようにする手法である.
$\pi_k$を混合係数として確率密度関数は
$$
p(w)=\prod_i p(w_i), \quad p(w_i)=\sum_{k=1}^M \pi_k \calN(w_i|\mu_k,\sigma_k^2).
$$
$p(w_i)$が確率分布なので混合係数は$\sum_k \pi_k=1$, $0 \le \pi_k \le 1$を満たす.
2乗ノルムの正規化項は平均0のガウス事前分布の負の対数尤度関数とみなせた.
ここでは複数個の重みに対応させるため混合ガウス分布を用いてみる.
$$
\Omega(w)=-\log p(w)=-\sum_i\log\left(\sum_{k=1}^M \pi_k \calN(w_i|\mu_k,\sigma_k^2)\right).
$$
最小化したい目的関数は誤差関数と正則化項の和で
$$
\tilde{E}(w)=E(w)+\Omega(w).
$$
$p(j)=\pi_j$とおいて負担率を導入する.
$$
\gamma_j(w_i)
 =p(j|w_i)
 =\frac{p(j)p(w_i|j)}{p(w_i)}=
 \frac{\pi_j \calN(w_i|\mu_j,\sigma_j^2)}{p(w_i)}.
$$
正規分布の微分
\begin{eqnarray*}
\dif{x}\calN(x|\mu,\sigma)&=&\calN(x|\mu,\sigma)\left(-\frac{x-\mu}{\sigma^2}\right),\\
\dif{\mu}\calN(x|\mu,\sigma)&=&\calN(x|\mu,\sigma)\left(\frac{x-\mu}{\sigma^2}\right)\\
\dif{\sigma}\calN(x|\mu,\sigma)&=&\calN(x|\mu,\sigma)\left(-\frac{1}{\sigma}+\frac{(x-\mu)^2}{\sigma^3}\right)
\end{eqnarray*}
を思い出しておく.
$\log p(w_i)$を$w_i$で微分すると
\begin{eqnarray*}
\dif{w_i}\log p(w_i)
 &=& \frac{1}{p(w_i)} \left(\sum_k \pi_k \dif{w_i} \calN(w_i|\mu_k,\sigma_k^2)\right)\\
 &=& \frac{1}{p(w_i)} \left(\sum_k \pi_k \calN(w_i|\mu_k,\sigma_k^2) \left(-\frac{w_i-\mu_k}{\sigma_k^2}\right)\right)\\
 &=&-\sum_k \frac{\pi_k \calN(w_i|\mu_k,\sigma_k)}{p(w_i)}\frac{w_i-\mu_k}{\sigma_k^2}\\
 &=&-\sum_k \gamma_k(w_i)\frac{w_i-\mu_k}{\sigma_k^2}.
\end{eqnarray*}
よって
$$
\diff{w_i}{\tilde{E}}=\diff{w_i}{E}+\sum_k \gamma_k(w_i)\frac{w_i-\mu_k}{\sigma_k^2}.
$$
同様に
\begin{eqnarray*}
\dif{\mu_k}\log p(w)
 &=&\sum_j \dif{\mu_k}\log p(w_j)
  =\sum_j \frac{\pi_k \calN(w_j|\mu_k,\sigma_k^2)}{p(w_j)}\frac{w_j-\mu_k}{\sigma_k^2}\\
 &=&\sum_j \gamma_k(w_j)\frac{w_j-\mu_k}{\sigma_k^2}.
\end{eqnarray*}
よって
$$
\diff{\mu_k}{\tilde{E}}=\sum_j \gamma_k(w_j)\frac{\mu_k-w_j}{\sigma_k^2}.
$$
\begin{eqnarray*}
\diff{\sigma_k}{\tilde{E}}
 &=&-\sum_j \dif{\sigma_k}\log p(w_j)
  =-\sum_j \frac{\pi_k \calN(w_j|\mu_k,\sigma_k^2)}{p(w_j)}
    \left(-\frac{1}{\sigma_k}+\frac{(w_j-\mu_k)^2}{\sigma_k^3}\right)\\
 &=&\sum_j \gamma_k(w_j)\left(\frac{1}{\sigma_k}-\frac{(w_j-\mu_k)^2}{\sigma_k^3}\right).
\end{eqnarray*}

$\pi_j$に関する制約より補助変数$\eta_j$を用いて
$$
\pi_j=\frac{\exp(\eta_j)}{\sum_k \exp(\eta_k)}
$$
と表すと\ref{takurasu}節式(\ref{ch4_mclass})より
$$
\diff{\eta_j}{\pi_k}=\pi_k(\delta_{kj}-\pi_j).
$$
よって
\begin{eqnarray*}
\diff{\eta_j}{\tilde{E}}
 &=& -\sum_i \dif{\eta_j} \log p(w_i)=-\sum_i\dif{\eta_j} \log \left(\sum_k \pi_k \calN(w_i|\mu_k,\sigma_k^2)\right)\\
 &=& -\sum_{i,k} \frac{\calN(w_i|\mu_k,\sigma_k^2)}{p(w_i)}\diff{\eta_j}{\pi_k}\\
 &=& -\sum_{i,k} \frac{\calN(w_i|\mu_k,\sigma_k^2)}{p(w_i)}\pi_k(\delta_{kj}-\pi_j)\\
 &=& -\sum_i \left(\frac{\pi_j \calN(w_i|\mu_j,\sigma_j^2)}{p(w_i)}
  - \frac{\pi_j \sum_k \pi_k \calN(w_i|\mu_k,\sigma_k^2)}{p(w_i)}\right)\\
 &=& -\sum_i (\gamma_j(w_i) - \pi_j) = \sum_i (\pi_j - \gamma_j(w_i)).
\end{eqnarray*}

\section{混合密度ネットワーク}
$$
p(t|x)=\sum_{k=1}^K \pi_k(x)\calN(t|\mu_k(x),\sigma_k^2(x)I)
$$
という分布のモデルを考える.
このモデルのパラメータを, $x$を入力としてえられるニューラルネットワークの出力となるようにとることで推論する.
前節と同様$\sum_k \pi_k(x)=1$, $0 \le \pi_k(x) \le 1$という制約があるので変数$a_l^\pi$を導入し
$$
\pi_k(x)=\frac{\exp(a_k^\pi)}{\sum_l \exp(a_l^\pi)}
$$
とする. 分散は0以上という制約があるので変数$a_k^\sigma$を導入し
$$
\sigma_k(x)=\exp(a_k^\sigma)
$$
とする. 平均は特に制約がないので
$$
\mu_{kj}(x)=a_{kj}^\mu
$$
とする.
$$
\calN_{nk}=\calN(t_n|\mu_k(x_n), \sigma_k^2(x_n)I)
$$
とおくとデータが独立の場合, 誤差関数は
$$
E(w)=-\sum_n \log \left(\sum_k \pi_k(x_n) \calN_{nk} \right).
$$
前節と同様$p(k|x)=\pi_k(x)$とおいて負担率を
$$
\gamma_{nk}(t_n|x_n)=p(k|t_n,x_n)=\frac{p(k|x_n)p(t_n|k)}{p(t_n|x_n)}=\frac{\pi_k \calN_{nk}}{\sum_l \pi_l \calN_{nl}}
$$
とする.
$$
\diff{a_k^\pi}{\pi_j}=\pi_j(\delta_{kj}-\pi_k)
$$
より
$$
\diff{a_k^\pi}{E_n}
 =-\frac{\sum_j \pi_j(\delta_{kj}-\pi_k)\calN_{nj}}{\sum_j \pi_j \calN_{nj}}
 =-(\gamma_{nk}-\pi_k)
 =\pi_k - \gamma_{nk}.
$$
$$
\calN(t|\mu,\sigma^2 I) =\frac{1}{(2\pi)^{L/2}}\frac{1}{\sigma^L}\exp\left(-\frac{1}{2\sigma^2}\sum_{l=1}^L(t_l-\mu_l)^2\right)
$$
なので
$$
\dif{\mu_l}\calN(t|\mu,\sigma^2 I)
 = \calN(t|\mu,\sigma^2 I)\left(-\frac{t_l-\mu_l}{\sigma^2}\right).
$$
よって
$$
\diff{a_{kl}^\mu}{E_n}
 =-\frac{\calN_{nk}}{\sum_j \pi_j \calN_{nj}}\frac{t_{nl}-\mu_{kl}}{\sigma_k^2}
 =\gamma_{nk}\left(\frac{\mu_{kl}-t_{nl}}{\sigma_k^2}\right).
$$
同様に
$$
\dif{\sigma}\calN(t|\mu,\sigma^2 I)
 = \calN(t|\mu,\sigma^2 I)\left(-\frac{L}{\sigma} + \frac{||t-\mu||^2}{\sigma^3}\right)
$$
より
$$
\dif{a_{kl}^\mu}\calN_{nj}=\delta_{jk}\calN_{nj}\left(\frac{t_{nl}-\mu_{kl}}{\sigma_k^2}\right).
$$
よって
$$
\diff{a_k^\sigma}{\calN_{nk}}
 =\diff{a_k^\sigma}{\sigma_k}\diff{\sigma_k}{\calN_{nk}}
 =\sigma_k \calN_{nk}\left(-\frac{L}{\sigma_k} + \frac{||t_n-\mu_k||^2}{\sigma_k^3}\right)
 =\calN_{nk}\left(-L + \frac{||t_n-\mu_k||^2}{\sigma_k^2}\right).
$$
よって
$$
\diff{a_k^\sigma}{E_n}=\gamma_{nk}\left(L-\frac{||t_n-\mu_k||^2}{\sigma_k^2}\right).
$$
条件付き平均についての密度関数の分散は
\begin{eqnarray*}
s^2(x)
 &=& E\left[||t-E[t|x]||^2\bigl|x\right]\\
 &=& \sum_k \pi_k \int \left(||t||^2-2\trans{t}E[t|x]+||E[t|x]||^2\right)\calN(t|\mu_k,\sigma_k^2 I)\,dt\\
 &=& \sum_k \pi_k \left(\sigma_k^2+||\mu_k||^2-2\trans{\mu_k}E[t|x]+||E[t|x]||^2\right)\\
 &=& \sum_k \pi_k(x)\left(\sigma_k(x)^2+||\mu_k-\sum_j \pi_j(x)\mu_j(x)||^2\right).
\end{eqnarray*}

\section{クラス分類のためのベイズニューラルネットワーク}
ロジスティックシグモイド出力を一つ持つネットワークによる2クラス分類問題を考える.
そのモデルの対数尤度関数は$t_n \in \{0, 1\}$, $y_n=y(x_n,w)$として
$$
\log p(\calD|w)=\sum_n (t_n \log y_n + (1-t_n)\log (1-y_n)).
$$
事前分布を
$$
p(w|\alpha)=\calN(w|0,\alpha^{-1}I)=\frac{1}{\sqrt{2\pi}^W|\alpha^{-1}|^{W/2}}\exp\left(-\half \alpha \inp{w}\right)
$$
とする（$W$は$w$に含まれるパラメータの総数）. ノイズがないので$\beta$を含まない.
$$
E(w)=\log p(\calD|w)+\frac{\alpha}{2}\inp{w}
$$
の最小化で$\wmap$を求め, $A=-\nabla^2 \log p(w|\calD)\bigl|_{w=\wmap}$を何らかの方法で求める.
ラプラス近似を使って事後分布をガウス近似すると
$$
q(w|\calD)=\calN(w|\wmap,A^{-1}).
$$
正規化項を求める\ref{ch4_bic}節式(\ref{ch4_approx})を使って
\begin{eqnarray*}
\log p(\calD|\alpha)
 &\approx& \log \left(p(\calD|\wmap)p(\wmap|\alpha)\sqrt{\frac{(2\pi)^W}{|A|}}\right)\\
 &=& \log p(\calD|\wmap)-\frac{W}{2}\log (2\pi) +\frac{W}{2}\log \alpha
 - \half \alpha \inp{\wmap}\\
 &+& \frac{W}{2}\log(2\pi)-\half \log |A|\\
 &=&-E(\wmap)-\half \log |A| + \frac{W}{2}\log \alpha.
\end{eqnarray*}
$$
E(\wmap)=-\sum_n (t_n \log y_n + (1-t_n)\log (1-y_n))+\half \alpha \inp{\wmap}.
$$
予測分布を考える.
出力ユニットの活性化関数を線形近似する.
\begin{eqnarray*}
a(x,w)
 &\approx& a(x,\wmap) + \nabla \trans{a(x,\wmap)}(w-\wmap)\\
 &&（\amap(x)=a(x,\wmap), \quad b = \nabla a(x,\wmap)として）\\
 &=& \amap(x) + \trans{b}(w-\wmap).
\end{eqnarray*}
\begin{eqnarray*}
p(a|x,\calD)
 &=& \int \delta(a-a(x,w))q(w|\calD)\,dw\\
 &=&\int \delta (a-\amap(x)-\trans{\wmap}b+\trans{w}b)q(w|\calD)\,dw.
\end{eqnarray*}
平均は
\begin{eqnarray*}
E[a]
 &=& \int ap(a|x,\calD)\,da = \int \delta(a-a(x,w))w(w)a\,dadw\\
 &=& \int a(x,w)q(w)\,dw=(\amap(x)-\trans{\wmap}b)+\int \trans{b}wq(w)\,dw\\
 &=& \amap(x)-\trans{\wmap}b+\trans{b}\wmap=\amap(x).
\end{eqnarray*}
分散は$\trans{w}b$が効くので
$$
\sigma_a^2(x)=\trans{b}A^{-1}b(x)
$$
予測分布は\ref{ch4_bayes}節式(\ref{ch4_probit})の近似式を使って
$$
p(t=1|x,\calD)=\int \sigma(a)p(a|x,\calD)\,da
 \approx \sigma(\kappa(\sigma_a^2)\amap(x)).
$$