HMM基本算法

隱馬爾科夫模型（Hidden Markov Model，以下簡稱HMM）作為語音信號的一種統計模型，在語音處理的各個領域中獲得了廣泛的應用。當然，隨着目前深度學習的崛起，尤其是RNN，LSTM等神經網絡序列模型的火熱，HMM的地位有所下降。但是作為一個經典的模型，學習HMM的模型和對應算法，對我們解決問題建模的能力提高以及算法思路的拓展還是很好的。本篇將介紹HMM應用到語音處理中經常會面臨的3大基本問題的解決方案，也就是給出三種基本算法：前向-后向算法、Baum_Welch算法以及Viterbi解碼算法。

前向-后向算法

​ 該算法用來計算給定一個觀察值序列 \(O=o_1,o_2,\ldots,o_T\) 以及一個模型 \(\lambda=(\pi,A,B)\) 時，由模型 \(\lambda\) 計算出O的概率 $ P(O| \lambda ) $ 。
HMM組成如下圖所示：

\(P(O|\lambda)\) 最直接的求取方式如下：

對一個固定的狀態序列 \(Q=q_1,q_2,\ldots,q_T\) ,有

\[P(O|Q,\lambda)=\mathop{\Pi}_{t=1}^TP(o_t|q_t,\lambda) = b_{q_1}(o_1)b_{q_2}(o_2) \ldots b_{q_T}(o_T)\tag{1-1} \]

其中，

\[b_{q_t}(o_t)=\left. b_{jk} \right|_{q_t = \theta_j,o_t=V_k}， 1 \leq t \leq T \tag{1-2} \]

而對給定\(\lambda\)，產生Q的概率為

\[P（Q|\lambda）= \pi_{q_1}a_{q_1q_2} \cdots a_{q_T-1q_T}\tag{1-3} \]

因此，所求概率為

\[\begin{eqnarray*} P(O|\lambda) & = & \sum_{all Q}P(O|Q,\lambda)P(Q|\lambda) \\ & = & \sum_{q_1,q_2,\ldots,q_T} \pi_{q_1}b_{q_1}(o_1)a_{q_T-1q_T}b_{q_T}(o_T) \end{eqnarray*}\tag{1-4} \]

​ 顯而易見，上式的計算量是驚人的，大約為\(2TN^T\)數量級，當\(N=5,T=100\)時，計算量達\(10^{72}\)，這是完全不能接受的。在此情況下，要想求出\(P(O|\lambda)\)，就必須尋求更為有效的算法，前向后向算法就是解決這一問題的一種有效算法。
​ 后面算法中為方便表示，對狀態\(\theta_i\)的形式簡記為\(i\)。

前向算法

定義前向變量為

\[\alpha_t(i) = P(o_1o_2 \cdots o_t,q_t=i|\lambda)\tag{1-5} \]

那么有，

• 初始化：對 \(1 \leq i \leq N\) ,有

\[\alpha_1(j) = \pi_ib_i(o_1)\tag{1-6} \]

• 遞推：對 \(1 \leq t \leq T-1,1 \leq j \leq N\),有

\[\alpha_t+1(j) = \left[\sum_{i=1}^N \alpha_t(i)a_{ij}\right]b_j(o_{t+1})\tag{1-7} \]

• 終止：

\[P(O|\lambda)=\sum_{i=1}^N \alpha_T(i)\tag{1-8} \]

式（1-7）中

\[b_j(o_t+1)=\left. b_{jk} \right|_{o_{t+1}=V_k}\tag{1-9} \]

這種算法計算量大為減少，只需要\(N^2T\)次運算。對\(N=5,T=100\)時，只需要大約3000次乘法計算。它是一種格形結構，圖給出前向算法示意圖：

后向算法

與前向算法類似，定義后向變量為

\[\beta_t(i)=P(o_{t+1}o_{t+2} \cdots o_T|q_t=i,\lambda)\tag{1-10} \]

類似前向算法，有

• 初始化：對\(1 \le i \le N\)，有

\[\beta_T(i)=1\tag{1-11} \]

• 遞推：對\(t=T-1,T-2,\cdots,1,1 \leq i \leq N\),有

\[\beta_t(i) = \sum_{j=1}^N a_{ij}b_j(o_{t+1}) \beta_{t+1}(j)\tag{1-12} \]

• 終止：

\[P(O|\lambda)=\pi_ib_i(o_1)\beta_1(i)\tag{1-13} \]

后向算法的計算量大約也在\(N^2T\)數量級，它也是一種格形結構，HMM后向算法中計算后向變量\(\beta_t(i)\)的遞推過程如下圖所示。

Viterbi算法

這個算法解決對給定一個觀察值序列\(O=o_1o_2 \cdots o_T\)和一個模型\(\lambda=(\pi,A,B)\)，如何確定一個最佳狀態序列\(Q^*=q_1^*,q_2^*,\cdots,q_T^*\)的問題。
“最佳”的意義有很多種，由不同的定義可得到不同的結論，這里討論的最佳意義上的狀態序列\(Q^*\)是指使\(P(Q,O|\lambda)\)最大時確定的狀態序列\(Q^*\)。這可用Viterbi算法來實現，其描述如下：
定義\(\delta_t(i)\)為時刻t時沿一條路徑\(q_1,q_2,\cdots,q_t\)，且q_t=i,產生出\(o_1o_2\cdots o_t\)的最大概率，即有

\[\delta_t(i)=\mathop{max}_{q_1,q_2,\cdots,q_{t-1}}P(q_1q_2 \cdots q_t,q_t=i,o_1o_2\cdots o_t|\lambda)\tag{2-1} \]

那么，求取最佳狀態序列$Q^* $的過程為

• 初始化：對$1 \le i \le N $,有

\[\delta_t(i)=\pi_i b_i(o_t)\tag{2-2} \]

\[\varphi_1(i)=0\tag{2-3} \]

• 遞推：對$2 \le t \le T，1 \le j \le N $，有

\[\delta_t(j)=\mathop{max}_{1 \le i \le N}[\delta_{t-1}(i)a_{ij}]b_j(o_t)\tag{2-4} \]

\[\varphi_t(j)=\mathop{argmax}_{1 \le i \le N}[\delta_{t-1}(i)a_{ij}]\tag{2-5} \]

• 終止：

\[P^*=\mathop{max}_{1\le i \le N}[\delta_T(i)]\tag{2-6} \]

\[q_T^*=\mathop{argmax}_{1 \le i \le N}[\delta_T(i)]\tag{2-7} \]

• 路徑回溯，確定最佳狀態序列：

\[q_t^*=\varphi_{t+1}(q_{t+1}^*),t=T-1,T-2,\cdots,1\tag{2-8} \]

對語音處理應用而言，\(P(Q,O|\lambda)\)動態范圍很大，或者說不同的Q使\(P(Q,O|\lambda)\)的值差別很大，而\(\mathop{\max}_q P(Q,O|\lambda)\)事實上是\({\sum}_qP(Q,O|\lambda)\)中舉足輕重的唯一成分，因此常常等價的使用\(max_qP(Q,O|\lambda)\)近似\(\sum_qP(Q,O|\lambda)\),那么，Viterbi 算法也就能用來計算\(P(O|\lambda)\)。
此外，上述的Viterbi算法也是一種格形結構，而且類似於前向算法。同樣，由后向算法的思想出發，亦可推導出Viterbi算法的另一種實現方式。

Baum-Welch算法

這個算法用於解決HMM訓練問題，即HMM參數估計問題。可描述為：給定一個觀察虛列\(O=o_1,o_2,\cdots,o_T\)，確定一個\(\lambda=(\pi,A,B)\)，使\(P(O|\lambda)\)最大。
顯然，由前向算法和后向算法定義的前向和后向變量，即式(1-5)、(1-10)，有

\[P(O|\lambda)=\sum_{i=1}^N \sum_{j=1}^N \alpha_t(i) a_{ij} b_j(o_{t+1}) \beta_{t+1}(j),1 \le t \le T-1 \tag{3-1} \]

求取\(\lambda\)使\(P(O|\lambda)\)最大是一個泛函極值問題。但由於給定的訓練序列有限，因而不存在一個估計\(\lambda\)的最佳方法。在這種情況下，Baum-Welch算法利用遞歸的思想，使\(P(O|\lambda)\)局部極大，最后得到模型參數\(\lambda=(\pi,A,B)\)。此外，用梯度方法也可以達到類似的目的。
定義\(\xi_t(i,j)\)為給定訓練序列O和模型\(\lambda\)時，HMM模型在t時刻處於i狀態，t+1時刻處於j狀態的概率，即

\[\xi_t(i,j)=P(q_t=i,q_{t+1}=j|O,\lambda)\tag{3-2} \]

可以推導出

\[\xi_t(i,j)=[\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)]\tag{3-3} \]

那么，HMM模型在時刻 t 處於 i 狀態的概率為

\[\gamma_t(i)=P(q_t=i|O,\lambda)=\sum_{j=1}^N \xi_t(i,j)=\alpha_t(i)\beta_t(i)/P(O|\lambda)\tag{3-4} \]

因此，\(\mathop{\sum}_{t=1}^{T-1} \gamma_t(i)\)表示從 i 狀態轉移出去次數的期望值 ，而\(\sum_{t=1}^{T-1}\xi_t(i,j)\)表示從 i 狀態轉移到狀態 j 次數的期望值。由此，導出Baum-Welch算法中著名的重估公式，即

\[\mathop{\pi_i}^-=\gamma_1(i)\tag{3-5} \]

\[\mathop{a_{ij}}^-=\mathop{\sum}_{t=1}^{T-1}\xi_t(i,j) / \mathop{\sum}_{t=1}^{T-1}\gamma_t(i)\tag{3-6} \]

\[\mathop{b_{jk}}^-=\mathop{\sum}_{t=1,o_k=V_k}^T\gamma_t(i) / \mathop{\sum}_{t=1}^T\gamma_t(i)\tag{3-7} \]

那么，HMM參數\(\lambda=(\pi,A,B)\)的求取過程為，根據觀察值序列O和選取的初始模型\(\lambda=(\pi,A,B)\)，由重估公式(3-5)、(3-6)和(3-7)求得一組新參數\(\mathop{\pi_i}^-\)，\(\mathop{a_{ij}}^-\)，\(\mathop{b_{jk}}^-\)，亦即得到了一個新模型\(\mathop{\lambda}^-=(\mathop{\pi}^-,\mathop{A}^-,\mathop{B}^-)\)。由下面的證明可以看出，\(P(O|\mathop{\lambda}^-)>P(O|\lambda)\)，即由重估公式得到的\(\mathop{\lambda}^-\)比\(\lambda\)在表示觀察值序列O方面要好，那么重復這個過程，逐步改進模型參數，直到\(P(O|\mathop{\lambda}^-)\)收斂，即不再明顯增大，此時的\(\mathop{\lambda}^-\)即為所求的模型。
應當指出，HMM訓練，或者稱為參數估計問題，是HMM在語音處理中應用的關鍵問題，與前面討論的兩個問題相比，這是最困難的問題。Baum-Welch算法只是得到廣泛應用的解決這一問題的經典方法，但並不是唯一的，也遠不是最完善的方法。

重估算法證明

Baum-Welch算法一般的證明方式是構造Q函數作為輔助函數，這和EM算法中的Q函數構造是一樣的。這里簡要介紹EM算法，這樣在介紹Q函數構造原則的同時，也將證明Baum-Welch算法的有效性。

EM算法是一種從“不完全數據”中求解模型分布參數的最大似然估計方法。所謂“不完全數據”一般是指兩種情況：一種情況是由於觀測過程本身的限制或者錯誤造成觀測數據成為有錯漏的“不完全”數據；另一種情況是對參數的似然函數直接優化十分困難，而引入額外的參數（隱含的或丟失的）后優化就會比較容易。於是定義原始觀測數據加上額外參數組成“完全數據”，原始觀測數據自然就成為“不完全數據”。在模式識別及相關領域中，后一種情況更常見。由於所要優化的似然函數是一個隨機變量，直接對其求最大化不好計算，但它的期望卻是一個確定性的函數，優化起來相對容易。這就是EM算法的基本思路。該算法中包含兩個主要方面：一是求期望，用E來表示，一是最大化，用M來表示；這也是這種方法名字的由來。EM算法在矢量優化和HMM模型的參數估計中有着很重要的作用。下面以HMM模型訓練為例介紹EM算法。

在HMM模型中，數據是由兩部分組成，一部分是可以觀測到的數據，如觀測特征序列O，稱為可觀測數據；另一部分數據無法觀測到，如狀態序列Q，稱為隱含序列。由這兩部分數據可以共同構成一個完全數據集（O,Q）。EM算法的目的是通過迭代的將完全數據集的對數似然度期望最大化，來實現對可觀測數據的對數似然度的最大化。

根據Bayes公式，完全數據集的似然度和可觀測數據集的似然度之間存在以下關系：

\[P(O,Q|\lambda)=P(Q|O,\lambda)P(O|\lambda)\tag{4-1} \]

那么，觀測數據的對數似然度可以表示為

\[logP(O|\lambda)=logP(O,Q|\lambda)-logP(Q|O,\lambda)\tag{4-2} \]

對於兩個參數集\(\lambda\)和\(\mathop{\lambda}^-\)，在已知O和\(\lambda\)的情況下，對公式（4-2）在完全數據集上求期望，則

\[E[logP(o|\mathop{\lambda}^-)|O,\lambda]=E[logP(O,Q|\mathop{\lambda}^-)|O,\lambda]-E[logP(Q|O,\mathop{\lambda}^-)|O,\lambda]\tag{4-3} \]

其中式（4-3）的左半部分為

\[\begin{eqnarray*} E[logP(o|\mathop{\lambda}^-)|O,\lambda] & = & \int logP(O|\mathop{\lambda}^-)P(Q|O,\lambda)dQ \\ & = & logP(O|\mathop{\lambda}^-) \\ & = & L(O,\mathop{\lambda}^-)\tag{4-4} \end{eqnarray*} \]

令

\[Q(\lambda,\mathop{\lambda}^-)=E[logP(O,Q|\mathop{\lambda}^-)|O,\lambda]=\int logP(O,Q|\mathop{\lambda}^-)P(Q|O,\lambda)dQ\tag{4-5} \]

而

\[H(\lambda,\mathop{\lambda}^-)=E[logP(Q|O,\mathop{\lambda}^-)|O,\lambda]=\int logP(Q|\mathop{\lambda}^-)P(Q|O,\lambda)dQ\tag{4-6} \]

則式（4-4）變為

\[L(O,\mathop{\lambda}^-)=Q(\lambda,\mathop{\lambda}^-)-H(\lambda,\mathop{\lambda}^-)\tag{4-7} \]

由Jensen不等式可以知道，\(H(\lambda,\mathop{\lambda}^-) \le H(\lambda,\lambda)\)。EM算法的基本原理在於，如果可以保證\(Q(\lambda,\mathop{\lambda}^-) \ge Q(\lambda,\lambda)\)不等式成立，那么不等式\(L(O,\mathop{\lambda}^-) \ge L(O,\lambda)\)一定成立。從上面分析可以看出，Q()函數實質上就是完全數據集的對數似然度的期望，這樣，通過將Q()函數最大化就可以實現觀測數據的對數似然度的最大化。
當隱藏數據為離散時，Q()函數一般表示為

\[Q(\lambda,\mathop{\lambda}^-)=\sum_{所有Q} \frac {P(O,Q|\lambda)}{P(O|\lambda)} logP(O,Q|\mathop{\lambda}^-)\tag{4-8} \]

EM算法一般描述——給定一個當前的參數集\(\lambda\)，可以通過如下方式獲得新的參數集\(\mathop{\lambda}^-\):

• 選擇初始參數$ \lambda $；
• 求期望，即在給定的參數集$ \lambda $上求Q()函數；
• 最大化，選擇$ \lambda=\mathop{argmax_{\mathop{\lambda}^-}}Q(\lambda,\mathop{\lambda}^-)$。