馬爾科夫模型(轉(zhuǎn)載)

1、隱馬爾可夫模型（一）馬爾可夫模型馬爾可夫模型（Markov Model)描述了一類隨機變量隨時間而變化的隨機函數(shù)。考察一個狀態(tài)序列（此時隨機變量為狀態(tài)值），這些狀態(tài)并不是相互獨立的，每個狀態(tài)的值依賴于序列中此狀態(tài)之前的狀態(tài)。數(shù)學描述：一個系統(tǒng)由N個狀態(tài)S= s1,s2,.sn,隨著時間的推移，該系統(tǒng)從一個狀態(tài)轉(zhuǎn)換成另一個狀態(tài)。Q= q1,q2,.qn為一個狀態(tài)序列，qiS,在t時刻的狀態(tài)為qt,對該系統(tǒng)的描述要給出當前時刻t所處的狀態(tài)st，和之前的狀態(tài)s1,s2,.st, 則t時刻位于狀態(tài)qt的概率為：P(qt=st|q1=s1,q2=s2,.qt-1=st-1)。這樣的模型叫馬爾可夫模型。特

2、殊狀態(tài)下，當前時刻的狀態(tài)只決定于前一時刻的狀態(tài)叫一階馬爾可夫模型，即P(qt=si|q1=s1,q2=s2,.qt-1=sj) =P(qt=si|qt-1=sj)。狀態(tài)之間的轉(zhuǎn)化表示為aij,aij=P(qt=sj|qt-1=si),其表示由狀態(tài)i轉(zhuǎn)移到狀態(tài)j的概率。其必須滿足兩個條件：   1.aij 0    2.=1對于有N個狀態(tài)的一階馬爾科夫模型，每個狀態(tài)可以轉(zhuǎn)移到另一個狀態(tài)（包括自己），則共有N2次狀態(tài)轉(zhuǎn)移，可以用狀態(tài)轉(zhuǎn)移矩陣表示。例如：一段文字中名詞、動詞、形容詞出現(xiàn)的情況可以用有3個狀態(tài)的y一階馬爾科夫模型M表示：

3、60;                 狀態(tài)s1:名詞         狀態(tài)s2:動詞       狀態(tài)s3:形容詞狀態(tài)轉(zhuǎn)移矩陣：     s1            s2     

4、;     s3                     A=  則狀態(tài)序列O=“名動形名”（假定第一個詞為名詞）的概率為：          P(O|M) = P(s1,s2,s3,s4 = P(s1)*p(s2|s1)p(s3|s2)p(s1|s

5、3)                                             =p(s1)*a12*a23*a31  

6、60;                                          =1*0.5*0.2*0.4      

7、                                       =0.04隱馬爾可夫模型（二）隱馬爾可夫模型的構(gòu)成     在馬爾可夫模型中，每一個狀態(tài)都是可觀察的序

8、列，是狀態(tài)關(guān)于時間的隨機過程，也成為可視馬爾可夫模型（Visible Markov Model,VMM）。隱馬爾科夫模型（Hidden Markov Model,HMM）中的狀態(tài)是不可見的，我們可以看到的是狀態(tài)表現(xiàn)出來的觀察值和狀態(tài)的概率函數(shù)。在隱馬模型中，觀察值是關(guān)于狀態(tài)的隨機過程，而狀態(tài)是關(guān)于時間的隨機過程，因此隱馬模型是一個雙重隨機過程。    當考慮潛在事件隨機生成表面事件時，可以用HMM解決。     舉個例子，說明隱馬模型：     有4個暗箱，放在暗處，每個箱子里有3種不

9、用顏色的球（紅、橙、藍），從箱子往外拿球是有一定規(guī)律的，現(xiàn)在工作人員從暗處的箱子中取球，去了5次，我們看到額觀察序列是：紅藍藍橙紅。這個過程就是一個隱馬模型。暗處的箱子表示狀態(tài)，箱子的個數(shù)表示狀態(tài)的個數(shù)，球的顏色表示狀態(tài)的輸出值，球的顏色個數(shù)表示狀態(tài)輸出觀察狀態(tài)的個數(shù)，從一個箱子轉(zhuǎn)換成另一個箱子表示狀態(tài)轉(zhuǎn)換，從暗處箱子中取出的球的觀察顏色表示狀態(tài)的輸出序列。      因此可以歸納隱馬模型的5個組成狀態(tài)：      （1）模型中的狀態(tài)個數(shù)N(例子中的箱子個數(shù))；   

10、;   （2）每個狀態(tài)的可以輸出的不同觀測值M（例子中的球的顏色數(shù)目）；      （3）狀態(tài)轉(zhuǎn)移矩陣A= aij(例子中aij表示從第i個箱子轉(zhuǎn)移到第j個箱子的概率），其中aij滿足條件：            I.  aij0， 1i,jN            II. aij= P(qt=sj

11、|qt-1=si),              III.   =1        （4）發(fā)射矩陣B=bj(k)，即從狀態(tài)sj觀察到符號vk的概率分布矩陣.(bj(k)在例子中表示從的j個箱子中拿出第k個顏色的概率），其中bj(k)滿足條件：          

12、  I.  bj(k)0， 1jN; 1kM            II. bj(k)= P(Ot=vk|qt=sj),             III. =1         （5）初始狀態(tài)概率分布 = j.(例子中一開始到第j個箱子的概率），其中j滿足條件

13、：            I.  j(k)0， 1jN            II. j= P(q1=sj),            III. =1        

14、60;  一般，一個HMM即為五元組=N,M,A,B,，為了簡便，常簡記為三元組=A,B,。          HMM有三個基本問題：         （1）評估問題：給定一個觀察序列O=O1O2.OT和模型=A,B,，如何快速地計算給定模型的條件下，觀察序列O=O1O2.OT的概率，即P(O|）？         （2）解碼問題：給定一

15、個觀察序列O=O1O2.OT和模型=A,B,，如何快速地選擇在給定模型的條件下在一定意義下”最優(yōu)“的狀態(tài)序列Q=Q1Q2.QT，是該狀態(tài)序列”最好地"解釋觀察序列？         （3）學習問題：給定一個觀察序列O=O1O2.OT，如何調(diào)整參數(shù)=A,B,，使得P(O|M)最大？          針對HMM的三個基本問題，相應的算法是：      

16、0;  （1）評估問題：前向后向算法         （2）解碼問題：維特比算法（Viterbi）         （3）學習問題：前向后向算法（BAUM-WELCH）。隱馬爾可夫模型（三）隱馬爾可夫模型的評估問題(前向算法）      隱馬模型的評估問題即，在已知一個觀察序列O=O1O2.OT，和模型=（A,B,的條件下，觀察序列O的概率，即P(O|  &

17、#160;                      如果窮盡所有的狀態(tài)組合，即S1S1.S1, S1S1.S2, S1S1.S3, ., S3S3.S3。這樣的話t1時刻有N個狀態(tài)，t2時刻有N個狀態(tài).tT時刻有N個狀態(tài)，這樣的話一共有N*N*.*N= NT種組合，時間復雜度為O(NT),計算時，就會出現(xiàn)“指數(shù)爆炸”，當T很大時，簡直無法計算這個值。為解決這一問題，Baum提

18、出了前向算法。      歸納過程      首先引入前向變量t(i):在時間t時刻，HMM輸出序列為O1O2.OT,在第t時刻位于狀態(tài)si的概率。      當T=1時，輸出序列為O1,此時計算概率為P(O1|）：假設有三個狀態(tài)（如下圖）1、2、3，輸出序列為O1，有三種可能一是狀態(tài)1發(fā)出，二是從狀態(tài)2發(fā)出，三是從狀態(tài)3發(fā)出。另外從狀態(tài)1發(fā)出觀察值O1得概率為b1(O1),從狀態(tài)2發(fā)出觀察值O1得概率為b2(O1),從狀態(tài)3發(fā)出觀察值

19、O1得概率為b3(O1)。因此可以算出     P(O1|）= 1*b1(O1)+2*b2(O1) +  3*b3(O1)= 1(1) + 1(2) + 1(3)                              

20、;             當T=2時，輸出序列為O1O2,此時計算概率為P(O1O2|）：假設有三個狀態(tài)（如下圖）1、2、3，輸出序列為O1，有三種可能一是狀態(tài)1發(fā)出，二是從狀態(tài)2發(fā)出，三是從狀態(tài)3發(fā)出。另外從狀態(tài)1發(fā)出觀察值O2得概率為b1(O2),從狀態(tài)2發(fā)出觀察值O2得概率為b2(O2),從狀態(tài)3發(fā)出觀察值O2得概率為b3(O2)。      要是從狀態(tài)1發(fā)出觀察值O2，可能從第一時刻的1、2或3狀態(tài)裝換過來，要是從狀

21、態(tài)1轉(zhuǎn)換過來，概率為1(1)*a11*b1(O2),要是從狀態(tài)2轉(zhuǎn)換過來，概率為1(2)*a21*b1(O2),要是從狀態(tài)3轉(zhuǎn)換過來，概率為1(3)*a31*b1(O2),因此     P(O1O2,q2=s1|）= 1(1)*a11*b1(O2)  + 1(2)*a21*b1(O2) + 1(3)*a31*b1(O2)=2(1)               

22、0;                           同理：P(O1O2,q2=s1|）= 1(1)*a12*b2(O2)  + 1(2)*a22*b2(O2) + 1(3)*a32*b2(O2)=2(2)           

23、;    P(O1O2,q2=s1|）= 1(1)*a13*b1(O2)  + 1(2)*a23*b3(O2) + 1(3)*a33*b3(O2)=2(3)     所以：P(O1O2|）=P(O1O2,q2=s1|）+ P(O1O2,q2=s1|）+ P(O1O2,q2=s1|）                    

24、;         =2(1) + 2(2) + 2(3)      以此類推。      前向算法       step1 初始化：1(i) = i*bi(O1), 1iN       step2 歸納計算:      

25、;                            step3 終結(jié)：                     

26、P(O|）=      時間復雜度      計算某時刻的某個狀態(tài)的前向變量需要看前一時刻的N個狀態(tài)，此時時間復雜度為O(N),每個時刻有N個狀態(tài)，此時時間復雜度為N*O(N)=O(N2),又有T個時刻，所以時間復雜度為T*O(N2)=O(N2T)。      程序例證             &#

27、160; 前向算法計算P(O|M)：        step1：1(1) =1*b1(red)=0.2*0.5=0.1          1(2)=2*b2(red)=0.4*0.4= 0.16          1(3)=3*b3(red)=0.4*0.7=0.21        ste

28、p2：2(1)=1(1)*a11*b1(white) + 1(2)*a21*b1(white) + 1(3)*a31*b1(white)                     .        step3:P(O|M) = 3(1)+3(2)+3(3)      &#

29、160; 程序代碼#include <stdio.h>#include <stdlib.h>#include <string.h>int main() float a33 = 0.5,0.2,0.3,0.3,0.5,0.2,0.2,0.3,0.5; float b32 = 0.5,0.5,0.4,0.6,0.7,0.3; float alpha43; int i,j,k, count = 1; /output list int list4 = 0,1,0,1; /step1:Initialization alpha00 = 0.2 * 0.5; a

30、lpha01 = 0.4 * 0.4; alpha02 = 0.4 * 0.7; /step2:iteration for (i=1; i<=3; i+) for(j=0; j<=2; j+) alphaij = 0; for(k=0; k<=2; k+) alphaij += alphai-1k * akj * bjlistcount; count += 1; for (i=0; i<=3; i+) for(j=0; j<=2; j+) printf("a%d%d=%fn",i+1,j+1,alphaij); /step3:end print

31、f("Forward:%fn", alpha30+alpha31+alpha32); return 0;     運行結(jié)果                  隱馬爾可夫模型（四）隱馬爾可夫模型的評估問題(后向算法）對于HMM的評估問題，利用動態(tài)規(guī)劃可以用前向算法，從前到后算出前向變量；也可以采用后向算法，從后到前算出后向變量。先介紹后向變量t(i)

32、:給定模型=（A,B,），并且在時間 時刻t 狀態(tài)為si 的前提下，輸出序列為Ot+1Ot+2.OT的概率，即                                    t(i)=P(Ot+1Ot+

33、2.OT|qt=si,)歸納過程    假設仍然有3個狀態(tài)                      當t=T時，按照定義：時間t  狀態(tài)qT 輸出為OT+1.的概率，從T+1開始的輸出是不存在的（因為T時刻是終止終止狀態(tài)），即T之后是空，是個必然事件，因此t(i)=1,11N    當t=T-1時，

34、0;                                          T-1(i)=P(OT|qT-1=si,) = ai1*b1（OT)*T(1)  

35、+  ai2*b2（OT)*T(2)  +  ai3*b3（OT)*T(3)      .    當t=1時，       1(1)=P(O2O3.OT|q2=s1,) = a11*b1（O2)*2(1) + a12*b2（O2)*2(2) + a13*b3（O2)*2(3)       1(2)=P(O2O3.OT|q2=s1,) = a21*b1（O2)*2(1)

36、 + a22*b2（O2)*2(2) + a23*b3（O2)*2(3)       1(3)=P(O2O3.OT|q2=s1,) = a31*b1（O2)*2(1) + a32*b2（O2)*2(2) + a33*b3（O2)*2(3)       P(O1O2.OT|) =                

37、                 =                                =  后向算法

38、60;   step1 初始化：T(i)=1, 11N    step2 歸納計算：                       1tT-1, 1iN    step3 求終結(jié)和：             &#

39、160;     P(O|）=  時間復雜度    計算某時刻在某個狀態(tài)下的后向變量需要看后一時刻的N個狀態(tài)，此時時間復雜度為O(N),每個時刻有N個狀態(tài)，此時時間復雜度為N*O(N)=O(N2),又有T個時刻，所以時間復雜度為T*O(N2)=O(N2T)。程序例證              后向算法    計算P(O|M)：    st

40、ep1：4(1) = 1          4(2) = 1          4(3) = 1    step2：3(1) = 4(1)*a11*b1(white) + 4(2)*a12*b2(white) + 4(3)*a13*b3(white)            

41、60;        .    step3:P(O|M) = 1*1(1)*b1(O1) + 2*1(2)*b2(O1) + 3*1(3)*b3(O1)程序代碼#include <stdio.h>#include <stdlib.h>#include <string.h>int main() float a33 = 0.5,0.2,0.3,0.3,0.5,0.2,0.2,0.3,0.5; float b32 = 0.5,0.5,0.4,0.6,0.7,0.3; floa

42、t result43; int list4 = 0,1,0,1; result30 = 1; result31 = 1; result32 = 1; int i,j,k, count = 3; for (i=2; i>=0; i-) for(j=0; j<=2; j+) resultij = 0; for(k=0; k<=2; k+) resultij += resulti+1k * ajk * bklistcount; count -= 1; for (i=0; i<=3; i+) for(j=0; j<=2; j+) printf("b%d%d =

43、%fn",i+1,j+1,resultij); printf("backward:%fn", result00*0.2*0.5+result01*0.4*0.4+result02*0.4*0.7); return 0;運行結(jié)果             隱馬爾可夫模型（五）隱馬爾可夫模型的解碼問題(維特比算法）閱讀目錄· HMM解碼問題· 維特比算法· 時間復雜度· 程序例證回到頂部HMM解碼問題&

44、#160;     給定一個觀察序列O=O1O2.OT,和模型=（A,B,），如何快速有效地選擇在一定意義下“最優(yōu)”的狀態(tài)序列Q=q1q2.qT，使該狀態(tài)最好地解釋觀察序列。                  一種想法是求出每個狀態(tài)的概率rt(i)最大(rt(i)=P(qt=si,O|）)，記q't(i)=argQmax(rt(i)，但是這樣做，忽略了狀態(tài)之間的關(guān)系，很

45、可能兩個狀態(tài)之間的概率為0，即aq't(i)q't+1(i)=0,這樣求得的“最優(yōu)”狀態(tài)序列是不合法的。      為防止狀態(tài)之間轉(zhuǎn)移概率為0（斷續(xù)問題），換一種思路，不是求單個狀態(tài)求得最大值，而是求得整個狀態(tài)序列最大值，即求                         

46、0;         Q'= argQmaxP(Q|O,）      此時用維特比算法，先定義下維特比變量t(i):在時間t，HMM沿著一條路徑到達狀態(tài)si，并輸出觀察序列O=O1O2.Ot的最大概率:      t(i)=max P(q1q2.qt=si,O1O2.Ot|)           &#

47、160;                           t                      t+1 &#

48、160;    上圖中，對于從t時刻三個到 t+1時刻的狀態(tài)1，到底取狀態(tài)1,2還是3，不是看單獨狀態(tài)1,2還是3的概率，而是看在狀態(tài)1,2,3各自的維特比變量值乘以相應的狀態(tài)轉(zhuǎn)換概率，從中選出最大值，假設2時最大，那么記下t+1時刻狀態(tài)1之前的路徑是t時刻的狀態(tài)2，以此類推。      t(i)的遞歸關(guān)系式: t+1(i)=maxj t(j)*aji*bi(Ot+1),為了記憶路徑，定義路徑變量t(i)，記錄該路徑上的狀態(tài)si的前一個狀態(tài)。回到頂部維特比算法    

49、;  step1 初始化：              t(i) = i*bi(O1), 1iN              t(i) = 0      step2 歸納計算：  t(i)=max1jN t-1(j)*aji*bi(Ot),2tT;1iN

50、             記憶路徑   t(i) = arg max1jNt-1(j)*aji*bi(Ot)      step3 終結(jié):            QT' = arg max1iN T(i)      

51、60;     P'(QT') = max1iN T(i) step4 路徑回溯:             qt'=t+1(qt+1') , t=T-1,T-2.1回到頂部時間復雜度      計算某時刻的某個狀態(tài)的前向變量需要比較前一時刻的N個狀態(tài)，此時時間復雜度為O(N),每個時刻有N個狀態(tài)，此時時間復雜度為N*O(N)=

52、O(N2),又有T個時刻，所以時間復雜度為T*O(N2)=O(N2T)?；氐巾敳砍绦蚶C                     step1 初始化：1(1) = 0.2*0.5=0.1 ，1(2) = 0.4*0.4=0.16， 1(3) = 0.4*0.7=0.21       step2 歸納計算：2(1) =m

53、ax0.1*0.5,0.16*0.3,0.21*0.2*0.6       .      step3 終結(jié)：最佳路徑是4(1)4(2)4(3)最大的一個對應的狀態(tài)      step4 回溯：從最后一個狀態(tài)往回返程序代碼 #include <stdio.h>#include <stdlib.h>#include <string.h>int main() float a33

54、 = 0.5,0.2,0.3,0.3,0.5,0.2,0.2,0.3,0.5; float b32 = 0.5,0.5,0.4,0.6,0.7,0.3; float result43; int list4 = 0,1,0,1; int max43; float tmp; /step1:Initialization result00 = 0.2*0.5; result01 = 0.4*0.4; result02 = 0.4*0.7; int i,j,k, count = 1, max_node; float max_v; /step2:歸納運算 for (i=1; i<=3; i+) fo

55、r(j=0; j<=2; j+) tmp = resulti-10 * a0j * bjlistcount; maxij = 0; for(k=1; k<=2; k+) if(resulti-1k * akj * bjlistcount > tmp) tmp = resulti-1k * akj* bjlistcount; maxij = k; resultij = tmp; max_v = result30; max_node = 0; for (k=1; k<=2; k+) if(result3k > max_v) max_v = result3k; max_

56、node = k; count += 1; /step3:終結(jié) for (i=0; i<=3; i+) for(j=0; j<=2; j+) printf("%d %d %fn",i+1,j+1,resultij); printf("Pmax=%fn", max_v); printf("step4:%d n", max_node+1); /step4:回溯 for(k=3; k>=1; k-) printf("step%d:%d n",k, maxkmax_node+1); max_node =

57、maxkmax_node; return 0;  運行結(jié)果        最終的序列是3 2 2 2隱馬爾可夫模型（六）隱馬爾可夫模型的評估問題(前向后向相結(jié)合算法）重新回顧：    前向變量t(i):在時刻t,在已知模型=（A,B,）的條件下，狀態(tài)處于si,輸出序列為O102.Ot,前向變量為t(i)    后向變量t(i):在時刻t,在已知模型=（A,B,）和狀態(tài)處于si的條件下，輸出序列為Ot+1Ot+2.OT,后向變量為t(i)公式推導：  

58、60; P(O,qt=si|） = P(O1O2.OT, qt=si|）                         =P(O1O2.Ot, qt=si,Ot+1Ot+2.OT|)              

59、;           =P(O1O2.Ot, qt=si|) * P(Ot+1Ot+2.OT|O1O2.Ot, qt=si,)                         =P(O1O2.Ot, qt=si|) * P(Ot

60、+1Ot+2.OT|qt=si,)                         =t(i) *  t(i)     P(O|）=案例分析：          分析：   

61、60;    P(q4=s3|O,M) =  P(q4=s3, O|M)/P(O|M)                            = P(O,q4=s3|M)/P(O|M)       &

62、#160;                    = 4(3) *  4(3)/       程序： #include <stdio.h>#include <stdlib.h>#include <string.h>int main() float a33 = 0.5,0.2,0.3,0.3,0

63、.5,0.2,0.2,0.3,0.5; float b32 = 0.5,0.5,0.4,0.6,0.7,0.3; float result_b83; float result_f83; float result, result_t; int list8 = 0,1,0,0,1,0,1,1; result_b70 = 1; result_b71 = 1; result_b72 = 1; result_f00 = 0.2 * 0.5; result_f01 = 0.4 * 0.4; result_f02 = 0.4 * 0.7; /Backward int i,j,k, count = 7; fo

64、r (i=6; i>=0; i-) for(j=0; j<=2; j+) result_bij = 0; for(k=0; k<=2; k+) result_bij += result_bi+1k * ajk * bklistcount; count -= 1; for (i=0; i<=7; i+) for(j=0; j<=2; j+) printf("b%d%d= %fn",i+1,j+1, result_bij); printf("Backward:%fn", result_b00*0.2*0.5+result_b01

65、*0.4*0.4+result_b02*0.4*0.7); /Forward count = 1; for (i=1; i<=7; i+) for(j=0; j<=2; j+) result_fij = 0; for(k=0; k<=2; k+) result_fij += result_fi-1k * akj * bjlistcount; count += 1; for (i=0; i<=7; i+) for(j=0; j<=2; j+) printf("a%d%d= %fn", i+1, j+1, result_fij); result =

66、 result_f70 + result_f71 + result_f72; printf("Forward: %fn", result); result_t = 0; for (i=0; i<=2; i+) result_t += result_f3i * result_b3i; printf("Result:%fn", result_f32*result_b32/result_t); return 0;        運行結(jié)果