隱馬爾可夫模型維特比算法詳解
關於隱馬爾可夫模型的維特比解碼算法網上已有一大批文章介紹，故本文不再介紹。
本文主要是在讀《自然語言處理簡明教程》和看HanLP 中文人名識別源碼過程中，對該算法的一次梳理，以防忘記。
隱馬模型有三個問題，其中二個是：
- 給定HMM模型 \(\lambda\) 和一個觀察序列O，確定觀察序列O出現的可能性\(P(O|\lambda)\)
- 給定HMM模型\(\lambda\) 和一個觀察序列O，確定產生O的最可能的隱藏序列Q
第一個問題用向前算法解決，可參考：隱馬爾科夫模型HMM（二）前向后向算法評估觀察序列概率
第二個問題用維特比算法解決，可參考：隱馬爾科夫模型HMM（四）維特比算法解碼隱藏狀態序列
下面記錄下我對這兩個算法的理解。
向前算法
將\(P(O|\lambda)\)簡寫成\(P(O)\)，觀察狀態是由隱藏狀態生成的，因此：任何一個可能的隱藏狀態\(Q=(q_1,q_2,...q_N)\)以一定的概率生成觀察狀態O。故：
根據貝葉斯公式：\(P(O,Q)=P(O|Q)*P(Q)\)
所以：\(P(O)=\Sigma_{Q}{P(O|Q)*P(Q)}\)
對於一個長度為T的觀察序列，書中第518頁指出：一共有\(N^T\)種可能的隱藏狀態。將這\(N^T\)個隱藏狀態生成序列O的概率求和，就得到了\(P(O)\)，但這種方法的時間復雜度是指數級的：\(O(N^T)\)
根據隱馬爾可夫模型的獨立輸出假設：
根據概率論中的鏈式法則：
再根據隱馬爾可夫模型的一階馬爾可夫鏈假設：
簡化一下，就是：
因此：對於一個特定的隱藏狀態\(Q^*\)
這里：\(p(q_i|q_{i-1})\)就是從隱藏狀態\(q_{i-1}\)轉移到隱藏狀態\(q_i\)的轉移概率；
\(p(o_i|q_i)\)就是隱藏狀態\(q_i\)生成觀察狀態\(o_i\) 的 發射概率。
正是根據轉移概率和發射概率 定義：\(\alpha_t(j)\)，從而采用動態規划的方法來求解：\(P(O|\lambda)\)
動態規划方法
定義：\(\alpha_t(j)=P(o_1,o_2,...o_t, q_t=j|\lambda)\)生成了\(t\)個觀察后，在時刻\(t\)隱藏狀態\(q\)取值為\(j\)的概率，用公式表示成：
這相當於動態規則的狀態方程。上面公式中有個\(\Sigma\)，在書中第586頁的圖9.8 解釋了 求和符號的意義：
- 在\(t-1\)時刻的每一個隱藏狀態\(q(i)\)的取值概率\(\alpha_{t-1}(i)\) 乘以 從狀態i 轉移 到 狀態j 的概率 再乘以 \(t\) 時刻狀態j 生成 觀察狀態\(o_t\)的概率
- 對上行中所說的： 每一個求和
對於動態規則，有兩個性質：最優子結構和重疊子問題。
\(\alpha_{t-1}(j)\)是\(\alpha_{t}(j)\)的子結構，\(\alpha_{t-2}(j)\)是\(\alpha_{t-1}(j)\)的子結構，因為問題的規模變小了。
重疊子問題：要求解\(\alpha_{t}(j)\)，需要求解\(\alpha_{t-1}(j)\)，要求解 \(\alpha_{t-1}(j)\)需要求解\(\alpha_{t-2}(j)\)……
那么\(\alpha_{t-2}(j)\)就是求解 \(\alpha_{t}(j)\) 和 \(\alpha_{t-1}(j)\) 的一個重疊子問題。
如果在自底向上求解過程中，把這些子問題記錄下來：將\(\alpha_{1}(j)\)、\(\alpha_{2}(j)\)……都保存起來，當使用到它們時，直接“查表”，那么計算起來會快很多。關於這種思想，可參考：動態規划之Fib數列類問題應用
書中有向前算法的詳細示例。這里不再介紹。
維特比算法
這里的維特比算法和上面的向前算法其實是非常相似的。向前算法是對 \(t-1\)時刻中的每個[隱藏狀態的概率 乘 轉移概率 乘 發射概率] 求和；
而維特比算法則是：根據 向前算法 \(t-1\)時的求和得到 \(t\)時刻的某個隱藏狀態概率，而 \(t\) 時刻 一共有N個隱藏狀態概率，算法選出這 N個隱藏狀態中 概率值 最大的那個 隱藏狀態。
而求解\(t\)時刻N個隱藏狀態概率的最大值有兩種方法：一種方法是暴力法，第587頁已經講了。另一種是動態規划方法，下面記錄一下動態規划方法的求解思路。
動態規划
\(v_t(j)\)表示：在 \(t-1\)時刻的 每個隱藏狀態(共有N個) 乘以 轉移概率 乘以 發射概率 得到一個結果，取這N個結果的最大值 （這也是最終求得的結果是最優的原因---不是貪心思路）作為 t 時刻 \(j\) 狀態的概率。由於一共有N個隱藏狀態，在\(t\)時刻，需要求解：\(v_t(1)\)、\(v_t(2)\)……\(v_t(N)\)
寫出動態規划的狀態方程如下：
\(t\)代表時間，范圍為1,2...T，\(a_{ij}\)代表隱藏狀態轉移矩陣的概率，即從隱藏狀態\(q_i\)轉移到\(q_j\)的概率。\(b_j(o_t)\)代表發射概率，即\(t\)時刻的隱藏狀態\(q_j\)生成觀察狀態\(o_t\)的概率。 具體詳細的示例解釋參考書上第587頁開始的講解。
下面用代碼驗證一下理論的正確性：
參考：通用維特比算法實現並針對《自然語言處理簡明教程》第9章隱馬爾可夫模型介紹，驗證了在觀察狀態 3 1 3 時最佳隱藏狀態為 H H H。具體驗證代碼如下，並加了一些注釋。
import static org.hapjin.hanlp.Viterbi.Activity.one;
import static org.hapjin.hanlp.Viterbi.Activity.two;
import static org.hapjin.hanlp.Viterbi.Activity.three;
import static org.hapjin.hanlp.Viterbi.Weather.hot;
import static org.hapjin.hanlp.Viterbi.Weather.cold;
public class Viterbi {
static enum Weather
{
cold,
hot,
}
static enum Activity
{
one,
two,
three,
}
static int[] states = new int[]{cold.ordinal(), hot.ordinal()};
// static int[] observations = new int[]{one.ordinal(), two.ordinal(),three.ordinal()};
static int[] observations = new int[]{three.ordinal(), one.ordinal(),three.ordinal()};
static double[] start_probability = new double[]{0.2, 0.8};
static double[][] transititon_probability = new double[][]{
{0.6, 0.4},//cold
{0.3, 0.7},//hot
};
static double[][] emission_probability = new double[][]{
{0.5, 0.4, 0.1},//cold
{0.2, 0.4, 0.4},//hot
};
public static void main(String[] args)
{
int[] result = Viterbi.compute(observations, states, start_probability, transititon_probability, emission_probability);
for (int r : result)
{
System.out.print(Weather.values()[r] + " ");
}
System.out.println();
}
public static int[] compute(int[] obs, int[] states, double[] start_p, double[][] trans_p, double[][] emit_p)
{
//動態規划中保存 當前最優結果, 供后續計算 直接 "查表"
double[][] V = new double[obs.length][states.length];//v_t(j)
//保存最優路徑
int[][] path = new int[states.length][obs.length];//[state][t]
for (int y : states)
{
V[0][y] = start_p[y] * emit_p[y][obs[0]];//t=0 (t=0代表起始隱藏狀態)
path[y][0] = y;
}
//時間復雜度: (T-1)*N*N=O(T*N^2)
// T-1
for (int t = 1; t < obs.length; ++t)
{
int[][] newpath = new int[states.length][obs.length];//應該是可以優化一下的.
//N 個隱藏狀態 即:{cold, hot}, N=2
for (int y : states)
{
double prob = -1;
int state;
//N
for (int y0 : states)
{
// v_{t-1}(i)*a_{ij}*b_j(o_t)
double nprob = V[t - 1][y0] * trans_p[y0][y] * emit_p[y][obs[t]];
//find max
if (nprob > prob)
{
prob = nprob;
state = y0;
// 記錄t時刻 隱藏狀態為y 時的最大概率
V[t][y] = prob;//t是第一個for循環參數, y 是第二個for循環參數
// 記錄路徑
System.arraycopy(path[state], 0, newpath[y], 0, t);//
newpath[y][t] = y;//將t時刻 最佳隱藏狀態 y 保存
}
}
}
path = newpath;
}//end outer for
double prob = -1;
int state = 0;
//找出最后那個時刻的 V_T(j) j=1,2...N 的最大值 對應的隱藏狀態y
for (int y : states)
{
if (V[obs.length - 1][y] > prob)
{
prob = V[obs.length - 1][y];
state = y;
}
}
return path[state];//根據上面 max{V_T(j)} 求得的y "回溯" 得到 最優路徑
}
}
《自然語言處理簡明教程》里面詳細介紹了HMM三個問題的求解過程，通俗易懂。
另外想學習一下概率圖模型，不知道有沒有好的書籍推薦？
參考鏈接：HanLP中人名識別分析
原文鏈接：http://www.cnblogs.com/hapjin/p/9033471.html