前面一篇介绍了隐马尔科夫模型的基本的一些概念,篇主要介绍三个问题的具体解决方法。如果对于概念不太理解的可以参考前一篇博客HMM模型基本概念,本篇博客主要介绍对于三个问题的主要推倒,内容主要基于统计学习方法这本书,但是在上面加上了一些自己的理解。下面一一介绍三个问题以及解决的办法。
概率问题
给定模型
λ=(A,B,π)
和观测序列
O=(o1,o2,...,oT)
,计算在模型
λ
下观测序列
O
出现的概率
P(O|λ)
。
直接计算法
直接计算法说白了就是暴力计算每一种情况的可能。对于所有可能的状态序列
I
求和,得到观测序列
O
的概率
P(O|λ)
,即:
P(O|λ)=∑IP(O|I,λ)P(I|λ)
=∑i1,i2,...,iTπi1bi1(o1)ai1i2bi2(o2)...aiT−1iTbiT(oT)
这种计算的缺点在于计算量很大,时间复杂度为
O(TNT)
。
前向后向算法
前向后向算法的核心是利用动态规划的思想减少计算的时间复杂度。
图1
前向算法
前向概率 给定隐马尔可夫模型
λ
,定义到时刻
t
部分观测序列为
o1,o2,...,ot
且状态为
qi
的概率为前向概率,记作
αt(i)=P(o1,o2,...,ot,it=qi|λ)(1)
然后可以递推求出前向概率
αt(i)
以及观测序列
P(O|λ)
盒子 |
1 |
2 |
3 |
4 |
红球数 |
5 |
3 |
6 |
8 |
白球数 |
5 |
7 |
4 |
2 |
那么这个前向概率到底是什么意思呢?估计好多人还是看的一脸懵逼。还是以之前的盒子与球模型,观测序列为
O={红,红,白,白,红}
,假设
t=3,i=1
,后面的序列我们不知道,那么
αt(i)=P(O={红,红,白},i=1|λ)
。即前面观测序列为[红,红,白],第三次丑的白色球是从盒子1中抽出的概率。
下面是对前向算法的形式化推导。
输入:隐马尔科夫模型
λ
,观测序列为
O
;
输出:观测序列概率
P(O|λ)
;
1.初值
α1(i)=πibi(o1)(2)
即求第一个观测值对应的状态为
i
的概率
2.递推,对t=1,2,…,T-1,
αt+1(i)=[∑j=1Nαt(j)aji]bi(ot+1),i=1,2,....,N(3)
3.终止
P(O|λ)=∑i=1NαT(i)(4)
上面这段公式推导的思想为,先求观测值为
o1
的概率,然后在前面的基础上求观测值为
o2
对应的概率,依次递推,最后求观测值为
ot
的概率。然后每一种概率对应了不同的状态,
αt(i)
则表示在状态为
i
时的概率,而
P(o1,o2,...,ot,iT=qi|λ)
,然后然后对
i
求和就能够得到上述(3)式。其递推过程图如下:
图2
前向算法就是利用这种方式将时间复杂度从
O(TNT)
降低到
O(TN2)
,至于降低的原因是减少直接利用了前面的计算结果,避免了每一次都需要重新计算。
还是以盒子与球模型为例,
λ=(A,B,π)
,状态集合
Q={1,2,3}
,观测集合
V={红,白}
.
设
T=3
,
O={红,白,红}
,求
P(O|λ)
1.计算初值
α1(1)=π1b1(o1)=0.2×0.5=0.10
α1(2)=π2b2(o1)=0.4×0.4=0.16
α1(3)=π3b3(o1)=0.4×0.7=0.28
2.递推计算
α2(1)=[∑j=13α1(j)aj1]b1(o2)
=(0.10×0.5+0.16×0.3+0.28×0.2)=0.154×0.5=0.77
α2(2)=[∑j=13α1(j)aj2]b2(o2)
α2(3)=[∑j=13α1(j)aj3]b3(o2)
同理可以得到
α3(1)=[∑j=13α2(j)aj1]b1(o3)=0.04187
α3(2)=[∑j=13α2(j)aj2]b2(o3)=0.03551
α3(3)=[∑j=13α2(j)aj3]b3(o3)=0.05284
3.终止
P(O|λ)=∑j=13α3(j)=0.13022
前向算法清楚了,其实后向概率也就清楚了。其本质就是和前向算法的思想是一样的,只不过方向相反,从后往前计算。以下图为例,先计算的是
qt
的概率,然后在计算
qt−1
的概率,以此类推。下面直接给出定义以及推导。
后向算法
后向概率
给定隐马尔科夫模型
λ
,定义在
t
时刻状态为
qt
,从
t+1
到
T
的部分观测序列为
ot+1,ot+2,...,oT
的概率为后向概率,记作
βt(i)=P(ot+1,ot+2,...,oT|it=qi,λ)
然后可以利用递推从后向前求解得到
P(O|λ)
,计算过程如下:
1.首先,对于最终时刻的所有状态
qi
规定
βT(i)=1
,即:
βT(i)=1,1,2,...,N(5)
本来应该对于这一步应该像前向算法一样存在初始状态概率,但是后向算法将初始概率放到最后一步计算,所有令
βT(i)
规定为=1
2.对于
t=T−1,T−2,...,1
βt(i)=[∑j=1Nβt+1(j)aji]bi(ot)(6)
这一步和书上不一样的,个人感觉这一步书上的写法是错误的。
3.终止
P(O|λ)=∑i=1Nπibi(o1)β1(i)(7)
并且可以将前向后向概率统一:
P(O|λ)=∑i=1N∑j=1Nαt(i)aijbj(ot+1)βt+1(j)(8)
最后化简能够得到
P(O|λ)=∑i=1Nαt+1(i)βt+1(i)(9)
老实说前向后向算法其实就是同一个东西,没太懂为什么要使用前向后向算法,明明前向算法就能得到结果。难道是前后同时计算速度更快?o(╯□╰)o
一些期望
这一部分主要是简化一些符号,为后面的计算做准备。
1.给定模型
λ
和观测
O
,在时刻
t
处于状态
qi
的概率,记为:
γi=P(it=qi|O,λ)
可以通过前向后向概率计算。
γi(i)=P(it=qi|O,λ)P(O|λ)
而
P(it=qi|O,λ)=αt(i)βt(i)
所以有
γt(i)=αt(i)βt(i)∑Ni=1αt(i)βt(i)(9)
2.给定模型
λ
和观测
O
,在时刻
t
处于状态
qi
且在时刻
t+1
处于状态
qj
的概率,记
ξt(i,j)=P(it=qi,it+1=qj|O,λ)(10)
所以有
ξt(i,j)=αt(i)aijbj(ot+1)βt+1(j)∑Ni=1∑Nj=1αt(i)aijbj(ot+1)βt+1(j)(11)
3.一些有用的期望
1)在观测
O
下状态
i
出现的期望值
∑t=1Tγt(i)(12)
2)在观测
O
下由状态
i
转移的期望值
如上图所示,从状态
t
出发,到其他的状态则称为状态转移,由于最后一个状态是不能转移到下一个状态因此,状态转移的期望值为:
∑t=1T−1γt(i)(13)
3).在观测为
O
下,由状态
i
转移到状态
j
的期望值为:
∑t=1T−1ξt(i,j)(14)
学习问题
学习问题是为了计算模型参数。在已经给定了观测序列
O
,根据是否给定状态序列
I
可以分为监督学习方法和非监督学习方法。监督学习的方法可以利用极大似然,非监督学习主要是利用Baum-Welch算法。
Baum-Welch算法
给定训练集数据为S个长度为
T
的观测序列
O={O1,O2,...,Os}
,而没有对应的状态序列,目标是学习隐马尔可夫模型
λ=(A,B,π)
,观测序列为
O
,状态序列为
I
,那么我们可以将
P(O|λ)
变成包含隐变量的概率模型:
P(O|λ)=∑IP(O|I,λ)P(I|λ)(15)
然后可以通过EM算法学习参数,步骤如下:
1.确定完全数据的对数似然函数
观测序列数据为
O=(o1,o2,...,oT)
,状态序列数据(隐数据)为
I=(i1,i2,...,it)
,完全数据为
(O,I)
,所以完全数据的对数似然函数为
logP(O,I|λ)
2.EM算法的E步,求
Q
函数
Q(λ,λ¯)
Q(λ,λ¯)=∑IP(O,I|λ¯)logP(O,I|λ)(16)
其中
P(O,I|λ)=πi1bi1(o1)ai1i2bi2(o2)....aiT−1iTbiT(oT)
所以有:
Q(λ,λ¯)=∑IP(O,I|λ¯){logπi1+∑t=1T−1logaitit+1+∑t=1Tbit(ot)}(17)
这里不得不吐槽一下李航统计学习方法的符号写法,看了半天才明白到底是怎么写的。
3.EM的M步,极大化
Q
函数
可以将上述
Q
函数拆分为三项,其中第一项为
∑IP(O,I|λ¯)logπi1
并且有
∑Ni=1πi=1
,然后利用拉格朗日乘子法,写出拉格朗日函数:
∑i=1NlogπiP(O,i1=i|λ¯)+γ(∑i=1Nπi−1)
对齐求偏导并且令结果为0,于是有
∂∂πi[∑i=1NπiP(O,i1=i|λ¯)+γ{∑i=1Nπi−1}]=0(17)
得
P(O,i1=i|λ¯)+γπi=0(18)
两边同时对
i
求和有
γ=−P(O|λ¯)
带入到上一步有
πi=P(O,i1=i|λ¯)P(O|λ¯)(19)
第二项可以写成
∑IP(O,I|λ¯)(∑t=1T−1logaitit+1)=∑i=1N∑j=1N∑t=1T−1P(O,it=i,it+1=j|λ¯)logaij(20)
和上面的类似,有约束条件
∑Nj=1aij=1
和
∑Ni=1aij=1
的拉格朗日乘子法可求出
∑i=1N∑j=1N∑t=1T−1P(O,it=i,it+1=j|λ¯)logaij+γ(∑j=1Naij−1)
两边同时对
aij
求导得
∑t=1T−1P(O,it=i,it+1=j|λ¯)+γaij=0(21)
两边同时对
j
求和有
γ=−∑t=1T−1P(O,it=i|λ¯)(22)
将式(22)带入(21)有
aij=∑T−1t=1P(O,it=i,it+1=j|λ¯)∑t=1T−1P(O,it=i,|λ¯)(23)
然后是对
bj(k)
的计算
∑IP(O,I|λ¯)(∑t=1Tlogbit(ot))=∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)(24)
并且约束条件为
∑Mk=1bj(k)=1
所以构造的拉格朗日函数为
∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)+γ(∑k=1Mbj(k)−1)=0
对
bj(k)
求导有
∂∂bj(k)[∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)+γ(∑k=1Mbj(k)−1)]=0
注意,只有当
ot=k
时,偏导才不为0所以有
∑t=1TP(O,it=j|λ¯)I(ot=k)+γbj(k)=0(25)
所以同时k求和,并且我们很容易得出
∑Mk=1I(ot=k)=1
γ=−∑t=1TP(O,ot=j|λ¯)(26)
所以将式子(26)带入(25)则有
bj(k)=∑Tt=1P(O,it=j|λ¯)I(ot=k)∑Tt=1P(O,it=j|λ¯)
Baum-Welch参数估计公式
aij=∑T−1t=1ξt(i,j)∑T−1t=1γt(i)
bj(k)=∑Tt=1,ot=kγt(j)∑Tt=1γt(j)
πi=γ1(i)
预测算法
回到上一篇博客的内容,预测问题是个什么问题呢?预测问题也叫做解码问题,即给定隐马尔科夫模型参数
λ=(A,B,π)
,以及观测序列
O=(o1,o2,o3,...,ot)
,求
P(I|O)
,即最有可能出现的状态
解决上面的问题主要有两种方法,一种是近似算法,另外一种是维比特算法。
近似算法
近似算法的思想其实很简单,在每一个
t
时刻选择最可能出现的状态
i∗t
,从而得到一个近似状态
I∗=(i∗1,i∗2,...,i∗T)
。
在t时刻处于
qi
的状态的概率为
γt(i)
所以在每一个t时刻最有可能的状态为
i∗t=argmax[γt(i)],t=1,2,3..T,1≤t≤N
从而得到状态序列
I∗=(i∗1,i∗2,...,i∗T)
.但是这种方法只保证了每一个t时刻最有可能的状态,不能保证整体,有点贪心的思想在里面
维比特算法
维比特算法其实就是利用动态规划的思想来求最大路径,类似于以前学的利用动态规划解决最短路径问题。文字描述什么的感觉不容易理解,直接上例子感觉跟容易理解。但是在说明例子之前,我们先定义两个符号
δ
和
ψ
方便后面的计算。
定义在时刻
t
状态为
i
的所有单个路径
(i1,i2,...,iT)
中概率最大的值为
δt(i)=maxi1,i2,...,iTP(it=i|λ),i=1,2,...,N,(28)
定义在时刻
t
状态为
i
的所有单个儿路径
(i1,i2,...,it−1,i)
中最大概率的路径的第
t−1
个节点为
ψt(i)=argmax1≤j≤Nψt1(j)aji
下面直接符号化维比特算法具体过程
输入:模型
λ=(A,B,π)
和观测
O=(o1,o2,...,oT)
输出:最优路径
I∗=(i∗1,i∗2,...,i∗T)
1)初始化
δ1(i)=πibi(o1),i=1,2,...,N
ψ1(i)=0
2)递推:对于
t=1,2,3...,T
δ1(i)=max1≤j≤N[δt−1(i)aji]bi(ot),i=1,2,...,N
ψi(t)=argmax1≤i≤N[δt−1aji],i=1,2,...,N
3).终止
P∗=max1≤i≤NδT(i)
i∗t=argmax1≤i≤N[δT(i)]
4).最优路径回溯。对于
t=T−1,T−2,...,1
I∗=(i∗1,i∗2,i∗3...,i∗T)
那么具体过程是什么样的呢?
还是以上盒子与球的模型为例
已知观测序列为
O=(红,白,红)
,试求最优状态序列
I∗=(i∗1,i∗2,...,i∗T)
(1)初始化,代入公式有
δ1(1)=0.10,δ1(2)=0.16,δ1(3)=0.28
记
ψ1(i)=0,i=1,2,3
最优路径求求截图如上
(2)t=2的时候如何计算呢,这里以
t=2,i=1
为例
δ2(1)=max1≤j≤3[δ1(j)aj1b1(o2)]
=maxj{0.1×0.5,0.16×0.3,0.28×0.2}×0.5
=0.028
根据计算得到当前路径从3到1的概率最大,所以有
ψ2(1)=3
根据上面的公式计算得到:
δ2(2)=0.0504,ψ2(2)=3
δ2(3)=0.042,ψ2(3)=3
δ3(1)=0.0.0756,ψ3(1)=2
δ3(2)=0.01008,ψ3(2)=2
δ3(3)=0.0147,ψ3(3)=3
所以最优路径概率为
P∗=0.0147
(3)倒推最优路径
最优路径最有一个状态对应的是3,而
ψ3(3)=3
,所以第二个状态为3,而
ψ2(3)=3
,因此第一个状态也为3
所以最优路径为
I∗=(i∗1,i∗2,i∗3)=(3,3,3)
HMM就这样水完了,感觉这篇博客写得像个草稿,有时间再将思想凝练下吧,暂时就这样窘o(╯□╰)o
参考文献
1.隐马尔可夫(HMM)、前/后向算法、Viterbi算法 再次总结
2.统计学习方法,李航