【机器学习】隐马尔可夫模型及其三个基本问题(四)状态序列预测算法及python实现
- 一、维特比算法
- 二、python实现
- 参考资料
隐马尔可夫模型状态序列预测问题是指给定模型 λ = [ A , B , ∏ ] \lambda = \left[ {A,B,\prod } \right] λ = [ A , B , ∏ ] 和观测序列 X = { x 1 , x 2 , ⋯   , x n } X = \left\{ {{x_1},{x_2}, \cdots ,{x_n}} \right\} X = { x 1 , x 2 , ⋯ , x n } ,求最可能出现的对应状态序列。本篇博文介绍状态序列预测算法中的维特比算法(Viterbi algorithm)[参考资料1]。
一、维特比算法
维特比算法实际是用动态规划解隐马尔科夫模型的预测问题,首先导入两个变量 δ \delta δ 和 ψ \psi ψ 。
定义在时刻
t t
t
状态为
q i q_i
q
i
的所有单个路径
( y 1 , y 2 , ⋯   , y t ) ({y_1},{y_2}, \cdots ,{y_t})
(
y
1
,
y
2
,
⋯
,
y
t
)
中概率最大值为:
δ t ( q i ) = max y 1 , y 2 , ⋯   , y t − 1 P ( y t = q i , y t − 1 , ⋯   , y 1 , x t , ⋯   , x 1 ∣ λ ) , i = 1 , 2 , ⋯   , N {\delta _t}\left( {{q_i}} \right) = \mathop {\max }\limits_{{y_1},{y_2}, \cdots ,{y_{t - 1}}} P\left( {{y_t} = {q_i},{y_{t - 1}}, \cdots ,{y_1},{x_t}, \cdots ,{x_1}\left| \lambda \right.} \right),i = 1,2, \cdots ,N
δ
t
(
q
i
)
=
y
1
,
y
2
,
⋯
,
y
t
−
1
max
P
(
y
t
=
q
i
,
y
t
−
1
,
⋯
,
y
1
,
x
t
,
⋯
,
x
1
∣
λ
)
,
i
=
1
,
2
,
⋯
,
N
则变量
δ \delta
δ
的递推公式为:
δ t + 1 ( q i ) = max y 1 , y 2 , ⋯   , y t P ( y t + 1 = q i , y t , ⋯   , y 1 , x t + 1 , ⋯   , x 1 ∣ λ ) = max 1 ≤ j ≤ N [ δ t ( j ) a j i ] b i ( x t + 1 ) , i = 1 , 2 , ⋯   , N ; t = 1 , 2 , ⋯   , n − 1 {\delta _{t + 1}}\left( {{q_i}} \right) = \mathop {\max }\limits_{{y_1},{y_2}, \cdots ,{y_t}} P\left( {{y_{t + 1}} = {q_i},{y_t}, \cdots ,{y_1},{x_{t + 1}}, \cdots ,{x_1}\left| \lambda \right.} \right) = \mathop {\max }\limits_{1 \le j \le N} \left[ {{\delta _t}\left( j \right){a_{ji}}} \right]{b_i}({x_{t + 1}}),i = 1,2, \cdots ,N;t = 1,2, \cdots ,n - 1
δ
t
+
1
(
q
i
)
=
y
1
,
y
2
,
⋯
,
y
t
max
P
(
y
t
+
1
=
q
i
,
y
t
,
⋯
,
y
1
,
x
t
+
1
,
⋯
,
x
1
∣
λ
)
=
1
≤
j
≤
N
max
[
δ
t
(
j
)
a
j
i
]
b
i
(
x
t
+
1
)
,
i
=
1
,
2
,
⋯
,
N
;
t
=
1
,
2
,
⋯
,
n
−
1
定义在时刻
t t
t
状态为
q i q_i
q
i
的所有单个路径
( y 1 , y 2 , ⋯   , y t − 1 , q i ) ({y_1},{y_2}, \cdots ,{y_{t - 1}},q_i)
(
y
1
,
y
2
,
⋯
,
y
t
−
1
,
q
i
)
中概率最大的路径的第
t − 1 t-1
t
−
1
个节点为:
ψ t ( q i ) = arg max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] , i = 1 , 2 , ⋯   , N {\psi _t}(q_i) = \arg \mathop {\max }\limits_{1 \le j \le N} \left[ {{\delta _{t - 1}}\left( j \right){a_{ji}}} \right],i = 1,2, \cdots ,N
ψ
t
(
q
i
)
=
ar
g
1
≤
j
≤
N
max
[
δ
t
−
1
(
j
)
a
j
i
]
,
i
=
1
,
2
,
⋯
,
N
维特比算法步骤
:
输入:模型
λ = [ A , B , ∏ ] \lambda = [A,B,\prod ]
λ
=
[
A
,
B
,
∏
]
和观测序列
X = { x 1 , x 2 , ⋯   , x n } X = \left\{ {{x_1},{x_2}, \cdots ,{x_n}} \right\}
X
=
{
x
1
,
x
2
,
⋯
,
x
n
}
输出:最优状态序列
Y = { y 1 , y 2 , ⋯   , y n } Y = \left\{ {{y_1},{y_2}, \cdots ,{y_n}} \right\}
Y
=
{
y
1
,
y
2
,
⋯
,
y
n
}
(1)初始化:
δ 1 ( q i ) = π i b i ( x 1 ) , i = 1 , 2 , ⋯   , N ψ 1 ( q i ) = 0 , i = 1 , 2 , ⋯   , N \begin{array}{l} {\delta _1}({q_i}) = {\pi _{_i}}{b_i}\left( {{x_1}} \right),i = 1,2, \cdots ,N\\ {\psi _1}\left( {{q_i}} \right) = 0,i = 1,2, \cdots ,N \end{array} δ 1 ( q i ) = π i b i ( x 1 ) , i = 1 , 2 , ⋯ , N ψ 1 ( q i ) = 0 , i = 1 , 2 , ⋯ , N
(2)递推:对于 t = 2 , 3 , ⋯   , n t = 2,3, \cdots ,n t = 2 , 3 , ⋯ , n
δ t ( q i ) = max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] b i ( x t ) , i = 1 , 2 , ⋯   , N ψ t ( q i ) = arg max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] , i = 1 , 2 , ⋯   , N \begin{array}{l} {\delta _t}\left( {{q_i}} \right) = \mathop {\max }\limits_{1 \le j \le N} \left[ {{\delta _{t - 1}}\left( j \right){a_{ji}}} \right]{b_i}({x_t}),i = 1,2, \cdots ,N\\ {\psi _t}\left( {{q_i}} \right) = \arg \mathop {\max }\limits_{1 \le j \le N} \left[ {{\delta _{t - 1}}\left( j \right){a_{ji}}} \right],i = 1,2, \cdots ,N \end{array} δ t ( q i ) = 1 ≤ j ≤ N max [ δ t − 1 ( j ) a j i ] b i ( x t ) , i = 1 , 2 , ⋯ , N ψ t ( q i ) = ar g 1 ≤ j ≤ N max [ δ t − 1 ( j ) a j i ] , i = 1 , 2 , ⋯ , N
(3)终止:
y n = arg max 1 ≤ i ≤ N [ δ n ( q i ) ] {y_n} = \arg \mathop {\max }\limits_{1 \le i \le N} \left[ {{\delta _n}\left( {{q_i}} \right)} \right] y n = ar g 1 ≤ i ≤ N max [ δ n ( q i ) ]
(4)最优路径回朔:对于 t = n − 1 , n − 2 , ⋯   , 1 t = n - 1,n - 2, \cdots ,1 t = n − 1 , n − 2 , ⋯ , 1
y t = ψ t + 1 ( y t + 1 ) {y_t} = {\psi _{t + 1}}\left( {{y_{t + 1}}} \right) y t = ψ t + 1 ( y t + 1 )
二、python实现
代码参考资料2
import
numpy as np
class
HMM
:
def
__init__
(
self
,
A
,
B
,
Pi
,
O
)
:
self
.
A
=
np
.
array
(
A
)
# 状态转移概率矩阵
self
.
B
=
np
.
array
(
B
)
# 观测概率矩阵
self
.
Pi
=
np
.
array
(
Pi
)
# 初始状态概率矩阵
self
.
O
=
np
.
array
(
O
)
# 观测序列
self
.
N
=
self
.
A
.
shape
[
0
]
# 状态取值个数
self
.
M
=
self
.
B
.
shape
[
1
]
# 观测值取值个数
def
viterbi
(
self
)
:
n
=
len
(
self
.
O
)
#观测序列和状态序列的长度
Y
=
np
.
zeros
(
n
)
#初始化状态序列
delta
=
np
.
zeros
(
(
n
,
self
.
N
)
)
psi
=
np
.
zeros
(
(
n
,
self
.
N
)
)
#初始化
for
i in
range
(
self
.
N
)
:
delta
[
0
,
i
]
=
self
.
Pi
[
i
]
*
self
.
B
[
i
,
self
.
O
[
0
]
]
psi
[
0
,
i
]
=
0
#递推
for
t in
range
(
1
,
n
)
:
for
i in
range
(
self
.
N
)
:
delta
[
t
,
i
]
=
self
.
B
[
i
,
self
.
O
[
t
]
]
*
np
.
array
(
[
delta
[
t
-
1
,
j
]
*
self
.
A
[
j
,
i
]
for
j in
range
(
self
.
N
)
]
)
.
max
(
)
psi
[
t
,
i
]
=
np
.
array
(
[
delta
[
t
-
1
,
j
]
*
self
.
A
[
j
,
i
]
for
j in
range
(
self
.
N
)
]
)
.
argmax
(
)
print
(
psi
)
#终止
Y
[
n
-
1
]
=
delta
[
n
-
1
,
:
]
.
argmax
(
)
print
(
Y
[
n
-
1
]
)
#回朔
for
t in
range
(
n
-
2
,
-
1
,
-
1
)
:
Y
[
t
]
=
psi
[
int
(
t
+
1
)
,
int
(
Y
[
t
+
1
]
)
]
return
Y
if
__name__
==
'__main__'
:
print
(
'----------------1.初始化HMM模型参数------------------'
)
A
=
[
[
0
,
1
,
0
,
0
]
,
[
0.4
,
0
,
0.6
,
0
]
,
[
0
,
0.4
,
0
,
0.6
]
,
[
0
,
0
,
0.5
,
0.5
]
]
B
=
[
[
0.5
,
0.5
]
,
[
0.3
,
0.7
]
,
[
0.6
,
0.4
]
,
[
0.8
,
0.2
]
]
Pi
=
[
0.25
,
0.25
,
0.25
,
0.25
]
print
(
'----------------2.观测序列---------------------------'
)
O
=
[
0
,
1
,
0
]
print
(
'----------------3.viterbi 算法-----------------------'
)
hmm
=
HMM
(
A
,
B
,
Pi
,
O
)
Y
=
hmm
.
viterbi
(
)
参考资料
1、《统计学习方法》李航
2、https://www.cnblogs.com/d-roger/articles/5719979.html