磷蝦群優化的改進粒子濾波算法

朱震曙，蔣長輝，薄煜明，吳盤龍

(南京理工大學 自動化學院，南京 210094)

在實際工程應用中，非線性系統無處不在.在這種應用背景下，基于卡爾曼濾波理論的濾波算法性能大大降低.

粒子濾波是一種基于蒙特卡羅方法和貝葉斯估計的統計濾波算法，利用大量帶權值的隨機樣本來逼近狀態.相對于傳統的狀態估計方法，粒子濾波沒有狀態函數、觀測函數非線性和非高斯的限制，因此在非線性和非高斯的系統估計中有著廣泛應用

標準的粒子濾波算法中，會出現權值退化現象[1-3]，表現為經過多次遞推以后，只有少量粒子的權值較大，其余粒子的權值幾乎為零.重采樣可以解決權值退化問題，但也會帶來樣本貧化現象，即高權值粒子被過度復制，有效粒子數減少，導致粒子的整體信息容量降低.針對上述問題，國內外許多學者進行了大量的研究.文獻[4]提出一種類電磁機制優化的粒子濾波，將采樣粒子看成具有相互吸引和排斥的電荷粒子，通過電磁吸引力引導粒子運動，但是本質上仍然無法完全消除由重采樣造成的樣本貧化問題.文獻[5]提出了基于引力場的粒子濾波算法，通過引力作用將粒子集中于高似然區域；文獻[6]提出了一種自適應智能重采樣粒子濾波算法，將粒子按權值大小分為兩個子集；文獻[7]采用了平方根嵌入式容積卡爾曼濾波產生的重要性密度函數，以此來提高濾波精度，這些算法在過程中都存在淘汰部分有效權值粒子的可能性，無法完全解決樣本貧化問題.

近年來基于智能算法優化的粒子濾波成為一個重要的研究方向.很多學者將遺傳算法、模擬退火算法、粒子群算法等引入到粒子濾波的重采樣過程中，利用智能算法尋優引導粒子向高似然區域運動，這樣可以提高粒子濾波的效率同時保持粒子的多樣性.文獻[8]提出一種利用遺傳算法中變異、交叉操作的粒子濾波算法，修正權值較小的粒子，仿真表明該算法可以有效地保持粒子的多樣性.文獻[9]將粒子群算法引入到粒子濾波，文獻[10]將蝙蝠算法引入到粒子濾波.文獻[11]將一種改進的人工螢火蟲算法和粒子濾波成功的結合，有效地提高了預測精度.群智能算法如果不針對粒子濾波的求解過程進行尋優策略的改進，就不能提高算法后期的尋優能力，在迭代過程容易陷入局部最優，如何快速尋找到全局最優解是亟待解決的問題.

磷蝦算法[12](krill herd algorithm，KH)是由Gandomi和Alavi在2012年提出的一種群體智能優化算法.該算法來自于模擬地球南極海洋中的磷蝦群體的行為，文獻[12]通過仿真實驗和常見的智能算法進行了對比.結果表明，磷蝦算法具有更好的魯棒性和全局搜索能力.但是直接將磷蝦算法和粒子濾波結合容易出現局部最優現象，影響濾波精度.文獻[13]提出了一種磷蝦群免疫粒子濾波算法，該算法使用人工免疫算法中的變異操作優化了磷蝦個體的分布情況，并將改進的算法用于粒子濾波中.但是根據文獻[12]的表述，只使用遺傳算法中的交叉操作效果更好.此外，該算法并未對磷蝦個體位置變化的過程進行改進，無法在粒子濾波的不同時刻實時調整，影響了算法的性能.文獻[14]在每次迭代中檢測磷蝦個體的運動狀態來動態調整慣性權重，減少算法的無效迭代次數，但這種調整并不完全適用于粒子濾波.文獻[15]將磷蝦種群分為不同的子種群，進行二次尋優，提高了磷蝦算法的局部收斂精度，但會降低收斂速度，不適合粒子濾波.

本文結合粒子濾波和磷蝦算法的特點，提出了一種改進的磷蝦算法優化粒子濾波 (improved krill herd algorithm optimized particle filter，IKH-PF)，動態地更新磷蝦算法中誘導和覓食行為中的權值，并將改進的遺傳算法交叉操作引入磷蝦個體的更新過程，保證前期快速全局尋優，后期高精度局部尋優，提高了粒子濾波的精度，同時在尋優過程中并未舍棄低權值粒子，保存了全部粒子的信息，有效解決了樣本貧化問題.

1 粒子濾波算法

粒子濾波基本思想是采用一組狀態空間中的隨機樣本對條件后驗概率密度函數進行近似，利用樣本均值代替傳統的積分運算，從而獲得對非線性系統的狀態估計如下：

xk=f(xk-1,vk-1)，

(1)

yk=h(xk,wk).

(2)

式中:xk為狀態值；f(·)為狀態函數；vk-1為系統噪聲；yk為觀測值；h(·)為觀測函數；wk為量測噪聲.

假設狀態初始概率密度p(x0|y0)=p(x0)，則預測方程為

狀態的更新方程為

其中

重要性函數q(x0:k|y1:k)可以改寫成如下形式:

則權值公式如下

(3)

式中δ(·)為狄克拉函數，重要性權值更新如下

當粒子分布不滿足貝葉斯濾波理論時，對粒子權值進行補償和更新為

(4)

最后進行權值歸一化，然后輸出估計的狀態如下：

(5)

(6)

2 磷蝦算法

磷蝦算法中，磷蝦個體的位置變化主要取決于以下3個方面: 1)磷蝦個體間的誘導運動；2)磷蝦個體的覓食運動；3)磷蝦個體的隨機擴散運動.

具體公式如下

ΔXi=Ni+Fi+Di.

(7)

式中：ΔXi為第i個磷蝦個體的位置變化;Ni為磷蝦個體受誘導運動引起的位置變化；Fi為磷蝦個體的覓食運動引起的位置變化；Di為個體的隨機擴散運動引起的位置變化.

2.1 磷蝦個體間的誘導運動

磷蝦個體受到其他個體影響的行為可以表示為:

(8)

(9)

在磷蝦群體中，相鄰個體之間的影響定義如下：

(10)

(11)

(12)

式中：rand為0和1之間的隨機數；I為當前迭代次數；Imax為最大迭代次數.

2.2 磷蝦個體的覓食運動

磷蝦個體的覓食行為引起的位置變化可以表示為

(13)

2.3 磷蝦個體隨機擴散運動

磷蝦個體的物理隨機擴散過程可以表示為

(14)

式中：Dmax為最大擴散速度；δ為隨機的矢量方向，δ∈[-1,1].

3 磷蝦算法優化的粒子濾波

在磷蝦算法中，將粒子濾波中粒子的狀態值作為磷蝦群的個體位置，則個體適應度值最好的磷蝦相當于權值最高的粒子.由此，粒子濾波的求解過程就可以轉化為磷蝦群的尋優過程.針對粒子濾波的特點，需要對磷蝦算法進行改進，來提高算法的尋優效率和保持粒子多樣性.

3.1 誘導權值和覓食權值的改進

根據文獻[16]的證明，磷蝦個體適應度值可以表示為：

磷蝦算法參數調整研究還很少，IKH-PF在求解過程中用到的適應度函數都為單峰值函數，本文主要參考文獻[16]對于磷蝦算法參數在單峰值函數下的測試結果，對權值wn和wf進行動態調整，設計了新的更新策略.當權值wn和wf取0.9時，在前幾次迭代中，磷蝦群就能迅速尋優；當權值wn和wf取0.1和0.3之間的值時，磷蝦群的尋優結果穩定且準確；當權值采用線性遞減時，算法尋優結果平穩，魯棒性強.綜上所述，本文提出線性函數和倒數函數相結合的非線性遞減策略，更新公式如下

(15)

式中：I、Imax分別為當前迭代次數和最大迭代次數；w1、w2為需要設置的參數.

當w1=0.2,w2=0.6,Imax=20時，得到的權值變化曲線如圖1所示.

從曲線圖可以看出，在前5次迭代中，權值迅速減少，后10次迭代權值近似線性遞減，符合之前的理論分析.

圖1 權值變化曲線

3.2 磷蝦群多樣性改進

為了保持磷蝦種群的多樣性，需要將遺傳算法引入磷蝦種群尋優中.適用于磷蝦算法的遺傳算法包括交叉和變異操作，只執行交叉操作的效果更好[12].經典的遺傳算法采用固定的交叉概率，但這種策略并不符合粒子濾波過程對于磷蝦算法尋優的要求.因此，本文提出一種改進的交叉概率變化公式，提高磷蝦群的多樣性，最終提高算法的性能并保持粒子的多樣性.

交叉概率可以控制磷蝦種群新個體產生的速度.針對粒子濾波的應用環境，尋優初期，交叉概率取較大值，提高磷蝦種群新個體產生的速度來保持種群的多樣性，從而提高算法前期的全局尋優速度；尋優后期，交叉概率取較小值，保存較優的磷蝦個體，從而提高局部收斂能力.由此設計的交叉概率更新公式如下

PC=PC0·e-2I/Imax,

(16)

式中PC0為需要設置的參數.

當PC0=0.9，Imax=20時，得到的交叉概率變化曲線如圖2所示.

圖2 交叉概率變化曲線

從圖2可以看出，隨著迭代次數的增加，交叉概率逐漸減小，符合粒子濾波求解過程的要求，同時在后期仍然能保持粒子的多樣性.

3.3 算法收斂性分析

在磷蝦算法中，令d(x(t),x*(t))表示當前磷蝦個體位置與全局最優位置之間的距離，簡化表示為d(x(t))，當x(t)≠x*(t)時，d(x(t))>0.對于一次迭代的解集X(t)={x(1),x(2),…,x(N)}，令d(x(t))=min{d(x(t)),x∈X}，可以得到每次迭代種群間的移動Δd(x(t))=d(x(t+1))-d(x(t)).以此定義迭代終止時間為

τ=min{t:d(x(t))=0}，

式中τ為第1個磷蝦個體到達全局最優位置的時刻.如果E[τ]存在有限值，則算法最終就是全局收斂的.

文獻[15]已經證明了標準磷蝦算法滿足以下兩個約束條件，并進一步證明在約束條件下算法是收斂的.

條件1存在一個規模為m的多項式λ0(m)>0，使得d(x(t))≤λ0(m).

該條件說明每代磷蝦個體和全局最優位置之間的距離在一個多項式范圍內.

條件2存在一個規模為m的多項式λ1(m)>0，使得：

該條件說明當磷蝦個體不斷接近全局最優位置時，種群的移動限定在一個多項式倒數的范圍內.

滿足以上兩個條件時，可以給出磷蝦算法的時間估計E[τ]為

E[τ|d(x(0))>0]≤λ(m)，

因此，算法是收斂的.

本文對于誘導和覓食運動中權值更新策略的改進，只是改變了不同時期磷蝦個體的運動速度，權值本身的取值仍在基本磷蝦算法的范圍之內，因此這項改進使算法仍然滿足條件1和2；磷蝦群多樣性改進中對于交叉概率更新的改進，同樣只是在不同時期采用不同的交叉概率，以增強群體的多樣性，對磷蝦群整體的移動沒有改變，因此沒有改變條件1和條件2.綜上所述，本文提出的IKH-PF在迭代過程中，磷蝦群是朝著全局最優位置移動的，存在可估的終止時間，算法是收斂的.

3.4 算法實現步驟

改進的磷蝦算法優化的粒子濾波具體實現如下：

步驟3更新粒子狀態. 根據提出的式(15)計算誘導和覓食權值，然后按照式(8)、(13)、(14)計算磷蝦個體的誘導、覓食和擴散運動，最后根據式(7)計算每次迭代中磷蝦個體的位置變化.

步驟4按照式(16)計算交叉概率，對磷蝦群個體進行交叉操作.

步驟5設置一個迭代次數，判斷是否達到設置的迭代次數，否則執行步驟3.

步驟6重新計算每個粒子的權值，按照式(4)對粒子權值進行補償和更新，按照式(5)對權值歸一化處理，最后按照式(6)輸出每個粒子當前時刻估計的狀態值.

步驟7重復執行步驟2～6，直到k=NK，完成NK個時刻的粒子狀態估計.

由于磷蝦算法的收斂性強，IKH-PF的終止條件為達到設置的迭代次數，使得磷蝦個體向最優區域移動，但避免最終收斂.這樣不僅可以保證粒子全部分布于高似然區域附近，保持粒子的多樣性，同時也可以降低算法的復雜度.

4 仿真實驗和分析

為了分析和驗證算法的性能，選取一種單靜態非增長模型，仿真環境為Matlab 2016a.

仿真實驗中采用的單靜態非增長模型的狀態方程和觀測方程如下：

8cos[1.2(t-1)]+ω(t),

式中:ω(t)、v(t)均為零均值高斯白噪聲，ω(t)的方差Q在不同仿真中分別設為1和5，v(t)的方差R=1.在磷蝦算法中，僅對粒子狀態值進行估計，所有參數均沒有單位.最大誘導速度Nmax設為0.2，個體覓食速度Vf設為0.1，最大隨機擴散速度Dmax設為0.05，最大迭代次數Imax=20，w1、w2分別設置為0.2和0.6.

為了進一步驗證算法的性能，在仿真實驗中同時和標準的粒子濾波、粒子群優化粒子濾波(PSO-PF)和蝙蝠算法優化粒子濾波(BA-PF)進行了對比.

4.1 算法的精度測試

算法的精度測試主要在以下3個條件下進行：

1)粒子數N=20，Q=1，仿真結果和估計誤差絕對值如圖3、4所示.

圖3 濾波狀態估計(N=20，Q=1)

圖4 估計誤差絕對值(N=20，Q=1)

2)粒子數N=50，Q=1，估計誤差絕對值如圖5所示.

圖5 估計誤差絕對值(N=50，Q=1)

3)粒子數N=100，Q=1，估計誤差絕對值如圖6所示.

圖6 估計誤差絕對值(N=100，Q=1)

不同條件下整體誤差采用均方根誤差為

表1 整體誤差對比

從圖3～6可以看出，在相同的粒子數條件下，相對于標準PF、PSO-PF和BA-PF，IKH-PF具有更好的狀態估計精度，這是因為磷蝦算法迭代尋優之后的粒子分布更加合理，從而提高了粒子濾波的精度.從圖4～6的誤差絕對值可以看出，相對于PSO-PF和BA-PF，IKH-PF由于改進了誘導和覓食權值更新策略，在一些容易陷入局部最優的位置，具有更好的收斂精度.從表1的數據可以看出，當實驗次數為50次時，本文提出的IKH-PF相對于其他3種算法，整體誤差是最小的；當噪聲方差變大時，IKH-PF誤差也沒有明顯變化，體現了該算法的穩定性.由于動態調整了權值和交叉概率，IKH-PF在粒子數為20時的精度已經高于其他算法在粒子數更多情況下的精度，體現了算法的運行效率.

4.2 粒子多樣性測試

為測試粒子的多樣性，將仿真的最大運行時間步長設置為50，粒子數設置為100.這里選取了運行步數k=8，27，47時的粒子分布情況.粒子分布情況如圖7～9所示.

圖7 k=8時粒子分布

圖8 k=27時粒子分布

圖9 k=47時粒子分布

從圖7～9可以看出，在標準的粒子濾波算法中，由于重采樣時只采用權值較大的粒子進行狀態估計，導致粒子分布于少數狀態值上，出現了樣本貧化現象，不利于整體狀態估計.對于進行了交叉操作的磷蝦算法，其優化的粒子大部分處于高似然區附近，也有少部分粒子分布于其他區域，整體保持了粒子的多樣性，有效地避免了樣本貧化現象.

5 結 論

1)為解決粒子濾波中的樣本貧化問題，本文結合粒子濾波和磷蝦算法的特點，采用改進的磷蝦算法優化粒子濾波.

2)針對粒子濾波，在磷蝦算法中改進了誘導和覓食運動的權值更新策略；同時為了保持粒子的多樣性，對磷蝦個體進行交叉操作，設計了新的交叉概率更新公式，以磷蝦群的尋優引導粒子向高似然區域移動.

3)仿真結果表明，提出的IKH-PF相對于標準的粒子濾波、PSO-PF和BA-PF具有更高的狀態估計精度，同時在運行后期能夠保持較好的粒子多樣性，避免了樣本貧化現象.