第19卷第9期計算機技術(shù)與發(fā)展Vol 19 No 9CoMPUTER TECHNOLOGY AND DEVELOPMENTSep.2009EM算法研究與應(yīng)用王愛平,張功營,劉方(安徽大學(xué)計算智能與信號處理教育部重點實驗室,安徽合肥230039)摘要:引入了可處理缺失數(shù)據(jù)的EM算法。EM算法是一種迭代算法,每一次迭代都能保證似然函數(shù)值增加,并且收斂到一個局部極大值。對EM算法的基本原理和實施步驟進行了分析。算法的命名,是因為算法的每一迭代包括兩步第步求期望( Expectation Step),稱為E步;第二步求極大值( Maximization Step),稱為M步。EM算法主要用來計算基于不完全數(shù)據(jù)的極大似然估計。在此基礎(chǔ)上把M算法融合到狀態(tài)空間模型的參數(shù)估計問題。給出了基于 Kalman平滑和算法的線性狀態(tài)空間模型參數(shù)估計方法。關(guān)鍵詞EM算法;狀態(tài)空間模型; Kalman中圖分類號:TP301.6獻標(biāo)識碼:A文章編號:1673-629X(2009)09-0108-03Research and application of EM algorithmWANG Ai-ping, ZHANG Gong-ying, LIU Fang( Ministry of Education Key Lab. of Intelligent Computing Signal ProcessingAnhui University, Hefei 230039, China)Abstrac: Following the description of traditional maximum likelihood estimation methods and the discussions on their disadvantages. EMalgorithm is an iterative algorithm, every iteration to ensure that the likelihood function can be increased, and the oonvergenoe to a localaximA. Presents an EM algorithm that can be used to deal with missing data problems, where the details of the EM algorithm and its realization procedure have been analyzed. Algorithm named because each iterative algorithm includes two steps: the first step in seeking expectations( Expectation Step), known as the e step; the second step for maxima( Maximization Step), kIown as step-by-step M.EM algorithm used to calculate the principel based on incomplete data, maximum likelihood estimation. This is then followed by applyingthe proposed em algorithm to the parameter estimation of state pace models. The paper also presents the Kalman smoothing basedrameter estimation methods for linear state space modelsKey words: EM algorithm; state- space model; Kalman1EM算法及其應(yīng)用間模型的參數(shù)估計問題EM算法是一種迭代算法,每一次迭代都能保證(1)似然函數(shù)值增加,并且收斂到一個局部極大值1,2。y, Hrr, Ur(2)算法的命名是因為算法的每一迭代包括兩個步:第一其中初始狀態(tài)x0,隨機誤差項a和v為獨立高斯分步求期望( Expectation Step),稱為E步;第二步求極大布,三者之間是相互獨立:xo~N(x00,P0),at計算基于不完全數(shù)據(jù)的極大似然估計要用來N0,Q),-iN0,R),E【u]=0,E[uxo]值( Maximization Step),稱為M步。EM算法=0,E[ax0]=0。假設(shè)模型噪聲2狀態(tài)的初始條件都服從高斯2EM算法的狀態(tài)空間模型參數(shù)估計分布,即:下面推導(dǎo)基于 Kalman平滑和算法的狀態(tài)空間模p(ro| 8)=(2x)2 I Poro I-to(x0型參數(shù)估計方法6-10。給出一般形式的線性狀態(tài)空)P((3)收稿日期:2008-12-26;修回日期:200-03-24基金項目:國家自然科學(xué)基金項目(60472065)dH中國煤化工∞p}-號(CNMHG(4)作者簡介:王愛平(1956-),女安徽合肥人,教授,主要從事計算機教學(xué)與研究。p(x!1x1,0)=(2)號1Qs-(x第9期王愛平等:EM算法研究與應(yīng)用109Fx-1)Q1(x-Fx1-1)(5) HrT)+ HPurh,])則完全數(shù)據(jù)的似然也服從高斯分布,完全數(shù)據(jù)的(10)似然可以表示為:p(x0r,y16)=p(x010)Ⅱp(x1x1,4算法實施步驟8) p(y, I xp, o)(6)現(xiàn)在把應(yīng)用EM算法進行狀態(tài)空間模型參數(shù)因此,完全數(shù)據(jù)的對數(shù)似然等于:估計的實施步驟歸納如下:LnP(Io:T, y1:T I8)=1)對參數(shù)0={R,Q,F,H,x00,Pool進行初始化In I Pono 1-4(xo-IoIo)Polo(xo-xouo)2)E步:根據(jù):-1,執(zhí)行卡爾曼濾波及平滑,計算吾Q=2x-Em)xB斯值3)M步:計算參數(shù)={R,Q,F,H,x0,P0新TInIr-1 2( -,x, yR-1( -, 3,的值;n+mN+ min(2r)4)重復(fù)參數(shù)直至收斂14。其中,N是樣本量,m是狀態(tài)變量的維度,n是觀測變5實驗仿真量的維度?？紤]如下的普通線性高斯?fàn)顟B(tài)1空間模型取3求解對數(shù)似然的期望的極大化根據(jù)EM算法基本原理,獲得以下式子y= HI t vr而狀態(tài)初值x0、噪聲干擾v2和v均被假設(shè)為高ELInp(Io: T, 1T I 0)]=-5hn I Pc斯分布,即:InIQ1-hnIR|(n+ m)N+mln(2x)xo -N(oo, Poro), w iid N(o, Q),v"iidtr(PoEL(xo-Ioio)(xo-Ioio)'J)N(0,R)其中參數(shù)實際取值為:F=0.9,H=0.2,Q=0.7-∑(QE[(x1-Fx-1)(x-Fx-1))R=0.2,狀態(tài)初始設(shè)置為x0~N(0,10)。1S(RE[(y-Hx)(y-Hx)1)(8)依據(jù) Kalman濾波和平滑公式,然后聯(lián)合EM算法對模型參數(shù)日=|F,H,Rn,Qu}進行估計EM算法迭代次數(shù)設(shè)置為400(當(dāng)然選代次數(shù)設(shè)置越大,能得到更A=2(ur a P)好的估計結(jié)果)。程序運行結(jié)果θ=10.9835,0.2049,0.2107,0B=∑(x1xt1r+P-1r)683},仿真結(jié)果如圖1~3所示C=∑(x1rgx11r+P:-1r)這樣,完全數(shù)據(jù)的對數(shù)似然期望可以表示為:ELInt(T, yi:TIni QI-lI R I(2x)-ftr(POOEIGCoT -Io)(Ic中國煤化工Por)-2∑(QCNMHGFBFJ)tr([( -Hrr)(y圖1似然函教的期望的Q(日,-1)收斂性110計算機技術(shù)與發(fā)展第19卷用D]廣州暨南大學(xué),2007[2]張杰,陽憲惠.多變量統(tǒng)計過程控制[M].北京:化工工業(yè)出版社,20[3]陳志國,劉文舉基于EM算法的極大似然參數(shù)估計探討J河南大學(xué)學(xué)報:自然科學(xué)版,2002,32(4):35-41[4]孫大飛, Dempster A P, Laird M,etl.MaxiSeries b,1997,39(1):1-38.[5] Meng X L, Rubin D B. Recent Extension to theEM algorithm[M]. Bayesian Statistics 4.Oord圖2參數(shù)F真實值(粗線)和估計值[6]胡壽松,王執(zhí)銓,胡維禮.最優(yōu)控制理論與系統(tǒng)[M].北京科學(xué)出版社,2005[7]郭雷,程代展,馮德興控制理論導(dǎo)論[M]北京:科學(xué)出版社,2005:1-1078 Andrieu C, Doucet A Online Expection-Maximization Type algorithms for Parameter Estimation inGeneral State Space Models[C]//in Proc. IEEE[s.L.]:[s.n.],2003:69-72.[9]賈沛璋朱征桃最優(yōu)估計及其應(yīng)用[M]北京科學(xué)出版社,1994[10] Parzen E On the estimation of a probability densityfunction andmode [J]. annals of Mathematical圖3參數(shù)H的收效過程Statistics,1962,33:1065-1076從上面的仿真圖形可以看出,EM算法用于狀態(tài)[1]waAP, Wang H. Minimising entropy and mean tracking空間模型參數(shù)估計,其估計結(jié)果還是不錯的,可以為時control for affine nonlinear and non Gaussian dynamic變不同的參數(shù)模型估計提供很好的幫助。stochastic system[J]. IEE Proceedings Control Theory & Ap-plication,2005,151(4):405-520[12] Wang A P, Wang H, Tan J. Optimal Filtering for Multivari6結(jié)束語ble Stochastic System via residual Probability Density Func-針對傳統(tǒng)的極大似然參數(shù)估計方法在解決實際問ion Shaping[c]//Proceedings of SICE 2005 Annual Confer-題時所存在的不足,引入了EM算法。討論了算法原ence.[sL]:[sn.],2005:215-219理和實施步驟,并把EM算法推廣到狀態(tài)空間模型的13]任光張均東時間序列狀態(tài)空間模型建模及應(yīng)用(M參數(shù)估計問題,然后推導(dǎo)了基于Kamn平滑和EM算大連:大連海事大學(xué)出版社,2000102-112法的狀態(tài)空間模型參數(shù)估計的過程,最后給出實驗仿14cuL. Wang H. Mininum entropy filtering for multivariatestochastic systems with non-Gaussian noises [J].IEEETransactions on Automatic Control, 2006, 51(4): 670-695[15]梁應(yīng)敞張賢達李衍達,等非高斯相關(guān)噪聲中高斯信號參考文獻:的時延估計[門]電子科學(xué)學(xué)刊,1997,19(5):607-612[1]陳學(xué)華狀態(tài)空間模型理論與算法及其在金融計量中的應(yīng)多史國計算機學(xué)會余我史《計算機技術(shù)與發(fā)展歡迎訂閱VLl中國煤化工CNMHG