エンジニアを目指す浪人のブログ

情報系に役立ちそうな応用数理をゆるめにメモします

EMアルゴリズムの基礎をまとめる

数学 統計学 情報学 機械学習

機械学習でよく用いられるEMアルゴリズム(expectation-maximization algorithm ; EM algorihm)を勉強していると，その目的あるいは用途として「観測変数と(観測できない)潜在変数がある確率モデルの尤度関数を最大化するパラメータを求める」と説明されている場合を目にします．よく用いられる応用としては正しそうですが，もう少し広いクラスにも適用可能かどうかという点にモヤモヤしてしまいました．また，対数尤度の更新毎の単調性の証明も知りたいと思いました．それらの点を含めた，確率モデルの具体例には踏み込まない抽象的なレイヤーにおけるEMアルゴリズムの基礎について，調べてまとめることにしました．

文献[1]がとてもわかりやすいのでそれをベースにしてまとめますが，類似の文献[2]も参考にしています．

=================================================================================

目次

{\displaystyle \; } 0. 導入
{\displaystyle \; } 1. 完全データと観測データ
{\displaystyle \; } 2. EMアルゴリズム
{\displaystyle \; } 3. 特別な場合の {\displaystyle Q } 関数
{\displaystyle \;\;\; } 3.1. 完全データに欠損データがある場合
{\displaystyle \;\;\; } 3.2. 完全データが独立同分布に従う確率変数の場合
{\displaystyle \; } 4. 単調性と収束
{\displaystyle \; } 5. EMアルゴリズムの別解釈
{\displaystyle \; } 6. MAP推定への適用


[ 0. 導入 ]

EMアルゴリズムは確率分布のパラメータ {\displaystyle \theta } の最尤推定値を探索するための手法です．

例を考えます．24時間の窓の外の気温 {\displaystyle x \in \mathbb{R}^{24} } の発生確率は季節 {\displaystyle \theta \in \{ \mathrm{summer},\mathrm{fall},\mathrm{winter},\mathrm{spring} \} } に依存するとし，その確率分布 {\displaystyle p(x | \theta) } を既知とします．{\displaystyle x } は観測できないが1日の平均気温 {\displaystyle y = \bar{x} } は観測可能で，{\displaystyle p(y | \theta) } を最大化する {\displaystyle \theta } の最尤推定値 {\displaystyle \hat{\theta} } を推定したいものとします．この問題が求解困難な場合にEMアルゴリズムは役立ちます．観測可能なデータ {\displaystyle y = \bar{x} } が与えられた下での平均的に起こりうる {\displaystyle \log p(x | \theta) } の計算とそれを最大化する {\displaystyle \theta } の探索を反復的に行います．このようにして {\displaystyle y} が与えられた下での {\displaystyle \theta } の最尤推定値 {\displaystyle \hat{\theta} } を探索します．

EMアルゴリズムは {\displaystyle p(y | \theta) } を最大化する {\displaystyle \theta } を得ることを保証しませんが，いくつかの理論保証があります．


記号を準備します．　

{\displaystyle \;\;\; Y \in \mathbb{R}^{d_1} \;\;\; } 観測データ確率変数
{\displaystyle \;\;\; X \in \mathbb{R}^{d_2} \;\;\; } 完全データ確率変数
{\displaystyle \;\;\; Z \in \mathbb{R}^{d_3} \;\;\; } 潜在データ確率変数
{\displaystyle \;\;\; y \in \mathbb{R}^{d_1} \;\;\; } 観測データ確率変数の実現値(観測できる)
{\displaystyle \;\;\; x \in \mathbb{R}^{d_2} \;\;\; } 完全データ確率変数の実現値(観測できない)
{\displaystyle \;\;\; z \in \mathbb{R}^{d_3} \;\;\; } 潜在データ確率変数の実現値(観測できない)
{\displaystyle \;\;\; \Theta \;\;\;\;\;\;\;\;\;\;\; } パラメータ空間
{\displaystyle \;\;\; \theta \in \Theta \;\;\;\;\;\; } 確率モデルのパラメータ
{\displaystyle \;\;\; \theta^{(m)} \in \Theta \;\; } {\displaystyle m } 回目更新時の {\displaystyle \theta }

{\displaystyle \;\;\; p(y \ | \theta) \;\;\;\;\; } 観測データの確率モデル
{\displaystyle \;\;\; p(x \ | \theta) \;\;\;\;\; } 完全データの確率モデル
{\displaystyle \;\;\; \mathcal{X} \;\;\;\;\;\;\;\;\;\;\;\; X } の台 ; {\displaystyle \ p(x \ | \theta) \gt 0 } となるような {\displaystyle x } の集合の閉包
{\displaystyle \;\;\; \mathcal{X}(y) \;\;\;\;\;\;\;\; X } の台 ; {\displaystyle \ p(x \ |y, \theta) \gt 0 } となるような {\displaystyle x } の集合の閉包

{\displaystyle \;\;\; E_{X|y,\theta} \left[ \cdot \right] \;\;\;\; = \int_{\mathcal{X}(y)} \cdot \; p(x|y,\theta) dx }
{\displaystyle \;\;\; D_{ \mathrm{KL} }( \cdot || \cdot ) \;\;\;\; } カルバック・ライブラー情報量

{\displaystyle \;\;\; p(\theta) \;\;\;\;\;\;\;\;\; } {\displaystyle \theta} の事前分布


閉包の定義は過去記事にあります．

 

[ 1. 完全データと観測データ ]

完全データ {\displaystyle X } と観測データ {\displaystyle Y } の関係を決定論的な関数 {\displaystyle T } で表現します．

(1.1){\displaystyle \;\;\; T : X \to Y, \;\;\; T : x \mapsto T(x) }

例えば，集合 {\displaystyle X } をその要素の平均に移す，ベクトル {\displaystyle X } をその {\displaystyle l_1 } ノルムに写す，などの使い方が考えられますが，最も有名なのは，完全データ {\displaystyle X } が観測データ {\displaystyle Y } と欠損データ(あるいは潜在データともいいます) {\displaystyle Z } の結合となる，すなわち，

(1.2.1){\displaystyle \;\;\; X= (Y,Z) }

(1.2.2){\displaystyle \;\;\; Y = T(X) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = T( \ (Y,Z) \ ) }

となる場合です．{\displaystyle T } は {\displaystyle X } から単に {\displaystyle Z } を取り除く関数です．混合正規分布(Gaussian mixture model ; GMM)や隠れマルコフモデル(hidden Markov model ; HMM)を表現する場合に相当します．


後の章で示すEMアルゴリズムの単調性が成り立たなくなるので， {\displaystyle T } が確率的な関数ではいけません．ただし，{\displaystyle N } を確率変動項として {\displaystyle Y=X+N } となるような場合を扱えないということではなく，その方法は {\displaystyle \tilde{X}=(X,N), \ Y=T(\tilde{X}) } とすることです．

 

[ 2. EMアルゴリズム ]

最尤推定とは，対数尤度関数を最大化する {\displaystyle \theta } の最尤推定値(maximum likelihood estimate ; MLE)を求めることです．

(2.1){\displaystyle \;\;\; L(\theta) = \log p(y \ | \theta) }

(2.2){\displaystyle \;\;\; \hat{\theta}_{\rm{MLE}} = \mathrm{arg} \max_{\theta \in \Theta} \ \log p(y \ | \theta) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ L(\theta) }


問題(2.2)を解くことが困難な場合に，EMアルゴリズムは有用かもしれません． {\displaystyle L(\theta) } の代わりに {\displaystyle \log p(x \ | \theta) } を用います．ただし{\displaystyle x } は観測できないので {\displaystyle \log p(x \ | \theta) } を直接的に最大化できません．条件付き確率分布 {\displaystyle p(x \ |y, \theta^{(m)}) } は計算できるので，その分布についての期待値を最大化するという思想です．そのために以下の {\displaystyle Q } 関数を導入します．条件付き確率の定義(文献[3]にあります)を用います．


(2.3){\displaystyle \;\;\; Q(\theta \ |\theta^{(m)}) = E_{X|y,\theta^{(m)}} \left[ \log p(X \ |\theta) \right] }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(x \ | \theta) \ p(x \ |y, \theta^{(m)}) \ dx }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(x \ | \theta) \ \frac{p(x ,y \ | \theta^{(m)}) }{p(y \ | \theta^{(m)}) } \ dx \;\;\;\;\;\;\;\;\;\;\; \because } 条件付き確率の定義
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(x \ | \theta) \ \frac{p(x , T(x) \ | \theta^{(m)}) }{p(y \ | \theta^{(m)}) } \ dx \;\;\;\;\;\; \because } (1.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(x \ | \theta) \ \frac{p(x \ | \theta^{(m)}) }{p(y \ | \theta^{(m)}) } \ dx }

 

EMアルゴリズムは以下の6つのステップで構成されています．
---------------------------------------------------------------------------------------------------------------------
Step 1: {\displaystyle \;\; }初期値 {\displaystyle \theta^{(m=0)} } を設定する．

Step 2: {\displaystyle \;\; }{\displaystyle \theta^{(m)} } と観測データ {\displaystyle y } を用いて {\displaystyle \ p(x \ |y, \theta^{(m)}) } を計算する．

Step 3: {\displaystyle \;\; }{\displaystyle \theta^{(m)} } を捨てる．

Step 4: {\displaystyle \;\; }Step 2で得た {\displaystyle p(x \ |y, \theta^{(m)}) } を用いて {\displaystyle Q(\theta \ |\theta^{(m)}) } を準備する．

Step 5: {\displaystyle \;\; }{\displaystyle \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ Q(\theta \ |\theta^{(m)}) } を探索する．

Step 6: {\displaystyle \;\; }{\displaystyle m=m+1 } としてStep 2に戻る．
---------------------------------------------------------------------------------------------------------------------

Step 6について補足します．EMアルゴリズム収束基準は明確ではありません．標準的なものとして以下の2つがあります．{\displaystyle \left\| \cdot \right\| } は適当なノルムとします．
{\displaystyle \;\;\; } ・ある {\displaystyle \epsilon \gt 0 } について {\displaystyle \left\| \theta^{(m+1)}- \theta^{(m)} \right\| \lt \epsilon } となるまで反復する．
{\displaystyle \;\;\; } ・ある {\displaystyle \epsilon \gt 0 } について {\displaystyle | L(\theta^{(m+1)})- L(\theta^{(m)})| \lt \epsilon } となるまで反復する．


伝統的な記述に従うと以下となります．
---------------------------------------------------------------------------------------------------------------------
E-Step: {\displaystyle \;\; p(x \ |y, \theta^{(m)}) } を計算して {\displaystyle Q(\theta \ |\theta^{(m)}) } を準備する．

M-Step: {\displaystyle \;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ Q(\theta \ |\theta^{(m)}) } を計算する．
---------------------------------------------------------------------------------------------------------------------

 

[ 3. 特別な場合の {\displaystyle Q } 関数 ]

特別な場合の {\displaystyle Q } 関数について論じます．本章の結果はEMアルゴリズムを混合正規分布や隠れマルコフモデルに適用する場合に用います．


[ 3.1. 完全データに欠損データがある場合 ]

完全データ {\displaystyle X } が(1.2.1)，{\displaystyle T } が(1.2.2)となるとき，{\displaystyle Q } 関数は以下のように完全データ {\displaystyle X } を消去して欠損データ {\displaystyle Z } についての条件付き期待値に変形できます．

(3.1){\displaystyle \;\;\; Q(\theta \ |\theta^{(m)}) = E_{X|y,\theta^{(m)}} \left[ \log p(X \ |\theta) \right] \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } (2.3)

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(x \ | \theta) \ p(x \ |y, \theta^{(m)}) \ dx \;\;\;\;\;\; \because } (2.3)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(y,z \ | \theta) \ p(y,z \ |y, \theta^{(m)}) \ dx \;\;\;\;\;\; \because } (1.2.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} \log p(y,z \ | \theta) \ p(z \ |y, \theta^{(m)}) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{Z}(y)} \log p(y,z \ | \theta) \ p(z \ |y, \theta^{(m)}) \ dz }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = E_{Z|y,\theta^{(m)}} \left[ \log p(y,Z \ |\theta) \right] }


[ 3.2. 完全データが独立同分布に従う確率変数の場合 ]

完全データ {\displaystyle X } が独立同分布({\displaystyle \mathrm{i.i.d.} })に従う {\displaystyle n } 個の確率変数 {\displaystyle X_i, \ i=1,\ldots,n } からなる場合を考えます．実現値 {\displaystyle x } は {\displaystyle n } 個の {\displaystyle \mathcal{X} } からなる直積集合 {\displaystyle \mathcal{X}^n = \mathcal{X} \times \cdots \times \mathcal{X} } の要素，すなわち {\displaystyle x \in \mathcal{X}^n , \; x_i \in \mathcal{X} } となります．確率モデルは {\displaystyle \forall x \in \mathcal{X}^n, \; \forall \theta \in \Theta } について以下となります．

(3.2.1){\displaystyle \;\;\; p(x \ | \theta) = \prod_{i=1}^n p(x_i \ | \theta) }

完全データと観測データの関係を以下とします．

(3.2.2){\displaystyle \;\;\; y_i = T(x_i), \; i=1,\ldots,n }


命題 6.1.
以下が成り立つ．

(3.2.3){\displaystyle \;\;\; Q(\theta \ |\theta^{(m)}) = \sum_{i=1}^n Q_{i}(\theta \ |\theta^{(m)}) }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \sum_{i=1}^n E_{X_i | y_i,\theta^{(m)}} \left[ \log p(X_i \ | \theta) \right] }

証明.

初めに以下を示す．

(3.2.4){\displaystyle \;\;\; p(x,y \ | \theta) = p(x,T(x) \ | \theta) \;\;\;\;\;\; \because } (1.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = p(x \ | \theta) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; =\prod_{i=1}^n p(x_i \ | \theta) \;\;\;\;\;\;\;\;\; \because } (3.2.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; =\prod_{i=1}^n p(x_i,T(x_i) \ | \theta) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; =\prod_{i=1}^n p(x_i,y_i \ | \theta) \;\;\;\;\;\; \because } (3.2.2)


これを用いて次に以下を示す．

(3.2.5){\displaystyle \;\;\; p(x_i \ | y, \theta) = \frac{ p(x_i,y \ | \theta) }{ p(y \ | \theta) } \;\;\;\;\;\;\;\; \because } 条件付き確率の定義

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \frac{ \int_{\mathcal{X}^{n-1} } p(x,y \ | \theta) \ dx_1 \ldots dx_{i-1}dx_{i+1}\ldots dx_n }{ \int_{\mathcal{X}^n} p(x,y \ | \theta) \ dx } }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \frac{ \int_{\mathcal{X}^{n-1} } \prod_{j=1}^n p(x_j,y_j \ | \theta) \ dx_1 \ldots dx_{i-1}dx_{i+1}\ldots dx_n }{ \int_{\mathcal{X}^n} \prod_{j=1}^n p(x_j,y_j \ | \theta) \ dx } \;\;\; \because } (3.2.4)

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \frac{ p(x_i,y_i \ | \theta) \prod_{j=1, \ j \neq i }^n \int_{\mathcal{X} } p(x_j,y_j \ | \theta) \ dx_j }{ \prod_{j=1}^n \int_{\mathcal{X}} p(x_j, y_j \ | \theta) \ dx_j } \;\;\; \because } {\displaystyle \; \mathrm{i.i.d.} }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \frac{ p(x_i,y_i \ | \theta) \prod_{j=1, \ j \neq i }^n p(y_j \ | \theta) }{ \prod_{j=1}^n p(y_j \ | \theta)} }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \frac{ p(x_i,y_i \ | \theta) }{ p(y_i \ | \theta)} }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = p(x_i \ | y_i, \theta) \;\;\;\;\;\; \because } 条件付き確率の定義


これを用いて最終的に以下を得る．

(3.2.6){\displaystyle \;\;\; Q(\theta \ |\theta^{(m)}) = E_{X|y,\theta^{(m)}} \left[ \log p(X \ |\theta) \right] \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } (2.3)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = E_{X|y,\theta^{(m)}} \left[ \log \prod_{i=1}^n p(X_i \ | \theta) \right] \;\;\;\;\;\;\;\;\;\; \because } (3.2.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = E_{X|y,\theta^{(m)}} \left[ \sum_{i=1}^n \log p(X_i \ | \theta) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \sum_{i=1}^n E_{X_i |y,\theta^{(m)}} \left[ \log p(X_i \ | \theta) \right] \;\;\;\;\;\;\;\;\;\;\;\; \because } {\displaystyle \; \mathrm{i.i.d.} }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \sum_{i=1}^n E_{X_i |y_i,\theta^{(m)}} \left[ \log p(X_i \ | \theta) \right] \;\;\;\;\;\;\;\;\;\;\; \because } (3.2.5)

(証明終わり)

 

[ 4. 単調性と収束 ]

はじめに以下の定理を示します．証明にはイェンゼンの不等式(Jensen's inequality)(過去記事にあります)と条件付き確率の定義を用います．途中の(※)では第1項のみが {\displaystyle \theta } に依存していることに注意します．

定理 4.1.
{\displaystyle \;\;\; \forall \theta \in \Theta, \;\; \left[ \ Q(\theta \ | \theta^{(m)}) \ge Q(\theta^{(m)} | \theta^{(m)}) \Longrightarrow L(\theta) \ge L(\theta^{(m)}) \ \right] }

証明.

{\displaystyle \;\;\; L(\theta)= \log p(y \ | \theta) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log \int_{\mathcal{X}(y)} p(x ,y \ | \theta) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log \int_{\mathcal{X}(y)} p(x ,T(x) \ | \theta) \ dx \;\;\;\;\;\;\;\;\; \because } (1.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log \int_{\mathcal{X}(y)} p(x \ | \theta) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log \int_{\mathcal{X}(y)} \frac{p(x \ | \theta)}{p(x \ | y, \theta^{(m)})} p(x \ | y, \theta^{(m)}) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log E_{X|y,\theta^{(m)}} \left[\frac{p(X \ | \theta)}{p(X \ | y, \theta^{(m)})} \right] }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\; \ge E_{X|y,\theta^{(m)}} \left[ \log \frac{p(X \ | \theta)}{p(X \ | y, \theta^{(m)})} \right] \;\;\;\;\; \because } イェンゼンの不等式
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = E_{X|y,\theta^{(m)}} \left[ \log p(X \ | \theta) \right] + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = Q(\theta \ |\theta^{(m)}) + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] \;\;\;\;\;\;\;\; \because } (2.3) {\displaystyle \;\;\;\;\;\; } (※)

{\displaystyle \;\;\;\;\;\;\;\;\;\;\; \ge Q(\theta^{(m)} \ |\theta^{(m)}) + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] \;\;\;\;\; \because } 定理の仮定
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = E_{X|y,\theta^{(m)}} \left[ \log p(X \ | \theta^{(m)}) \right] + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = E_{X|y,\theta^{(m)}} \left[ \log \frac{p(X \ | \theta^{(m)})}{p(X \ | y, \theta^{(m)})} \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} p(x \ | y, \theta^{(m)}) \log \frac{p(x \ | \theta^{(m)})}{p(x \ | y, \theta^{(m)})} \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} p(x \ | y, \theta^{(m)}) \log \frac{p(x,T(x) \ | \theta^{(m)})}{p(x\ | y, \theta^{(m)})} \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} p(x \ | y, \theta^{(m)}) \log \frac{p(x,y \ | \theta^{(m)})}{p(x \ | y, \theta^{(m)})} \ dx \;\;\;\;\; \because } (1.1)
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \int_{\mathcal{X}(y)} p(x \ | y, \theta^{(m)}) \log p(y \ | \theta^{(m)}) \ dx \;\;\;\;\;\;\;\;\;\; \because } 条件付き確率の定義
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log p(y \ | \theta^{(m)}) \int_{\mathcal{X}(y)} p(x \ | y, \theta^{(m)}) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = \log p(y \ | \theta^{(m)}) }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\; = L(\theta^{(m)}) \;\;\;\;\;\; } (※※)

(証明終わり)


EMアルゴリズムの単調性(monotonicity)は上の定理を用いて容易に示すことができます．{\displaystyle L(\theta^{(m)}) } は更新毎に単調増加し少なくとも最悪(最小)にはならないことを意味します．

{\displaystyle \;\;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ Q(\theta \ |\theta^{(m)}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } Step 5

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\; \Rightarrow \;\;\; Q(\theta^{(m+1)} \ | \theta^{(m)}) \ge Q(\theta^{(m)} | \theta^{(m)}) }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\; \Rightarrow \;\;\; L(\theta^{(m+1)}) \ge L(\theta^{(m)}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } 定理 4.1.


上に有界で単調増加な数列は収束する(文献[4]にあります)ので，{\displaystyle L(\theta) } が {\displaystyle \forall \theta \in \Theta } で上に有界であれば列 {\displaystyle \{ L(\theta^{(m)}) \}_{m=0,1,\ldots} } の収束は保証されます(このとき {\displaystyle \theta } の有界性については何もいっていないので列 {\displaystyle \{ \theta^{(m)} \}_{m=0,1,\ldots} } の収束は保証しません)．


EMアルゴリズムにおける {\displaystyle \{ \theta^{(m)} \}_{m=0,1,\ldots} } の一般的な収束定理は存在しません．{\displaystyle L(\theta) } と {\displaystyle Q(\theta \ | \theta') } と初期値 {\displaystyle \theta^{(m=0)} } に依存します．

 

[ 5. EMアルゴリズムの別解釈 ]

{\displaystyle \mathcal{X}(y) } を台にもち密度が {\displaystyle \tilde{p}(x) } となる完全データの確率分布 {\displaystyle \tilde{P} } と {\displaystyle \theta } の関数 {\displaystyle F } を以下のように定義します．密度 {\displaystyle p(x \ |y, \theta) } の確率分布を {\displaystyle P_{\theta} } とします．

(5.1){\displaystyle \;\;\; F( \tilde{P},\theta ) = L(\theta)- D_{ \mathrm{KL} }( \tilde{P} || P_{\theta} ) }


EMアルゴリズムは以下のように2つの最大化問題の反復として表現することもできます．カルバック・ライブラー情報量を最小化する確率分布 {\displaystyle \tilde{P} } を探索し，その分布を用いて関数 {\displaystyle F } を最大化する {\displaystyle \theta } を探索します．
---------------------------------------------------------------------------------------------------------------------
Max Step 1: {\displaystyle \;\; \tilde{P}^{(m+1)} = \mathrm{arg} \max_{ \tilde{P} } \ F( \tilde{P},\theta^{(m)} ) } を計算する．

Max Step 2: {\displaystyle \;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ F( \tilde{P}^{(m+1)} ,\theta ) } を計算する．
---------------------------------------------------------------------------------------------------------------------


Max Step 1は2章のE-Stepと(ほぼ)同じであることを示します．

{\displaystyle \;\;\; \tilde{P}^{(m+1)} = \mathrm{arg} \max_{ \tilde{P} } \ F( \tilde{P},\theta^{(m)} ) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{ \tilde{P} } \ \left[ L ( \theta^{(m)} )- D_{ \mathrm{KL} } ( \tilde{P} || P_{ \theta^{(m)} } ) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \min_{ \tilde{P} } \ D_{ \mathrm{KL} } ( \tilde{P} || P_{ \theta^{(m)} } ) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = P_{ \theta^{(m)} } }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = p(x \ |y, \theta^{(m)}) }


Max Step 2は2章のM-Stepと同じであることを示します．

{\displaystyle \;\;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ F( \tilde{P}^{(m+1)} ,\theta ) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ L ( \theta )- D_{ \mathrm{KL} } ( \tilde{P}^{(m+1)} || P_{ \theta } ) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ \log p(y \ | \theta) - D_{ \mathrm{KL} } ( \tilde{P}^{(m+1)} || P_{ \theta } ) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log p(y \ | \theta) \ dx - D_{ \mathrm{KL} } ( \tilde{P}^{(m+1)} || P_{ \theta } ) \right] }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log \frac{p(x \ | \theta)}{p(x \ |y, \theta)} \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; - \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log \frac{ p(x \ |y, \theta^{(m)}) }{ p(x \ |y, \theta) }\ dx \;\;\; \because } Max Step 1
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log p(x \ | \theta) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; - \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log p(x \ |y, \theta^{(m)}) \ dx }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \int_{\mathcal{X}(y)} p(x \ |y, \theta^{(m)}) \log p(x \ | \theta) \ dx }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ E_{X|y,\theta^{(m)}} \left[ \log p(X \ |\theta) \right] }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ Q(\theta \ |\theta^{(m)}) \;\;\; \because } (2.3)


[ 6. MAP推定への適用 ]

EM アルゴリズムは最尤推定を行うアルゴリズムですが，MAP推定にも適用可能です．MAP推定とは，事後分布を最大化する {\displaystyle \theta } のMAP推定値(maximum a posteriori estimate)を求めることです．ベイズの定理(文献[5]にあります)を用います．

(6.1){\displaystyle \;\;\; \hat{\theta}_{\rm{MAP}} = \mathrm{arg} \max_{\theta \in \Theta} \ \log p(\theta \ |y ) }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \log \frac{ p(y \ |\theta ) p(\theta)}{p(y)} \;\;\;\;\;\;\;\;\;\;\;\; \because } ベイズの定理
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ \log p(y \ |\theta ) + \log p(\theta) \right] }
{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ L(\theta) + \log p(\theta) \right] \;\;\;\;\;\;\;\;\;\; \because } (2.1)


MAP推定のためのEMアルゴリズムは以下となります．MAP E-Stepは2章のE-Stepと同じです．
---------------------------------------------------------------------------------------------------------------------
MAP E-Step: {\displaystyle \;\; p(x \ |y, \theta^{(m)}) } を計算して {\displaystyle Q(\theta \ |\theta^{(m)}) } を準備する．

MAP M-Step: {\displaystyle \;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ Q(\theta \ |\theta^{(m)}) + \log p(\theta) \right] } を計算する．
---------------------------------------------------------------------------------------------------------------------


以下の定理を示します．証明は定理 4.1.の証明とほぼ同じです．

定理 7.1.

{\displaystyle \forall \theta \in \Theta } で以下が成り立つ．

{\displaystyle \;\;\; Q(\theta \ | \theta^{(m)}) + \log p(\theta) \ge Q(\theta^{(m)} | \theta^{(m)}) + \log p(\theta^{(m)}) }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \Longrightarrow L(\theta) + \log p(\theta) \ge L(\theta^{(m)}) + \log p(\theta^{(m)}) }

証明.

{\displaystyle \;\;\; L(\theta) + \log p(\theta) }

{\displaystyle \;\;\;\;\;\;\;\;\; \ge Q(\theta \ |\theta^{(m)}) + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] + \log p(\theta) \;\;\;\;\;\;\;\;\;\;\; \because } (※) に {\displaystyle + \log p(\theta) }
{\displaystyle \;\;\;\;\;\;\;\;\; \ge Q(\theta^{(m)} \ |\theta^{(m)}) + E_{X|y,\theta^{(m)}} \left[ - \log p(X \ | y, \theta^{(m)}) \right] + \log p(\theta^{(m)}) \;\;\; \because } 定理の仮定
{\displaystyle \;\;\;\;\;\;\;\;\; \ge L(\theta^{(m)}) + \log p(\theta^{(m)}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } (※※)

(証明終わり)


MAP推定のためのEMアルゴリズムの単調性は上の定理を用いて容易に示すことができます．{\displaystyle L(\theta^{(m)}) + \log p(\theta^{(m)}) } は更新毎に単調増加し少なくとも最悪(最小)にはならないことを意味します．

{\displaystyle \;\;\; \theta^{(m+1)} = \mathrm{arg} \max_{\theta \in \Theta} \ \left[ Q(\theta \ |\theta^{(m)}) + \log p(\theta) \right] \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } MAP M-Step

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\; \Rightarrow \;\;\; Q(\theta^{(m+1)} \ | \theta^{(m)}) + \log p(\theta^{(m+1)}) \ge Q(\theta^{(m)} | \theta^{(m)}) + \log p(\theta^{(m)}) }

{\displaystyle \;\;\;\;\;\;\;\;\;\;\;\;\; \Rightarrow \;\;\; L(\theta^{(m+1)}) + \log p(\theta^{(m+1)}) \ge L(\theta^{(m)}) + \log p(\theta^{(m)}) \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \because } 定理 7.1.

=================================================================================

以上，EMアルゴリズムの基礎をまとめました．

 


参考文献
[1] University of Washington Yihua Chen先生らによるノート https://vannevar.ece.uw.edu/techsite/papers/documents/UWEETR-2010-0002.pdf
[2] Google Research Maya R. Gupta先生らによるノート http://mayagupta.org/publications/EMbookGuptaChen2010.pdf
[3] Wikipedia Conditional probability のページ https://en.wikipedia.org/wiki/Conditional_probability
[4] Wikipedia Monotone convergence theorem のページ https://en.wikipedia.org/wiki/Monotone_convergence_theorem
[5] Wikipedia Bayes' theorem のページ https://en.wikipedia.org/wiki/Bayes%27_theorem