兩種馬爾可夫鏈狀態(tài)轉(zhuǎn)移概率矩陣的估計與比較

陳 雯，呂王勇，2，李思奇，代 娟，鄧 柙

(1.四川師范大學(xué) 數(shù)學(xué)科學(xué)學(xué)院， 成都 610068； 2.可視化計算與虛擬現(xiàn)實四川省重點實驗室， 成都 610068)

馬爾可夫鏈模型是以俄國數(shù)學(xué)家A.A.Markov命名的一種動態(tài)隨機模型，通過分析隨機變量現(xiàn)實的運動情況來預(yù)見這些變量未來的運動情況[1]。目前，馬爾可夫鏈模型在自然科學(xué)、工程技術(shù)、社會科學(xué)、經(jīng)濟研究等領(lǐng)域有著廣泛的應(yīng)用[2-5]。馬爾可夫模型通過研究系統(tǒng)對象不同狀態(tài)的初始概率和狀態(tài)之間的轉(zhuǎn)移概率來進行預(yù)測，因此，狀態(tài)轉(zhuǎn)移概率的確定稱為馬爾可夫模型預(yù)測的關(guān)鍵。對于馬爾可夫鏈的轉(zhuǎn)移概率矩陣，張二艷等[6]提出了用應(yīng)用統(tǒng)計方法進行估計；李成燮[7]和李永立[8]分別利用線性方程組和非線性優(yōu)化法進行估計；文士發(fā)等[9]通過構(gòu)造優(yōu)化模型并將其轉(zhuǎn)化為線性規(guī)劃模型來進行估計；聶篤忠等[10]介紹了一種迭代求解方法來確定狀態(tài)概率轉(zhuǎn)移矩陣；唐小我等[11]在最小二乘指標下導(dǎo)出了馬爾科夫鏈狀態(tài)轉(zhuǎn)移概率矩陣的估計公式；Yue Dequan等[12]通過迭代公式得到轉(zhuǎn)移概率矩陣的顯式解析表達式；賀思輝[13]等以加權(quán)平均的方式計算經(jīng)驗轉(zhuǎn)移概率矩陣。

在馬爾可夫鏈預(yù)測方法中，重要的是確定狀態(tài)轉(zhuǎn)移概率矩陣。通過驗證馬爾可夫鏈任意等間隔子序列也是馬爾可夫鏈，提出了2種間隔為k的轉(zhuǎn)移概率矩陣的估計方法：第1種方法是利用間隔為1的一步轉(zhuǎn)移概率矩陣和C-K方程計算間隔為k的轉(zhuǎn)移概率矩陣；第2種方法是通過頻率代替概率來得到間隔為k的轉(zhuǎn)移概率矩陣。在現(xiàn)有文獻中對轉(zhuǎn)移概率矩陣估計的比較大多都是通過預(yù)測精度來實現(xiàn)，而忽略了對轉(zhuǎn)移概率矩陣自身的比較。針對這一問題，提出用偏差和均方誤差來比較轉(zhuǎn)移概率矩陣的估計，并利用計算機進行仿真實驗。

1 馬爾可夫鏈任意等間隔子序列也是馬爾可夫鏈

設(shè)隨機序列{Xn，n=0，1，2，…}可能取值為E={1，2，…，m}，初始分布為q=(q1，q2，…，qm)，且序列滿足馬爾可夫性，即對任意的正整數(shù)n，i0，i1，…，in-1，i，j∈E有

P(Xn+1=j|Xn=i，Xn-1=in-1，…，X1=i1，

X0=i0)=P(Xn+1=j|Xn=i)=pij

故認為序列{Xn，n=0，1，2，…}是齊次馬爾可夫鏈，pij是從狀態(tài)i到狀態(tài)j的一步轉(zhuǎn)移概率，則其一步轉(zhuǎn)移概率矩陣為

序列的l步轉(zhuǎn)移概率為

則其l步轉(zhuǎn)移概率矩陣為

由C-K方程可知

P(l)=Pl

對上述序列{Xn，n=0，1，2，…}進行間隔為k的采樣，得到子序列{X0，Xk，X2k，…，XNk，…}，其可能取值仍為E={1，2，…，m}。對任意的正整數(shù)n，間隔為k的子序列的一步轉(zhuǎn)移概率pkij為

pkij=P{X(n+1)k=j|Xnk=i，X(n-1)k=jn-1，…，

Xk=j1，X0=j0}=

P{X(n+1)k=j|Xnk=i}

(1)

則子序列滿足馬爾可夫性，其一步轉(zhuǎn)移概率矩陣Pk為

(2)

由上述分析可知，原序列間隔為k的子序列滿足馬爾可夫性，且其一步轉(zhuǎn)移概率就是原序列的k步轉(zhuǎn)移概率，其一步轉(zhuǎn)移概率矩陣就是原序列的k步轉(zhuǎn)移概率矩陣。綜上，齊次馬爾可夫鏈的任意等間隔子序列也是齊次馬爾可夫鏈。

2 轉(zhuǎn)移概率矩陣的估計

齊次馬爾可夫鏈通常用于預(yù)測，因此確定齊次馬爾可夫鏈的關(guān)鍵是要確定其轉(zhuǎn)移概率矩陣。設(shè)有一個長度為N的樣本序列{xn，n=1，2，…，N}，其為馬氏鏈{Xn，n=0，1，2，…}的一次觀測。介紹2種估計間隔為k的齊次馬爾可夫鏈{Xn，n=0，1，2，…}轉(zhuǎn)移概率矩陣的方法。第1種方法是利用間隔為1的一步轉(zhuǎn)移概率矩陣和C-K方程計算間隔為k的轉(zhuǎn)移概率矩陣；第2種方法是通過頻率代替概率得到間隔為k的轉(zhuǎn)移概率矩陣。

1) 方法1：基于C-K方程計算轉(zhuǎn)移概率矩陣

用f1ij表示樣本序列x1，x2，…，xN中間隔為1的從狀態(tài)i到狀態(tài)j的頻數(shù)，i，j∈E，從而得到間隔為1的一步轉(zhuǎn)移頻數(shù)矩陣

由式(2)和C-K方程可知

2) 方法2：基于轉(zhuǎn)移頻率計算轉(zhuǎn)移概率矩陣

對子序列{x0，xk，x2k，…，x?N/k」k}，用fkij表示子序列中從狀態(tài)i轉(zhuǎn)移到達狀態(tài)j的頻數(shù)，i，j∈E，從而得到間隔為k的轉(zhuǎn)移頻數(shù)矩陣

則間隔為k的轉(zhuǎn)移概率矩陣的估計為

3 轉(zhuǎn)移概率矩陣估計優(yōu)劣的比較

采用偏差(Bias)和均方誤差(MSE)作為評價指標[14]比較以上2種方法的優(yōu)劣。由于2種方法的比較指標計算方法相同，下面以方法1為例。

偏差(Bias)為

均方誤差(MSE)為

估計的偏差和均方誤差都是一個矩陣，利用矩陣的最大值、最小值和極差來比較上述兩個指標。在方法1中，偏差的比較指標為偏差最大值(Α.B.Maxk)、偏差最小值(Α.B.Mink)、偏差極差(Α.B.Rangek)，其中k表示間隔，且

同理，均方誤差的比較指標為均方誤差最大值(Α.M.Maxk)、均方誤差最小值(Α.M.Mink)、均方誤差極差(Α.M.Rangek)。

4 轉(zhuǎn)移概率矩陣估計的仿真

設(shè)齊次馬爾可夫鏈{Xn}的初始分布為

其一步轉(zhuǎn)移概率矩陣為

則齊次馬爾可夫鏈間隔為k的真實轉(zhuǎn)移概率矩陣可以通過其一步轉(zhuǎn)移概率矩陣和C-K方程得到

和

均方誤差為

1) 偏差的比較。偏差的比較指標有最大值(B.Maxk)、最小值(B.Mink)、極差(B.Rangek)，k表示間隔，則有

通過表1可知：方法1的偏差的最大值和極差相對更??；方法2的偏差的最小值相對更小。因此，以偏差為評判標準時，方法1優(yōu)于方法2。

表1 轉(zhuǎn)移概率矩陣偏差

2) 均方誤差的比較。均方誤差的比較指標有最大值(M.Maxk)、最小值(M.Mink)、極差(M.Rangek)，k表示間隔，通過表2可知：方法1的均方誤差的最大值、最小值、極差都相對更小。

表2 轉(zhuǎn)移概率矩陣均方誤差

根據(jù)圖1可知：隨著間隔k的增大，方法1的均方誤差最大值先下降再趨于平穩(wěn)，方法2的均方誤差最大值先上升后下降再小幅度的上下起伏；方法1的均方誤差最大值的折線圖遠遠低于方法2，說明當k≥2時，方法1的均方誤差最大值遠遠小于方法2的均方誤差最大值。根據(jù)圖2可知：隨著間隔k的增大方法1的均方誤差最小值先上升再趨于平穩(wěn)，方法2的均方誤差最小值先上升再小幅度的上下起伏；方法1的均方誤差最小值的折線圖遠遠低于方法2，說明當k≥2時，方法1的均方誤差最小值遠遠小于方法2的均方誤差最小值。

圖1 MSE最大值折線圖 圖2 MSE最小值折線圖

根據(jù)圖3可知：隨著間隔k的增大，方法1的均方誤差極值先下降再趨于平穩(wěn)，方法2的均方誤差極值上升后下降再小幅度的上下起伏；方法1的均方誤差極值的折線圖遠低于方法2，說明當k≥2時，方法1的均方誤差極值遠小于方法2的均方誤差極值。因此，以均方誤差為評判標準時，方法1優(yōu)于方法2。

圖3 MSE極值折線圖

綜上，將偏差和均方誤差這2個評價指標結(jié)合起來，估計間隔為k(k=2，3，4，5，6，7，8)的轉(zhuǎn)移概率矩陣時方法1優(yōu)于方法2。且多次任意改變初始分布和一步轉(zhuǎn)移概率矩陣進行仿真，所得結(jié)果都驗證了上述結(jié)論的正確性。故在估計齊次馬爾可夫鏈的間隔為k的轉(zhuǎn)移概率矩陣時應(yīng)當選擇第1種方法，即利用間隔為1的一步轉(zhuǎn)移概率矩陣和C-K方程計算間隔為k的轉(zhuǎn)移概率矩陣。

5 結(jié)論

對齊次有限馬爾可夫過程的轉(zhuǎn)移概率矩陣進行估計，首先證明了馬爾可夫鏈任意等間隔的子序列仍是馬爾可夫鏈，提出了2種估計齊次有限馬爾可夫過程的間隔為k的轉(zhuǎn)移概率矩陣的方法，通過仿真利用偏差和均方誤差對2種估計方法進行比較。結(jié)果表明第1種估計方法，即利用間隔為1的一步轉(zhuǎn)移概率矩陣和C-K方程計算間隔為k的轉(zhuǎn)移概率矩陣相對較好。