1前言

這裡介紹維特比算法，主要是其在解決HMM模型中預測問題中起到了很大的作用，之前也粗略介紹過維特比算法，但是不是很詳細，這裡再詳細介紹一下。HMM預測問題也稱為解碼(decoding)問題。已知模型

和觀測序列

，求給定的觀測序列條件概率

最大的狀態序

即給定觀測序列，求最有可能的對應的狀態序列。對於該問題，有兩種算法：近似算法與維特比算法（Viterbi algorithm），我們主要介紹維特比算法。

維特比算法實際是用動態規劃解隱馬爾科夫模型預測問題，即用動態規劃(dynamic programming)求概率最大路徑（最優路徑）。這時一條路徑對應著一個狀態序列。

2 維特比算法

輸入：模型

和觀測序列

;

輸出：最優路徑

。

(1) 初始化：

(2)遞推.對t=2,3,…,T

需要注意的是

面向的是t-1時刻得到當前的轉移率，並沒有與狀態概率相乘。

(3) 終止

例：HMM模型

，題目再述：有三個盒子，每個盒子中有紅、白兩種球，其中專業概率相關參數如下：。

已知觀測序列

，試求最優狀態序列，即最優路徑

.

解：如下圖所示，要在所有可能的路徑中選擇一條最優路徑，求狀態i觀測O1為紅的概率，記此概率為

，則

求最優路徑

(1)初始化.在t=1時，對每一個狀態i，i=1,2,3,求狀態為i觀測O1為紅球的概率，記此概率為

，上例已知：

則：

帶入實際數據：

記

.

(2)在t=2時，對每個狀態i, i=1,2,3,求在t=1時狀態j觀測為紅並在t=2時狀態為i觀測O2為白的路徑的最大概率，記此最大概率為

，則：

同時，對每個狀態i,i=1,2,3,記錄概率最大路徑的的前一個狀態j:

計算有：

需要注意的是：在t=2也有三個狀態的可能性。其中

表示從t=1的各個可能的狀態到t=2 狀態為1觀測為白的最大值。

是因為，當前路徑中，上一個狀態中狀態3的概率最大，這裡記住上一個狀態，以便回溯。這裡計算

不是很清楚，以

計算為例，如下：

可得當i=3時取最大，即有

同樣，在t=3時，

(3)以P*表示最優路徑的概率，則：

最優路徑的終點是：

(4)由最優路徑的終

,逆向查找,

:

在t=2時，

在t=1時，

於是求得最優路徑，即最優狀態序列

3 Python 實現

以下代碼是個人根據李航老師那本書進行書寫的，也難免有些bug，如果有的話，也希望各位友友提出，共同學學習和進步。

 1def viterbi(A, B, Pi, Obser, state):
 2 """
 3 計算預測狀態
 4 :para:A 狀態轉移矩陣
 5 :para:B 發射矩陣
 6 :para:Pi 初始化矩陣
 7 :para:Obser 觀測序列
 8 :parar:state 狀態集合
 9 :return: 返回兩個值，第一個值是整個過程的維特比計算矩陣，第二個是預測序列的索引
10 """
11 import numpy as np
12 row, col = len(Obser), len(state)
13 res = np.zeros((row, col)) # 豎向矩陣
14 res2 = np.zeros_like(res)
15 # print(res2)
16 # 轉換為矩陣計算
17 A, B, Pi = np.array(A), np.array(B), np.array(Pi)
18 # 初始化
19 res[0, :] = B.T[0]*Pi
20 # 後續循環狀態(2-t狀態)
21 for i in range(1, row):
22 # 循環隱藏狀態數，計算當前狀態每個隱藏狀態的概率
23 ob = Obser[i] # 當前觀察值
24 tempres, tempres2 = [], []
25 for j in range(col):
26 # 以盒子1為例， 其他盒子轉移到盒子 

27 # print(A[:, j]) # 表示A中的第j列數據, 即其他盒子轉移到盒子j的概率
28 # print(res[i - 1]) # res中第i-1行數的值 即delta(i-1)
29 # print(B[:, ob]) # 發射矩陣中的第ob列（由觀測值確定）
30 # delta j的計算
31 delta = A[:, j]*res[i - 1]*B[j][ob]
32 # Psi # 獲取最大值的索引
33 tempres2.append(np.argmax(A[:, j]*res[i - 1]))
34 tempres.append(np.max(delta))
35 res[i, :] = np.array(tempres) # 結果矩陣賦值
36 res2[i, :] = np.array(tempres2)
37 # 通過res和res2回溯
38 result = []
39 # 最後一行直接計算
40 result.append(np.argmax(res[row-1, :]))
41 i = row - 1
42 while i > 0:
43 result.append(res2[i][np.argmax(res[i, :])])
44 i -= 1
45 result.reverse() # 我們是逆向添加的
46 return res, result
47
48
49if __name__ == "__main__":
50 # 隱藏狀態, 為方便計算這裡將隱層
51 invisiable = {0: '盒子1', 1: '盒子2', 2: '盒子3'}
52 invisiable_ls = [0, 1, 2]
53 # 初始狀態 pi
54 pi = [0.2, 0.4, 0.4]
55 # 轉移矩陣 A
56 trainsion_probility = [
57 [0.5, 0.2, 0.3],
58 [0.3, 0.5, 0.2],
59 [0.2, 0.3, 0.5]
60 ]
61 # 發射矩陣B
62 emission_probility = [
63 [0.5, 0.5],
64 [0.4, 0.6],
65 [0.7, 0.3] 

66 ]
67 # 觀測序列
68 obs_dic = {0: "紅", 1: "白"}
69 # obs_seq = [0, 1, 0, 1, 0, 0, 0, 1, 1, 1]
70 obs_seq = [0, 1, 0]
71 print("觀測序列為：")
72 for i in obs_seq:
73 print(obs_dic[i], end=" ")
74 print("")
75 # 結果
76 res, result = viterbi(trainsion_probility, emission_probility, pi, obs_seq, invisiable_ls)
77 print("res:\\n", res)
78 print("預測序列為：")
79 for i in result:
80 print(invisiable[i], end=" ")

輸出結果：

1觀測序列為：
2紅 白 紅
3res:
4 [[0.1 0.16 0.28 ]
5 [0.028 0.0504 0.042 ]
6 [0.00756 0.01008 0.0147 ]]
7預測序列為：
8盒子3 盒子3 盒子3

Reference

李航的《統計機器學習》

4寫在後面

在Typora寫完公式，又去寫代碼，頭暈腦脹，這個微信還不支持公式，這個排版也是醉醉的，對各個感興趣的可以看看原文！寫作不易，如果對您有用的話，點贊轉發吧！

點擊“瞭解更多”查看原文，獲得更佳的閱讀體驗。