鄭作虎　歐陽東升　李黎明

(95806部隊(duì)　北京　100076)

0　引言

在雷達(dá)信號處理領(lǐng)域，如何實(shí)現(xiàn)對目標(biāo)運(yùn)動狀態(tài)的準(zhǔn)確估計是研究的重點(diǎn)內(nèi)容，由于目標(biāo)運(yùn)動模型的非線性非高斯特性，非線性濾波算法在雷達(dá)目標(biāo)跟蹤中有著廣泛的應(yīng)用，目前，主要包括擴(kuò)展卡爾曼濾波[1](Extend Kalman Filter，EKF)、不敏卡爾曼濾波[2](Unscented Kalman Filter, UKF)和粒子濾波算法[3](Particle Filter, PF)。其中，EKF是一種次優(yōu)濾波算法，通過一階泰勒級數(shù)展開近似非線性系統(tǒng)，進(jìn)而得到線性系統(tǒng)的最優(yōu)解，但一階泰勒級數(shù)近似的線性化誤差較大，導(dǎo)致算法濾波精度較差；UKF利用不敏變換較好解決了非線性問題，通過一系列采樣點(diǎn)逼近概率密度函數(shù)，不需要進(jìn)行線性化處理，但在非高斯噪聲環(huán)境下性能下降；PF是一種基于蒙特卡羅方法和遞推貝葉斯估計的統(tǒng)計濾波方法，其基本思想是通過序貫重要性采樣(Sequential Importance Sample, SIS)算法采樣一組隨機(jī)加權(quán)粒子來近似后驗(yàn)概率密度函數(shù)(Probability Density Function, PDF)，以均值運(yùn)算代替積分運(yùn)算，最終獲得狀態(tài)的最小均方估計，當(dāng)粒子數(shù)足夠時，可以認(rèn)為這些粒子的統(tǒng)計特性近似于狀態(tài)的后驗(yàn)概率分布的統(tǒng)計特性[4]。

與EKF和UKF算法相比，PF算法適用于任何非高斯非線性系統(tǒng)環(huán)境，但存在粒子退化問題[5]。解決粒子退化問題的一種方法就是構(gòu)造好的重要性密度函數(shù)[6]，選擇標(biāo)準(zhǔn)是逼近于真實(shí)的后驗(yàn)PDF，且便于采樣；解決粒子退化問題的另一種方法就是重采樣[7]，重采樣方法在一定程度解決了粒子退化問題，但大量復(fù)制權(quán)值大的粒子也降低了粒子的多樣性，導(dǎo)致了粒子枯竭問題。

1　算法基本原理和實(shí)現(xiàn)步驟

1.1　EKF算法

設(shè)非線性非高斯?fàn)顟B(tài)空間模型為

狀態(tài)方程：

xk=f(xk-1)+wk

(1)

量測方程：

zk=h(xk)+vk

(2)

式中：xk為k時刻N(yùn)維狀態(tài)矢量；yk為k時刻M維量測矢量；wk、vk分別為k時刻系統(tǒng)噪聲和量測噪聲，相互獨(dú)立，協(xié)方差矩陣分別為Qk和Rk；f(·)，h(·)為有界非線性映射。

(3)

(4)

(5)

(6)

由式(4)、(6)可知，非線性系統(tǒng)通過泰勒級數(shù)展開近似為線性系統(tǒng)，根據(jù)KF算法原理，可得EKF算法的基本公式為：

狀態(tài)預(yù)測：

(7)

協(xié)方差預(yù)測：

(8)

量測預(yù)測：

(9)

相應(yīng)的協(xié)方差

(10)

增益：

(11)

狀態(tài)估計：

(12)

協(xié)方差估計：

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

(13)

EKF算法的實(shí)現(xiàn)步驟為：

(1)系統(tǒng)初始化

(14)

(3)利用式(8)計算狀態(tài)預(yù)測協(xié)方差Pk|k-1；

(4) 利用式(11)計算濾波增益Kk；

擴(kuò)展卡爾曼濾波算法由于其針對非線性模型的簡單直觀性，近年來已經(jīng)得到了廣泛的應(yīng)用，但它也有以下不足：(1)由于采用一階泰勒級數(shù)近似線性化，線性化誤差較大，導(dǎo)致EKF濾波算法精度變差，甚至發(fā)生濾波發(fā)散問題；(2)EKF濾波算法需要計算非線性系統(tǒng)的雅可比矩陣，但通常不易實(shí)現(xiàn)；(3)EKF算法要求非線性系統(tǒng)必須是連續(xù)可微的，限制了其應(yīng)用范圍[8]。

1.2　UKF算法

不敏卡爾曼濾波基于不敏變換(Unscented Transformation，UT)，在UT過程中，重點(diǎn)是采樣Sigma采樣點(diǎn)，實(shí)際中一般采用對稱采樣，計算簡單，濾波精度可達(dá)到2階泰勒級數(shù)展開。UT主要實(shí)現(xiàn)步驟為：

Step1首先計算2nx+1個Sigma采樣點(diǎn)和相對應(yīng)的權(quán)值

(15)

(16)

Step2每個Sigma采樣點(diǎn)通過非線性函數(shù)傳播

yj=f(χj)j=0,…,2nx

(17)

Step3ny的估計均值和協(xié)方差為

(18)

(19)

基于UT，UKF算法的實(shí)現(xiàn)步驟為：(1)利用式(14)進(jìn)行系統(tǒng)初始化。

(2)利用式(15)、(16)計算狀態(tài)Sigma采樣點(diǎn)ξj(k-1|k-1)和其相應(yīng)的權(quán)值Wj，根據(jù)式(1)計算Sigma采樣點(diǎn)的一步預(yù)測值：

ξj(k|k-1)=f(ξj(k-1|k-1))

(20)

(3)利用ξj(k|k-1)及權(quán)值Wj，根據(jù)式(18)、(19)計算狀態(tài)預(yù)測值和相應(yīng)的協(xié)方差預(yù)測值

(21)

(22)

(4)根據(jù)式(2)，計算量測采樣點(diǎn)的一步預(yù)測值

ζj(k|k-1)=h(ξj(k|k-1))

(23)

(5)根據(jù)式(18)、(19)計算量測預(yù)測值和相應(yīng)的協(xié)方差預(yù)測值

(24)

(25)

(6)計算濾波增益值

(26)

相對EKF算法，UKF算法濾波精度有了很大的提高，應(yīng)用范圍也寬廣了很多，但算法僅具有二階濾波精度，且僅適用于高斯系統(tǒng)。

1.3　PF算法

根據(jù)蒙特卡羅方法，對于k時刻目標(biāo)狀態(tài)的后驗(yàn)PDF，可用一系列加權(quán)粒子表示為

(27)

(28)

(29)

假定系統(tǒng)當(dāng)前狀態(tài)和未來量測無關(guān)，且系統(tǒng)狀態(tài)服從馬爾科夫過程，則π(x0:k|z1:k)可分解為

π(x0:k|z1:k)=π(xk|xk-1,zk)π(x0:k-1|z1:k-1)

(30)

同理，狀態(tài)后驗(yàn)PDF可表示為

∝p(zk|xk)p(xk|xk-1)p(x0:k-1|z1:k-1)

(31)

將式(30)、式(31)代入式(29)，則重要性權(quán)值可表示為

(32)

在PF中，重要性密度函數(shù)與真實(shí)的后驗(yàn)PDF的逼近程度直接影響著算法濾波性能，逼近程度越高，粒子權(quán)值方差越小，濾波性能越好。研究表明，使粒子權(quán)值方差最小的最優(yōu)重要性密度函數(shù)為

(33)

將式(33)代入式(32)可得

(34)

(35)

將式(35)代入式(32)可得

(36)

將權(quán)值歸一化

(37)

則系統(tǒng)狀態(tài)的最小均方估計為

(38)

直接根據(jù)狀態(tài)先驗(yàn)PDF采樣粒子，容易實(shí)現(xiàn)，權(quán)值計算也僅與前一時刻的粒子權(quán)值和似然函數(shù)有關(guān)，計算簡單。但是，在量測精度較高的場合(量測噪聲較小)，狀態(tài)先驗(yàn)PDF與狀態(tài)后驗(yàn)PDF存在較大偏差，這樣根據(jù)狀態(tài)先驗(yàn)PDF采樣的粒子與真實(shí)狀態(tài)偏差較大，隨著迭代過程的進(jìn)行，權(quán)值小的粒子越來越多，權(quán)值大的粒子越來越少，導(dǎo)致權(quán)值方差越來越大，加劇了粒子退化問題。

為了解決粒子退化問題，標(biāo)準(zhǔn)PF通常在SIS之后進(jìn)行重采樣。通過復(fù)制大權(quán)值粒子，去掉小權(quán)值粒子，重采樣算法從估計的后驗(yàn)PDF上重新抽取N個獨(dú)立同分布的粒子集，其相應(yīng)權(quán)值為1/N，能夠降低粒子權(quán)值方差，提高有效粒子數(shù)，改善濾波精度，但重采樣算法也導(dǎo)致了樣本枯竭問題，損失了粒子的多樣性。

標(biāo)準(zhǔn)PF算法的實(shí)現(xiàn)步驟為：

(1)初始化；

(2)計算k時刻粒子權(quán)值；

(3)重采樣；

(4)k時刻系統(tǒng)狀態(tài)估計。

2　實(shí)驗(yàn)與分析

采用文獻(xiàn)[9]中典型的一維非線性非高斯模型，對EKF、UKF、PF三種算法的濾波性能進(jìn)行比較分析。狀態(tài)方程和量測方程分別為：

k=0,1,…,k-1

(39)

(40)

參數(shù)設(shè)置如下：狀態(tài)噪聲wk～Gamma(3,2)，量測噪聲vk～N(0,1)，初始狀態(tài)x0=1，初始狀態(tài)方差p(x0)=1，時間序列長度K=60，UKF的參數(shù)取κ=2，PF粒子數(shù)N=100。本文對各濾波方法進(jìn)行M=100次的Monte Carlo仿真，各濾波算法每一時刻的均方根誤差為

(41)

為了比較EKF、UKF、PF算法濾波性能，圖1給出了三種算法濾波狀態(tài)估計值與真實(shí)值的比較曲線。從圖中可以看出，PF算法的濾波狀態(tài)最接近于真實(shí)狀態(tài)，濾波性能優(yōu)于EKF、UKF算法；因較好地解決了非線性問題，UKF算法濾波性能優(yōu)于EKF；EKF存在線性化誤差，濾波性能最差。

為了更精確比較三種算法的濾波性能，圖2給出了Monte Carlo仿真100次時三種算法均方根誤差曲線的比較，表1給出了均方根誤差的均值和方差比較。從中可以看出，PF算法的均方根誤差最小，表明最接近于真實(shí)的狀態(tài)值，且最平穩(wěn)，濾波性能最優(yōu)，UKF算法次之，EKF算法最差。

3　結(jié)束語

本文詳細(xì)介紹了基本PF算法的基本原理和實(shí)現(xiàn)步驟，并在典型的一維非線性非高斯模型下綜合比較了PF與EKF、UKF三種算法的濾波性能。實(shí)驗(yàn)結(jié)果表明，在非線性非高斯條件下，PF算法的濾波性能優(yōu)于其他兩種算法，且隨著粒子數(shù)的增多，算法濾波性能更優(yōu)。

表1EKF、UKF、PF算法均方根誤差比較

算法均方根誤差均值方差EKF36．32102291UKF16．615929．8329PF1．12300．3190

參考文獻(xiàn):

[1]夏楠,邱天爽,李景春,等.一種卡爾曼濾波與粒子濾波相結(jié)合的非線性濾波算法[J].電子學(xué)報,2013,41(1):148-152.

[2]申正義,王晴晴,王洪林,等. 基于修正SUKF的固定單站無源定位跟蹤算法[J].火控雷達(dá)技術(shù),2015,44(1):38-42.

[3]王首勇,萬洋,劉俊凱,等.現(xiàn)代雷達(dá)目標(biāo)檢測理論與方法(第二版)[M].北京:科學(xué)出版社,2015.

[4]萬洋,王首勇.非線性濾波算法的性能分析[J]. 空軍雷達(dá)學(xué)院學(xué)報,2010,24(2):111-114.

[5]鄭作虎,王首勇.一種基于GHF的高斯粒子濾波算法[J].控制與決策,2014,29(9):1698-1702.

[6]張俊根.粒子濾波及其在目標(biāo)跟蹤中的應(yīng)用研究[D].西安:西安電子科技大學(xué),2011.

[7]趙玲玲.目標(biāo)跟蹤中的粒子濾波與概率假設(shè)密度濾波研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2011.

[8]陳亮.機(jī)動目標(biāo)跟蹤關(guān)鍵技術(shù)研究[D].哈爾濱:哈爾濱工程大學(xué),2012.

[9]gordon n, salmond d, smith a.Novel Approach to Nonlinear/Non-Gaussian Bayesian State Estimation[J].IEE Proceeding on Radar,Sonar and Navigation,1993,140(2):107-113.