隱馬爾科夫模型HMM（二）前向後向演算法評估觀察序列概率
在隱馬爾科夫模型HMM（一）HMM模型中，我們講到了HMM模型的基礎知識和HMM的三個基本問題，本篇我們就關注於HMM第一個基本問題的解決方法，即已知模型和觀測序列，求觀測序列出現的概率。
1. 回顧HMM問題一：求觀測序列的概率
首先我們回顧下HMM模型的問題一。這個問題是這樣的。我們已知HMM模型的引數$\lambda = (A, B, \Pi)$。其中$A$是隱藏狀態轉移概率的矩陣，$B$是觀測狀態生成概率的矩陣， $\Pi$是隱藏狀態的初始概率分佈。同時我們也已經得到了觀測序列$O =\{o_1,o_2,...o_T\}$,現在我們要求觀測序列$O$在模型$\lambda$下出現的條件概率$P(O|\lambda)$。
乍一看，這個問題很簡單。因為我們知道所有的隱藏狀態之間的轉移概率和所有從隱藏狀態到觀測狀態生成概率，那麼我們是可以暴力求解的。
我們可以列舉出所有可能出現的長度為$T$的隱藏序列$I = \{i_1,i_2,...,i_T\}$,分佈求出這些隱藏序列與觀測序列$O =\{o_1,o_2,...o_T\}$的聯合概率分佈$P(O,I|\lambda)$，這樣我們就可以很容易的求出邊緣分佈$P(O|\lambda)$了。
具體暴力求解的方法是這樣的：首先，任意一個隱藏序列$I = \{i_1,i_2,...,i_T\}$出現的概率是：$$P(I|\lambda) = \pi_{i_1} a_{i_1i_2} a_{i_2i_3}... a_{i_{T-1}\;\;i_T}$$
對於固定的狀態序列$I = \{i_1,i_2,...,i_T\}$，我們要求的觀察序列$O =\{o_1,o_2,...o_T\}$出現的概率是：$$P(O|I, \lambda) = b_{i_1}(o_1)b_{i_2}(o_2)...b_{i_T}(o_T)$$
則$O$和$I$聯合出現的概率是：$$P(O,I|\lambda) = P(I|\lambda)P(O|I, \lambda) = \pi_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)...a_{i_{T-1}\;\;i_T}b_{i_T}(o_T)$$
然後求邊緣概率分佈，即可得到觀測序列$O$在模型$\lambda$下出現的條件概率$P(O|\lambda)$：$$P(O|\lambda) = \sum\limits_{I}P(O,I|\lambda) = \sum\limits_{i_1,i_2,...i_T}\pi_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)...a_{i_{T-1}\;\;i_T}b_{i_T}(o_T)$$
雖然上述方法有效，但是如果我們的隱藏狀態數$N$非常多的那就麻煩了，此時我們預測狀態有$N^T$種組合，演算法的時間複雜度是$O(TN^T)$階的。因此對於一些隱藏狀態數極少的模型，我們可以用暴力求解法來得到觀測序列出現的概率，但是如果隱藏狀態多，則上述演算法太耗時，我們需要尋找其他簡潔的演算法。
前向後向演算法就是來幫助我們在較低的時間複雜度情況下求解這個問題的。
2. 用前向演算法求HMM觀測序列的概率
前向後向演算法是前向演算法和後向演算法的統稱，這兩個演算法都可以用來求HMM觀測序列的概率。我們先來看看前向演算法是如何求解這個問題的。
前向演算法本質上屬於動態規劃的演算法，也就是我們要通過找到區域性狀態遞推的公式，這樣一步步的從子問題的最優解擴充到整個問題的最優解。
在前向演算法中，通過定義“前向概率”來定義動態規劃的這個區域性狀態。什麼是前向概率呢, 其實定義很簡單：定義時刻$t$時隱藏狀態為$q_i$, 觀測狀態的序列為$o_1,o_2,...o_t$的概率為前向概率。記為：$$\alpha_t(i) = P(o_1,o_2,...o_t, i_t =q_i | \lambda)$$
既然是動態規劃，我們就要遞推了，現在我們假設我們已經找到了在時刻$t$時各個隱藏狀態的前向概率，現在我們需要遞推出時刻$t+1$時各個隱藏狀態的前向概率。
從下圖可以看出，我們可以基於時刻$t$時各個隱藏狀態的前向概率，再乘以對應的狀態轉移概率，即$\alpha_t(j)a_{ji}$就是在時刻$t$觀測到$o_1,o_2,...o_t$，並且時刻$t$隱藏狀態$q_j$, 時刻$t+1$隱藏狀態$q_i$的概率。如果將想下面所有的線對應的概率求和，即$\sum\limits_{j=1}^N\alpha_t(j)a_{ji}$就是在時刻$t$觀測到$o_1,o_2,...o_t$，並且時刻$t+1$隱藏狀態$q_i$的概率。繼續一步，由於觀測狀態$o_{t+1}$只依賴於$t+1$時刻隱藏狀態$q_i$, 這樣$[\sum\limits_{j=1}^N\alpha_t(j)a_{ji}]b_i(o_{t+1})$就是在在時刻$t+1$觀測到$o_1,o_2,...o_t，o_{t+1}$，並且時刻$t+1$隱藏狀態$q_i$的概率。而這個概率，恰恰就是時刻$t+1$對應的隱藏狀態$i$的前向概率，這樣我們得到了前向概率的遞推關係式如下：$$\alpha_{t+1}(i) = \Big[\sum\limits_{j=1}^N\alpha_t(j)a_{ji}\Big]b_i(o_{t+1})$$
我們的動態規劃從時刻1開始，到時刻$T$結束，由於$\alpha_T(i)$表示在時刻$T$觀測序列為$o_1,o_2,...o_T$，並且時刻$T$隱藏狀態$q_i$的概率，我們只要將所有隱藏狀態對應的概率相加，即$\sum\limits_{i=1}^N\alpha_T(i)$就得到了在時刻$T$觀測序列為$o_1,o_2,...o_T$的概率。
下面總結下前向演算法。
輸入：HMM模型$\lambda = (A, B, \Pi)$，觀測序列$O=(o_1,o_2,...o_T)$
輸出：觀測序列概率$P(O|\lambda)$
1) 計算時刻1的各個隱藏狀態前向概率：$$\alpha_1(i) = \pi_ib_i(o_1),\; i=1,2,...N$$
2) 遞推時刻$2,3,...T$時刻的前向概率：$$\alpha_{t+1}(i) = \Big[\sum\limits_{j=1}^N\alpha_t(j)a_{ji}\Big]b_i(o_{t+1}),\; i=1,2,...N$$
3) 計算最終結果：$$P(O|\lambda) = \sum\limits_{i=1}^N\alpha_T(i)$$
從遞推公式可以看出，我們的演算法時間複雜度是$O(TN^2)$，比暴力解法的時間複雜度$O(TN^T)$少了幾個數量級。
3. HMM前向演算法求解例項
這裡我們用隱馬爾科夫模型HMM（一）HMM模型中盒子與球的例子來顯示前向概率的計算。
我們的觀察集合是:$$V=\{紅，白\}，M=2$$
我們的狀態集合是：$$Q =\{盒子1，盒子2，盒子3\}， N=3 $$
而觀察序列和狀態序列的長度為3.
初始狀態分佈為：$$\Pi = (0.2,0.4,0.4)^T$$
狀態轉移概率分佈矩陣為：
$$A = \left( \begin{array} {ccc} 0.5 & 0.2 & 0.3 \\ 0.3 & 0.5 & 0.2 \\ 0.2 & 0.3 &0.5 \end{array} \right) $$
觀測狀態概率矩陣為：
$$B = \left( \begin{array} {ccc} 0.5 & 0.5 \\ 0.4 & 0.6 \\ 0.7 & 0.3 \end{array} \right) $$
球的顏色的觀測序列:$$O=\{紅，白，紅\}$$
按照我們上一節的前向演算法。首先計算時刻1三個狀態的前向概率：
時刻1是紅色球，隱藏狀態是盒子1的概率為：$$\alpha_1(1) = \pi_1b_1(o_1) = 0.2 \times 0.5 = 0.1$$
隱藏狀態是盒子2的概率為：$$\alpha_1(2) = \pi_2b_2(o_1) = 0.4 \times 0.4 = 0.16$$
隱藏狀態是盒子3的概率為：$$\alpha_1(3) = \pi_3b_3(o_1) = 0.4 \times 0.7 = 0.28$$
現在我們可以開始遞推了，首先遞推時刻2三個狀態的前向概率：
時刻2是白色球，隱藏狀態是盒子1的概率為：$$\alpha_2(1) = \Big[\sum\limits_{i=1}^3\alpha_1(i)a_{i1}\Big]b_1(o_2) = [0.1*0.5+0.16*0.3+0.28*0.2 ] \times 0.5 = 0.077$$
隱藏狀態是盒子2的概率為：$$\alpha_2(2) = \Big[\sum\limits_{i=1}^3\alpha_1(i)a_{i2}\Big]b_2(o_2) = [0.1*0.2+0.16*0.5+0.28*0.3 ] \times 0.6 = 0.1104$$
隱藏狀態是盒子3的概率為：$$\alpha_2(3) = \Big[\sum\limits_{i=1}^3\alpha_1(i)a_{i3}\Big]b_3(o_2) = [0.1*0.3+0.16*0.2+0.28*0.5 ] \times 0.3 = 0.0606$$
繼續遞推，現在我們遞推時刻3三個狀態的前向概率：
時刻3是紅色球，隱藏狀態是盒子1的概率為：$$\alpha_3(1) = \Big[\sum\limits_{i=1}^3\alpha_2(i)a_{i1}\Big]b_1(o_3) = [0.077*0.5+0.1104*0.3+0.0606*0.2 ] \times 0.5 = 0.04187$$
隱藏狀態是盒子2的概率為：$$\alpha_3(2) = \Big[\sum\limits_{i=1}^3\alpha_2(i)a_{i2}\Big]b_2(o_3) = [0.077*0.2+0.1104*0.5+0.0606*0.3 ] \times 0.4 = 0.03551$$
隱藏狀態是盒子3的概率為：$$\alpha_3(3) = \Big[\sum\limits_{i=1}^3\alpha_3(i)a_{i3}\Big]b_3(o_3) = [0.077*0.3+0.1104*0.2+0.0606*0.5 ] \times 0.7 = 0.05284$$
最終我們求出觀測序列:$O=\{紅，白，紅\}$的概率為：$$P(O|\lambda) = \sum\limits_{i=1}^3\alpha_3(i) = 0.13022 $$
4. 用後向演算法求HMM觀測序列的概率
熟悉了用前向演算法求HMM觀測序列的概率，現在我們再來看看怎麼用後向演算法求HMM觀測序列的概率。
後向演算法和前向演算法非常類似，都是用的動態規劃，唯一的區別是選擇的區域性狀態不同，後向演算法用的是“後向概率”，那麼後向概率是如何定義的呢？
定義時刻$t$時隱藏狀態為$q_i$, 從時刻$t+1$到最後時刻$T$的觀測狀態的序列為$o_{t+1},o_{t+2},...o_T$的概率為後向概率。記為：$$\beta_t(i) = P(o_{t+1},o_{t+2},...o_T| i_t =q_i , \lambda)$$
後向概率的動態規劃遞推公式和前向概率是相反的。現在我們假設我們已經找到了在時刻$t+1$時各個隱藏狀態的後向概率$\beta_{t+1}(j)$，現在我們需要遞推出時刻$t$時各個隱藏狀態的後向概率。如下圖，我們可以計算出觀測狀態的序列為$o_{t+2},o_{t+3},...o_T$， $t$時隱藏狀態為$q_i$, 時刻$t+1$隱藏狀態為$q_j$的概率為$a_{ij}\beta_{t+1}(j)$, 接著可以得到觀測狀態的序列為$o_{t+1},o_{t+2},...o_T$， $t$時隱藏狀態為$q_i$, 時刻$t+1$隱藏狀態為$q_j$的概率為$a_{ij}b_j(o_{t+1})\beta_{t+1}(j)$, 則把下面所有線對應的概率加起來，我們可以得到觀測狀態的序列為$o_{t+1},o_{t+2},...o_T$， $t$時隱藏狀態為$q_i$的概率為$\sum\limits_{j=1}^{N}a_{ij}b_j(o_{t+1})\beta_{t+1}(j)$，這個概率即為時刻$t$的後向概率。
這樣我們得到了後向概率的遞推關係式如下：$$\beta_{t}(i) = \sum\limits_{j=1}^{N}a_{ij}b_j(o_{t+1})\beta_{t+1}(j)$$
現在我們總結下後向演算法的流程,注意下和前向演算法的相同點和不同點：
輸入：HMM模型$\lambda = (A, B, \Pi)$，觀測序列$O=(o_1,o_2,...o_T)$
輸出：觀測序列概率$P(O|\lambda)$
1) 初始化時刻$T$的各個隱藏狀態後向概率：$$\beta_T(i) = 1,\; i=1,2,...N$$
2) 遞推時刻$T-1,T-2,...1$時刻的後向概率：$$\beta_{t}(i) = \sum\limits_{j=1}^{N}a_{ij}b_j(o_{t+1})\beta_{t+1}(j),\; i=1,2,...N$$
3) 計算最終結果：$$P(O|\lambda) = \sum\limits_{i=1}^N\pi_ib_i(o_1)\beta_1(i)$$
此時我們的演算法時間複雜度仍然是$O(TN^2)$。
5. HMM常用概率的計算
利用前向概率和後向概率，我們可以計算出HMM中單個狀態和兩個狀態的概率公式。
1）給定模型$\lambda$和觀測序列$O$,在時刻$t$處於狀態$q_i$的概率記為:$$\gamma_t(i) = P(i_t = q_i | O,\lambda) = \frac{P(i_t = q_i ,O|\lambda)}{P(O|\lambda)} $$
利用前向概率和後向概率的定義可知：$$P(i_t = q_i ,O|\lambda) = \alpha_t(i)\beta_t(i)$$
於是我們得到：$$\gamma_t(i) = \frac{ \alpha_t(i)\beta_t(i)}{\sum\limits_{j=1}^N \alpha_t(j)\beta_t(j)}$$
2）給定模型$\lambda$和觀測序列$O$,在時刻$t$處於狀態$q_i$，且時刻$t+1$處於狀態$q_j$的概率記為:$$\xi_t(i,j) = P(i_t = q_i, i_{t+1}=q_j | O,\lambda) = \frac{ P(i_t = q_i, i_{t+1}=q_j , O|\lambda)}{P(O|\lambda)} $$
而$P(i_t = q_i, i_{t+1}=q_j , O|\lambda)$可以由前向後向概率來表示為:$$P(i_t = q_i, i_{t+1}=q_j , O|\lambda) = \alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)$$
從而最終我們得到$\xi_t(i,j)$的表示式如下：$$\xi_t(i,j) = \frac{\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)}{\sum\limits_{r=1}^N\sum\limits_{s=1}^N\alpha_t(r)a_{rs}b_s(o_{t+1})\beta_{t+1}(s)}$$
3) 將$\gamma_t(i)$和$\xi_t(i,j)$在各個時刻$t$求和，可以得到：
在觀測序列$O$下狀態$i$出現的期望值$\sum\limits_{t=1}^T\gamma_t(i)$
在觀測序列$O$下由狀態$i$轉移的期望值$\sum\limits_{t=1}^{T-1}\gamma_t(i)$
在觀測序列$O$下由狀態$i$轉移到狀態$j$的期望值$\sum\limits_{t=1}^{T-1}\xi_t(i,j)$
上面這些常用的概率值在求解HMM問題二，即求解HMM模型引數的時候需要用到。我們在這個系列的第三篇來討論求解HMM引數的問題和解法。
（歡迎轉載，轉載請註明出處。歡迎溝通交流： liujianping-ok@163.com）