NLP --- 隱馬爾可夫HMM（EM演算法（期望最大化演算法））

阿新 • • 發佈：2018-12-17

期望最大化 (Expectation Maximization) 演算法最初是由 Ceppellini[2] 等人 1950 年在討論基因頻率的估計的時候提出的。後來又被 Hartley[3] 和Baum[4] 等人發展的更加廣泛。目前引用的較多的是 1977 年 Dempster[5]等人的工作。它主要用於從不完整的資料中計算最大似然估計。後來經過其他學者的發展，這個演算法也被用於聚類等應用。

因此如果大家理解了K-means聚類，那麼這個演算法就容易理解了，這裡我們先一起來回歸一下K-means的聚類思想，其實他的聚類思想很簡單，我們來看一下演算法過程：

k-means的k就是最終聚集的簇數，這個要你事先自己指定。k-means在常見的機器學習演算法中算是相當簡單的，基本過程如下：

首先任取k個樣本點作為k個簇的初始中心；
對每一個樣本點，計算它們與k箇中心的距離，把它歸入距離最小的中心所在的簇；
等到所有的樣本點歸類完畢，重新計算k個簇的中心；
重複以上過程直至樣本點歸入的簇不再變動。

通過上面大家可以看到，為了找到合適的聚類中心點，k-means剛開始給的初始中心是任意的，然後不斷的更新迭代後按照上面演算法步驟不停的更改k箇中心的位置，直到他們都幾乎不再變化位置，也就是說當k個位置幾乎不再變化時基本上就說明了到達了真實的資料中心了。這裡有一個思想就是我剛開始不知道資料的最佳K箇中心點，但是我可以通過資料不停的迭代，進而達到理想的中心點，而這個思想也是EM演算法的核心思想。
下面先通過一個簡單的例子引入EM演算法的核心思想，然後在使用數學進行證明這樣做的合理性。

EM

考慮一個投擲硬幣的實驗：現在我們有兩枚硬幣 A 和 B，這兩枚硬幣和普通的硬幣不一樣，他們投擲出正面的概率和投擲出反面的概率不一定相同。我們將 A 和 B 投擲出正面的概率分別記為 $\theta _A$ 和 $\theta _B$ 。我們現在獨立地做 5 次試驗：隨機的從這兩枚硬幣中抽取 1 枚，投擲 10 次，統計出現正面的次數。那麼我們就得到了如表格1的實驗結果。

試驗代號

投擲的硬幣

出現正面的次數

在這個實驗中，我們記錄兩組隨機變數：

$X = (X_1,X_2,X_3,X_4,X_5)$

$Z = (Z_1,Z_2,Z_3,Z_4,Z_5)$

其中 $X_i\in \left \{ 1,2,3,4,5,6,7,8,9,10 \right \}$ ，代表試驗 $i$ 中出現的正面的次數， $Z_i\in \left \{ A,B \right \}$ 代表這次試驗投擲的是硬幣A還是B。

我們的目標是通過這個實驗來估計 $\theta =(\theta _A,\theta _B)$ 的數值。這個實驗中的引數估計就是有完整資料的引數估計，這個是因為我們不僅僅知道每次試驗中投擲出正面的次數，我們還知道每次試驗中投擲的是硬幣 A 還是 B。
一個很簡單也很直接的估計 θ 的方法如公式(1)所示。

實際上這樣的估計就是統計上的極大似然估計 (maximum likelihood estimation) 的結果。用 $P(X,Z|\theta )$ 來表示 $X,Z$ 的聯合概率分佈（其中帶有引數 θ），那麼對於上面的實驗，我們可以計算出他們出現我們觀察到的結果即 $x^0 = (5,9,8,4,7),z^0 = (B,A,A,B,A)$ 的概率。函式 $P(X=x^{(0)},Z = z^{0}|\theta )$ 就叫做 θ 的似然函式。我們將它對 θ 求偏導並令偏導數為 0，就可以得到如(1) 的結果。

對其求偏導，這裡只對 $\theta _A$ 求偏導， $\theta _B$ 也是一樣的求，這裡我們把係數統一使用C來表示了，關於 $\theta _B$ 的式子使用 $f(\theta _B)$ 進行代替，然後整理一下上式得：

$P(X=x^0,Z = z^0|\theta ) = Cf(\theta _B)\theta _A^{24}(1-\theta _A)^6$

對其求偏導，並令其為0的：

$\frac{\partial P(X=x^0,Z = z^0|\theta )}{\partial \theta _A} = Cf(\theta _B)[24\theta _A^{23}(1-\theta _A)^6-6\theta ^{24}(1-\theta _A)^5]=0$

此時化簡解得 $\theta _A$ 如下：

$24(1-\theta _A) - 6\theta _A = 0$

$30\theta _A = 24$

$\theta _A = 0.8$

同理可求 $\theta _B =0.45$

由此可見上面利用均值和這裡使用最大釋然函式的結果是一樣的，因此這裡的（1）式稱為極大釋然估計。這裡大家應該明白了吧。下面我們把難度加大：

我們將這個問題稍微改變一下，我們將我們所觀察到的結果修改一下：我們現在只知道每次試驗有幾次投擲出正面，但是不知道每次試驗投擲的是哪個硬幣，也就是說我們只知道表1中第一列和第三列。這個時候我們就稱 $Z$ 為隱藏變數 (Hidden Variable)， $X$ 稱為觀察變數 (Observed Variable)。這個時候再來估計引數 $\theta _A$ 和 $\theta _B$ ，就沒有那麼多資料可供使用了，這個時候的估計叫做不完整資料的引數估計（這裡X就可認為是HMM中的可見變量了，Z就是隱藏的狀態變量了）。
如果我們這個時候有某種方法（比如，正確的猜到每次投擲硬幣是 A還是 B），這樣的話我們就可以將這個不完整的資料估計變為完整資料估計。
當然我們如果沒有方法來獲得更多的資料的話，那麼下面提供了一種在這種不完整資料的情況下來估計引數 θ 的方法。我們用迭代的方式來進行：

(1) 我們先賦給 θ 一個初始值，這個值不管是經驗也好猜的也好，反正我們給它一個初始值。在實際使用中往往這個初始值是有其他演算法的結果給出的，當然隨機給他分配一個符合定義域的值也可以。這裡我們就給定 $\theta _A$ = 0.7, $\theta _B$ = 0.4
(2) 然後我們根據這個來判斷或者猜測每次投擲更像是哪枚硬幣投擲的結果。比如對於試驗 1，如果投擲的是 A，那麼出現 5 個正面的概率為 $C_{10}^5\times 0.7^5\times (1-0.7^5)\approx 0.1029$ ；如果投擲的是B，出現5個正面的概率為 $C_{10}^5\times 0.4^5\times (1-0.4)^5\approx 0.2007$ ；基於試驗1的試驗結果，可以判斷這個試驗投擲的是硬幣 A 的概率為 0.1029/(0.1029+0.2007)=0.3389，是 B 的概率為 0.2007/(0.1029+0.2007)= 0.6611。因此這個結果更可能是投擲 B 出現的結果。
(3) 假設上一步猜測的結果為 B,A,A,B,A，那麼根據這個猜測，可以像完整資料的引數估計一樣 (公式2) 重新計算 θ 的值。

這樣一次一次的迭代 2-3 步驟直到收斂，我們就得到了 θ 的估計。這裡大家應該能感受到這和K-means很像對吧，你可能有疑問，這個方法靠譜麼？事實證明，它確實是靠譜的。下面我們就從數學角度來說明他的合理性。

數學證明

首先來明確一下我們的目標：我們的目標是在觀察變數 X 和給定觀察樣本 $x_1,x_2,...,x_n$ 的情況下，極大化對數似然函式

這個函式大家應該都能看懂吧，看不懂的請返回看我的上一篇部落格。

其中只包含觀察變數的概率密度函式（其實是求邊緣密度，即對另一個變數進行求和就可以了得到邊緣密度了）：

其中 Z 為隱藏隨機變數, $\left \{ Z_j \right \}$ 是 Z 的所有可能的取值。那麼把（6）式代入（5）式：

這裡我們引入了一組引數（不要怕多，我們後面會處理掉它的）α，它滿足∀可能的j, $\alpha _j$ ∈ (0,1]和
$\sum _j\alpha _j=1$ 到這裡，先介紹一個凸函式的性質，或者叫做凸函式的定義。f(x) 為凸函式， $\forall i = 1,2,3,..,n,\lambda _i\in [0,1],\sum_{i=1}^{n}\lambda _i = 1$ ,對 f(x) 定義域中的任意 n 個 $x_1,x_2,...,x_n$ 有 :

$f(\sum_{i=1}^{n}\lambda _ix_i)\leq \sum_{i=1}^{n}\lambda _if(x_i)$ $\left ( 8 \right )$

上式稱為Jensen不等式

對於嚴格凸函式，上面的等號只有在 $x_1= x_2=...=x_n$ 的時候成立。關於凸函式的其他性質不再贅述。對數函式是一個嚴格凸函式。因而我們可以有下面這個結果：

上式是 $\left ( 7 \right )$ 式根據 $\left ( 8 \right )$ 式得到的，只是不同的是ln是凹函式，但是同樣適用Jensen不等式，只是不等號方向不一樣罷了。

現在我們根據等號成立的條件來確定 $\alpha _j$ 即

這裡需仔細推一下，我們知道 $\left ( 9 \right )$ 的成立條件是 $x_1= x_2=...=x_n$ ，那麼聯合概率密度 $P(X_i,Z_j)$ 是不變的，等於一個常數，和j無關，和j無關的意思是當j確定時，無論i等於几上式都是一個常數，大家可以這樣理解，因此可以寫成 $\left ( 10 \right )$ 式，對上式的j求和同時把 $\alpha _j$ 乘到右邊：

$\sum_{j}P(X_i=x_i,Z=z_j) = \sum _j\alpha _j*c = c$

我們仔細看看等式的左邊，他是對j求和，而且P是聯合概率密度，對某一個維度求和不就是求邊緣密度嗎？所以上式可以寫為如下：

$\sum_{j}P(X_i=x_i,Z=z_j) = P(X_i=x_i)=\sum _j\alpha _j*c = c$

也就是說：

$P(X_i=x_i) = c$

此時把上式代入（10）式，可得（11）式

其中 c 是一個與 j 無關的常數。因為 $\sum _j\alpha _j=1$ 稍作變換就可以得到

上式是什麼意思呢？大家還記得條件概率的定義嗎？

$P(B|A) = \frac{P(AB)}{P(A)}$

因此（11）式就是條件概率的表示式即為 $P(Z=z_j|X=x_i) = \frac{P(X_i=x_i,Z=z_j)}{P(X=x_i)}$ .，這個式子代表什麼意思呢？其實就是在我們知道了扔的硬幣為正的情況下這個硬幣為A或者為B的概率，就是再求隱藏的條件概率。

現在來解釋一下我們得到了什麼。 $\alpha _j$ 就是 $Z=z_j$ 在 $X=x_i$ 下的條件概率或者後驗概率。求 α 就是求隱藏隨機變數 Z 的條件分佈。總結一下目前得到的公式就是

好，得到上式以後，我們如何求 $\theta$ 呢，這裡需要和大家說明的是，這裡的 $\theta$ 在上式沒體現處理，他是包含在P和 $\alpha$ 中的，這裡大家需要理解，那我們如何求解引數 $\theta$ 使的 $l(\theta )$ 達到最大呢？這裡我們需要明確上式是很複雜的，直接求導不好求，因此需要想辦法進行間接求極大值，怎麼求呢？這裡給出一個簡單的例子：

$l(\theta )=d(\theta )q\left ( \theta \right )$

假設有一個待求極值的函式 $l(\theta )$ ，如上式，其中 $d(\theta )$ 很複雜無法求導，而 $q(\theta )$ 相對簡單一點，那麼針對這樣的函式我們如何求極值呢？這裡提供一種迭代的方法，首先我給 $d(\theta )$ 任意賦值一個 $\theta _0$ 的值，此時 $d(\theta _0)$ 就是一個常數了，那麼上式就可以寫成 $l(\theta )=d(\theta_0 )q\left ( \theta \right )$ ，此時在對 $l(\theta )$ 求導得到 $\theta _1$ ,再把 $\theta _1$ 代入 $d(\theta )$ ，在類似的計算，直到 $\theta$ 不再發生變化此時的 $\theta$ 就是所求的值，過程如下：