はじめに
PRML解答例まとめを参照
演習 10.28
10.2節で導入したベイズ混合ガウスモデルを10.4節で議論した指数分布族とその共役事前分布のモデルとして書き換えよ．すなわち，一般的な結果
\begin{aligned} \ln q^{\star}(\mathbf{Z}) &=\mathbb{E}_{\eta}[\ln p(\mathbf{X}, \mathbf{Z} \mid \eta)]+\text { const } \\ &=\sum_{n=1}^{N}\left\{\ln h\left(\mathbf{x}_{n}, \mathbf{z}_{n}\right)+\mathbb{E}\left[\boldsymbol{\eta}^{\mathrm{T}}\right] \mathbf{u}\left(\mathbf{x}_{n}, \mathbf{z}_{n}\right)\right\}+\text { const } \end{aligned} \tag{10.115}
q^{\star}(\eta)=f\left(\nu_{N}, \boldsymbol{\chi}_{N}\right) g(\eta)^{\nu_{N}} \exp \left\{\nu_{N} \eta^{\mathrm{T}} \boldsymbol{\chi}_{N}\right\} \tag{10.119}
を用いて，特定の場合の結果
q^{\star}(\mathbf{Z})=\prod_{n=1}^{N} \prod_{k=1}^{K} r_{n k}^{z_{n k}} \tag{10.48}
q^{\star}(\boldsymbol{\pi})=\operatorname{Dir}(\boldsymbol{\pi} \mid \boldsymbol{\alpha}) \tag{10.57}
q^{\star}\left(\boldsymbol{\mu}_{k}, \mathbf{\Lambda}_{k}\right)=\mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \mathbf{\Lambda}_{k}\right)^{-1}\right) \mathcal{W}\left(\mathbf{\Lambda}_{k} \mid \mathbf{W}_{k}, \nu_{k}\right) \tag{10.59}
を導け．
ベイズ混合ガウス分布における、指数型分布族のそれぞれの関数形を導出して、一般的な結果(10.115),(10.119)に代入していく。
１変数ガウス分布と、指数型分布族の標準形との対応関係は、演習2.57により(2.220)〜(2.223)のとおり導出済み。
多変数ガウス分布（混合分布ではない）との対応関係は、これを拡張して、p(\mathbf{x}|\boldsymbol{\eta})=h(\mathbf{x})g(\boldsymbol{\eta})\exp [\boldsymbol{\eta}^T \mathbf{u}(\mathbf{x})]において、
\begin{aligned}
\boldsymbol{\eta}
:=\left[\begin{array}{c}
\boldsymbol{\eta}_1 \\
\boldsymbol{\eta}_2
\end{array}\right]
&\leftrightarrow
\left[\begin{array}{c}
\Lambda \boldsymbol{\mu} \\
-\frac{1}{2}\vec{\Lambda}
\end{array}\right] \\
\mathbf{u}(\mathbf{x}) &\leftrightarrow \left[\begin{array}{c}
\mathbf{x} \\
\mathbf{x}\mathbf{x}^T
\end{array}\right] \\
h(\mathbf{x})&\leftrightarrow \frac{1}{(2\pi)^{D/2}}\\
g(\boldsymbol{\eta})&\leftrightarrow |-2\boldsymbol\eta _2|^{1/2} \exp \left( \frac{1}{4}\boldsymbol{\eta}_1^T \boldsymbol\eta _2 ^{-1}\boldsymbol\eta _1\right)
\end{aligned}
と書ける。ただし、行列の上に矢印（\rightarrow）が書かれているのは、行列の各要素を並べたD \times D次元のベクトルを意味する。後の式変形の都合で、g(\boldsymbol{\eta})の要素を\boldsymbol{\eta}と\mathbf{u}(\mathbf{x})に押し込めて、
\begin{aligned}
\boldsymbol{\eta}
&\leftrightarrow
\left[\begin{array}{c}
\Lambda \boldsymbol{\mu} \\
-\frac{1}{2}\vec{\Lambda}\\
\boldsymbol{\mu}^T \boldsymbol\Lambda \boldsymbol\mu \\
\ln |\boldsymbol{\Lambda} |
\end{array}\right] \\
\mathbf{u}(\mathbf{x}) &\leftrightarrow \left[\begin{array}{c}
\mathbf{x} \\
\mathbf{x}\mathbf{x}^T\\
-\frac{1}{2}\\
\frac{1}{2}
\end{array}\right] \\
h(\mathbf{x})&\leftrightarrow \frac{1}{(2\pi)^{D/2}}\\
g(\boldsymbol{\eta})&\leftrightarrow 1
\end{aligned}
と書き直す。
今考えているベイズ混合ガウスモデル:
\begin{aligned}
p(\mathbf{X}, \mathbf{Z}, \boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol\Lambda )
&=p(\mathbf{X}| \mathbf{Z}, \boldsymbol{\mu}, \boldsymbol\Lambda )
p(\mathbf{Z}| \boldsymbol{\pi} )
p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol\Lambda )\\
&=\left(\prod_{n=1}^N \prod_{k=1}^K \pi_k^{z_{nk}}\right)
\left(\prod_{n=1}^N \prod_{k=1}^K \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})^{z_{nk}} \right)
p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol\Lambda )\\
&=\left(\prod_{n=1}^N \prod_{k=1}^K \left\{ \pi_k
\mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})\right\}^{z_{nk}} \right)
p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol\Lambda )\\
\end{aligned}
の形に徐々に近づけていく。まずは\pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})を指数型分布族の標準形に対応づけるには、
\begin{aligned}
\boldsymbol{\eta}
&\leftrightarrow
\left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right] \\
\mathbf{u}(\mathbf{x}_n) &\leftrightarrow \left[\begin{array}{c}
\mathbf{x}_n \\
\overrightarrow{\mathbf{x}_n\mathbf{x}_n^T}\\
-\frac{1}{2}\\
\frac{1}{2}\\
1
\end{array}\right] \\
h(\mathbf{x}_n)&\leftrightarrow \frac{1}{(2\pi)^{D/2}}\\
g(\boldsymbol{\eta})&\leftrightarrow 1
\end{aligned}
次に、\left\{ \pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})\right\}^{z_{nk}}を指数型分布族の標準形に対応づけるには、
\begin{aligned}
\boldsymbol{\eta}
&\leftrightarrow
\left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right] \\
\mathbf{u} (\mathbf{x}_n,z_{nk})&\leftrightarrow z_{nk} \left[\begin{array}{c}
\mathbf{x}_n \\
\overrightarrow{\mathbf{x}_n\mathbf{x}_n^T}\\
-\frac{1}{2}\\
\frac{1}{2}\\
1
\end{array}\right] \\
h(\mathbf{x}_n,z_{nk})&\leftrightarrow \left( \frac{1}{(2\pi)^{D/2}} \right)^{z_{nk}}\\
g(\boldsymbol{\eta})&\leftrightarrow 1
\end{aligned}
最後に、\prod_{k=1}^K \left\{ \pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})\right\}^{z_{nk}}を指数型分布族の標準形に対応づけるには、
\begin{aligned}
\boldsymbol{\eta}
&\leftrightarrow
\left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right] _{k=1,\cdots,K}\\
\mathbf{u} (\mathbf{x}_n,\mathbf{z}_{n})&\leftrightarrow
\left[ z_{nk} \left[\begin{array}{c}
\mathbf{x}_n \\
\overrightarrow{\mathbf{x}_n\mathbf{x}_n^T}\\
-\frac{1}{2}\\
\frac{1}{2}\\
1
\end{array}
\right] \right] _{k=1,\cdots,K}\\
h(\mathbf{x}_n,\mathbf{z}_n)&\leftrightarrow \prod_{k=1}^K \left( \frac{1}{(2\pi)^{D/2}} \right)^{z_{nk}}\\
g(\boldsymbol{\eta})&\leftrightarrow 1
\end{aligned}
ここで、[\ \ \ ]_{k=1,\cdots,K}とは、k=1に対応するベクトル、k=2に対応するベクトル\cdotsと順に並べてできる長いベクトルを表す。
今得られた対応関係を(10.116)式に代入して、
\begin{aligned}
q^\star (\mathbf{z}_n)
&= h(\mathbf{x}_n, \mathbf{z}_n) g(\mathbb{E}[\boldsymbol\eta])\exp \{\mathbb{E} [\boldsymbol{\eta}^T] \mathbf{u}(\mathbf{x}_n , \mathbf{z}_n )\}\\
&= \left\{ \prod_{k=1}^K \left(\frac{1}{(2\pi)^{D/2}}\right)^{z_{nk}}\right\}
\cdot 1 \cdot
\exp \left\{ \sum_{k=1}^K
\left( \mathbb{E}\left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right]
\right)^T \left( z_{nk} \left[\begin{array}{c}
\mathbf{x}_n \\
\overrightarrow{\mathbf{x}_n\mathbf{x}_n^T}\\
-\frac{1}{2}\\
\frac{1}{2}\\
1
\end{array}
\right]
\right)
\right\}\\
&= \left\{ \prod_{k=1}^K \left(\frac{1}{(2\pi)^{D/2}}\right)^{z_{nk}}\right\}
\exp \left\{ \sum_{k=1}^K z_{nk}
\cdot\mathbb{E}_{\boldsymbol\mu_k , \boldsymbol\Lambda_k}\left[
\boldsymbol{\mu}_k^T \Lambda_k
\mathbf{x}_n
-\frac{1}{2}\mathbf{x}_n^T\Lambda_k
\mathbf{x}_n
-\frac{1}{2}
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k
+\frac{1}{2}
\ln |\boldsymbol{\Lambda}_k|
+\ln\pi_k
\right]
\right\}\\
&= \prod_{k=1}^K \left(\frac{1}{(2\pi)^{D/2}}
\exp \left\{ -\frac{1}{2}
\mathbb{E}_{\boldsymbol\mu_k , \boldsymbol\Lambda_k}\left[
(\mathbf{x}_n-\boldsymbol{\mu}_k)^T \Lambda_k (\mathbf{x}_n-\boldsymbol{\mu}_k)\right]
+\frac{1}{2}
\mathbb{E}[\ln |\boldsymbol{\Lambda}_k|]
+\mathbb{E}[\ln\pi_k]
\right\}
\right)^{z_{nk}}\\
&= \prod_{k=1}^K \left(
\exp \left\{ - \frac{D}{2}\ln (2\pi)-\frac{1}{2}
\mathbb{E}_{\boldsymbol\mu_k , \boldsymbol\Lambda_k}\left[
(\mathbf{x}_n-\boldsymbol{\mu}_k)^T \Lambda_k (\mathbf{x}_n-\boldsymbol{\mu}_k)\right]
+\frac{1}{2}
\mathbb{E}[\ln |\boldsymbol{\Lambda}_k|]
+\mathbb{E}[\ln\pi_k]
\right\}
\right)^{z_{nk}}\\
&=\prod_{k=1}^K \rho _{nk}^{z_{nk}}
\end{aligned}
を得る。（最後の式変形は、(10.46)式の\rho _{nk}の定義より。）以上より、
\begin{aligned}
q^\star (\mathbf{Z})
&= \prod_{n=1}^{N} q^\star (\mathbf{z}_n)\\
&= \prod_{n=1}^{N}\prod_{k=1}^K \rho _{nk}^{z_{nk}}
\end{aligned}
と(10.48)式を得る。
次に、(10.119)式を用いて(10.57)式と(10.59)式を導く。
\begin{aligned}
q^\star (\boldsymbol\eta )
&\propto g(\boldsymbol\eta)^{\nu_N} \exp \left[ \nu _N \boldsymbol \eta^T \boldsymbol{\chi} _N \right]\\
&= 1 \cdot \exp \left[ \boldsymbol\eta^T \left( \nu_0\boldsymbol\chi_0+\sum_{n=1}^N \mathbb{E}_{z_n}[\mathbf{u}(\mathbf{x}_n,\mathbf{z}_n)]\right)\right]\\
&= \exp \left(
\nu_0\boldsymbol\eta^T \boldsymbol\chi_0
+\sum_{n=1}^N\sum_{k=1}^K \mathbb{E}[z_{nk}] \left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right]
^T \left[\begin{array}{c}
\mathbf{x}_n \\
\overrightarrow{\mathbf{x}_n\mathbf{x}_n^T}\\
-\frac{1}{2}\\
\frac{1}{2}\\
1
\end{array}
\right]
\right)\\
&= \exp \left[
\nu_0\boldsymbol\eta^T \boldsymbol\chi_0
+\sum_{n=1}^N\sum_{k=1}^K r_{nk} \left(
\boldsymbol{\mu}_k^T \Lambda_k
\mathbf{x}_n
-\frac{1}{2}\mathbf{x}_n^T\Lambda_k
\mathbf{x}_n
-\frac{1}{2}
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k
+\frac{1}{2}
\ln |\boldsymbol{\Lambda}_k|
+\ln\pi_k
\right)
\right]\\
&= \exp \left[
\nu_0\boldsymbol\eta^T \boldsymbol\chi_0
+\sum_{n=1}^N\sum_{k=1}^K r_{nk} \left( -\frac{1}{2}
(\mathbf{x}_n-\boldsymbol{\mu}_k)^T \Lambda_k (\mathbf{x}_n-\boldsymbol{\mu}_k)
+\frac{1}{2}
\ln |\boldsymbol{\Lambda}_k|
+\ln\pi_k
\right)
\right]
\end{aligned}
事前確率分布の項\exp\left[\nu_0\boldsymbol\eta^T\boldsymbol\chi_0\right]は、図10.5の事前確率分布
\begin{aligned}
p(\boldsymbol\pi,\boldsymbol\mu,\boldsymbol\Lambda)
&= p(\boldsymbol\pi)p(\boldsymbol\mu|\boldsymbol\Lambda)p(\boldsymbol\Lambda)\\
&\propto \prod_{k=1}^K \pi_k ^{\alpha_0-1}
|\boldsymbol\Lambda_k |^{1/2}\exp \left[-\frac{1}{2}
(\boldsymbol\mu_k-\mathbf{m}_0)^T (\beta_0\boldsymbol\Lambda_k)
(\boldsymbol\mu_k-\mathbf{m}_0)\right]
|\boldsymbol\Lambda_k|^{(\nu_0-D-1)/2}\exp\left(-\frac{1}{2}{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k) \right)\\
&=
\exp \left[\sum_{k=1}^K\left\{
(\alpha_0-1)\ln \pi_k
-\frac{1}{2}
(\boldsymbol\mu_k-\mathbf{m}_0)^T (\beta_0\boldsymbol\Lambda_k)
(\boldsymbol\mu_k-\mathbf{m}_0)
+\frac{\nu_0-D}{2}
\ln
|\boldsymbol\Lambda_k|
-\frac{1}{2}{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k) \right\}\right]
\end{aligned}
に一致するように\nu_0\boldsymbol\chi_0を選ぶことは可能（指数型分布族の事前共役分布の形を所与とすれば、この事実は証明不要な気もするが、念のため本問の解答の最後で\nu_0\boldsymbol\chi_0の具体的な形を構築する）。
\begin{aligned}
q^
\star (\boldsymbol\eta )
\propto &
\exp \left[\sum_{k=1}^K\left\{
(\alpha_0-1)\ln \pi_k
-\frac{1}{2}
(\boldsymbol\mu_k-\mathbf{m}_0)^T (\beta_0\boldsymbol\Lambda_k)
(\boldsymbol\mu_k-\mathbf{m}_0)
+\frac{\nu_0-D}{2}
\ln
|\boldsymbol\Lambda_k|
-\frac{1}{2}{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k) \right\}\right]\\
&\cdot \exp \left[
\sum_{n=1}^N\sum_{k=1}^K r_{nk} \left( -\frac{1}{2}
(\mathbf{x}_n-\boldsymbol{\mu}_k)^T \Lambda_k (\mathbf{x}_n-\boldsymbol{\mu}_k)
+\frac{1}{2}
\ln |\boldsymbol{\Lambda}_k|
+\ln\pi_k
\right)
\right]
\end{aligned}
この同時確率分布は\boldsymbol\piに依存する項だけ積の形でくくり出せて、
\begin{aligned}
q^\star (\boldsymbol\eta )
&\propto
\exp \left[\sum_{k=1}^K\left\{
(\alpha_0-1)\ln \pi_k
+\sum_{n=1}^N r_{nk}
\ln\pi_k
\right\}
\right]\\
&=\exp \left[\sum_{k=1}^K\left\{
(\alpha_0+N_k-1)\ln \pi_k
\right\}
\right]\\
&=\prod_{k=1}^K \pi_k^{\alpha_0+N_k-1}
\end{aligned}
であり、(10.57)式のディリクレ分布が導かれた。
残りの項のうち、\boldsymbol\mu_kに依存する項のみをくくり出すと、
\begin{aligned}
q^\star (\boldsymbol\eta )
\propto \ &
\exp \left[
-\frac{1}{2}\left\{ \sum_{k=1}^K \boldsymbol\mu_k^T \left(\beta_0 \boldsymbol\Lambda_k + \sum_{n=1}^N r_{nk}\boldsymbol\Lambda_k\right)\boldsymbol\mu_k
-2\left(\boldsymbol\mu_k^T\beta_0\boldsymbol\Lambda_k \mathbf{m}_0 + \sum_{n=1}^N r_{nk} \boldsymbol\mu_k^T \Lambda_k \mathbf{x}_n
\right)
\right\}
\right]\\
=:\ &\exp \left[
-\frac{1}{2}\left\{ \sum_{k=1}^K \boldsymbol\mu_k^T \left((\beta_0+N_k) \boldsymbol\Lambda_k \right)\boldsymbol\mu_k
-2\boldsymbol\mu_k^T\left(\boldsymbol\Lambda_k (\beta_0 \mathbf{m}_0 + N_{k} \overline{\mathbf{x}_k}
) \right)
\right\}
\right]\\
= \ &|\boldsymbol\Lambda_k|^{1/2} \exp \left[
-\frac{1}{2} \sum_{k=1}^K \left( \boldsymbol\mu_k - \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}\right)^T \left((\beta_0+N_k) \boldsymbol\Lambda_k \right)\left( \boldsymbol\mu_k - \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}\right)
\right]
\\
& \cdot |\boldsymbol\Lambda_k|^{-1/2}\exp \left[ \frac{1}{2} \frac{1}{\beta_0+N_k} (\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^T\Lambda_k(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})\right]\\
\end{aligned}
１行目は、\boldsymbol\mu_kが(10.59)式のガウス分布に従うことを表す。
最後に、その他の項をくくり出すと、
\begin{aligned}
q^\star (\boldsymbol\eta )
\propto &
\exp \left[\sum_{k=1}^K\left\{
-\frac{1}{2}
\mathbf{m}_0^T \beta_0\boldsymbol\Lambda_k
\mathbf{m}_0
+\frac{\nu_0-D}{2}
\ln
|\boldsymbol\Lambda_k|
-\frac{1}{2}{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k)
+\sum_{n=1}^N r_{nk} \left( -\frac{1}{2}
\mathbf{x}_n^T \Lambda_k \mathbf{x}_n
+\frac{1}{2}
\ln |\boldsymbol{\Lambda}_k|
\right)
\right\}
\right]
\\ \cdot & |\boldsymbol\Lambda_k|^{-1/2}\exp \left[ \frac{1}{2} \frac{1}{\beta_0+N_k} (\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^T\Lambda_k(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})\right]\\
=& \prod_{k=1}^K
|\boldsymbol\Lambda_k|^{\frac{\nu_0+N_k-D-1}{2}}
\exp \left[
-\frac{1}{2}
\left(
\mathbf{m}_0^T \beta_0\boldsymbol\Lambda_k
\mathbf{m}_0
+{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k)
+\sum_{n=1}^N r_{nk}
\mathbf{x}_n^T \Lambda_k \mathbf{x}_n
-\frac{1}{\beta_0+N_k} (\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^T\Lambda_k(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})
\right)
\right]\\
=& \prod_{k=1}^K
|\boldsymbol\Lambda_k|^{\frac{\nu_0+N_k-D-1}{2}}
\exp \left[
-\frac{1}{2}
\left\{
{\rm Tr} \left( \mathbf{m}_0 \mathbf{m}_0^T \beta_0\boldsymbol\Lambda_k
\right)
+{\rm Tr}\left(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k\right)
+{\rm Tr} \left( \sum_{n=1}^N r_{nk}
\mathbf{x}_n \mathbf{x}_n^T \Lambda_k \right)
-\frac{1}{\beta_0+N_k} {\rm Tr}\left(
(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})
(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^T\Lambda_k
\right)
\right\}
\right]\\
=& \prod_{k=1}^K
|\boldsymbol\Lambda_k|^{\frac{\nu_0+N_k-D-1}{2}}
\exp \left[
-\frac{1}{2}
{\rm Tr}\left\{ \left( \beta_0 \mathbf{m}_0 \mathbf{m}_0^T
+\mathbf{W}_0^{-1} + \sum_{n=1}^N r_{nk}
\mathbf{x}_n \mathbf{x}_n^T -\frac{1}{\beta_0+N_k}
(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})
(\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^T
\right)\Lambda_k
\right\}
\right]
\\
\end{aligned}
となる。最後の式は、|\Lambda|^{*}\exp [-\frac{1}{2}Tr(*\Lambda)]という形をしており、(10.59)式のウィシャート分布を得る。（一見すると(10.62)式の形と異なるように見えるが、\mathbf{S}_kの定義に戻って計算すると、上式と一致することが示せる。
おまけ：
{\nu}_0\boldsymbol\chi_0の具体的な形状は、以下の通り。
\begin{aligned}
{\nu}_0\boldsymbol\chi_0
&=
\left[\begin{array}{c}
\beta_0 \mathbf{m}_0 \\
\beta_0 \overrightarrow{\mathbf{m}_0 \mathbf{m}_0^T}+\overrightarrow{\mathbf{W}_0^{-1}}\\
-\frac{1}{2}\beta_0\\
\frac{\nu_0-D-1}{2}\\
\alpha_0-1
\end{array}\right] _{k=1,\cdots,K}
\end{aligned}
念のため、上の式を再現できることを確認しておく。
\begin{aligned}
{\nu}_0 \boldsymbol{\eta} ^T
\boldsymbol\chi_0
&=
\left[\begin{array}{c}
\Lambda_k \boldsymbol{\mu}_k \\
-\frac{1}{2}\vec\Lambda_k\\
\boldsymbol{\mu}_k^T \boldsymbol\Lambda_k \boldsymbol\mu_k \\
\ln |\boldsymbol{\Lambda}_k|\\
\ln\pi_k
\end{array}\right] _{k=1,\cdots,K}^T
\left[\begin{array}{c}
\beta_0 \mathbf{m}_0 \\
\beta_0 \overrightarrow{\mathbf{m}_0 \mathbf{m}_0^T}+\overrightarrow{\mathbf{W}_0^{-1}}\\
-\frac{1}{2}\beta_0\\
\frac{\nu_0-D-1}{2}\\
\alpha_0-1
\end{array}\right] _{k=1,\cdots,K}\\
&=
\sum_{k=1}^K\left\{
(\alpha_0-1)\ln \pi_k
-\frac{1}{2}
(\boldsymbol\mu_k-\mathbf{m}_0)^T (\beta_0\boldsymbol\Lambda_k)
(\boldsymbol\mu_k-\mathbf{m}_0)
+\frac{\nu_0-D-1}{2}
\ln
|\boldsymbol\Lambda_k|
-\frac{1}{2}{\rm Tr}(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k) \right\}
\end{aligned}
となる。なお、最後のTr(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k)の項への変形にあたって、一般の正方行列A,Bについて
\begin{aligned}
{\rm Tr}[AB]={\rm Tr}\left[\left( \sum_j a_{ij}b_{jk}\right)_{ik} \right]
=\sum_{i, j}a_{ij}b_{ji}
\end{aligned}
が成り立つこと、およびBが対称行列であれば最右辺が\sum_{i,j} a_{ij}b_{ij}に一致することを用いた。
演習 10.29
二階微分を計算することで，関数f(x) = \ln (x)は0 \lt x \lt \inftyで上に凸であることを示せ．
g(\eta)=\min _{x}\{\eta x-f(x)\} \tag{10.133}
で定義される双対な関数
g(\eta)の形を求め，
f(x)=\min _{\eta}\{\eta x-g(\eta)\} \tag{10.132}
に従って
\eta x - g(\eta)を
\etaについて最小化すると，実際に関数
\ln(x)が得られることを確かめよ．
※P.209のようにf(x)とg(\eta)が双対の働きとなっていることを示す。
f(x) = \ln(x)の2階微分は\displaystyle f''(x) = - \frac{1}{x^{2}}である。これは0 \lt x \lt \inftyでf''(x) < 0となるので上に凸である。
(10.133)式に代入して
g(\eta)=\min _{x}\{\eta x-\ln(x)\}
最小値を求めるためxについて微分すると
\frac{dg}{dx} = \eta - \frac{1}{x}
\displaystyle x=\frac{1}{\eta}のとき最小となり、このときg(\eta) = 1+\ln(\eta)となる。
今度はこれを(10.132)式に代入して\displaystyle f(x) = \min _{\eta}\{\eta x-1-\ln(\eta)\}となる。これも同様に\etaについて微分すると
\frac{df}{d\eta} = x - \frac{1}{\eta}
これは\displaystyle \eta = \frac{1}{x}のとき最小となり、このとき
f(x) = 1-\left( 1+\ln\left(\frac{1}{x}\right) \right) = \ln(x)
となり題意通り\ln(x)が得られた。
演習 10.30
二階微分を計算することで，対数ロジスティック関数f(x)=-\ln \left(1+e^{-x}\right)が上に凸であることを示せ．対数ロジスティック関数を点x = \xiのまわりで一次までテイラー展開することで，変分上界
\sigma(x) \leqslant \exp (\eta x-g(\eta)) \tag{10.137}
を直接導出せよ．
後の問題設定のため、対数ロジスティック関数であるf(x)を
\sigma(x) = \frac{1}{1+e^{-x}} \tag{10.134}
を用いて
f(x)=-\ln \left(1+e^{-x}\right)=\ln \left(\frac{1}{1+e^{-x}}\right)=\ln \sigma(x)
と表しておく。これより
\begin{aligned}
\frac{d \sigma}{d x} &=\frac{e^{-x}}{\left(1+e^{-x}\right)^{2}} = (\sigma(x))^2e^{-x} \\
&=\frac{1}{1+e^{-x}}-\frac{1}{\left(1+e^{-x}\right)^2} \\
&=\sigma(x)-(\sigma(x))^{2} \\
&=\sigma(x)(1-\sigma(x))
\end{aligned}
である。これからf(x)が上に凸であること、すなわちf''(x) \lt 0を示す。
f^{\prime}(x)=\frac{d}{d x} \ln \sigma(x)=\frac{1}{\sigma(x)} \frac{d}{d x} \sigma(x)=\frac{1}{\sigma(x)}(\sigma(x))^{2} e^{-x}=\sigma(x) e^{-x}
\begin{aligned} f^{\prime \prime}(x) &=\frac{d}{d x}\left(\sigma(x) e^{-x}\right) \\
&=\frac{d}{d x} \sigma(x) \cdot e^{-x}-\sigma(x) e^{-x} \\
&=\{\sigma(x)(1-\sigma(x))-\sigma(x)\} e^{-x} \\
&=-\{\sigma(x)\}^{2} e^{-x} \lt 0 \end{aligned}
これよりf(x)は上に凸であることが示された。
次に\ln \sigma(x)をx=\xiの周りで一次のテイラー展開を行うと
\begin{aligned} \ln \sigma(x) &=\ln \sigma(\xi)+\left.\frac{d\ln \sigma(x)}{d x}\right|_{x=\xi}(x-\xi)+O\left(\xi^{2}\right) \\
&\simeq \ln \sigma(\xi)+\sigma(\xi) e^{-\xi}(x-\xi) \end{aligned}
\ln \sigma(x)は上に凸（concave）な関数なので、このテイラー展開は
\ln \sigma(x) \leq \ln \sigma(\xi)+\sigma(\xi) e^{-\xi}(x-\xi) \tag{A}
となる。
変分上界を求めるために、10.5節の議論のように（(10.127)あたりの変形）\eta = \sigma(\xi)e^{-\xi}とおいたときの\sigma(\xi)と\xiを\etaで表すことを目指す。
\begin{aligned} \eta &=\left(1+e^{-\xi}\right)^{-1} e^{-\xi} \\
&=\frac{e^{-\xi}}{1+e^{-\xi}} \\
&=1-\frac{1}{1+e^{-\xi}} \\
&=1-\sigma(\xi) \end{aligned}
これより\sigma(\xi) = 1 - \etaとなる。また、\eta = \sigma(\xi)e^{-\xi}の両辺の対数をとって
\begin{aligned} \ln \eta &=\ln \sigma(\xi)-\xi \\ \xi &=\ln \sigma(\xi)-\ln \eta \\
&=\ln (1-\eta)-\ln \eta \end{aligned}
となる。これらを(\textrm{A})に代入して
\begin{aligned} \ln \sigma(x) & \leq \ln (1-\eta)+\eta(x-\ln (1-\eta)+\ln \eta) \\
&=\eta x+(1-n) \ln (1-\eta)+\eta \ln \eta \\
&=\eta x-g(\eta)\hspace{1em}(\because(10.136))
\end{aligned}
これは変分上界
\sigma(x) \leqslant \exp (\eta x-g(\eta)) \tag{10.137}
と等価である。すなわち、題意の通りテイラー展開から直接(10.137)式が導出された。
演習 10.31
xについて二階微分を求めることで，関数f(x) = -\ln(e^{x/2}+e^{-x/2})はxの上に凸な関数であることを示せ．次に，変数x^2についての二階微分を考え，これがx^2については下に凸な関数で、あることを示せ．f(x)のグラフをxおよびx^2について描いてみよ．関数f(x)を変数x^2について\xi^2を中心として一次までテイラー展開することにより，ロジスティックシグモイド関数の下界
\sigma(x) \geqslant \sigma(\xi) \exp \left\{(x-\xi) / 2-\lambda(\xi)\left(x^{2}-\xi^{2}\right)\right\} \tag{10.144}
を直接導出せよ．
微分を計算する時の係数が煩わしいので、x/2をxと定義し直す。（証明すべき凸性に影響はない。最後の計算結果で元に戻す。）
\begin{aligned}
f(x) &= -\ln(e^{x}+e^{-x})\\
f'(x) &= -\frac{e^{x}-e^{-x}}{e^x+e^{-x}}\\
f''(x) &= -\frac{(e^{x}+e^{-x})^2-(e^{x}-e^{-x})^2}{(e^{x}+e^{-x})^2}\\
&= -\frac{(2e^x)\cdot (2e^{-x})}{(e^{x}+e^{-x})^2}\\
&= -\frac{4}{(e^{x}+e^{-x})^2} < 0
\end{aligned}
y=x^2とおく。f(y)=-\ln (e^{\sqrt{y}}+e^{-\sqrt{y}})を真面目にyで2階微分するのは面倒なので、以下のように解く。
\begin{aligned}
\frac{d}{dy}f(x) &= \frac{dx}{dy}\frac{df(x)}{dx}\\
\frac{d^2}{dy^2}f(x) &= \frac{d}{dy}\left( \frac{dx}{dy} \frac{df(x)}{dx} \right)\\
&= \frac{d^2x}{dy^2}\frac{df(x)}{dx} + \frac{dx}{dy} \frac{d}{dy} \left( \frac{df(x)}{dx} \right)\\
&= \frac{d^2x}{dy^2}\frac{df(x)}{dx} + \left( \frac{dx}{dy} \right)^2 \cdot \frac{d^2 f(x)}{dx^2}\\
\end{aligned}
ここで、
\begin{aligned}
\frac{dx}{dy} &= \left( \frac{dy}{dx} \right)^{-1} = \left( \frac{d}{dx} x^2 \right)^{-1} = (2x)^{-1} = \frac{1}{2x} \\
\frac{d^2x}{dy^2} &= \frac{d}{dy} \left( \frac{dx}{dy} \right) = \frac{d}{dy} \left( \frac{1}{2x} \right)= \frac{dx}{dy}\frac{d}{dx} \left( \frac{1}{2x} \right)
= \frac{1}{2x} \left( -\frac{1}{2x^2}\right) = -\frac{1}{4x^3}\\
\end{aligned}
であるから、
\begin{aligned}
\frac{d^2}{dy^2}f(x) &= \frac{d^2x}{dy^2}\frac{df(x)}{dx} + \left( \frac{dx}{dy} \right)^2 \cdot \frac{d^2 f(x)}{dx^2}\\
&= \left(-\frac{1}{4x^3} \right)\left(-\frac{e^x-e^{-x}}{e^x+e^{-x}} \right) + \left( \frac{1}{2x}\right)^2 \left( - \frac{4}{(e^x+e^{-x})^2}\right)\\
&=\frac{(e^x-e^{-x})(e^x+e^{-x})-4x}{4x^3 (e^x+e^{-x})^2}\\
&=\frac{e^{2x}-e^{-2x}-4x}{4x^3 (e^x+e^{-x})^2}
\end{aligned}
ここで、分子がx>0のとき正、x<0のとき負であることは、微分して普通に示しても良いが、
\begin{aligned}
e^{2x} &= 1+(2x)+ \frac{1}{2!}(2x)^2+\frac{1}{3!}(2x)^3+\frac{1}{4!}(2x)^4 +\frac{1}{5!}(2x)^5 +\cdots\\
e^{-2x} &= 1-(2x)+ \frac{1}{2!}(2x)^2-\frac{1}{3!}(2x)^3+\frac{1}{4!}(2x)^4 -\frac{1}{5!}(2x)^5 +\cdots\\
e^{2x}-e^{-2x} &= 4x+ 2\cdot \left( \frac{1}{3!}(2x)^3+\frac{1}{5!}(2x)^5\cdots \right)\\
\end{aligned}
に注目すれば明らか。よって、右辺全体はx>0でもx<0でも正なので、\frac{d^2}{dy^2}f(x)>0と言える。
\begin{aligned}
f(x)&=f(x)|_{x=\xi} + (x^2-\xi^2) \left( \left. \frac{df(x)}{dy} \right|_{x=\xi} \right)+ \cdots\\
&\geq f(x)|_{x=\xi} + (x^2-\xi^2) \left(\left. \frac{df(x)}{dy} \right|_{x=\xi}\right)\\
&= f(\xi) + (x^2-\xi^2)\left( \frac{dx}{dy} \left.\frac{f(x)}{dx}\right|_{x=\xi}\right)\\
&= f(\xi) + (x^2-\xi^2)\left( -\frac{1}{4\xi} \tanh(\xi) \right)\\
&\equiv f(\xi) - \lambda(\xi) (x^2-\xi^2)\\
\end{aligned}
となる（不等号は、x^2に関する凸性より）。従って、
\begin{aligned}
\ln \sigma(x) &= \frac{x}{2} + f(x)\\
&\geq \frac{x}{2} + f(\xi) - \lambda(\xi) (x^2-\xi^2)\\
&= \frac{x-\xi}{2} + \left( \frac{\xi}{2}+f(\xi)\right) - \lambda(\xi) (x^2-\xi^2)\\
&= \frac{x-\xi}{2} + \ln \sigma (\xi) - \lambda(\xi) (x^2-\xi^2)\\
\end{aligned}
両辺のexpをとって、(10.144)式を得る。
演習 10.32
ロジスティック回帰の変分ベイズ推定を時系列に従って学習し，データ点が一回に一つだけ到着し，次のデータ点が到着するまでに処理して廃棄しなければならない場合について考える．このとき，事後分布のガウス近似は下界
p(t \mid \mathbf{w})=e^{a t} \sigma(-a) \geqslant e^{a t} \sigma(\xi) \exp \left\{-(a+\xi) / 2-\lambda(\xi)\left(a^{2}-\xi^{2}\right)\right\} \tag{10.151}
を用いることで常に保持できることを示せ．この場合，分布は事前分布で初期化され，データ点が取り込まれるごとに。対応する変分パラメータ
\xi_nが最適化される．
※(方針)
(10.151)式を利用しガウス変分事後分布の逐次的な更新式が得られることを確認する
データN個のときの変分事後分布を
q_N(\mathbf{w})=\mathcal{N}(\mathbf{w}\mid\mathbf{m}_N,\mathbf{S}_N^{-1})
とする．N+1個めのデータが与えられたとき(10.151)より
p(t_{N+1}|,\mathbf{w})=e^{\mathbf{w}^{\top} \phi_{N+1} t_{N+1}} \sigma\left(-\mathbf{w}^{\top} \phi_{N+1}\right) \geq e^{\mathbf{w}^{\top} \phi_{N+1} t_{N+1}} \sigma\left(\xi_{N+1}\right) \exp \left[-\left(\mathbf{w}^{\top} \phi_{N+1}+\xi_{N+1}\right) / 2-\lambda\left(\xi_{N+1}\right)\left\{\left(\mathbf{w}^{\top} \phi_{N+1}\right)^{2}-\xi_{N+1}^{2}\right\}\right].
この式から
p(t_{N+1},\mathbf{w})=p(t_{N+1}\mid\mathbf{w})p(\mathbf{w})\geq p(\mathbf{w})e^{\mathbf{w}^{\top} \phi_{N+1} t_{N+1}} \sigma\left(\xi_{N+1}\right) \exp \left[-\left(\mathbf{w}^{\top} \phi_{N+1}+\xi_{N+1}\right) / 2-\lambda\left(\xi_{N+1}\right)\left\{\left(\mathbf{w}^{\top} \phi_{N+1}\right)^{2}-\xi_{N+1}^{2}\right\}\right].
を得る両辺の対数をとって
\begin{aligned}
\ln p(t_{N+1}\mid\mathbf{w})p(\mathbf{w})&\geq\ln p(\mathbf{w})+\mathbf{w}^{\top} \phi_{N+1} t_{N+1}+\ln \sigma\left(\xi_{N+1}\right)-\left(\mathbf{w}^{\top} \phi_{N+1}+\xi_{N+1}\right) / 2-\lambda\left(\xi_{N+1}\right)\left\{\left(\mathbf{w}^{\top} \phi_{N+1}\right)^{2}-\xi_{N+1}^{2}\right\}\\
&\simeq -\frac{1}{2}(\mathbf{w}-\mathbf{m}_N)^{\top}\mathbf{S}_N^{-1}(\mathbf{w}-\mathbf{m}_N)+\mathbf{w}^{\top} \phi_{N+1}\left(t_{N+1}-\frac{1}{2}\right)-\lambda\left(\xi_{N+1}\right) \mathbf{w}^{\top} \phi_{N+1} \phi_{N+1}^{\top} \mathbf{w}+Const
\end{aligned}
これは\mathbf{w}の二次形式になっているため，N+1個目のデータが与えられたときの変分事後分布をq_{N+1}(\mathbf{w})とするとガウス分布になり，平方完成することで
q_{N+1}(\mathbf{w})=\mathcal{N}(\mathbf{w}\mid\mathbf{m}_{N+1},\mathbf{S}_{N+1}^{-1})
\mathbf{m}_{N+1}＝\mathbf{S}_{N+1}\left[\mathbf{S}_{N}^{-1} \mathbf{m}_{N}+\left(t_{N+1}-1 / 2\right) \boldsymbol{\phi}_{N+1}\right]
\mathbf{S}_{N+1}^{-1}＝2 \lambda\left(\xi_{N+1}\right) \boldsymbol{\phi}_{N+1} \boldsymbol{\phi}_{N+1}^{\mathrm{T}}+\mathbf{S}_{N}^{-1}
が得られる．
演習 10.33
Q\left(\boldsymbol{\xi}, \boldsymbol{\xi}^{\text {old }}\right)=\sum_{n=1}^{N}\left\{\ln \sigma\left(\xi_{n}\right)-\xi_{n} / 2-\lambda\left(\xi_{n}\right)\left(\boldsymbol{\phi}_{n}^{\mathrm{T}} \mathbb{E}\left[\mathbf{w} \mathbf{w}^{\mathrm{T}}\right] \boldsymbol{\phi}_{n}-\xi_{n}^{2}\right)\right\}+\mathrm{const} . \tag{10.161}
で定義した量Q\left(\boldsymbol{\xi}, \boldsymbol{\xi}^{\text {old }}\right)を変分パラメータ\xi_nについて微分することで，
\left(\xi_{n}^{\text {new }}\right)^{2}=\boldsymbol{\phi}_{n}^{\mathrm{T}} \mathbb{E}\left[\mathbf{w} \mathbf{w}^{\mathrm{T}}\right] \boldsymbol{\phi}_{n}=\boldsymbol{\phi}_{n}^{\mathrm{T}}\left(\mathbf{S}_{N}+\mathbf{m}_{N} \mathbf{m}_{N}^{\mathrm{T}}\right) \boldsymbol{\phi}_{n} \tag{10.163}
で与えられるベイズロジスティック回帰モデルの
\xi_nの更新式を示せ．
Preparation
\begin{aligned}
\frac{\partial \sigma(x)}{\partial x} &= \frac{e^{-x}}{(1 + e^{-x})^2} \\
&= \sigma(x)(1 - \sigma(x))
\end{aligned}
\begin{aligned}
\frac{\partial \log(\sigma(x))}{\partial x} &= (1 + e^{-x})\frac{e^{-x}}{(1 + e^{-x})^2} \\
&= \frac{e^{-x}}{1 + e^{-x}} = 1 - \sigma(x)
\end{aligned}
\lambda(\xi) = \frac{1}{2\xi}\left[\sigma(\xi) - \frac{1}{2}\right] \tag{10.150}
\begin{aligned}
\mathbf{S}_N &= \mathbb{E}[\mathbf{ww}^{\mathrm{T}}] - \mathbb{E}[\mathbf{w}]\mathbb{E}[\mathbf{w}^{\mathrm{T}}] \\
&= \mathbb{E}[\mathbf{ww}^{\mathrm{T}}] - \mathbf{m}_N\mathbf{m}^{\mathrm{T}}_N
\end{aligned}
Solution
\begin{aligned}
\frac{\partial Q}{\partial \xi_n} &= (1 - \sigma(\xi_n)) - \frac{1}{2} + 2\xi_n\lambda(\xi_n) - \lambda'(\xi_n)(\mathbf{\phi}^{\mathrm{T}}_n\mathbb{E}\left[\mathbf{ww}^{\mathrm{T}}\right]\mathbf{\phi}_n - \xi^2_n) \\
&= -2\xi_n\lambda(\xi_n) + 2\xi_n\lambda(\xi_n) - \lambda'(\xi_n)(\mathbf{\phi}^{\mathrm{T}}_n\mathbb{E}\left[\mathbf{ww}^{\mathrm{T}}\right]\mathbf{\phi}_n - \xi^2_n) \\
&= 0
\end{aligned}
よって、
(\xi^{\mathrm{new}}_n)^2 = \mathbf{\phi}^{\mathrm{T}}_n\mathbb{E}\left[\mathbf{ww}^{\mathrm{T}}\right]\mathbf{\phi}_n = \mathbf{\phi}^{\mathrm{T}}_n\left(\mathbf{S}_N + \mathbf{m}_N\mathbf{m}^{\mathrm{T}}_N\right)\mathbf{\phi}_n
演習 10.34
この演習問題では，4.5節のベイスロジスティック回帰モデルの変分パラメータ\boldsymbol{\xi}の更新式を，
\begin{aligned} \mathcal{L}(\boldsymbol{\xi})&= \frac{1}{2} \ln \frac{\left|\mathbf{S}_{N}\right|}{\left|\mathbf{S}_{0}\right|}+\frac{1}{2} \mathbf{m}_{N}^{\mathrm{T}} \mathbf{S}_{N}^{-1} \mathbf{m}_{N}-\frac{1}{2} \mathbf{m}_{0}^{\mathrm{T}} \mathbf{S}_{0}^{-1} \mathbf{m}_{0} \\ &+\sum_{n=1}^{N}\left\{\ln \sigma\left(\xi_{n}\right)-\frac{1}{2} \xi_{n}+\lambda\left(\xi_{n}\right) \xi_{n}^{2}\right\} \end{aligned} \tag{10.164}
で与えられる下界を直接最大化することで導出する．これには
\mathcal{L}(\boldsymbol{\xi})の
\xi_nに附する微分を0としてみよ．行列式の対数の微分には結果
\frac{d}{d \alpha} \ln |\mathbf{A}|=\operatorname{Tr}\left(\mathbf{A}^{-1} \frac{d}{d \alpha} \mathbf{A}\right) \tag{3.117}
を，変分事後分布
q(\mathbf{w})の平均と分散には
\mathbf{m}_{N}=\mathbf{S}_{N}\left(\mathbf{S}_{0}^{-1} \mathbf{m}_{0}+\sum_{n=1}^{N}\left(t_{n}-1 / 2\right) \phi_{n}\right) \tag{10.157}
\mathbf{S}_{N}^{-1}=\mathbf{S}_{0}^{-1}+2 \sum_{n=1}^{N} \lambda\left(\xi_{n}\right) \phi_{n} \phi_{n}^{\mathrm{T}} \tag{10.158}
の式を利用せよ．
Preparation
\frac{\partial}{\partial x}(\mathbf{A}^{-1}) = -\mathbf{A}^{-1}\frac{\partial\mathbf{A}}{\partial x}\mathbf{A}^{-1} \tag{C.21}
\begin{aligned}
\mathbf{m}^{\mathrm{T}}_N\mathbf{S}^{-1}_{N}\mathbf{m}_N &= [\mathbf{S}_{N}\mathbf{S}^{-1}_{N}\mathbf{m}_N]^{\mathrm{T}}\mathbf{S}^{-1}_{N}[\mathbf{S}_{N}\mathbf{S}^{-1}_{N}\mathbf{m}_N] \\
&= [\mathbf{S}_{N}\mathbf{a}_N]^{\mathrm{T}}\mathbf{S}^{-1}_{N}[\mathbf{S}_{N}\mathbf{a}_N] \\
&= \mathbf{a}^{\mathrm{T}}_N\mathbf{S}^{-1}_{N}\mathbf{a}_N.
\end{aligned}
ただし、
\mathbf{a}_N = \mathbf{S}^{-1}_{N}\mathbf{m}_N.
Solution
\begin{aligned}
\frac{\partial \mathcal{L}(\mathbf{\xi})}{\partial \xi_n} &= \frac{1}{2}\mathrm{Tr}\left(\mathbf{S}^{-1}_{N}\frac{\partial \mathbf{S}_{N}}{\partial \xi_n}\right) + \frac{1}{2} \mathrm{Tr}\left(\mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N \frac{\partial \mathbf{S}_{N}}{\partial \xi_n}\right) + \lambda'(\xi_n)\xi^2_n \\
&= \frac{1}{2}\mathrm{Tr}\left((\mathbf{S}^{-1}_{N} + \mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N)\frac{\partial \mathbf{S}_{N}}{\partial \xi_n}\right) + \lambda'(\xi_n)\xi^2_n \\
&= -\frac{1}{2}\mathrm{Tr}\left((\mathbf{S}^{-1}_{N} + \mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N)\mathbf{S}_N[2\lambda'(\xi_n)\mathbf{\phi}_N\mathbf{\phi}^{\mathrm{T}}_N]\mathbf{S}_N\right) + \lambda'(\xi_n)\xi^2_n \\
&= -\lambda'(\xi_n)\left\{\mathrm{Tr}\left((\mathbf{S}^{-1}_{N} + \mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N)\mathbf{S}_N\mathbf{\phi}_N\mathbf{\phi}^{\mathrm{T}}_N\mathbf{S}_N\right) - \xi^2_n\right\} \\
&= 0.
\end{aligned}
よって、
\begin{aligned}
\xi^2_n &= \mathrm{Tr}\left((\mathbf{S}^{-1}_{N} + \mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N)\mathbf{S}_N\mathbf{\phi}_N\mathbf{\phi}^{\mathrm{T}}_N\mathbf{S}_N\right) \\
&= (\mathbf{S}_N\mathbf{\phi}_n)^{\mathrm{T}}(\mathbf{S}^{-1}_{N} + \mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N)(\mathbf{S}_N\mathbf{\phi}_n) \\
&= \mathbf{\phi}^{\mathrm{T}}_n(\mathbf{S}_N + \mathbf{S}_N\mathbf{a}_N\mathbf{a}^{\mathrm{T}}_N\mathbf{S}_N)\mathbf{\phi}_n \\
&= \mathbf{\phi}^{\mathrm{T}}_n(\mathbf{S}_N + \mathbf{m}_N\mathbf{m}^{\mathrm{T}}_N)\mathbf{\phi}_n
\end{aligned}
演習 10.35
変分ベイズロジスティック回帰モデルの下界\mathcal{L}(\boldsymbol{\xi})についての結果
\begin{aligned} \mathcal{L}(\boldsymbol{\xi})&= \frac{1}{2} \ln \frac{\left|\mathbf{S}_{N}\right|}{\left|\mathbf{S}_{0}\right|}+\frac{1}{2} \mathbf{m}_{N}^{\mathrm{T}} \mathbf{S}_{N}^{-1} \mathbf{m}_{N}-\frac{1}{2} \mathbf{m}_{0}^{\mathrm{T}} \mathbf{S}_{0}^{-1} \mathbf{m}_{0} \\ &+\sum_{n=1}^{N}\left\{\ln \sigma\left(\xi_{n}\right)-\frac{1}{2} \xi_{n}+\lambda\left(\xi_{n}\right) \xi_{n}^{2}\right\} \end{aligned} \tag{10.164}
を導出せよ．これには
\mathcal{L}(\boldsymbol{\xi})を定める積分
\ln p(\mathbf{t})=\ln \int p(\mathbf{t} \mid \mathbf{w}) p(\mathbf{w}) \mathrm{d} \mathbf{w} \geqslant \ln \int h(\mathbf{w}, \boldsymbol{\xi}) p(\mathbf{w}) \mathrm{d} \mathbf{w}=\mathcal{L}(\boldsymbol{\xi}) \tag{10.159}
の中で，ガウス事前分布の式に
q(\mathbf{w}) = \mathcal{N}(\mathbf{w}\mid \mathbf{m}_0, \mathbf{\Sigma}_0)を，尤度関数に下界
h(\mathbf{w}, \boldsymbol{\xi})を代入するのが最も簡単である．次に指数関数の中で
\mathbf{w}に依存する項をまとめ，平方完成することでガウス分布の積分を導出せよ．多次元ガウス分布の正規化項についての標準的な結果を使ってこれを計算し，最後に対数をとって
(10.164)を求めよ．
-
h(\mathbf{w}, \xi)p(\mathbf{w})を計算する
\begin{aligned}
h(\mathbf{w}, \xi)p(\mathbf{w}) &= \mathrm{N}(\mathbf{w}|\mathbf{m}_0, \mathbf{S}_0)\prod^N \sigma(\xi_n)\exp\{\mathbf{w}^{\mathrm{T}}\mathbf{\phi}_nt_n - (\mathbf{w}^{\mathrm{T}}\mathbf{\phi}_n + \xi_n)/2 - \lambda(\xi_n)([\mathbf{w}^{\mathrm{T}}\mathbf{\phi}]^2 - \xi^2_n)\} \\
&= \{(2\pi)^{-W/2} \cdot |\mathbf{S}_0|^{-1/2} \cdot \prod^N \sigma(\xi_n)\} \cdot \exp\{-\frac{1}{2}(\mathbf{w} - \mathbf{m}_0)^{\mathrm{T}}\mathbf{S}^{-1}_0(\mathbf{w} - \mathbf{m}_0)\} \\
& \ \ \ \cdot \prod^N \exp\{\mathbf{w}^{\mathrm{T}}\mathbf{\phi}_nt_n - (\mathbf{w}^{\mathrm{T}}\mathbf{\phi}_n + \xi_n)/2 - \lambda(\xi_n)([\mathbf{w}^{\mathrm{T}}\mathbf{\phi}]^2 - \xi^2_n)\} \\
&= \left\{(2\pi)^{-W/2} \cdot |\mathbf{S}_0|^{-1/2} \cdot \prod^N \sigma(\xi_n) \cdot \exp\left(-\frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 - \sum^N \frac{\xi_n}{2} + \sum^N \lambda(\xi_n)\xi^2_n\right)\right\} \\
& \ \ \ \cdot \exp\left\{-\frac{1}{2}\mathbf{w}^{\mathrm{T}}\left(\mathbf{S}^{-1}_0 + 2\sum^N\lambda(\xi_n)\mathbf{\phi}_n\mathbf{\phi}^{\mathrm{T}}_n\right)\mathbf{w} + \mathbf{w}^{\mathrm{T}}\left(\mathbf{S}^{-1}_0\mathbf{m}_0 + \sum^N \mathbf{\phi}_n\left(t_n - \frac{1}{2}\right)\right) \right\} \\
&= \left\{(2\pi)^{-W/2} \cdot |\mathbf{S}_0|^{-1/2} \cdot \prod^N \sigma(\xi_n) \cdot \exp\left(-\frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 - \sum^N \frac{\xi_n}{2} + \sum^N \lambda(\xi_n)\xi^2_n\right)\right\} \\
& \ \ \ \cdot \exp\left\{\frac{1}{2}\mathbf{w}^{\mathrm{T}}\mathrm{S}^{-1}_N\mathbf{w} + \mathbf{w}^{\mathrm{T}}\mathrm{S}^{-1}_N\mathbf{m}_N\right\} \ (\because (10.157)-(10.158)) \\
&= \left\{(2\pi)^{-W/2} \cdot |\mathbf{S}_0|^{-1/2} \cdot \prod^N \sigma(\xi_n) \cdot \exp\left(-\frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 - \sum^N \frac{\xi_n}{2} + \sum^N \lambda(\xi_n)\xi^2_n\right) + \frac{1}{2}\mathbf{m}^{\mathrm{T}}_N\mathbf{S}^{-1}_N\mathbf{m}_N\right\} \\
& \ \ \ \cdot \exp\left\{-\frac{1}{2}(\mathbf{w} - \mathbf{m}_N)^{\mathrm{T}}\mathrm{S}^{-1}_N(\mathbf{w} - \mathbf{m}_N)\right\}
\end{aligned}
\begin{aligned}
\int h(\mathbf{w}, \xi)p(\mathbf{w})d\mathbf{w} &= (2\pi)^{-W/2} \cdot |\mathbf{S}_0|^{-1/2} \\
& \ \ \ \cdot \prod^N \sigma(\xi_n) \cdot \exp\left(-\frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 - \sum^N \frac{\xi_n}{2} + \sum^N \lambda(\xi_n)\xi^2_n\right) + \frac{1}{2}\mathbf{m}^{\mathrm{T}}_N\mathbf{S}^{-1}_N\mathbf{m}_N \\
& \ \ \ \cdot (2\pi)^{W/2} \cdot |\mathbf{S}_N|^{1/2} \\
&= \left(\frac{|\mathbf{S}_N|}{|\mathbf{S}_0|}\right)^{1/2} \prod^N \sigma(\xi_n) \cdot \exp\left(-\frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 - \sum^N \frac{\xi_n}{2} + \sum^N \lambda(\xi_n)\xi^2_n\right) + \frac{1}{2}\mathbf{m}^{\mathrm{T}}_N\mathbf{S}^{-1}_N\mathbf{m}_N
\end{aligned}
\mathcal{L}(\xi) = \frac{1}{2}\log\frac{|\mathbf{S}_N|}{|\mathbf{S}_0|} - \frac{1}{2}\mathbf{m}^{\mathrm{T}}_0\mathbf{S}^{-1}_0\mathbf{m}_0 + \frac{1}{2} \mathbf{m}^{\mathrm{T}}_N\mathbf{S}^{-1}_N\mathbf{m}_N + \sum^N \left\{\log\sigma(\xi_n) - \frac{\xi_n}{2} + \lambda(\xi_n)\xi^2_n \right\}
演習 10.36
10.7節で議論したADFの近似法を考える．このとき.因子f_j(\boldsymbol{\theta})を含めることで，モデルエビデンスを
p_{j}(\mathcal{D}) \simeq p_{j-1}(\mathcal{D}) Z_{j} \tag{10.242}
のように更新できることを示せ．ここでZ_jは正規化定数であり
Z_{j}=\int f_{j}(\theta) q^{\backslash j}(\boldsymbol{\theta}) \mathrm{d} \boldsymbol{\theta} \tag{10.197}
で与えられる
p_0({\mathcal{D}} ) = 1と初期化して上の結果を再帰的に適用することで，
p(\mathcal{D}) \simeq \prod_{j} Z_{j} \tag{10.243}
を導け．
※
演習 10.37
10.7節のEP法のアルゴリスムについて考え，定義
p(\mathcal{D}, \theta)=\prod_{i} f_{i}(\boldsymbol{\theta}) \tag{10.188}
における因子の一つ
f_0(\boldsymbol{\theta})が，近似分布
q(\boldsymbol{\theta})と同じ形の指数分布族になっていると仮定する．このとき因子
\tilde{f}_0(\boldsymbol{\theta})を
f_0(\boldsymbol{\theta})で初期化すれば，EP法での
\tilde{f}_0(\boldsymbol{\theta})の更新は
\tilde{f}_0(\boldsymbol{\theta})を変えないことを示せ．典型的には，これは因子の一つが事前分布
p(\boldsymbol{\theta})の場合に起こり，事前分布に対応する因子は一度だけ厳密な形で取り込まれ，以降は更新する必要はないことがわかる．
q(\boldsymbol{\theta})の初期値は
q_{\mathrm{init}}(\boldsymbol{\theta}) = \tilde{f_0}(\boldsymbol{\theta})\prod_{i \neq 0}\tilde{f_i}(\boldsymbol{\theta})
と表される。ここで、
q^{\setminus 0}(\boldsymbol{\theta}) = \prod_{i \neq 0}\tilde{f_i}(\boldsymbol{\theta}).
q^{new}(\boldsymbol{\theta})は、モーメント一致法を用いて、
q^{new}(\boldsymbol{\theta}) = q^{\setminus 0}(\boldsymbol{\theta})f_0(\boldsymbol{\theta}).
問題の設定より、\tilde{f_0}(\boldsymbol{\theta})の初期値はf_0(\boldsymbol{\theta})なので、
q^{new}(\boldsymbol{\theta}) = q_{\mathrm{init}}(\boldsymbol{\theta}).
また、
Z_0 = \int q^{\setminus 0}(\boldsymbol{\theta})f_0(\boldsymbol{\theta})d\boldsymbol{\theta} = \int q_{\mathrm{init}}(\boldsymbol{\theta})d\boldsymbol{\theta} = 1.
従って、(10.207)の更新式は、
\tilde{f_0}(\boldsymbol{\theta}) = Z_0\frac{q^{new}(\boldsymbol{\theta})}{q^{\setminus 0}(\boldsymbol{\theta})} = f_0(\boldsymbol{\theta})
となる。
演習 10.38
この演習問題と次の演習問題では，雑音データ問題でのEP法の結果(10.214)-(10.224)を確かめる．除算を行う公式
q^{\backslash j}(\boldsymbol{\theta})=\frac{q(\boldsymbol{\theta})}{\widetilde{f}_{j}(\boldsymbol{\theta})} \tag{10.205}
から始めて，
\mathbf{m}^{\backslash n}=\mathbf{m}+v^{\backslash n} v_{n}^{-1}\left(\mathbf{m}-\mathbf{m}_{n}\right) \tag{10.214}
\left(v^{\backslash n}\right)^{-1}=v^{-1}-v_{n}^{-1}\tag{10.215}
を，指数関数の中を平方完成して平均と分散を見出すことで導け．また
Z_{j}=\int q^{\backslash j}(\boldsymbol{\theta}) f_{j}(\boldsymbol{\theta}) \mathrm{d} \boldsymbol{\theta} \tag{10.206}
で定義される正規化定数Z_nは，雑音データ問題については
Z_{n}=(1-w) \mathcal{N}\left(\mathbf{x}_{n} \mid \mathbf{m}^{\backslash n},\left(v^{\backslash n}+1\right) \mathbf{I}\right)+w \mathcal{N}\left(\mathbf{x}_{n} \mid \mathbf{0}, a \mathbf{I}\right) \tag{10.216}
で与えられることを示せ．これには一般的な結果
p(\mathbf{y})=\mathcal{N}\left(\mathbf{y} \mid \mathbf{A} \boldsymbol{\mu}+\mathbf{b}, \mathbf{L}^{-1}+\mathbf{A} \mathbf{\Lambda}^{-1} \mathbf{A}^{\mathrm{T}}\right) \tag{2.115}
を利用すればよい．
(10.205)、(10.212)、(10.213)式より、
\begin{aligned}
q^{\setminus j}(\boldsymbol{\theta}) &= \frac{q(\boldsymbol{\theta})}{\tilde{f_j}(\boldsymbol{\theta})} = \frac{\mathcal{N}(\boldsymbol{\theta}|\mathbf{m}, v\mathbf{I})}{s_n \mathcal{N}(\boldsymbol{\theta}|\mathbf{m}_n, v_n\mathbf{I})} \\
&\propto \frac{\exp\left\{-\frac{1}{2}(\boldsymbol{\theta} - \mathrm{m})^{\mathrm{T}}(v\mathbf{I})^{-1}(\boldsymbol{\theta} - \mathrm{m})\right\}}{\exp\left\{-\frac{1}{2}(\boldsymbol{\theta} - \mathrm{m}_n)^{\mathrm{T}}(v_n\mathbf{I})^{-1}(\boldsymbol{\theta} - \mathrm{m}_n)\right\}} \\
&\propto \exp\left\{-\frac{1}{2}(\boldsymbol{\theta} - \mathrm{m})^{\mathrm{T}}(v\mathbf{I})^{-1}(\boldsymbol{\theta} - \mathrm{m}) + \frac{1}{2}(\boldsymbol{\theta} - \mathrm{m}_n)^{\mathrm{T}}(v_n\mathbf{I})^{-1}(\boldsymbol{\theta} - \mathrm{m}_n)\right\} \\
&= \exp\left\{-\frac{1}{2}(\boldsymbol{\theta}^{\mathrm{T}}\mathbf{A}\boldsymbol{\theta} + \boldsymbol{\theta}^{\mathrm{T}}\mathbf{B})\right\} + \mathrm{const}
\end{aligned}
まず、
\begin{aligned}
(\boldsymbol{\Sigma}^{\setminus n})^{-1} &= \mathbf{A} = (v\mathbf{I})^{-1} - (v_n\mathbf{I})^{-1} \\
&= (v^{-1} - v^{-1}_n)\mathbf{I}^{-1}
\end{aligned}
また、
-2(\boldsymbol{\Sigma}^{\setminus n})^{-1}\mathbf{m}^{\setminus n} = \mathbf{B} = 2\left(-(v\mathbf{I})^{-1}\mathbf{m} + (v_n\mathbf{I})^{-1}\mathbf{m}_n\right)
より、
\begin{aligned}
\mathbf{m}^{\setminus n} &= -(\boldsymbol{\Sigma}^{\setminus n})^{-1}\left(-(v\mathbf{I})^{-1}\mathbf{m} + (v_n\mathbf{I})^{-1}\mathbf{m}_n\right) \\
&= -(v^{\setminus n})\left(-v^{-1}\mathbf{m} + (v_n^{-1}\mathbf{m}_n\right) \\
&= v^{\setminus n}v^{-1}\mathbf{m} - \frac{v^{\setminus n}}{v_n}\mathbf{m}_n \\
&= v^{\setminus n}((v^{\setminus n})^{-1} - v^{-1}_n)\mathbf{m} - \frac{v^{\setminus n}}{v_n}\mathbf{m}_n \\
&= \mathbf{m} + \frac{v^{\setminus n}}{v_n}(\mathbf{m} - \mathbf{m}_n)
\end{aligned}
次に、(10.206)、(10.209)式より、
\begin{aligned}
Z_n &= \int q^{\setminus n}(\boldsymbol{\theta})f_n(\boldsymbol{\theta})d\boldsymbol{\theta} \\
&= \int \mathcal{N}(\boldsymbol{\theta}|\mathbf{m}^{\setminus n}, v^{\setminus n}\mathbf{I})\left\{(1 - w)\mathcal{N}(\mathbf{x}_n|\boldsymbol{\theta}, \mathbf{I}) + w\mathcal{N}(\mathbf{x}_n|\mathbf{0}, a\mathbf{I})\right\} \\
&= (1 - w)\int \mathcal{N}(\boldsymbol{\theta}|\mathbf{m}^{\setminus n}, v^{\setminus n}\mathbf{I})\mathcal{N}(\mathbf{x}_n|\boldsymbol{\theta}, \mathbf{I})d\boldsymbol{\theta} + w\int \mathcal{N}(\boldsymbol{\theta}|\mathbf{m}^{\setminus n}, v^{\setminus n}\mathbf{I})\mathcal{N}(\mathbf{x}_n|\mathbf{0}, a\mathbf{I})d\boldsymbol{\theta} \\
&= (1 - w) \mathcal{N}(\mathbf{x}_n|\mathbf{m}^{\setminus n}, (v^{\setminus n} + 1)\mathbf{I}) + w\mathcal{N}(\mathbf{x}_n|\mathbf{0}, a\mathbf{I}). (\because (2.113)-(2.115))
\end{aligned}
演習 10.39
雑音データ問題でのEP法におけるq^{\text {new }}(\boldsymbol{\theta})の平均と分散は
\mathbf{m}^{\text {new }}=\mathbf{m}^{\backslash n}+\rho_{n} \frac{v^{\backslash n}}{v^{\backslash n}+1}\left(\mathbf{x}_{n}-\mathbf{m}^{\backslash n}\right) \tag{10.217}
v^{\text {new }}=v^{\backslash n}-\rho_{n} \frac{\left(v^{\backslash n}\right)^{2}}{v^{\backslash n}+1}+\rho_{n}\left(1-\rho_{n}\right) \frac{\left(v^{\backslash n}\right)^{2}\left\|\mathbf{x}_{n}-\mathbf{m}^{\backslash n}\right\|^{2}}{D\left(v^{\backslash n}+1\right)^{2}} \tag{10.218}
で与えられることを示せ．このためには，最初にq^{\text {new }}(\boldsymbol{\theta})の下での\boldsymbol{\theta}と\boldsymbol{\theta}\boldsymbol{\theta}^{\mathrm T}の期待値が
\begin{aligned} \mathbb{E}[\boldsymbol{\theta}] &=\mathbf{m}^{\backslash n}+v^{\backslash n} \nabla_{\mathbf{m} \backslash n} \ln Z_{n} \\ \mathbb{E}\left[\boldsymbol{\theta}^{\mathrm{T}} \boldsymbol{\theta}\right] &=2\left(v^{\backslash n}\right)^{2} \nabla_{v^{\backslash n}} \ln Z_{n}+2 \mathbb{E}[\boldsymbol{\theta}]^{\mathrm{T}} \mathbf{m}^{\backslash n}-\left\|\mathrm{m}^{\backslash n}\right\|^{2}+v^{\backslash n} D \end{aligned}
となることを証明し，Z_nについての結果
Z_{n}=(1-w) \mathcal{N}\left(\mathbf{x}_{n} \mid \mathbf{m}^{\backslash n},\left(v^{\backslash n}+1\right) \mathbf{I}\right)+w \mathcal{N}\left(\mathbf{x}_{n} \mid \mathbf{0}, a \mathbf{I}\right) \tag{10.216}
を用いる．次に，
\widetilde{f}_{j}(\boldsymbol{\theta})=Z_{j} \frac{q^{\text {new }}(\boldsymbol{\theta})}{q \backslash j(\theta)} \tag{10.207}
を用いて指数関数の中を平方完成することで(10.220)-(10.222)の結果を証明せよ．最後に，(10.208)を用いて(10.223)の結果を導け．
Preparation
\begin{aligned}
\frac{\partial \mathcal{N}(x|\mu, \sigma^2)}{\partial \mu} &= \mathcal{N}(x|\mu, \sigma^2)\cdot\frac{x - \mu}{\sigma^2} \\
\frac{\partial \mathcal{N}(x|\mu, \sigma^2)}{\partial \sigma^2} &= \mathcal{N}(x|\mu, \sigma^2)\cdot\left(\frac{(x - \mu)^2}{(\sigma^2)^2} - \frac{1}{\sigma^2}\right) \\
\end{aligned}
\begin{aligned}
\rho_n &= \frac{1}{Z_n}(1 - w)\mathcal{N}(\mathbf{x}_n|\mathbf{m}^{\setminus n}, (v^{\setminus n} + 1)\mathbf{I}) \\
&= \frac{1}{Z_n}(1 - w)\frac{Z_n - w\mathcal{N}(\mathbf{x}_n|\mathbf{0}, a\mathbf{I})}{1 - w} \ (\because (10.216)) \\
&= 1 - \frac{w}{Z_n}\mathcal{N}(\mathbf{x}_n|\mathbf{0}, a\mathbf{I}).
\end{aligned}
Solution
まず、
\begin{aligned}
\nabla_{\mathbf{m}^{\setminus n}} \ln Z_n &= \frac{1}{Z_n}\nabla_{\mathbf{m}^{\setminus n}} \int q^{\setminus n}(\boldsymbol{\theta}) f_n(\boldsymbol{\theta}) d\boldsymbol{\theta} \\
&= \frac{1}{Z_n}\int q^{\setminus n}(\boldsymbol{\theta}) f_n(\boldsymbol{\theta})\left\{-\frac{1}{v^{\setminus n}}(\mathbf{m}^{\setminus n} - \boldsymbol{\theta})\right\} d\boldsymbol{\theta} \\
&= -\frac{\mathbf{m}^{\setminus n}}{v^{\setminus n}} + \frac{\mathbb{E}(\boldsymbol{\theta})}{v^{\setminus n}} \ (\because (10.203)???)
\end{aligned}
上式を整理すると、
\begin{aligned}
\mathbb{E}(\boldsymbol{\theta}) &= \mathbf{m}^{\setminus n} + v^{\setminus n}\nabla_{\mathbf{m}^{\setminus n}} \ln Z_n \\
&= \mathbf{m}^{\setminus n} + v^{\setminus n} \frac{1}{Z_n} (1 - w)\mathcal{N}(\mathbf{x}_n|\mathbf{m}^{\setminus n}, (v^{\setminus n} + 1)\mathbf{I})\cdot\frac{\mathbf{x}_n - \mathbf{m}^{\setminus n}}{v^{\setminus n} + 1} \\
&= \mathbf{m}^{\setminus n} + v^{\setminus n}\rho_n\frac{\mathbf{x}_n - \mathbf{m}^{\setminus n}}{v^{\setminus n} + 1}
\end{aligned}
同様に、
\begin{aligned}
\nabla_{v^{\setminus n}} \ln Z_n &= \frac{1}{Z_n}\nabla_{v^{\setminus n}} \int q^{\setminus n}(\boldsymbol{\theta}) f_n(\boldsymbol{\theta}) d\boldsymbol{\theta} \\
&= \frac{1}{Z_n} \int q^{\setminus n}(\boldsymbol{\theta}) f_n(\boldsymbol{\theta})\left\{\frac{1}{2(v^{\setminus n})^2}(\mathbf{m}^{\setminus n} - \boldsymbol{\theta})^{\mathrm{T}}(\mathbf{m}^{\setminus n} - \boldsymbol{\theta}) - \frac{D}{2v^{\setminus n}}\right\} \\
&= \frac{1}{2(v^{\setminus n})^2}\left\{\mathbb{E}(\boldsymbol{\theta}^{\mathrm{T}}\boldsymbol{\theta}) - 2\mathbb{E}(\boldsymbol{\theta}^{\mathrm{T}})\mathbf{m}^{\setminus n} + \|\mathbf{m}^{\setminus n}\|^2 \right\} - \frac{D}{2v^{\setminus n}}.
\end{aligned}
整理すると、
\mathbb{E}(\boldsymbol{\theta}^{\mathrm{T}}\boldsymbol{\theta}) = 2(v^{\setminus n})^2 \nabla_{v^{\setminus n}} \ln Z_n + 2\mathbb{E}(\boldsymbol{\theta}^{\mathrm{T}})\mathbf{m}^{\setminus n} - \|\mathbf{m}^{\setminus n}\|^2 + Dv^{\setminus n}
となる。また、
\begin{aligned}
\nabla_{v^{\setminus n}} \ln Z_n &= \frac{1}{Z_n} (1 - w)\mathcal{N}(\mathbf{x}_n|\mathbf{m}^{\setminus n}, (v^{\setminus n} + 1)\mathbf{I}) \left(\frac{1}{2(v^{\setminus n} + 1)^2}\|\mathbf{x}_n - \mathbf{m}^{\setminus n}\|^2 - \frac{D}{2(v^{\setminus n} + 1)}\right) \\
&= \rho_n \left(\frac{1}{2(v^{\setminus n} + 1)^2}\|\mathbf{x}_n - \mathbf{m}^{\setminus n}\|^2 - \frac{D}{2(v^{\setminus n} + 1)}\right) \\
\end{aligned}
以上の結果と、
v\mathbf{I} = \mathbb{E}(\boldsymbol{\theta}\boldsymbol{\theta}^{\mathrm{T}}) - \mathbb{E}(\boldsymbol{\theta})\mathbb{E}(\boldsymbol{\theta^{\mathrm{T}}})
を組み合わせることで、(10.218)式を得る。
Discussion