# PRML第10章演習問題解答（10.19〜）
<head>
<style>
div.panel-primary {
border: 1px solid #000;
margin: 10px 5px;
padding: 16px 10px 0px;
}
</style>
</head>
## 演習 10.19
<div class="panel-primary">
ベイズ混合ガウスモデルの変分ベイズ法における予測分布
$$p(\widehat{\mathrm{x}} \mid \mathbf{X}) \simeq \frac{1}{\widehat{\alpha}} \sum_{k=1}^{K} \alpha_{k} \operatorname{St}\left(\widehat{\mathbf{x}} \mid \mathbf{m}_{k}, \mathbf{L}_{k}, \nu_{k}+1-D\right) \tag{10.81}$$
を導出せよ．
</div>
P.197の$(10.81)$式の導出を行うためにこの節での手順を一から踏む。
$$
p(\mathbf{Z} \mid \pi)=\prod_{n=1}^{N} \prod_{k=1}^{K} \pi_{k}^{z_{n k}} \tag{10.37}
$$
$$
p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol{\mu}, \boldsymbol{\Lambda})=\prod_{n=1}^{N} \prod_{k=1}^{K} \mathcal{N}\left(\mathbf{x}_{n} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right)^{z_{n k}} \tag{10.38}
$$
データセット$\mathbf{X}$に対する新しい観測値$\mathbf{\hat{x}}$の予測分布についてこれに対応する潜在分布$\mathbf{\hat{z}}$が存在し、よって予測分布はしたがって以下で与えられる（純粋なベイズの定理・周辺化を用いた）。
$$
p(\widehat{\mathbf{x}} \mid \mathbf{X})=\sum_{\widehat{\mathbf{z}}} \iiint p(\widehat{\mathbf{x}} \mid \widehat{\mathbf{z}}, \boldsymbol{\mu}, \boldsymbol{\Lambda}) p(\widehat{\mathbf{z}} \mid \pi) p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Lambda} \mid \mathbf{X}) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu} \mathrm{d} \boldsymbol{\Lambda} \tag{10.78}
$$
$(10.37)$と$(10.38)$を代入して、
$$
p(\widehat{\mathbf{x}} \mid \mathbf{X})=\sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Lambda} \mid \mathbf{X}) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu} \mathrm{d} \boldsymbol{\Lambda} \tag{10.79}
$$
となる。真の事後分布$p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Lambda} \mid \mathbf{X})$を変分近似で置き換える。このとき
$$
q(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Lambda})=q(\boldsymbol{\pi}) \prod_{j=1}^{K} q\left(\boldsymbol{\mu}_{j}, \boldsymbol{\Lambda}_{j}\right) \tag{10.55}
$$
を用いて、和$\displaystyle \sum_{k=1}^{K}$のうち1つの項に注目する（＝$k$を固定する）。$j\neq k$であるような$j$についての$\displaystyle \int d\mu_{j}\int d\Lambda_{j}$を考えると、積分の中身は$\displaystyle \int q(\boldsymbol{\mu}_{j}, \boldsymbol{\Lambda}_{j})d\boldsymbol{\mu}_{j}d\boldsymbol{\Lambda}_{j} = 1$（確率の定義より）となるので、$k$番目の積分$\displaystyle \int d\mu_{k}\int d\Lambda_{k}$しか残らない。これより
$$
\begin{aligned}
p(\widehat{\mathbf{x}} \mid \mathbf{X}) &=\sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) p(\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\Lambda} \mid \mathbf{X}) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu} \mathrm{d} \boldsymbol{\Lambda} \\
&\simeq \sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) q(\boldsymbol{\pi}) q\left(\boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}\right) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu}_{k} \mathrm{~d} \boldsymbol{\Lambda}_{k} \\
&=\sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) \underbrace{\operatorname{Dir}(\boldsymbol{\pi}\mid \boldsymbol{\alpha})}_{(10.57)} \underbrace{\mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \boldsymbol{\Lambda}_{k}\right)^{-1}\right) \mathcal{W}\left(\boldsymbol{\Lambda}_{k} \mid \mathbf{W}_{k}, \boldsymbol{\nu}_{k}\right)}_{(10.59)} \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu}_{k} \mathrm{d} \boldsymbol{\Lambda}_{k} \\
&=\sum_{k=1}^{K} \int\pi_{k} \operatorname{Dir}(\boldsymbol{\pi}\mid \boldsymbol{\alpha}) \mathrm{d} \boldsymbol{\pi} \int \left[ \int \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) \mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \boldsymbol{\Lambda}_{k}\right)^{-1}\right) \mathrm{d} \boldsymbol{\mu}_{k} \right]\mathcal{W}\left(\boldsymbol{\Lambda}_{k} \mid \mathbf{W}_{k}, \boldsymbol{\nu}_{k}\right) \mathrm{d} \boldsymbol{\Lambda}_{k} \cdots (\textrm{A})\\
\end{aligned}
$$
$\boldsymbol{\pi}$の積分に関係するのは$\pi_{k}\operatorname{Dir}(\boldsymbol{\pi}\mid \boldsymbol{\alpha})$のみで、$\int\pi_{k} \operatorname{Dir}(\boldsymbol{\pi}\mid \boldsymbol{\alpha}) \mathrm{d} \boldsymbol{\pi}$はディリクレ分布以下での$\pi_{k}$の期待値であるから
$$
\int\pi_{k} \operatorname{Dir}(\boldsymbol{\pi}\mid \boldsymbol{\alpha}) \mathrm{d} \boldsymbol{\pi} = \frac{\alpha_{k}}{\widehat{\alpha}}\ (\because{(B.17)})
$$
次に$\displaystyle \int \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) \mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \boldsymbol{\Lambda}_{k}\right)^{-1}\right) \mathrm{d} \boldsymbol{\mu}_{k}$について、これを
$$
\begin{aligned}
p\left(\boldsymbol{\mu}_{k}\right)&=\mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \boldsymbol{\Lambda}_{k}\right)^{-1}\right) \\
p\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}\right)&=\mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right)
\end{aligned}
$$
とみなして$(2.115)$の公式を用いると
$$
\begin{aligned}
& \int \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}^{-1}\right) \mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k},\left(\beta_{k} \boldsymbol{\Lambda}_{k}\right)^{-1}\right) d \boldsymbol{\mu}_{k} \\
=&\ \mathcal{N}\left(\widehat{\mathbf{x}} \mid \mathbf{m}_{k},\left(\boldsymbol{\Lambda}_{k}^{-1}+\beta_{k}^{-1} \boldsymbol{\Lambda}_{k}^{-1}\right)\right) \\
=&\ \mathcal{N}\left(\widehat{\mathbf{x}} \mid \mathbf{m}_{k},\left(1+\beta_{k}^{-1}\right) \boldsymbol{\Lambda}_{k}^{-1}\right) \end{aligned}
$$
となる。以上から$(\textrm{A})$式に戻ると
$$
\begin{aligned} p(\widehat{\mathbf{x}} \mid \mathbf{X})\simeq & \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \int \mathcal{N}\left(\widehat{\mathbf{x}} \mid \mathbf{m}_{k},\left(1+\beta_{k}^{-1}\right) \boldsymbol{\Lambda}_{k}^{-1}\right) \mathcal{W}\left(\boldsymbol{\Lambda}_{k} \mid \mathbf{W}_{k}, \nu_{k}\right) \mathrm{d} \boldsymbol{\Lambda}_{k} \\
=&\ \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \int \frac{1}{(2 \pi)^{D / 2}} \frac{\left|\boldsymbol{\Lambda}_{k}\right|^{1 / 2}}{\left(1+\beta_{k}^{-1}\right)^{D / 2}} \exp \left\{-\frac{\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{k}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)}{2\left(1+\beta_{k}^{-1}\right)}\right\} \\ & B\left(\mathbf{W}_{k}, \nu_{k}\right)\left|\boldsymbol{\Lambda}_{k}\right|^{\left(\nu_{k}-D-1\right) / 2} \exp \left\{-\frac{1}{2} \operatorname{Tr}\left[\mathbf{W}_{k}^{-1} \boldsymbol{\Lambda}_{k}\right]\right\} \mathrm{d} \boldsymbol{\Lambda}_{k} \\
=&\ \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \int \frac{B\left(\mathbf{W}_{k}, \nu_{k}\right)}{(2 \pi)^{D / 2}} \frac{\left|\boldsymbol{\Lambda}_{k}\right|^{\left((\nu_{k}+1)-D-1\right) / 2}}{\left(1+\beta_{k}^{-1}\right)^{D / 2}} \exp \left\{-\frac{\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{k}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)}{2\left(1+\beta_{k}^{-1}\right)}\right\} \exp \left\{-\frac{1}{2} \operatorname{Tr}\left[\mathbf{W}_{k}^{-1} \boldsymbol{\Lambda}_{k}\right]\right\} \mathrm{d} \boldsymbol{\Lambda}_{k} \\
=&\ \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \int \frac{B\left(\mathbf{W}_{k}, \nu_{k}\right)}{(2 \pi)^{D / 2}} \frac{\left|\boldsymbol{\Lambda}_{k}\right|^{\left((\nu_{k}+1)-D-1\right) / 2}}{\left(1+\beta_{k}^{-1}\right)^{D / 2}} \exp \left\{-\frac{1}{2} \operatorname{Tr}\left[\left(\frac{\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}}}{1+\beta_{k}^{-1}}+\mathbf{W}_{k}^{-1}\right) \boldsymbol{\Lambda}_{k}\right]\right\} \mathrm{d} \boldsymbol{\Lambda}_{k} \\
=&\ \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \frac{B\left(\mathbf{W}_{k}, \nu_{k}\right)}{(2 \pi)^{D / 2} \left(1+\beta_{k}^{-1}\right)^{D / 2}}
\int \left|\boldsymbol{\Lambda}_{k}\right|^{\left((\nu_{k}+1)-D-1\right) / 2} \exp \left\{-\frac{1}{2} \operatorname{Tr}\left[\left(\frac{\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}}}{1+\beta_{k}^{-1}}+\mathbf{W}_{k}^{-1}\right) \boldsymbol{\Lambda}_{k}\right]\right\} \mathrm{d} \boldsymbol{\Lambda}_{k}
\end{aligned}
$$
ここで、$\int$の中身は
$$
\begin{aligned}
\mathbf{W^{\prime}}_{k}^{-1} &=\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} + \mathbf{W}_{k}^{-1} \\
{\nu^{\prime}}_{k} &= \nu_{k}+1
\end{aligned}
$$
としたときのウィシャート分布$\mathcal{W}(\boldsymbol{\Lambda}_{k}\mid \mathbf{W^{\prime}}_{k}, {\nu^{\prime}}_{k})$となっているので、この積分結果は正規化定数である$B(\mathbf{W^{\prime}}_{k}, {\nu^{\prime}}_{k})$の逆数になることがわかる。すなわち
$$
p(\widehat{\mathbf{x}} \mid \mathbf{X}) \simeq \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \frac{1}{(2 \pi)^{D / 2} \left(1+\beta_{k}^{-1}\right)^{D / 2}}\frac{B\left(\mathbf{W}_{k}, \nu_{k}\right)}{B(\mathbf{W^{\prime}}_{k}, {\nu^{\prime}}_{k})} \tag{B}
$$
となる。この正規化定数部分をさらに展開していく。
$$
\begin{aligned}
\frac{B(\mathbf{W}_k,\nu_k)}{B(\mathbf{W^{\prime}}_{k},\nu_k+1)}
&=\frac{\left|\mathbf{W}_{k}\right|^{-\frac{\nu_{k}}{2}}\left(2^{\frac{\nu_{k} D}{2}} \pi^{\frac{D(D-1)}{4}} \prod_{i=1}^{D} \Gamma\left(\frac{\nu_{k}+1-i}{2}\right)\right)^{-1}}{\left|\mathbf{W}^{\prime}_{k}\right|^{-\frac{\nu_{k}+1}{2}}\left(2^{\frac{\left(\nu_{k}+1\right) D}{2}} \pi^{\frac{D(D-1)}{4}} \prod_{i=1}^{D} \Gamma\left(\frac{\nu_{k}+2-i}{2}\right)\right)^{-1}} ~~~(\because\ (B.79)) \\
&=\frac{\left|\mathbf{W}_{k}\right|^{-\frac{\nu_{k}}{2}}}{\left|\mathbf{W}^{\prime}_{k}\right|^{-\frac{\nu_{k}+1}{2}}} 2^{\frac{D}{2}} \frac{\prod_{i=1}^{D} \Gamma\left(\frac{\nu_{k}+2-i}{2}\right)}{\prod_{i=1}^{D} \Gamma\left(\frac{\nu_{k}+1-i}{2}\right)} \\
&=2^{D/2}\frac{\left|\mathbf{W}_{k}\right|^{-\frac{\nu_{k}}{2}}}{\left|\left\{\mathbf{W}_{k}^{-1}+\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}}\right\}^{-1}\right|^{-\frac{\nu_{k}}{2}}} \\
&~~~~\frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right) \Gamma\left(\frac{\nu k}{2}\right) \Gamma\left(\frac{\nu_{k}-1}{2}\right) \cdots \Gamma\left(\frac{\nu_{k}+2-D}{2}\right)}{\Gamma\left(\frac{\nu k}{2}\right) \Gamma\left(\frac{\nu_{k}-1}{2}\right) \cdots \Gamma\left(\frac{\nu_{k}+2-D}{2}\right) \Gamma\left(\frac{\nu_{k}+1-D}{2}\right)} \\
&=2^{D/2}\left|\mathbf{W}_{k}\right|^{-\frac{\nu_{k}}{2}}\left|\mathbf{W}_{k}^{-1}\left\{\mathbf{I}+\mathbf{W}_{k}\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}}\right\}\right|^{-\frac{\nu_{k}+1}{2}}\frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)} \\
&=2^{D/2}\left|\mathbf{W}_{k}\right|^{1/2}\left|\mathbf{I}+\mathbf{W}_{k}\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}}\right|^{-\frac{\nu_{k}+1}{2}}\frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)} \\
&=2^{D/2}\left|\mathbf{W}_{k}\right|^{1/2}\left[1+\left\{\mathbf{W}_{k}\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\right\}^{\mathrm{T}}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\right]^{-\frac{\nu_{k}+1}{2}}\frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)} ~~ (\because (\textrm{C}.15))\\
&=2^{D/2}\left|\mathbf{W}_{k}\right|^{1/2}\left\{1+\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \mathbf{W}_{k}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\right\}^{-\frac{\nu_{k}+1}{2}}\frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)}~~ (\because \mathbf{W}_{k}^{\mathrm T} = \mathbf{W}_{k})
\end{aligned}
$$
これを$(\textrm{B})$に代入して
$$
\begin{aligned}
p(\widehat{\mathbf{x}} \mid \mathbf{X}) &\simeq \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \frac{\Gamma\left(\frac{\nu_{k}+1}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)}\frac{\left|\mathbf{W}_{k}\right|^{1/2}}{\pi^{D / 2} \left(1+\beta_{k}^{-1}\right)^{D / 2}}\left\{1+\left(1+\beta_{k}^{-1}\right)^{-1}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \mathbf{W}_{k}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\right\}^{-\frac{\nu_{k}+1}{2}} \\
&= \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \frac{\Gamma\left(\frac{\nu_{k}+1-D}{2} + \frac{D}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)}
\frac{\left|\frac{\nu_{k}+1-D}{1+\beta_{k}^{-1}}\mathbf{W}_{k}\right|^{1/2}}{\pi^{D / 2} \left(\nu_{k}+1-D\right)^{D / 2}} \\
&~~~~\left\{1+\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \left( \frac{1}{\nu_{k}+1-D}\frac{\nu_{k}+1-D}{1+\beta_{k}^{-1}}\mathbf{W}_{k} \right)\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)\right\}^{-\frac{\nu_{k}+1-D}{2} - \frac{D}{2}} \\
&= \sum_{k=1}^{K} \frac{\alpha_{k}}{\widehat{\alpha}} \frac{\Gamma\left(\frac{\nu_{k}+1-D}{2} + \frac{D}{2}\right)}{\Gamma\left(\frac{\nu_{k}+1-D}{2}\right)} \frac{\left|\mathbf{L}_{k}\right|^{1/2}}{\left\{\pi (\nu_{k}+1-D)\right\}^{D/2}}\left( 1 + \frac{\Delta^{2}}{\nu_{k}+1-D}\right)^{-\frac{\nu_{k}+1-D}{2} - \frac{D}{2}} \\
&= \frac{1}{\widehat{\alpha}}\sum_{k=1}^{K}\alpha_{k}\operatorname{St} \left( \widehat{\mathbf{x}} \mid \mathbf{m}_{k}, \mathbf{L}_{k}, \nu_{k}+1-D \right) ~~ (\because (\textrm{B}.68))
\end{aligned}
$$
となる。ここで、
$$
\mathbf{L}_{k} =\frac{\nu_{k}+1-D}{1+\beta_{k}^{-1}} \mathbf{W}_{k} = \frac{(\nu_{k}+1-D)\beta_{k}}{(1+\beta_{k})} \mathbf{W}_{k} \tag{10.82}
$$
$$
\Delta^{2} =\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)^{\mathrm{T}} \mathbf{L}_{k}\left(\widehat{\mathbf{x}}-\mathbf{m}_{k}\right)
$$
とした。これより$(10.81)$を得た。
## 演習 10.20
<div class="panel-primary">
この演習問題では，データ集合のサイズ$N$が大きくなった場合の混合ガウスモデルの変分ベイズ法による解を考え，これが(期待通り)9章のEMアルゴリズムに基づく最尤推定の解に近づくことを示す．この演習問題を解くには，付録Bの結果が有用であろう．最初に，精度の事後分布$q^{\star}(\mathbf{\Lambda}_k)$が最尤推定値の周囲に鋭い分布を持つことを示せ．平均の事後分布$q^{\star}(\boldsymbol{\mu}_k \mid \mathbf{\Lambda}_k)$についても同様のことを示せ．次に，混合比の事後分布$q^{\star}(\boldsymbol{\pi})$について考え，これも最尤推定値の周囲に鋭く分布することを示せ．同様に，大きな$N$については負担率は対応する最尤推定値と等しくなることを，大きな$x$についてのディガンマ関数の次の漸近的な結果
$$
\psi(x)=\ln x+O(1 / x)
$$
を利用して示せ．最後に
$$
p(\widehat{\mathbf{x}} \mid \mathbf{X}) \simeq \sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \mathbf{\Lambda}_{k}^{-1}\right) q(\boldsymbol{\pi}) q\left(\boldsymbol{\mu}_{k}, \mathbf{\Lambda}_{k}\right) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu}_{k} \mathrm{~d} \mathbf{\Lambda}_{k} \tag{10.80}
$$
を用いて，大きな$N$については予測分布は混合ガウス分布になることを示せ．
</div>
(10.59)式の導出を考えると、(演習10.13より)
$$
q^{\star}\left(\mathbf{\Lambda}_{k}\right)=\mathcal{W}\left(\mathbf{\Lambda}_{k} \mid \mathbf{W}_{k}, \nu_{k}\right)\\
q^{\star}\left(\boldsymbol{\mu}_{k} \mid \boldsymbol{\Lambda}_{k}\right)=\mathcal{N}\left(\boldsymbol{\mu}_{k} \mid \mathbf{m}_{k}, \beta_{k} \boldsymbol{\Lambda}_{k}\right)
$$
となる。
これらの分布について、$N \rightarrow \infty$のとき、
$N_{k} \rightarrow \infty$であり、
(10.60)~(10.63)式より、
$$
\beta_{k} \rightarrow N_{k}\\
\mathbf{m}_{k} \rightarrow \overline{\mathrm{x}}_{k}\\
\mathbf{W}_{k} \rightarrow N_{k}^{-1} \mathbf{S}_{k}^{-1}\\
\nu_{k} \rightarrow N_{k}
$$
である。
これらと(B.79)~(B.81)式より,
$$
\mathrm{E}\left[\boldsymbol{\Lambda}_{k}\right]=\nu_{k} \mathbf{W}_{k} \rightarrow \mathbf{S}_{k}^{-1}
$$
$$
\begin{aligned}
-\ln B\left(\mathbf{W}_{k}, \nu_{k}\right)&=-\ln (|\mathbf{W}_{k}|^{-\nu_{k} / 2}\left(2^{\nu_{k} D / 2} \pi^{D(D-1) / 4} \prod_{i=1}^{D} \Gamma\left(\frac{\nu_{k}+1-i}{2}\right)\right)^{-1})\\
&\rightarrow-\frac{N_{k}}{2}\left(D \ln N_{k}+\ln \left|\mathbf{S}_{k}\right|-D \ln 2\right)+\sum_{i=1}^{D} \ln \Gamma\left(\frac{N_{k}+1-i}{2}\right)\\
&\rightarrow-\frac{N_{k}}{2}\left(D \ln N_{k}+\ln \left|\mathbf{S}_{k}\right|-D \ln 2\right)+\sum_{i=1}^{D} \frac{N_{k}}{2}\left(\ln N_{k}-\ln 2-1\right)~~ (\because (\textrm1.146))\\
& \rightarrow-\frac{N_{k} D}{2}\left(\ln N_{k}-\ln 2-\ln N_{k}+\ln 2+1\right)-\frac{N_{k}}{2} \ln \left|\mathbf{S}_{k}\right| \\
&=-\frac{N_{k}}{2}\left(\ln \left|\mathbf{S}_{k}\right|+D\right)
\end{aligned}
$$
$$
\begin{aligned}
\mathbb{E}[\ln |\boldsymbol{\Lambda}_{k}|] &=\sum_{i=1}^{D} \psi\left(\frac{\nu_{k}+1-i}{2}\right)+D \ln 2+\ln |\mathbf{W}_{k}|\\
& \rightarrow D \ln \frac{N_{k}}{2}+D \ln 2-D \ln N_{k}-\ln \left|\mathbf{S}_{k}\right| \\
&=-\ln \left|\mathbf{S}_{k}\right|
\end{aligned}
$$
ただし、$\psi(\cdot)$は(10.241)式:
$$\psi(x)=\ln x+O(1 / x)$$
のディガンマ分布。
よって(B.82)式より
$$
\begin{aligned}
\mathrm{H}[\boldsymbol{\Lambda}_{k}]&=-\ln B(\mathbf{W}_{k}, \nu_{k})-\frac{(\nu_{k}-D-1)}{2} \mathbb{E}[\ln |\boldsymbol{\Lambda}_{k}|]+\frac{\nu_{k} D}{2}\\
& \rightarrow 0
\end{aligned}
$$
これにより $q^{\star}\left(\boldsymbol{\Lambda}_{k}\right)$ については示された。
また、
$$
\mathbf{m}_{k} \rightarrow \overline{\mathrm{x}}_{k}\\
\beta_{k} \mathbf{\Lambda}_{k} \rightarrow \beta_{k} \nu_{k} \mathbf{W}_{k} \rightarrow N_{k} \mathbf{S}_{k}^{-1}
$$
より $q^{\star}\left(\boldsymbol{\mu}_{k} \mid \boldsymbol{\Lambda}_{k}\right)$ についても示された。
$q^{\star}(\pi)$については、
(10.56),(10.57)式にある通り
$$
q^{\star}(\pi)=\operatorname{Dir}(\pi \mid \alpha)\\
\alpha_{k}=\alpha_{0}+N_{k}
$$
であり、$\alpha_{k} \rightarrow N_{k}$である。
(B.17),(B.19)式より、
$$
\begin{aligned}
\mathbb{E}\left[\pi_{k}\right]&=\frac{\alpha_{k}}{\overline{\alpha}}\\
&\rightarrow \frac{N_{k}}{N}
\end{aligned}
$$
$$
\begin{aligned}
\operatorname{cov}\left[\mu_{j} \mu_{k}\right]&=-\frac{\alpha_{j} \alpha_{k}}{\widehat{\alpha}^{2}(\widehat{\alpha}+1)}\\
&\rightarrow 0
\end{aligned}
$$
よって$q^{\star}(\pi)$についても示された。
最後に(10.80)式より、
$$
\begin{aligned}
p(\widehat{\mathbf{x}} \mid \mathbf{X}) &\simeq \sum_{k=1}^{K} \iiint \pi_{k} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \mathbf{\Lambda}_{k}^{-1}\right) q(\boldsymbol{\pi}) q\left(\boldsymbol{\mu}_{k}, \mathbf{\Lambda}_{k}\right) \mathrm{d} \boldsymbol{\pi} \mathrm{d} \boldsymbol{\mu}_{k} \mathrm{~d} \mathbf{\Lambda}_{k}\\
&\rightarrow \sum_{k=1}^{K} \frac{\alpha_{k}}{\bar{\alpha}} \iint \mathcal{N}\left(\widehat{\mathbf{x}} \mid \boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}\right) q\left(\boldsymbol{\mu}_{k}, \boldsymbol{\Lambda}_{k}\right) \mathrm{d} \boldsymbol{\mu}_{k} \mathrm{~d} \boldsymbol{\Lambda}_{k}\\
&\rightarrow \sum_{k=1}^{K} \frac{N_{k}}{N} \mathcal{N}\left(\widehat{\mathbf{x}} \mid \overline{\mathbf{x}}_{k}, \mathbf{W}_{k}\right)
\end{aligned}
$$
ただし最後の行は$q^{\star}\left(\boldsymbol{\Lambda}_{k}\right)$と$q^{\star}\left(\boldsymbol{\mu}_{k} \mid \boldsymbol{\Lambda}_{k}\right)$が特定の位置についてのデルタ関数と近似した。
これにより示された。
## 演習 10.21
<div class="panel-primary">
$K$個の混合要素を持つ混合モデルにおいて，混合要素の入れ替えについての対称性から得られる，同値なパラメータ設定の数は$K!$であることを示せ．
</div>
P.197によれば
> 例として，一つの観測値$x$についての二混合のガウス混合分布を考えよう．パラメータの値は$\pi_{1} = a$，$\pi_{2} = b$，$\pi_{3} = c$，$\pi_{4} = d$，$\pi_{5} = e$，$\pi_{6} = f$とする．このとき，二つの混合要素を入れ替えた別の設定$\pi_{1} = b$，$\pi_{2} = a$，$\pi_{3} = d$，$\pi_{4} = c$，$\pi_{5} = f$，$\pi_{6} = e$も，対称性から同じ$p(x)$を与える．
とあるように、もし$K$個の混合要素が存在する場合は、それらを入れ替えることで同値なパラメータ設定が可能なので、一般に$K!$個存在することは明らかである。
## 演習 10.22
<div class="panel-primary">
これまでにガウス混合モデルの事後分布の持つそれぞれの峰は，$K!$個ある同値な峰の一つであることを見てきた．変分ベイズ推論のアルゴリズムを実行した結果，近似事後分布$q$がどれかの峰の周りに局所化して得られたとしよう．このとき，完全な事後分布はこうした分布$q$の$K!$個の混合分布となり，各混合要素が峰となって同じ混合係数を持つ．この混合分布$q$の混合要素の間の重なりが無視できる程度だと仮定すると，結果として得られる全体の下界は，$q$の一つの混合要素の下界に項$\ln K!$を加えたものになることを示せ．
</div>
今、p166並びに、演習9.24より
\begin{align}
\ln p(\mathbf{X}) = L(q) + KL(q||p)
\end{align}
が成り立つ。
この時、$KL(q||p)$はKLダイバージェンスであり、1.6.1の議論から、$KL(q||p) \geq 0$である。よって、$\ln p(\mathbf{X}) \geq L(q)$であり、$L(q)$は、$\ln p(\mathbf{X})$の下界である。よって、本題は、求めたい真の分布を$p(\mathbf{Z}|\mathbf{X})$として、$L(p(\mathbf{Z}|\mathbf{X}))$を求めれば良い。
まず、各峰は、pの真の各峰を$r_i$ $(i \in \{1, 2... K!\})$とおくと、pは単純に各$r_i$の平均で表すことができる。
\begin{align}
p(\mathbf{Z}|\mathbf{X}) \simeq \sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})
\end{align}
また、各峰の重なりが無視できるという仮定から、$r_k \neq 0 \rightarrow r_{i \neq k = 0}=0$が成り立つ。
ここで、ある真の峰$r_k$の近似の峰をqとおく。すなわち、この問題では、$p$を混合要素$r_i$を重ね合わせたものと見なし、その近似であるqによって下界を表すことを目指す。
すると、10.3式から、
\begin{align}
L(p(\mathbf{Z}|\mathbf{X})) &= \int p(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{p(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} \\
&= \int\sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} \\
&= \frac{1}{K!} \int r_1(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} +
\frac{1}{K!} \int r_2(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} + \cdots
\frac{1}{K!} \int r_{K!}(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\sum_i^{K!} \frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} \\
&= \frac{1}{K!} \int r_1(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\frac{1}{K!} r_1(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} +
\frac{1}{K!} \int r_2(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{ \frac{1}{K!} r_2(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} + \cdots
\frac{1}{K!} \int r_{K!}(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\frac{1}{K!} r_i{K!}\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} &\because r_k \neq 0 \rightarrow r_{i \neq k = 0}=0 \\
&= \frac{1}{K!} \sum_i^{K!} \int r_i(\mathbf{Z}|\mathbf{X}) \ln \{\frac{p(\mathbf{Z}, \mathbf{X})}{\frac{1}{K!} r_i(\mathbf{Z}|\mathbf{X})} \} d\mathbf{Z} \\
&= \frac{1}{K!} \sum_i^{K!} \int r_i(\mathbf{Z}|\mathbf{X}) \{ \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_i(\mathbf{Z}|\mathbf{X})} + \ln K!\} d\mathbf{Z} \\
&= \frac{1}{K!} \{ \sum_i^{K!} \int r_i(\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_i(\mathbf{Z}|\mathbf{X})} d\mathbf{Z} +
\ln K! \sum_i^{K!} \int r_i(\mathbf{Z}|\mathbf{X}) d\mathbf{Z}\}
\\
&= \frac{1}{K!} \sum_i^{K!} \int r_i(\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_i(\mathbf{Z}|\mathbf{X})} d\mathbf{Z} +
\ln K! &\because \int r_i(\mathbf{Z}|\mathbf{X}) d\mathbf{Z} = r_i(\mathbf{X}|\mathbf{X}) = 1\\
&=\frac{1}{K!} \{ \int r_1 (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_1(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}+
\int r_2 (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_2(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}
\}+\cdots +
\int r_k (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_k(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}
\}+\cdots +
\int r_{K!} (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_{K!}(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}
\}+
\ln K! \\
&=\frac{1}{K!} K! \int q(\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ q(\mathbf{Z}|\mathbf{X})} d\mathbf{Z} +
\ln K! &\because r_i\text{はそれぞれ同値であり、積分は同じ。詳細は最後 }\\
&= L(q) + \ln K! &\because (10.3)
\end{align}
今、$L(q)$は一つの混合要素の下界なので、題意は満たされた。
最後から2番目の式変形について、まず、自明に、$\int r_k (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_k(\mathbf{Z}|\mathbf{X})} d\mathbf{Z} =\int q(\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ q(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}$が成り立つ。
そして、$i \neq k$について、
\begin{align}
\int r_i (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ r_i(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}
&= \int r_i (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{Z}| \mathbf{X})d\mathbf{Z}
+\int r_i (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{X}) d\mathbf{Z}
-\int r_i (\mathbf{Z}|\mathbf{X}) \ln r_i (\mathbf{Z}|\mathbf{X}) d\mathbf{Z} \\
&= \int r_i (\mathbf{Z}|\mathbf{X}) \ln \frac{r_i(\mathbf{Z}| \mathbf{X})}{K!} d\mathbf{Z}
+\int r_i (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{X}) d\mathbf{Z}
-\int r_i (\mathbf{Z}|\mathbf{X}) \ln r_i (\mathbf{Z}|\mathbf{X}) d\mathbf{Z} &\because r_i(\mathbf{Z} \notin \mathbf{Z}_i |\mathbf{X}) = 0\\
&= \int r_k (\mathbf{Z}|\mathbf{X}) \ln \frac{r_k(\mathbf{Z}| \mathbf{X})}{K!} d\mathbf{Z}
+\int r_k (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{X}) d\mathbf{Z}
-\int r_k (\mathbf{Z}|\mathbf{X}) \ln r_k (\mathbf{Z}|\mathbf{X}) d\mathbf{Z} &\because r_i\text{はそれぞれ同値であり、積分は同じ}\\
&= \int r_k (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{Z}| \mathbf{X})d\mathbf{Z}
+\int r_k (\mathbf{Z}|\mathbf{X}) \ln p(\mathbf{X}) d\mathbf{Z}
-\int r_k (\mathbf{Z}|\mathbf{X}) \ln r_k (\mathbf{Z}|\mathbf{X}) d\mathbf{Z} &\because r_k(\mathbf{Z} \notin \mathbf{Z}_k |\mathbf{X}) = 0\\
&= \int r_k (\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{r_k (\mathbf{Z}|\mathbf{X})} d\mathbf{Z}\\
&= \int q(\mathbf{Z}|\mathbf{X}) \ln \frac{p(\mathbf{Z}, \mathbf{X})}{ q(\mathbf{Z}|\mathbf{X})} d\mathbf{Z}
\end{align}
最後は冗長かもしれないのでその時はご教示ください。
## 演習 10.23
<div class="panel-primary">
混合係数$\{ \pi_k \}$に事前分布を与えない変分ベイズガウス混合モデルを考えよう．代わりに混合係数はパラメータとして扱い，対数周辺尤度の下界を最大化する際に値を求める．ラグランジュ乗数法を用いて，混合係数の和が$1$になる制約条件の下でこの下界を混合係数について最大化すると，再推定式
$$\pi_{k}=\frac{1}{N} \sum_{n=1}^{N} r_{n k} \tag{10.83}$$
の結果が得られることを示せ．この際，下界のすべての項を考える必要はなく，$\{\pi_k\}$に依存する項だけを考えればよいことに注意せよ．
</div>
変分ベイズガウス混合モデルでは、下界は(10.70)式で与えられる。
本問では、**混合係数$\{ \pi_k \}$に事前分布を与えない**パラメータとして扱うため、対数周辺尤度として第２項のみを考えれば良い。つまり、
$$
\mathscr{L} \propto \mathbb{E}[\ln p(\mathbf{Z} \mid \pi)]=\sum_{n=1}^{N} \sum_{k=1}^{K} r_{n k} \ln \pi_{k} \tag{10.72}
$$
$$
L=\sum_{n=1}^{N} \sum_{k=1}^{K} r_{n k} \ln \pi_{k}+\lambda \cdot\left(\sum_{k=1}^{K} \pi_{k}-1\right)
$$
上式のLagrangianについて、$\pi_k$について微分し、=0とおくと、
$$
\frac{\partial L}{\partial \pi_{k}}=\frac{\sum_{n=1}^{N}r_{nk}}{\pi_{k}}+\lambda=\frac{N_{k}}{\pi_{k}}+\lambda=0 \tag{A}
$$
(A)の両辺に$\pi_k$をかけ、$\sum_{k=1}^{K}$をとると、
$$
\sum_{k=1}^{K} N_{k}+\lambda \sum_{k=1}^{K}\pi_{k}=0
$$
$\sum_{k=1}^{K} N_{k}=N$, $\sum_{k=1}^{K}\pi_{k}=1$より、$$\lambda=-N$$(A)に代入し、$\pi_{k}$について解くと、$$\pi_{k}=\frac{N_{k}}{N}=\underline{\frac{1}{N} \sum_{n=1}^{N} r_{n k}}$$
## 演習 10.24
<div class="panel-primary">
10.2節でガウス混合モデルを最尤推定で扱う際に現れる特異性は，ベイズ的な解では現れないことを見た．こうした特異性は，ベイズモデルを最大事後機率(MAP)推定を使って解く際には現れるかどうか議論せよ．
</div>
最尤推定で現れる特異性とは、9.2.1節で議論した$\left|\boldsymbol{\Lambda}_{k}\right| \rightarrow \infty$に発散してしまうことを意味している。ベイズモデルではこのようなことが起きないことを示す。
混合ガウス分布の事後確率は、(10.9),(10.38),(10.40),(10.50)を利用すれば以下となる。
$$
\begin{aligned}
\mathbb{E}_{q(\mathbf{Z})} &[\ln p(\mathbf{X} \mid \mathbf{Z}, \boldsymbol{\mu}, \boldsymbol{\Lambda}) p(\boldsymbol{\mu}, \mathbf{\Lambda})] \\
=& \frac{1}{2} \sum_{n=1}^{N} r_{k n}\left(\ln \left|\boldsymbol{\Lambda}_{k}\right|-\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{k}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{k}\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{k}\right)\right) \\
&+\ln \left|\boldsymbol{\Lambda}_{k}\right|-\beta_{0}\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right)^{\mathrm{T}} \boldsymbol{\Lambda}_{k}\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right) \\
&+\left(\nu_{0}-D-1\right) \ln \left|\boldsymbol{\Lambda}_{k}\right|-\operatorname{Tr}\left[\mathbf{W}_{0}{ }^{1} \boldsymbol{\Lambda}_{k}\right]+\text { const. }
\end{aligned}
$$
これを（10.51)-(10.53)を利用して$\mathbf{\Lambda}_{k}$について整理すると（$\mathbf{\Lambda}_{k}$と無関係な項は無視）
$$
\left(\nu_{0}+N_{k}-D\right) \ln \left|\boldsymbol{\Lambda}_{k}\right|-\operatorname{Tr}\left[\left(\mathbf{W}_{0}^{-1}+\beta_{0}\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right)\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right)^{\mathrm{T}}+N_{k} \mathbf{S}_{k}\right) \boldsymbol{\Lambda}_{k}\right]
$$
(C.24),(C.28)を使用して$\mathbf{\Lambda}_{k}$で微分してゼロとおくと以下になる。
$$
\mathbf{\Lambda}_{k}^{-1}=\frac{1}{\nu_{0}+N_{k}-D}\left(\mathbf{W}_{0}^{-1}+\beta_{0}\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right)\left(\boldsymbol{\mu}_{k}-\mathbf{m}_{0}\right)^{\mathrm{T}}+N_{k} \mathbf{S}_{k}\right)
$$
ウィシャート分布の正式より、$\left|\boldsymbol{\Lambda}_{k}^{-1}\right|$はゼロになることはない。
よって、ベイズモデルでは、$\left|\boldsymbol{\Lambda}_{k}\right| \rightarrow \infty$に発散することはないことが示された。
## 演習 10.25
<div class="panel-primary">
10.2節で議論したベイズ混合ガウス分布の変分ベイス法による解では，事後分布について分解した近似
$$q(\mathbf{Z})=\prod_{i=1}^{M} q_{i}\left(\mathbf{Z}_{i}\right) \tag{10.5}$$
を用いた図10.2で見たように，こうした分解の仮定はパラメータ空間で、の事後分布の特定の方向の分散を過小評価してしまう．この影響がモデルエピデンスの変分近似に及ぼす影響について質的に議論せよ．さらに，この影響が混合モデルの混合要素数に関してどう変わるか述べよ．これから，変分ガウス混合モデルが最適な混合要素数を過小評価しがちか，過大評価しがちか説明せよ．
</div>
混合成分の数が増えると、相関している可能性のある変数の数も増える一方、平均場近似の式(10.5)を用いるとそれらの相関を表現することができない(図10.2,3)。その結果、KLダイバージェンスの最小化を行うときに複数の山をつぶして近似してしまうことが考えられるため、過小評価する。
## 演習 10.26
<div class="panel-primary">
ベイズ線形回帰モデルの変分ベイズ法による解法を拡張し，$\beta$についてガンマ超事前分布$\textrm{Gam}(\beta\mid c_0, d_0)$を導入して，分解された変分事後分布$q(\mathbf{w}) q(\alpha) q(\beta)$を仮定して変分ベイズ法によって解け．変分事後分布の三つの因子の更新式を導出し，さらに下界および予測分布の式を求めよ．
</div>
※
$\beta$を含めた全ての変数の同時分布は
$$
p(\mathbf{t}, \mathbf{w}, \alpha, \beta)=p(\mathbf{t}|\mathbf{w}, \beta)p(\mathbf{w}|\alpha)p(\alpha)p(\beta)
$$
と書くことができる．本文中の議論をなぞって，$\mathbf{w}, \alpha, \beta$の尤度関数と事前分布を
$$
\begin{aligned}
&p(\mathbf{t} \mid \mathbf{w}, \beta, \mathbf{X})=\prod_{n=1}^{N} N\left(\mathbf{t}_{n} \mid \mathbf{w}^{\mathrm{T}} \phi_{n}, \beta^{-1}\right) \\
&p(\mathbf{w} \mid \alpha)=N\left(\mathbf{w} \mid 0, \alpha^{-1} \mathbf{I}\right) \\
&p(\alpha)=\operatorname{Gam}\left(\alpha \mid a_{0}, b_{0}\right) \\
&p(\beta)=\operatorname{Gam}\left(\beta \mid c_{0}, d_{0}\right)
\end{aligned}
$$
と書くことができる．
ここで変分推論の枠組みで考え，問題中の設定から変分事後分布は
$$
q(\mathbf{w}, \alpha, \beta) = q(\mathbf{w})q(\alpha)q(\beta)
$$
と分解できるとする．
$q(\mathbf{w}),q(\alpha),q(\beta)$の更新式を求める．まず$q(\alpha)$から10.1節で導出した一般的な結果(10.9)を用いて
$$
\begin{aligned}
\ln q^*(\alpha)&=\mathbb{E}_{\mathbf{w}, \beta}[\ln p(\mathbf{t}, \mathbf{w}, \alpha, \beta \mid \mathbf{X})]\\
&=\mathbb{E}_{\mathbf{w}, \beta}[\ln p(\mathbf{t} \mid \mathbf{w}, \beta, \mathbf{X}) p(\mathbf{w} \mid \alpha) p(\alpha) p(\beta)]\\
&=\mathbb{E}_{\mathbf{w}, \beta}[\ln p(\mathbf{t} \mid \mathbf{w}, \beta, \mathbf{X})]+\mathbb{E}_{\mathbf{w}}[\ln \beta(\mathbf{w} \mid \alpha)]+\ln p(\alpha)+\mathbb{E}_{\beta}[\ln p(\beta)]\\
&=\mathbb{E}_{\mathbf{w}}\left[\ln N\left(\mathbf{w} \mid 0, \alpha^{-1} \mathbf{I}\right)\right]+\ln \operatorname{Gam}\left(\alpha \mid a_{0}, b_{0}\right)+\textrm{const}\\
&=\mathbb{E}_{\mathbf{w}}\left[\ln \frac{1}{(2 \pi)^{\frac{M}{2}}} \frac{1}{\left(\alpha^{-1}\right)^{\frac{1}{2}}} \operatorname{exp}\left(-\frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} \mathbf{w}\right)\right]+\ln \frac{1}{\Gamma\left(a_{0}\right)} b_{0}^{a_{0}} \alpha^{a_{0}-1} e^{-b_{0} \alpha}+ \textrm{const.}\\
&=\frac{M}{2} \ln \alpha-\frac{\alpha}{2} \mathbb{E}\left[\mathbf{w}^{\mathrm{T}} \mathbf{w}\right]+\left(a_{0}-1\right) \ln \alpha-b_{0} \alpha+ \textrm{const.}\\
&=\left(\frac{M}{2}+a_{0}-1\right) \ln \alpha-\left(\frac{1}{2} \mathbb{E}\left[\mathbf{w}^{\mathrm{T}} \mathbf{w}\right]+b_{0}\right) \alpha+ \textrm{const.}
\end{aligned}
$$
ここで$\beta$を導入した場合にも$\alpha$に依存しない項はconstに押し込んで計算することができるため，(10.92)-(10.95)式までの議論をそのまま用いることができる．
$$
q^*(\alpha)=\operatorname{Gam}\left(\alpha \mid a_{N}, b_{N}\right) , a_N=\frac{M}{2} + a_0, b_N=\frac{1}{2}\mathbb{E}\left[\mathbf{w}^{\mathrm{T}} \mathbf{w}\right]+b_0
$$
を得る．次に$q(\mathbf{w})$について(10.9)より
$$
\begin{aligned}
\ln q^*(\mathbf{w})&=\mathbb{E}_{\alpha, \beta}[\ln \beta(\mathbf{t}, \mathbf{w}, \alpha, \beta \mid x)]\\
&=\mathbb{E}_{\alpha, \beta}[\ln p(\mathbf{t} \mid \mathbf{w}, \beta, X) p(\mathbf{w} \mid \alpha) \gamma(\alpha) p(\beta)]\\
&=\mathbb{E}_{\beta}\left[\ln \prod_{n=1}^{N} N\left(\mathbf{t}_{n} \mid \mathbf{w}^{\mathrm{T}} \phi_{n}, \beta^{-1}\right)\right]+\mathbb{E}_{\alpha}\left[\ln N\left(w \mid 0, \alpha^{-1} I\right)\right]+ \textrm{const.}\\
&=\mathbb{E}_{\beta}\left[\sum_{n=1}^N\ln \frac{1}{(2 \pi)^{\frac{M}{2}}} \frac{1}{\left(\beta^{-1}\right)^{\frac{1}{2}}} \operatorname{exp}\left\{-\frac{\beta}{2}(\mathbf{t}_n-\mathbf{w}^{\mathrm{T}}\phi_{n})^2 \right\}\right]+\mathbb{E}_{\alpha}\left[\ln \frac{1}{(2 \pi)^{\frac{M}{2}}} \frac{1}{\left(\alpha^{-1}\right)^{\frac{1}{2}}} \operatorname{exp}\left(-\frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} \mathbf{w}\right)\right]+ \textrm{const.}\\
&=\mathbb{E}_{\beta}\left[\beta\right]\left(\mathbf{w}^{\mathrm{T}}\Phi^{\mathrm{T}}\mathbf{t}-\frac{1}{2}\mathbf{w}^{\mathrm{T}}\Phi^{\mathrm{T}}\Phi\mathbf{w}\right)-\frac{1}{2}\mathbb{E}_{\alpha}\left[\alpha\right]\mathbf{w}^{\mathrm{T}}\mathbf{w}+ \textrm{const.}\\
&=-\frac{1}{2}\mathbf{w}^{\mathrm{T}}\left(\mathbb{E}_{\beta}[\beta] \Phi^{\mathrm{T}} \Phi+\mathbb{E}_{\alpha}[\alpha]\mathbf{I}\right) \mathbf{w}+\mathbb{E}_{\beta}[\beta] \mathbf{w}^{\mathrm{T}} \Phi^{\mathrm{T}} \mathbf{t}+ \textrm{const.}\\
\end{aligned}
$$
これは$\mathbf{w}$に関して2次形式なのでガウス分布になり，平方完成すると
$$
q^*(\mathbf{w})=\mathcal{N}(\mathbf{w}\mid \mathbf{m}_N, \mathbf{S}_N)
$$
$$
\mathbf{m}_N=\mathbb{E}_{\beta}[\beta]\mathbf{S}_N\mathbf{\Phi}^{\mathrm{T}}\mathbf{t}
$$
$$
\mathbf{S}_N=\mathbb{E}_{\alpha}[\alpha]\mathbf{I}+\mathbb{E}_{\beta}[\beta]\mathbf{\Phi}^{\mathrm{T}}\mathbf{\Phi}
$$
を得る．最後に$q(\beta)$について(10.9)より
$$
\begin{aligned}
\ln q^{\star}(\beta) &=\mathbb{E}_{\mathbf{w}, \alpha}[\ln p(\mathbf{t}, \mathbf{w}, \alpha, \beta \mid \mathbf{X})]\\
&=\mathbb{E}_{\mathbf{w}, \alpha}[\ln p(\mathbf{t} \mid \mathbf{w}, \beta, \mathbf{X}) p(\mathbf{w} \mid \alpha) p(\alpha) p(\beta)]\\
&=\mathbb{E}_{\mathbf{w}}[\ln p(\mathbf{t} \mid \mathbf{w}, \beta)]+\ln p(\beta)+\text { const } \\
&= \frac{N}{2} \cdot \ln \beta-\frac{\beta}{2} \cdot \mathbb{E}\left[\sum_{n=1}^{N}\left(t_{n}-\mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}_{n}\right)^{2}\right]+\left(c_{0}-1\right) \ln \beta-d_{0} \beta +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\frac{\beta}{2} \cdot \mathbb{E}\left[\left\|\mathbf{\Phi}_{\mathbf{w}}-\mathbf{t}\right\|^{2}\right]-d_{0} \beta +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot\left\{\frac{1}{2} \cdot \mathbb{E}\left[\|\mathbf{\Phi} \mathbf{w}-\mathbf{t}\|^{2}\right]+d_{0}\right\} +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot\left\{\frac{1}{2} \cdot \mathbb{E}\left[\mathbf{w}^{\mathrm{T}} \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}_{\mathbf{w}}-2 \mathbf{t}^{\mathrm{T}} \mathbf{\Phi}_{\mathbf{w}}+\mathbf{t}^{\mathrm{T}} \mathbf{t}\right]+d_{0}\right\} +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot\left\{\frac{1}{2} \cdot \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbb{E}\left[\mathbf{w} \mathbf{w}^{\mathrm{T}}\right]\right]-\mathbf{t}^{\mathrm{T}} \mathbf{\Phi} \mathbb{E}[\mathbf{w}]+\frac{1}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+d_{0}\right\} +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot\left\{\frac{1}{2} \cdot \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}\left(\mathbf{m}_{N} \mathbf{m}_{N}^{\mathrm{T}}+\mathbf{S}_{N}\right)\right]-\mathbf{t}^{\mathrm{T}} \mathbf{\Phi} \mathbf{m}_{N}+\frac{1}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+d_{0}\right\} +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot\left\{\frac{1}{2} \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbf{S}_{N}\right]+\frac{1}{2} \mathbf{m}_{N}^{\mathrm{T}} \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbf{m}_{N}-\mathbf{t}^{\mathrm{T}} \mathbf{\Phi} \mathbf{m}_{N}+\frac{1}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+d_{0}\right\} +\text { const }\\
&=\left(\frac{N}{2}+c_{0}-1\right) \cdot \ln \beta-\beta \cdot \frac{1}{2}\left\{\operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbf{S}_{N}\right]+\left\|\mathbf{\Phi} \mathbf{m}_{N}-\mathbf{t}\right\|^{2}+2 d_{0}\right\}+\text { const }
\end{aligned}
$$
これより
$$
q^{\star}(\beta)=\operatorname{Gam}\left(\beta \mid c_{N}, d_{N}\right)
$$
$$
c_N=\frac{N}{2}+c_0
$$
$$
d_N=d_{0}+\frac{1}{2}\left\{\operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbf{S}_{N}\right]+\left\|\mathbf{\Phi} \mathbf{m}_{N}-\mathbf{t}\right\|^{2}\right\}
$$
以上から各因子の更新式が得られた．
次に変分下界を求める．変分下界は本文中の式(10.107)を$\beta$を考慮した形に修正すれば得られ，考えるべき項は$\mathbb{E}\left[\ln p(\beta)\right], -\mathbb{E}\left[\ln q^*(\beta)\right]$の二つであるので，それぞれの計算をして，ディガンマ関数$\varphi(a)=\frac{d}{da}\ln\Gamma(a)$として
$$
\begin{aligned}
\mathbb{E}[\ln p(\beta)] &=\left(c_{0}-1\right) \mathbb{E}[\ln \beta]-d_{0} \mathbb{E}[\beta]+c_{0} \ln d_{0}-\ln \Gamma\left(c_{0}\right) \\
&=\left(c_{0}-1\right) \cdot\left(\varphi\left(c_{N}\right)-\ln d_{N}\right)-d_{0} \frac{c_{N}}{d_{N}}+c_{0} \ln d_{0}-\ln \Gamma\left(c_{0}\right)
\end{aligned}
$$
ここで(B.26)(ガンマ分布の関数形についての定義式),(B.30)(ガンマ分布に従う確率変数の自然対数の期待値がディガンマ関数に紐づけられる式)をそれぞれ用いた．また
$$
-\mathbb{E}\left[\ln q^{\star}(\beta)\right]=\left(c_{N}-1\right) \cdot \varphi\left(c_{N}\right)-c_{N}+\ln d_{N}-\ln \Gamma\left(c_{N}\right)
$$
ここでガンマ分布に従う確率変数のエントロピーについての式(B.31)を用いた．(10.107)-(10.112)の式を修正することで$\beta$を考慮に入れた変分下界を得る．
最後に予測分布を考える．
これも本文中の議論を$\beta$を考慮したものに修正して得ることができて(10.105),(10.106)から
$$
\begin{aligned}
p(t \mid \mathbf{x}, \mathbf{t}) &=\int p(t \mid \mathbf{x}, \mathbf{w}) p(\mathbf{w} \mid \mathbf{t}) \mathrm{d} \mathbf{w} \\
& \simeq \int p(t \mid \mathbf{x}, \mathbf{w}) q(\mathbf{w}) \mathrm{d} \mathbf{w} \\
&=\int \mathcal{N}\left(t \mid \mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}(\mathbf{x}), \beta^{-1}\right) \mathcal{N}\left(\mathbf{w} \mid \mathbf{m}_{N}, \mathbf{S}_{N}\right) \mathrm{d} \mathbf{w} \\
&=\mathcal{N}\left(t \mid \mathbf{m}_{N}^{\mathrm{T}} \boldsymbol{\phi}(\mathbf{x}), \sigma^{2}(\mathbf{x})\right)
\end{aligned}
$$
ここで分散は
$$
\sigma^{2}(\mathbf{x})=\frac{1}{\mathbb{E}\left[\beta\right]}+\boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \mathbf{S}_{N} \boldsymbol{\phi}(\mathbf{x})
$$
である．
## 演習 10.27
<div class="panel-primary">
付録Bで与えられている公式を用いて， 線形基底関数回帰モデルの変分下界は
$$
\begin{aligned} \mathcal{L}(q)&= \mathbb{E}[\ln p(\mathbf{w}, \alpha, \mathbf{t})]-\mathbb{E}[\ln q(\mathbf{w}, \alpha)] \\
&= \mathbb{E}_{\mathbf{w}}[\ln p(\mathbf{t} \mid \mathbf{w})]+\mathbb{E}_{\mathbf{w}, \alpha}[\ln p(\mathbf{w} \mid \alpha)]+\mathbb{E}_{\alpha}[\ln p(\alpha)] \\
&-\mathbb{E}_{\alpha}[\ln q(\mathbf{w})]_{\mathbf{w}}-\mathbb{E}[\ln q(\alpha)] \end{aligned} \tag{10.107}
$$
の形で書け，その各項は
$$\begin{aligned} \mathbb{E}[\ln p(\mathbf{t} \mid \mathbf{w})]_{\mathbf{w}}=& \frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+\beta \mathbf{m}_{N}^{\mathrm{T}} \mathbf{\Phi}^{\mathrm{T}} \mathbf{t} \\
&-\frac{\beta}{2} \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}\left(\mathbf{m}_{N} \mathbf{m}_{N}^{\mathrm{T}}+\mathbf{S}_{N}\right)\right] \end{aligned} \tag{10.108}$$
$$\begin{aligned} \mathbb{E}[\ln p(\mathbf{w} \mid \alpha)]_{\mathbf{w}, \alpha}=&-\frac{M}{2} \ln (2 \pi)+\frac{M}{2}\left(\psi\left(a_{N}\right)-\ln b_{N}\right) \\
&-\frac{a_{N}}{2 b_{N}}\left[\mathbf{m}_{N}^{\mathrm{T}} \mathbf{m}_{N}+\operatorname{Tr}\left(\mathbf{S}_{N}\right)\right] \end{aligned}\tag{10.109}$$
$$\begin{aligned} \mathbb{E}[\ln p(\alpha)]_{\alpha}=&\ a_{0} \ln b_{0}+\left(a_{0}-1\right)\left[\psi\left(a_{N}\right)-\ln b_{N}\right] \\
&-b_{0} \frac{a_{N}}{b_{N}}-\ln \Gamma\left(a_{0}\right) \end{aligned} \tag{10.110}$$
$$-\mathbb{E}[\ln q(\mathbf{w})]_{\mathbf{w}}=\frac{1}{2} \ln \left|\mathbf{S}_{N}\right|+\frac{M}{2}[1+\ln (2 \pi)] \tag{10.111}$$
$$-\mathbb{E}[\ln q(\alpha)]_{\alpha}=\ln \Gamma\left(a_{N}\right)-\left(a_{N}-1\right) \psi\left(a_{N}\right)-\ln b_{N}+a_{N} \tag{10.112}$$
となることを示せ．
</div>
※演習10.16, 10.17のように各項の確率分布に適切なものを当てはめて計算していくだけ。
$$
\begin{aligned}
\mathbb{E}_{\mathbf{w}}[\ln p(\mathbf{t} \mid \mathbf{w})] &=\mathbb{E}_{\mathbf{w}}\left[\ln \prod_{n=1}^{N} \mathcal{N}\left(t_{n} \mid \mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}_{n}, \beta^{-1}\right)\right]\hspace{1em}(\because(B.87))\\
&=\mathbb{E}_{\mathbf{w}}\left[\sum_{n=1}^{N} \ln \mathcal{N}\left(t_{n} \mid \mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}_{n}, \beta^{-1}\right)\right] \\
&=\mathbb{E}_{\mathbf{w}}\left[\sum_{n=1}^{N} \ln \left\{\left(\frac{\beta}{2 \pi}\right)^{\frac{1}{2}} \exp \left\{-\frac{\beta}{2}\left(t_{n}-\mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}_{n}\right)^{2}\right\}\right.\right.\\
&=\frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbb{E}_{\mathbf{w}}\left[\sum_{n=1}^{N}\left(t_{n}-\mathbf{w}^{\mathrm{T}} \boldsymbol{\phi}_{n}\right)^{2}\right] \\
&=\frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbb{E}_{\mathbf{w}}\left[(\mathbf{t}-\mathbf{\Phi} \mathbf{w})^{\mathrm{T}}(\mathbf{t}-\mathbf{\Phi} \mathbf{w})\right] \\
&=\frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+\beta \mathbb{E}_{\mathbf{w}}\left[\mathbf{w}^{\mathrm{T}}\right] \mathbf{\Phi}^{\mathrm{T}} \mathbf{t}-\frac{\beta}{2} \mathbb{E}_{\mathbf{w}}\left[\mathbf{w}^{\mathrm{T}} \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbf{w}\right] \\
&=\frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+\beta \mathbf{m}_{N}^{\mathrm{T}} \mathbf{\Phi}^{\mathrm{T}} \mathbf{t}-\frac{\beta}{2} \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi} \mathbb{E}_{\mathbf{w}}\left[\mathbf{ww}^{\mathrm{T}}\right]\right] \\
&=\frac{N}{2} \ln \left(\frac{\beta}{2 \pi}\right)-\frac{\beta}{2} \mathbf{t}^{\mathrm{T}} \mathbf{t}+\beta_{m N}^{\mathrm{T}} \mathbf{\Phi} \mathbf{t}-\frac{\beta}{2} \operatorname{Tr}\left[\mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}\left(\mathbf{m}_{N} \mathbf{m}_{N}^{\mathrm{T}}+\mathbf{S}_{N}\right)\right]
\end{aligned}
$$
----
$$
\begin{aligned} \mathbb{E}_{\mathbf{w}, \alpha}[\ln p(\mathbf{w} \mid \alpha)] &=\mathbb{E}_{\mathbf{w}, \alpha}\left[\ln \mathcal{N}\left(\mathbf{w} \mid \mathbf{0}, \alpha^{-1} \mathbf{I}\right)\right] \\
&=\mathbb{E}_{\mathbf{w}, \alpha}\left[\ln \left\{\left(\frac{\alpha}{2 \pi}\right)^{M / 2} \exp \left\{-\frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} w\right\}\right]\right.\\
&=\mathbb{E}_{\mathbf{w}, \alpha}\left[\frac{M}{2} \ln \left(\frac{\alpha}{2 \pi}\right)\right]-\frac{1}{2} \mathbb{E}_{\mathbf{w}, \alpha}\left[\alpha \mathbf{w}^{\mathrm{T}} \mathbf{w}\right] \\
&=-\frac{M}{2} \ln (2 \pi)+\frac{M}{2} \mathbb{E}_{\alpha}[\ln \alpha]-\frac{\mathbb{E}_{\alpha}[\alpha]}{2} \mathbb{E}_{\mathbf{w}}\left[\mathbf{w}^{\mathrm{T}} \mathbf{w}\right] \\
&=-\frac{M}{2} \ln (2 \pi)+\frac{M}{2} \underbrace{\left(\psi\left(a_{N}\right)-\ln b_{N}\right)}_{(B.30)} - \underbrace{\frac{a_{N}}{2 b_{N}}}_{(B.27)} \mathbb{E}_{\mathbf{w}}\left[\operatorname{Tr}\left(\mathbf{w} \mathbf{w}^{\mathrm{T}}\right)\right] \\
&=-\frac{M}{2} \ln (2 \pi)+\frac{M}{2} \left(\psi\left(a_{N}\right)-\ln b_{N}\right) - \frac{a_{N}}{2 b_{N}} \left[\mathbf{m}_{N}^{\mathrm{T}} \mathbf{m}_{N}+\operatorname{Tr}\left(\mathbf{S}_{N}\right)\right]
\end{aligned}
$$
ここで$\mathbb{E}_{\mathbf{w}}\left[\operatorname{Tr}\left(\mathbf{w} \mathbf{w}^{\mathrm{T}}\right)\right]$の変形についてはトレースと期待値の交換性と
$$
\begin{aligned} & \operatorname{Tr}\left[\mathbb{E}_{\mathbf{w}}\left[\mathbf{w} \mathbf{w}^{\mathrm{T}}\right]\right] \\
= & \operatorname{Tr}\left[\operatorname{cov}[\mathbf{w}]+\mathbb{E}_{\mathbf{w}}[\mathbf{w}] \mathbb{E}_{\mathbf{w}}\left[\mathbf{w}^{\mathrm{T}}\right]\right] \quad(\because(1.42)) \\
= & \operatorname{Tr}\left[\mathbf{S}_{N}+\mathbf{m}_{N} \mathbf{m}_{N}^{\mathrm{T}}\right] \\
= &\ \mathbf{m}_{N}^{\mathrm{T}} \mathbf{m}_{N}+\operatorname{Tr}\left(\mathbf{S}_{N}\right)
\end{aligned}
$$
を用いた。
----
$$
\begin{aligned} \mathbb{E}_{\alpha}[\ln p(\alpha)] &=\mathbb{E}_{\alpha \sim q(\alpha)}\left[\ln \operatorname{Gam}\left(\alpha \mid a_{0}, b_{0}\right)\right] \\
&=\mathbb{E}_{\alpha \sim q(\alpha)}\left[\ln \left\{\frac{1}{\Gamma\left(a_{0}\right)} b_{0}^{a_{0}} \alpha^{a_{0}-1} e^{-b_{0} \alpha}\right\}\right] \\
&=\mathbb{E}_{\alpha \sim q(\alpha)}\left[-\ln \Gamma\left(a_{0}\right)+a_{0} \ln b_{0}+\left(a_{0}-1\right) \ln \alpha-b_{0} \alpha\right] \\
&=a_{0} \ln b_{0}+\left(a_{0}-1\right) \mathbb{E}_{\alpha}[\ln \alpha]-b_{0} \mathbb{E}_{\alpha}[\alpha]-\ln \Gamma\left(a_{0}\right) \\
&=a_{0} \ln b_{0}+\left(a_{0}-1\right)\left(\psi\left(a_{N}\right)-b_{N}\right)-b_{0} \frac{a_{N}}{b_{N}}-\ln \Gamma\left(a_{0}\right) \end{aligned}
$$
----
$$
\begin{aligned}
-\mathbb{E}_{\mathbf{w}}\left[\ln q(\mathbf{w})\right] &=-\mathbb{E}_{\mathbf{w} \sim q(\mathbf{w})}\left[\ln \mathcal{N}\left(\mathbf{w} \mid \mathbf{m}_{N}, \mathbf{S}_{N}\right)\right] \\
&=-\mathbb{E}_{\mathbf{w} \sim q(\mathbf{w})}\left[\ln \left\{\left(\frac{1}{2 \pi}\right)^{\frac{M}{2}} \frac{1}{\left|\mathbf{S}_{N}\right|^{\frac{1}{2}}} \exp \left\{-\frac{1}{2}\left(\mathbf{w}-\mathbf{m}_{N}\right)^{\mathrm{T}} \mathbf{S}_{N}^{-1}\left(\mathbf{w}-\mathbf{m}_{N}\right)\right\}\right]\right.\\
&=\frac{M}{2} \ln (2 \pi)+\frac{1}{2} \ln \left|\mathbf{S}_{N}\right|+\frac{1}{2} \operatorname{Tr}\left[\mathbb{E}_{\mathbf{w}}\left[\left(\mathbf{w}-\mathbf{m}_{N}\right)\left(\mathbf{w}-\mathbf{m}_{N}\right)^{\mathrm{T}}\right] \mathbf{S}_{N}^{-1}\right] \\
&=\frac{M}{2} \ln (2 \pi)+\frac{1}{2} \ln \left|\mathbf{S}_{N}\right|+\frac{1}{2} \operatorname{Tr}\left[\operatorname{cov}[\mathbf{w}] \mathbf{S}_{N}^{-1}\right] \\
&=\frac{M}{2} \ln (2 \pi)+\frac{1}{2} \ln \left|\mathbf{S}_{N}\right|+\frac{1}{2} M \\
&=\frac{1}{2} \ln \left|\mathbf{S}_{N}\right|+\frac{M}{2}[1+\ln (2 \pi)]
\end{aligned}
$$
----
$$
\begin{aligned}-\mathbb{E}_{\alpha}[\ln q(\alpha)] &=-\mathbb{E}_{\alpha \sim q(\alpha)}\left[\ln \operatorname{Gam}\left(\alpha \mid a_{N}, b_{N}\right)\right] \\
&=-\mathbb{E}_{\alpha \sim q(\alpha)}\left[-\ln \Gamma\left(a_{N}\right)+a_{N} \ln b_{N}+\left(a_{N}-1\right) \ln \alpha-b_{N} \alpha\right] \\
&=\ln \Gamma\left(a_{N}\right)-a_{N} \ln b_{N}-(a_N - 1)\mathbb{E}_{\alpha \sim q(\alpha)}[\ln \alpha]+b_{N} \mathbb{E}_{\alpha \sim q(\alpha)}[\alpha] \\
&=\ln \Gamma\left(a_{N}\right)-a_{N} \ln b_{N}-\left(a_{N}-1\right)\left(\psi\left(a_{N}\right)-\ln b_{N}\right)+b_{N} \frac{a_{N}}{b_{N}} \\
&=\ln \Gamma\left(a_{N}\right)-\left(a_{N}-1\right) \psi\left(a_{N}\right)-\ln b_{N}+a_{N} \end{aligned}
$$
## 演習 10.28
<div class="panel-primary">
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}$$
を導け．
</div>
ベイズ混合ガウス分布における、指数型分布族のそれぞれの関数形を導出して、一般的な結果(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{align}
\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{align}
と書ける。ただし、行列の上に矢印（$\rightarrow$）が書かれているのは、行列の各要素を並べた$D \times D$次元のベクトルを意味する。後の式変形の都合で、$g(\boldsymbol{\eta})$の要素を$\boldsymbol{\eta}$と$\mathbf{u}(\mathbf{x})$に押し込めて、
\begin{align}
\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{align}
と書き直す。
今考えているベイズ混合ガウスモデル:
\begin{align}
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{align}
の形に徐々に近づけていく。まずは$\pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})$を指数型分布族の標準形に対応づけるには、
\begin{align}
\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{align}
次に、$\left\{ \pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})\right\}^{z_{nk}}$を指数型分布族の標準形に対応づけるには、
\begin{align}
\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{align}
最後に、$\prod_{k=1}^K \left\{ \pi_k \mathcal{N} (\mathbf{x}_n|\boldsymbol\mu _k,\Lambda_k^{-1})\right\}^{z_{nk}}$を指数型分布族の標準形に対応づけるには、
\begin{align}
\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{align}
ここで、$[\ \ \ ]_{k=1,\cdots,K}$とは、$k=1$に対応するベクトル、$k=2$に対応するベクトル$\cdots$と順に並べてできる長いベクトルを表す。
今得られた対応関係を(10.116)式に代入して、
\begin{align}
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{align}
を得る。（最後の式変形は、(10.46)式の$\rho _{nk}$の定義より。）以上より、
\begin{align}
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{align}
と(10.48)式を得る。
----
次に、(10.119)式を用いて(10.57)式と(10.59)式を導く。
\begin{align}
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{align}
事前確率分布の項$\exp\left[\nu_0\boldsymbol\eta^T\boldsymbol\chi_0\right]$は、図10.5の事前確率分布
\begin{align}
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{align}
に一致するように$\nu_0\boldsymbol\chi_0$を選ぶことは可能（指数型分布族の事前共役分布の形を所与とすれば、この事実は証明不要な気もするが、念のため本問の解答の最後で$\nu_0\boldsymbol\chi_0$の具体的な形を構築する）。
\begin{align}
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{align}
この同時確率分布は$\boldsymbol\pi$に依存する項だけ積の形でくくり出せて、
\begin{align}
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{align}
であり、(10.57)式のディリクレ分布が導かれた。
残りの項のうち、$\boldsymbol\mu_k$に依存する項のみをくくり出すと、
\begin{align}
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{align}
１行目は、$\boldsymbol\mu_k$が(10.59)式のガウス分布に従うことを表す。
最後に、その他の項をくくり出すと、
\begin{align}
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{align}
となる。最後の式は、$|\Lambda|^{*}\exp [-\frac{1}{2}$Tr$(*\Lambda)]$という形をしており、(10.59)式のウィシャート分布を得る。（一見すると(10.62)式の形と異なるように見えるが、$\mathbf{S}_k$の定義に戻って計算すると、上式と一致することが示せる。
----
**おまけ：**
${\nu}_0\boldsymbol\chi_0$の具体的な形状は、以下の通り。
\begin{align}
{\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{align}
念のため、上の式を再現できることを確認しておく。
\begin{align}
{\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{align}
となる。なお、最後のTr$(\mathbf{W}_0^{-1}\boldsymbol\Lambda_k)$の項への変形にあたって、一般の正方行列$A,B$について
\begin{align}
{\rm Tr}[AB]={\rm Tr}\left[\left( \sum_j a_{ij}b_{jk}\right)_{ik} \right]
=\sum_{i, j}a_{ij}b_{ji}
\end{align}
が成り立つこと、および$B$が対称行列であれば最右辺が$\sum_{i,j} a_{ij}b_{ij}$に一致することを用いた。
## 演習 10.29
<div class="panel-primary">
二階微分を計算することで，関数$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)$が得られることを確かめよ．
</div>
※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
<div class="panel-primary">
二階微分を計算することで，対数ロジスティック関数$f(x)=-\ln \left(1+e^{-x}\right)$が上に凸であることを示せ．対数ロジスティック関数を点$x = \xi$のまわりで一次までテイラー展開することで，変分上界
$$\sigma(x) \leqslant \exp (\eta x-g(\eta)) \tag{10.137}$$
を直接導出せよ．
</div>
後の問題設定のため、対数ロジスティック関数である$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
<div class="panel-primary">
$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}$$
を直接導出せよ．
</div>
微分を計算する時の係数が煩わしいので、$x/2$を$x$と定義し直す。（証明すべき凸性に影響はない。最後の計算結果で元に戻す。）
+ $x$に関して上に凸であることの証明
\begin{align}
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{align}
+ $x^2$に関して下に凸であることの証明
$y=x^2$とおく。$f(y)=-\ln (e^{\sqrt{y}}+e^{-\sqrt{y}})$を真面目に$y$で2階微分するのは面倒なので、以下のように解く。
\begin{align}
\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{align}
ここで、
\begin{align}
\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{align}
であるから、
\begin{align}
\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{align}
ここで、分子が$x>0$のとき正、$x<0$のとき負であることは、微分して普通に示しても良いが、
\begin{align}
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{align}
に注目すれば明らか。よって、右辺全体は$x>0$でも$x<0$でも正なので、$\frac{d^2}{dy^2}f(x)>0$と言える。
+ グラフを書いてみよ
$f(x)$のグラフを$x$について描く
$f(x)$のグラフを$x^2$について描く
+ 下界(10.44)を直接導出する（$x$の定義は元に戻す）
\begin{align}
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{align}
となる（不等号は、$x^2$に関する凸性より）。従って、
\begin{align}
\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{align}
両辺のexpをとって、(10.144)式を得る。
## 演習 10.32
<div class="panel-primary">
ロジスティック回帰の変分ベイズ推定を時系列に従って学習し，データ点が一回に一つだけ到着し，次のデータ点が到着するまでに処理して廃棄しなければならない場合について考える．このとき，事後分布のガウス近似は下界
$$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$が最適化される．
</div>
※(方針)
(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}^{T}+\mathbf{S}_{N}^{-1}
$$
が得られる．
## 演習 10.33
<div class="panel-primary">
$$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$の更新式を示せ．
</div>
### Preparation
- $\sigma(x)$の微分
\begin{aligned}
\frac{\partial \sigma(x)}{\partial x} &= \frac{e^{-x}}{(1 + e^{-x})^2} \\
&= \sigma(x)(1 - \sigma(x)) \notag
\end{aligned}
- $\log(\sigma(x))$の微分
\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) \notag
\end{aligned}
- $\lambda(\xi)$
$$
\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
<div class="panel-primary">
この演習問題では，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}$$
の式を利用せよ．
</div>
### 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
<div class="panel-primary">
変分ベイズロジスティック回帰モデルの下界$\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)$を求めよ．
</div>
- $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}
- $\mathbf{w}$で積分する
\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
<div class="panel-primary">
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}$$
を導け．
</div>
※
## 演習 10.37
<div class="panel-primary">
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})$の場合に起こり，事前分布に対応する因子は一度だけ厳密な形で取り込まれ，以降は更新する必要はないことがわかる．
</div>
$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
<div class="panel-primary">
この演習問題と次の演習問題では，雑音データ問題での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}$$
を利用すればよい．
</div>
(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
<div class="panel-primary">
雑音データ問題での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{align} \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{align}$$
となることを証明し，$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)$の結果を導け．
</div>
### 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}
- $\rho_n$
\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)式を得る。
Sign in with Wallet
Connect another wallet