許江湖，黃 亮，劉 忠

(海軍工程大學(xué) 電子工程學(xué)院，湖北 武漢430033)

0 引 言

粒子濾波又稱序貫蒙特卡羅方法，是一種基于蒙特卡羅方法和遞推貝葉斯估計的統(tǒng)計濾波方法。粒子濾波完全突破了Kalman 濾波理論框架，對系統(tǒng)的過程噪聲和量測噪聲沒有任何限制。粒子濾波(PF)作為解決非線性、非高斯動態(tài)系統(tǒng)的參數(shù)估計和狀態(tài)濾波問題的有效方法，近年來隨著計算機運算能力的急劇增長和計算成本的不斷降低，已成為研究熱點［1－3］。

粒子退化是標準粒子濾波算法的主要缺陷。采用重采樣方法在一定程度上可以抑制粒子退化現(xiàn)象的發(fā)生。然而重采樣帶來的新問題是:權(quán)值越大的粒子的子代越來越多，而權(quán)值較小的粒子被剔除，最糟糕的情形是新的粒子集實際都是一個權(quán)值最大的粒子的子代，即所謂“樣本枯竭”現(xiàn)象，從而導(dǎo)致粒子集的多樣性變差，不足以用來近似表征后驗密度，難以保證估計精度，特別是在樣本受限條件下，這種粒子多樣性減弱對濾波精度的影響更為突出，甚至導(dǎo)致濾波發(fā)散現(xiàn)象［4］。為此，不少學(xué)者嘗試利用進化算法(包括遺傳算法、進化策略、進化規(guī)劃)來解決粒子濾波的退化問題［4－8］，并取得了一定的理論成果。2006 年，Mehrabian 和Lucas 提出了一種新穎的數(shù)值優(yōu)化算法，稱為入侵野草(雜草)優(yōu)化(IWO)算法［9］。同進化算法一樣，IWO算法也是一種隨機搜索仿生學(xué)優(yōu)化算法，該算法模仿了野草入侵的種子空間擴散、生長、繁殖和競爭性消亡的基本過程，具有很強的魯棒性和自適應(yīng)性。與進化算法相比，IWO 算法簡單易于實現(xiàn)，不需要遺傳操作算子，能簡單而有效地收斂于問題的最優(yōu)解，是一種強有力的智能優(yōu)化算法。研究表明，IWO 算法在性能上優(yōu)于進化算法和人工蜂群克隆算法，并且在高維問題上優(yōu)于粒子群算法［10］。因此，IWO 算法自提出以來，一直受到學(xué)者們的關(guān)注和運用［10－15］。

本文將IWO 算法與粒子濾波算法相結(jié)合，提出基于IWO 算法的粒子濾波算法(IWOPF)。該算法將IWO 應(yīng)用于粒子重采樣，以保證粒子的有效性和多樣性。最后通過計算機仿真比較該算法與基于進化算法的粒子濾波算法以及傳統(tǒng)粒子濾波算法的濾波性能和計算時。

1 粒子濾波簡介

離散時間非線性動態(tài)系統(tǒng)的狀態(tài)方程和量測方程分別為:

式中:k 為采樣時刻;xk∈Rnx為k 時刻的系統(tǒng)狀態(tài)向量;fk:Rnx×Rnv→Rnx為系統(tǒng)狀態(tài)xk-1的非線性函數(shù);vk為過程噪聲，其協(xié)方差為Qk，均值為0;nx和nv分別為狀態(tài)和過程噪聲的維數(shù);zk∈Rnz為k 時刻系統(tǒng)狀態(tài)的量測向量;hk:Rnx×Rnn→Rnz為系統(tǒng)狀態(tài)xk的非線性函數(shù);nk為量測噪聲，其協(xié)方差為Rk，均值為0，且與vk相互獨立;nz和nn分別為量測值和量測噪聲的維數(shù)。

濾波過程的任務(wù)是通過可獲得的系統(tǒng)觀測值zk估計出系統(tǒng)狀態(tài)xk，也即要求得到系統(tǒng)狀態(tài)的后驗概率分布p(xk| z1:k)。粒子濾波的核心思想是利用一系列隨機樣本的加權(quán)和表示狀態(tài)后驗概率密度:

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

標準粒子濾波算法選擇最易于實現(xiàn)的先驗概率密度作為重要密度函數(shù)，即

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

當(dāng)前時刻粒子的權(quán)重被評估后，通過引入重采樣技術(shù)來改善粒子退化問題。

2 IWO 算法

IWO 算法中，“野草”表示隨機產(chǎn)生的可行解，“種子”表示野草的后代，“種群”表示所有野草的集合。一般而言，IWO 算法的執(zhí)行過程要經(jīng)歷以下4 個步驟:

1)初始化種群

在這個步驟中，需要確定最大種群數(shù)目Pmax、最大迭代次數(shù)itermax、問題維數(shù)D、最大和最小可生成種子數(shù)Smax和Smin、非線性調(diào)節(jié)指數(shù)n、標準差初始值σinit和最終值σfinal以及初始搜索空間X，并隨機生成N0個解。

2)繁殖

各野草產(chǎn)生的種子個數(shù)為:

式中:f 為當(dāng)前野草的適應(yīng)度值;fmax和fmin分別為當(dāng)前種群中野草對應(yīng)的最大和最小適應(yīng)度值。

從式(7)可以看出，種群中的成員能夠散播的種子數(shù)是根據(jù)該成員的適應(yīng)度值以及種群所有成員的最大和最小適應(yīng)度值來決定的。在用優(yōu)化算法求解問題過程中，可能會直覺地認為可行解比不可行解具有更好的適應(yīng)度值，因此會放棄不可行解。然而，這種做法忽視了進化算法具有概率性，其產(chǎn)生的解是在不可行和可行之間不斷轉(zhuǎn)化的。Mehrabian和Lucas 認為，在進化過程中一些不可行解可能帶有更多的有用信息，如果它可能穿過一個不可行區(qū)域(特別是非凸的可行搜索空間)，則系統(tǒng)更容易到達最優(yōu)點。因此，在IWO 算法中，繁殖過程按照自然界中的繁殖法則，給予不可行的個體生存和繁殖的機會，只是這種機會相對較小。

3)空間散布

空間散布體現(xiàn)了算法的隨機性和適應(yīng)性。種群產(chǎn)生的種子個體按正態(tài)分布N(0，σ)在其父代野草個體附近的D 維空間進行散布。每次迭代對應(yīng)不同的標準差，隨著迭代的進行，標準差從初始值σinit減小到最終值σfinal，當(dāng)前迭代的標準差為:

式(8)確保了在較遠區(qū)域產(chǎn)生種子的概率以非線性的方式逐漸降低，以使適應(yīng)度值好的個體聚集在一起。

4)競爭性生存

經(jīng)過數(shù)代繁殖后，產(chǎn)生的后代數(shù)目將超過環(huán)境資源的承受能力，通過預(yù)先設(shè)定的最大種群數(shù)目Pmax對種群數(shù)量進行控制。在算法迭代過程中，種群中的所有野草及其后代按適應(yīng)度值從大到小進行排序，只有適應(yīng)度值最好的前Pmax個個體可以生存下來，其余的個體被清除。

綜合上述步驟可知:在算法運行過程中，一個野草個體可以有多個種子后代，種子按正態(tài)分布在父代野草周圍。只有相對較優(yōu)秀的種子才有機會進入下一輪迭代，一旦種子存活下來進入一輪迭代，所有存活下來的野草和種子共同構(gòu)成一輪迭代的種群。這種機制給予那些適應(yīng)度值低的個體繁殖的機會，如果它們的后代適應(yīng)度值更好，這些后代就可以生存下來。當(dāng)算法運行到預(yù)先設(shè)定的最大迭代次數(shù)itermax時，算法終止。

3 IWOPF 算法

前面已經(jīng)述及，粒子退化是標準粒子濾波算法的主要缺陷。而粒子濾波與IWO 算法的相似之處在于它們都有一個初始化的群體，每個個體代表系統(tǒng)的一個可能解，這些個體都根據(jù)一定的規(guī)則進行狀態(tài)轉(zhuǎn)變，它們都對適應(yīng)度較高的個體進行復(fù)制。因此，用IWO 算法來解決粒子濾波的退化問題具有可行性。

IWOPF 算法將IWO 算法應(yīng)用于粒子重采樣以解粒子濾波的退化問題。即通過IWO 算法中的繁殖、空間散布、競爭性生存等操作來保證粒子的多樣性和有效性。具體而言，IWOPF 算法的執(zhí)行步驟如下:

1)初始化

設(shè)置粒子數(shù)目N，從先驗分布p(x0)中采集粒子令確定最大迭代次數(shù)itermax、問題維數(shù)D、最大和最小可生成種子數(shù)Smax和Smin、非線性調(diào)節(jié)指數(shù)n、標準差初始值σinit和最終值σfinal。而最大種群數(shù)目Pmax= N。

for k = 1，2，…。

2)重要性采樣

3)更新權(quán)值

根據(jù)當(dāng)前觀測值zk，計算每個粒子的權(quán)值，并歸一化:

4)基于IWO 算法的重采樣

②繁殖:根據(jù)式(7)計算每個粒子產(chǎn)生的種子數(shù);

③空間散布:根據(jù)式(8)計算當(dāng)前迭代的標準差;

④競爭性生存:把所有粒子及其后代的適應(yīng)度值從大到小進行排序，只有適應(yīng)度值最好的前Pmax= N 個粒子可以生存下來，并重新定義這些粒子為| i = 1，2，…，N};

⑤判斷算法運行是否達到預(yù)先設(shè)定的最大迭代次數(shù)itermax，如果是，則轉(zhuǎn)步驟5，并令否則轉(zhuǎn)步驟4 中的①;

5)狀態(tài)估計

4 仿真結(jié)果與分析

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

運動模型為:

其中:

采樣周期T = 1 s，仿真時間為80 s，轉(zhuǎn)彎角速度為目標初始位置為(7，8)km，初始速度(－0.08， －0.06)km/s，系統(tǒng)過程噪聲是均值為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有關(guān)參數(shù)如表1 所示。圖1 和圖2 分別為粒子數(shù)N =400 時，蒙特卡羅仿真次數(shù)為50 次，實例1 基于X 軸方向和Y 軸方向目標位置的均方根誤差(RMSE)曲線。表2 和表3 分別給出了不同粒子數(shù)目情況下，4 種濾波算法對實例1 的濾波平均RMSE 和50次仿真所用的時間。其中，平均RMSE 的定義如下［16］:平均

式中:MC 為仿真次數(shù);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 參數(shù)表Tab.1 IWOPF parameter values

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

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

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

表5 不同粒子數(shù)目情況下仿真實例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 算法明顯優(yōu)于另外3 種粒子濾波算法。而從表3 和表5 可看出，隨著粒子數(shù)目的增加，標準PF計算時間基本上是線性增加，而ESPF、GPF 和IWOPF 則是非線性增加。另外，IWOPF 的計算時間明顯少于ESPF 和GPF，這是由于IWOPF 沒有交叉和變異這2 種遺傳操作算子，算法相對簡單。綜上所述，IWOPF 繼承了IWO 算法簡單而有效的優(yōu)點，相對于基于進化算法的粒子濾波，可以在大幅度縮短計算時間的情況下，有效提高濾波精度。

5 結(jié) 語

IWO 是近年來提出一種簡單、有效的基于種群的新穎數(shù)值優(yōu)化算法。自提出以來，一直受到國內(nèi)外學(xué)術(shù)界的關(guān)注，并在許多領(lǐng)域得到成功運用。本文將IWO 與PF 相結(jié)合，提出了一種基于IWO 算法的PF算法，該算法將IWO 運用于粒子濾波的重采樣環(huán)節(jié)，以保證粒子有效性和多樣性。計算機仿真結(jié)果表明，基于IWO 算法的PF 算法繼承了IWO 算法簡單而有效的優(yōu)點，相對于基于進化算法的粒子濾波，可以在大幅度縮短計算時間的情況下，有效提高濾波精度。由于IWO 算法是近年來才提出一種智能優(yōu)化算法，雖然其在許多應(yīng)用研究領(lǐng)域取得了長足的發(fā)展，便遠未達到成熟的階段，其與粒子濾波相結(jié)合也有許多問題值得深入地分析與探討。例如:算法中的有一組參數(shù)需要預(yù)先設(shè)置，如何根據(jù)不同的情況設(shè)置合理的參數(shù)組合是下一步深入研究的問題。

［1］李遠征，盧朝陽，高全學(xué)，等．基于多特征融合的均值遷移粒子濾波跟蹤算法［J］． 電子與信息學(xué)報，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］．控制理論與應(yīng)用，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］朱志宇．粒子濾波算法及其應(yīng)用［M］． 北京:科學(xué)出版社，2010:73 －76．

［6］李翠蕓，姬紅兵．快速Metropolis－Hastings 變異的遺傳重采樣粒子濾波器［J］．系統(tǒng)工程與電子技術(shù)，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］． 自動化學(xué)報，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］．華中科技大學(xué)學(xué)報(自然科學(xué)版)，2011，39(10):55－60．

［11］韓毅，蔡建湖，李延來，等．野草算法及其研究進展［J］．計算機科學(xué)，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］胡士強，敬忠良．粒子濾波原理及其應(yīng)用［M］．北京:科學(xué)出版社，2010:70 －71．