🧠

PRML 第10章(10.28から10.39まで)解答例

2022/05/15に公開
8

はじめに

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)に代入していく。

1変数ガウス分布と、指数型分布族の標準形との対応関係は、演習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}\sum_{k=1}^K \left\{ \boldsymbol\mu_k^{\mathrm{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^{\mathrm{T}}\beta_0\boldsymbol\Lambda_k \mathbf{m}_0 + \sum_{n=1}^N r_{nk} \boldsymbol\mu_k^{\mathrm{T}} \boldsymbol\Lambda_k \mathbf{x}_n \right) \right\} \right]\\ =\ &\exp \left[ -\frac{1}{2}\sum_{k=1}^K \left\{ \boldsymbol\mu_k^{\mathrm{T}} \left((\beta_0+N_k) \boldsymbol\Lambda_k \right)\boldsymbol\mu_k -2\boldsymbol\mu_k^{\mathrm{T}}\left(\boldsymbol\Lambda_k (\beta_0 \mathbf{m}_0 + \underbrace{N_{k} \overline{\mathbf{x}_k}}_{\sum_{n=1}^N r_{nk}\mathbf{x}_n} ) \right) \right\} \right]\\ = & \exp \left[ -\frac{1}{2}\sum_{k=1}^K \left\{ \left( \boldsymbol\mu_k - \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}\right)^{\mathrm{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.\right. \\ \\ &\hspace{2em} \left.\left.-\frac{1}{\beta_0+N_k}\left(\beta_0 \mathbf{m}_0+N_k \overline{\mathbf{x}_k}\right)^{\mathrm{T}} \boldsymbol\Lambda_k\left(\beta_0 \mathbf{m}_0+N_k \overline{\mathbf{x}_k}\right)\right\}\right] \\ = & \prod_{k=1}^K |\boldsymbol\Lambda_{k}|^{\frac{1}{2}}\exp \left[-\frac{1}{2} \left\{ \left( \boldsymbol\mu_k - \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}\right)^{\mathrm{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\} \right] \\ \cdot & \prod_{k=1}^{K} |\boldsymbol\Lambda_{k}|^{-\frac{1}{2}}\exp \left[\frac{1}{2} \frac{1}{\beta_0+N_k}\left(\beta_0 \mathbf{m}_0+N_k \overline{\mathbf{x}_k}\right)^{\mathrm{T}} \boldsymbol\Lambda_k\left(\beta_0 \mathbf{m}_0+N_k \overline{\mathbf{x}_k}\right)\right] \\ \propto & \prod_{k=1}^K |\boldsymbol\Lambda_{k}|^{\frac{1}{2}}\exp \left[-\frac{1}{2} \left\{ \left( \boldsymbol\mu_k - \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}\right)^{\mathrm{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\} \right] \\ = & \prod_{k=1}^{K} \mathcal{N}\left(\boldsymbol\mu_k \mid \frac{\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k}}{\beta_0+N_k}, \left((\beta_0+N_k)\boldsymbol\Lambda_k\right)^{-1}\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^{\mathrm{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^{\mathrm{T}} \boldsymbol\Lambda_k \mathbf{x}_n +\frac{1}{2} \ln |\boldsymbol{\Lambda}_k| \right) \right\} \right] \\ & \cdot \prod_{k=1}^{K}|\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})^{\mathrm{T}}\boldsymbol\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^{\mathrm{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^{\mathrm{T}} \Lambda_k \mathbf{x}_n -\frac{1}{\beta_0+N_k} (\beta_0\mathbf{m}_0+N_k\overline{\mathbf{x}_k})^{\mathrm{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^{\mathrm{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^{\mathrm{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})^{\mathrm{T}}\boldsymbol\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^{\mathrm{T}} +\mathbf{W}_0^{-1} + \sum_{n=1}^N r_{nk} \mathbf{x}_n \mathbf{x}_n^{\mathrm{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})^{\mathrm{T}} \right)\boldsymbol\Lambda_k \right\} \right] \\ \end{aligned}

となる。最後の式は、|\boldsymbol\Lambda|^{*}\exp [-\frac{1}{2}Tr(*\boldsymbol\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 \inftyf''(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/2xと定義し直す。(証明すべき凸性に影響はない。最後の計算結果で元に戻す。)

  • 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}
  • x^2に関して下に凸であることの証明

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と言える。

  • グラフを書いてみよ
    f(x)のグラフをxについて描く

    f(x)のグラフをx^2について描く

  • 下界(10.44)を直接導出する(xの定義は元に戻す)

\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

  • \sigma(x)の微分
\begin{aligned} \frac{\partial \sigma(x)}{\partial x} &= \frac{e^{-x}}{(1 + e^{-x})^2} \\ &= \sigma(x)(1 - \sigma(x)) \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) \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

この演習問題では,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}
  • \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

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}
  • \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)式を得る。

Discussion

ChoikoChoiko

演習問題10.28について、μkのくくり出しの部分とその他の項のくくり出しのところですが、平方完成により外出しされている部分にΠが必要かと思います。私の勘違いでしたらスルーして下さいませ。

ChoikoChoiko

まず、(10.57)のディリクレ分布を導いた後の部分ですが、
「残りの項のうち、μkに依存する項のみをくくり出すと、」の4行下の式部分(”おまけ”の上へ11行目)はexpの前に(”1〜Kをかける”という記号としての)Πをつけるか、もしくはexp()内に(”1〜Kを足す”という記号としての)Σをつける必要があると思います。当該部分の2行上と3行上にΣ記号がありますが、その外側を{}で囲むのではなく内側を{}で囲むべきかと(この点も先に触れるべきでした)。くくり出す元の式がそうなっております。
同じく
(10.59)のガウス分布を導いた後の部分ですが、
「最後にその他の項を、くくり出すと」の2行下の式部分(”おまけ”の上へ7行目)も同様と思います。

ChoikoChoiko

ご修正いただき、大変恐縮です。
2点だけ、ございます。
・(10.57)のディリクレ分布を導いた後の、「残りの項のうち、μkに依存する項のみをくくり出すと」のすぐ下に1行目から3行目までの部分ですが、-1/2の右にある”{”はΣのすぐ右に記載されるのが良いと考えます。
※記載の仕方として、Σ{a+b}≠Σa+bという前提での話です。
・「残りの項のうち、μkに依存する項のみをくくり出すと」の下に5行目と7行目ですが、ΣをΠで外出ししていますのでΣを除外すべきところ、残ったままとなっているのではないかと思います。

ChoikoChoiko

お忙しい中、ご修正下さり、誠にありがとうございました。

ChoikoChoiko

演習問題10.34のPreparationの”後で使う式変形”の3行目ですが、Sに付されている”−1”(逆行列)は不要と思いますがいかがでしょうか。