EMアルゴリズム


※上記の広告は60日以上更新のないWIKIに表示されています。更新することで広告が下部へ移動します。

EMアルゴリズム(Expectation-Maximization Algorithm)

確率モデルが観測できない変数(潜在変数/隠れ変数)に依存する場合に最尤法を実施するためのアルゴリズム。非常に多くの応用に使え、変分ベイズ法などの基礎を成すのでとっても大事なアルゴリズム。

推定の流れ

このアルゴリズムでは、現在のパラメータ(\theta)と観測変数(X)から得られる情報を使って、潜在変数(Z)の条件付き確率(p(Z|X,\theta))を求めるステップ(E-Step)と、観測変数の事後確率の期待値を最大化するフェーズ(M-Step)の二段階の処理を繰り返しながらパラメータを最適化する。
つまり、尤度関数(P(X,Z|\theta))の最大化を一発で計算したいところなのだが、それは難しいので、今の\thetaを使って計算されるP(Z|X,\theta)を利用して、尤度関数の期待値(\sum_{Z}P(Z|X,\theta)P(X,Z|\theta))の最大化という問題に置き換えている。実際には、最後の式を2つ目の\thetaを構成する各変数について偏微分して0になる値を探すなどを行って更新することになる。

詳細な定式化(間違ってるかも)

パラメータを\theta、観測変数をX、潜在変数をZとする。
目的は、観測変数のパラメータに対する確率p(X|\theta)の最大化。ところが、潜在変数があるのでp(X,Z|\theta)の周辺化で置き換えて話を進められるようにする。そこでベイズの公式による式変形を行う。

p(X|\theta)=\frac{p(X,Z|\theta)}{p(Z|X,\theta)}

さて、ここで天下りながら分布q(Z)を導入する。これは潜在変数Zの事後分布の近似分布である。式変形のミソは、この分布q(Z)p(Z|X,\theta)の間の違い(距離)が現れるように前式を変形していくことである。そこで、まず右辺の分子分母をq(Z)で割る。

p(X|\theta)=\frac{\frac{p(X,Z|\theta)}{q(Z)}}{\frac{p(Z|X,\theta)}{q(Z)}}

次に両辺に対数をとる。

\ln p(X|\theta)=\ln \frac{p(X,Z|\theta)}{q(Z)} - \ln \frac{p(Z|X,\theta)}{q(Z)}

さらに両辺にq(Z)をかけてZに関し周辺化を行う。

\sum_{Z}q(Z)\ln P(X|\theta) = \sum_{Z}q(Z)\ln \frac{p(X,Z|\theta)}{q(Z)} - \sum_{Z}q(Z)\ln \frac{p(Z|X,\theta)}{q(Z)}

左辺において、\ln P(X|\theta)Zに関係せず、\sum_{Z}q(Z)は1となる。
右辺の第一項をL(q, \theta)、第二項をKL(q||p)とかくと、下記式を得ることができる。

\ln P(X|\theta) = L(q, \theta) + KL(q||p)

この式において最適化で変更できるのはq\thetaである。EMアルゴリズムは、この2つの視点で交互に最大化する。
まずEステップでは、KL(q||p)は0以上なので、KL(q||p)を最小化することでL(q, \theta) を最大化する。KL距離の定義よりKL(q||p)q=pのとき最小化となる。すなわち、\thetaを固定しqp(Z|X,\theta)を設定すればよい。

q(Z)=p(Z|X,\theta)

次いでMステップでは、qを固定し\thetaに関して最大化する。すなわち、下記式を\thetaに関して最大化する。ここで、Eステップ時点で利用した\thetaの値は、\theta^{old}として固定されている。

\sum_{Z}p(Z|X,\theta^{old})\ln \frac{p(X,Z|\theta)}{p(Z|X,\theta^{old})}
 = \sum_{Z}p(Z|X,\theta^{old})\ln(X,Z|\theta) - \sum_{Z}p(Z|X,\theta^{old})\ln(Z|X, \theta^{old})
 = Q(\theta, \theta^{old}) + const.

上述において、L(q, \theta) \ln p(X|\theta)の下界をなしており、Eステップではこの下界を最大化するqを手に入れている。一方、Mステップでは、L(q, \theta) \thetaに関して最大化されるわけだが、このとき新しい\thetaを使ったp(Z|X,\theta)と今までのq(Z)との間に新たな違い(距離)が生まれる。そのため、EステップとMステップを繰り返すことで徐々に\ln p(X|\theta)を最大化していく必要がある。

実装例

言語

Python 2.6 + scipy + matplotlib

問題設定

2つの正規分布からなる2次元の混合正規分布に対して平均、分散、負担率を推定する。
詳細はソース参照。

ソース


結果

  • コンソール最後の部分
<< 29 >>
[ 0 ]
Mu=
[[ 21.05175707]
 [ 21.65353366]]
Sigma=
[[  95.11897414    9.67571604]
 [   9.67571604  116.52378045]]
Pi=
0.674809728052
[ 1 ]
Mu=
[[-14.92553657]
 [-10.79587385]]
Sigma=
[[ 137.40967457   -2.5813963 ]
 [  -2.5813963    95.53296063]]
Pi=
0.325190271948
 

  • グラフ

補足

この例ではうまくいくが、もとの正規分布の分散が大きすぎると一つの正規分布でほとんどのデータを説明し、もう一つの正規分布がごくわずか(1点)のデータを説明しようとしてしまう。このあたりは、分散の監視などが必要。
また、そもそもの収束判定には尤度関数の値の変化率を見るのが一般的。
|新しいページ|検索|ページ一覧|RSS|@ウィキご利用ガイド | 管理者にお問合せ
|ログイン|
添付ファイル


統計・機械学習


計算機


プログラミング



Raspberry Pi


その他



工事中


関係ない話題


管理者用



プロフィール

但馬慶行
(TAJIMA Yoshiyuki)

機械学習とロボットを趣味で勉強中。
連絡先はトップページ