1. 隐马尔科夫模型的定义
隐马尔可夫模型是关于时序的概率模型，描述一个有隐藏的马尔科夫链随机生成不可观测的状态随机序列，再由各个状态生成一个观测从而产生观测随机序列的过程。
1.1 组成
- 所有可能状态的集合：
Q = { q 1 , q 2 , … … q N } ( N 是 可 能 的 状 态 的 集 合 ) Q=\{ q_1,q_2,……q_N \} (N是可能的状态的集合) Q={q1,q2,……qN}(N是可能的状态的集合) - 所有可能观的集合
V = { v 1 , v 2 , … … v N } ( M 是 可 能 的 观 测 数 ) V = \{ v_1,v_2,……v_N \} (M是可能的观测数) V={v1,v2,……vN}(M是可能的观测数) - I: 长度为T的状态序列
- O: 对应的观测序列
I = ( i 1 , i 2 , … … , i T ) ; O = ( o 1 , o 2 , … … ， o T ) I=(i_1,i_2,……,i_T);O=(o_1,o_2,……，o_T) I=(i1,i2,……,iT);O=(o1,o2,……，oT)
- 状态转移概率分布
A : 状 态 转 移 矩 阵 ： A = [ a i j ] N ∗ N a i j = P ( i T + 1 = q j ∣ i T = q i ) , i = 1 , 2 , … … N ； j = 1 , 2 ， … … N A:状态转移矩阵：A=[a_{ij}]_{N*N} \\ a_ij = P(i_{T+1} = q_j | i_T = q_i ), i=1,2,……N； j= 1,2，……N A:状态转移矩阵：A=[aij]N∗Naij=P(iT+1=qj∣iT=qi),i=1,2,……N；j=1,2，……N
t时刻处于咋混给他qi的条件下在时刻他+转移到状态qj的概率
- 观测概率分布
B = [ b j ( k ) ] N ∗ N b j ( k ) = P ( o t = v k ∣ i t = q j ) , k = 1 , 2 , … … M ; j = 1 , 2 , … … N B= [b_j(k)]_{N*N} \\ b_j(k) = P(o_t = v_k| i_t=q_j),k=1,2,……M;j=1,2,……N B=[bj(k)]N∗Nbj(k)=P(ot=vk∣it=qj),k=1,2,……M;j=1,2,……N
在t时刻处于状态qj的条件下生成观测vk的概率
- 初始概率分布：
π = ( π i ) π i = P ( i 1 = q i ) , i = 1 , 2 , … … N \pi = (\pi_i) \\ \pi_i = P(i_1 = q_i),i=1,2,……N π=(πi)πi=P(i1=qi),i=1,2,……N
t=1时刻处于状态qi的状态
隐马尔可夫模型的三元符号表示：
λ
=
(
A
,
B
,
π
)
\lambda = (A,B,\pi)
λ=(A,B,π)
隐马尔可夫的两个假设条件：
- 隐马尔可夫链t时刻的状态只与t-1时刻的状态有关，与其他时刻的状态及观测无关；
- 任意时刻的观测只依赖于该时刻的马尔可夫链的状态，与其他观测及状态无关。
2.隐马尔科夫要解决的问题：
- 识别问题（在给定参数下计算样本出现的概率）：
给 定 λ = ( A , B , π ) ， 计 算 其 产 生 观 测 序 列 O = [ o 1 , o 2 … … o n ] 的 概 率 P ( O ∣ λ ) 给定\lambda = (A,B,\pi) ，计算其产生观测序列O= [o_1,o_2……o_n]的概率P(O|\lambda) 给定λ=(A,B,π)，计算其产生观测序列O=[o1,o2……on]的概率P(O∣λ) - 学习问题（由观测样本得到模型参数）：
给 定 观 测 序 列 O = [ o 1 , o 2 … … o n ] ， 调 整 模 型 λ = ( A , B , π ) 参 数 ， 使 得 该 序 列 出 现 的 概 率 P ( O ∣ λ ) 最 大 给定 观测序列O= [o_1,o_2……o_n] ，调整模型\lambda = (A,B,\pi) 参数， 使得该序列出现的概率P(O|\lambda)最大 给定观测序列O=[o1,o2……on]，调整模型λ=(A,B,π)参数，使得该序列出现的概率P(O∣λ)最大 - 解码问题（由观测样本得到隐状态）：
给 定 观 测 序 列 O = [ o 1 , o 2 … … o n ] 和 模 型 λ = ( A , B , π ) ， 找 到 与 此 观 测 序 列 最 匹 配 的 状 态 序 列 给定 观测序列O= [o_1,o_2……o_n] 和模型\lambda = (A,B,\pi) ，找到与此观测序列最匹配的状态序列 给定观测序列O=[o1,o2……on]和模型λ=(A,B,π)，找到与此观测序列最匹配的状态序列
3.隐马尔科夫模型的例子
盒子和球模型： 假设有4个盒子，每个盒子里都装有红、白两种颜色的球，盒子里的红、白球数由下表给出：
| 盒子编号 | 1 | 2 | 3 | 4 |
| 红球数 | 5 | 3 | 6 | 4 |
| 白球数 | 5 | 7 | 4 | 2 |
按照下列方法抽球，产生一个球的颜色的观测序列
- 以等概率随机取一个盒子，再从这个盒子随机取一个球，记录其颜色，然后放回；
- 然后从当前盒子随机转移到下一个盒子，规则是：如果当前盒子是盒子1，那么下一盒子一定是盒子2；如果当前是盒子2或3，那么分别以概率0.4和0.6转移到左边或右边的盒子；如果当前是盒子4，那么各以0.5的概率停留在盒子4或转移到盒子3；
- 确定转移的盒子后，再从这个盒子里随机抽出1个球，记录其颜色，放回；
- 如此重复5次，得到一个球的颜色的观测序列： O ={红，红，白，白 ，红 }
在这个过程中，观察者只能观测到球的颜色的序列，观测不到球是从哪个盒子取出的，即观测不到盒子的序列。在这个例子中有两个随机序列，一个是盒子的序列（状态序列），一个是球的颜色的观测序列（观测序列）。前者是隐藏的，只有后者是可观测的。这是一个隐马尔可夫模型的例子。根据所给条件，可以明确状态集合、观测集合、序列长度以及模型的三要素。
盒子对应状态，状态集合是：
Q
=
盒
子
1
,
盒
子
2
,
盒
子
3
,
盒
子
4
，
N
=
4
Q= {盒子1,盒子2,盒子3,盒子4}，N=4
Q=盒子1,盒子2,盒子3,盒子4，N=4
球的颜色对应观测，观测的集合是：
V
=
红
，
白
，
M
=
2
V= {红，白}，M=2
V=红，白，M=2
状态序列和观测序列长度T= 5
初始概率分布为：
π
=
(
0.25
,
0.25
,
0.25
,
0.25
)
T
\pi = (0.25,0.25,0.25,0.25)^T
π=(0.25,0.25,0.25,0.25)T
状态转移概率分布为：
[
0
1
0
0
0.4
0
0.6
0
0
0.4
0
0.6
0
0
0.5
0.5
]
\begin{bmatrix} 0&1&0&0 \\ 0.4&0&0.6&0 \\ 0&0.4&0&0.6 \\ 0&0&0.5&0.5 \\ \end{bmatrix}
⎣⎢⎢⎡00.400100.4000.600.5000.60.5⎦⎥⎥⎤
观测概率分布为：
[
0.5
0.5
0.3
0.7
0.6
0.4
0.8
0.2
]
\begin{bmatrix} 0.5&0.5 \\ 0.3&0.7 \\ 0.6&0.4 \\ 0.8&0.2\\ \end{bmatrix}
⎣⎢⎢⎡0.50.30.60.80.50.70.40.2⎦⎥⎥⎤
4.代码实现
4.1 隐马尔可夫数据生成
import breeze.linalg.{DenseMatrix, DenseVector, sum}
import scala.collection.mutable.ListBuffer
/**
*
* @param pi 隐状态初始概率分布
* @param stateTransitionMatrix 状态转移矩阵
* @param confusionMatrix 观测状态生成矩阵
*/
case class hmm(
pi: DenseVector[Double],
stateTransitionMatrix: DenseMatrix[Double],
confusionMatrix: DenseMatrix[Double]
) {
// 状态集合个数
val n: Int = stateTransitionMatrix.cols
// 根据给定的概率分布随机返回数据
def getDistData(dist: DenseVector[Double]): Int = {
var initState: Int = 0
for (i <- 0 until dist.length) {
if (math.random <= dist.slice(0, i).toArray.sum) {
initState = i
return initState
}
}
initState
}
// 根据给定的参数生成观测序列
def generate(t: Int): Array[Int] = {
// require(true)
// 根据初试概览向量随机生成初始状态
val initState: Int = getDistData(pi)
// 生成第一个观测
val inner = confusionMatrix(initState, ::).inner
val initData: Int = getDistData(inner)
//生成余下的状态和序列
val datas = new ListBuffer[Int]
datas.append(initData)
for (i <- 1 until t) {
val st = getDistData(stateTransitionMatrix(initState, ::).inner)
datas append getDistData(confusionMatrix(st, ::).inner)
}
}
}
4.2 前向算法
步骤：
- 计算初值
α 1 ( i ) = π i b i ( o 1 ) , i = 1 , 2 , … … N \alpha_1(i) = \pi_ib_i(o_1) ,i =1,2,……N α1(i)=πibi(o1),i=1,2,……N - 递推对t=1，2，……T-1
α t + 1 ( i ) = [ ∑ j = 1 N α t ( j ) a j i ] b i ( o t + 1 ) , i = 1 , 2 , … … N \alpha_{t+1}(i) = [\sum_{j=1}^N \alpha_t(j)a_{ji}]b_i(o_{t+1}) ,i=1,2,……N αt+1(i)=[j=1∑Nαt(j)aji]bi(ot+1),i=1,2,……N - 终止计算
P ( O ∣ λ ) = ∑ i = 1 N α T ( i ) P(O|\lambda) = \sum_{i=1}^N \alpha_T(i) P(O∣λ)=i=1∑NαT(i)
/**
* 前向算法
*
* @param o 观测序列 观测序列
* @return 产生观测序列的概率
*/
def forwardAlgorithm(o: DenseVector[Int]) = {
// 计算初值 α1(i)=πi∗b(i)(𝑂(1))
var alphati = confusionMatrix(::, o(0)) * pi
//α(t+1)(i)= [∑ α(t)(j) a(j)(i)]*b(i)(o(t+1))
for (t <- 1 until o.length) {
val bit = confusionMatrix(::, o(t))
val dlt = new ListBuffer[Double]()
for (j <- 0 until n) {
val aji = stateTransitionMatrix(::, j)
dlt.append(sum(alphati * aji))
}
alphati = DenseVector(dlt.toArray) * bit
}
// step3 终止计算(概率)：𝑃(𝑂|𝜆)=∑ a(t)(i)
sum(alphati)
}
4.3 后向算法
步骤：
1.
β
T
(
i
)
=
1
,
i
=
1
,
2
,
…
…
N
\beta_T(i) =1 ,i=1,2,……N
βT(i)=1,i=1,2,……N
2. 对t=T-1,T-2,……1
β
t
(
i
)
=
∑
j
=
1
N
a
i
j
b
j
(
o
t
+
1
)
β
t
+
1
(
j
)
,
i
=
1
,
2
,
…
…
N
\beta_t(i)=\sum_{j=1}^N a_{ij}b_j(o_{t+1})\beta_{t+1}(j) ,i=1,2,……N
βt(i)=j=1∑Naijbj(ot+1)βt+1(j),i=1,2,……N
3. 终止计算
P
(
O
∣
λ
)
=
∑
i
=
1
N
π
i
b
i
(
o
1
)
β
1
(
i
)
P(O|\lambda) = \sum_{i=1}^N \pi_ib_i(o_1)\beta_1(i)
P(O∣λ)=i=1∑Nπibi(o1)β1(i)
/**
* 后向算法
*
* @param o 观测序列 观测序列
* @return 产生观测序列的概率
*/
def backwardAlgorithm(o: DenseVector[Int]) = {
// βT(i) = 1, i=1,2,……N
var betati: DenseVector[Double] = DenseVector(Array.fill(n)(1.0))
// βt(i) = ∑ a(i)(j)b(j)(o(t+1))β(t+1)(j)
for (t <- 1 until o.length reverse) {
println("betati:" + betati)
val beat = new ListBuffer[Double]()
for (i <- 0 until n) {
val aij = stateTransitionMatrix(i, ::).inner
val bjT = confusionMatrix(::, o(t))
sum (aij * bjT * betati )
beat.append( sum (aij * bjT * betati ))
}
betati = DenseVector(beat.toArray)
}
sum (pi * confusionMatrix(::, o(0)) * betati)
}
def main(args: Array[String]): Unit = {
val pi = DenseVector(Array(0.2, 0.4, 0.4)) //pi
val stateTransitionMatrix = DenseMatrix((0.5, 0.2, 0.3), (0.3, 0.5, 0.2), (0.2, 0.3, 0.5)) //A
val confusionMatrix = DenseMatrix((0.5, 0.5), (0.4, 0.6), (0.7, 0.3)) //B
val o = DenseVector(Array(0, 1, 0, 1)) //O
val hmmModel = hmm(pi, stateTransitionMatrix, confusionMatrix)
val forward = hmmModel.forwardAlgorithm(o)
val backward = hmmModel.backwardAlgorithm(o)
println(forward) //0.026862016
println(backward) //0.026862016
}
参考资料：
https://zhuanlan.zhihu.com/p/85454896
《统计学习方法》
https://www.zhihu.com/question/20962240/answer/1290190531
6928
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
