許江湖,黃 亮,劉 忠
(海軍工程大學 電子工程學院,湖北 武漢430033)
粒子濾波又稱序貫蒙特卡羅方法,是一種基于蒙特卡羅方法和遞推貝葉斯估計的統計濾波方法。粒子濾波完全突破了Kalman 濾波理論框架,對系統的過程噪聲和量測噪聲沒有任何限制。粒子濾波(PF)作為解決非線性、非高斯動態系統的參數估計和狀態濾波問題的有效方法,近年來隨著計算機運算能力的急劇增長和計算成本的不斷降低,已成為研究熱點[1-3]。
粒子退化是標準粒子濾波算法的主要缺陷。采用重采樣方法在一定程度上可以抑制粒子退化現象的發生。然而重采樣帶來的新問題是:權值越大的粒子的子代越來越多,而權值較小的粒子被剔除,最糟糕的情形是新的粒子集實際都是一個權值最大的粒子的子代,即所謂“樣本枯竭”現象,從而導致粒子集的多樣性變差,不足以用來近似表征后驗密度,難以保證估計精度,特別是在樣本受限條件下,這種粒子多樣性減弱對濾波精度的影響更為突出,甚至導致濾波發散現象[4]。為此,不少學者嘗試利用進化算法(包括遺傳算法、進化策略、進化規劃)來解決粒子濾波的退化問題[4-8],并取得了一定的理論成果。2006 年,Mehrabian 和Lucas 提出了一種新穎的數值優化算法,稱為入侵野草(雜草)優化(IWO)算法[9]。同進化算法一樣,IWO算法也是一種隨機搜索仿生學優化算法,該算法模仿了野草入侵的種子空間擴散、生長、繁殖和競爭性消亡的基本過程,具有很強的魯棒性和自適應性。與進化算法相比,IWO 算法簡單易于實現,不需要遺傳操作算子,能簡單而有效地收斂于問題的最優解,是一種強有力的智能優化算法。研究表明,IWO 算法在性能上優于進化算法和人工蜂群克隆算法,并且在高維問題上優于粒子群算法[10]。因此,IWO 算法自提出以來,一直受到學者們的關注和運用[10-15]。
本文將IWO 算法與粒子濾波算法相結合,提出基于IWO 算法的粒子濾波算法(IWOPF)。該算法將IWO 應用于粒子重采樣,以保證粒子的有效性和多樣性。最后通過計算機仿真比較該算法與基于進化算法的粒子濾波算法以及傳統粒子濾波算法的濾波性能和計算時。
離散時間非線性動態系統的狀態方程和量測方程分別為:

式中:k 為采樣時刻;xk∈Rnx為k 時刻的系統狀態向量;fk:Rnx×Rnv→Rnx為系統狀態xk-1的非線性函數;vk為過程噪聲,其協方差為Qk,均值為0;nx和nv分別為狀態和過程噪聲的維數;zk∈Rnz為k 時刻系統狀態的量測向量;hk:Rnx×Rnn→Rnz為系統狀態xk的非線性函數;nk為量測噪聲,其協方差為Rk,均值為0,且與vk相互獨立;nz和nn分別為量測值和量測噪聲的維數。
濾波過程的任務是通過可獲得的系統觀測值zk估計出系統狀態xk,也即要求得到系統狀態的后驗概率分布p(xk| z1:k)。粒子濾波的核心思想是利用一系列隨機樣本的加權和表示狀態后驗概率密度:

通常很難從p(x0:k| z1:k)中直接采樣,解決方法是先從一個已知的、容易采樣且盡量接近后驗概率分布的重要性概率密度函數q(x0:k| z1:k)中采樣,這時權值的遞推公式為:

標準粒子濾波算法選擇最易于實現的先驗概率密度作為重要密度函數,即

這時,式(4)可以進一步簡化為:

當前時刻粒子的權重被評估后,通過引入重采樣技術來改善粒子退化問題。
IWO 算法中,“野草”表示隨機產生的可行解,“種子”表示野草的后代,“種群”表示所有野草的集合。一般而言,IWO 算法的執行過程要經歷以下4 個步驟:
1)初始化種群
在這個步驟中,需要確定最大種群數目Pmax、最大迭代次數itermax、問題維數D、最大和最小可生成種子數Smax和Smin、非線性調節指數n、標準差初始值σinit和最終值σfinal以及初始搜索空間X,并隨機生成N0個解。
2)繁殖
各野草產生的種子個數為:

式中:f 為當前野草的適應度值;fmax和fmin分別為當前種群中野草對應的最大和最小適應度值。
從式(7)可以看出,種群中的成員能夠散播的種子數是根據該成員的適應度值以及種群所有成員的最大和最小適應度值來決定的。在用優化算法求解問題過程中,可能會直覺地認為可行解比不可行解具有更好的適應度值,因此會放棄不可行解。然而,這種做法忽視了進化算法具有概率性,其產生的解是在不可行和可行之間不斷轉化的。Mehrabian和Lucas 認為,在進化過程中一些不可行解可能帶有更多的有用信息,如果它可能穿過一個不可行區域(特別是非凸的可行搜索空間),則系統更容易到達最優點。因此,在IWO 算法中,繁殖過程按照自然界中的繁殖法則,給予不可行的個體生存和繁殖的機會,只是這種機會相對較小。
3)空間散布
空間散布體現了算法的隨機性和適應性。種群產生的種子個體按正態分布N(0,σ)在其父代野草個體附近的D 維空間進行散布。每次迭代對應不同的標準差,隨著迭代的進行,標準差從初始值σinit減小到最終值σfinal,當前迭代的標準差為:

式(8)確保了在較遠區域產生種子的概率以非線性的方式逐漸降低,以使適應度值好的個體聚集在一起。
4)競爭性生存
經過數代繁殖后,產生的后代數目將超過環境資源的承受能力,通過預先設定的最大種群數目Pmax對種群數量進行控制。在算法迭代過程中,種群中的所有野草及其后代按適應度值從大到小進行排序,只有適應度值最好的前Pmax個個體可以生存下來,其余的個體被清除。
綜合上述步驟可知:在算法運行過程中,一個野草個體可以有多個種子后代,種子按正態分布在父代野草周圍。只有相對較優秀的種子才有機會進入下一輪迭代,一旦種子存活下來進入一輪迭代,所有存活下來的野草和種子共同構成一輪迭代的種群。這種機制給予那些適應度值低的個體繁殖的機會,如果它們的后代適應度值更好,這些后代就可以生存下來。當算法運行到預先設定的最大迭代次數itermax時,算法終止。
前面已經述及,粒子退化是標準粒子濾波算法的主要缺陷。而粒子濾波與IWO 算法的相似之處在于它們都有一個初始化的群體,每個個體代表系統的一個可能解,這些個體都根據一定的規則進行狀態轉變,它們都對適應度較高的個體進行復制。因此,用IWO 算法來解決粒子濾波的退化問題具有可行性。
IWOPF 算法將IWO 算法應用于粒子重采樣以解粒子濾波的退化問題。即通過IWO 算法中的繁殖、空間散布、競爭性生存等操作來保證粒子的多樣性和有效性。具體而言,IWOPF 算法的執行步驟如下:
1)初始化
設置粒子數目N,從先驗分布p(x0)中采集粒子令確定最大迭代次數itermax、問題維數D、最大和最小可生成種子數Smax和Smin、非線性調節指數n、標準差初始值σinit和最終值σfinal。而最大種群數目Pmax= N。
for k = 1,2,…。
2)重要性采樣
3)更新權值
根據當前觀測值zk,計算每個粒子的權值,并歸一化:

4)基于IWO 算法的重采樣

②繁殖:根據式(7)計算每個粒子產生的種子數;
③空間散布:根據式(8)計算當前迭代的標準差;
④競爭性生存:把所有粒子及其后代的適應度值從大到小進行排序,只有適應度值最好的前Pmax= N 個粒子可以生存下來,并重新定義這些粒子為| i = 1,2,…,N};
⑤判斷算法運行是否達到預先設定的最大迭代次數itermax,如果是,則轉步驟5,并令否則轉步驟4 中的①;
5)狀態估計

為了驗證本文算法的有效性,通過2 個仿真實例比較本文提出的IWOPF、標準PF、文獻[4]提出的基于進化采樣的粒子濾波 (ESPF)以及文獻[5]提出的遺傳粒子濾波(GPF)的濾波性能以及仿真時間。2 個仿真實例采用的運動模型以及初始條件均相同,但量測模型不同。
運動模型為:

其中:


采樣周期T = 1 s,仿真時間為80 s,轉彎角速度為目標初始位置為(7,8)km,初始速度(-0.08, -0.06)km/s,系統過程噪聲是均值為0,方差為Q = diag {[3,3]}的加性高斯白噪聲;
仿真實例1 的量測模型為:

其中,量測噪聲是均值為0,方差為R = diag{[200,0.005]}的加性高斯白噪聲;
仿真實例2 的量測模型為:

其中,量測噪聲是均值為0,方差為R = 0.005 的加性高斯白噪聲。
進行50 次蒙特卡羅仿真。仿真中,粒子初始分布為高斯分布,其均值為[7.3,- 0.06,7.7,-0.08]× 103,方差為diag[100,4,100,4]。IWOPF有關參數如表1 所示。圖1 和圖2 分別為粒子數N =400 時,蒙特卡羅仿真次數為50 次,實例1 基于X 軸方向和Y 軸方向目標位置的均方根誤差(RMSE)曲線。表2 和表3 分別給出了不同粒子數目情況下,4 種濾波算法對實例1 的濾波平均RMSE 和50次仿真所用的時間。其中,平均RMSE 的定義如下[16]:平均

式中:MC 為仿真次數;L 為仿真長度;x(k)為真實值;)為估計值。
表4 和表5 則是針對實例2 的濾波情況。

圖1 實例1 的X 軸方向位置RMSE 曲線Fig.1 RMSE in position on X axis for example 1

圖2 實例1 的Y 軸方向位置RMSE 曲線Fig.2 RMSE in position on Y axis for example 1

表1 IWOPF 參數表Tab.1 IWOPF parameter values

表2 不同粒子數目情況下仿真實例1 的位置估計平均RMSE 比較Tab.2 The mean comparison of RMSE in position of example 1 for different sample size 單位:m

表3 不同粒子數目情況下仿真實例1 的50 次仿真所用時間Tab.3 The computation time of 50 Monte Carlo simulations of example 1 for different sample size 單位:s

表4 不同粒子數目情況下仿真實例2 的位置估計平均RMSE 比較Tab.4 The mean comparison of RMSE in position of example 2 for different sample size 單位:m

表5 不同粒子數目情況下仿真實例2 的50 次仿真所用時間Tab.5 The computation time of 50 Monte Carlo simulations of example 2 for different sample size 單位:s
從圖1 ~圖2、表2 及表4 可看出,本文提出的IWOPF 算法明顯優于另外3 種粒子濾波算法。而從表3 和表5 可看出,隨著粒子數目的增加,標準PF計算時間基本上是線性增加,而ESPF、GPF 和IWOPF 則是非線性增加。另外,IWOPF 的計算時間明顯少于ESPF 和GPF,這是由于IWOPF 沒有交叉和變異這2 種遺傳操作算子,算法相對簡單。綜上所述,IWOPF 繼承了IWO 算法簡單而有效的優點,相對于基于進化算法的粒子濾波,可以在大幅度縮短計算時間的情況下,有效提高濾波精度。
IWO 是近年來提出一種簡單、有效的基于種群的新穎數值優化算法。自提出以來,一直受到國內外學術界的關注,并在許多領域得到成功運用。本文將IWO 與PF 相結合,提出了一種基于IWO 算法的PF算法,該算法將IWO 運用于粒子濾波的重采樣環節,以保證粒子有效性和多樣性。計算機仿真結果表明,基于IWO 算法的PF 算法繼承了IWO 算法簡單而有效的優點,相對于基于進化算法的粒子濾波,可以在大幅度縮短計算時間的情況下,有效提高濾波精度。由于IWO 算法是近年來才提出一種智能優化算法,雖然其在許多應用研究領域取得了長足的發展,便遠未達到成熟的階段,其與粒子濾波相結合也有許多問題值得深入地分析與探討。例如:算法中的有一組參數需要預先設置,如何根據不同的情況設置合理的參數組合是下一步深入研究的問題。
[1]李遠征,盧朝陽,高全學,等.基于多特征融合的均值遷移粒子濾波跟蹤算法[J]. 電子與信息學報,2010,32(2):411 -415.LI Yuan-zheng,LU Zhao-yang,GAO Quan-xue,et al.Particle filter and mean shift tracking method based on multi-feature fusion [J]. Journal of Electronics &Information Technology,2010,32(2):411 -415.
[2]VILA J P.Enhanced consistency of the resampled convolution particle filter[J].Statistics and Probability Letters,2012,82:786-797.
[3]PARK C B,LEE S W.Real-time 3D pointing gesture recognition for mobile robots with cascade HMM and particle filter[J].Image and Vision Computing,2011(29):51-63.
[4]胡振濤,潘泉,梁彥,等.基于進化采樣的粒子濾波算法[J].控制理論與應用,2009,26(3):269 -273.HU Zhen-tao,PAN Quan,LIANG Yan,et al. The particle filter algorithm based on evolution sampling[J]. Control Theory & Applications,2009,26(3):269 -273.
[5]朱志宇.粒子濾波算法及其應用[M]. 北京:科學出版社,2010:73 -76.
[6]李翠蕓,姬紅兵.快速Metropolis-Hastings 變異的遺傳重采樣粒子濾波器[J].系統工程與電子技術,2009,31(8):1968-1972.LI Cui-yun,JI Hong-bing.Genetic resampling particle filter based on fast Metropolis-Hastings mutation[J]. Systems Engineering and Electronics,2009,31(8):1968 -1972.
[7]葉龍,王京玲,張勤. 遺傳重采樣粒子濾波器[J]. 自動化學報,2007,33(8):885 -887.YE Long,WANG Jing-ling,ZHANG Qin. Genetic resampling particle filter[J]. Acta Automatica Sinica,2007,33(8):885 -887.
[8]PARK S,HWANG J,ROU K,et al. A new particle filter inspired by biological evolution: genetic filter [C]//Proceedings of World Academy of Science,Engineering and Technology,Bangkok,Thailand:IEEE,2007(21):459-463.
[9]MEHRABIAN A R,LUCAS C.A novel numerical optimization algorithm inspired from weed colonization[J].Ecological Informatics,2006,1(4):355 -366.
[10]張帥,王營冠,夏凌楠.離散二進制入侵雜草算法[J].華中科技大學學報(自然科學版),2011,39(10):55-60.
[11]韓毅,蔡建湖,李延來,等.野草算法及其研究進展[J].計算機科學,2011,38(3):20 -23.HAN Yi,CAI Jian-hu,LI Yan-lai,et al. Invasive weed optimization and its advances[J].Computer Science,2011,39(10):55 -60.
[12]MOHAMADI M F,KOMJANI N,MOUSAVI P.Application of invasive weed optimization to design a broadband patch antenna with symmetric radiation pattern[J].IEEE Antennas and Wirwless Propagation Lwtters,2011(10):1369-1372.
[13]GOURAB G R,SWAGATAM D,PRITHWISH C,et al.Design of non-uniform circular antenna arrays using a modified invasive weed optimization algorithm[J]. IEEE Transactions on Antenns and Propagation,2011,59(1):110-118.
[14]MEHRABIAN A R,YOUSEFI-Koma A.A novel technique for optimal placement of piezoelectric actuators on smart structures[J].Journal of the Franklin Institute,2011(348):12-23.
[15]KUNDU D,SURESH K,GHOSH S,et al. Multi-objective optimization with artificial weed colonies[J]. Information Sciences,2011(181):2441 -2454.
[16]胡士強,敬忠良.粒子濾波原理及其應用[M].北京:科學出版社,2010:70 -71.