非T噪聲條件下自適應(yīng)目標(biāo)跟蹤方法研究

張 揚(yáng),馬天力,高 嵩,陳超波

(西安工業(yè)大學(xué)電子信息工程學(xué)院,陜西 西安 710021)

1 引言

近年來,隨著航空、航天、航海事業(yè)的不斷發(fā)展,目標(biāo)跟蹤技術(shù)在對(duì)海底、海面、陸地和空中等領(lǐng)域取得了極大的應(yīng)用。在目標(biāo)跟蹤過程中,傳感器測量信息不僅包含目標(biāo)信號(hào),同時(shí)也包含隨機(jī)量測噪聲。目標(biāo)跟蹤要解決的基本問題就是從被噪聲污染的量測信息中盡可能充分地濾除干擾噪聲,使系統(tǒng)的狀態(tài)估計(jì)更加準(zhǔn)確。傳統(tǒng)的濾波方法多采用卡爾曼濾波器[1],其基于最小均方誤差準(zhǔn)則,對(duì)于準(zhǔn)確建模且噪聲服從高斯分布的動(dòng)態(tài)系統(tǒng),KF能夠遞歸得到狀態(tài)的一致最小方差線性無偏估計(jì),但要求量測噪聲滿足高斯分布的假設(shè)。通常,量測噪聲的統(tǒng)計(jì)特性可以事先由傳感器的物理特性得到。但對(duì)于跟蹤系統(tǒng)而言,由于外界干擾、加速度的物理特性和環(huán)境等因素影響,很難用一個(gè)準(zhǔn)確的統(tǒng)計(jì)特性來表示量測噪聲,所以其統(tǒng)計(jì)特性往往是未知且時(shí)變的。針對(duì)量測噪聲時(shí)變的目標(biāo)跟蹤問題,一般采用自適應(yīng)卡爾曼濾波器(AKF),例如Sage-Husa自適應(yīng)濾波器(SHAKF)[2-4],基于新息的自適應(yīng)濾波器(IAKF)[5-8]和交互式多模型自適應(yīng)濾波器(IMMAKF)[9,10]。SHAKF可以用最大后驗(yàn)準(zhǔn)則估計(jì)量測噪聲的統(tǒng)計(jì)量,但是它不能保證量測噪聲協(xié)方差的收斂性,易造成濾波發(fā)散[11]。IAKF可以通過新息協(xié)方差對(duì)量測噪聲進(jìn)行估計(jì)和調(diào)整,使濾波算法能更好地適應(yīng)變化的噪聲統(tǒng)計(jì)特性[12]。但是,其只適用于緩慢變化的量測噪聲協(xié)方差,因?yàn)閷?shí)現(xiàn)量測噪聲協(xié)方差的可靠估計(jì)需要很大的數(shù)據(jù)窗口[13]。IMMAKF可以使用多個(gè)模型來描述目標(biāo)的運(yùn)動(dòng)狀態(tài),通過有效的加權(quán)融合進(jìn)行目標(biāo)運(yùn)動(dòng)狀態(tài)估計(jì),很好地克服了單模型濾波算法估計(jì)誤差較大的問題,但其計(jì)算比較復(fù)雜[14]。

在實(shí)際的工程應(yīng)用中,跟蹤目標(biāo)容易受到來自外界的干擾,量測噪聲存在野值,導(dǎo)致噪聲分布出現(xiàn)重尾特征。例如在雷達(dá)跟蹤系統(tǒng)中,由于目標(biāo)位置的隨機(jī)擺動(dòng),使得對(duì)目標(biāo)位置的量測伴隨著閃爍噪聲的出現(xiàn),這種閃爍噪聲導(dǎo)致系統(tǒng)的量測呈現(xiàn)非高斯特性[15]。直接使用上述方法會(huì)導(dǎo)致算法性能降低,濾波精度下降。為解決重尾量測噪聲的濾波問題,Guo等[16]提出了一種混合高斯估計(jì)算法(GM-CDMKF)解決非高斯噪聲的目標(biāo)跟蹤問題,在不增加過多計(jì)算的情況下,改善了非高斯噪聲的濾波性能。文獻(xiàn)[17]研究了用混合拉普拉斯模型處理重尾噪聲問題,該方法使用梯度下降法求解負(fù)似然對(duì)數(shù)最小值來估計(jì)參數(shù)。Roth等[18]將量測噪聲的統(tǒng)計(jì)特性近似為T分布,推導(dǎo)出了一種魯棒T分布濾波算法,該方法借助T分布對(duì)量測異常值建模,最終通過T分布來近似系統(tǒng)的后驗(yàn)概率密度函數(shù)。Zhu等[19]提出基于T分布量測噪聲的VBT濾波方法,利用變分貝葉斯學(xué)習(xí)算法進(jìn)行模型的參數(shù)估計(jì),可以解決重尾噪聲模型的濾波問題,提高了估計(jì)精度。但是基于T分布的濾波算法性能可能會(huì)因自由度參數(shù)的選擇不當(dāng)而降低。文獻(xiàn)[20]將粒子濾波算法用于T分布的非高斯噪聲模型,該算法基于蒙特卡洛思想,將粒子進(jìn)行重采樣以獲得足夠的統(tǒng)計(jì)數(shù)據(jù),但其在改善跟蹤性能的同時(shí)也需要大量的粒子,使計(jì)算量大大增加[21]。在目標(biāo)跟蹤的實(shí)際應(yīng)用中,為了保持最優(yōu)濾波系統(tǒng)的有效性,需要比較精確的表征傳感數(shù)據(jù)的統(tǒng)計(jì)特性。而可靠性較差的傳感器誘導(dǎo)出的異常量測值,傳感器發(fā)生故障或受到環(huán)境中的電磁干擾,量測模型往往伴隨著不可避免的重尾偏斜噪聲,這會(huì)嚴(yán)重降低以高斯分布或T分布為前提假設(shè)的目標(biāo)跟蹤系統(tǒng)的濾波性能。

針對(duì)量測噪聲重尾不對(duì)稱條件下的線性目標(biāo)跟蹤系統(tǒng)的狀態(tài)估計(jì)問題,本文提出了一種解決重尾不對(duì)稱量測噪聲的VBST算法。每次時(shí)間更新后,用變分貝葉斯學(xué)習(xí)估計(jì)得到近似后驗(yàn)分布,再根據(jù)量測噪聲對(duì)目標(biāo)的速度和位置狀態(tài)進(jìn)行更新,從而達(dá)到自適應(yīng)濾波的目的。通過與傳統(tǒng)的卡爾曼濾波算法、VBT算法以及粒子濾波算法對(duì)比,驗(yàn)證本文所提算法的濾波性能和魯棒性。

2 問題描述

建立如(1)式的偏斜T分布量測噪聲的系統(tǒng)模型

xk=Axk|k-1+wk

zk=Cxk+ek

(1)

其中,xk∈n是狀態(tài)向量。A∈n×n是狀態(tài)轉(zhuǎn)移矩陣,wk∈n是過程噪聲,假設(shè)其服從均值為零,方差為Q的高斯白噪聲。zk∈m是在k時(shí)刻的量測值,C∈m×n是量測模型矩陣。ek∈m是獨(dú)立于wk的重尾不對(duì)稱的非高斯量測噪聲

(2)

當(dāng)i=1時(shí),式(2)為一維偏斜T分布,其概率密度函數(shù)可以表示為[22]

(3)

其中

(4)

(5)

T(·;0,1,v)是T分布的累積分布函數(shù)。Γ(·)表示伽馬函數(shù)。圖1給出了不同偏斜系數(shù)下偏斜T分布的概率密度函數(shù)圖[23],其中,偏斜T分布的其它參數(shù)設(shè)置為:μ=0,σ=1,v=4。

圖1 不同偏斜系數(shù)下的偏斜T分布概率密度函數(shù)

如圖1所示,若偏斜系數(shù)大于0,偏斜T的概率密度函數(shù)左偏;若偏斜系數(shù)小于0,偏斜T分布的概率密度函數(shù)右偏。

當(dāng)δ趨于0時(shí),偏斜T分布近似為學(xué)生T分布;當(dāng)δ趨于0且v趨于∞時(shí),偏斜T分布近似為正態(tài)分布[20]。

根據(jù)狀態(tài)空間模型,p(xk|z1:k-1)的概率密度函數(shù)為

p(xk|z1:k-1)=N(xk;xk|k-1,Pk|k-1)

(6)

z1,…zN是來自ST(zk;Cxk,R,Δ,v)的N個(gè)獨(dú)立的觀測值,其似然函數(shù)為[24]:

(7)

(8)

(9)

其中,N(·)表示高斯分布,HN(·)表示半正態(tài)分布,G(·)表示伽馬分布。

基于貝葉斯規(guī)則,其聯(lián)合概率分布為

(10)

圖2表示偏斜T分布量測噪聲系統(tǒng)的圖模型,圖中方框表示先驗(yàn)的超參數(shù),深色的圓形表示系統(tǒng)的量測。

圖2 偏斜T分布的圖模型

3 基于偏斜T噪聲的魯棒狀態(tài)估計(jì)

3.1 變分貝葉斯原理分析

基于貝葉斯定理

(11)

其中,z表示觀測數(shù)據(jù)的集合,模型的未知參數(shù)為φ={X,u,Λ},隱變量為X={x1,x2,…,xN},先驗(yàn)的超參數(shù)ψ={Δ,v}。則

p(z|φ,ψ)=p(z|φ)

(12)

根據(jù)式(12)可以將式(11)寫為

(13)

邊緣似然函數(shù)的對(duì)數(shù)logp(z|φ)表示為

logp(z|φ)=F(q(φ))+KL(q(φ)‖p(φ|z,ψ))

(14)

其中,F(q(φ))為自由能量函數(shù)[25-26]。

(15)

KL(q(φ)‖p(φ|z,ψ))是近似后驗(yàn)q(φ)與真實(shí)后驗(yàn)p(φ|z,ψ)之間的Kullback-Liebler散度[27]。

(16)

變分的思想是尋找一個(gè)簡單易求的概率分布q(φ)去逼近真實(shí)后驗(yàn)分布p(φ|z,ψ),即讓KL散度值達(dá)到最小。由于KL散度是非負(fù)的,當(dāng)且僅當(dāng)q(φ)=p(φ|z,ψ)時(shí)取值為零,最小化KL散度相當(dāng)于最大化自由能量函數(shù),因此可以通過使自由能量達(dá)到最大來實(shí)現(xiàn)目標(biāo)。

根據(jù)均值域逼近理論[28],q(φ)可以分解為

(17)

在VB步驟中,分解的真實(shí)后驗(yàn)分布與近似后驗(yàn)分布的最小KL散度見下式

(18)

其中,Q=qx(xk)qu(uk)qΛ(Λk)。

根據(jù)式(17),近似后驗(yàn)分布表示為

p(xk,uk,Λk|z1:k)≈qx(xk)qu(uk)qΛ(Λk)

(19)

則z,φ聯(lián)合分布的對(duì)數(shù)表示如下式所示

log(p(z|φ)p(φ|ψ))

(20)

3.2 遞歸VB卡爾曼濾波的偏斜T算法

通過求解logqx(xk),logqu(uk)和logqΛ(Λk)可以得到其估計(jì)的后驗(yàn)概率密度函數(shù)。

對(duì)于目標(biāo)狀態(tài)

logqx(xk)

(21)

其后驗(yàn)概率密度函數(shù)的形式如式(22)所示

qx(xk)=N(xk;xk|k,Pk|k)

(22)

其中,狀態(tài)xk|k,協(xié)方差矩陣Pk|k和卡爾曼濾波增益Kx的更新方程如下式所示。

(23)

Pk|k=(I-KkC)Pk|k-1

(24)

(25)

同理,對(duì)于模型的未知參數(shù)uk

logqu(uk)

(26)

其后驗(yàn)概率密度函數(shù)的形式如式(27)所示。

qu(uk)=TN(uk;uk|k,Uk|k)

(27)

其中,TN(·)表示截?cái)嗾龖B(tài)分布。qu(uk)的后驗(yàn)概率密度參數(shù)的更新方程[29]如下式所示

uk|k=Ku(zk-Cxk|k)

(28)

(29)

=Δ(Δ2+R)-1

(30)

對(duì)于模型的未知參數(shù)Λk

logqΛ(Λk)

(31)

其后驗(yàn)概率密度函數(shù)的形式如下式所示

(32)

其中

φk=R-1((zk-Cxk|k)(zk-Cxk|k)T+CPk|kCT)

(33)

其中,式(21)(26)(31)中的cx、cu和cΛ是和變量xk、uk和Λk相關(guān)的常數(shù)。

4 仿真研究

為了驗(yàn)證所提算法的有效性,分別采用卡爾曼濾波算法(KF),VBT算法,粒子濾波算法(PF)以及VBST算法對(duì)目標(biāo)跟蹤系統(tǒng)進(jìn)行仿真。

假設(shè)目標(biāo)的初始位置在x方向500m處,速度為1.5m/s,在y方向500m處,速度為1m/s。觀測站的位置在原點(diǎn)(0,0),粒子數(shù)設(shè)為1000。協(xié)方差矩陣P0=diag([2,10-2,2,10-2]),采樣時(shí)間為0.1s,仿真時(shí)間200s,過程噪聲為零均值的高斯白噪聲, 其協(xié)方差矩陣Q=diag([10-2,10-3,10-2,10-3])。

4.1 實(shí)驗(yàn)一:非時(shí)變偏斜T噪聲模型

假設(shè)量測噪聲的偏斜系數(shù)為0.9,自由度為10,R=2×I2×2。本文主要通過四種算法的均方根誤差(RMSE)來說明算法的濾波性能。

圖3為真實(shí)軌跡和KF算法、VBT算法、粒子濾波算法以及VBST算法的濾波軌跡。從圖中可以直觀看出VBST算法的濾波性能要明顯優(yōu)于另外三種算法,粒子濾波算法次之,KF和VBT的算法相對(duì)而言較差。

圖3 真實(shí)軌跡及四種算法濾波后軌跡

圖4和圖5是四種算法的位置均方根誤差和速度均方根誤差,從圖中可以看出,VBST算法具有更小的位置和速度均方根誤差,其濾波性能更好。主要原因是VBST算法不斷的更新隱變量和參數(shù),直到算法收斂,使得到的近似后驗(yàn)分布逼近真實(shí)后驗(yàn)分布。雖然VBT算法也用到了變分貝葉斯方法,但是T分布不能很好地?cái)M合重尾非對(duì)稱噪聲,所以在本系統(tǒng)的量測噪聲模型中,其估計(jì)的后驗(yàn)分布與真實(shí)后驗(yàn)分布相差太大,造成濾波性能下降。

圖4 四種算法的位置均方根誤差

圖5 四種算法的速度均方根誤差

4.2 實(shí)驗(yàn)二:時(shí)變偏斜T噪聲模型

為了驗(yàn)證所提算法的魯棒性,考慮時(shí)變的偏斜T分布噪聲模型,實(shí)驗(yàn)參數(shù)設(shè)置如下:t<150s設(shè)置偏斜系數(shù)為0.9,自由度為10,R=2×I2×2。t∈[150,200]s改變噪聲分布,增大噪聲協(xié)方差,設(shè)置偏斜系數(shù)為-0.6,自由度為40,R=5×I2×2,其仿真結(jié)果如圖6、7和8所示。

圖6 時(shí)變噪聲下的四種算法濾波軌跡

從圖6可以看出,改變噪聲分布后,系統(tǒng)的量測噪聲在x方向525m處突然增大,KF和VBT算法的濾波效果明顯變差,粒子濾波算法也有輕微的波動(dòng),但是VBST的濾波效果并未受到很大影響。

圖7和圖8為四種算法在不同量測噪聲下的位置和速度均方根誤差,從圖中可以看出,在150s改變噪聲分布后,KF和VBT算法的均方根誤差明顯增大,粒子濾波算法相比較KF和VBT算法而言,濾波效果較好,位置和速度的RMSE值也小,而VBST算法的位置均方根誤差和速度均方根誤差都比較小,說明在時(shí)變噪聲條件下本文所提算法和粒子濾波算法較于其它兩種算法而言魯棒性更強(qiáng)。

圖7 時(shí)變噪聲下的四種算法位置均方根誤差

圖8 時(shí)變噪聲下的四種算法速度均方根誤差

以實(shí)驗(yàn)一的仿真模型為例,主要通過求解變分下界函數(shù)L,也表示自由能量函數(shù)F,分析VBST算法的收斂性。設(shè)置迭代次數(shù)為30,進(jìn)行100次蒙特卡羅實(shí)驗(yàn),仿真結(jié)果如圖9所示,由下圖可以看出,L最終穩(wěn)定在-0.63±0.03,說明算法收斂。

圖9 變分自由能量變化過程

表1為單次蒙特卡羅實(shí)驗(yàn)四種算法的相對(duì)運(yùn)行時(shí)間對(duì)比,假設(shè)卡爾曼濾波算法運(yùn)算時(shí)間為單位1。

表1 不同算法的運(yùn)行時(shí)間

從表1可以看出,粒子濾波算法的計(jì)算量較大,用時(shí)較長。KF和VBT算法的運(yùn)算時(shí)間短,但濾波效果差。VBST算法濾波性能較好,運(yùn)行時(shí)間較粒子濾波而言小的多。從圖8圖9可以看出粒子濾波和VBST算法的濾波性能相似,但從算法的復(fù)雜度來看,粒子濾波的時(shí)間復(fù)雜度為O(n+Nn+2Nn2),VBST算法的時(shí)間復(fù)雜度為O(Mn),其中,n=200,N=1000,M=30,所以粒子濾波的計(jì)算要比VBST算法復(fù)雜得多。綜上分析,本文所提的VBST算法濾波精度更高,性能更加優(yōu)越。

5 結(jié)論

本文考慮量測噪聲為偏斜T分布時(shí)的目標(biāo)跟蹤問題,提出了一種基于變分貝葉斯理論的VBST濾波算法。實(shí)驗(yàn)結(jié)果表明,在量測噪聲為偏斜T分布時(shí),本文所提算法具有較高的估計(jì)精度和較好的濾波效果。在時(shí)變噪聲條件下,VBST算法也有很好的魯棒性。鑒于VBST算法的優(yōu)越性能,未來考慮非線性系統(tǒng)模型和混合不確定量測噪聲模型。