EMアルゴリズム(混合正規分布)

潜在変数が含まれているモデルのハイパパラメータの最尤推定量を求めるアルゴリズムをEMアルゴリズムと言う. これを混合正規分布のパラメータ推定に使ってみたくなったので調べてみた.

まず, K個の正規分布の線型結合である混合正規分布は\bf xを確率変数とすれば以下のように書ける.

p\({\bf x}\) = \sum_{k=1}^{K} \pi_k p\({\bf x} | {\bf \mu_k}, \Sigma_k\)

ここで潜在変数z=\(z_1, z_2, \cdots, z_k, \cdots, z_K\), z_k \in {0,1}を導入する. zxが第kの正規分布に属していたらz_k=1でその他の要素がすべて0となるベクトルである. つまり, xがどの正規分布に属するかを表す.
実際のサンプルはxだけなので, zは直接得られない. これは, ある混合正規分布から生起したデータxは観測できるけど, そのサンプルはどの正規分布から生起したかわかんないよねという状況.

zを導入したことで, zを観測した場合のx事後分布は

p\({\bf x} | z_k = 1\)=p\({\bf x}| {\bf \mu}_k, \Sigma_k\)
と書ける. これは, 第kの正規分布から生起したのがわかったんだから, xはその正規分布に従うでしょという状況. あたりまえだけど念のため.

そして, zの確率分布p\(z\)

p\(z_k=1\)=\pi_k
と定義する. これは, 第kの正規分布から\pi_kの確率でxが生起するという事前分布を与えている.
そうすると. 混合正規分布は以下のようにzについて周辺化したものとして書ける.

p\(\bf x\) = \sum_{\bf z} p\({\bf z}\) p\({\bf x}|{\bf z}\) \\
= \sum_{k=1}^{K} \pi_k p\({\bf x} | {\bf \mu}_k, \Sigma_k\)

ここからが本題. 観測したサンプルから混合正規分布のパラメータを推定する. 最尤推定の時の同じように対数尤度関数を考える.
観測したN個のサンプルの行列をX=\left {{\bf x}_1, \cdots, {\bf x}_N \right}^Tとし, 全て同じ混合正規分布から独立に生起したとすれば, 対数尤度関数は

\ln p\(X|{\bf \pi},{\bf \mu},\Sigma\)=\sum_{n=1}^{N} \ln\left{ \sum_{k=1}{^K}\pi_k p\({\bf x}|{\bf \mu}_k, \Sigma_k\) \right}
と書ける. ただし,

{\bf \pi}=\(\pi_1, \cdots ,\pi_K\) \\
{\bf \mu}=\({\bf \mu}_1, \cdots {\bf \mu}_K \) \\
{\bf \Sigma}=\(\Sigma_1, \cdots, \Sigma_K \)
である.

ここで, x_nを観測した場合のzの条件付き確率を\gamma\(z_{nk}\)で定義する.

\gamma\(z_{nk}\) = p\(z_k=1|{\bf x}_n\) \\
= \frac{p\(z_k=1\) p\({\bf x}_n\ | z_k=1)}{p\({\bf x}_n\)} \\
= \frac{\pi_k p\({\bf x}_n | {\bf \mu}_k, \Sigma_k\)}{\sum_{j=1}^{K} \pi_j p\({\bf x}_n | {\bf \mu}_j, \Sigma_j\)}
そして, 対数尤度関数を最大化することを考える. そのような{\bf \mu},\Sigmaでは, 各パラメータについての偏微分係数が0になってなければならない. そこで, 各パラメータで偏微分して0と置くことで

{\bf \mu}_k = \frac{1}{N_k} \sum_{n=1}^N \gamma \( z_{nk}\) {\bf x}_n\\
\Sigma_k = \frac{1}{N_k} \sum_{n=1}^N \gamma \( z_{nk}\) \({\bf x}_n -{\bf \mu}_k \)\({\bf x}_n -{\bf \mu}_k \)^T \\
を得る. ただし,
N_k = \sum_{n=1}^{N}\gamma \( z_{nk}\)
である. \gamma\(z_{nk}\)は"{\bf x}_nを観測した場合に第kの正規分布に属する確率"を表しているので, 平均と分散共分散行列は重み付き平均を取っていることになる. 例えば, \gamma\(z_{nk}\)が大きいほどn番目のサンプルが第kの正規分布に与える影響が大きくなる.
また, {\bf \pi}については, \sum_{k=1}^{K} \pi_k = 1という制約条件があるため, ラグランジュの未定乗数法を用いると
\pi_k =\frac{N_k}{N} \\
を得る.


上記の式で求めたパラメータで対数尤度関数が最大化されるはずだが, 各パラメータが式の右辺に含まれていてるので直接計算できない. そこでEMアルゴリズムを用いる. EMアルゴリズムは以下に示すEステップとMステップを繰り返すことで, 最初に設定した適当なパラメータから, 対数尤度を最大化するパラメータに収束させてゆく.

  1. {\bf \mu}_k^t,\Sigma_k^t,\pi_k^tを適当に設定
  2. Eステップ. 現在のパラメータ{\bf \mu}_k^t, \Sigma_k^t, \pi_k^t\gamma \(z_{nk}\)を計算
  3. Mステップ. パラメータの更新
     {\bf \mu}_k^{t+1} = \frac{1}{N_k} \sum_{n=1}^{N} \gamma \( z_{nk} \) {\bf x}_n \\ {\bf \Sigma}_k^{t+1} = \frac{1}{N_k} \sum_{n=1}^{N} \gamma\(z_{nk}\) \({\bf x}_n - {\bf \mu}_k^{t+1}\) ({\bf x}_n - {\bf \mu}_k^{t+1}\)^T \\ \pi_k^{t+1} = \frac{N_k}{N}
  4. 対数尤度を計算し, 前回のステップからの変化率が十分小さければ収束したとして終了する. そうでない場合は2へ

対数尤度関数には極大解が多数存在することが普通なので, EMアルゴリズムは対数尤度関数の大域的に最大の解に収束するとは限らない.
また, パラメータの初期値のはK-meanクラスタリングで得られた値を用いるのが一般的. 得られたクラスタの平均と共分散行列をそれぞれの初期値として設定し, 混合係数の初期値には, 各クラスタに所属するサンプル数を全体のサンプル数で割った値を用いる.

SVM

SOM

AdaBoost

ボルツマンマシン



カウンタ -
最終更新:2010年06月27日 11:35