隨機場(Random Field)是在同一向量空間上一組隨機變量組成的集合。通常情況下我們將存在相關關系的一組隨機變量作為隨機場進行研究。通常我們使用λ來表示一個參數給定的隨機場。
概率圖模型采用圖來表示隨機變量之間的相關關系， 最常見的概率圖模型是采用有向無環圖的貝葉斯網絡和采用無向圖的馬爾科夫隨機場。
馬爾科夫隨機場遵循馬爾科夫假設， 即在隨機過程中下一個狀態僅與當前狀態有關， 與歷史狀態無關。
最簡單的馬爾科夫隨機場的概率圖是鏈式的， 除了兩端外其它節點上的隨機變量只與相鄰的兩個節點有關， 我們將這種模型稱為馬爾科夫鏈。
馬爾科夫隨機場常用於分析隨機過程， 隨機場上每一個節點代表隨機系統在某個時刻可能處於的一個狀態， 節點對應的隨機變量代表了某個時刻系統處於該狀態的概率。
給定系統初始時刻處於各個狀態的概率，每個節點依據與它相鄰的節點的先驗概率計算下一時刻系統處於該狀態的概率。 逐次迭代，觀測變量的概率分布將會生成一個可以描述隨機過程的時間序列， 我們稱為觀測序列。
隱馬爾科夫模型(Hiden Markov Model, HMM)增加了一組觀測隨機變量， 一條隱藏的馬爾科夫鏈將產生一個不可觀測的隨機過程。鏈上每個節點先驗結果將在觀測變量上產生一個可觀測的概率分布，我們可以借助觀測變量分析隱馬爾科夫鏈上的隨機過程。
除了馬爾科夫假設外， HMM顯然依賴一個嚴格的觀測獨立性假設，即觀測結果僅依賴馬爾科夫鏈當前的狀態與其它時刻的狀態和觀測無關。
我們可以使用三個要素來描述一個含有N個狀態M個觀測值的HMM模型:
-
初始狀態概率向量π：N維向量， 表示初始時刻隨機過程處於對應狀態的概率。
-
狀態轉移概率矩陣A: N * N 矩陣， a[i,j]表示從狀態i轉移到j的概率.
-
觀測轉移概率矩陣B: N * M矩陣， b[j, k]表示處於狀態j時產生觀測值k的概率。
HMM有三類基本問題：
-
概率計算問題: 計算給定模型下出現特定觀測序列的概率。
-
學習問題, 又稱參數估計問題： 給定觀測序列， 調整模型參數使得給定觀測序列出現的概率P(O|λ)最大。
-
預測問題， 又稱解碼問題： 給定模型， 求使觀測序列出現概率最大的狀態序列。
示例： 這里以天氣與地面濕度為例來說明HMM中三類問題求解的方法。
我們將天氣分為3個等級, 不同天氣下地面濕度不同， 觀察者只能觀察地面濕度無法直接觀察天氣。
假設地面濕度僅與當天天氣有關，每種天氣導致地面濕度的概率為:
| 天氣 |
晴 |
多雲 |
陰 |
雨 |
| 干燥 |
0.8 |
0.6 |
0.5 |
0.2 |
| 濕潤 |
0.2 |
0.4 |
0.5 |
0.8 |
為了便於表示，下文使用wet， dry代表濕潤和干燥兩種狀態。
假設每天天氣僅與前一天天氣有關， 之間的概率關系為：
-
若今天為晴，0.5的概率明天為晴, 0.5的概率明天多雲
-
若今天為多雲，0.4的概率明天為晴, 0.6的概率明天陰
-
若今天為陰，0.7的概率明天為雨, 0.3的概率明天為多雲
-
若今天為雨，0.2的概率明天為多雲， 0.3的概率明天為陰, 0.5的概率明天為雨
根據統計，一年中四種天氣占得的比例為(0.3, 0.1, 0.2, 0.4)
該示例中的隨機過程可以用HMM模型描述， 示例中存在兩個假設：
-
地面濕度僅與當天天氣有關: 觀測獨立性假設
-
每天的天氣僅與前一天有關: 馬爾科夫假設
從這個示例中可以看出， HMM模型的假設相當嚴格在使用時應慎重。
該示例中存在兩個序列:
根據天氣的統計數據得到初始狀態矩陣:
\[\pi = (0.3, 0.1, 0.2, 0.4) \]
根據盒子轉移的規則，寫出狀態轉移矩陣:
\[A = \left[ \begin{array}{ccc} 0.5 & 0.5 & 0 & 0 \\ 0.4 & 0 & 0.6 & 0 \\ 0 & 0.7 & 0 & 0.3 \\ 0 & 0.2 & 0.3 & 0.5 \\ \end{array} \right] \]
使用行向量\(pi_i\)表示某一天各天氣的概率, p * A即為下一天各天氣的概率:
\[pi_1 = pi * A = (0.19, 0.37, 0.18, 0.26) \]
根據每種天氣下地面濕度的概率，寫出觀測轉移矩陣:
\[B = \left[ \begin{array}{ccc} 0.8 & 0.2 \\ 0.6 & 0.4 \\ 0.5 & 0.5 \\ 0.2 & 0.8 \\ \end{array} \right] \]
使用行向量\(pi_i\)表示某一天各天氣的概率, p * B即為當天地面各濕度狀態的概率:
\[o = \pi * B = (0.48, 0.52) \]
常見資料中常以小寫英文字母表示矩陣的元素，本問題中下標僅表示時間
用A[i,j]表示矩陣A中第i行j列的元素， 行號與列號均從1開始
概率計算問題
HMM的定義給出了計算給定模型下出現特定觀測序列概率的方法，狀態序列概率:
\[P(I\mid\lambda)=P(i_1)P(i_2\mid i_1)\cdots P(i_T\mid i_{T-1})= \pi[i_1] * A[i_1,i_2] \cdots A[i_{T-1},i_T] \]
給定狀態序列下， 觀測序列概率:
\[P(O\mid I, \lambda)=P(o_1\mid i_1)\cdots P(o_T\mid i_T)=B[i_1,o_1]* B[i_2,o_2] \cdots B[i_T,(o_T] \]
上式中，\(o_i\)表示觀測序列中第i個元素。給定模型下觀測序列條件概率:
\[P(O\mid\lambda)=\sum_I P(O\mid I, \lambda)P(I\mid\lambda)=\sum_{\underbrace{i_1,i_2, \cdots ,i_T}_{N\times N\times \cdots \times N}}\pi[i_1] * A[i_1,i_2] * B[i_1,o_1] \cdots A[i_{T-1},i_T]*B[i_T,o_T] \]
但是該方法幾乎窮舉了所有可能的觀測序列, 時間復雜度高達\(O(TN^T)\)。
前向-后向算法通過部分概率和減少了計算量。
給定模型λ， 定義到時刻t出現部分觀測序列\(O_t\)且當前狀態為\(p_i\)的概率為前向概率， 記作：
\[\alpha_t(i) = P(O_{1,2, \cdots t},p_i | \lambda) \]
α在t=1時的初值為：
\[\alpha_1(i)=\pi_i * B[i,o_1],\quad i=1,2,\cdots,N \]
通過遞推得到α的后續值：
\[\alpha_{t+1}(i)= \sum_{j=1}^N(\alpha_t(j) * a[j,i]) * b[i,o_{t+1}] \quad i=1,2,\cdots,N \]
最終得到給定觀測序列概率P(O|λ):
\[P(O\mid \lambda)=\sum_{i=1}^N\alpha_T(i) \]
時間復雜度大為降低， 其原因在於α遞推過程中避免了計算不需要的部分概率和。
類似地可以定義后向概率:
\[\beta_t(i) = P(O_{t+1， \cdots T}, p_i | \lambda) \]
初值, 注意后向法是從后向前遞推的:
\[\beta_T(i)=1,\quad i=1,2,\cdots,N \]
得到給定觀測序列概率P(O|λ):
\[P(O\mid \lambda)=\sum_{i=1}^N\pi_i b[i,o_1]\beta_1(i) \]
給定HMM模型:
\[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] \]
\[\pi = (0.2, 0.4, 0.4) O = (1, 2, 1) \]
(1)求初值：
\[\alpha_1(1) = \pi_0[1] * B[1,1] = 0.2 * 0.5 = 0.1 \\ \alpha_1(2) = \pi_0[2] * B[2,1] = 0.4 * 0.4 = 0.16 \\ \alpha_1(3) = \pi_0[3] * B[3,1] = 0.4 * 0.7 = 0.28 \]
(2)遞推
\[\alpha_2(1) = \sum_{j=1}^3(\alpha_1(j) * A[j, 1] )* B[1,2] = (0.1 * 0.5 + 0.16 * 0.3 + 0.28 * 0.2 ) * 0.5 = 0.077 \\ \alpha_2(2) = \sum_{j=1}^3(\alpha_1(j) * A[j, 2] )* B[2,2] = (0.1 * 0.2 + 0.16 * 0.5 + 0.28 * 0.3 ) * 0.6 = 0.1104 \\ \alpha_2(3) = \sum_{j=1}^3(\alpha_1(j) * A[j, 3] )* B[3,2] = (0.1 * 0.3+ 0.16* 0.2 + 0.28 * 0 .5 ) * 0.3 = 0.0606 \]
\[\alpha_3(1) = \sum_{j=1}^4(\alpha_2(j) * A[j, 1] )* B[1,1] = 0.04187\\ \alpha_3(2) = \sum_{j=1}^4(\alpha_2(j) * A[j, 2] )* B[2,1] = 0.03551 \\ \alpha_3(3) = \sum_{j=1}^4(\alpha_2(j) * A[j, 3] )* B[3,1] = 0.05284 \\ \]
(3)求概率
\[P(O | \lambda) = \sum_{j=1}^3[\alpha_3(j)] = 0.13022 \]
參數估計問題
參數估計問題或稱為學習問題， 是指如何調整模型參數使得給定觀測序列出現的概率P(O|λ)最大的問題。根據訓練中是否給出了與觀測序列對應的狀態序列分為監督學習方法和非監督學習方法。
監督學習需要訓練數據中提供隨機過程中的狀態序列和觀測序列，非監督學習只需提供觀測序列。 在擁有狀態序列時可以直接使用訓練樣本中狀態轉移的頻率估計概率。
最大期望算法(Expectation-Maximization algorithm, EM)是一種極大似然參數估計法， 適用於依賴於隱藏變量的概率模型， 隱馬爾科夫模型是其典型的應用場景。
EM算法進行多次迭代對參數值進行估計，每次迭代分為計算期望和最大化兩步。
HMM的非監督學習方法即是著名的Baum-Welch算法， 它是EM算法在HMM中的特例。監督學習依賴人工標注的數據過多， 實際應用中我們通常選擇Baum-Welch算法進行訓練。
Baum-Welch算法只需要S個長度為T的觀測序列就可進行訓練。首先隱馬爾科夫模型是一個含有隱變量的概率模型:
\[P(O|\lambda) = \sum P(O|I,\lambda)P(I|\lambda) \]
確定模型全數據的對數似然函數\(log P(O,I|\lambda)\),根據EM算法寫出Q函數:
\[Q(\lambda|\lambda \overline) = \sum log P(O,I|\lambda)P(O,I|\lambda \overline) \]
其中\(\lambda*\)是參數的當前估計值， \(\lambda\)是要極大化的模型參數。即EM算法中的計算期望。
代入概率表達式：
\[P(O,I\mid\lambda)=P(O\mid I, \lambda)P(I\mid\lambda)=\pi_{i_1} * a[i_1,i_2] * b[i_1,o_1] \cdots a[i_{T-1},i_T] * b[i_T,o_T] \]
將Q改寫為:
\[Q(\lambda,\overline{\lambda})=\sum_Ilog\pi_{i_1}P(O,I|,\overline{\lambda})+\sum_I\left(\sum_{t=1}^{T-1}log\ a[i_t,i_{t+1}]\right)P(O,I|,\overline{\lambda})+\sum_I\left(\sum_{t=1}^{T}log\ b[i_t,o_t]\right)P(O,I|,\overline{\lambda}) \]
展開上式的第一項：
\[\sum_Ilog\pi_{i_1}P(O,I|,\overline{\lambda}) = \sum_{i=1}^Nlog\pi_{i}P(O,i_1=i|,\overline{\lambda}) \]
注意到\(\sum \pi_{i} = 1\), 對各項使用拉格朗日乘子法極大化， 求出最優參數:
\[\pi_i = \frac{P(O,i_1=i|,\overline{\lambda})}{P(O|\overline{\lambda})} = \gamma_1(i) \]
\[a_{ij} = \frac{\sum_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|,\overline{\lambda})}{\sum_{t=1}^{T-1}P(O,i_t=i|,\overline{\lambda})} = \frac{\sum_{t=1}^{T-1}\xi_t(i,j)}{\sum_{t=1}^{T-1}\gamma_t(j)} \]
\[b_{jk} = \frac{\sum_{t=1,o_t=v_k}^{T}P(O,i_t=j|,\overline{\lambda})}{\sum_{t=1}^{T}P(O,i_t=j|,\overline{\lambda})}= \frac{\sum_{t=1,o_t=v_k}^{T}\gamma_t(j)}{\sum_{t=1}^{T}\gamma_t(j)} \]
注意到上述表達式中出現了一些特殊的概率值, 在t時刻處於i狀態的概率:
\[\gamma_t(i)=P(i_t=q_i|O,\lambda)=\frac{P(i_t=q_i,O|\lambda)}{P(O|\lambda)}=\frac{\alpha_t(i)\beta_t(i)}{\sum_{j=1}^N\alpha_t(j)\beta_t(j)} \]
在時刻t處於i狀態， 時刻t+1處於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)}=\frac{\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)}{\sum_{i=1}^N\sum_{j=1}^N\alpha_t(i) * a[i,j] * b[j,o_{t+1}] * \beta_{t+1}(j)} \]
將它們代入參數估計式中， 即可得到最終表達式。
解碼問題
解碼問題是指，給定模型參數 求使給定觀測序列出現概率最大的狀態序列的問題，即最大化P(I | O, λ)的問題。很容易想到的估計方法是選擇使得各時刻系統處於i狀態的概率最大的狀態作為該時刻狀態:
\[I = arg\quad max\quad \gamma_t(i) \]
然而這種方法容易陷入局部最優解， 更有效的方法是采用Viterbi算法估計全局最優解。
每個狀態序列都對應着一條狀態轉移路徑，最可能的狀態序列對應的路徑出現概率也最高。如果我們把出現的概率定義為路徑的代價話， 我們就可以使用最優路徑算法求得出現概率最高的專業路徑。
Viterbi算法即采用動態規划算法求得這條最優路徑(概率最大路徑)。引入變量\(\delta\)表示在時刻t狀態為i的所有路徑中最大概率值:
\[\delta_t(i)=\max_{i_1,i_2,...,i_{t-1}}P(i_t=i,i_{t-1},...,i_1,o_t,...,o_1|\lambda) \]
定義最優路徑中t-1個節點為:
\[\psi_t(i) = arg \quad max \quad \big[ \delta_{t-1}(j) a[j,i] \big] \]
初始化:
\[\delta_1(i) = \pi_ib[i,o_1] \]
\[\psi_1(i) = 0 \]
遞推：
\[\delta_t(i) = max(\delta_{t-1}(j)a[j,i]) * b[i,o_t] \]
\[\psi_t(i) = arg \quad max(\delta_{t-1}(j)*a[j,i]) \]
回溯:
\[i_t^* = \psi_{t+1}(i_{t+1}^*) \]
求得的\(i_t\)序列即為所求狀態序列
