沈鋒， 徐廣輝， 桑靖

(哈爾濱工程大學自動化學院，黑龍江哈爾濱150001)

一種自適應變分貝葉斯容積卡爾曼濾波方法

沈鋒， 徐廣輝， 桑靖

(哈爾濱工程大學自動化學院，黑龍江哈爾濱150001)

針對應用于非線性系統(tǒng)模型的容積卡爾曼濾波工作性能會受觀測噪聲參數(shù)變化的影響而降低的問題，提出一種自適應的變分貝葉斯容積卡爾曼濾波算法。在每一次更新步驟中，將系統(tǒng)狀態(tài)與變化的觀測噪聲統(tǒng)計信息一起作為隨機變量，并用變分貝葉斯方法進行估計，在迭代逼近得到噪聲方差后，再利用容積卡爾曼濾波對系統(tǒng)狀態(tài)進行更新。仿真實驗證明變分貝葉斯容積卡爾曼濾波算法在非線性系統(tǒng)的濾波問題中能夠較好跟蹤變化的觀測噪聲方差，相比容積卡爾曼濾波擁有較好的估計性能。

變分貝葉斯;容積卡爾曼濾波;自適應;非線性系統(tǒng)

0 引 言

作為卡爾曼濾波的衍生，擴展卡爾曼濾波(EKF)、無跡卡爾曼濾波(UKF)、容積卡爾曼(c KF)等成熟的非線性濾波算法自提出以來已經(jīng)受到了廣泛而深入的研究［1-2］。但此類非線性濾波算法在實際非線性系統(tǒng)的應用中會有一定的局限性，原因在于這些算法本身的建立依賴于準確的模型、確定的系統(tǒng)參數(shù)、已知的噪聲統(tǒng)計特性，而在實際工程領域，由于人們對工況的認知有限，建?？赡軙霈F(xiàn)誤差，即使模型參數(shù)、噪聲統(tǒng)計特性得到確定，在系統(tǒng)實際運行中，系統(tǒng)本身存在攝動，同時容易受到外界的干擾，系統(tǒng)噪聲或觀測噪聲也可能隨之發(fā)生變化，這種情形下，此類非線性濾波器的工作性能會發(fā)生退化，極端情形下甚至無法正常工作。

針對這種情形，自適應的濾波算法開始受到了關注。就目前的研究方法來講，自適應濾波算法主要有貝葉斯法、最大似然法、相關法以及協(xié)方差匹配法［3］。其中，相關法由于它在計算上的相對便捷，受到了更多關注與研究。而事實上，后3種方法都可以視為自適應濾波算法在貝葉斯框架下的特例，它們在廣義上都屬于貝葉斯法。然而，在貝葉斯法中，由于涉及概率密度函數(shù)的計算，積分項都過于復雜，一般情況下難以得到確定的解析解。這使得貝葉斯法通常只能在理論上得到解釋，而在實際應用上則顯得非常有限。

通常情況下，在貝葉斯準則的基礎上，采樣方法被引入用以隨機逼近得到近似的精確的參數(shù)或模型。一個比較典型的方法便是蒙特卡洛馬爾可夫鏈McMc(Monte carlo Markov chain)方法。它通過構造馬爾可夫鏈(Markov chain)的極限不變分布來模擬高維數(shù)的積分運算［4］。然而McMc等隨機采樣方法的估計精確度是以犧牲計算量為代價的，在許多實時性要求較高的工程應用領域有很大的局限性，而且如何判定馬爾可夫鏈的收斂也是一個難題。

除了隨機逼近，還有一類被稱作確定性逼近的方法，變分貝葉斯方法便是其中的一種。Attias H在文獻［5］中的詳細論述被認為是變分貝葉斯方法最早的理論體系。該方法提議引入一個新的形式簡單的分布，通過迭代更新變分參數(shù)，不斷最大化待估計參數(shù)的邊緣似然函數(shù)的下界來逼近參數(shù)的真實后驗分布，直至算法收斂。

雖然變分貝葉斯方法在提出之初的幾年里普遍是用于系統(tǒng)模型的參數(shù)估計，但是近幾年，因為變分貝葉斯方法相比McMc等采樣方法在估計上的快速性，使得它的應用領域已經(jīng)從圖像處理［6］、盲源分離［7］、語音增強［8］、信道估計［9］等參數(shù)推斷領域延伸到了狀態(tài)估計領域，如今，在濾波問題方面，變分貝葉斯方法也受到了許多學者的關注和研究: Vrettas M D將變分貝葉斯方法用于隨機動態(tài)系統(tǒng)的參數(shù)估計問題當中［10］;Boujemma A將變分貝葉斯方法用于動態(tài)斷面X射線照相術的卡爾曼濾波之中［11］;Sarkka S提出了將變分貝葉斯方法用于估計經(jīng)典卡爾曼濾波算法中的觀測噪聲［12］;Sun J L則通過變分貝葉斯，設計了針對未知系統(tǒng)輸入下的兩步卡爾曼濾波器［13］;孫世軍則用基于時間序列的變分貝葉斯方法用于盲源分離問題之中去估計源信號與混合矩陣［14］。

本文在文獻［12-13］的研究成果的基礎上，針對工程領域中普遍存在的非線性系統(tǒng)模型，用概率的方法，從廣義貝葉斯濾波角度著重設計了基于變分貝葉斯的噪聲自適應容積卡爾曼濾波器。該濾波器可以自適應地運用變分貝葉斯方法，動態(tài)地跟蹤觀測噪聲，從而有效地提高估計精確度。文中首先給出變分貝葉斯學習的基本原理，并簡單的闡述了容積卡爾曼濾波方法;接著給出了本文所提出的變分貝葉斯容積卡爾曼濾波方法原理及實現(xiàn)步驟，并對與傳統(tǒng)的容積卡爾曼濾波方法進行了計算機仿真對比;最后給出了結論。

1 變分貝葉斯

在參數(shù)估計的問題中，在獲得了觀測樣本的數(shù)據(jù)集Z后，根據(jù)貝葉斯準則，核心是參數(shù)集θ的后驗概率密度函數(shù)的計算

而式(1)計算的一個難點在于分母，邊緣似然概率密度函數(shù)p(Z)的計算。正如引言所述，通常情況下p(Z)的計算難以得到精確的解析解，引入變分貝葉斯方法，引入一個簡單的近似分布函數(shù)q(θ)，并取p(Z)的對數(shù)形式:

另一方面，又有

通過求偏導，可以得到F(q(θ))的極值，其通解表達式為［4］

可以發(fā)現(xiàn)，每個參數(shù)θi的近似分布可以通過求對數(shù)聯(lián)合概率密度函數(shù)關于其他參數(shù)分布q(θj≠i)的期望求得，所以每一個參數(shù)分布的計算都依賴于其他參數(shù)的分布。這就形成了迭代的機制:在給定先驗知識的情況下，初始化參數(shù)值，通過變分貝葉斯算法循環(huán)迭代計算，進行參數(shù)更新，直至自由能量(對數(shù)邊緣似然函數(shù)的下界)達到最大值，判定算法收斂并結束，此時可以得到系統(tǒng)模型參數(shù)的估計值。

2 容積卡爾曼濾波

針對所有離散的線性與非線性狀態(tài)空間模型，它們的概率表示形式如下:與代表了一般貝葉斯意義下的系統(tǒng)狀態(tài)方程與系統(tǒng)觀測方程。對于非線性系統(tǒng)，考慮加性的系統(tǒng)噪聲與觀測噪聲，分別有其中與h(xk)為非線性函數(shù)，wk～N(0，qk)，vk～N(0，rk)均服從高斯分布。所以，對于狀態(tài)空間模型有

而狀態(tài)變量xk服從均值為mk，協(xié)方差為Pk的高斯分布，即于是，在k-1時刻觀測后，先驗概率而在k時刻觀測后，后驗概率所以，在遞歸的容積卡爾曼濾波算法中［17］，有

預測方程:

更新方程:

3 變分貝葉斯容積卡爾曼濾波

在容積卡爾曼濾波中，觀測噪聲的統(tǒng)計矩信息被認為是確定的，且保持不變?，F(xiàn)在在變分貝葉斯容積卡爾曼濾波算法中，把動態(tài)觀測噪聲的方差rk和狀態(tài)變量xk當作隨機變量，作為待估計的參數(shù)。在k-1觀測時刻后，兩者的聯(lián)合概率密度函數(shù)的先驗分布為

在k時刻觀測后，聯(lián)合概率密度函數(shù)的后驗分布為

這樣，式(10)與式(11)就構成了廣義貝葉斯濾波理論的預測方程和更新方程。但是除了一些比較簡單的情形，式(10)與式(11)當中的積分運算很難得到解析解，所以引入變分貝葉斯方法求其次優(yōu)近似解。

隨機變量xk與rk被認為相互獨立，根據(jù)先驗知識，認為它們分別服從高斯分布與逆Gamma分布。因此在k-1時刻后，

經(jīng)歷了第k時刻的觀測后，在變分貝葉斯方法中，引入一個新的分布來代替真實的后驗分布為了推導的簡潔，文章后面在寫法上均略去了新的分布對于觀測量z1:k的依賴，同時考慮到兩個變量相互獨立，所以同樣對于聯(lián)合概率密度函數(shù)的近似分布有q(xk，rk)=q(xk)q(rk)。

式(5)給出了變分貝葉斯方法近似解的通解表達式，采用容積卡爾曼濾波算法中一樣的處理方法，對hk(.)的線性化處理，對于狀態(tài)變量:

而對于觀測噪聲的方差:

其中，c是一個與分布形式無關的常數(shù)?？梢园l(fā)現(xiàn)，式(13)與式(14)在形式上分別是高斯分布與逆Gamma分布，只是參數(shù)與先驗分布的不同。這是因為選取的高斯分布與逆Gamma分布都屬于共軛指數(shù)分布域，而正是共軛性保證了變量在先驗分布與后驗分布在形式上的一致性［15］。

近似的高斯分布的參數(shù)擁有如下表達式:

其中

近似的逆Gamma分布表示為

注意到在(18)中，方差陣對角線上元素的期望的計算需要已知的逆Gamma分布的信息，根據(jù)逆Gamma分布的性質，即

而在(19)中，根據(jù)容積卡爾曼的采樣策略，期望部分可以繼續(xù)展開如下:

于是，式(13)～式(19)就構成了全部對觀測噪聲rk與系統(tǒng)狀態(tài)變量xk進行聯(lián)合估計的變分貝葉斯算法的步驟:首先給定先驗知識并初始化各分布參數(shù)，然后開始對近似分布和各個期望進行迭代計算，得到更新的分布式(13)與式(14)后，分別計算參數(shù)在分布下的新的期望式(18)與式(19)，進而再利用式(15)與式(17)中的更新方程對分布參數(shù)進行更新。如此循環(huán)迭代計算，直至算法收斂，得到逼近的估計值。

在用變分貝葉斯方法得到估計的觀測噪聲方差后，可以將其與容積卡爾曼濾波算法相融合，對非線性系統(tǒng)模型進行自適應的容積卡爾曼濾波處理，濾波算法總結如下。

初始化:

預測步驟:

式(8);

更新步驟:

式(15)～式(19);

需要指出的是在預測步驟中，對于超參數(shù)αk-1與βk-1的預測，借鑒了文獻［12］中觀測噪聲方差的建模方法，這樣使得觀測噪聲的方差的分布參數(shù)能夠平穩(wěn)地發(fā)生變化。其中ρ是一個在(0，1］內的變化因子。而對比容積卡爾曼濾波器，該濾波器在預測步驟中對噪聲方差進行動態(tài)建模，在更新步驟中先用變分貝葉斯方法，迭代地估計觀測噪聲的方差，并在得到噪聲的統(tǒng)計信息后對系統(tǒng)的狀態(tài)進行更新，如此反復直至濾波算法結束。

4 仿真分析

為了驗證變分貝葉斯容積卡爾曼濾波器(VB-c KF)的性能，用設計的濾波器去考察一個存在變化噪聲方差的非線性模型:

其中rk為待估計的、不確定的、服從高斯分布的白噪聲的方差。本實驗中，在［1:200］、［201:400］、［401:600］時刻，依次取rk為1、5、3，并分別用容積卡爾曼濾波器(c KF)與設計的濾波器進行性能比較。

圖1列出了在實驗中，容積卡爾曼濾波算法和變分貝葉斯容積卡爾曼濾波算法對狀態(tài)跟蹤的對比，圖中用狀態(tài)的均方根誤差作為縱軸，仿真次數(shù)作為橫軸。從中可以看出，由于人為地取觀測噪聲的方差值使之發(fā)生波動，而容積卡爾曼濾波缺乏對噪聲的實時估計與信息更新，在精確度上要遜色于變分貝葉斯容積卡爾曼濾波。

圖1 濾波器性能對比Fig.1 Performance com parison between two different filters

圖2給出了變分貝葉斯容積卡爾曼濾波器對噪聲方差的跟蹤情況。從中可以看出，在兩次噪聲發(fā)生劇烈的跳變后，利用變分貝葉斯方法迭代地逼近真實的噪聲方差，不僅有著很快的收斂速度，能夠較快地給出方差的估計值，而且在估計精確度上也有著較好的表現(xiàn)。對應體現(xiàn)在圖1中，相比容積卡爾曼濾波，變分貝葉斯容積卡爾曼濾波器有著更小的誤差。

圖2 噪聲方差跟蹤圖Fig.2 Tracking of noise variance

5 結 論

針對典型的非線性系統(tǒng)模型，在容積卡爾曼濾波器的基礎上結合變分貝葉斯方法，本文提出了一種自適應變分貝葉斯容積卡爾曼濾波方法，該方法對狀態(tài)量進行更新之前，實時地跟蹤觀測噪聲的方差，使濾波器的參數(shù)設計能夠盡可能和工作狀況相匹配，從而達到自適應濾波的目的，仿真結果驗證了本文所提方法的有效性。

［1］ S.J.Julier，J.K.Uhlmann，H.F.Durrant-Whyte，A new approach for filtering nonlinear systems［c］//In Proc.1995 American control，conference，Seattle，Washington，1995:1628 -1632.

［2］ D.Simon，Optimal state estimation［M］.Hoboken，New Jersey: John Wiley&Sons，Inc.，2006:397-407，447-452.

［3］ R.Mehra，Approaches to adaptive filtering［J］.IEEE Transactions on Automatic control，1972，17(5):693-698.

［4］ c.P.Robert and G.casella，Introducing Monte carlo Methods with R［c］.Springer，2010:167-175.

［5］ D.J.c.MacKay，Developments in probabilistic modelling with neural networks-ensemble learning:in Neural Networks:Artificial Intelligence and Industrial Applications［c］//Proc.3rd Annu. Symposium on Neural Networks，Nijmegen，Netherlands，14-15 Sept.1995.Berlin:Springer，1995:191-198.

［6］ H.Attias，Inferring parameters and structure of latent variable models by variational Bayes［c］//In Proc.15th conf.Annu. conf.on Uncertainty in Artificial Intelligence，San Francisco，c A，1999:21-30.

［7］ G.chantas，N.Galatsanos，A.Likas.Variational Bayesian image restoration based on a product of t-distributions image prior［J］. IEEE Transactions on Image Processing，2008，17(10):1795 -1805.

［8］ A.Honkela，H.Valpola，A.Llin.Blind separation of nonlinear mixtures by variational Bayesian learning［J］.Digital Signal Processing，2007，17(5):914-934.

［9］ Q.Huang，J.Yang，Y.Zhou，Variational Bayesian method forspeech enhancement［J］.Neurocomputing，2007，70(16-18): 3063-3067.

［10］ K.Harada，H.Sakai，Variational Bayesian blind estimation of SIMO channels［c］//In Proc.IEEE Int.conf.Acoustics，Speech，and Signal Processing，Dallas，TX，Mar.14-19，2010:3218-3221.

［11］ M.D.Vrettas，D.cornford，M.Opper，Estimating parameters in stochastic systems:A variational Bayesian approach［J］.Physica D，2011，240(23):1877-1900.

［12］ B.Ait-El-Fquih T.Rodet，Variational Bayesian Kalman filtering in dynamical tomography［c］//In Proc.IEEE International conference on Acoustics，Speech，and Signal Processing，Prague，c zech Republic，May 22-27，2011:4004-4007.

［13］ S.Sarkka，A.Nummenmaa，Recursive noise adaptive Kalman filtering by variational Bayesian approximations［J］.IEEE Transactions on Automatic control，2009，54(3):596-600.

［14］ J.Sun，J.Zhou X.Gu，Variational Bayesian two-stage Kalman filter for systemswith unknown inputs［J］.Procedia Engineering，2012，29(1):2265-2273.

［15］ S.Sun，c.Peng，W.Hou，et al，Blind source separation with time series variational Bayes expectation maximization algorithm［J］.Digital Signal Processing，2012，22(1):17-33.

［16］ M.Opper and O.Winther，F(xiàn)rom na?vemean field theory to the tap equations［M］.In Advanced Mean Field Methods，M.Opper and D.Saad，Eds.cambridge，MA:MIT Press，2001:7-20.

［17］ M.J.Beal，Variational algorithms for approximate Bayesian inference［D］.Ph.D.dissertation，Univ.college London，London，2003.

(編輯:劉素菊)

Adaptive variational Bayesian cubature Kalman filtering

SHEN Feng， XU Guang-hui， SANG Jing
(college of Automation，Harbin Engineering University，Harbin 150001，china)

Focusing on the performance of cubature Kalman filteringmay be degraded due to the fact that in practical situations the statistics ofmeasurement noisemight change.An adaptive variational Bayesian cubature Kalman filtering algorithm was proposed which can be used in non-linear system models.In each update step of proposed method，both system state and time-variantmeasurement noise were recognized as random variables to estimate.Measurements noise variances were approximated by variational Bayes，thereafter，system states were updated by cubature Kalman filtering.Simulation results demonstrate the proposed filter can well track measurement noise for a non-linear system and outperforms cubature Kalman filter.

variational Bayes;cubature Kalman filtering;adaptive;non-linear system

10.15938/j.emc.2015.04.015

TP 202

A

1007-449X(2015)04-0094-06

2013-05-29

國家自然科學基金(61102107，61374208);中央高?；究蒲袠I(yè)務費專項資金(HEUc FX41310)

沈 鋒(1981—)，男，副教授，研究方向為非線性濾波技術等;徐廣輝(1987—)，男，博士研究生，研究方向為自適應信號處理技術等;桑 靖(1988—)，女，碩士研究生，研究方向為多傳感器信息融合技術等。

沈 鋒