🧠

PRML 第12章(12.16から12.29まで)解答例

2022/05/15に公開

はじめに

PRML解答例まとめを参照

演習 12.17

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)が出ることを示せ.


(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

12.2.4節で説明された因子分析モデルについて,独立なパラメータの数の表現を与える式を導け.


因子分析における共分散行列 (式12.65) より、独立なパラメータの数は

  • D \times M行列\mathbf{W}

及び

  • M \times M次元の対角行列\mathbf{\Psi}

に依存する。また、確率的主成分分析における議論 (pp.293-294) と同様に、M(M - 1)/2個のパラメータは回転に関して冗長。つまり、因子分析モデルにおける独立なパラメータの数は

\begin{aligned} & DM + D - M(M - 1)/2 \\ =& D(M + 1) - M(M - 1)/2 \end{aligned}

となる。

演習 12.19

12.2.4節で説明した因子分析モデルが潜在変数空間の座標の回転に対して不変であることを示せ.


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

2次導関数を考えることによって12.2.4節で説明した因子分析モデルの対数尤度関数の,パラメータ\boldsymbol{\mu}に対する唯一の停留点が,

\overline{\mathbf{x}}=\frac{1}{N} \sum_{n=1}^{N} \mathbf{x}_{n} \tag{12.1}

で定義されたサンプル平均で与えられることを示せ.さらに,この停留点が最大値を与えることを示せ.


※対数尤度関数\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

因子分析についての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}}で置き換えられることに注意せよ.


事前分布は、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

因子分析モデルの完全データ対数尤度関数の期待値の式を書き下せ.また,対応する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}

を導け.


※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

確率的主成分分析モデルの離散個の混合を考え,その確率的有向グラフィカルモデルを描け.個々の主成分分析モデルは\mathbf{W},\boldsymbol{\mu},\sigma^2という独自のパラメータ値を持つ.次に,これらのパラメータ値が混合モデルの各要素に共有される場合のグラフを描け.


離散個(K個)の主成分分析モデルが存在し、混合係数\pi_{k}によって1つの確率的主成分分析モデルが構成されていることを考える。

前半は「個々の主成分分析モデルは\mathbf{W},\boldsymbol{\mu},\sigma^2という独自のパラメータ値を持つ」場合、後半はこれらが共通している場合を考える。

確率的有向グラフィカルモデルは以下の通り。左が独立な場合。右が共通の場合。

(公式解答によれば\mathbf{s}はパラメータ\boldsymbol{\pi}によって規定される、各主成分分析モデルの混合比を規定するK値の潜在変数……と言っているが、これまで扱ってきた混合係数\pi_{k}とは違うんだろうか?)

演習 12.24

2.3.7節において,スチューデントのt分布が,連続的潜在変数について周辺化されたガウス分布の無限個の混合とみなされることを学んだ.この表現を利用して,多変数のスチューデントt分布の対数尤度関数を,与えられた観測データ集合に対して最大化し, EステップとMステップの形を導け.


(2.161)式の離散変数バージョンで\etaを潜在変数とみなすと、完全データ対数尤度関数は、

\begin{aligned} \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{aligned}

となる。(\boldsymbol{\eta}は、N個の潜在変数\eta_nを並べたベクトル)

\boldsymbol{\eta}は潜在変数なので\boldsymbol{\eta}で期待値を取ると、

\begin{aligned} \mathbb{E}_{\eta}[\ln p(\mathbf{X}, \boldsymbol{\eta} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)] &=-\frac{1}{2} \sum_{n=1}^{N}\left\{D\left(\ln (2 \pi)-\mathbb{E}\left[\ln \eta_{n}\right]\right)-\ln |\Lambda|+\mathbb{E}\left[\eta_{n}\right]\left(\mathbf{x}_{n}^{\mathrm{T}} \boldsymbol{\Lambda} \mathbf{x}_{n}-2 \mathbf{x}_{n}^{\mathrm{T}} \boldsymbol{\Lambda} \boldsymbol{\mu}+\boldsymbol{\mu}^{\mathrm{T}} \boldsymbol{\Lambda} \boldsymbol{\mu}\right)(*)\right.\\ &\left.+2 \ln \Gamma\left(\frac{\nu}{2}\right)-\nu(\ln \nu-\ln 2)-(\nu-2) \mathbb{E}\left[\ln \eta_{n}\right]+\nu \mathbb{E}\left[\eta_{n}\right]\right\} \end{aligned}

と、(12.53)式に対応する式が導かれる。(第1行はガウス分布の展開から、第2行はガンマ分布の展開から出てくる。)

次に、Eステップでの\mathbb{E} [\eta_n]\mathbb{E} [\ln\eta_n]の更新式を導く。\boldsymbol\etaの確率分布は、

\begin{aligned} 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{aligned}

となる。ガウス分布の精度\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{aligned} 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{aligned}

である。このパラメータa_n, b_nを用いると、ガンマ分布の公式(B.27)と(B.30)より、

\begin{aligned} \mathbb{E} [\eta_n] &= \frac{a_n}{b_n}\\ \mathbb{E} [\ln\eta_n] &= \psi (a_n)-\ln b_n \end{aligned}

と書ける。(\psi( \cdot )は、(B.25)で定義されるディガンマ関数)

Mステップの更新式(1/3)
上記の(*)式を\boldsymbol{\mu}で微分した値をゼロとおくと、

\begin{aligned} \frac{\partial}{\partial \boldsymbol{\mu}} \mathbb{E}_{\boldsymbol{\eta}}[\ln p(\mathbf{X}, \boldsymbol{\eta} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)] &=-\frac{1}{2} \sum_{n=1}^{N}\left\{\mathbb{E}\left[\eta_{n}\right]\left(-2 \mathbf{\Lambda} \mathbf{x}_{n}+2 \boldsymbol{\Lambda} \boldsymbol{\mu}\right)\right\} \\ &=\mathbf{\Lambda}\left(\sum_{n=1}^{N} \mathbb{E}\left[\eta_{n}\right] \mathbf{x}_{n}-\boldsymbol{\mu} \sum_{n=1}^{N} \mathbb{E}\left[\eta_{n}\right]\right)=0 \end{aligned}

\boldsymbol{\Lambda}は精度行列なので、逆行列(=分散行列)が存在する。その逆行列を左から乗じて整理すると、

\begin{aligned} \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{aligned}

Mステップの更新式(2/3)
上記の(*)式を\boldsymbol{\Lambda}で微分した値をゼロと置くのだが、準備として任意のベクトル\boldsymbol\xi,\boldsymbol\zetaに対して、

\begin{aligned} \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{aligned}

であることに注目すると、

\begin{aligned} \frac{\partial}{\partial \boldsymbol{\Lambda}} \mathbb{E}_{\boldsymbol{\eta}}[\ln p(\mathbf{X}, \boldsymbol{\eta} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)] &=-\frac{1}{2} \sum_{n=1}^{N}\left\{-\Lambda^{-1}+\mathbb{E}\left[\eta_{n}\right]\left(\mathbf{x}_{n} \mathbf{x}_{n}^{\mathrm{T}}-2 \boldsymbol{\eta} \mathbf{x}_{n}^{\mathrm{T}}+\boldsymbol{\eta} \boldsymbol{\eta}^{\mathrm{T}}\right)\right\}=0 \\ & \Leftrightarrow \Lambda_{\mathrm{ML}}=\left(\frac{1}{N} \sum_{n=1}^{N} \mathbb{E}\left[\eta_{n}\right]\left(\mathbf{x}_{n} \mathbf{x}_{n}^{\mathrm{T}}-2 \boldsymbol{\eta} \mathbf{x}_{n}^{\mathrm{T}}+\boldsymbol{\eta} \boldsymbol{\eta}^{\mathrm{T}}\right)\right)^{-1} \end{aligned}

Mステップの更新式(3/3)
上記の(*)式を\nuで微分した値をゼロと置くと、

\begin{aligned} \frac{\partial}{\partial \nu} \mathbb{E}_{\eta}[\ln p(\mathbf{X}, \boldsymbol{\eta} \mid \boldsymbol{\mu}, \boldsymbol{\Lambda}, \nu)] &=-\frac{1}{2} \sum_{n=1}^{N}\left\{\frac{2}{\Gamma\left(\frac{\nu}{2}\right)} \frac{\partial}{\partial \nu} \Gamma(\nu / 2)-(\ln \nu-\ln 2)-1-\mathbb{E}\left[\ln \eta_{n}\right]+\mathbb{E}\left[\eta_{n}\right]\right\}=0 \\ & \Leftrightarrow \frac{\psi(\nu / 2)}{\Gamma(\nu / 2)}-(\ln \nu-\ln 2)-1-\frac{1}{N} \sum_{n=1}^{N}\left(\mathbb{E}\left[\ln \eta_{n}\right]-\mathbb{E}\left[\eta_{n}\right]\right)=0 \end{aligned}

となる。これは\nuについて解析的に解くことはできない。

演習 12.25

潜在変数の空間分布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}行列は単位行列に比例したままとなり,それゆえ確率的主成分分析は通常の主成分分析がそうであるように,データ空間の軸の回転の下で共変的である.


確率的主成分分析の対数尤度関数は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}は座標系の回転に相当する

演習 12.26

\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}

で与えられる主成分による射影に影響しないことを示せ.


※ (復習)\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

k(\mathbf{x}, \mathbf{x}^{\prime}) = \mathbf{x}^{\mathrm T}\mathbf{x}^{\prime}で与えられる線形のカーネル関数を選択すれば,カーネル主成分分析の特別な場合として普通の主成分分析アルゴリズムが再現できることを示せ.


※カーネル関数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

関数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}

を使え.


固定された密度(確率密度)は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のときの図を使用)

(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

2つの変数z_1z_2が独立で,p(z_1,z_2) = p(z_1)p(z_2)となると仮定する.これらの2つの変数に対する共分散行列が対角的であることを示せ.これは変数の独立性が,2つの変数が無相関となることの十分条件であることを示している.次に,2つの変数y_1y_2を考え,y_10を中心に対称に分布しており,また,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であることが,独立性の十分条件にはならないことを示している.


※前半は演習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_10を中心に対称に分布し(\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ならば確率0y_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_1y_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_1y_2は独立ではない場合でも共分散は0となることがあることが示された。

Discussion