基于近似貝葉斯計算的HMM隱狀態估計

2014-03-20 12:03胡錫健

山東理工大學學報（自然科學版） 2014年1期

陳琦, 胡錫健

(新疆大學數學與系統科學學院, 新疆烏魯木齊 830046)

隱馬爾可夫模型(Hidden Markov Model，HMM)是由Baum L E在20世紀60年代末提出的一個非常重要的統計分析模型，其在語音識別、故障診斷、生物統計等領域中有廣泛的應用．評估問題、隱狀態估計(解碼)問題、參數估計(學習)問題是HMM的三個基本問題.其中，隱狀態估計問題在某種程度上其實就是一個貝葉斯估計問題.當HMM的隱狀態取連續值時，通常將其稱為狀態空間模型.對于線性高斯系統，可以利用卡爾曼濾波(Kalman Filter)得到解析解.對于狀態空間是離散的、有限的馬爾可夫鏈的系統，運用基于網格的濾波方法(Grid-based methods)可以得到后驗概率分布的解析解.但實際上這些條件在很多動態隨機系統中很難滿足，大部分是非線性、非高斯的，無法求出解析解，因此，只能求助于次優或逼近方法.文獻[1]提出擴展卡爾曼濾波(Extended Kalman Filter,EKF),但EKF存在線性化誤差，只適用于弱非線性系統.文獻[2]提出UT的變換，并在EKF濾波框架里得到廣泛應用，被稱為UKF，這種方法不必對非線性方程進行線性化和計算雅克比矩陣，均值和方差均可精確到二階，比標準的EKF精度更高.但由于這兩種非線性估計方法都是把pk(xk|y1:k)用高斯逼近，如果真實的分布是非高斯的，那么就不能精確表示狀態的真實分布.然而，對于狀態空間是連續的但是可以劃分成有限的Ns個單元，那么仍然可以用基于網格的逼近方法來近似后驗概率密度[3].近年來，將粒子濾波用于處理非線性、非高斯系統的估計，為貝葉斯估計提供了有效的新解法.但當隱狀態空間維數很高時離子濾波方法有以下兩個不足之處：計算量特別大；不能精確近似目標分布.此外，離子濾波還要求似然函數是逐點可知的，但在實際情況中通常不能滿足，而用近似貝葉斯計算方法(Approximate Bayesian Compute,ABC)可以避免這些問題.本文將HMM隱態估計問題看作一個貝葉斯濾波問題，首先詳細介紹粒子濾波用于HMM隱狀態估計，然后將近似貝葉斯計算思想引入HMM隱狀態估計中.

1 問題描述

HMM是一對離散時間隨機過程{Xn}n≥1、{Yn}n≥1，其中xn∈dx是隱狀態，yn∈dy是觀測值。HMM隱狀態估計，就是給定觀測序列y0:n及模型參數θ，對隱狀態作出恰當的判斷，即在某種意義下求與此觀測序列對應的最優隱狀態序列x0:n.這一問題的回答取決于對模型狀態物理意義的理解，根據不同的準則，可能會選出不同的最優狀態序列．

單點最優原則下的隱狀態估計問題其實是一個最優濾波問題[4]，HMM模型的隱狀態是由φk|n(xk|y0:n)而得到的，即時刻k(0≤k≤n)下最有可能的隱狀態為

(1)

(2)

2 基于離子濾波的HMM隱狀態估計

粒子濾波算法是通過非參數化的蒙特卡洛(MonteCarlo)模擬方法來實現遞推的貝葉斯濾波，適用于任何能用狀態空間模型描述的非線性系統，其核心思想是得到一批有權重的離散隨機采樣點(粒子)來近似狀態變量的后驗概率密度函數，然后利用他們得到估計值，當得到的粒子很多時，其精度可以逼近最優估計.

2.1 離子濾波估計過程描述

HMM隱狀態估計的核心是后驗分布φ0:n|n(y0:n，·)，假設其密度函數存在，則

(3)

(4)

L(y0:n)φ0:n|n(x0:n|y0:n)=L(y0:n-1)φ0:n-1|n-1(x0:n-1|y0:n-1)q(xn-1,xn)g(xn,yn)

(5)

(6)

得φk|n(xk|y0:n)的估計為

(7)

2.2 建議分布選取

故(4)式的最優建議分布為

(8)

2.3 重抽樣

在進行序貫重要性采樣時，隨著 k的增大，將會出現權值退化問題，當n較大時，就將大量的時間浪費在一些對估計結果不起任何作用的粒子更新上，針對權值退化問題，Gordon等提出采用重抽樣(Resampling)方法，即評估粒子權值后重抽樣粒子以減少小權值的粒子而復制大權值的粒子.

2.4 粒子濾波算法

(a)當k=0時

(b)當k=1,2,…，n時

第k步抽樣：

第k步重抽樣：

由于

(9)

由(7)式知

(10)

(11)

設

(12)

(13)

3 基于ABC(Approximate Bayesian Compute)的HMM隱狀態估計

3.1 ABC的引入

粒子濾波算法是基于重點抽樣來近似目標分布的序列，通常目標分布的維數隨著時間參數的增加而增加.在維數很高時，離子濾波有計算量大、不能準確地近似目標分布等缺點.第一個缺點歸因于進行模擬和估計權值時的計算量都很大，第二個缺點是由于目標分布本身的復雜性.此外，離子濾波通常要求g(xk,yk)是逐點可得的，但在一些應用中不能滿足該要求.采用基于近似貝葉斯計算的方法可以解決此類問題.

近似貝葉斯計算已經成為貝葉斯范例的一個重要方面，在人們感興趣的許多實際問題中，貝葉斯模型不可得，或者由于計算量太大而不能被擬合，或者所需的模擬方法不能處理問題的復雜性.文獻[5]介紹了后驗密度π(x|y)∝g(y|x)μ(x)的如下近似：

(14)

其中：Aε,y={u:ρ(s(u),s(y))<ε}；u是一個輔助數據點，∏A是集合A的示性函數，g是似然函數，μ是先驗密度，s:dy→ds是一個數據的匯總統計量，ρ是一個距離測度.通常ds比dy小很多.如果s是充分的，可以看到當ε趨于0時，πε(x|y)趨于π(x|y).可以采用拒絕抽樣[5]、馬爾可夫蒙特卡洛[6]或序貫蒙特卡洛[7-8])等方法對(14)式進行統計推斷，其中的關鍵點是，若不能對g(u|x)進行數值估計，但對任一x，如果可以從g(·|x)進行抽樣，仍就可以從πε(x|y)生成樣本.

3.2 ABC目標的構建

πn(x1:n|y1:n)的邊際ABC目標為

(15)

3.3 ABC目標的SMC算法

第1步：

4 結束語

文中將HMM的隱狀態估計問題看成一個貝葉斯估計問題，在單點原則下就是一個最優濾波問題，而解決濾波問題的方法很多，但適用的條件及范圍各不相同，其得到的估計結果精度也不相同.離子濾波可以很好地對大部分HMM的隱狀態的后驗密度進行估計，但其也有自身的一些缺點，并且當遇到g(xk,yk)逐點不可得時，無法對后驗密度進行估計.而基于ABC的方法可以對此類情況進行估計，從而使HMM隱狀態估計問題得到了全面解決，而且也可以考慮將ABC方法用到其它更廣泛的領域中.

[1] Anderson B D,Moore J B.Optimal filtering[M].New Jersey:Prentice-Hall,1979.

[2] Julier S J, Uhlmann J K.A new extension of the Kalman filter to nonlinear systems[C]// The 11th Int Symp on Aerospace:Simulation and Controls. London: 1997:10-12.

[3] 張詩桂，朱立新，趙義正.離子濾波算法研究進展與展望[J].自動化技術與應用 ,2010,29(6):1-16.

[4] 朱成文，李兵，胡奎，等.HMM隱狀態的粒子濾波估計[J].計算機工程與應用,2012,48(8):161-167.

[5] Pritchard J K, Seielstad M T, Perez-Lezaun A,etal. Population growth of human Y chromosome microsatellites[J]. Mol Biol Evol.,1999,16, 1 791-1 798.

[6] Majoram P, Molitor J , Plagnol V,etal. Markov chain Monte Carlo without likelihoods[J]. Proc Natl Acad Sci,2003,100, 15 324-15 328.

[7] Beaumont M A, Cornuet J M, Marin J M, et al.Adaptivity for ABC algorithms: the ABC- PMC scheme[J]. Biometrika, 2009,96, 983-990.

[8] Doucet A, Johansen A M. A tutorial on particle filtering and smoothing: fifteen years later[J]. Oxford Handbook of Nonlinear Filtering,2009,7：54-89.

[9] Jasra A ,Singh S S,Martin J S,etal. Filtering via approximate Bayesian computation[J]. Statistics and Computing, 2012,22(6):1 223-1 237.

[10] Moral P D, Doucet A, Jasra A. An adaptive sequential Monte Carlo method for approximate Bayesian computation[J]. Statistics and Computing, 2012,22(5):1 009-1 020.