统计学习方法读书笔记第十章:隐马尔科夫模型
统计学习方法读书笔记第十章:隐马尔科夫模型
隐马尔科夫模型是可用于标注问题的统计学模型,描述由隐藏的马尔科夫链随机生成观测序列的过程,属于生成模型。
隐马尔科夫模型的基本概念
-
隐马尔科夫模型的定义
隐马尔科夫模型 隐马尔科夫模型是关于时序的概率模型,描述有一个隐藏的马尔科夫链随机生成不可观测的状态随机序列,再由各个状态生成一个观测而产生观测随机序列的过程。隐藏的马尔科夫链随机生成的状态的序列,称为状态序列;每个状态生成一个观测,而由此产生的观测的随机序列,称为观测序列。序列的每一个位置又可以看作是一个时刻。
隐马尔科夫模型由初始概率分布、状态转移概率分布以及观测概率分布确定。应马尔科夫模型的形式定义如下:
设 Q Q Q是所有可能的状态的集合, V V V是所有可能的观测的集合。
Q = { q 1 , q 2 , ⋯   , q N } , V = { v 1 , v 2 , ⋯   , v M } Q=\{q_1,q_2,\cdots,q_N\}, V=\{v_1,v_2,\cdots,v_M\} Q={q1,q2,⋯,qN},V={v1,v2,⋯,vM}
其中, N N N是可能的状态数, M M M是可能的观测数。
I I I是长度为 T T T的状态序列, O O O是对应的观测序列。
I = ( i 1 , i 2 , ⋯   , i T ) , O = ( o 1 , o 2 , ⋯   , o T ) I=(i_1,i_2,\cdots,i_T), O=(o_1,o_2,\cdots,o_T) I=(i1,i2,⋯,iT),O=(o1,o2,⋯,oT)
A A A是状态转移概率矩阵:
(1) A = [ a i j ] N × N A=[a_{ij}]_{N\times N} \tag{1} A=[aij]N×N(1)
其中,
(2) a i j = P ( i t + 1 = q j ∣ i t = q i ) , i = 1 , 2 , ⋯   , N ; j = 1 , 2 , ⋯   , N a_{ij}=P(i_{t+1}=q_j|i_t=q_i), i=1,2,\cdots,N;j=1,2,\cdots,N \tag{2} aij=P(it+1=qj∣it=qi),i=1,2,⋯,N;j=1,2,⋯,N(2)
是在时刻 t t t处于状态 q i q_i qi的条件下在时刻 t + 1 t+1 t+1转移到状态 q j q_j qj的概率。
B B B是观测概率矩阵:
(3) B = [ b j ( k ) ] N × M B=[b_j(k)]_{N\times M} \tag{3} B=[bj(k)]N×M(3)
其中,
(4) b j ( k ) = P ( o t = v k ∣ i t = q j ) , k = 1 , 2 , ⋯   , M ; j = 1 , 2 , ⋯   , N b_j(k)=P(o_t=v_k|i_t=q_j), k=1,2,\cdots,M; j=1,2,\cdots,N \tag{4} bj(k)=P(ot=vk∣it=qj),k=1,2,⋯,M;j=1,2,⋯,N(4)
是在时刻 t t t处于状态 q j q_j qj的条件下生成观测 v k v_k vk的概率。
π \pi π是初始状态概率向量:
(5) π = ( π i ) \pi=(\pi_i) \tag{5} π=(πi)(5)
其中,
(6) π i = P ( i 1 = q i ) , i = 1 , 2 , ⋯   , N \pi_i=P(i_1=q_i), i=1,2,\cdots,N \tag{6} πi=P(i1=qi),i=1,2,⋯,N(6)
是时刻 t = 1 t=1 t=1处于状态 q i q_i qi的概率。
隐马尔科夫模型由初始状态概率向量 π \pi π、状态转移概率矩阵 A A A和观测概率矩阵 B B B决定。 π \pi π和 A A A决定状态序列, B B B决定观测序列。因此,隐马尔科夫模型 λ \lambda λ可以用三原符号表示,即
(7) λ = ( A , B , π ) \lambda=(A,B,\pi) \tag{7} λ=(A,B,π)(7)
A A A, B B B, π \pi π称为隐马尔科夫模型的三要素。
状态转移概率矩阵 A A A与初始状态概率向量 π \pi π确定了隐藏的马尔科夫链,生成不可观测的状态序列。观测概率矩阵 B B B确定了如何从状态生成观测,与状态序列综合确定了如何产生观测序列。
从定义可知,隐马尔科夫模型作了两个基本假设:
(1) 齐次马尔科夫性假设,即假设隐藏的马尔科夫链在任意时刻 t t t的状态只依赖于其前一时刻的状态,与其他时刻的状态及观测无关,也与时刻 t t t无关。
(8) P ( i t ∣ i t − 1 , o t − 1 , ⋯   , i 1 , o 1 ) = P ( i t ∣ i t − 1 ) , t = 1 , 2 , ⋯   , T P(i_t|i_{t-1},o_{t-1},\cdots,i_1,o_1)=P(i_t|i_{t-1}), t=1,2,\cdots,T \tag{8} P(it∣it−1,ot−1,⋯,i1,o1)=P(it∣it−1),t=1,2,⋯,T(8)
(2) 观测独立性假设,即假设任意时刻的观测只依赖于该时刻的马尔科夫链的状态,与其他观测及状态无关。
(9) P ( o t ∣ i T , o T , i T − 1 , o T − 1 , ⋯   , i t + 1 , o t + 1 , i t , i t − 1 , o t − 1 , ⋯   , i 1 , o 1 ) = P ( o t ∣ i t ) P(o_t|i_T,o_T,i_{T-1},o_{T-1},\cdots,i_{t+1},o_{t+1},i_t,i_{t-1},o_{t-1},\cdots,i_1,o_1)=P(o_t|i_t) \tag{9} P(ot∣iT,oT,iT−1,oT−1,⋯,it+1,ot+1,it,it−1,ot−1,⋯,i1,o1)=P(ot∣it)(9)
隐马尔科夫模型可以用于标注,这时状态对应着标记。标注问题是给定观测的序列预测其对应的标记序列。可以假设标注问题的数据是由隐马尔科夫模型生成的。这样我们可以利用隐马尔科夫模型的学习与预测算法进行标注。 -
观测序列的生成过程
根据马尔科夫模型的定义,可以将一个长度为 T T T的观测序列 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT)的生成过程描述如下:
算法1(观测序列的生成)
输入:隐马尔科夫模型 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π),观测序列长度 T T T;
输出:观测序列 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT)。
(1) 按照初始状态分布 π \pi π产生状态 i 1 i_1 i1;
(2) 令 t = 1 t=1 t=1;
(3) 按照状态 i t i_t it的观测概率分布 b i t ( k ) b_{i_t}(k) bit(k)生成 o t o_t ot;
(4) 按照状态 i t i_t it的状态转移概率分布 { a i t i t + 1 } \{a_{i_ti_{t+1}}\} {aitit+1}产生状态 i t + 1 , i t + 1 = 1 , 2 , ⋯   , N i_{t+1},i_{t+1}=1,2,\cdots,N it+1,it+1=1,2,⋯,N;
(5) 令 t = t + 1 t=t+1 t=t+1;如果 t < T t<T t<T,转步(3);否则,终止。 -
隐马尔科夫模型的3个基本问题
隐马尔科夫模型有3个基本为题:
(1) 概率模型问题。给定模型 λ = ( a , N , π ) \lambda=(a,N,\pi) λ=(a,N,π)和观测序列 O = ( o 1 , o , 2 ⋯   , o T ) O=(o_1,o_,2\cdots,o_T) O=(o1,o,2⋯,oT),计算在模型KaTeX parse error: Expected 'EOF', got '\ambda' at position 1: \̲a̲m̲b̲d̲a̲下观测序列 O O O出现的概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
(2) 学习问题。已知观测序列 O = ( o 1 , o , 2 ⋯   , o T ) O=(o_1,o_,2\cdots,o_T) O=(o1,o,2⋯,oT),估计模型参数 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π)的参数,使得在该模型下观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)最大。即用极大似然估计的方法估计参数。
(3) 预测问题,也称为解码问题。已知模型 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π)和序列 O = ( o 1 , o , 2 ⋯   , o T ) O=(o_1,o_,2\cdots,o_T) O=(o1,o,2⋯,oT),求给定观测序列条件概率 P ( I ∣ O ) P(I|O) P(I∣O)最大的状态序列 I = ( i 1 , i 2 , ⋯   , i T ) I=(i_1,i_2,\cdots,i_T) I=(i1,i2,⋯,iT)。即给定观测序列,求最有可能的对应的状态序列。
概率计算算法
本节介绍计算观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的前向与后向算法。先介绍概念上可行但计算上不可行的直接计算法。
- 直接计算法
给定模型 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π)和观测序列 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT),计算观测序列 O O O出现的概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。最直接的方法是按概率公式直接计算。通过列举所有可能的长度为 T T T的状态序列 I = ( i 1 , i 2 , ⋯   , i T ) I=(i_1,i_2,\cdots,i_T) I=(i1,i2,⋯,iT),求各个状态序列与观测序列 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT)的联合概率 P ( O , I ∣ λ ) P(O,I|\lambda) P(O,I∣λ),然后对所有可能的状态序列求和,得到 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
状态序列 I = ( i 1 , i 2 , ⋯   , i T ) I=(i_1,i_2,\cdots,i_T) I=(i1,i2,⋯,iT)的概率是
(10) P ( I ∣ λ ) = π i 1 a i 1 i 2 a i 2 i 3 ⋯ a i T − 1 i T P(I|\lambda)=\pi_{i_1}a_{i_1i_2}a_{i_2i_3}\cdots a_{i_{T-1}i_T} \tag{10} P(I∣λ)=πi1ai1i2ai2i3⋯aiT−1iT(10)
对固定的状态序列 I = ( i 1 , i 2 , ⋯   , i T ) I=(i_1,i_2,\cdots,i_T) I=(i1,i2,⋯,iT),观测序列 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT)的概率是 P ( O ∣ I , λ ) P(O|I,\lambda) P(O∣I,λ)
(11) P ( O ∣ I , λ ) = b i 1 ( o 1 ) b i 2 ( o 2 ) ⋯ b i T ( o T ) P(O|I,\lambda)=b_{i_1}(o_1)b_{i_2}(o_2)\cdots b_{i_T}(o_T) \tag{11} P(O∣I,λ)=bi1(o1)bi2(o2)⋯biT(oT)(11)
O O O和 I I I同时出现的联合概率为
(12) P ( O , I ∣ λ ) = P ( O ∣ I , λ ) P ( I ∣ λ ) = π i 1 b i 1 ( o 1 ) a i 1 i 2 b i 2 ( o 2 ) ⋯ a i T − 1 i T b i T ( o T ) \begin{aligned} P(O,I|\lambda)&=P(O|I,\lambda)P(I|\lambda) \\ &=\pi_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)\cdots a_{i_{T-1}i_T}b_{i_T}(o_T) \tag{12} \end{aligned} P(O,I∣λ)=P(O∣I,λ)P(I∣λ)=πi1bi1(o1)ai1i2bi2(o2)⋯aiT−1iTbiT(oT)(12)
然后,对所有可能的状态序列 I I I求和,得到观测序列 O O O的概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ),即
(13) P ( O ∣ λ ) = ∑ I P ( O ∣ I , λ ) P ( I ∣ λ ) = ∑ i 1 , i 2 , ⋯   , i T π i 1 b i 1 ( o 1 ) a i 1 i 2 b i 2 ( o 2 ) ⋯ a i T − 1 i T b i T ( o T ) \begin{aligned} P(O|\lambda)&=\sum_{I}P(O|I,\lambda)P(I|\lambda) \\ &=\sum_{i_1,i_2,\cdots,i_T}\pi_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)\cdots a_{i_{T-1}i_T}b_{i_T}(o_T) \tag{13} \end{aligned} P(O∣λ)=I∑P(O∣I,λ)P(I∣λ)=i1,i2,⋯,iT∑πi1bi1(o1)ai1i2bi2(o2)⋯aiT−1iTbiT(oT)(13)
但是,利用上式计算量很大,是 O ( T N T ) O(TN^T) O(TNT)阶的,这种算法不可行。
下面介绍计算观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的有效算法:前向-后向算法。 - 前向算法
前向概率 给定隐马尔科夫模型 λ \lambda λ,定义到时刻 t t t部分观测序列为 o 1 , o 2 , ⋯   , o t o_1,o_2,\cdots,o_t o1,o2,⋯,ot且状态为 q i q_i qi的概率为前向概率,记作
(14) α t ( i ) = P ( o 1 , o 2 , ⋯   , o t , i t = q i ∣ λ ) \alpha_t(i)=P(o_1,o_2,\cdots,o_t,i_t=q_i|\lambda) \tag{14} αt(i)=P(o1,o2,⋯,ot,it=qi∣λ)(14)
可以递推地求得前向概率 α t ( i ) \alpha_t(i) αt(i)及观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
观测序列概率的前向算法
输入:隐马尔科夫模型 λ \lambda λ,观测序列 O O O;
输出:观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
(1) 初值
(15) α 1 ( i ) = π i b i ( o 1 ) , i = 1 , 2 , ⋯   , N \alpha_1(i)=\pi_ib_i(o_1),i=1,2,\cdots,N \tag{15} α1(i)=πibi(o1),i=1,2,⋯,N(15)
(2) 递推 对 t = 1 , 2 , ⋯   , T − 1 t=1,2,\cdots,T-1 t=1,2,⋯,T−1,
(116) α t + 1 ( i ) = [ ∑ j = 1 N α t ( j ) a j i ] b i ( o t + 1 ) , i = 1 , 2 , ⋯   , N \alpha_{t+1}(i)=\bigg[\sum_{j=1}^{N}\alpha_t(j)a_{ji}\bigg]b_i(o_{t+1}), i=1,2,\cdots,N \tag{116} αt+1(i)=[j=1∑Nαt(j)aji]bi(ot+1),i=1,2,⋯,N(116)
(3) 终止
(17) P ( O ∣ λ ) = ∑ i = 1 N α T ( i ) P(O|\lambda)=\sum_{i=1}^N\alpha_T(i) \tag{17} P(O∣λ)=i=1∑NαT(i)(17)
前向算法,步骤(1)初始化前向概率,是初始时刻的状态 i 1 = q i i_1=q_i i1=qi和观测 o 1 o_1 o1的联合概率。步骤(2)是前向概率的递推公式,计算到时刻 t + 1 t+1 t+1部分观测序列为 o 1 , o 2 , ⋯   , o t , o t + 1 o_1,o_2,\cdots,o_t,o_{t+1} o1,o2,⋯,ot,ot+1且在时刻 t + 1 t+1 t+1处于状态 q i q_i qi的前向概率,如下图所示。在(2)式的方括弧里,既然 α t ( j ) \alpha_t(j) αt(j)是到时刻 t t t观测到 o 1 , o 2 , ⋯   , o t o_1,o_2,\cdots,o_t o1,o2,⋯,ot并在时刻 t t t处于状态 q j q_j qj的前向概率,那么乘积 α t ( j ) a j i \alpha_t(j)a_{ji} αt(j)aji就是到时刻 t t t观测到 o 1 , o 2 , ⋯   , o t o_1,o_2,\cdots,o_t o1,o2,⋯,ot并在时刻 t t t处于状态 q j q_j qj而在时刻 t + 1 t+1 t+1到达状态 q i q_i qi的联合概率。对这个乘积在时刻 t t t的所有可能的 N N N个状态 q j q_j qj求和,其结果就是到时刻 t t t观测为 o 1 , o 2 , ⋯   , o t o_1,o_2,\cdots,o_t o1,o2,⋯,ot并在时刻 t + 1 t+1 t+1处于状态 q i q_i qi的联合概率。方括弧里的值与观测概率 b i ( o t + 1 ) b_i(o_{t+1}) bi(ot+1)的乘积恰好是时刻 t + 1 t+1 t+1观测到 o 1 , o 2 , ⋯   , o t , o t + 1 o_1,o_2,\cdots,o_t,o_{t+1} o1,o2,⋯,ot,ot+1并在时刻 t + 1 t+1 t+1处于状态 q i q_i qi的前向概率 α t + 1 ( i ) \alpha_{t+1}(i) αt+1(i)。步骤(3)给出 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的计算公式。因为
α T ( i ) = P ( o 1 , o 2 , ⋯   , o T , i T = q i ∣ λ ) \alpha_T(i)=P(o_1,o_2,\cdots,o_T,i_T=q_i|\lambda) αT(i)=P(o1,o2,⋯,oT,iT=qi∣λ)
所以
P ( O ∣ λ ) = ∑ i = 1 N α T ( i ) P(O|\lambda)=\sum_{i=1}^{N}\alpha_T(i) P(O∣λ)=i=1∑NαT(i)
如下图所示,前向算法实际是基于“状态序列的路径结构”递推计算 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的算法。前向算法高效的关键是其局部计算前向概率,然后利用路径结构将前向概率“递推”到全局,得到 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。具体地,在时刻 t = 1 t=1 t=1,计算 α 1 ( i ) \alpha_1(i) α1(i)的 N N N个值 ( i = 1 , 2 , ⋯   , N ) (i=1,2,\cdots,N) (i=1,2,⋯,N);在各个时刻 t = 1 , 2 , ⋯   , T − 1 t=1,2,\cdots,T-1 t=1,2,⋯,T−1,计算 α t + 1 ( i ) \alpha_{t+1}(i) αt+1(i)的 N N N个值 ( i = 1 , 2 , ⋯   , N ) (i=1,2,\cdots,N) (i=1,2,⋯,N),并且每个 α t + 1 ( i ) \alpha_{t+1}(i) αt+1(i)的计算利用前一时刻 N N N个 α t ( j ) \alpha_t(j) αt(j)。减少计算量的原因在于每一次计算直接引用前一个时刻的计算结果,避免重复计算。这样,利用前向概率计算 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的计算量是 O ( N 2 T ) O(N^2T) O(N2T)阶的,而不是直接计算的 O ( T N T ) O(TN^T) O(TNT)阶。
- 后向算法
后向概率 给定隐马尔科夫模型 λ \lambda λ,定义在时刻 t t t状态为 q i q_i qi的条件下,从 t + 1 t+1 t+1到 T T T的部分观测序列为 o t + 1 , o t + 2 , ⋯   , o T o_{t+1},o_{t+2},\cdots,o_T ot+1,ot+2,⋯,oT的概率为后向概率,记作
(18) β t ( i ) = P ( o t + 1 , o t + 2 , ⋯   , o T ∣ i t = q i , λ ) \beta_t(i)=P(o_{t+1},o_{t+2},\cdots,o_T|i_t=q_i,\lambda) \tag{18} βt(i)=P(ot+1,ot+2,⋯,oT∣it=qi,λ)(18)
可以利用递推的方法求得后向概率 β t ( i ) \beta_t(i) βt(i)及观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
观测序列概率的后向算法
输入:隐马尔科夫模型 λ \lambda λ,观测序列 O O O;
输出:观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)。
(1) (19) β T ( i ) = 1 , i = 1 , 2 , ⋯   , N \beta_T(i)=1,i=1,2,\cdots,N \tag{19} βT(i)=1,i=1,2,⋯,N(19)
(2) 对 t = T − 1 , T − 2 , ⋯   , 1 t=T-1,T-2,\cdots,1 t=T−1,T−2,⋯,1
(20) β 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,\cdots,N \tag{20} βt(i)=j=1∑Naijbj(ot+1)βt+1(j),i=1,2,⋯,N(20)
(3) (21) P ( O ∣ λ ) = ∑ i = 1 N π i b i ( o 1 ) β 1 ( i ) P(O|\lambda)=\sum_{i=1}^N\pi_{i}b_i(o_1)\beta_1(i) \tag{21} P(O∣λ)=i=1∑Nπibi(o1)β1(i)(21)
步骤(1)初始化后向概率,对最终时刻的所有状态 q i q_i qi规定 β T ( i ) = 1 \beta_T(i)=1 βT(i)=1。步骤 ( 2 ) (2) (2)是后向概率的递推公式。如下图所示,为了计算在时刻 t t t状态为 q i q_i qi条件下时刻 t + 1 t+1 t+1之后的观测序列 o t + 1 , o t + 2 , ⋯   , o T o_{t+1},o_{t+2},\cdots,o_T ot+1,ot+2,⋯,oT的后向概率 β t ( i ) \beta_t(i) βt(i),只需要考虑在时刻 t + 1 t+1 t+1所有可能的 N N N个状态 q j q_j qj的转移概率(即 a i j a_{ij} aij项),以及在此状态下的观测 o t + 1 o_{t+1} ot+1的观测概率(即 b j ( o t + 1 ) b_j(o_{t+1}) bj(ot+1)项),然后考虑状态 q j q_j qj之后的观测序列的后向概率(即 β t + 1 ( j ) \beta_{t+1}(j) βt+1(j)项)。步骤(3)求 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)的思路与步骤(2)一致,只是初始概率 π i \pi_i πi代替转移概率。
利用前向概率和后向概率的定义可以将观测序列概率 P ( O ∣ λ ) P(O|\lambda) P(O∣λ)统一写成
(22) P ( O ∣ λ ) ) = ∑ i = 1 N ∑ j = 1 N α t ( i ) a i j b j ( o t + 1 ) β t + 1 ( j ) , t = 1 , 2 , ⋯   , T − 1 P(O|\lambda))=\sum_{i=1}^N\sum_{j=1}^N\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j),t=1,2,\cdots,T-1 \tag{22} P(O∣λ))=i=1∑Nj=1∑Nαt(i)aijbj(ot+1)βt+1(j),t=1,2,⋯,T−1(22)
此式当 t = 1 t=1 t=1和 t = T − 1 t=T-1 t=T−1时分别为前向概率计算公式和后向概率计算公式。
- 一些概率与期望值的计算
利用前向概率和后向概率,可以得到关于单个状态和两个状态概率的计算公式。
- 给定模型
λ
\lambda
λ和观测
O
O
O,在时刻
t
t
t处于状态
q
i
q_i
qi的概率。记
(23) γ t ( i ) = P ( i t = q i ∣ O , λ ) \gamma_t(i)=P(i_t=q_i|O,\lambda) \tag{23} γt(i)=P(it=qi∣O,λ)(23)
可以通过前向后向概率计算。事实上,
γ t ( i ) = P ( i t = q i ∣ O , λ ) = P ( i t = q i , O ∣ λ ) P ( O ∣ λ ) \gamma_t(i)=P(i_t=q_i|O,\lambda)=\frac{P(i_t=q_i,O|\lambda)}{P(O|\lambda)} γt(i)=P(it=qi∣O,λ)=P(O∣λ)P(it=qi,O∣λ)
由前向概率 α t ( i ) \alpha_t(i) αt(i)和后向概率 β t ( i ) \beta_t(i) βt(i)定义可知:
α t ( i ) β t ( i ) = P ( i t = q i , O ∣ λ ) \alpha_t(i)\beta_t(i)=P(i_t=q_i,O|\lambda) αt(i)βt(i)=P(it=qi,O∣λ)
于是得到:
(24) γ t ( i ) = α t ( i ) β t ( i ) P ( O ∣ λ ) = α t ( i ) β t ( i ) ∑ j = 1 N α t ( j ) β t ( j ) \gamma_t(i)=\frac{\alpha_t(i)\beta_t(i)}{P(O|\lambda)}=\frac{\alpha_t(i)\beta_t(i)} \tag{24}{\sum_{j=1}^N\alpha_t(j)\beta_t(j)} γt(i)=P(O∣λ)αt(i)βt(i)=∑j=1Nαt(j)βt(j)αt(i)βt(i)(24) - 给定模型
λ
\lambda
λ和观测
O
O
O,在时刻
t
t
t处于状态
q
i
q_i
qi且在时刻
t
+
1
t+1
t+1处于状态
q
j
q_j
qj的概率。记
(25) ξ t ( i , j ) = P ( i t = q i , i t + 1 = q j , O ∣ λ ) P ( O ∣ λ ) \xi_t(i,j)=\frac{P(i_t=q_i,i_{t+1}=q_j,O|\lambda)}{P(O|\lambda)} \tag{25} ξt(i,j)=P(O∣λ)P(it=qi,it+1=qj,O∣λ)(25)
可以通过前向后向概率计算:
ξ t ( i , j ) = P ( i t = q i , i t + 1 = q j , O ∣ λ ) P ( O ∣ λ ) = P ( i t = q i , i t + 1 = q j , O ∣ λ ) ∑ i = 1 N ∑ j = 1 N P ( i t = q i , i t + 1 = q j , O ∣ λ ) \xi_t(i,j)=\frac{P(i_t=q_i,i_{t+1}=q_j,O|\lambda)}{P(O|\lambda)}=\frac{P(i_t=q_i,i_{t+1}=q_j,O|\lambda)}{\sum_{i=1}^N\sum_{j=1}^NP(i_t=q_i,i_{t+1}=q_j,O|\lambda)} ξt(i,j)=P(O∣λ)P(it=qi,it+1=qj,O∣λ)=∑i=1N∑j=1NP(it=qi,it+1=qj,O∣λ)P(it=qi,it+1=qj,O∣λ)
而
P ( i t = q i , i t + 1 = q j , O ∣ λ ) = α t ( i ) a i j b j ( o t + 1 ) β t + 1 ( j ) P(i_t=q_i,i_{t+1}=q_j,O|\lambda)=\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j) P(it=qi,it+1=qj,O∣λ)=αt(i)aijbj(ot+1)βt+1(j)
所以
(26) ξ t ( i , j ) = α t ( i ) a i j b j ( o t + 1 ) β t + 1 ( j ) ∑ i = 1 N ∑ j = 1 N α t ( i ) a i j b j ( o t + 1 ) β t + 1 ( j ) \xi_t(i,j)=\frac{\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)}{\sum_{i=1}^N\sum_{j=1}^N\alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)} \tag{26} ξt(i,j)=∑i=1N∑j=1Nαt(i)aijbj(ot+1)βt+1(j)αt(i)aijbj(ot+1)βt+1(j)(26) - 将
γ
t
(
i
)
\gamma_t(i)
γt(i)和
ξ
t
(
i
,
j
)
\xi_t(i,j)
ξt(i,j)对各个时刻
t
t
t求和,可以得到一些有用的期望值:
(1) 在观测 O O O下状态 i i i出现的期望值
(27) ∑ t = 1 T γ t ( i ) \sum_{t=1}^T\gamma_t(i) \tag{27} t=1∑Tγt(i)(27)
(2) 在观测 O O O下由状态 i i i转移的期望值
(28) ∑ t = 1 T − 1 γ t ( i ) \sum_{t=1}^{T-1}\gamma_t(i) \tag{28} t=1∑T−1γt(i)(28)
(3) 在观测 O O O下由状态 i i i转移到状态 j j j的期望值
(29) ∑ t = 1 T − 1 ξ t ( i , j ) \sum_{t=1}^{T-1}\xi_t(i,j) \tag{29} t=1∑T−1ξt(i,j)(29)
学习算法
隐马尔科夫模型的学习,根据训练数据是包括观测序列和对应的状态序列还是只有观测序列,可以分别由监督学习与非监督学习实现。本节首先介绍监督学习算法,而后介绍非监督学习算法--Baum-Welch算法(也就是EM算法)。
- 监督学习方法
假设已给训练数据包含 S S S个长度相同的观测序列和对应的状态序列 { ( O 1 , I 1 ) , ( O 2 , I 2 ) , ⋯   , ( O S , I S ) } \{(O_1,I_1),(O_2,I_2),\cdots,(O_S,I_S)\} {(O1,I1),(O2,I2),⋯,(OS,IS)},那么可以用极大似然估计法来估计隐马尔科夫模型的参数。具体方法如下:
- 转移概率
a
i
j
a_{ij}
aij的估计
设样本中时刻 t t t处于状态 i i i时刻 t + 1 t+1 t+1转移到状态 j j j的频数为 A i j A_{ij} Aij,那么状态转移概率 a i j a_{ij} aij的估计是
(30) a i j ^ = A i j ∑ j = 1 N A i j , i = 1 , 2 , ⋯   , N ; j = 1 , 2 , ⋯   , N \hat{a_{ij}}=\frac{A_{ij}}{\sum_{j=1}^NA_{ij}}, i=1,2,\cdots,N; j=1,2,\cdots,N \tag{30} aij^=∑j=1NAijAij,i=1,2,⋯,N;j=1,2,⋯,N(30) - 观测概率
b
j
(
k
)
b_j(k)
bj(k)的估计
设样本中状态为 j j j并观测为 k k k的频数是 B j k B_{jk} Bjk,那么状态为 j j j观测为 k k k的概率 b j ( k ) b_j(k) bj(k)的估计是
(31) b j ( k ) ^ = B i j ∑ k = 1 M B i j , j = 1 , 2 , ⋯   , N ; k = 1 , 2 , ⋯   , M \hat{b_j(k)}=\frac{B_{ij}}{\sum_{k=1}^MB_{ij}}, j=1,2,\cdots,N; k=1,2,\cdots,M \tag{31} bj(k)^=∑k=1MBijBij,j=1,2,⋯,N;k=1,2,⋯,M(31) - 初始状态概率
π
i
\pi_i
πi的估计
π
i
^
\hat{\pi_i}
πi^为
S
S
S个样本中初始状态为
q
i
q_i
qi的频率。
由于监督学习需要使用训练数据,而人工标准训练数据往往代价很高,有时就会利用非监督学习的方法。
- Baum-Welch算法
假设给定训练数据只包含 S S S个长度为 T T T的观测序列 { O 1 , O 2 , ⋯   , O S } \{O_1,O_2,\cdots,O_S\} {O1,O2,⋯,OS}而没有对应的状态序列,目标是学习隐马尔科夫模型 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π)的参数。我们将观测序列数据看作观测数据 O O O,状态序列数据看作不可观测的隐数据 I I I,那么隐马尔科夫模型事实上是一个含有隐变量的概率模型
(32) P ( P ∣ λ ) = ∑ I P ( O ∣ I , λ ) P ( I ∣ λ ) P(P|\lambda)=\sum_IP(O|I,\lambda)P(I|\lambda) \tag{32} P(P∣λ)=I∑P(O∣I,λ)P(I∣λ)(32)
它的参数学习可以由EM算法实现。
- 确定完全数据的对数似然函数
所有观测数据写成 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT),所有隐数据写成 I = ( i 1 , i 2 , ⋯   , i T ) I=(i_1,i_2,\cdots,i_T) I=(i1,i2,⋯,iT),完全数据是 ( O , I ) = ( o 1 , o 2 , ⋯   , o T , i 1 , i 2 , ⋯   , i T ) (O,I)=(o_1,o_2,\cdots,o_T,i_1,i_2,\cdots,i_T) (O,I)=(o1,o2,⋯,oT,i1,i2,⋯,iT)。完全数据的对数似然函数是 l o g P ( O , I ∣ λ ) logP(O,I|\lambda) logP(O,I∣λ)。 - EM算法的E步:求Q函数
Q
(
λ
,
λ
ˉ
)
Q(\lambda,\bar\lambda)
Q(λ,λˉ)
(33) Q ( λ , λ ˉ ) = ∑ I l o g P ( O , I ∣ λ ) P ( O , I ∣ λ ˉ ) Q(\lambda,\bar\lambda)=\sum_IlogP(O,I|\lambda)P(O,I|\bar\lambda) \tag{33} Q(λ,λˉ)=I∑logP(O,I∣λ)P(O,I∣λˉ)(33)
其中, λ ˉ \bar\lambda λˉ是隐马尔科夫模型参数的当前估计值, λ \lambda λ是要极大化的隐马尔科夫模型参数。
P ( O , I ∣ λ ) = π i 1 b i 1 ( o 1 ) a i 1 i 2 b i 2 ( o 2 ) ⋯ a i T − 1 i T b i T ( o T ) P(O,I|\lambda)=\pi_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)\cdots a_{i_{T-1}i_T}b_{i_T}(o_T) P(O,I∣λ)=πi1bi1(o1)ai1i2bi2(o2)⋯aiT−1iTbiT(oT)
于是Q函数可以写成:
(34) Q ( λ , λ ˉ ) = ∑ I l o g π i 1 P ( O , I ∣ λ ˉ ) + ∑ I ( ∑ t = 1 T − 1 l o g a i t i t + 1 ) P ( O , I ∣ λ ˉ ) + ∑ I ( ∑ t = 1 T l o g b i t ( o t ) ) P ( O , I ∣ λ ˉ ) \begin{aligned} Q(\lambda,\bar\lambda)&=\sum_Ilog\pi_{i_1}P(O,I|\bar\lambda) \\ &+\sum_I\bigg(\sum_{t=1}^{T-1}loga_{i_ti_{t+1}}\bigg)P(O,I|\bar\lambda)+\sum_I\bigg(\sum_{t=1}^Tlogb_{i_t}(o_t)\bigg)P(O,I|\bar\lambda) \tag{34} \end{aligned} Q(λ,λˉ)=I∑logπi1P(O,I∣λˉ)+I∑(t=1∑T−1logaitit+1)P(O,I∣λˉ)+I∑(t=1∑Tlogbit(ot))P(O,I∣λˉ)(34)
式中求和都是对所有训练数据的序列长度 T T T进行的。 - EM算法的M步:极大化Q函数
Q
(
λ
,
λ
ˉ
)
Q(\lambda,\bar\lambda)
Q(λ,λˉ)求模型参数
A
A
A,
B
B
B,
π
\pi
π
由于要极大化的参数在式(34)中单独地出现在3个项中,所以只需对各项分别极大化。
(1) 式(34)的第1项可以写成:
∑ I l o g π i 0 P ( O , I ∣ λ ˉ ) = ∑ i = 1 N l o g π i P ( O , i 1 = i ∣ λ ˉ ) \sum_Ilog\pi_{i_0}P(O,I|\bar\lambda)=\sum_{i=1}^Nlog\pi_iP(O,i_1=i|\bar\lambda) I∑logπi0P(O,I∣λˉ)=i=1∑NlogπiP(O,i1=i∣λˉ)
注意到 π i \pi_i πi满足约束条件 ∑ i = 1 N π i = 1 \sum_{i=1}^N\pi_i=1 ∑i=1Nπi=1,利用拉格朗日乘子法,写出拉格朗日函数:
∑ i = 1 N l o g π i P ( O , i 1 = i ∣ λ ˉ ) + γ ( ∑ i = 1 N π i − 1 ) \sum_{i=1}^Nlog\pi_iP(O,i_1=i|\bar\lambda)+\gamma\bigg(\sum_{i=1}^N\pi_i-1\bigg) i=1∑NlogπiP(O,i1=i∣λˉ)+γ(i=1∑Nπi−1)
对其求偏导数并令结果为0
(35) ∂ ∂ π i [ ∑ i = 1 N l o g π i P ( O , i 1 = i ∣ λ ˉ ) + γ ( ∑ i = 1 N π i − 1 ) ] \frac{\partial}{\partial\pi_i}\bigg[\sum_{i=1}^Nlog\pi_iP(O,i_1=i|\bar\lambda)+\gamma\bigg(\sum_{i=1}^N\pi_i-1\bigg)\bigg] \tag{35} ∂πi∂[i=1∑NlogπiP(O,i1=i∣λˉ)+γ(i=1∑Nπi−1)](35)
得
P ( O , i 1 = i ∣ λ ˉ ) + γ π i = 0 P(O,i_1=i|\bar\lambda)+\gamma\pi_i=0 P(O,i1=i∣λˉ)+γπi=0
对 i i i求和得到 γ \gamma γ
γ = − P ( O ∣ λ ˉ ) \gamma=-P(O|\bar\lambda) γ=−P(O∣λˉ)
带入式(35)即得
(36) π i = P ( O , i 1 = i ∣ λ ˉ ) P ( O ∣ λ ˉ ) \pi_i=\frac{P(O,i_1=i|\bar\lambda)}{P(O|\bar\lambda)} \tag{36} πi=P(O∣λˉ)P(O,i1=i∣λˉ)(36)
(2) 式(34)的第2项可以写成
∑ I ( ∑ t = 1 T − 1 l o g a i t i t + 1 ) P ( O , I ∣ λ ˉ ) = ∑ i = 1 N ∑ j = 1 N ∑ t = 1 T − 1 l o g a i j P ( O , i t = i , i t + 1 = j ∣ λ ˉ ) \sum_I\bigg(\sum_{t=1}^{T-1}loga_{i_ti_{t+1}}\bigg)P(O,I|\bar\lambda)=\sum_{i=1}^N\sum_{j=1}^N \sum_{t=1}^{T-1}loga_{ij}P(O,i_t=i,i_{t+1}=j|\bar\lambda) I∑(t=1∑T−1logaitit+1)P(O,I∣λˉ)=i=1∑Nj=1∑Nt=1∑T−1logaijP(O,it=i,it+1=j∣λˉ)
类似第1项,应用具有约束条件 ∑ j = 1 N a i j = 1 \sum_{j=1}^Na_{ij}=1 ∑j=1Naij=1的拉格朗日乘子法可以求出
(37) a i j = ∑ t = 1 T − 1 P ( O , i t = i , i t + 1 = j ∣ λ ˉ ) ∑ t = 1 T − 1 P ( O , i t = i ∣ λ ˉ ) a_{ij}=\frac{\sum_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\bar\lambda)}{\sum_{t=1}^{T-1}P(O,i_t=i|\bar\lambda)} \tag{37} aij=∑t=1T−1P(O,it=i∣λˉ)∑t=1T−1P(O,it=i,it+1=j∣λˉ)(37)
(3) 式(34)的第3项为
∑ I ( ∑ t = 1 T l o g b i t ( o t ) ) P ( O , I ∣ λ ˉ ) = ∑ j = 1 N ∑ t = 1 T l o g b j ( o t ) P ( O , i t = j ∣ λ ˉ ) \sum_I\bigg(\sum_{t=1}^Tlogb_{i_t}(o_t)\bigg)P(O,I|\bar\lambda)=\sum_{j=1}^N\sum_{t=1}^Tlogb_j(o_t)P(O,i_t=j|\bar\lambda) I∑(t=1∑Tlogbit(ot))P(O,I∣λˉ)=j=1∑Nt=1∑Tlogbj(ot)P(O,it=j∣λˉ)
同样用拉格朗日乘子法,约束条件是 ∑ k = 1 M b j ( k ) = 1 \sum_{k=1}^Mb_j(k)=1 ∑k=1Mbj(k)=1。注意,只有在 o t = v k o_t=v_k ot=vk时 b j ( o t ) b_j(o_t) bj(ot)对 b j ( k ) b_j(k) bj(k)的偏导数才不为0,以 I ( o t = v k ) I(o_t=v_k) I(ot=vk)表示。求得
(38) b j ( k ) = ∑ t = 1 T P ( O , i t = j ∣ λ ˉ ) I ( o t = v k ) ∑ t = 1 T P ( O , i t = j ∣ λ ˉ ) b_j(k)=\frac{\sum_{t=1}^TP(O,i_t=j|\bar\lambda)I(o_t=v_k)}{\sum_{t=1}^TP(O,i_t=j|\bar\lambda)} \tag{38} bj(k)=∑t=1TP(O,it=j∣λˉ)∑t=1TP(O,it=j∣λˉ)I(ot=vk)(38)
- Baum-Welch模型参数估计公式
将式(36)~式(38)中的各概率分别用 γ t ( i ) \gamma_t(i) γt(i), ξ t ( i , j ) \xi_t(i,j) ξt(i,j)表示,则可将相应的公式写成:
(39) a i j = ∑ t = 1 T − 1 ξ t ( i , j ) ∑ t = 1 T − 1 γ t ( i ) a_{ij}=\frac{\sum_{t=1}^{T-1}\xi_t(i,j)}{\sum_{t=1}^{T-1}\gamma_t(i)} \tag{39} aij=∑t=1T−1γt(i)∑t=1T−1ξt(i,j)(39)
(40) b j ( k ) = ∑ t = 1 , o t = v k T γ t ( j ) ∑ t = 1 T γ t ( j ) b_j(k)=\frac{\sum_{t=1,o_t=v_k}^T\gamma_t(j)}{\sum_{t=1}^T\gamma_t(j)} \tag{40} bj(k)=∑t=1Tγt(j)∑t=1,ot=vkTγt(j)(40)
(41) π i = γ 1 ( i ) \pi_i=\gamma_1(i) \tag{41} πi=γ1(i)(41)
其中, γ t ( i ) \gamma_t(i) γt(i), ξ t ( i , j ) \xi_t(i,j) ξt(i,j)分别由式(24)及式(26)给出。式(39)~式(41)就是Baum-Welch算法,它是EM算法在隐马尔科夫模型学习中的具体实现,由Baum和Welch提出。
算法4(Baum-Welch算法)
输入:观测数据 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT);
输出:隐马尔科夫模型参数。
(1) 初始化
对 n = 0 n=0 n=0,选取 a i j ( 0 ) a_{ij}^{(0)} aij(0), b j ( k ) ( 0 ) b_j(k)^{(0)} bj(k)(0), π i ( 0 ) \pi_i^{(0)} πi(0),得到模型 λ ( 0 ) = ( A ( 0 ) , B ( 0 ) , π ( 0 ) ) \lambda^{(0)}=(A^{(0)},B^{(0)},\pi^{(0)}) λ(0)=(A(0),B(0),π(0))。
(2) 递推。对 n = 1 , 2 , ⋯ n=1,2,\cdots n=1,2,⋯,
a i j = ∑ t = 1 T − 1 ξ t ( i , j ) ∑ t = 1 T − 1 γ t ( i ) b j ( k ) = ∑ t = 1 , o t = v k T γ t ( j ) ∑ t = 1 T γ t ( j ) π i = γ 1 ( i ) a_{ij}=\frac{\sum_{t=1}^{T-1}\xi_t(i,j)}{\sum_{t=1}^{T-1}\gamma_t(i)} \\ b_j(k)=\frac{\sum_{t=1,o_t=v_k}^T\gamma_t(j)}{\sum_{t=1}^T\gamma_t(j)} \\ \pi_i=\gamma_1(i) aij=∑t=1T−1γt(i)∑t=1T−1ξt(i,j)bj(k)=∑t=1Tγt(j)∑t=1,ot=vkTγt(j)πi=γ1(i)
右端各值按观测 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT)和模型 λ ( n ) = ( A ( n ) , B ( n ) , π ( n ) ) \lambda^{(n)}=(A^{(n)},B^{(n)},\pi^{(n)}) λ(n)=(A(n),B(n),π(n))计算。式中 γ t ( i ) \gamma_t(i) γt(i), ξ t ( i , j ) \xi_t(i,j) ξt(i,j)由式(24)和式(26)给出。
(3) 终止。得到模型参数 λ ( n + 1 ) = ( A ( n + 1 ) , B ( n + 1 ) , π ( n + 1 ) ) \lambda^{(n+1)}=(A^{(n+1)},B^{(n+1)},\pi^{(n+1)}) λ(n+1)=(A(n+1),B(n+1),π(n+1))。
预测算法
下面介绍隐马尔科夫模型预测的两种算法:近似算法与维特比算法。
-
近似算法
近似算法的想法是,在每个时刻 t t t选择在该时刻最有可能出现的状态 i t ∗ i_t^* it∗,从而得到一个状态序列 I ∗ = ( i 1 ∗ , i 2 ∗ , ⋯   , i T ∗ ) I^*=(i_1^*,i_2^*,\cdots,i_T^*) I∗=(i1∗,i2∗,⋯,iT∗),将它作为预测的结果。
给定隐马尔科夫模型 λ \lambda λ和观测序列 O O O,在时刻 t t t处于状态 q i q_i qi的概率 γ t ( i ) \gamma_t(i) γt(i)是
(42) γ t ( i ) = α t ( i ) β t ( i ) P ( O ∣ λ ) = α t ( i ) β t ( i ) ∑ j = 1 N α t ( j ) β t ( j ) \gamma_t(i)=\frac{\alpha_t(i)\beta_t(i)}{P(O|\lambda)}=\frac{\alpha_t(i)\beta_t(i)}{\sum_{j=1}^N\alpha_t(j)\beta_t(j)} \tag{42} γt(i)=P(O∣λ)αt(i)βt(i)=∑j=1Nαt(j)βt(j)αt(i)βt(i)(42)
在每一时刻 t t t最有可能的状态 i t ∗ i_t^* it∗是
(43) i t ∗ = a r g max 1 ≤ i ≤ N [ γ t ( i ) ] , t = 1 , 2 , ⋯   , T i_t^*=arg\max_{1\leq i\leq N}[\gamma_t(i)], t=1,2,\cdots,T \tag{43} it∗=arg1≤i≤Nmax[γt(i)],t=1,2,⋯,T(43)
从而得到状态序列 I ∗ = ( i 1 ∗ , i 2 ∗ , ⋯   , i T ∗ ) I^*=(i_1^*,i_2^*,\cdots,i_T^*) I∗=(i1∗,i2∗,⋯,iT∗)。
近似算法的优点是计算简单,其缺点是不能保证预测的状态序列整体是最有可能的状态序列,因为预测的状态序列可能有实际不发生的部分。事实上,上述方法得到的状态序列中有可能存在转移概率为0的相邻专题,即对某些 i i i, j j j, a i j = 0 a_{ij}=0 aij=0时。尽管如此,近似算法仍然是有用的。 -
维特比算法
维特比算法实际是动态规划解隐马尔科夫模型预测问题,即用动态规划求概率最大路径(最优路径)。这时一条路径对应着一个状态序列。
根据动态规划原理,最优路径具有这样的特性:如果最优路径在时刻 t t t通过结点 i t ∗ i_t^* it∗,那么这一路径从结点 i t ∗ i_t^* it∗到终点 i T ∗ i_T^* iT∗的部分路径,对于从 i t ∗ i_t^* it∗到 i T ∗ i_T^* iT∗的所有可能的部分路径来说,必须是最优的。因为加入不是这样,那么从 i t ∗ i_t^* it∗到 i T ∗ i_T^* iT∗就有另一条更好的部分路径存在,如果把它和从 i 1 ∗ i_1^* i1∗到 i t ∗ i_t^* it∗的部分路径连接起来,就会形成一条比原来的路径更优的路径,这是矛盾的,依据这一原理,我们只需从时刻 t = 1 t=1 t=1开始,递推地计算在时刻 t t t状态为 i i i的各条部分路径的最大概率,直到得到时刻 t = T t=T t=T状态为 i i i的各条路径的最大概率。时刻 t = T t=T t=T的最大概率即为最优路径的概率 P ∗ P^* P∗,最优路径的终结点 i T ∗ i_T^* iT∗也同时得到。之后,为了找出最优路径的各个结点,从终结点 i T ∗ i_T^* iT∗开始,由后向前逐步求得结点 i T − 1 ∗ , ⋯   , i 1 ∗ i_{T-1}^*,\cdots,i_1^* iT−1∗,⋯,i1∗,得到最优路径 I ∗ = ( i 1 ∗ , i 2 ∗ , ⋯   , i T ∗ ) I^*=(i_1^*,i_2^*,\cdots,i_T^*) I∗=(i1∗,i2∗,⋯,iT∗)。这就是维特比算法。
首先导入两个变量 δ \delta δ和 ψ \psi ψ。定义在时刻 t t t状态为 i i i的所有单个路径 ( i 1 , i 2 , ⋯   , i t ) (i_1,i_2,\cdots,i_t) (i1,i2,⋯,it)中概率最大值为
(44) δ t ( i ) = max i 1 , i 2 , ⋯   , i t − 1 P ( i t = i , i t − 1 , ⋯   , i 1 , o t , ⋯   , o 1 ∣ λ ) , i = 1 , 2 , ⋯   , N \delta_t(i)=\max_{i_1,i_2,\cdots,i_{t-1}}P(i_t=i,i_{t-1},\cdots,i_1,o_t,\cdots,o_1|\lambda), i=1,2,\cdots,N \tag{44} δt(i)=i1,i2,⋯,it−1maxP(it=i,it−1,⋯,i1,ot,⋯,o1∣λ),i=1,2,⋯,N(44)
由定义可得变量 δ \delta δ的递推公式:
(45) δ t + 1 ( i ) = max i 1 , i 2 , ⋯   , i t P ( i t + 1 = i , i t , ⋯   , i 1 , o t + 1 , ⋯   , o 1 ∣ λ ) = max 1 ≤ j ≤ N [ δ t ( j ) a j i ] b i ( o t + 1 ) , i = 1 , 2 , ⋯   , N ; t = 1 , 2 , ⋯   , T − 1 \begin{aligned} \delta_{t+1}(i)&=\max_{i_1,i_2,\cdots,i_t}P(i_{t+1}=i,i_t,\cdots,i_1,o_{t+1},\cdots,o_1|\lambda) \\ &=\max_{1\leq j\leq N}[\delta_t(j)a_{ji}]b_i(o_{t+1}), i=1,2,\cdots,N; t=1,2,\cdots,T-1 \tag{45} \end{aligned} δt+1(i)=i1,i2,⋯,itmaxP(it+1=i,it,⋯,i1,ot+1,⋯,o1∣λ)=1≤j≤Nmax[δt(j)aji]bi(ot+1),i=1,2,⋯,N;t=1,2,⋯,T−1(45)
定义在时刻 t t t状态为 i i i的所有单个路径 ( i 1 , i 2 , ⋯   , i t − 1 , i ) (i_1,i_2,\cdots,i_{t-1},i) (i1,i2,⋯,it−1,i)中概率最大的路径的第 t − 1 t-1 t−1个结点为
(46) ψ t ( i ) = a r g max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] , i = 1 , 2 , ⋯   , N \psi_t(i)=arg\max_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}], i=1,2,\cdots,N \tag{46} ψt(i)=arg1≤j≤Nmax[δt−1(j)aji],i=1,2,⋯,N(46)
下面介绍维特比算法。
算法5(维特比算法)
输入:模型 λ = ( A , B , π ) \lambda=(A,B,\pi) λ=(A,B,π)和观测 O = ( o 1 , o 2 , ⋯   , o T ) O=(o_1,o_2,\cdots,o_T) O=(o1,o2,⋯,oT);
输出:最优路径 I ∗ = ( i 1 ∗ , i 2 ∗ , ⋯   , i T ∗ ) I^*=(i_1^*,i_2^*,\cdots,i_T^*) I∗=(i1∗,i2∗,⋯,iT∗)。
(1) 初始化
δ 1 ( i ) = π i b i ( o 1 ) , i = 1 , 2 , ⋯   , N ψ 1 ( i ) = 0 , i = 1 , 2 , ⋯   , N \delta_1(i)=\pi_ib_i(o_1), i=1,2,\cdots,N \\ \psi_1(i)=0, i=1,2,\cdots,N δ1(i)=πibi(o1),i=1,2,⋯,Nψ1(i)=0,i=1,2,⋯,N
(2) 递推。对 t = 2 , 3 , ⋯   , T t=2,3,\cdots,T t=2,3,⋯,T
δ t ( i ) = max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] b i ( o t ) , i = 1 , 2 , ⋯   , N ψ t ( i ) = a r g max 1 ≤ j ≤ N [ δ t − 1 ( j ) a j i ] , i = 1 , 2 , ⋯   , N \delta_t(i)=\max_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}]b_i(o_t), i=1,2,\cdots,N \\ \psi_t(i)=arg\max_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}], i=1,2,\cdots,N δt(i)=1≤j≤Nmax[δt−1(j)aji]bi(ot),i=1,2,⋯,Nψt(i)=arg1≤j≤Nmax[δt−1(j)aji],i=1,2,⋯,N
(3) 终止
P ∗ = max 1 ≤ i ≤ N δ T ( i ) i T ∗ = a r g max 1 ≤ i ≤ N [ δ T ( i ) ] P^*=\max_{1\leq i\leq N}\delta_T(i) \\ i_T^*=arg\max_{1\leq i\leq N}[\delta_T(i)] P∗=1≤i≤NmaxδT(i)iT∗=arg1≤i≤Nmax[δT(i)]
(4) 最优路径回溯。对 t = T − 1 , T − 2 , ⋯   , 1 t=T-1,T-2,\cdots,1 t=T−1,T−2,⋯,1
i t ∗ = ψ t + 1 ( i t + 1 ∗ ) i_t^*=\psi_{t+1}(i_{t+1}^*) it∗=ψt+1(it+1∗)
求得最优路径 I ∗ = ( i 1 ∗ , i 2 ∗ , ⋯   , i T ∗ ) I^*=(i_1^*,i_2^*,\cdots,i_T^*) I∗=(i1∗,i2∗,⋯,iT∗)。