白曉波，邵景峰，田建剛

1.西安工程大學(xué) 管理學(xué)院，西安 710048

2.陸軍邊海防學(xué)院 工程基礎(chǔ)系，西安 710108

1 引言

粒子濾波（particle filter）[1]基于蒙特卡洛思想，對非線性非高斯系統(tǒng)進行狀態(tài)估計。但是粒子濾波存在權(quán)值退化問題，針對該問題，主要采用重采樣的方法解決[2]，但是重采樣時只保留權(quán)值較大粒子，而舍棄小權(quán)值粒子，這又導(dǎo)致粒子貧化問題。這就必定影響到系統(tǒng)狀態(tài)估計的精度。針對粒子濾波中粒子權(quán)值退化和粒子貧化問題，很多學(xué)者進行了大量研究，如文獻[3-5]提出了相應(yīng)的重采樣或改進算法。陳世明等人[6]提出一種基于引力場的粒子濾波算法以避免傳統(tǒng)粒子濾波算法中粒子貧化與退化現(xiàn)象。王向前等人[7]基于粒子濾波和高斯-鄰域搜索濾波的兩個實際特征，提出鄰域迭代重采樣粒子濾波。陳波等人[8]提出了改進的高斯粒子濾波，從而狀態(tài)變量服從線性變化，觀測方程為非線性變化的系統(tǒng)模型有顯著效果。這些改進方法都在很大程度上提高了粒子濾波的精度，但是，上述方法仍然基于重采樣方法，未能從根本上解決粒子貧化問題。王法勝等人[9]在粒子濾波的綜述中指出，樣本匱乏或粒子貧化問題依然是該研究領(lǐng)域的熱點。

為了徹底解決權(quán)值退化和粒子貧化問題，有學(xué)者在粒子濾波中融入了群體智能優(yōu)化算法。同時，文獻[10]指出基于群體智能算法優(yōu)化粒子濾波，是粒子濾波新的發(fā)展方向。核心思想是，利用群體尋優(yōu)思想替代重采樣。如文獻[11]根據(jù)蟻群算法的基本思想，將樣本粒子隨時間傳遞的過程等效為螞蟻行進路線，避免了重采樣方法引起的樣本匱乏問題。文獻[12]對基于傳統(tǒng)優(yōu)化算法以及機器學(xué)習(xí)方法的樣本匱乏解決技術(shù)進行了綜述。而比較有代表性的是田夢楚等人[13-14]和朱文超等人[15]利用螢火蟲算法對粒子濾波進行智能優(yōu)化，以及陳志敏等人[16]利用蝙蝠算法優(yōu)化粒子濾波。利用螢火蟲和蝙蝠算法的迭代尋優(yōu)，代替標準粒子濾波的重采樣，使粒子在迭代過程中合理地分布于高似然區(qū)和低似然區(qū)，保證了粒子的多樣性，且提高了粒子濾波的精度。

然而，螢火蟲算法和蝙蝠算法有其局限性，蝙蝠算法的尋優(yōu)能力主要依靠蝙蝠個體之間的相互作用和影響，但個體本身缺乏變異機制，一旦受到某個局部極值約束后自身很難擺脫；而且在進化過程中，種群中的超級蝙蝠可能會吸引其他個體迅速向其周圍聚集，使得種群多樣性大幅下降，這必然會在一定程度上影響到粒子多樣性。螢火蟲算法雖然考慮到了個體間的相互影響，但是個體向最優(yōu)值趨近的過程中，可能會出現(xiàn)螢火蟲移動距離大于個體與最優(yōu)值間距的情況，導(dǎo)致個體更新自己位置時跳過了最優(yōu)值，進而最優(yōu)值發(fā)現(xiàn)率降低，影響算法的收斂速度。

綜上，鑒于螢火蟲算法和蝙蝠算法的局限性，結(jié)合煙花算法（fireworks algorithm，F(xiàn)WA）[17-18]的優(yōu)勢，并對其進行改進，以更好地對粒子濾波進行優(yōu)化有著重要的理論和工程意義。

2 粒子分布多樣性度量

在文獻[13-14]和文獻[16]中，均體現(xiàn)了群體智能啟發(fā)式算法對粒子多樣性的促進作用。但是，各文獻并未對算法中粒子多樣性的程度進行度量。因此，這里提出在不同時刻粒子群多樣性的定義。

定義1設(shè)粒子濾波過程中，在k時有粒子集合為I={x1k,x2k,…,xNk}，粒子集合I的中心粒子為xck，各粒子xik到中心粒子xck的距離，則粒子多樣性指標表示如下。

以上粒子濾波中粒子多樣性度量的定義是后續(xù)研究工作的基礎(chǔ)。

3 改進煙花算法優(yōu)化粒子濾波

基于螢火蟲算法和蝙蝠算法的局限性，考慮粒子群個體之間的相互作用與影響、個體本身的變異機制，以及算法在迭代過程中出現(xiàn)的粒子個體跳出最優(yōu)值，出現(xiàn)震蕩，降低濾波精度的問題，先對煙花算法進行改進，再對粒子濾波進行優(yōu)化。

3.1 標準煙花算法

標準煙花算法[18]由譚營教授等人于2010年提出，算法過程如下。

步驟1隨機初始化N個煙花種群。

步驟2計算煙花Xi的爆炸半徑Ri，產(chǎn)生的火花數(shù)Si，計算公式如下。

步驟3產(chǎn)生爆炸火花。計算煙花Xi生成的第j（1≤j≤Si）個火花的第k維坐標公式如下。

步驟4變異算子，以增加種群多樣性。

步驟5映射規(guī)則。修正火花坐標，公式如下。

步驟6選擇進入下一次迭代的煙花種群。

3.2 標準FWA不能直接優(yōu)化PF的原因

標準的煙花算法不能直接用于優(yōu)化粒子濾波，主要體現(xiàn)在如下兩方面。

（1）火花變異操作。標準FWA的火花變異，采用均值和方差為1的高斯變異，如式（5）所示，雖然能夠在一定程度上促進爆炸火花的隨機性，若直接應(yīng)用于FWA，促進粒子多樣性的程度較低，對PF的濾波精度提高也極其有限。而PF針對的非線性、非高斯問題，需要粒子具有較好的多樣性，即在狀態(tài)值的高、低似然區(qū)均分布合理。那么，單一的高斯變異就不能最大程度地滿足實際系統(tǒng)。

（2）火花選擇策略。在文獻[18]中標準FWA的煙花選擇策略需要計算當前個體到候選集合其他個體之間的距離，若候選集合的待選火花數(shù)為N，有K個候選集合，那么，選擇策略的時間復(fù)雜度就為O(N2)。這樣的時間復(fù)雜度若直接用于優(yōu)化PF，必然在一定程度上降低PF實時性。

3.3 改進煙花算法

對煙花算法進行改進的具體過程如下。

步驟1隨機初始化N個煙花種群。

步驟2利用式（2）、式（3）對每個煙花分別計算Ri和Si。其中為常數(shù)，控制爆炸半徑，fmin=min(f(Xi))為當前火花集合中適應(yīng)度最小值；m為常數(shù)；fmax=max(f(Xi))為適應(yīng)度最大值，f(Xi)為煙花Xi的適應(yīng)度值，ξ表示計算機中的極小正常數(shù)。

為了防止Si過多或過少，對Si進行修正。具體方法如下：

rou(?)為取整函數(shù)，a、b為給定常數(shù)。

步驟3利用式（4）產(chǎn)生Xi的第j個、第k維坐標，同標準FWA的步驟3。

步驟4混合高斯變異火花。為了增加種群多樣性，通過變異算子產(chǎn)生新的火花，標準FWA中采用高斯函數(shù)生成高斯隨機數(shù)，但是在PF中，由于系統(tǒng)狀態(tài)非線性、噪聲數(shù)據(jù)非高斯，以致利用高斯函數(shù)產(chǎn)生的變異算子不能很好地適應(yīng)粒子濾波。因此，將標準FWA的高斯函數(shù)進一步改進為混合高斯變異算子，提高生成更優(yōu)良煙花的概率。

設(shè)煙花Xi按照預(yù)設(shè)參數(shù)產(chǎn)生M個火花，k=1,2,…,z表示火花xj的第k維坐標。則混合高斯分布下的坐標如下。

式中，f(T)為概率密度函數(shù)，c為階段分布系數(shù)，N為設(shè)定的高斯分布個數(shù)，βi、μi和σi分別為混合分布中第i個分布的比例系數(shù)、均值和標準差。T（Tmin≤T≤Tmax）為設(shè)定煙花隨機爆炸半徑，其分布函數(shù)表示如下。

那么，根據(jù)分布函數(shù)性質(zhì)得F(Tmax)=1，c(F(Tmax)-F(Tmin))=1，得截斷參數(shù)c=1/(F(Tmax)-F(Tmin))，轉(zhuǎn)換為標準正態(tài)分布函數(shù)，如式（11）所示。

Φ(?)為標準正態(tài)分布函數(shù)。

步驟5利用式（6）修正火花坐標。同標準FWA步驟5的映射規(guī)則。

步驟6選擇策略。標準煙花算法采用精英選擇策略將適應(yīng)度最好的煙花留作下一次迭代，具體計算方法詳見文獻[18]，其較復(fù)雜的計算過程不適用于粒子濾波，否則必然降低算法的實時性。因此，對選擇方法進行簡化。具體過程如下。

（1）第g代爆炸煙花為，煙花Xi,g的第g+1代火花集合為表示第k維坐標。則第g+1代的候選火花集合為C={xl,c|1≤l≤Si}?；谖墨I[19]中的CR公式，改進自適應(yīng)交叉因子CR的計算方法，對每個執(zhí)行如下操作。

式中，jrand為[1,2,…,Si]之間的隨機整數(shù)，T為迭代次數(shù)。CRi為交叉概率，是一條開口向上的拋物線。在種群向下變異的初始階段，CR值較小，且增速較慢，有利于煙花算法的全局搜索能力和多樣性；隨著迭代的進行，CR不斷增大，有利于煙花算法快速收斂，并提高局部搜索能力。進而，有利于PF的粒子在系統(tǒng)狀態(tài)的高、低似然區(qū)合理分布。

（2）對于候選火花集合C={xl,c|1≤l≤Si}，采用如下方式，用ui,g+1取代Xi,g向下迭代。

pα為隨機選取的火花向下傳遞的概率，rndF(?)表示從集合C中隨機選取一個xl,c?；诟倪M的選擇策略能夠在一定程度上將大權(quán)值粒子直接用于下一次迭代，也考慮了父代的遺傳信息傳遞給子代。且選擇策略的時間復(fù)雜度從標準FWA的O(N2)降為O(N)，從而不影響優(yōu)化后粒子濾波的實時性。

3.4 改進煙花算法的收斂性分析

從改進的煙花算法（improved fireworks algorithm，imp-FWA）整個過程分析，主要分為3個步驟：一是爆炸產(chǎn)生火花；二是煙花混合高斯變異；三是選擇下一代種群。因此，為了分析imp-FWA的收斂性，按這3個步驟依次進行。首先，引入以下定義及性質(zhì)。

定義2[20]對于方陣A∈Rn×n：

（1）若1≤ ?i,j≤n,aij≥ 0，則稱A為非負的(A≥ 0)。

（2）若1≤ ?i,j≤n,aij>0，則稱A為全正的(A>0)。

（3）若存在一個整數(shù)m≥1，使Am>0，則稱A為本原的。

（4）若A非負，且在相同的行、列初等變換后得到如下形式：

且C、T是方陣，則稱A是可約的。

（5）若A非負且，則稱A為隨機的。

（6）若A為隨機陣，且同一列中元素全相同，則稱A是穩(wěn)定的。

（7）若A為隨機陣，且每列中至少一個正數(shù)，則稱A是列容的。

引理1[21]如果矩陣Y是本原的隨機陣，那么Yk收斂到一個全正穩(wěn)定的隨機陣：

(y1,y2,…,yn)T唯一滿足：

為了研究imp-FWA的馬氏鏈的表示，設(shè)所有煙花的坐標空間為S，煙花位置為X，X∈S，馬氏鏈狀態(tài)空間為G，則狀態(tài)空間維數(shù)d=|G|=|S|N。

設(shè)imp-FWA每次迭代對應(yīng)的隨機變量為：R(t)={Xi(t)|1≤i≤n,1≤t＜∞}，表示煙花的位置集合，每次迭代都基于上一步的位置信息，但與原狀態(tài)沒有關(guān)系。因此，imp-FWA是馬爾可夫過程。為了簡單起見，根據(jù)imp-FWA主要的3個步驟，對imp-FWA的收斂性進行分步證明。

（1）爆炸產(chǎn)生的火花。將爆炸火花的爆炸算子的概率矩陣表示為F=(fij)L×L，fij表示爆炸后狀態(tài)由i變?yōu)閖。根據(jù)全概率公式得，F(xiàn)為隨機陣。

（2）煙花混合高斯變異。imp-FWA利用混合高斯變異，使得在每一次迭代時，每個煙花都具有一定概率向當前最優(yōu)位置移動。因此，引入?yún)?shù)ρ＞0，表示個體向當前最優(yōu)位置移動的概率。則一個群體的狀態(tài)由i到j(luò)的概率就定義為M=(mi,j)L×L。那么，如果兩個群體狀態(tài)i和j有相等的煙花數(shù)Ζij，則mij=，N為煙花種群數(shù)。從而，M為全正矩陣。

（3）選擇下一代種群。通過選擇算子，煙花狀態(tài)發(fā)生改變，設(shè)sij表示狀態(tài)i變?yōu)闋顟B(tài)j，則選擇運算對應(yīng)的概率矩陣為一隨機矩陣，表示為S=(sij)L×L，根據(jù)全概率公式得。

定理1通過混合高斯變異改進的imp-FWA能夠以概率1全局收斂到最優(yōu)解。

證明imp-FWA每一步進化過程的馬氏鏈轉(zhuǎn)移概率矩陣表示為P=FMS。

設(shè)U=FM，V=US。F是隨機矩陣，那么F各行中至少存在一個大于零的元素，因此：

從而U>0。而S亦為隨機陣，同理可得：

由定義2可知P為本原矩陣，再利用引理1可知最優(yōu)爆炸火花的位置不在馬氏鏈極限分布中的概率為0，且包括最優(yōu)位置的火花的極限分布之和為1。因此，定理1得證。

3.5 改進的FWA優(yōu)化PF

基于3.3節(jié)改進的煙花算法和標準PF（詳見文獻[9]），提出 FWA-PF（fireworks algorithm particle filter）。核心思想是利用imp-FWA對標準粒子濾波中的每一個粒子進行迭代爆炸尋優(yōu)，取代PF的重采樣，使粒子向目標狀態(tài)后驗概率值高的方向移動，且保持良好的多樣性，以提高狀態(tài)估計的準確性。但是，經(jīng)過煙花算法的處理后，各粒子在狀態(tài)空間的分布被改變，其分布密度函數(shù)也不能用p(xk|y1,k-1)表示，從而喪失了貝葉斯濾波理論基礎(chǔ)。因此，還需改進粒子權(quán)值的表示方法。根據(jù)文獻[22]的思想，改進后的粒子權(quán)值表示如下：

基于改進煙花算法，優(yōu)化后的粒子濾波算法過程如圖1所示。

FWA-PF算法的詳細實施過程如下。

步驟1初始化FWA參數(shù)，交叉因子CRmin、CRmax，煙花爆炸半徑控制參數(shù)，火花數(shù)控制參數(shù)m，迭代次數(shù)初值T=0，ξ為計算機中能夠表示的最小正值，爆炸火花數(shù)修正參數(shù)a、b，火花被選概率pα。

步驟2初始化粒子集合，并初始化粒子權(quán)值ωi，以權(quán)值代替式（2）、式（3）中的適應(yīng)度值。

步驟3對每個粒子Xi進行煙花爆炸變異。

步驟4利用式（15）修正粒子權(quán)值。

步驟5式（16）權(quán)值歸一化，計算方法如下：

步驟6系統(tǒng)狀態(tài)輸出：

3.6 FWA-PF時間復(fù)雜度

從FWA-PF的詳細過程來分析，影響算法效率的因素主要有兩個：一是粒子數(shù)N；二是每個粒子爆炸所產(chǎn)生的火花數(shù)Si。但是Si的取值具有隨機性，由式（3）計算、式（7）修正得出，若式（7）中的控制常數(shù)b無限接近1，則Si?m，那么FWA-PF在最壞情況下的時間復(fù)雜度為O(3m×N)。而在實際應(yīng)用時，式（3）中的fmax和f(Xi)為歸一化的粒子權(quán)值，則：

式（18）的值向0靠近，那么，若a取值較大，Si＜am，Si的取值大多為rou(am)，F(xiàn)WA-PF的時間復(fù)雜度為O(3rou(am)×N)，0＜a＜b＜1。從而在大多數(shù)情況下，F(xiàn)WA-PF具有較低的時間復(fù)雜度。

3.7 FWA-PF收斂性分析

為了分析FWA-PF的收斂性，首先對FWA-PF過程以及后續(xù)分析過程中用到的部分符號進行說明。

（1）狀態(tài)過程。FWA-PF依然遵循標準PF的基本特性，即狀態(tài)過程X為一階Markov過程，將初始分布和狀態(tài)轉(zhuǎn)移概率分布分別表述為p(dx0)和K(dxt|xt-1)，則動態(tài)系統(tǒng)的表述如下。

假設(shè)p(dyt|xt)概率密度存在，表示為p(yt|xt)[23]，則：

且p(yt|x0:t)=p(yt|xt)。

（3）FWA-PF新粒子的迭代生成，其實質(zhì)是利用爆炸算子、混合高斯變異和選擇算子共同作用的結(jié)果，也就是每次迭代的粒子按照如下分布進行的分步生成策略：Gh(dxt|xt-1,y1:t)為煙花爆炸混合高斯分布函數(shù)，假設(shè)存在其概率密度，則用Gh(xt|xt-1,y1:t)表示。

以下是分析過程中用到的符號及說明?！現(xiàn)‖表示有界函數(shù)的最大值，表示大于等于ε的最小整數(shù)。Xk:t稱作擴展狀態(tài)，表示k時到t時的狀態(tài)軌跡，Xk:t={Xk,Xk+1,…,Xt}，小寫字母的xk:t={xk,xk+1,…,xt}或zk:t={zk,zk+1,…,zt}表示其實現(xiàn)。Yk:t為k時到t時的量測軌跡，其實現(xiàn)為yk:t={yk,yk+1,…,yt}。π0:t|m(dx0:t)表示量測為y1:m時，x0:t的后驗概率分布p(dx0:t|y1:m)，π0:t|m(dx0:t)=p(X0:t∈dx0:t|Y1:m=y1:m)。(υ,τ)=∫τυ表示函數(shù)υ和τ內(nèi)積。

表示狀態(tài)過程X由0到k時的轉(zhuǎn)移概率分布，又因為X為一階Markov過程[24]，則：

δ為Dirac-Delta函數(shù)。

3.7.1 定義與假設(shè)

定義3R(t+1)nx上的函數(shù)h(x0:t)，若 (π0:t|t(dx0:t),|h(x0:t)|p)＜∞,p≥1成立，則h(x0:t)屬于Lt,p空間。由測度論[25]定義Lt,p空間的模如下：

假設(shè)1若量測軌跡y1:t已知，且t>0時(π0:t|t-1(dx0:t),p(yt|x0:t))>0成立。

假設(shè)2假設(shè)在FWA-PF中的任意參數(shù)μj,t（1≤j≤n,n為參數(shù)個數(shù)），在t>0時滿足條件。(π0:t|t-1(dx0:t),p(yt|x0:t))≥μj,t>0。實際應(yīng)用時，雖然難以提前找到符合條件的μi,t，但是若假設(shè)1成立，那么總能找到{μi,t|1≤i≤n,t>0}滿足假設(shè)2。

假設(shè)3t>0時，的模：

引理2假設(shè)1和假設(shè)3同時滿足，若h(x0:t)∈Lt.p，那么：

引理3若函數(shù)h(x0:t)有界，那么h(x0:t)∈Lt,p。

3.7.2 FWA-PF收斂性證明

基于3.7.1節(jié)的3個假設(shè)，得出如下結(jié)論。

引理4若假設(shè)1、假設(shè)2和假設(shè)3都滿足，那么對于?h(x0)∈Lt,4均有與N無關(guān)的一組數(shù)和，使以下兩式成立：

fw(?)表示獲得N個初始煙花樣本的分布函數(shù)，fw(dx0)=π*0N:0|0(dx0:0)。

引理5當假設(shè)1、假設(shè)2和假設(shè)3都滿足，若?H(x0:t-1)∈Lt-1,4均有與N無關(guān)的一組數(shù)b*t-1|t-1(H(x0:t-1))和η*t-1|t-1(H(x0:t-1))，使以下兩式成立：

那么，對于?h(x0:t)∈Lt,4就具有一組與N無關(guān)的數(shù)qt|t-1(h(x0:t))和ηt|t-1(h(x0:t))，使以下兩式成立：

引理6假設(shè)1、假設(shè)2和假設(shè)3滿足時，若?H(x0:t)∈Lt,4，均存在數(shù)qt|t-1(H(x0:t))和ηt|t-1(H(x0:t))（與N無關(guān)），使式（29）、式（30）成立，那么對于?h(x0:t)∈Lt,4均存在與N無關(guān)的一組數(shù)，使以下兩式成立：

引理7若假設(shè)1、假設(shè)2和假設(shè)3都滿足，且對于?H(x0:t)∈Lt,4均存在與N無關(guān)的一組數(shù)和，使式（31）和式（32）成立，那么就有數(shù)σt（σt與N無關(guān)），使得不成立的概率為，并且當。

引理8假設(shè)1、假設(shè)2和假設(shè)3都滿足，且對于均存在與N無關(guān)的一組數(shù)和，使式（31）和式（32）成立，且當時，對于均存在一組與N無關(guān)的數(shù)，使得以下兩式成立：

引理9假設(shè)1、假設(shè)2和假設(shè)3都滿足，若，均存在一組與N無關(guān)的數(shù)和使式（33）和式（34）成立，那么對于均存在與N無關(guān)的一組數(shù)和，使以下兩式成立：

引理10假設(shè)1、假設(shè)2和假設(shè)3都滿足，若均存在與N無關(guān)的一組數(shù)和，使式（35）、式（36）成立，則對于Lt,4，均存在一組與N無關(guān)的數(shù)，使下式成立。

假設(shè)4對于FWA-PF，假設(shè)其迭代次數(shù)有限。

基于引理4至引理10，可得：若假設(shè)1、假設(shè)2、假設(shè)3和假設(shè)4都滿足，則對于成立。表示在經(jīng)驗分布FWA-PF的估計值；表示其最優(yōu)估計值。下，h(x0:t)的均值作為

3.8 FWA-PF運行機制

FWA-PF的基本思想是利用FWA的煙花爆炸產(chǎn)生的火花，進行混合高斯變異，再通過選擇策略選擇新的火花用作下一次迭代，以取代標準PF的粒子重采樣過程，以實現(xiàn)粒子的多樣性。其運行機制如圖2所示。

在文獻[16]中描述了標準PF在經(jīng)過多次重采樣后的粒子分布情況，即很多粒子都集中于高似然區(qū)，這就喪失了粒子的多樣性。因此，在圖2中，重點表述了FWA-PF中粒子的初始分布、優(yōu)化后的運動趨勢和優(yōu)化后的分布情況。以下是對圖2的進一步說明。

（1）圖2（a）是FWA-PF的粒子初始分布，粒子在高低似然區(qū)均有分布。

（2）圖2（b）是FWA-PF的各粒子經(jīng)爆炸、混合高斯變異，以及選擇以后，除兩個方形藍色的粒子在一定概率下直接用作下一代，位置不變外，其他粒子均向各自的方向移動。爆炸的隨機性保證了粒子的多樣性，又由于使用了式（6）的火花坐標修正方法，因此各粒子的分布并不會超過火花爆炸邊界。

（3）圖2（c）是FWA-PF各粒子在圖2（b）所示的趨勢運動后的位置，由于其他粒子的運動，原來紅色的最優(yōu)粒子離真實值較遠，黑色實心粒子成為最優(yōu)粒子，而方形藍色粒子位置保持不變。黑色菱形粒子為上一次迭代時最優(yōu)粒子的新位置。

4 實驗仿真

基于改進煙花算法imp-FWA，在第3章設(shè)計了FWA-PF算法，重點從以下四方面說明算法的可行性。

（1）改進煙花算法imp-FWA的BenchMark測試。主要通過進化算法的6個BenchMark標準函數(shù)（Sphere、Quadric、Ackley、Rosenbrock、Griewangk、Rastrigin）測試imp-FWA的收斂性和性能。

（2）FWA-PF粒子多樣性評價。FWA-PF基于改進的煙花算法，其性能受初始參數(shù)的影響，煙花爆炸半徑控制參數(shù)，火花數(shù)控制參數(shù)m，火花數(shù)修正參數(shù)a、b。若變異的火花數(shù)較少，可保留較多父代的信息，但降低了種群的多樣性，以致尋優(yōu)效果不好；若變異的火花數(shù)較多，可增加種群多樣性，但最優(yōu)煙花對變異個體產(chǎn)生的貢獻少，導(dǎo)致算法的收斂速度緩慢，以致實時性較低。因此，根據(jù)定義1，先對兩個重要參數(shù)對粒子多樣性的影響進行分析。

Fig.2 Particle mechanism optimization of FWA-PF圖2 FWA-PF粒子優(yōu)化機制

M=50為實驗仿真次數(shù)；k為濾波次數(shù)；j為狀態(tài)估計值；xj為量測值。

（4）FWA-PF與ADE-PF[3]、FA-PF[14]和BA-PF[16]做性能和效率的對比分析。

實驗環(huán)境、狀態(tài)方程和量測方程。

（1）實驗仿真環(huán)境。CPU：Intel?CoreTMi5-4200U@160 GHz 2.30 GHz；內(nèi)存：4 GB，Windows 7 64位操作系統(tǒng)，500 GB硬盤，Matlab R2012a。

（2）狀態(tài)方程和量測方程。

Q=10為系統(tǒng)過程噪聲方差，R=10為量測噪聲方差，randn為(0,1)之間的隨機數(shù)，粒子數(shù)N=200，修正參數(shù)a=0.2，b=0.8，即煙花爆炸產(chǎn)生的火花數(shù)在[20,80]之間，式（9）中的參數(shù)N=3。

4.1 imp-FWA的BenchMark測試

對于進化搜索算法，通常利用BenchMark標準函數(shù)測試其性能。這里主要利用6個常用BenchMark函數(shù)測試imp-FWA的性能，并與標準FWA進行對比。用于測試的6個BenchMark函數(shù)如表1所示。

Table 1 6 BenchMark test functions表1 6個BenchMark測試函數(shù)

由圖3的6個BenchMark測試函數(shù)的對比得出以下結(jié)論：imp-FWA的收斂速度較FWA慢，但是由圖3（b）、圖3（c）和圖3（f）可知，在Ackley、Quadric和Rastrigin函數(shù)上，即使到搜索后期，imp-FWA要優(yōu)于FWA，這得益于改進的選擇策略。而其他圖上的結(jié)果在搜索后期FWA和imp-FWA具有相似的收斂效果。尤其在圖3（d）中，由于Rosenbrock函數(shù)的信息較少，其搜索方向不明，故imp-FWA和FWA都具有相對穩(wěn)定的階段。綜上，imp-FWA可以取得較好的收斂性，且由于選擇策略的簡化，使其整個搜索效率提高，能夠用來優(yōu)化粒子濾波。

4.2 簡化的火花選擇策略對PF精度的影響

Fig.3 Comparison of convergence of imp-FWAand FWA圖3 imp-FWA與FWA收斂性對比

標準FWA的選擇策略是將適應(yīng)度最好的火花傳遞給下一代，其目的在于使得迭代過程中火花種群能夠快速地向最優(yōu)值靠近，而改進煙花算法重點目的在于解決PF中粒子權(quán)值退化和粒子在迭代過程中喪失多樣性的問題。因此，在改進火花選擇策略時，重點是為了保證PF粒子的多樣性，火花的選擇具有一定概率下的隨機性。同時，在一定概率下將父代的優(yōu)良火花可以直接傳遞給下一代。如式（12）～式（14）所示。與標準FWA相比，imp-FWA只是在一定概率下將適應(yīng)度最好的火花傳遞給下一代，導(dǎo)致了一部分適應(yīng)度好的火花的丟失，在尋優(yōu)過程的前期收斂性低于標準FWA，但全局收斂性略好于FWA，符合4.1節(jié)BenchMark測試結(jié)果。若優(yōu)化PF，粒子多樣性的增加能夠提高濾波精度。另外，F(xiàn)WA-PF通過混合高斯變異算子來增加火花的多樣性。

以下為實驗仿真的方式，對比imp-FWA優(yōu)化PF（記為PF-1）和標準FWA優(yōu)化PF（記為PF-2）的濾波效果。如圖4所示。

在圖4中不難發(fā)現(xiàn)PF-1在離散點上更接近真實值。imp-FWA和標準FWA優(yōu)化PF，在不同粒子數(shù)和相同方差下（過程方差和量測方差Q=R=20）進行50次獨立測試的濾波性能運行時間對比。如表2所示。

Fig.4 Filter effect comparison of PF-1 and PF-2圖4 PF-1與PF-2濾波效果對比

從表2的結(jié)果分析，在相同的粒子數(shù)、過程方差和量測方差下，imp-FWA優(yōu)化后的PF的精度略高于標準FWA優(yōu)化的PF，且運行效率也高。綜合4.1節(jié)和4.2節(jié)的實驗結(jié)果，得出以下結(jié)論?；鸹ㄟx擇策略的簡化，雖然在一定程度上降低了FWA的收斂速度，但是應(yīng)用于粒子濾波的優(yōu)化時，卻能保證粒子的多樣性以及父代優(yōu)良粒子在一定概率下直接傳給下一代。相比將標準FWA直接用來優(yōu)化PF，imp-FWA優(yōu)化PF更能提高PF的性能和運行效率。

Table 2 Performance comparison of optimized PF by imp-FWAand by FWA表2 imp-FWA和FWA優(yōu)化PF性能對比

4.3 FWA-PF粒子多樣性指標

Fig.5 Particle distribution(k=100)圖5 粒子分布(k=100)

在圖5中，k時的真值為15.858 9。從粒子的分布來看，PF的粒子大多集中于高似然區(qū)，而FWA-PF的分布在高低似然區(qū)，從而FWA-PF比PF的粒子多樣性更好。由式（1）得100個粒子多樣性指標ρ=25.145 3。FWA-PF的粒子多樣性指標ρ=25.400 6。從而也說明了粒子多樣性對濾波精度的促進作用。理論上，煙花爆炸半徑的增大，增加了爆炸產(chǎn)生火花的搜尋范圍，爆炸產(chǎn)生的火花數(shù)的增加，有利于增加粒子的多樣性。同時，鑒于實驗的隨機性，對主要兩個參數(shù)和m在不同取值時，重復(fù)進行50次實驗，計算粒子多樣性指標ρ的均值，分析煙花半徑和爆炸火花數(shù)對粒子多樣性指標的影響。具體實驗結(jié)果如圖6所示。

Fig.6 Relation of parameter,mand particle diversity index圖6 參數(shù)、m與粒子多樣性指標關(guān)系

4.4 主要參數(shù)取值對FWA-PF性能影響

從理論上分析，粒子多樣性指標越大，粒子的分布就越為合理，且能提高濾波精度，在4.1節(jié)中，F(xiàn)WA-PF的估計值為15.553 0，標準PF的估計值為15.470 8，真值為15.858 9。那么，F(xiàn)WA-PF的濾波值更接近真值，濾波精度更高。同樣基于實驗存在的隨機性，利用式（20）計算各參數(shù)下FWA-PF的RMSE均值。

Fig.7 Filtering comparison of PF and FWA-PF(m=110,=5)圖7 PF與FWA-PF濾波對比(=5,m=110)

根據(jù)圖7的對比結(jié)果分析，F(xiàn)WA-PF比標準PF的估計值更接近真值，但是其效果并不明顯。FWA-PF的均方根誤差RMSE=1.826 3，標準PF的均方根誤差RMSE=2.177 7。從而說明煙花爆炸半徑參數(shù)不變，增加火花爆炸數(shù)，能夠提高粒子濾波的濾波精度，這一結(jié)論與圖4所示的實驗結(jié)論相吻合，即火花爆炸數(shù)的增大能夠增加粒子多樣性指標，而粒子多樣性指標大，則濾波精度也高。

Fig.8 Filtering comparison of FWA-PF and PF(m=120,=5)圖8 FWA-PF與PF濾波對比(=5,m=120)

圖8中，通過對火花數(shù)的增加，F(xiàn)WA-PF的估計值比m=110時更接近真值，原因與圖7的實驗相同。標準PF的RMSE為2.245 8，F(xiàn)WA-PF的RMSE為1.629 8。

由于實驗仿真過程中存在隨機性，針對每個火花數(shù)控制參數(shù)m進行50次獨立實驗，并求得RMSE均值，分析火花數(shù)控制參數(shù)m對FWA-PF濾波精度的影響，如表3所示。

Table 3 Effect ofmon precision of FWA-PF filtering(=5)表3 m對FWA-PF濾波精度影響(=5)

Table 3 Effect ofmon precision of FWA-PF filtering(=5)表3 m對FWA-PF濾波精度影響(=5)

序號1 2 3 4 5火花爆炸數(shù)m100 110 120 130 140 RMSE均值1.928 8 1.877 9 1.750 4 1.725 6 1.706 1

在表3中，在煙花爆炸半徑不變，增加火花爆炸數(shù)能夠提高FWA-PF的濾波精度，這和圖6所示的結(jié)果一致，即火花爆炸數(shù)的增加能夠提高粒子多樣性，而粒子多樣性指標的增加又促進了濾波精度的提高。

其他參數(shù)不變，分析煙花爆炸半徑R的增加對FWA-PF精度的影響。

（3）參數(shù)m=100，=10。實驗結(jié)果如圖9所示。

Fig.9 Filtering comparison of FWA-PF and PF(m=100,=10)圖9 FWA-PF與PF濾波效果對比圖(m=100,=10)

（4）參數(shù)m=100，=15。實驗結(jié)果如圖10所示。

Fig.10 Filtering comparison of FWA-PF and PF(m=100,=15)圖10 FWA-PF與PF濾波效果對比圖(m=100,=15)

對于FWA-PF，在圖9和圖10的仿真效果上，并未能很好地體現(xiàn)爆炸半徑參數(shù)對濾波精度的影響，但是從均方根誤差RMSE來分析，RMSE(R=10)=1.867 0，RMSE(R=15)=1.852 9。即增大煙花爆炸半徑時，F(xiàn)WA-PF的濾波精度有所提高，均高于標準粒子濾波的濾波精度，在圖9和圖10中，RMSE(PF)=1.946 0。鑒于實驗仿真時結(jié)果的隨機性，針對不同的爆炸半徑各進行50次仿真實驗，得出實驗結(jié)果如表4所示。

在表4中，隨著火花爆炸半徑的增加，RMSE逐漸減小，說明濾波精度增加，即火花爆炸半徑的增加能夠提高FWA-PF濾波精度。

Table 4 Effect ofon precision of FWA-PF filtering(m=100)表4 對FWA-PF濾波精度影響(m=100)

Table 4 Effect ofon precision of FWA-PF filtering(m=100)表4 對FWA-PF濾波精度影響(m=100)

序號1 2 3 4 5火花爆炸半徑 10 15 20 25 30 RMSE均值1.920 3 1.896 8 1.775 2 1.742 8 1.726 7

4.5 FWA-PF與其他算法的對比

為了進一步說明FWA-PF的有效性，利用相同的真實數(shù)據(jù)，將 FWA-PF（m=120,=20）與 ADE-PF（adaptive differential evolution particle filtering）[3]、FAPF（firefly algorithm optimized particle filter）[14]和 BAPF（bat algorithm optimized particle filter）[16]進行以下對比。

（1）粒子數(shù)N=100,K=100,Q=10,R=10。以上4種濾波算法濾波效果對比如圖11所示。

Fig.11 4 algorithms of filtering effect contrast diagram圖11 4種算法濾波效果對比圖

在圖11實驗結(jié)果中，4種濾波算法的RMSE均值分別為：RMSE(FWA-PF)=1.650 7，RMSE(BA-PF)=2.180 1，RMSE(FA-PF)=2.190 1，RMSE(ADE-PF)=2.219 0。從RMSE均值可得出以下結(jié)論：FWA-PF的濾波精度優(yōu)于其他3種濾波算法。

（2）不同粒子數(shù)算法性能與效率對比。不同粒子數(shù)時各進行50次實驗，對比其4種濾波算法的性能和效率。K=200，實驗結(jié)果如表5所示。

在表5的實驗結(jié)果中，F(xiàn)WA-PF在相同參數(shù)條件下，濾波精度均高于其他3種算法。在算法運行效率上，略低于BA-PF，與FA-PF基本相當，但優(yōu)于ADEPF。通過以上實驗，充分說明了改進的FWA優(yōu)化粒子濾波的可行性。

4.6 FWA-PF的實施方法

根據(jù)4.3節(jié)和4.4節(jié)的實驗結(jié)論，F(xiàn)WA-PF中的粒子多樣性和濾波精度均受火花爆炸數(shù)和爆炸半徑的影響，即粒子多樣性和濾波精度與煙花爆炸半徑和火花爆炸數(shù)都成正相關(guān)關(guān)系，但是，火花爆炸數(shù)的增加必然引起算法時間復(fù)雜度增加，并降低FWA-PF的實時性。那么，在解決實際問題時，取得合適的參數(shù)就顯得尤為重要。而對于FWA-PF，對其性能的要求，主要體現(xiàn)在兩方面，一是濾波精度，二是運行效率。因此，下面結(jié)合這兩個因素給出FWA-PF中爆炸半徑R和火花爆炸數(shù)S的設(shè)置過程如下。

步驟1在特定軟硬件環(huán)境下（操作系統(tǒng)、處理器性能、內(nèi)存和硬盤容量），以及基本的粒子數(shù)、狀態(tài)過程方差和量測過程方差，設(shè)定合理的濾波精度閾值p（常以RMSE均值表示）和運行時間閾值e。

步驟2設(shè)置初始火花爆炸半徑R和爆炸數(shù)S。

步驟3運行程序，經(jīng)過多次采樣（一般不少于50次），計算均方根誤差均值RMSEavg和運行時間均值tavg。在具體調(diào)整參數(shù)時又分為以下4種情況。

Table 5 Results of experimental simulation表5 實驗仿真結(jié)果

（1）若RMSEavg≤p，且tavg≤e，說明系統(tǒng)達到了濾波精度和實時性要求，不需要再進行調(diào)整。進入步驟4，結(jié)束調(diào)整。

（2）若RMSEavg>p，但tavg≤e，說明實時性滿足要求，而精度不足。那么就優(yōu)先考慮以步長l增加煙花爆炸半徑，其次以步長s增加爆炸粒子數(shù)，兩者同時增加，重復(fù)步驟3。

（3）若RMSEavg≤p，但tavg>e，說明濾波精度達到要求，而運行效率不足。對于這種情況，優(yōu)先考慮以步長s減少爆炸粒子數(shù)，重復(fù)步驟3。

（4）若RMSEavg>p，但tavg>e，說明濾波精度和運行效率均未達到要求，那么就優(yōu)先考慮以步長l增加煙花爆炸半徑，并以步長s減少火花爆炸數(shù)，重復(fù)步驟3。

步驟4結(jié)束調(diào)整，設(shè)定參數(shù)火花爆炸半徑R和火花爆炸數(shù)S。

通過以上4個步驟，實現(xiàn)FWA-PF兩個重要參數(shù)的設(shè)定。

5 結(jié)束語

針對粒子濾波中的粒子退化和喪失多樣性問題，將煙花算法（FWA）的火花變異算法改為混合變異算法，以此提高FWA在粒子濾波中對不同分布的適應(yīng)性，從而保證了粒子的多樣性。通過實驗仿真的方法，從濾波的精度和程序運行效率兩方面，驗證了利用改進的FWA優(yōu)化標準粒子濾波的可行性。然而，F(xiàn)WA算法的性能受制于兩個主要參數(shù)：煙花爆炸半徑和爆炸火花數(shù)。這兩個參數(shù)影響到粒子的多樣性和濾波性能。因此，針對實際應(yīng)用時的狀態(tài)方程和量測方程，如何自適應(yīng)地選擇最合適的FWA參數(shù)是下一步主要的研究內(nèi)容。