變結構交互式多模型濾波和平滑算法

陳維義, 何 凡,*, 劉國強, 毛偉偉

(1. 海軍工程大學兵器工程學院, 湖北 武漢 430030; 2. 海軍士官學校兵器系, 安徽 蚌埠 233000)

0 引 言

目前,卡爾曼濾波、粒子濾波、基于機器學習的濾波及其改進形式已被廣泛應用在目標跟蹤、參數估計和狀態預測等領域[1-4]。對于基于單模型的濾波方法而言,當目標機動十分復雜時,會出現濾波模型與目標機動模型不匹配的情況,從而導致濾波精度大幅降低[5]。

針對這一問題,交互式多模型(interacting multiple model, IMM)[6-8]展現出更加優越的性能。近些年,IMM算法得到了很多學者的關注,很多學者在結構和參數等層面對IMM算法進行了改進。文獻[9]基于恒速和當前統計模型設計了IMM算法,首先利用最小二乘法估計當前統計模型的平均速度,然后將當前統計模型應用于IMM算法。該方法提高了模型精度,從而提高了濾波精度。文獻[10]從多個方面對IMM算法進行了改進,包括采用改進的卡爾曼濾波器作為子濾波器、不同模型之間的非對稱狀態估計和基于熵的模型概率更新公式。文獻[7,11]提出了IMM的一種替代方法,模型集中的模型由勻加速模型構成,降低了模型集的復雜程度。在此基礎上,文獻[12]提出了一種自適應IMM算法,首先利用濾波器對目標的加速度進行估計,然后選取估計加速度附近的值構建模型。該方法可以減少模型集中模型的數量,在降低計算量的同時提高了模型精度。文獻[13]基于二階馬爾可夫鏈提出了一種二階IMM算法,該算法利用了更多的先驗信息,提高了濾波精度。由于上述IMM模型集中的模型種類和數量不變,所以又稱其為固定結構IMM(fixed structure IMM, FSIMM)。

為了避免因模型不匹配而造成的精度誤差問題,在使用IMM算法時,應使用盡可能多的模型覆蓋目標機動模型。但值得注意的是,單個模型集中模型過多,同樣會降低濾波精度[14]。在此背景下,變結構IMM(variable structure IMM,VSIMM)應運而生。VSIMM經過不斷的發展和改進,大致可以分為4類:模型組切換(model group switching, MGS)、可能模型集(likely mode set, LMS)、期望模型增強(expected mode augmentation, EMA)和自適應網格(adaptive grid,AG)[15]。其中,MGS將模型集分為模型子集,一個時間步長只選擇一個模型子集進行估計,模型子集之間根據模型子集轉移概率進行切換[16]。LMS將模型分為3種類型:不可能的、重要的和主要的模型,在每個時間步長里,用于估計的模型子集由主要的模型和接近主要的模型構成[17-18]。與MGS類似,EMA將一個大的模型集分成小的模型子集,然后計算下一個時間步長所有模型子集的概率,選擇概率最大的模型子集用于估計[19]。AG算法結合圖理論,將所有的模型構成一個網格,利用先驗信息和當前數據得到一個局部細化網格,細化網格中的模型構成候選模型子集,然后根據規則選擇模型,構成下一時刻用于濾波的模型子集[20]。

值得注意的是,文獻[19]所提到的EMA方法在計算似然函數和模型子集概率時,對公式進行了近似。針對這一問題,本文首先設計了一種新的VSIMM(novel VSIMM, NVSIMM)濾波算法,給出了精確的數學模型。另外,在前向NVSIMM濾波的基礎上,對數據進行平滑處理,設計了VSIMM平滑(VSIMM smoothing, VSIMMS)算法,前向NVSIMM和后向VSIMMS相結合即為VSIMM濾波和平滑(VSIMM filtering and smoothing, VSIMMFS)算法。

本文安排如下:第1節對濾波和平滑問題進行了描述,并給出了IMM濾波算法和IMM平滑算法的數學模型;第2節設計了VSIMMFS算法,包括前向NVSIMM濾波和后向VSIMM平滑兩個部分;第3節對本文所提算法的有效性進行仿真驗證和討論;最后,第4節給出了本文的結論。

1 IMM濾波和平滑算法

1.1 問題描述

假設目標可能有r個運動模型,模型集記為Ω={M1,M2,…,Mr},模型之間的轉換概率矩陣為

(1)

式中:pij(1≤i≤r,1≤j≤r)為模型i到模型j的轉移概率。

系統離散化狀態方程如下:

(2)

模型j的量測方程為

(3)

1.2 前向IMM濾波

IMM算法框圖如圖1所示。

圖1 IMM算法框圖Fig.1 Block diagram of the IMM algorithm

IMM算法一般可以分為以下4個步驟[6-7]。

步驟 1前向輸入交互

(4)

(5)

(6)

步驟 2濾波

分別基于每一個模型進行卡爾曼濾波:

步驟 2.1預測狀態

(7)

步驟 2.2預測協方差矩陣

(8)

步驟 2.3計算卡爾曼增益

(9)

步驟 2.4濾波估計

(10)

步驟 2.5濾波值的協方差

(11)

步驟 3模型概率更新

模型概率更新公式為

(12)

(13)

(14)

步驟 4輸出

IMM濾波的狀態估計和協方差估計分別為

(15)

(16)

1.3 后向IMM平滑

在前向IMM濾波的基礎上,對數據進一步進行平滑處理,得到后向IMM平滑算法,其算法框圖如圖2所示,包括以下幾個步驟[21]。

步驟 1后向狀態交互

后向混合狀態和對應的協方差為

(17)

(18)

(19)

(20)

步驟 2卡爾曼平滑[22]

平滑值和對應的協方差為

(21)

(22)

步驟 3平滑模型概率更新

平滑后的模型概率計算公式為

(23)

步驟 4IMM平滑狀態交互

后向IMM平滑的狀態輸出和對應的協方差為

(24)

(25)

圖2 IMM平滑算法框圖Fig.2 Block diagram of the IMM smoothing algorithm

2 VSIMMFS算法

VSIMMFS算法包括前向NVSIMM濾波和后向VSIMM平滑兩個部分,這兩個部分分別在后續兩個小節中展開。

2.1 VSIMM濾波

為了避免模型集中模型數量過多,減少模型誤差,本節設計了NVSIMM算法。NVSIMM濾波算法包括多個模型子集,每個模型子集獨立運行IMM,選取概率最高的模型子集的估計結果作為最終的估計狀態輸出。NVSIMM算法框圖如圖3所示。

圖3 NVSIMM算法框圖Fig.3 Block diagram of NVSIMM algorithm

算法包括以下4個步驟。

步驟 1并行獨立IMM狀態估計

步驟 2計算每個模型子集的似然函數

模型子集的似然函數為

(26)

(27)

步驟 3計算模型子集概率

記模型子集概率p(Πk(m)|Zk)=κk(n),則有:

(28)

式中:p(zk|Zk-1,Πk(n))在步驟2中已經求出,分母為歸一化因子,p(Πk(n)|Zk-1)根據全概率公式有:

(29)

式中:κk-1(m)=p(Πk-1(m)|Zk-1)為上一時刻的模型子集概率,應用馬爾可夫性質可知模型子集轉移概率與觀測值無關,即:

(30)

步驟 4NVSIMM濾波狀態估計

選擇概率最大的模型子集的IMM估計結果作為最終的狀態估計輸出,首先求出概率最高的模型子集編號:

(31)

從而可以得到最終的狀態估計和對應的協方差為

(32)

Pk|k=Pk|k(nm)

(33)

2.2 VSIMMS算法

本節在NVSIMM的基礎上,對數據進行進一步平滑,得到VSIMMS算法。VSIMMS算法框圖如圖4所示。

圖4 VSIMMS算法框圖Fig.4 Block diagram of VSIMMS algorithm

算法包括以下4個步驟。

步驟 1并行獨立IMM狀態平滑

步驟 2計算平滑后的模型子集的似然函數

應用全概率公式和馬爾可夫性質,可以得到模型子集的似然函數為

(34)

式中:p(zk|Zk-1,Πk(m))為前向NVSIMM濾波過程中模型子集的似然函數,應用全概率公式和馬爾可夫性質,p(Πk(m)|Πt(n))可以轉化為

(35)

(36)

(37)

步驟 3計算平滑后的模型子集概率

根據貝葉斯定理可以得到模型子集概率為

(38)

(39)

步驟 4VSIMMS狀態估計

概率最高的模型子集編號為

(40)

判斷概率最高的模型子集在Tswitch時間段內是否發生變化。若發生變化,VSIMMS的狀態估計和對應的協方差為

(41)

Pt|k=Pt|k(Ω)

(42)

若不發生變化,VSIMMS的狀態估計和對應的協方差為

(43)

(44)

3 仿真與討論

為了驗證VSIMMFS算法的有效性,本節首先針對三維空間中的目標運動軌跡展開分析討論。目標的初始位置、速度和加速度分別為(10 km, 40 km, 30 km)、(300 m/s, 0, 0)和(0, 0, 0)。目標的機動參數如下:

(1) 30 s的勻速直線運動;

(2) 30 s的勻加速運動,加速度大小為(-10 m/s2, -10 m/s2, -10 m/s2);

(3) 30 s的勻速轉彎運動,轉彎角速率為(0.2 rad/s, 0.2 rad/s, 0);

(4) 30 s的勻速直線運動。

目標運動軌跡如圖5所示。

圖5 目標運動軌跡Fig.5 Target motion trajectory

利用蒙特卡羅法分析對比IMM、NVSIMM和VSIMMFS的性能,其中IMM算法包括3個模型:勻速直線運動(constant velocity, CV)模型、勻加速運動(constant acceleration, CA)模型和勻速轉彎(constant turn, CT)模型[23-25]。VSIMMFS算法包括兩個模型子集,第一個模型子集包括CV、CA,第二個模型子集包括CV、CT;NVSIMM濾波算法包括的模型子集和VSIMM相同。可以注意到,VSIMMFS和NVSIMM僅需兩個原始模型子集即可包括所有的模型,但這兩種算法的原始模型子集并不包括模型子集CA、CT。設模型轉移概率矩陣為

假設系統過程噪聲的標準差為1 m/s2,測量噪聲的標準差為50 m,Tswitch=20 s,仿真結果如圖6～圖9所示。圖6給出了IMM算法中模型概率變化曲線,其中,概率最大的曲線所對應的模型與目標實際機動模型相同。圖7給出了NVSIMM和VSIMMFS中模型子集概率曲線,由圖7 可知,平滑后(VSIMMFS)的模型子集概率大于NVSIMM的模型子集概率,說明平滑后的模型子集與目標機動模型匹配程度更高。另外,對比圖6和圖7可知,圖6中概率最大的曲線所對應的模型屬于圖7中概率高的模型子集,也就是說,概率最高的模型子集始終包括目標實際的機動模型。圖8給出了IMM、NVSIMM和VSIMMFS對目標運動軌跡的跟蹤效果,這3種方法都能夠對目標軌跡進行較好的估計,三者之間的差異具體體現在圖9中。圖9給出了IMM、NVSIMM和VSIMMFS對目標位置估計的均方根誤差(root mean square error, RMSE)曲線。從圖9中可以看出,VSIMMFS的RMSE最小,VSIMMFS次之,IMM的RMSE最大。值得注意的是,概率最大的模型子集發生變化時,NVSIMM算法的性能會短暫下降,而VSIMMFS算法對位置估計的RMSE與IMM相同,原因在于在設計VSIMMFS算法時,當最大概率的模型子集短暫發生時,利用基于模型全集的IMM估計值作為VSIMMFS的估計值,從而避免了VSIMMFS算法短暫的性能下降。

圖6 IMM算法中的模型概率Fig.6 Model probability in IMM algorithm

圖7 NVSIMM和VSIMMFS中模型子集概率Fig.7 Model subset probabilities in NVSIMM and VSIMMFS

圖8 目標運動軌跡跟蹤效果Fig.8 Effect of target trajectory tracking

圖9 位置跟蹤的RMSEFig.9 RMSE of position tracking

為了進一步體現VSIMMFS算法的優越性,將VSIMMFS與基于深度強化學習的目標跟蹤算法(maneuvering target tracking based on deep reinforcement learning, MTTDRL)[27]、粒子濾波(particle filter, PF)[28]、基于粒子群優化的PF(particle swarm optimization algorithm-based PF, PSO-PF)[29]、基于混沌PSO-PF(chaos PSO-PF, CPSO-PF)[26,30]進行對比。為了保證對比的有效性,仿真參數參考文獻[26]和文獻[27]設計,目標的運動軌跡如圖10所示。圖11分析對比了上述算法對目標位置的跟蹤性能。由圖11(a)可知,相比于PF、PSO-PF和CPSO-PF,VSIMMFS的位置跟蹤的RMSE較小;從圖11(b)可以看出,利用VSIMMFS對目標跟蹤時的位置偏差更小。所以綜上可知,VSIMMFS的跟蹤性能優于上述其他算法。

圖10 目標運動軌跡Fig.10 Target motion trajectory

圖11 算法性能對比Fig.11 Performance comparison of algorithms

4 結 論

為了進一步提高對機動目標的跟蹤精度,本文設計了一種VSIMMFS算法。該方法包括前向VSIMM濾波和后向VSIMMS兩個部分,通過并行獨立運行IMM濾波和IMM平滑,選擇概率最大的模型子集的估計結果作為最終的狀態估計值。通過與IMM、NVSIMM、PF、PSO-PF、CPSO-PF和MTTDRL算法比較,VSIMMFS具有更好的跟蹤性能。