大白話講解EM算法

001、一個非常簡單的例子

假設現在有兩枚硬币1和2，,随機抛擲後正面朝上機率分别為P1，P2。為了估計這兩個機率，做實驗，每次取一枚硬币，連擲5下，記錄下結果，如下：

硬币	結果	統計
1	正正反正反	3正-2反
2	反反正正反	2正-3反
1	正反反反反	1正-4反
2	正反反正正	3正-2反
1	反正正反反	2正-3反

可以很容易地估計出P1和P2，如下：

P1 = （3+1+2）/ 15 = 0.4

P2= （2+3）/10 = 0.5

到這裡，一切似乎很美好，下面我們加大難度。

010、加入隐變量z

還是上面的問題，現在我們抹去每輪投擲時使用的硬币标記，如下：

硬币	結果	統計
Unknown	正正反正反	3正-2反
Unknown	反反正正反	2正-3反
Unknown	正反反反反	1正-4反
Unknown	正反反正正	3正-2反
Unknown	反正正反反	2正-3反

好了，現在我們的目标沒變，還是估計P1和P2，要怎麼做呢？

顯然，此時我們多了一個隐變量z，可以把它認為是一個5維的向量（z1,z2,z3,z4,z5)，代表每次投擲時所使用的硬币，比如z1，就代表第一輪投擲時使用的硬币是1還是2。但是，這個變量z不知道，就無法去估計P1和P2，是以，我們必須先估計出z，然後才能進一步估計P1和P2。

但要估計z，我們又得知道P1和P2，這樣我們才能用最大似然機率法則去估計z，這不是雞生蛋和蛋生雞的問題嗎，如何破？

答案就是先随機初始化一個P1和P2，用它來估計z，然後基于z，還是按照最大似然機率法則去估計新的P1和P2，如果新的P1和P2和我們初始化的P1和P2一樣，請問這說明了什麼？（此處思考1分鐘）

這說明我們初始化的P1和P2是一個相當靠譜的估計！

就是說，我們初始化的P1和P2，按照最大似然機率就可以估計出z，然後基于z，按照最大似然機率可以反過來估計出P1和P2，當與我們初始化的P1和P2一樣時，說明是P1和P2很有可能就是真實的值。這裡面包含了兩個互動的最大似然估計。

如果新估計出來的P1和P2和我們初始化的值差别很大，怎麼辦呢？就是繼續用新的P1和P2疊代，直至收斂。

這就是下面的EM初級版。

011、EM初級版

我們不妨這樣，先随便給P1和P2賦一個值，比如：

P1 = 0.2

P2 = 0.7

然後，我們看看第一輪抛擲最可能是哪個硬币。

如果是硬币1，得出3正2反的機率為 0.2*0.2*0.2*0.8*0.8 = 0.00512

如果是硬币2，得出3正2反的機率為0.7*0.7*0.7*0.3*0.3=0.03087

然後依次求出其他4輪中的相應機率。做成表格如下：

輪數	若是硬币1	若是硬币2
1	0.00512	0.03087
2	0.02048	0.01323
3	0.08192	0.00567
4	0.00512	0.03087
5	0.02048	0.01323

按照最大似然法則：

第1輪中最有可能的是硬币2

第2輪中最有可能的是硬币1

第3輪中最有可能的是硬币1

第4輪中最有可能的是硬币2

第5輪中最有可能的是硬币1

我們就把上面的值作為z的估計值。然後按照最大似然機率法則來估計新的P1和P2。

P1 = （2+1+2）/15 = 0.33

P2=（3+3）/10 = 0.6

設想我們是全知的神，知道每輪抛擲時的硬币就是如本文第001部分标示的那樣，那麼，P1和P2的最大似然估計就是0.4和0.5（下文中将這兩個值稱為P1和P2的真實值）。那麼對比下我們初始化的P1和P2和新估計出的P1和P2：

初始化的P1	估計出的P1	真實的P1	初始化的P2	估計出的P2	真實的P2
0.2	0.33	0.4	0.7	0.6	0.5

看到沒？我們估計的P1和P2相比于它們的初始值，更接近它們的真實值了！

可以期待，我們繼續按照上面的思路，用估計出的P1和P2再來估計z，再用z來估計新的P1和P2，反複疊代下去，就可以最終得到P1 = 0.4，P2=0.5，此時無論怎樣疊代，P1和P2的值都會保持0.4和0.5不變，于是乎，我們就找到了P1和P2的最大似然估計。

這裡有兩個問題：

1、新估計出的P1和P2一定會更接近真實的P1和P2？

答案是：沒錯，一定會更接近真實的P1和P2，數學可以證明，但這超出了本文的主題，請參閱其他書籍或文章。

2、疊代一定會收斂到真實的P1和P2嗎？

答案是：不一定，取決于P1和P2的初始化值，上面我們之是以能收斂到P1和P2，是因為我們幸運地找到了好的初始化值。

100、EM進階版

下面，我們思考下，上面的方法還有沒有改進的餘地？

我們是用最大似然機率法則估計出的z值，然後再用z值按照最大似然機率法則估計新的P1和P2。也就是說，我們使用了一個最可能的z值，而不是所有可能的z值。

如果考慮所有可能的z值，對每一個z值都估計出一個新的P1和P2，将每一個z值機率大小作為權重，将所有新的P1和P2分别權重相加，這樣的P1和P2應該會更好一些。

所有的z值有多少個呢？顯然，有2^5=32種，需要我們進行32次估值？？

不需要，我們可以用期望來簡化運算。

輪數	若是硬币1	若是硬币2
1	0.00512	0.03087
2	0.02048	0.01323
3	0.08192	0.00567
4	0.00512	0.03087
5	0.02048	0.01323

利用上面這個表，我們可以算出每輪抛擲中使用硬币1或者使用硬币2的機率。比如第1輪，使用硬币1的機率是：

0.00512/(0.00512+0.03087)=0.14

使用硬币2的機率是1-0.14=0.86

依次可以算出其他4輪的機率，如下：

輪數	z_i=硬币1	z_i=硬币2
1	0.14	0.86
2	0.61	0.39
3	0.94	0.06
4	0.14	0.86
5	0.61	0.39

上表中的右兩清單示期望值。看第一行，0.86表示，從期望的角度看，這輪抛擲使用硬币2的機率是0.86。相比于前面的方法，我們按照最大似然機率，直接将第1輪估計為用的硬币2，此時的我們更加謹慎，我們隻說，有0.14的機率是硬币1，有0.86的機率是硬币2，不再是非此即彼。這樣我們在估計P1或者P2時，就可以用上全部的資料，而不是部分的資料，顯然這樣會更好一些。

這一步，我們實際上是估計出了z的機率分布，這步被稱作E步。

結合下表：

硬币	結果	統計
Unknown	正正反正反	3正-2反
Unknown	反反正正反	2正-3反
Unknown	正反反反反	1正-4反
Unknown	正反反正正	3正-2反
Unknown	反正正反反	2正-3反

我們按照期望最大似然機率的法則來估計新的P1和P2：

以P1估計為例，第1輪的3正2反相當于

0.14*3=0.42正

0.14*2=0.28反

依次算出其他四輪，清單如下：

輪數	正面	反面
1	0.42	0.28
2	1.22	1.83
3	0.94	3.76
4	0.42	0.28
5	1.22	1.83
總計	4.22	7.98

P1=4.22/(4.22+7.98)=0.35

可以看到，改變了z值的估計方法後，新估計出的P1要更加接近0.4。原因就是我們使用了所有抛擲的資料，而不是之前隻使用了部分的資料。

這步中，我們根據E步中求出的z的機率分布，依據最大似然機率法則去估計P1和P2，被稱作M步。

101、總結

以上，我們用一個實際的小例子，來實際示範了EM算法背後的idea，共性存于個性之中，通過這個例子，我們可以對EM算法究竟在幹什麼有一個深刻感性的認識，掌握EM算法的思想精髓。

Reference:

http://pan.baidu.com/s/1i4NfvP7