# PRML第12章演習問題解答（12.16〜12.29） <head> <style> div.panel-primary { border: 1px solid #000; margin: 10px 5px; padding: 16px 10px 0px; } </style> </head> [演習12.15まで](https://hackmd.io/GrOUx8qbSEmmEw_1WgPX8A?view) ## 演習 12.16 <div class="panel-primary"> ![](https://i.imgur.com/kXPibhJ.png) 図12.11においていくつかのデータの値がランダムに欠損しているデータ集合に確率的主成分分析を適用する例を説明した．この状況において，確率的主成分分析モデルの尤度関数を最大化するEMアルゴリズムを導け．$\{\mathbf{z}_n\}$に加えて，ベクトル$\{\mathbf{x}_n\}$の要素である欠損データ値も潜在変数として扱う必要があることに注意せよ．すべてのデータ値が観測される特別な場合に，このアルゴリズムが12.2.2節で導かれた確率的主成分分析のEMアルゴリズムに帰着することを示せ． </div> まず、尤度を各データの各要素の和の形に分解する。 \begin{align} p(\mathbf{X}|\mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) &= \int p(\mathbf{Z}) p(\mathbf{X}| \mathbf{Z}, \mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) d \mathbf{Z} \\ &= \prod_n^N \int p(\mathbf{z}_{n}) p(\mathbf{x}_n| \mathbf{z}_n, \mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) d \mathbf{z}_n \\ &= \prod_n^N \int p(\mathbf{z}_{n})\prod_i^D \mathcal{N}(x_{ni}| \mathbf{w}_i \mathbf{z}_n+ \mathbf{\mu}_i, \mathbf{\sigma^2}) d \mathbf{z}_n \\ \end{align} 今、$x_{ni}$は、$\mathbf{x}_n$の$i$番目の要素を示す。$\mathbf{\mu}_i$は$\mathbf{\mu}$の$i$番目の要素を、$\mathbf{w}_i$は、$\mathbf{W}$の$i$行目を指す。ランダム欠損であると仮定すると、上の式から欠損を積分消去できることによって尤度を求められる。よって、観察された変数、欠損した変数をそれぞれ$\mathbf{x}_n^o$、$\mathbf{x}_n^m$として書き直すと、 \begin{align} p(\mathbf{X}|\mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) &= \prod_n^N \int p(\mathbf{z}_{n})\prod_{x_{ni} \in \mathbf{x}_n^o} \mathcal{N}(x_{ni}| \mathbf{w}_i \mathbf{z}_n+ \mathbf{\mu}_i, \mathbf{\sigma^2}) \int \prod_{x_{ni} \in \mathbf{x}_n^m} \mathcal{N}(x_{ni}| \mathbf{w}_i \mathbf{z}_n+ \mathbf{\mu}_i, \mathbf{\sigma^2}) d \mathbf{x}_n^m d \mathbf{z}_n\\ &= \prod_n^N \int p(\mathbf{z}_{n})\prod_{x_{ni} \in \mathbf{x}_n^o} \mathcal{N}(x_{ni}| \mathbf{w}_i \mathbf{z}_n+ \mathbf{\mu}_i, \mathbf{\sigma^2}) d \mathbf{z}_n \\ &= \prod_n^N \int p(\mathbf{z}_{n}) \mathcal{N}(\mathbf{x_{n}^o}| \mathbf{z}_n, \mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) d \mathbf{z}_n \\ &= \prod_n p(\mathbf{x}_n^o|\mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) \end{align} すなわち、観察された$\mathbf{x}_n^o$から求められる。 ここで、観察/欠損を区別するために$\iota_{ni}$を導入する。これは$x_{ni}$が観察されたときに1を取り、欠損の時に0を取る変数である。よって、(12.32)式を書き直すと、 \begin{align} p(\mathbf{x}|\mathbf{z}) = \prod_i^D \mathcal{N}(x_{ni}|\mathbf{w}_i\mathbf{z}+ \mu_i, \sigma^2)^{\iota_{ni}} \space \space (Ex12.16.1) \end{align} になり、完全データ対数尤度は、これを用いて、 \begin{align} \ln p(\mathbf{X}, \mathbf{Z}|\mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) &= \sum_n^N \{ \ln p(\mathbf{z}_n) + \sum_i^D \ln \iota_{ni} \mathcal{N}(x_{ni}| \mathbf{w}_i\mathbf{z}+ \mu_i, \sigma^2) \} \end{align} と書ける。あとは、欠損がない場合(p. 295)の通り、$p(\mathbf{z}) = \mathcal{N}(\mathbf{z}|\mathbf{0}, \mathbf{I})(12.31)$ と(Ex.12.16.1)を用いて、完全データ対数尤度の期待値をとると、 \begin{align} E[\ln p(\mathbf{X}, \mathbf{Z}|\mathbf{\mu}, \mathbf{W}, \mathbf{\sigma^2}) ] = -\sum_n\{\frac{M}{2} \ln (2\pi) + \frac{1}{2}Tr(E[\mathbf{z}_n\mathbf{z}_n^T]) + \sum_i^D \iota_{ni}\{\ln (2\pi \sigma) + \frac{1}{2\sigma^2}(x_{ni}-\mu_{ni})^2 - \frac{1}{\sigma^2} E[\mathbf{z}_n]^T\mathbf{w}_i^T(x_{ni}-\mu_{ni})+ \frac{1}{2\sigma^2}Tr(E[\mathbf{z}\mathbf{z}^T]\mathbf{w}_i^T\mathbf{w}_i) \}\} \end{align} あとはこの式を更新するEMアルゴリズムを求めていけば良い。教科書と同様に$\mu$で微分して0を取る$\mu$を求めると、 \begin{align} \mu_i^{new} = \frac{1}{\sum_m^M \iota_{mi}} \sum_n \iota_{ni} x_{ni} \end{align} そしてEステップにおいて(12.54)、(12.55)相当を求めると、 \begin{align} E[\mathbf{z}_n] = \mathbf{M}^{-1}\mathbf{W}^T\mathbf{y}_n \end{align} なお、$y_n$は、$\mathbf{x}_n - \mathbf{\bar{x}}$の観察された部分である。$\mathbf{W}_n$も$\mathbf{W}$の観察された$\mathbf{x}_n$に対応する。(12.41)に対応するのは、 \begin{align} \mathbf{M}_n = \mathbf{W}_n^T\mathbf{W}_n + \sigma^2 \mathbf{I} \end{align} である。そして、(12.55)式も \begin{align} E[\mathbf{z}_n \mathbf{z}_n^T ] = \sigma^2 \mathbf{M}_n^{-1} + E[\mathbf{z}_n]E[ \mathbf{z}_n]^T \end{align} そして、Mステップも同様に、 \begin{align} \mathbf{W}_{new} = [\sum_n^N \mathbf{y}_n E[\mathbf{z}_n]^T][\sum_n^N E[\mathbf{z}_n \mathbf{z}_n^T] ]^{-1} \end{align} \begin{align} \sigma_{new}^2 = \frac{1}{\sum_n^N \sum_i^D \iota_{ni}} \sum_n^N \sum_i^D \iota_{ni} \{ (x_{ni}-\mu_{ni}^{new})^2 - 2 E[\mathbf{z}_n]^T(\mathbf{w}_i^{new})^T(x_{ni}-\mu_{ni}^{new})+ Tr(E[\mathbf{z}\mathbf{z}^T](\mathbf{w}_i^{new})^T\mathbf{w}_i^{new}) \}\} \end{align} なお、$\mathbf{w}_i^{bew}$は$\mathbf{W}_{new}$の$i$番目の行である。 そして、もし全て観察された場合には、$\iota_{ni}$は全て1であり、$\mathbf{y}_n = \mathbf{x}_n$であり、$\mathbf{W}_n = \mathbf{W}$であり、$\mu^{new} = \bar{\mathbf{x}}$である。この時、教科書12.54-12.57になる。 ## 演習 12.17 <div class="panel-primary"> $D\times M$行列を$\mathbf{W}$とおき，その列ベクトルが，$D$次元のデータ空間に埋め込まれた$M$次元の線型部分空間を定義するものとする．$\boldsymbol{\mu}$を$D$次元ベクトルとおく．$n=1,\ldots,N$に対しデータ集合$\mathbf{x}_n$が与えられたものとする．このとき，各$\mathbf{x}_n$を，$M$次元ベクトルの集合$\mathbf{z}_n$からの線形写像を用いて$\mathbf{W}\mathbf{z}_n + \boldsymbol{\mu}$で近似することを考える．再構成のコスト関数は，二乗和からなる $$ J=\sum_{n=1}^{N}\left\|\mathbf{x}_{n}-\boldsymbol{\mu}-\mathbf{W}{\mathbf{z}_{n}}\right\|^{2} \tag{12.95} $$ を使う．まず，$\boldsymbol{\mu}$についての$J$の最小化が，$\mathbf{x}_n$と$\mathbf{z}_n$をそれぞれ平均$0$の変数$\mathbf{x}_n - \overline{\mathbf{x}}$と$\mathbf{z}_n - \overline{\mathbf{z}}$で置き換えたのと同様な表現に導くことを示せ．ただし$\overline{\mathbf{x}}$と$\overline{\mathbf{z}}$はサンプル平均を表す．次に，$\mathbf{W}$を固定したまま$J$を$\mathbf{z}_n$について最小化すると，主成分分析のEステップ(12.58)が出ることを示せ．また{$\mathbf{z}_n$}を固定したまま$J$を$\mathbf{W}$について最小化すると，主成分分析のM ステップ(12.59)が出ることを示せ． </div> (i) $\boldsymbol{\mu}$による$J$の最小化 $$ \begin{aligned} \frac{\partial J}{\partial \boldsymbol{\mu}} &=\sum_{n=1}^{N} 2\left(\boldsymbol{\mu}-\left(\mathbf{x}_{n}-\mathbf{W}\mathbf{z}_{n}\right)\right\} \\ &=2 N\{\boldsymbol{\mu}-(\overline{\mathbf{x}}-\mathbf{W} \overline{\mathbf{z}})\} \end{aligned} $$ より$J$を最小化する$\boldsymbol{\mu}$は $$ \boldsymbol{\mu}=\mathbf{x}-\mathbf{W} \overline{\mathbf{z}} $$ であり、このとき、 $$ J=\sum_{n=1}^{N}\left\|\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)-\mathbf{W}\left(\mathbf{z}_{n}-\overline{\mathbf{z}}\right)\right\|^{2} $$ となる。 ここで $$ \mathbf{u}_{n}=\mathbf{x}_{n}-\overline{\mathbf{x}}\\ \mathbf{v}_{n}=\mathbf{z}_{n}-\overline{\mathbf{z}} $$ とすると $$ J=\sum_{n=1}^{N}\left\|\mathbf{u}_{n}-\mathbf{W}\mathbf{v}_{n}\right\|^{2} $$ と表せる。 (ii) $\mathbf{z}_{n}$による$J$の最小化 $$ \begin{aligned} \frac{\partial J}{\partial \mathbf{z}_{n}} &=\frac{\partial J}{\partial \mathbf{v}_{n}} \\&=\frac{\partial}{\partial \mathbf{v}_{n}}\left\|\mathbf{u}_{n}-\mathbf{W}\mathbf{v}_{n}\right\|^{2} \\ &=-2 \mathbf{W}^{T}\left(\mathbf{u}_{n}-\mathbf{W}\mathbf{v}_{n}\right) \end{aligned} $$ より、 $$ \mathbf{u}_n=\mathbf{W}\mathbf{v}_{n} $$ すなわち $$ \mathbf{v}_n=\left(\mathbf{W}^{T}\mathbf{W}\right)\mathbf{W}^T\mathbf{u}_{n} $$ が$J$を最小化する$\mathbf{v}_{n}$であり、これは $$ \tilde{\mathbf{X}}^{T}=\left(\mathbf{u}_{1} \mathbf{u}_{2} \cdots \mathbf{u}_{N}\right)\\ \mathbf{\Omega}=\left(\mathbf{v}_{1} \mathbf{v}_{2} \cdots \mathbf{v}_{N}\right) $$ とすると(12.58)式である $$ \mathbf{\Omega}=\left(\mathbf{W}^{T}\mathbf{W}\right)\mathbf{W}^T\tilde{\mathbf{X}}^{T} $$ と等しい。 (iii) $\mathbf{W}$による$J$の最小化 $$ \frac{\partial J}{\partial \mathbf{W}}=-2 \sum_{n=1}^{N}\left(\mathbf{u}_{n}-\mathbf{W} \mathbf{v}_{n}\right) \mathbf{v}_{n}^{T} $$ より、$J$を最小化する$\mathbf{W}$は $$ W \sum_{n=1}^{N} \mathbf{v}_{n}\mathbf{v}_{n}^{T}=\sum_{n=1}^{N} \mathbf{u}_{n} \mathbf{v}_{n}^{T}\\ \mathbf{W} \mathbf{\Omega} \mathbf{\Omega}^{T}=\tilde{\mathbf{X}}^{T} \mathbf{\Omega}^{T}\\ \mathbf{W}=\tilde{\mathbf{X}}^{T} \mathbf{\Omega}^{T}\left(\mathbf{\Omega} \mathbf{\Omega}^{T}\right)^{-1} $$ となりこれは(12.59)式に一致する ## 演習 12.18 <div class="panel-primary"> 12.2.4節で説明された因子分析モデルについて，独立なパラメータの数の表現を与える式を導け． </div> ※ 因子分析における共分散行列 (式12.65)　より、独立なパラメータの数は - $D \times M$行列$\mathbf{W}$ 及び - $M \times M$次元の対角行列$\mathbf{\Psi}$ に依存する。また、確率的主成分分析における議論 (pp.293-294) と同様に、$M(M - 1)/2$個のパラメータは回転に関して冗長。つまり、因子分析モデルにおける独立なパラメータの数は \begin{align*} & DM + D - M(M - 1)/2 \\ =& D(M + 1) - M(M - 1)/2 \end{align*} となる。 ## 演習 12.19 <div class="panel-primary"> 12.2.4節で説明した因子分析モデルが潜在変数空間の座標の回転に対して不変であることを示せ． </div> ※P.302の観測変数の周辺分布$p(\mathbf{x}) = \mathcal{N}(\mathbf{x}\mid \boldsymbol{\mu}, \mathbf{WW}^{\mathrm T}+\mathbf{\Psi})$が、潜在変数空間$\mathbf{z}$での回転$\mathbf{z} \to \mathbf{\tilde{z}}$としても不変、すなわち同じ式となることを示す。 P.289のような$M \times M$の直交行列（回転行列）$\mathbf{R}$を定義し（$\mathbf{RR}^{\mathrm T} = \mathbf{I}$）、潜在変数空間上で$\mathbf{z}$を$\mathbf{R}$で回転させた$\mathbf{\tilde{z}}$を考える。すなわち$\mathbf{\tilde{z}} = \mathbf{Rz}$となる。潜在変数上の事前分布は $$ p(\mathbf{\tilde{z}}) = \mathcal{N}(\mathbf{\tilde{z}} \mid \mathbf{0}, \mathbf{I}) \tag{12.31} $$ であり、これについて展開していくと $$ \begin{aligned} p(\widetilde{\mathbf{z}}) &=\mathcal{N}\left(\mathbf{R}\mathbf{z} \mid \mathbf{0}, \mathbf{I}\right) \\ &=\left(\frac{1}{2 \pi}\right)^{M / 2} \exp \left\{-\frac{1}{2}\left(\mathbf{R}\mathbf{z}\right)^{\mathrm T}\left(\mathbf{R}\mathbf{z}\right)\right\} \\ &=\left(\frac{1}{2 \pi}\right)^{M / 2} \exp \left\{-\frac{1}{2} \mathbf{z}^{\mathrm T} \mathbf{z}\right\} \\ &=\mathcal{N}(\mathbf{z} \mid \mathbf{0}, \mathbf{I}) \\ &=p(\mathbf{z}) \end{aligned} $$ であり、事前分布は回転不変である。また、潜在変数$\mathbf{z}$を与えられたときの観測変数$\mathbf{x}$の条件付き分布$p(\mathbf{x}\mid \mathbf{z})$について、P.289で用いた$\mathbf{\widetilde{W}} = \mathbf{WR}$を用いると $$ \begin{aligned} p(\mathbf{x} \mid \widetilde{\mathbf{z}}) &=N(\mathbf{x} \mid \mathbf{W} \tilde{\mathbf{z}}+\boldsymbol{\mu}, \mathbf{\Psi}) \\ &=\mathcal{N}(\mathbf{x} \mid \mathbf{WRz}+\boldsymbol{\mu}, \mathbf{\Psi}) \\ &=\mathcal{N}(\mathbf{x} \mid \widetilde{\mathbf{W}}\mathbf{z}+\boldsymbol{\mu}, \mathbf{\Psi}) \end{aligned} $$ $(2.115)$を用いて周辺分布$p(\mathbf{x})$を計算すると $$ \begin{aligned} p(\mathbf{x}) &=\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Psi}+\widetilde{\mathbf{W}} \mathbf{I} \widetilde{\mathbf{W}}^{\mathrm T}\right) \\ &=\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Psi}+\mathbf{WRR}^{\mathrm T} \mathbf{W}^{\mathrm T}\right) \\ &=\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Psi}+\mathbf{WW}^{\mathrm T}\right) \end{aligned} $$ これは回転しない場合でも同式である。したがって潜在変数空間における回転不変性が示された。 ## 演習 12.20 <div class="panel-primary"> 2次導関数を考えることによって12.2.4節で説明した因子分析モデルの対数尤度関数の，パラメータ$\boldsymbol{\mu}$に対する唯一の停留点が， $$ \overline{\mathbf{x}}=\frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_{n} \tag{12.1} $$ で定義されたサンプル平均で与えられることを示せ．さらに，この停留点が最大値を与えることを示せ． </div> ※対数尤度関数$\ln L$の$\boldsymbol{\mu}$に対する微分を計算すれば良い。 $$ \begin{aligned} \ln L &=\ln p\left(\mathbf{X} \mid \boldsymbol{\mu}, \mathbf{W}, \sigma^{2}\right) \\ &=\sum_{n=1}^{N} \ln p\left(\mathbf{x}_{n} \mid \boldsymbol{\mu}, \mathbf{W}, \sigma^{2}\right) \end{aligned} $$ 今、因子分析モデルでは$p(\mathbf{x}) = \mathcal{N}(\mathbf{x}\mid \boldsymbol{\mu}, \mathbf{C})$と置いているので $$ \begin{aligned} \ln L &=\sum_{n=1}^{N} \ln \mathcal{N}\left(\mathbf{x}_{n} \mid \boldsymbol{\mu}, \mathbf{C} \right) \\ &=\sum_{n=1}^{N} \ln \left(\left(\frac{1}{2 \pi}\right)^{D / 2}\left\|\mathbf{C}^{-1}\right\|^{1 / 2} \exp \left(-\frac{1}{2}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)^{\mathrm T} \mathbf{C}^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)\right)\right) \\ &=\sum_{n=1}^{N}\left\{-\frac{D}{2} \ln (2 \pi)-\frac{1}{2} \ln \|\mathbf{C}\|-\frac{1}{2}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)^{\mathrm T} \mathbf{C}^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)\right) \end{aligned} $$ $\boldsymbol{\mu}$について微分すると $$ \frac{\partial \ln L}{\partial \boldsymbol{\mu}} = \sum_{n=1}^{N}\mathbf{C}^{-1}(\mathbf{x}_{n}-\boldsymbol{\mu}) = 0 \\ N\boldsymbol{\mu} = \sum_{n=1}^{N}\mathbf{x}_{n} $$ よって停留点はサンプル平均である$\displaystyle \overline{\mathbf{x}}=\frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_{n}$となる。 また対数尤度関数のヘッセ行列$\mathbf{H}$を計算すると $$ \mathbf{H} = \nabla \otimes \nabla \ln L = -N\mathbf{C}^{-1} $$ $\mathbf{C} = \mathbf{WW}^{\mathrm T} + \mathbf{\Psi}$は正定置行列なので、ヘッセ行列は負定置になる。これより停留点は極大値となり、唯一の停留点なので極大かつ最大となる。 ## 演習 12.21 <div class="panel-primary"> 因子分析についてのEMアルゴリズムのEステップに対する公式 $$ \mathbb{E}\left[\mathbf{z}_{n}\right] =\mathbf{G} \mathbf{W}^{\mathrm{T}} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \tag{12.66} $$ $$\mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm{T}}\right] =\mathbf{G}+\mathbb{E}\left[\mathbf{z}_{n}\right] \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm{T}} \tag{12.67} $$ を導け．演習問題12.20の結果から，パラメータ$\boldsymbol{\boldsymbol{\mu}}$はサンプル平均$\overline{\mathbf{x}}$で置き換えられることに注意せよ． </div> 事前分布は、$p(\mathbf{z}) = \mathcal{N}(\mathbf{z}\mid \mathbf{0}, \mathbf{I})$ 因子分析モデルでは$p(\mathbf{x} \mid \mathbf{z})=\mathcal{N}(\mathbf{x} \mid \mathbf{W z}+\boldsymbol{\mu}, \Psi)$ またここでは$\boldsymbol{\mu} = \mathbf{\bar{x}}$と置き換えられる。$(2.116)$式を利用して事後分布を求めると $$ p(\mathbf{z} \mid \mathbf{x})=\mathcal{N}\left(\mathbf{z} \mid \mathbf{G}\left\{\mathbf{w}^{\mathrm T} \mathbf{\Psi}^{-1}(\mathbf{x}-\bar{\mathbf{x}})\right\}, \mathbf{G}\right) $$ ただし、ここで $$ \mathbf{G}=\left(\mathbf{I}+\mathbf{W}^{\mathrm{T}} \mathbf{\Psi}^{-1} \mathbf{W}\right)^{-1} \tag{12.68} $$ である。 この事後分布の期待値を計算すると $$ \begin{aligned} \mathbb{E}_{\mathbf{z}_{n} \sim p\left(\mathbf{z}_{n} \mid \mathbf{x}_{n}\right)}\left[\mathbf{z}_{n}\right] &=\mathbb{E}\left[\mathcal{N}\left(\mathbf{z}_{n} \mid \mathbf{G}\left(\mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1}(\mathbf{x}_{n}-\bar{\mathbf{x}})\right), \mathbf{G}\right)\right] \\ &=\mathbf{G}\mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1}(\mathbf{x}_{n}-\bar{\mathbf{x}}) \end{aligned} $$ また $$ \begin{aligned} \mathbb{E}_{\mathbf{z}_{n} \sim p\left(\mathbf{z}_{n} \mid \mathbf{x}_{n}\right)}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right] &=\operatorname{cov}\left[\mathbf{z}_{n}\right]+\mathbb{E}\left[\mathbf{z}_{n}\right] \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} \\ &=\mathbf{G}+\mathbb{E}\left[\mathbf{z}_{n}\right] \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} \end{aligned} $$ 以上で導かれた。 ## 演習 12.22 <div class="panel-primary"> 因子分析モデルの完全データ対数尤度関数の期待値の式を書き下せ．また，対応するMステップの式 $$ \mathbf{W}_{\text {new }}=\left[\sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm{T}}\right]\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm{T}}\right]\right]^{-1} \tag{12.69} $$ と $$ \mathbf{\Psi}_{\text {new }}=\operatorname{diag}\left\{\mathbf{S}-\mathbf{W}_{\text {new }} \frac{1}{N} \sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm{T}}\right\} \tag{12.70} $$ を導け． </div> ※P.294のような形でまず完全データの対数尤度関数$\ln L$を書くと（途中で演習12.20の結果である$\boldsymbol{\mu} = \mathbf{\overline{x}}$を用いている） $$ \begin{aligned} \ln L=& \sum_{n=1}^{N}\left\{\ln p\left(\mathbf{x}_{n} \mid \mathbf{z}_{n}\right) + \ln p\left(\mathbf{z}_{n}\right)\right\} \\ =& \sum_{n=1}^{N}\left\{\ln \mathcal{N}\left(\mathbf{x}_{n} \mid \mathbf{Wz}_{n} + \boldsymbol{\mu}, \mathbf{\Psi}\right) + \ln \mathcal{N}\left( \mathbf{z}_{n}\mid \mathbf{0}, \mathbf{I}\right) \right\} \\ =& \sum_{n=1}^{N}\left[-\frac{D}{2} \ln (2 \pi)-\frac{1}{2} \ln |\mathbf{\Psi}|-\frac{1}{2}\left\{\left(\mathbf{x}_{n}-\mathbf{Wz}_{n}-\boldsymbol{\mu}\right)^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\mathbf{Wz}_{n}-\boldsymbol{\mu}\right)\right\}\right. \\ &\left.-\frac{M}{2} \ln (2 \pi)-\frac{1}{2} \mathbf{z}_{n}^{\mathrm T} \mathbf{z}_{n}\right] \\ =&\ \frac{1}{2} \sum_{n=1}^{N}\left\{-M \ln (2 \pi)-\mathbf{z}_{n}^{\mathrm{T}} \mathbf{z}_{n}-D \ln (2 \pi)-\ln |\mathbf{\Psi}|\right.\\ &\left.-\left(\mathbf{x}_{n}-\overline{\mathbf{x}}-\mathbf{W} \mathbf{z}_{n}\right)^{\mathrm{T}} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}-\mathbf{W} \mathbf{z}_{n}\right)\right\} \end{aligned} $$ この対数尤度関数について期待値をとると $$ \begin{aligned} \mathbb{E}[\ln L]&=\frac{1}{2} \sum_{n=1}^{N}\left\{-\ln |\mathbf{\Psi}|-\mathbb{E}\left[\left(\mathbf{Wz}_{n}\right)^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{Wz}_{n}\right)\right] \right.\\ &\left.+\ 2 \mathbb{E}\left[\left(\mathbf{Wz}_{n}\right)^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\right]-\mathbb{E}\left[\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\right]\right\}+\text { const } \\ &=\frac{1}{2} \sum_{n=1}^{N}\left\{-\ln |\mathbf{\Psi}|-\operatorname{Tr}\left[\mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right] \mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1} \mathbf{W}\right] \right.\\ &\left.+\ 2 \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} \mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\right\} -\frac{1}{2} N \operatorname{Tr}\left[\mathbf{S} \mathbf{\Psi}^{-1}\right]+\text { const } \end{aligned} $$ 途中でトレースと期待値演算子の交換、またデータ共分散行列$\mathbf{S}$ $$ \mathbf{S}=\frac{1}{N} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm{T}} \tag{12.3} $$ を利用した。 この$\mathbb{E}[\ln L]$について$\mathbf{W}$の更新を行うために微分を取る。用意として $$ \begin{aligned} \frac{\partial}{\partial \mathbf{W}} \operatorname{Tr}\left[\mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right] \mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1} \mathbf{W}\right] &=\ \frac{\partial}{\partial \mathbf{W}} \operatorname{Tr}\left[\mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]\right] \\ &=\ \mathbf{\Psi}^{-1} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]+\left(\mathbf{\Psi}^{-1}\right)^{\mathrm T} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n}\mathbf{z}_{n}^{\mathrm T}\right]^{\mathrm T} (\because \textrm{Matrix Cookbook (117)})\\ &=\ 2 \mathbf{\Psi}^{-1} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n} ^{\mathrm T}\right] (\because \mathbf{\Psi}と\mathbb{E}[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}]は対称行列) \end{aligned} $$ $$ \frac{\partial}{\partial \mathbf{W}}\mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} \mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) = 2\mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} (\because \textrm{Matrix Cookbook (71)}) $$ これらを用いると $$ \frac{\partial}{\partial \mathbf{W}}\mathbb{E}[\ln L] = \frac{1}{2}\sum_{n=1}^{N} \left\{ -2 \mathbf{\Psi}^{-1} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n} ^{\mathrm T}\right] + 2\mathbf{\Psi}^{-1}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} \right\} $$ となり、停留点を求めるため$\displaystyle \frac{\partial}{\partial \mathbf{W}}\mathbb{E}[\ln L] = 0$として$\mathbf{W} \to \mathbf{W}_{\textrm{new}}$とすると $$ \mathbf{\Psi}^{-1} \mathbf{W}_{\textrm{new}}\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]=\mathbf{\Psi}^{-1} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T} $$ $$ \mathbf{W}_{\text {new }}=\left[\sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm{T}}\right]\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm{T}}\right]\right]^{-1} \tag{12.69} $$ を得られる。 次に$\mathbf{\Psi}$の更新について、同様に$\mathbf{\Psi}$の微分を計算すると $$ \begin{aligned} \frac{\partial}{\partial \mathbf{\Psi}} \mathbb{E}\left[\ln L\right]=\frac{1}{2} \sum_{n=1}^{N}\{&-\mathbf{\Psi}^{-\mathrm T}+\left(\mathbf{\Psi}^{-1} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right] \mathbf{W}^{\mathrm T} \mathbf{\Psi}^{-1}\right)^{\mathrm T} \\ &\left.-2 \mathbf{\Psi}^{-\mathrm T} \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T} \mathbf{\Psi}^{-\mathrm T}\right\} \\ &-\frac{N}{2}\left(-\mathbf{\Psi}^{-\mathrm T} \mathbf{S}^{\mathrm T} \mathbf{\Psi}^{-\mathrm T}\right) \end{aligned} $$ 上の変形にはMatrix Cookbookの(57), (61), (63)の公式を利用した。左と右からそれぞれ$\mathbf{\Psi}^{\mathrm T}$をかけ、$0$に等しいとすると $$ \frac{1}{2} \sum_{n=1}^{N}\left\{-\mathbf{\Psi}^{\mathrm T}+\left(\mathbf{W} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]^{\mathrm T} \mathbf{W}^{\mathrm T}\right)-2 \mathbf{W} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T}\right\}+\frac{N}{2} \mathbf{S}^{\mathrm T}=0 $$ $$ N \mathbf{\Psi}^{\mathrm T}-\mathbf{W}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]^{\mathrm T}\right] \mathbf{W}^{\mathrm T}+2 \mathbf{W}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T}\right] -N \mathbf{S}^{\mathrm T}=0 $$ 転置をとり、$\mathbf{\Psi}$を分離しつつ、$\mathbf{W} \to \mathbf{W}_{\text{new}}, \mathbf{\Psi} \to \mathbf{\Psi}_{\text{new}}$とすると(対称行列なので$\mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right] = \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]^{\mathrm T}$である) $$ \mathbf{\Psi}_{\text {new }}=\mathbf{S}-\frac{2}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T}\right]+\frac{1}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]\right] \mathbf{W}_{\text{new}}^{\mathrm T} $$ $(12.69)$を上式の最後の$\mathbf{W}_{\text{new}}^{\mathrm T}$にのみ適用させると $$ \begin{aligned} \mathbf{\Psi}_{\text {new }} &=\mathbf{S}-\frac{2}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T}\right] \\ &+\frac{1}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]\right]\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n} \mathbf{z}_{n}^{\mathrm T}\right]\right]^{-\mathrm T}\left[\sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right) \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T}\right]^{\mathrm T} \\ &=\mathbf{S}-\frac{2}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n}\right]^{\mathrm T}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\right]+\frac{1}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T}\right] \\ &=\mathbf{S}-\frac{1}{N} \mathbf{W}_{\text{new}}\left[\sum_{n=1}^{N} \mathbb{E}\left[\mathbf{z}_{n}\right]\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm T}\right] \end{aligned} $$ $\mathbf{\Psi}$はもとより対角行列なので、$\textrm{diag}$演算子をつければ$(12.70)$式を得る。 ## 演習 12.23 <div class="panel-primary"> 確率的主成分分析モデルの離散個の混合を考え，その確率的有向グラフィカルモデルを描け．個々の主成分分析モデルは$\mathbf{W},\boldsymbol{\mu},\sigma^2$という独自のパラメータ値を持つ．次に，これらのパラメータ値が混合モデルの各要素に共有される場合のグラフを描け． </div> 離散個（$K$個）の主成分分析モデルが存在し、混合係数$\pi_{k}$によって1つの確率的主成分分析モデルが構成されていることを考える。 前半は「個々の主成分分析モデルは$\mathbf{W},\boldsymbol{\mu},\sigma^2$という独自のパラメータ値を持つ」場合、後半はこれらが共通している場合を考える。 確率的有向グラフィカルモデルは以下の通り。左が独立な場合。右が共通の場合。 ![](https://i.imgur.com/QqRtN8I.png) （公式解答によれば$\mathbf{s}$はパラメータ$\boldsymbol{\pi}$によって規定される、各主成分分析モデルの混合比を規定する$K$値の潜在変数……と言っているが、これまで扱ってきた混合係数$\pi_{k}$とは違うんだろうか？） ## 演習 12.24 <div class="panel-primary"> 2.3.7節において，スチューデントのt分布が，連続的潜在変数について周辺化されたガウス分布の無限個の混合とみなされることを学んだ．この表現を利用して，多変数のスチューデントt分布の対数尤度関数を，与えられた観測データ集合に対して最大化し， EステップとMステップの形を導け． </div> (2.161)式の離散変数バージョンで$\eta$を潜在変数とみなすと、完全データ対数尤度関数は、 \begin{align} \ln p(\mathbf{X}, \boldsymbol{\eta} | \boldsymbol\mu , \boldsymbol\Lambda , \nu ) = \sum _{n=1}^N \left\{ \ln \mathcal{N} \big( \mathbf{x}_n | \boldsymbol\mu , (\eta_n \boldsymbol\Lambda ) ^{-1} \big) + \ln {\rm Gam} (\eta_n |\frac{\nu}{2},\frac{\nu}{2}) \right\} \end{align} となる。（$\boldsymbol{\eta}$は、$N$個の潜在変数$\eta_n$を並べたベクトル） $\boldsymbol{\eta}$は潜在変数なので$\boldsymbol{\eta}$で期待値を取ると、 \begin{align} \mathbb{E}_\boldsymbol\eta [\ln p (\mathbf{X}, \boldsymbol\eta | \boldsymbol\mu , \boldsymbol\Lambda , \nu )] = -\frac{1}{2} \sum _{n=1}^N \left\{ D (\ln(2\pi)-\mathbb{E} [\ln\eta_n])-\ln|\Lambda| +\mathbb{E} [\eta_n] (\mathbf{x}_n^{\rm T}\boldsymbol\Lambda\mathbf{x}_n-2\mathbf{x}_n^{\rm T}\boldsymbol\Lambda\boldsymbol\mu+\boldsymbol\mu^{\rm T}\boldsymbol\Lambda\boldsymbol\mu)\\ +2\ln\Gamma \left(\frac{\nu}{2}\right) -\nu (\ln \nu - \ln 2 ) - (\nu -2) \mathbb{E} [\ln\eta_n] + \nu\mathbb{E} [\eta_n] \right\} \tag{*} \end{align} と、(12.53)式に対応する式が導かれる。（第１行はガウス分布の展開から、第２行はガンマ分布の展開から出てくる。） 次に、Eステップでの$\mathbb{E} [\eta_n]$、$\mathbb{E} [\ln\eta_n]$の更新式を導く。$\boldsymbol\eta$の確率分布は、 \begin{align} p(\boldsymbol{\eta} | \mathbf{X}, \boldsymbol\mu , \boldsymbol\Lambda , \nu ) &= \prod _{n=1}^N p(\eta_n | \mathbf{x}_n, \boldsymbol\mu , \boldsymbol\Lambda , \nu ) \\ &\propto \prod _{n=1}^N p(\mathbf{x}_n | \eta_n, \boldsymbol\mu , \boldsymbol\Lambda , \nu ) p( \eta_n | \boldsymbol\mu , \boldsymbol\Lambda , \nu )\\ &= \prod _{n=1}^N \mathcal{N} (\mathbf{x}_n | \boldsymbol\mu , (\eta_n \boldsymbol\Lambda ^{-1})) {\rm Gam} (\eta_n |\frac{\nu}{2},\frac{\nu}{2}) \end{align} となる。ガウス分布の精度$\eta_n$の事前分布がGam$(\eta_n |\nu/2,\nu/2)$であるときの事後分布を表すから、積の中身はガンマ分布Gam$(\eta_n | a_n, b_n)$に従う。パラメータ$a_n, b_n$の値は、(2.150)式と(2.151)式より、 \begin{align} a_n &= \frac{\nu+D}{2}\\ b_n &= \frac{\nu + (\mathbf{x}_n - \boldsymbol\mu)^{\rm T}\boldsymbol\Lambda (\mathbf{x}_n - \boldsymbol\mu) }{2} \end{align} である。このパラメータ$a_n, b_n$を用いると、ガンマ分布の公式(B.27)と(B.30)より、 \begin{align} \mathbb{E} [\eta_n] &= \frac{a_n}{b_n}\\ \mathbb{E} [\ln\eta_n] &= \psi (a_n)-\ln b_n \end{align} と書ける。（$\psi( \cdot )$は、(B.25)で定義されるディガンマ関数） **Mステップの更新式(1/3)** 上記の(*)式を$\boldsymbol{\mu}$で微分した値をゼロとおくと、 \begin{align} \frac{\partial}{\partial\mu} \mathbb{E} _\boldsymbol\eta [\ln p (\mathbf{X}, \boldsymbol\eta | \boldsymbol\mu , \boldsymbol\Lambda , \nu )] &= -\frac{1}{2} \sum_{n=1}^N \{ \mathbb{E} [\eta_n] (-2\boldsymbol\Lambda \mathbf{x}_n + 2 \boldsymbol\Lambda \boldsymbol\mu )\}\\ &= \boldsymbol\Lambda \left( \sum_{n=1}^N \mathbb{E} [\eta_n] \mathbf{x}_n - \boldsymbol\mu \sum_{n=1}^N \mathbb{E}[\eta_n] \right)=0 \end{align} $\boldsymbol{\Lambda}$は精度行列なので、逆行列（＝分散行列）が存在する。その逆行列を左から乗じて整理すると、 \begin{align} \boldsymbol\mu_{\rm ML} = \frac{\sum_{n=1}^N \mathbb{E} [\eta_n] \mathbf{x}_n }{\sum_{n=1}^N \mathbb{E}[\eta_n]} \end{align} **Mステップの更新式(2/3)** 上記の(*)式を$\boldsymbol{\Lambda}$で微分した値をゼロと置くのだが、準備として任意のベクトル$\boldsymbol\xi,\boldsymbol\zeta$に対して、 \begin{align} \frac{\partial}{\partial \boldsymbol\Lambda} (\boldsymbol\xi^{\rm T}\boldsymbol\Lambda\boldsymbol\zeta ) &= \frac{\partial}{\partial \boldsymbol\Lambda} {\rm Tr} (\boldsymbol\xi^{\rm T}\boldsymbol\Lambda\boldsymbol\zeta )\\ &= \frac{\partial}{\partial \boldsymbol\Lambda} {\rm Tr} (\boldsymbol\Lambda\boldsymbol\zeta \boldsymbol\xi^{\rm T})\\ &= \boldsymbol\zeta \boldsymbol\xi^{\rm T} \tag{C.25より} \end{align} であることに注目すると、 \begin{align} \frac{\partial}{\partial\boldsymbol\Lambda} \mathbb{E} _\boldsymbol\eta [\ln p (\mathbf{X}, \boldsymbol\eta | \boldsymbol\mu , \boldsymbol\Lambda , \nu )] &= -\frac{1}{2} \sum_{n=1}^N \left\{ -\Lambda^{-1} + \mathbb{E} [\eta_n] (\mathbf{x}_n \mathbf{x}_n^{\rm T} -2 \boldsymbol\eta \mathbf{x}_n^{\rm T} + \boldsymbol\eta\boldsymbol\eta^{\rm T}) \right\} =0\\ &\Leftrightarrow \Lambda_{\rm ML}= \left( \frac{1}{N}\sum_{n=1}^N \mathbb{E} [\eta_n] (\mathbf{x}_n \mathbf{x}_n^{\rm T} -2 \boldsymbol\eta \mathbf{x}_n^{\rm T} + \boldsymbol\eta\boldsymbol\eta^{\rm T}) \right) ^{-1}\\ \end{align} **Mステップの更新式(3/3)** 上記の(*)式を$\nu$で微分した値をゼロと置くと、 \begin{align} \frac{\partial}{\partial\nu} \mathbb{E} _\boldsymbol\eta [\ln p (\mathbf{X}, \boldsymbol\eta | \boldsymbol\mu , \boldsymbol\Lambda , \nu )] &= -\frac{1}{2} \sum_{n=1}^N \left\{ \frac{2}{\Gamma (\frac{\nu}{2})}\frac{\partial}{\partial\nu}\Gamma(\nu/2) - (\ln\nu - \ln 2) -1 - \mathbb{E} [\ln \eta_n] + \mathbb{E} [\eta_n] \right\} =0\\ &\Leftrightarrow \frac{\psi(\nu/2)}{\Gamma (\nu/2)} - (\ln\nu - \ln 2) -1 - \frac{1}{N}\sum_{n=1}^N (\mathbb{E} [\ln \eta_n] - \mathbb{E} [\eta_n] )=0 \end{align} となる。これは$\nu$について解析的に解くことはできない。 ## 演習 12.25 <div class="panel-primary"> 潜在変数の空間分布$p(\mathbf{z}) = \mathcal{N}(\mathbf{x}\mid \mathbf{0}, \mathbf{I})$を持つ線形ガウス潜在変数モデルを考える．観測変数に対する条件付き分布を$p(\mathbf{x} \mid \mathbf{z})=\mathcal{N}(\mathbf{x} \mid \mathbf{W z}+\boldsymbol{\mu}, \mathbf{\Phi})$とおく．ただし$\mathbf{\Phi}$はノイズ項の対称かつ正定値の任意の共分散行列である．$D \times D$行列$\mathbf{A}$を使ってデータ変数に非特異的な線形変換$\mathbf{x} \to \mathbf{Ax}$を行うものとする．もし$\boldsymbol{\mu}_{\mathrm{ML}}$と$\mathbf{W}_{\mathrm{ML}}$と$\mathbf{\Phi}_{\mathrm{ML}}$がもともとの非変換データに対応した最尤解を表すものとすれば，$\mathbf{A} \mu_{\mathrm{ML}}$と$\mathbf{AW}_{\mathrm{ML}}$と$\mathbf{A} \mathbf{\Phi}_{\mathrm{ML}} \mathbf{A}^{\mathrm{T}}$が。変換されたデータ集合の最尤解に対応していることを示せ．最後に，以下の2つの場合にモデルの形が保存されることを示せ．(i) $\mathbf{A}$が対角行列で，$\mathbf{\Phi}$も対角行列．これは因子分析に対応する．変換された$\mathbf{\Phi}$は対角的なままであり，したがって，因子分析モデルは，データ変数の要素ごとの尺度変更に**共変的(covariant)** である．(ii) $\mathbf{A}$が対角的で，$\mathbf{\Phi}$が単位行列に比例する場合，つまり$\mathbf{\Phi} = \sigma^2\mathbf{I}$である場合，これは確率的主成分分析に対応している．変換された$\mathbf{\Phi}$行列は単位行列に比例したままとなり，それゆえ確率的主成分分析は通常の主成分分析がそうであるように，データ空間の軸の回転の下で共変的である． </div> ※ 確率的主成分分析の対数尤度関数は12.2節の議論から以下のようにかける． $$ \begin{aligned} L(\boldsymbol{\mu}, \mathbf{W}, \mathbf{\Phi})=-\frac{N D}{2} \ln (2 \pi)-\frac{N}{2} \ln \left|\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right| -\frac{1}{2} \sum_{n=1}^{N}\left\{\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)^{\mathrm{T}}\left(\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right)^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}\right)\right\} \end{aligned} $$ いま$\mathbf{x} \to \mathbf{Ax}$の変換を考えているので $$ \begin{aligned} &L_{\mathbf{A}}(\boldsymbol{\mu}, \mathbf{W}, \mathbf{\Phi})=-\frac{N D}{2} \ln (2 \pi)-\frac{N}{2} \ln \left|\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right|-\frac{1}{2} \sum_{n=1}^{N}\left\{\left(\mathbf{A} \mathbf{x}_{n}-\boldsymbol{\mu}\right)^{\mathrm{T}}\left(\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right)^{-1}\left(\mathbf{A} \mathbf{x}_{n}-\boldsymbol{\mu}\right)\right\} \end{aligned} $$ 以上のように書き換えることができる．ここで$\boldsymbol{\mu}$の最尤解を求めるので$\boldsymbol{\mu}$についての微分を0とおいて$\boldsymbol{\mu}$について解くと $$ \boldsymbol{\mu}_{\mathrm{A}}=\frac{1}{N} \sum_{n=1}^{N} \mathbf{A} \mathbf{x}_{n}=\mathbf{A} \overline{\mathbf{x}}=\mathbf{A} \boldsymbol{\mu}_{\mathrm{ML}} $$ が得られる．これを対数尤度関数に代入して$\displaystyle \mathbf{S}=\frac{1}{N} \sum_{n=1}^{N}\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)\left(\mathbf{x}_{n}-\overline{\mathbf{x}}\right)^{\mathrm{T}}$を用いると $$ \begin{aligned} L_{\mathbf{A}}(\boldsymbol{\mu}, \mathbf{W}, \mathbf{\Phi})=-\frac{N D}{2} \ln (2 \pi)-\frac{N}{2} \ln \left|\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right|-\frac{1}{2} \sum_{n=1}^{N} \operatorname{Tr}\left\{\left(\mathbf{W} \mathbf{W}^{\mathrm{T}}+\mathbf{\Phi}\right)^{-1} \mathbf{A} \mathbf{S} \mathbf{A}^{\mathrm{T}}\right\} \end{aligned} $$ 以上のように書ける．ここで以下の定義を用いることで $$ \mathbf{\Phi}_{\mathbf{A}}=\mathbf{A} \mathbf{\Phi}^{-1} \mathbf{A}^{\mathrm{T}}, \quad \mathbf{W}_{\mathbf{A}}=\mathbf{A} \mathbf{W} $$ 対数尤度関数の最終項について変換$A$をあらわに用いない形で書き換えることができる $$ \begin{aligned} &L_{\mathbf{A}}\left(\boldsymbol{\mu}_{\mathbf{A}}, \mathbf{W}_{\mathbf{A}}, \mathbf{\Phi}_{\mathbf{A}}\right)=-\frac{N D}{2} \ln (2 \pi)-\frac{N}{2} \ln \left|\mathbf{W}_{\mathbf{A}} \mathbf{W}_{\mathbf{A}}^{\mathrm{T}}+\mathbf{\Phi}_{\mathbf{A}}\right|-\frac{1}{2} \sum_{n=1}^{N}\left\{\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathbf{A}}\right)^{\mathrm{T}}\left(\mathbf{W}_{\mathbf{A}} \mathbf{W}_{\mathbf{A}}^{\mathrm{T}}+\mathbf{\Phi}_{\mathbf{A}}\right)^{-1}\left(\mathbf{x}_{n}-\boldsymbol{\mu}_{\mathbf{A}}\right)\right\}-N \ln |\mathbf{A}| \end{aligned} $$ いまこの式は最初の対数尤度関数と比べて$-\ln A$の項だけが異なっている．従って変換$A$をによるデータセットの$\mathbf{W}, \mathbf{\Phi}$の最尤解は変換前の最尤解と同じ形で書け， $$ \mathbf{W}_{\mathbf{A}}=\mathbf{AW}_{\mathrm{ML}} $$ $$ \mathbf{\Phi}_{\mathbf{A}}=\mathbf{A} \mathbf{\Phi}_{\mathrm{ML}} \mathbf{A}^{\mathrm{T}} $$ となることがわかる． ． 次にこの変換の前後で$\mathbf{\Phi}$についての制約が保たれるか確認する． (i) $\mathbf{A}$が対角行列で$\mathbf{\Phi}$ も対角行列のとき(因子分析に対応) このとき $$ \mathbf{\Phi}_{A}=\mathbf{A}\mathbf{\Phi}^{-1}\mathbf{A}^{\mathrm{T}} $$ と書くことができて，これは対角行列同士の積であるため，$\mathbf{\Phi}_A$も対角行列であることがわかる． 従ってこの場合データベクトルの要素ごとの尺度を変更すると対応する$\mathbf{\Phi}_A$の要素の尺度変更にデータの尺度変更が吸収されることがわかる． (ii) $\mathbf{A}$が対角行列で$\mathbf{\Phi}$ が単位行列に比例する場合(確率的主成分分析に対応) このとき $$ \mathbf{\Phi}_{A}=\mathbf{A}\mathbf{\Phi}^{-1}\mathbf{A}^{\mathrm{T}} $$ であることから変換$\mathbf{A}$の前後ともに$\mathbf{\Phi}=\sigma^2\mathbf{I}$を満たすとき $$ \mathbf{A}\mathbf{A}^{\mathrm{T}}=\mathbf{I} $$ となるので$\mathbf{A}$は直交行列である必要がある．このとき変換$\mathbf{A}$は座標系の回転に相当する ![](https://i.imgur.com/dNe3lFQ.png) ## 演習 12.26 <div class="panel-primary"> $$ \mathbf{K} \mathbf{a}_{i}=\lambda_{i} N \mathbf{a}_{i} \tag{12.80} $$ を満たす任意のベクトル$\mathbf{a}_i$が $$ \mathbf{K}^{2} \mathbf{a}_{i}=\lambda_{i} N \mathbf{K} \mathbf{a}_{i} \tag{12.79} $$ も満たすことを示せ．また固有値$\lambda$を持つ$(12.80)$の任意の解に，$\mathbf{K}$のゼロ固有値に属する固有ベクトルのいかなる定数倍を加えたものもまた固有値$\lambda$を持つ$(12.79)$の解であることを示せ．最後にゼロ固有値に属する固有ベクトルの成分を加える変更が， $$ y_{i}(\mathbf{x})=\boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \mathbf{v}_{i}=\sum_{n=1}^{N} a_{i n} \boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \boldsymbol{\phi}\left(\mathbf{x}_{n}\right)=\sum_{n=1}^{N} a_{i n} k\left(\mathbf{x}, \mathbf{x}_{n}\right) \tag{12.82} $$ で与えられる主成分による射影に影響しないことを示せ． </div> ※ （復習）$\mathbf{K}$は6章のP.3で定義されているグラム行列であり、$n$番目の行が$\boldsymbol{\phi}(\mathbf{x}_{n}^{\mathrm T})$で表される計画行列$\mathbf{\Phi}$に対して$\mathbf{K} = \mathbf{\Phi}\mathbf{\Phi}^{\mathrm T}$で定義される。これは$N \times N$の対称行列であり、その要素は $$ K_{n m}=\boldsymbol{\phi}\left(\mathbf{x}_{n}\right)^{\mathrm{T}} \boldsymbol{\phi}\left(\mathbf{x}_{m}\right)=k\left(\mathbf{x}_{n}, \mathbf{x}_{m}\right) \tag{6.6} $$ である。 $(1)$ $(12.80)$が成立する前提で$(12.79)$が成立することを示す。$(12.79)$の左辺に$(12.80)$を代入すると $$ \mathbf{K}^2\mathbf{a}_{i} = \mathbf{K}(\lambda_{i}N\mathbf{a}_{i}) = \lambda_{i}N\mathbf{Ka}_{i} $$ となるので成立する。 $(2)$ $(12.80)$の解のうち、固有値$\lambda_{i}$に対応する固有ベクトルを$\mathbf{a}_{i}$とし、$\mathbf{K}$のゼロ固有値に属する固有ベクトルを$\mathbf{b}_{i}$とおく。このとき定数$m$を用いた線形結合のベクトル$\widetilde{\mathbf{a}_{i}} = \mathbf{a}_{i} + m\mathbf{b}_{i}$を$(12.79)$の左辺に代入したとき、右辺が得られることを示す。 設定から$\mathbf{Kb}_{i} = \lambda_{i}N\mathbf{b}_{i}( = \mathbf{0})$となるので $$ \begin{aligned} \mathbf{K}^{2}\widetilde{\mathbf{a}_{i}} &= \mathbf{K}^{2}(\mathbf{a}_{i} + m\mathbf{b}_{i}) \\ &=\mathbf{K}^{2}\mathbf{a}_{i} + m\mathbf{K}^{2}\mathbf{b}_{i} \\ &=\lambda_{i}N\mathbf{K}\mathbf{a}_{i} + m\mathbf{K}(\lambda_{i}N\mathbf{b}_{i}) \\ &=\lambda_{i}N\mathbf{K}(\mathbf{a}_{i} + m\mathbf{b}_{i}) \\ &=\lambda_{i}N\mathbf{K}\widetilde{\mathbf{a}_{i}} \end{aligned} $$ よって$\widetilde{\mathbf{a}_{i}}$は$(12.79)$の解であることが示された。 $(3)$ $(12.82)$の$a_{in}$を$\widetilde{\mathbf{a}_{i}}$の成分$\tilde{a}_{in}$に置き換えても$(12.82)$の結果が変化しないことを示す。$(2)$で用いた$\mathbf{K}\mathbf{b}_{i}=0$より、$\sum_{n=1}^{N} b_{in}k(\mathbf{x}, \mathbf{x}_{n}) = 0$であるから、 $$ \begin{aligned} \boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \tilde{\mathbf{v}}_{i}&=\sum_{n=1}^{N} \tilde{a}_{i n} \boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \boldsymbol{\phi}\left(\mathbf{x}_{n}\right) \\ &=\sum_{n=1}^{N} (a_{in}+mb_{in}) \boldsymbol{\phi}(\mathbf{x})^{\mathrm{T}} \boldsymbol{\phi}\left(\mathbf{x}_{n}\right) \\ &=\sum_{n=1}^{N} a_{in}k(\mathbf{x}, \mathbf{x}_{n}) + m\sum_{n=1}^{N} b_{in}k(\mathbf{x}, \mathbf{x}_{n})\\ &=\sum_{n=1}^{N} a_{in}k(\mathbf{x}, \mathbf{x}_{n}) \end{aligned} $$ よって変化しないことが示された。 ## 演習 12.27 <div class="panel-primary"> $k(\mathbf{x}, \mathbf{x}^{\prime}) = \mathbf{x}^{\mathrm T}\mathbf{x}^{\prime}$で与えられる線形のカーネル関数を選択すれば，カーネル主成分分析の特別な場合として普通の主成分分析アルゴリズムが再現できることを示せ． </div> ※カーネル関数$k(\mathbf{x},\mathbf{x}^{\prime}) = \mathbf{x}^{\mathrm T}\mathbf{x}^{\prime}$を用いて、あるデータ点$\mathbf{x}_{m}$について$(12.82)$を用いた固有ベクトル$i$への射影 $$ y_{i}(\mathbf{x}) = \sum_{n=1}^{N} a_{i n} k\left(\mathbf{x}, \mathbf{x}_{n}\right)\tag{12.82} $$ が、通常の主成分分析での射影と同じ結果、すなわち、あるデータ点$\mathbf{x}_{m}$が固有ベクトル$\mathbf{u}_{i}$へ射影された$\mathbf{u}_{i}^{\mathrm T}\mathbf{x}_{m}$に一致することを示せば良い。これは $$ \begin{aligned} \sum_{n=1}^{N} a_{i n} k\left(\mathbf{x}_{m}, \mathbf{x}_{n}\right) &=\sum_{n=1}^{N} a_{m} \mathbf{x}_{m}^{\mathrm T} \mathbf{x}_{n} \\ &=\left(\sum_{n=1}^{N} a_{i n} \mathbf{x}_{n}^{\mathrm T}\right) \mathbf{x}_{m} \\ &\equiv {\mathbf{u}_{i}^{\star}}^{\mathrm T}\mathbf{x}_{m} \end{aligned} $$ この$\displaystyle {\mathbf{u}_{i}^{\star}}\equiv \sum_{n=1}^{N} a_{i n} \mathbf{x}_{n}$が、通常の主成分分析の主成分の定義である$\mathbf{Su}_{i} = \lambda_{i}\mathbf{u}_{i}$によって定義される$\mathbf{u}_{i}$と一致することを示す。左辺から展開していって $$ \begin{aligned} \mathbf{Su}_{i} &=\left(\frac{1}{N} \sum_{m=1}^{N} \mathbf{x}_{m} \mathbf{x}_{m}^{\mathrm T}\right)\left(\sum_{n=1}^{N} a_{i n} \mathbf{x}_{n}\right) \\ &=\frac{1}{N} \sum_{m=1}^{N} \mathbf{x}_{m} \sum_{n=1}^{N} a_{i n} \mathbf{x}_{m}^{\mathrm T} \mathbf{x}_{n} \\ &=\frac{1}{N} \sum_{m=1}^{N} \mathbf{x}_{m} \sum_{n=1}^{N} a_{i n} k\left(\mathbf{x}_{m}, \mathbf{x}_{n}\right) \\ &=\frac{1}{N} \sum_{m=1}^{N} \mathbf{x}_{m} \mathbf{K} \mathbf{a}_{i} \\ &=\frac{1}{N} \sum_{m=1}^{N} \mathbf{x}_{m} \lambda_{i}N \mathbf{a}_{i} \quad(\because(12.80)) \\ &=\lambda_{i} \sum_{m=1}^{N} a_{i m} \mathbf{x}_{m} \\ &=\lambda_{i} \mathbf{u}_{i}^{\star} \end{aligned} $$ 以上から、題意が示された。 ## 演習 12.28 <div class="panel-primary"> 関数$f(x)$による非線形変数変換$y=f(x)$を施すことにより，ある固定された密度$q(x)$から，任意の密度関数$p(y)$が得られることを示せ．また$f(x)$が満足する微分方程式を書き下し，密度の変換を説明する図を描け．ただし$q(x)$はいたるところで非ゼロとし，$f(x)$は単調と仮定する．したがって$0\leqslant f^{\prime}(x) \lt \infty$である．上記を示す際，変数変換による確率密度の変換に関する性質 $$ \begin{aligned} p_{y}(y) &=p_{x}(x)\left|\frac{\mathrm{d} x}{\mathrm{~d} y}\right| \\ &=p_{x}(g(y))\left|g^{\prime}(y)\right| \end{aligned} \tag{1.27} $$ を使え． </div> 固定された密度（確率密度）は$q(x) \gt 0$（非ゼロ）で、非線形変数変換$f(x)$は$0\leqslant f^{\prime}(x) \lt \infty$より単調増加関数。単調増加なので、逆関数が存在し、$x = f^{-1}(y)$が存在する。このとき1.2.1 確率密度や演習1.4でも見たように、$(1.27)$の変数変換を用いて $$ p(y) = q(f^{-1}(y))\left|\frac{\mathrm{d} f^{-1}(y)}{\mathrm{~d} y}\right| $$ と書くことができる。 今、$f$の制約は単調であることのみなので、$x$に対する確率質量を$y$に任意に配分することができる。密度の変換はこんなイメージ（演習問題1.4のときの図を使用） ![](https://i.imgur.com/FC9IHEz.png) $(1.27)$について赤いガウス分布$p_{x}(x)$について$x=g(y)$（この場合はシグモイド曲線）で変換すると、緑の曲線$p_{x}(g(y))$が得られる。しかし$(1.27)$でやっているのはそれにさらに$|g^{\prime}(y)|$をかけており、これで微小区間$dx, dy$を用いて$p_y(y)dy = p_{x}(x)dx$の微小面積を保存するように変換している。これにより、紫の曲線が得られる。この問題では$p_{x}(x)$が$q(x)$、$g^{-1}(x) = f(x)$, $p_{y}(y) = p(y)$である。$f(x)$の傾きを$0\leqslant f^{\prime}(x) \lt \infty$の範囲で任意に調節すれば、自在に$p(y)$の形状を設定できる。 $y=f(x)$が満足する微分方程式は、上式を変形して $$ \left| \frac{dy}{dx}\right| = \left| \frac{dy}{df^{-1}(y)}\right| = \frac{q(f^{-1}(y))}{p(y)} = \frac{q(x)}{p(f(x))} $$ である。 ## 演習 12.29 <div class="panel-primary"> 2つの変数$z_1$と$z_2$が独立で，$p(z_1,z_2) = p(z_1)p(z_2)$となると仮定する．これらの2つの変数に対する共分散行列が対角的であることを示せ．これは変数の独立性が，2つの変数が無相関となることの十分条件であることを示している．次に，2つの変数$y_1$と$y_2$を考え，$y_1$が$0$を中心に対称に分布しており，また，$y_2 = {y_1}^2$を満たすと仮定する．条件付き分布$p(y_2\mid y_1)$を書き下し，これが$y_1$に依存しており，2つの変数が独立とはならないことを示せ．次に，これら2つの変数に対する共分散行列が対角的であることを示せ．これを示すために，関係式$p\left(y_{1}, y_{2}\right)=p\left(y_{1}\right) p\left(y_{2} \mid y_{1}\right)$を用いて，非対角成分が$0$となることを示せ．この反例は，相関が$0$であることが，独立性の十分条件にはならないことを示している． </div> ※前半は演習1.6とほぼ同じ $z_1, z_2$の共分散$\operatorname{cov}[z_1,z_2]$を求める。 $$ \begin{aligned} \operatorname{cov}\left[z_{1}, z_{2}\right] &=\mathbb{E}\left[\left(z_{1}-\mathbb{E}\left[z_{1}\right]\right)\left(z_{2}-\mathbb{E}\left[z_{2}\right]\right)\right] \\ &=\mathbb{E}\left[z_{1} z_{2}\right]-\mathbb{E}\left[z_{1}\right] \mathbb{E}\left[z_{2}\right] \\ &=\iint z_{1} z_{2} p\left(z_{1}, z_{2}\right) d z_{1} d z_{2}-\mathbb{E}\left[z_{1}\right] \mathbb{E}\left[z_{2}\right] \\ &=\iint z_{1} z_{2} p\left(z_{1}\right) p\left(z_{2}\right) d z_{1} d z_{2}\quad (\because p(z_1,z_2) = p(z_1)p(z_2))\\ &=\int z_{1} p\left(z_{1}\right) d z_{1} \int z_{2} p\left(z_{2}\right) d z_{2}-\mathbb{E}\left[z_{1}\right] \mathbb{E}\left[z_{2}\right] \\ &=\mathbb{E}\left[z_{1}\right] \mathbb{E}\left[z_{2}\right]-\mathbb{E}\left[z_{1}\right] \mathbb{E}\left[z_{2}\right] \\ &=0 \end{aligned} $$ $z_1, z_2$の共分散行列は $$ \begin{aligned} \Sigma &=\left(\begin{array}{cc}\operatorname{var}\left[z_{1}\right] & \operatorname{cov}\left[z_{1}, z_{2}\right] \\ \operatorname{cov}\left[z_{1}, z_{2}\right] & \operatorname{var}\left[z_{2}\right]\end{array}\right) \\ &=\left(\begin{array}{cc}\operatorname{var}\left[z_{1}\right] & 0 \\ 0 & \operatorname{var}\left[z_{2}\right]\end{array}\right) \end{aligned} $$ となるので、共分散行列が対角的であることになる。 後半は、$\operatorname{cov}[y_1,y_2]=0$であっても$y_1, y_2$が独立でないことがあることを示す問題。 問題設定から、$y_1$が$0$を中心に対称に分布し（$\mathbb{E}[y_1]=0$）、かつ$y_2 = {y_1}^{2}$の関係が成立する2変数について考えると、$y_1$が与えられた下での$y_2$の条件つき確率$p(y_2\mid y_1)$は、性質上$y_2 \neq y_1^2$ならば確率$0$で$y_2 = y_1^2$に常に存在するので、数式上 $$ p(y_2 \mid y_1) = \delta(y_2-y_1^{2}) $$ となる。 これは明らかに$p(y_2\mid y_1)$が$y_1$に依存しているので、$y_1$と$y_2$は独立ではない。この共分散を求めると $$ \begin{aligned} \operatorname{cov}\left[y_{1}, y_{2}\right] &=\mathbb{E}\left[y_{1}, y_{2}\right]-\mathbb{E}\left[y_{1}\right] \mathbb{E}\left[y_{2}\right] \\ &=\iint y_{1} y_{2} p\left(y_{1}, y_{2}\right) d y_{1} d y_{2}-0 \quad (\because \mathbb{E}[y_1]=0)\\ &=\iint \delta\left(y_{2}-y_{1}^{2}\right) p\left(y_{1}\right) y_{1} y_{2} d y_{1} d y_{2} \\ &=\int y_{1} p\left(y_{1}\right) \underbrace{\int \delta\left(y_{2}-y_{1}^{2}\right) y_{2} d y_{2}}_{y_1^2} d y_{1} \\ &=\int y_{1} p\left(y_{1}\right) y_{1}^{2} d y_{1} \\ &=\int y_{1}^{3} p\left(y_{1}\right) d y_{1} \end{aligned} $$ 最後に、$p(y_1)$は偶関数で$y_1^3$は奇関数であるから、この積分は$0$となる。すなわち、$y_1$と$y_2$は独立ではない場合でも共分散は$0$となることがあることが示された。