隐马尔科夫模型数学理论推导(HMM)

一、隐马尔科夫模型基本概念

隐马尔可夫模型由初始状态概率向量 π π 、状态转移概率矩阵 A 和观测概率矩阵 B 决定。A,B, π π 称为隐马尔可夫模型的三要素。

λ=(A,B,π) λ = ( A , B , π )

上图中是一个简单的描述股票市场的马尔科夫模型:其中隐藏状态为 {Bull,Bear,Even} { B u l l , B e a r , E v e n } ,观测状态为 {up,down,unchanged} { u p , d o w n , u n c h a n g e d }

根据图模型,我们可以轻易地写出状态转移概率矩阵 A (Let Bull = 1, Bear = 2, Even = 3):

A=0.60.50.40.20.30.10.20.20.5 A = [ 0.6 0.2 0.2 0.5 0.3 0.2 0.4 0.1 0.5 ]

观测概率矩阵 B (Let up = 1, down = 2, unchanged= 3):

B=0.70.10.30.10.60.30.20.304 B = [ 0.7 0.1 0.2 0.1 0.6 0.3 0.3 0.3 04 ]

隐马尔可夫模型作了两个 基本假设 yt y t 为t时刻的观测, qt q t 为t时刻的状态):

  1. p(qt|q1,,qt1,y1,,yt1)=p(qt|qt1) p ( q t | q 1 , … , q t − 1 , y 1 , … , y t − 1 ) = p ( q t | q t − 1 )
  2. p(yt|q1,,qt1,qt,y1,,yt1)=p(yt|qt) p ( y t | q 1 , … , q t − 1 , q t , y 1 , … , y t − 1 ) = p ( y t | q t )

对第一条假设的解释:隐藏的马尔可夫链在任意时刻t的状态只依赖于其 前一时刻 的状态, 与其他时刻的状态及观测无关,也与时刻t无关 ;

对第二条假设的解释:任意时刻的观测只依赖于 该时刻 的马尔可夫链的状态,与其他观测及状态无关 。

下面提出 HMM 的三个 基本问题

(记: λ=(A,B,π) λ = ( A , B , π ) Q={q1,,qT} Q = { q 1 , … , q T } Y={y1,,yT} Y = { y 1 , … , y T }

  • 概率计算问题(计算其产生观测序列的概率 ) 计算 p(Y|λ) p ( Y | λ )
  • 参数学习问题(训练模型使选取的参数能最好的描述观测数据 )  λMLE=argmaxλ p(Y|λ) ⇒   λ M L E = a r g m a x λ   p ( Y | λ )
  • 解码(decoding)问题(找到与此观测序列最匹配的隐状态序列 )  argmaxQ p(Y|Q,λ) ⇒   a r g m a x Q   p ( Y | Q , λ )

下面,我们主要讨论如何解决最常用到的前两个基本问题。

二、概率计算算法

1、直接计算法

p(Y|λ)=q1,,qT=q1,,qT(=Q p(Y,Q|λ)=Q p(Y|Q,λ) p(Q|λ)=q1,,qTp(y1,,yT|q1,,qT;λ)p(q1,,qT|λ)bq1(y1)bqT(yT)p(q1|λ)p(q2|q1,λ)p(q3|q1,q2,λ)p(qT|q1,,qT1,λ)bq1(y1)bqT(yT))(p(q1|λ)p(q2|q1,λ)p(q3|q2,λ)p(qT|qT1,λ))=q1,,qT(bq1(y1)bqT(yT))(πq1aq1,q2aqT1,T)=q1,,qTπq1t=2Taqt1,tbqt(yt)(1)(2)(3)(4)(5)(6)(7)(8)(9) (1) p ( Y | λ ) = ∑ Q   p ( Y , Q | λ ) (2) = ∑ Q   p ( Y | Q , λ )   p ( Q | λ ) (3) = ∑ q 1 , … , q T p ( y 1 , … , y T | q 1 , … , q T ; λ ) p ( q 1 , … , q T | λ ) (4) = ∑ q 1 , … , q T b q 1 ( y 1 ) … b q T ( y T ) ⋅ p ( q 1 | λ ) p ( q 2 | q 1 , λ ) p ( q 3 | q 1 , q 2 , λ ) (5) … p ( q T | q 1 , … , q T − 1 , λ ) (6) = ∑ q 1 , … , q T ( b q 1 ( y 1 ) … b q T ( y T ) ) ⋅ ( p ( q 1 | λ ) p ( q 2 | q 1 , λ ) p ( q 3 | q 2 , λ ) (7) … p ( q T | q T − 1 , λ ) ) (8) = ∑ q 1 , … , q T ( b q 1 ( y 1 ) … b q T ( y T ) ) ⋅ ( π q 1 a q 1 , q 2 … a q T − 1 , T ) (9) = ∑ q 1 , … , q T π q 1 ∏ t = 2 T a q t − 1 , t b q t ( y t )

但是这种方法的计算量很大,是 O(TNT) O ( T N T ) ,因此这种算法在实际不可行。

2、前向和后向算法


Note: Picture source

由上图定义了前向概率(左)和后向概率(右),前向概率描述了 y1 y 1 yt y t 和 t 时刻为第 i i 个状态时的联合分布,后向概率在已知t 时刻为第 i 个状态时描述了 yt+1 y t + 1 yT y T 的联合分布;

前向算法:

t=1 时,
αi(1)=p(y1,q1=i|λ)=p(q1=i|λ) p(y1|q1=i,λ)=πibi(y1) α i ( 1 ) = p ( y 1 , q 1 = i | λ ) = p ( q 1 = i | λ )   p ( y 1 | q 1 = i , λ ) = π i b i ( y 1 )

t=2 时,

αj(1)=p(y1,y2,q2=j|λ)=i=1kp(y1,y2,q1=i,q2=j|λ)=i=1kp(q1=i) p(y1|q1=i)p(q2=j|q1=i) p(q2=j|q1=i)p(y2|q2=j)=[i=1kαi(1)ai,j]bj(y2)(33)(34)(35)(36)(37) (33) α j ( 1 ) = p ( y 1 , y 2 , q 2 = j | λ ) (34) = ∑ i = 1 k p ( y 1 , y 2 , q 1 = i , q 2 = j | λ ) (35) = ∑ i = 1 k p ( q 1 = i )   p ( y 1 | q 1 = i ) p ( q 2 = j | q 1 = i ) (36)   ⋅ p ( q 2 = j | q 1 = i ) p ( y 2 | q 2 = j ) (37) = [ ∑ i = 1 k α i ( 1 ) a i , j ] b j ( y 2 )

……

因此,t 2 时, αj(t+1)=[ki=1αi(t)ai,j]bj(yt+1) α j ( t + 1 ) = [ ∑ i = 1 k α i ( t ) a i , j ] b j ( y t + 1 )

又显然有 p(Y|λ)=ki=1αi(T) p ( Y | λ ) = ∑ i = 1 k α i ( T ) ……(1)

下图可以直观地理解这个过程:

这里写图片描述

前向算法:

输入:模型参数 λ λ ,观测 Y Y

输出:p(Y|λ)

(1) 初值

αi(1)=πibi(y1) α i ( 1 ) = π i b i ( y 1 )

(2) 递推 对 t=1,……,T-1

αj(t+1)=[i=1kαi(t)ai,j]bj(yt+1) α j ( t + 1 ) = [ ∑ i = 1 k α i ( t ) a i , j ] b j ( y t + 1 )

(3) 终止

p(Y|λ)=i=1kαi(T) p ( Y | λ ) = ∑ i = 1 k α i ( T )

后向算法:

t=T 时, βi(T)=1 β i ( T ) = 1

t=T-1 时,

βi(T1)=p(yT|qT1=i)=j=1kp(yT,qT=j|qT1=i)=j=1k p(qT=j|qT1=i) p(yT|qT=j,qT1=i)=j=1k p(qT=j|qT1=i) p(yT|qT=j)=j=1kai,jbj(yT)(15)(16)(17)(18)(19) (15) β i ( T − 1 ) = p ( y T | q T − 1 = i ) (16) = ∑ j = 1 k p ( y T , q T = j | q T − 1 = i ) (17) = ∑ j = 1 k   p ( q T = j | q T − 1 = i )   p ( y T | q T = j , q T − 1 = i ) (18) = ∑ j = 1 k   p ( q T = j | q T − 1 = i )   p ( y T | q T = j ) (19) = ∑ j = 1 k a i , j b j ( y T )

t=T-2 时,

βi(T2)=p(yT,yT1|qT2=i)=j=1kp(yT,yT1,qT1=j,qT=l|qT2=i)=j=1kl=1k p(qT=j|qT1=i) p(yT|qT=j)p(qT1=j|qT2=i) p(yT1|qT1=j)=j=1k p(qT=j|qT1=i) p(yT|qT=j)=j=1kai,jbj(yT1)βj(T1)(20)(21)(22)(23)(24)(25) (20) β i ( T − 2 ) = p ( y T , y T − 1 | q T − 2 = i ) (21) = ∑ j = 1 k p ( y T , y T − 1 , q T − 1 = j , q T = l | q T − 2 = i ) (22) = ∑ j = 1 k ∑ l = 1 k   p ( q T = j | q T − 1 = i )   p ( y T | q T = j ) ⋅ (23) p ( q T − 1 = j | q T − 2 = i )   p ( y T − 1 | q T − 1 = j ) (24) = ∑ j = 1 k   p ( q T = j | q T − 1 = i )   p ( y T | q T = j ) (25) = ∑ j = 1 k a i , j b j ( y T − 1 ) β j ( T − 1 )

……

因此, tT1 t ≤ T − 1 时:

βi(t)=j=1kai,jbj(yt+1)βj(t+1) β i ( t ) = ∑ j = 1 k a i , j b j ( y t + 1 ) β j ( t + 1 )

又显然有 p(Y|λ)=ki=1πibi(y1)βi(1) p ( Y | λ ) = ∑ i = 1 k π i b i ( y 1 ) β i ( 1 ) ……(2)

后向算法:

输入: 模型参数 λ λ ,观测 Y Y

输出:p(Y|λ)

(1) 初值

βi(T)=1 β i ( T ) = 1

(2) 递推 对 t=T-1,……,1

βi(t)=j=1kai,jbj(yt+1)βj(t+1) β i ( t ) = ∑ j = 1 k a i , j b j ( y t + 1 ) β j ( t + 1 )

(3)终止

p(Y|λ)=i=1kπibi(y1)βi(1) p ( Y | λ ) = ∑ i = 1 k π i b i ( y 1 ) β i ( 1 )

前向算法 和 后向算法 的统一

利用前面的定义可以将观测序列概率 p(Y|λ) p ( Y | λ ) 统一:

p(Y|λ)=i=1kj=1kαi(t)ai,jbj(yt+1)βj(t+1)t=1,,T1 p ( Y | λ ) = ∑ i = 1 k ∑ j = 1 k α i ( t ) a i , j b j ( y t + 1 ) β j ( t + 1 ) , t = 1 , … , T − 1

此式当 t=1 和 t=T-1 时分别为式 (1) 和(2)。

特定状态的概率

给定模型 λ λ 和 观测 Y Y ,在时刻 t 处于状态 i 的概率如下;

p(qt=i|Y,λ)=p(Y,qt=i|λ)p(Y|λ)=p(Y,qt=i|λ)kj=1p(Y,qt=j|λ)=αi(t)βi(t)kj=1αi(t)βi(t) p ( q t = i | Y , λ ) = p ( Y , q t = i | λ ) p ( Y | λ ) = p ( Y , q t = i | λ ) ∑ j = 1 k p ( Y , q t = j | λ ) = α i ( t ) β i ( t ) ∑ j = 1 k α i ( t ) β i ( t )

其中,

p(Y,qt=i|λ)=p(Y|qt=i,λ)p(qt=i|λ)=p(y1,,yt|qt=i)p(yt+1,,yT|qt=i)p(qt=i|λ)=p(y1,,yt,qt=i|λ)p(yt+1,,yT|qt=i)=αi(t)βi(t)(26)(27)(28)(29) (26) p ( Y , q t = i | λ ) = p ( Y | q t = i , λ ) p ( q t = i | λ ) (27) = p ( y 1 , … , y t | q t = i ) p ( y t + 1 , … , y T | q t = i ) p ( q t = i | λ ) (28) = p ( y 1 , … , y t , q t = i | λ ) p ( y t + 1 , … , y T | q t = i ) (29) = α i ( t ) β i ( t )

三、参数学习算法

根据训练数据是包括观测序列和对应的状态序列还是只有观测序列,可以分别由监督学习与非监督学习实现。但由于监督学习需要使用训练数据,而人工标注数据代价也往往很高,因此我们会利用非监督的学习方法来学习参数。

将观测序列数据看作观测数据 Y,状态序列数据看作不可观测的隐数据 Q,那么隐马尔可夫模型事实上是一个含有隐变量的概率模型 :

p(Y|λ)=Q p(Y|Q,λ) p(Q|λ) p ( Y | λ ) = ∑ Q   p ( Y | Q , λ )   p ( Q | λ )

我们先回顾一下EM算法:EM算法推导(收敛性证明和在GMM中的应用)

在 HMM 中,我们可以写成如下:

λ(g+1)=argmaxλ(qln(p(Y,q|λ))p(q|Y,λ(g)))=argmaxλ(qln(p(Y,q|λ))p(q,Y|λ(g))p(Y|λ(g)))=argmaxλ(qln(p(Y,q|λ))p(q,Y|λ(g)))(30)(31)(32) (30) λ ( g + 1 ) = a r g m a x λ ( ∫ q l n ( p ( Y , q | λ ) ) p ( q | Y , λ ( g ) ) ) (31) = a r g m a x λ ( ∫ q l n ( p ( Y , q | λ ) ) p ( q , Y | λ ( g ) ) p ( Y | λ ( g ) ) ) (32) = a r g m a x λ ( ∫ q l n ( p ( Y , q | λ ) ) p ( q , Y | λ ( g ) ) )

p(Y|λ(g)) p ( Y | λ ( g ) ) 为常数,因此可以省略)

E step:

求 Q 函数,见如下:

Q(λ,λ(g))=qln(p(Y,q|λ))p(q,Y|λ(g))=q1,,qT(lnπq1)p(q,Y|λ(g))+q1,,qT(t=1Tlnaqt1,qt)p(q,Y|λ(g))+q1,,qT(t=1Tlnbqt(yt))p(q,Y|λ(g)) Q ( λ , λ ( g ) ) = ∫ q l n ( p ( Y , q | λ ) ) p ( q , Y | λ ( g ) ) = ∑ q 1 , … , q T ( l n π q 1 ) p ( q , Y | λ ( g ) ) + ∑ q 1 , … , q T ( ∑ t = 1 T l n a q t − 1 , q t ) p ( q , Y | λ ( g ) ) + ∑ q 1 , … , q T ( ∑ t = 1 T l n b q t ( y t ) ) p ( q , Y | λ ( g ) )

M step:

极大化 Q 函数,求模型参数 A,B, π π

观察上述 Q 函数,要极大化的参数分别单独地出现在3个项中,所以只需对各项分别极大化;

(1)第一项可以写出:

q1,,qT(lnπq1)p(q,Y|λ(g))=i=1k(lnπi)p(q1=i,Y|λ(g)) ∑ q 1 , … , q T ( l n π q 1 ) p ( q , Y | λ ( g ) ) = ∑ i = 1 k ( l n π i ) p ( q 1 = i , Y | λ ( g ) )

s.t.i=1kπi=1 s . t . ∑ i = 1 k π i = 1

利用拉格朗日乘子法,即可求解;

πi=p(q1=i,Y|λ(g))ki=1p(q1=i,Y|λ(g)) π i = p ( q 1 = i , Y | λ ( g ) ) ∑ i = 1 k p ( q 1 = i , Y | λ ( g ) )

(1)第二项可以写出:

q1,,qT(t=1Tlnaqt1,qt)p(q,Y|λ(g))=i=1kj=1kt=1Tlnai,j p(qt1=i,qt=j,Y|λ(g)) ∑ q 1 , … , q T ( ∑ t = 1 T l n a q t − 1 , q t ) p ( q , Y | λ ( g ) ) = ∑ i = 1 k ∑ j = 1 k ∑ t = 1 T l n a i , j   p ( q t − 1 = i , q t = j , Y | λ ( g ) )

s.t.i=1kaij=1 s . t . ∑ i = 1 k a i j = 1

同样利用拉格朗日乘子法,即可求解:

ai,j=Tt=1p(qt1=i,qt=j,Y|λ(g))Tt=1p(qt1=i,Y|λ(g)) a i , j = ∑ t = 1 T p ( q t − 1 = i , q t = j , Y | λ ( g ) ) ∑ t = 1 T p ( q t − 1 = i , Y | λ ( g ) )

(3)第二项可以写出:

q1,,qT(t=1Tlnbqt(yt))p(q,Y|λ(g))=j=1kt=1T(lnbj(yt))p(qt=j,Y|λ(g)) ∑ q 1 , … , q T ( ∑ t = 1 T l n b q t ( y t ) ) p ( q , Y | λ ( g ) ) = ∑ j = 1 k ∑ t = 1 T ( l n b j ( y t ) ) p ( q t = j , Y | λ ( g ) )

s.t.i=1kbj(yt)=1 s . t . ∑ i = 1 k b j ( y t ) = 1

同样利用拉格朗日乘子法,即可求解(注意:只有在 yt=vl y t = v l 时,偏导数才不为0):

bj(yt=vl)=Tt=1p(qt=j,Y=vl|λ(g))Tt=1p(qt=j,Y|λ(g)) b j ( y t = v l ) = ∑ t = 1 T p ( q t = j , Y = v l | λ ( g ) ) ∑ t = 1 T p ( q t = j , Y | λ ( g ) )

上面的 EM 算法又称为 Baum-Welch 算法。

Baum-Welch 算法

输入:观测序列Y

(1) 初始化 λ0=(A(0),B(0),π(0)) λ 0 = ( A ( 0 ) , B ( 0 ) , π ( 0 ) )

(2) 递推直至EM算法收敛

π(g+1)i=p(q1=i,Y|λ(g))ki=1p(q1=i,Y|λ(g)) π i ( g + 1 ) = p ( q 1 = i , Y | λ ( g ) ) ∑ i = 1 k p ( q 1 = i , Y | λ ( g ) )

a(g+1)i,j=Tt=1p(qt1=i,qt=j,Y|λ(g))Tt=1p(qt1=i,Y|λ(g)) a i , j ( g + 1 ) = ∑ t = 1 T p ( q t − 1 = i , q t = j , Y | λ ( g ) ) ∑ t = 1 T p ( q t − 1 = i , Y | λ ( g ) )

bj(l)(g+1)=Tt=1p(qt=j,Y=vl|λ(g))Tt=1p(qt=j,Y|λ(g)) b j ( l ) ( g + 1 ) = ∑ t = 1 T p ( q t = j , Y = v l | λ ( g ) ) ∑ t = 1 T p ( q t = j , Y | λ ( g ) )

(3) 终止。得到参数 λ(n+1)=(A(n+1),B(n+1),π(n+1)) λ ( n + 1 ) = ( A ( n + 1 ) , B ( n + 1 ) , π ( n + 1 ) )

四、参考资料


[1] 李航《统计学习方法》
[2] 徐亦达教授的自视频
[3] machine-learning-notes.Professor Richard Xu .

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值