EMアルゴリズムは、混合正規分布のパラメータ推定でよく登場します。ところが、難解な数式で説明されることが多く、理解が難しいです。
そこで、混合正規分布よりも簡単な例から、徐々にステップアップして、最終的には、混合正規分布のEMアルゴリズムが理解できるように解説しました。
尚、本記事は、ベイズの定理の知識を前提としています。
例１：２つの袋とボール
問題設定
袋Aと袋Bがあり、それぞれの袋には、赤いボールと青いボールが入っています。ただし、
- 袋Aは赤いボールが多い
- 袋Bは青いボールが多い
ということだけ分かっています。
袋をランダムに選んで、ボールを1個取り出す実験を何回か行います。（取り出したボールは実験のたびに、元に戻すと仮定します。）
この実験では、ボールをどちらの袋から取り出したのかは観測できません。 観測できるのは、次の結果ようなだけです。
| 実験 | 結果 |
|---|---|
| 1 | 赤 |
| 2 | 青 |
| 3 | 赤 |
| 4 | 赤 |
| 5 | 青 |
| 6 | 青 |
| 7 | 青 |
ここで知りたいのは、袋Aと袋B、それぞれに赤いボールがどのくらいの割合で入っているのか？ です。
つまり、下の$\theta_A$と$\theta_B$の値です。
\begin{align}
\theta_A&=P(\text{赤}|A)\\
\theta_B&=P(\text{赤}|B)
\end{align}
難しいポイント： どの袋から出たのか分からない
もし、どの袋から取り出したのか分かっていたら、$\theta_A$と$\theta_B$を推定するのは簡単です。例えば、
| 実験 | 結果 | 袋 |
|---|---|---|
| 1 | 赤 | A |
| 2 | 青 | B |
| 3 | 赤 | B |
| 4 | 赤 | A |
| 5 | 青 | B |
| 6 | 青 | B |
| 7 | 青 | A |
と分かっていれば、袋Aと袋Bを分けて、それぞれ赤の割合を計算できます。
\begin{align}
\theta_A&=P(\text{赤}|A)=\frac{P(\text{赤},A)}{P(A)}=\frac{2/7}{3/7}=\frac{2}{3}\\
\theta_B&=P(\text{赤}|B)=\frac{P(\text{赤},B)}{P(B)}=\frac{1/7}{4/7}=\frac{1}{4}
\end{align}
しかし、実際には、
| 実験 | 結果 | 袋 |
|---|---|---|
| 1 | 赤 | ? |
| 2 | 青 | ? |
| 3 | 赤 | ? |
| 4 | 赤 | ? |
| 5 | 青 | ? |
| 6 | 青 | ? |
| 7 | 青 | ? |
です。
つまり、「どのデータがどの袋に属しているのか分からない」 という問題になっています。
これがEMアルゴリズムが扱う典型的な問題です。
とりあえず袋の中身を仮定する
ボールがどの袋から得られたのか分からない状況では、$\theta_A$と$\theta_B$を求めるのはとても難しいです。
ですので、とりあえず、適当な値を仮定します。
- 袋Aは赤いボールが多い
- 袋Bは青いボールが多い
という情報がありますので、例えば、
袋Aから赤が出る確率を
\theta_A = P(\text{赤}|A)=0.8
袋Bから赤が出る確率を
\theta_B = P(\text{赤}|B)=0.3
とおいてしまいましょう。（何も情報がなければ0.5でもよいです）
このボールはどちらの袋？
ある実験で赤が出たとします。この赤いボールは、袋Aと袋Bのどちらから出たのでしょうか？
袋Aなら赤が出る確率は、
P(\text{赤}|A)=0.8
袋Bなら、
P(\text{赤}|B)=0.3
です。
そこで、ベイズの定理から、「観測した赤いボールが袋Aから出た確率」を下の通り求めることができます。
\begin{align}
P(A|\text{赤})
& =\frac{P(\text{赤}|A)P(A)}{P(\text{赤}|A)P(A)+P(\text{赤}|B)P(B)} \cdots(1)\\
& =\frac{P(\text{赤}|A)}{P(\text{赤}|A)+P(\text{赤}|B)} \cdots(2)\\
& =
\frac{0.8}{0.8+0.3}\\
& = \frac{8}{11}
\approx0.727
\end{align}
となります。
(1)式から(2)式への変形は、問題設定で、「袋Aを選択する確率と袋Bを選択する確率は等しい」という条件を置いていたため成り立ちます。
つまり、赤いボールが出たなら、約72.7%の確率で袋Aから出たと考えられるということです。
同様に、「観測した青いボールが袋Aから出た確率」は次のように求めることができます。
\begin{align}
P(A|\text{青})
& =\frac{P(\text{青}|A)P(A)}{P(\text{青}|A)P(A)+P(\text{青}|B)P(B)} \\
& =\frac{P(\text{青}|A)}{P(\text{青}|A)+P(\text{青}|B)} \\
& =\frac{1-P(\text{赤}|A)}{(1-P(\text{赤}|A))+(1-P(\text{赤}|B))} \\
&=\frac{0.2}{0.2+0.7}=\frac{2}{9}\approx0.222
\end{align}
つまり、青いボールが出たなら、約22.2%の確率で袋Aから出たと考えられるということです。
「袋がAかBか」を決めない
ここがEMアルゴリズムの重要なポイントです。
例えば、赤いボールについて、「これは袋Aから出た！」と決めつけるのではありません。
この決めつけは、
P(A|\text{赤})=1
と同じですが、そうではなく、先ほど求めた値である
P(A|\text{赤})=0.727
を使うのです。
では、袋の中身を更新しよう
次に、この確率を使って袋の中身を推定します。ベイズの定理より、
\begin{align}
\theta_A & = P(\text{赤}|A)\\
& =\frac{P(A|\text{赤})P(\text{赤})}{P(A|\text{赤})P(\text{赤})+P(A|\text{青})P(\text{青})} \\
& =\frac{0.727\times\frac{3}{7}}{0.727\times\frac{3}{7}+0.222\times\frac{4}{7}} \\
& \approx0.711
\end{align}
一方、Bについては、
\begin{align}
\theta_B & = P(\text{赤}|B)\\
& =\frac{P(B|\text{赤})P(\text{赤})}{P(B|\text{赤})P(\text{赤})+P(B|\text{青})P(\text{青})} \\
& =\frac{(1-0.727)\times\frac{3}{7}}{(1-0.727)\times\frac{3}{7}+(1-0.222)\times\frac{4}{7}} \\
& \approx0.208
\end{align}
です。
最初は、
\begin{align}
\theta_A&=0.8 \\
\theta_B&=0.3
\end{align}
と仮定していましたが、更新後は、
\begin{align}
\theta_A&\approx0.711\\
\theta_B&\approx0.208
\end{align}
となりました。
「このボールはどちらの袋？」をもう一度
ここで終わりではありません。
新しい値
\begin{align}
P(\text{赤}|A)&=0.711\\
P(\text{赤}|B)&=0.208
\end{align}
を使って、もう一度、このボールはAとBのどちらから出たのだろう？ を計算します。
赤の場合、
P(A|\text{赤})
=
\frac{0.711}{0.711+0.208}
\approx0.774
となります。
最初は、
P(A|\text{赤})\approx0.727
だったので、
0.727 \rightarrow 0.774
と変化しました。
同じようにして、$P(B|\text{赤})$も計算すると、初期値から変化します。
EMアルゴリズムの流れ
EMアルゴリズムでは、EステップとMステップを繰り返し行います。
繰り返しによって、最初適当に決めた値が、正しい値に収束します。
Eステップのポイント
-
袋の中身（$\theta_A, \theta_B$）を仮定して、ボールがどちらの袋から得られるかを推定する
- このとき、「赤は袋Aから得られたはずだ！」と推定するのではなく、「赤が袋Aから得られた確率（期待値）」を使う。
- いいかえると、求めたかった未知のパラメータを仮置きして、観測できなかったパラメータを期待値（Expectation）により決定する。
-
これがEステップである。式にすると下の通りである。
P(A∣\text{赤})=\frac{\theta_A}{\theta_A+\theta_B}
P(A∣\text{青})=\frac{1-\theta_A}{(1-\theta_A)+(1-\theta_B)}
Mステップのポイント
- ボールがどちらの袋から得られるかの情報を使って、袋の中身（$\theta_A, \theta_B$）を推定する。
- いいかえると、「Eステップで求めた値から、パラメータ$\theta_A, \theta_B$の尤もらしさ（もっともらしさ）を最大化（Maximization）」する。
- これがMステップである。式にすると下の通りである。
\theta_A=\frac{\displaystyle\sum_{i:\,x_i=\text{赤}}P(A\mid x_i)}{\displaystyle\sum_{i}P(A\mid x_i)}
\theta_B=\frac{\displaystyle\sum_{i:\,x_i=\text{赤}}P(B\mid x_i)}{\displaystyle\sum_{i}P(B\mid x_i)}
例２：袋の選択確率が未知
問題設定
例1では、「袋Aと袋Bを選択する確率は等しい」という条件を置いていました。つまり、
P(A)=P(B)=0.5
でした。
今回は、袋A、Bが選ばれる確率が0.5とは限らず、かつ、その値が未知であるとします。
それ以外の条件は例１と同じです。袋の選択確率も含めて、EMアルゴリズムで推定してみましょう。
袋の選択確率を表すパラメータを追加する
まず、袋Aが選ばれる確率を
\pi=P(A)
とします。
すると、袋Bが選ばれる確率は、$\pi$ を使って表現できます。
P(B)=1-\pi
したがって、今回推定したいパラメータは3つになります。
\theta_A=P(\text{赤}|A)
\theta_B=P(\text{赤}|B)
\pi=P(A)
このボールはどちらの袋？(Eステップ)
まずは例1と同じように、適当な値から始めてみましょう。
例えば、
\begin{align}
\theta_A&=0.8\\
\theta_B&=0.3\\
\pi&=0.5
\end{align}
とします。
このとき、赤いボールが観測されたときに、このボールが袋Aから出た確率はどのくらいか？ を考えます。
ベイズの定理を使うと、
\begin{align}
P(A|\text{赤})
&=
\frac{P(\text{赤}|A)P(A)}
{P(\text{赤}|A)P(A)+P(\text{赤}|B)P(B)}
\\
&=
\frac{\theta_A\pi}
{\theta_A\pi+\theta_B(1-\pi)}
\end{align}
となります。
今回の初期値では、
\theta_A=0.8
\theta_B=0.3
\pi=0.5
なので、
\begin{align}
P(A|\text{赤})
&=
\frac{0.8\times0.5}
{0.8\times0.5+0.3\times0.5}
\\
&=
\frac{0.8}{0.8+0.3}
\\
&\approx0.727
\end{align}
となります。これは例1と同じ値です。
今回はたまたま、
\pi=0.5
と仮定しているからです。
青の場合も同様にすると、
\begin{align}
P(A|\text{青})
&=
\frac{(1-\theta_A)\pi}
{(1-\theta_A)\pi+(1-\theta_B)(1-\pi)}\\
&\approx0.222
\end{align}
となります。
袋の中身を更新（Mステップ）
次に、この確率を使って袋の中身を推定します。ベイズの定理より、
\begin{align}
\theta_A & = P(\text{赤}|A)\\
& =\frac{P(A|\text{赤})P(\text{赤})}{P(A|\text{赤})P(\text{赤})+P(A|\text{青})P(\text{青})} \\
& =\frac{0.727\times\frac{3}{7}}{0.727\times\frac{3}{7}+0.222\times\frac{4}{7}} \\
& \approx0.711
\end{align}
一方、Bについては、
\begin{align}
\theta_B & = P(\text{赤}|B)\\
& =\frac{P(B|\text{赤})P(\text{赤})}{P(B|\text{赤})P(\text{赤})+P(B|\text{青})P(\text{青})} \\
& =\frac{(1-0.727)\times\frac{3}{7}}{(1-0.727)\times\frac{3}{7}+(1-0.222)\times\frac{4}{7}} \\
& \approx0.208
\end{align}
です。ここまでは例1と同じです。
例2では、それに加えて$\pi$（袋Aの選択確率）も更新します。
\begin{align}
\pi & = P(A)\\
& = \frac{N回の実験のうち、袋Aから行われた回数の期待値}{実験回数N} \\
&= \frac{\sum_{i=1}^N P(A|x_i)}{N} \\
&=\frac{0.727\times3+0.222\times4}{7}\\
&\approx0.439
\end{align}
となります。なお、全確率の公式を使って、計算することも可能です。
\begin{align}
\pi &=p(A)\\
&=P(A|赤)P(赤)+P(A|青)P(青)\\
&=0.727\times\frac{3}{7}+0.222\times\frac{4}{7}\\
&\approx0.439
\end{align}
つまり、7回の実験のうち、約43.9%が袋Aから行われたと考えるのが妥当そうだということです。
最初は、
\pi=0.5
と仮定していましたが、データから推定することで更新されました。
EMアルゴリズムの流れ
例2でもEステップとMステップを繰り返します。
これにより、適当に決めた初期値が正しい値へ収束します。
Eステップ
P(A|\text{赤})
=
\frac{
\pi\theta_A
}{
\pi\theta_A+(1-\pi)\theta_B
}
P(A|\text{青})
=
\frac{
\pi(1-\theta_A)
}{
\pi(1-\theta_A)+(1-\pi)(1-\theta_B)
}
Mステップ
\theta_A
=
\frac{
\displaystyle\sum_{i:x_i=\text{赤}}P(A|x_i)
}{
\displaystyle\sum_iP(A|x_i)
}
\theta_B
=
\frac{
\displaystyle\sum_{i:x_i=\text{赤}}P(B|x_i)
}{
\displaystyle\sum_iP(B|x_i)
}
\pi
=
\frac{1}{N}
\sum_{i=1}^{N}P(A|x_i)
例3：いよいよ混合正規分布へ
ここまでの例では、観測値は「赤」「青」の2種類しかありませんでした。
しかし、これを連続値に変えてみます。
例えば、テストの点数が、
| 人に振られたID | 点数 |
|---|---|
| 1 | 52 |
| 2 | 48 |
| 3 | 91 |
| 4 | 87 |
| 5 | 55 |
| 6 | 93 |
| ... | ... |
のように観測されたとします。
そして、このデータは、実は「2つの異なる集団」が混ざってできていると考えます。例えば、
- 集団A：平均が50点くらい
- 集団B：平均が90点くらい
のような感じです。
ところが、各点数がどちらの集団から来たのかは分かりません。 つまり、下の状況です。
| ID | 点数 | 集団A or B |
|---|---|---|
| 1 | 52 | ? |
| 2 | 48 | ? |
| 3 | 91 | ? |
| 4 | 87 | ? |
| 5 | 55 | ? |
| 6 | 93 | ? |
| ... | ... | ... |
これは、例1・例2の
「ボールは袋Aと袋Bのどちらから出たのか分からない」
という問題とまったく同じ構造です。ただし、赤・青の2種類だった観測値が、今度は連続値になっています。
例2では、袋A、袋Bから赤が出る確率と、袋Aの選択確率
\begin{align}
\theta_A&=P(\text{赤}|A)\\
\theta_B&=P(\text{赤}|B)\\
\pi&=P(A)
\end{align}
として、これを推定しました。
今回は連続値なので、正規分布を仮定し、その平均、分散を推定することにします。
つまり、
\begin{align}
\mu_A&=\text{集団Aのテストの平均点}\\
\sigma^2_A&=\text{集団Aのテストの点の分散}\\
\mu_B&=\text{集団Bのテストの平均点}\\
\sigma^2_B&=\text{集団Bのテストの点の分散}\\
\pi&=P(A)=\text{集団Aの人数の割合}\\
\end{align}
とりあえず、求めるパラメータを仮定する
例1・例2と同じように、最初から正しいパラメータを知っているわけではありません。
そこで、とりあえず適当な値を仮定します。例えば、
\begin{align}
\mu_A&=50\\
\sigma_A^2&=100\\
\mu_B&=90\\
\sigma_B^2&=100\\
\pi&=0.5
\end{align}
とします。
つまり、
- 集団Aの正規分布の平均=50、分散=100
- 集団Bの正規分布の平均=90、分散=100
- 集団AとBの人数比=1:1
と仮定したわけです。
データは集団AとBのどちらから来た？(Eステップ)
例えば、テストの点数$x=52$というデータが観測されたとします。
この52という値は、
- 正規分布Aから出たのか？
- 正規分布Bから出たのか？
分かりません。
しかし、現在のパラメータを使えば、ベイズの定理で、推定することができます。
\begin{align}
P(A|x)
&=
\frac{P(x|A)P(A)}
{P(x|A)P(A)+P(x|B)P(B)}\\
&=
\frac{\pi P(x|A)}
{\pi P(x|A)+(1-\pi)P(x|B)}
\end{align}
(連続分布の場合、$P$ではなく$f$を使う場合が多いですが、ここでは、例1，2と統一するため、$P$を用います。)
次の表は、全てのデータについて、$P(A|x)$を求めたときのイメージです。
| ID | 点数 (x) | データが集団Aのものである確率 P(A|x) |
|---|---|---|
| 1 | 52 | 0.92 |
| 2 | 48 | 0.94 |
| 3 | 91 | 0.06 |
| 4 | 87 | 0.1 |
| 5 | 55 | 0.91 |
| 6 | 93 | 0.04 |
| ... | ... | ... |
ここまでは例1とほとんど同じです。
正規分布のパラメータと集団の割合を更新する。（Mステップ）
Eステップで求めた値をもとに、パラメータを推定します。
ここからは例1と大きく異なります。なぜなら、推定するパラメータがまるで違うからです。
まず、$\mu_A, \sigma_A^2$からです。
\begin{align}
\mu_A
&=
\frac{
\sum_{i=1}^{N} x_i P(A|x_i)
}{
\sum_{i=1}^{N}P(A|x_i)
}\\
\sigma_A^2
&=
\frac{
\sum_{i=1}^{N}
P(A|x_i)(x_i-\mu_A)^2
}{
\sum_{i=1}^{N}P(A|x_i)
}\\
\end{align}
つまり、平均については、テストの点数$x$の重み付き平均になっています。重みは「Aから来た可能性が高いデータ」ほど重くします。分散も同じです。
同様に、
\begin{align}
\mu_B
&=
\frac{
\sum_{i=1}^{N}(1-P(A|x_i))x_i
}{
\sum_{i=1}^{N}(1-P(A|x_i))
}
\\
\sigma_B^2
&=
\frac{
\sum_{i=1}^{N}
(1-P(A|x_i))(x_i-\mu_B)^2
}{
\sum_{i=1}^{N}(1-P(A|x_i))
}
\end{align}
最後に$\pi$ですが、これは例2と同じです。
\pi
=
\frac{1}{N}
\sum_{i=1}^{N}P(A|x_i)
となります。
EMアルゴリズムの流れ
例3でもEステップとMステップを繰り返します。
これにより、適当に決めた初期値が正しい値へ収束します。
EMアルゴリズムの本質
ここまでの例を通して、EMアルゴリズムの本質を整理すると、次のようになります。
観測データだけでは、直接パラメータを求めるのが難しい。
そこで、
① 隠れている情報を現在のパラメータから確率的に推定する
↓
② その確率を使ってパラメータを更新する
↓
③ 更新したパラメータでもう一度、隠れている情報を推定する
↓
④ これを繰り返す
という方法です。
これが、EMアルゴリズムの基本構造です。