基于PHD的粒子濾波檢測前跟蹤改進算法

(1.海軍航空大學， 山東煙臺 264001;2.中國民航科學技術研究院機場研究所， 北京 100028)

0 引言

現代戰場環境下隨著以無人機為代表的航空器的廣泛應用，其雷達散射截面積較小，回波信號功率、信噪比會大大降低，從而導致雷達的檢測性能急劇下降[1-2]。檢測前跟蹤(TBD)方法利用多幀原始觀測數據聯合實現檢測和跟蹤過程，與傳統檢測方法相比，針對低信噪比目標的檢測性能更加優異，受到越來越多的關注[3]。典型的TBD方法包括Hough變換[4]、動態規劃[5]、粒子濾波[6]等實現形式。由于當今多目標跟蹤[7]問題中，應用隨機集理論(RFS)能夠避免數據關聯操作，將隨機集理論遷移到TBD領域也已成為了一大研究趨勢。

為減少RFS理論中多目標最優貝葉斯遞推的算法復雜度，Mahler[8]率先提出次最優的概率假設密度(PHD)濾波算法。該濾波迭代多目標后驗概率密度的一階矩，完成多目標的數目估計；其次根據目標估計數，在概率假設密度函數的峰值點處提取目標狀態[9]。作為RFS理論框架下的典型濾波器，PHD濾波器自始至終迭代傳遞多目標整體的后驗概率密度，因此在跟蹤多個數量變化的目標同時，能有效地規避復雜的數據關聯與配對，在保證算法實時性的同時，有效地提高跟蹤精度[10]。將PHD濾波器應用至TBD中，可以更好地解決多目標TBD問題，但該方面研究仍處于起步階段，性能還需不斷完善，為此不少文獻[11-12]提出了諸多改進方法。

本文針對現有基于PHD的粒子濾波檢測前跟蹤的位置估計精度較低的缺陷進行改進，受文獻[13]的啟發，結合當前量測值，在濾波的預測和更新步驟之間加入多目標粒子群優化，使粒子經過量測驅動向后驗概率密度較大方向移動；然后在目標狀態提取階段，使用基于密度的DBSCAN數據聚類方法，將經由粒子群優化算法改良分布形態的粒子進行聚類提取目標的估計狀態。經過仿真實驗表明，改進后的算法相比原有的PHD-TBD在多目標位置估計上性能改善明顯。

1 系統運動模型和觀測模型

1.1 目標運動模型

目標運動模型設計為CV運動模型。多目標群體xk,1,…,xk,Num中第t個目標的運動方程為

xk,t=Fkxk-1,t+wk

(1)

1.2 雷達系統觀測模型

借鑒文獻[14]在關于多普勒雷達的觀測模型，雷達觀測范圍內共有Ur個距離分辨單元，Ud個多普勒單元和Ub個方位單元。在雷達觀測的Ur×Ud×Ub個單元中，得到的功率觀測值為

(2)

(3)

式中:R，D和B表示與距離、多普勒和方位單元尺寸分辨率，分別根據帶寬、積累時間和波束寬度決定；Lr，Ld和Lb分別表示3個觀測維度上的損耗系數；rk,t，dk,t和bk,t表示k時刻目標t所處的距離、多普勒和方位單元;ri，dj和bl表示(i,j,l)單元的距離、多普勒和方位信息。

依據量測模型，目標似然函數可表示為

p(zk|xk,1,xk,2,…,xk,Num)=

(4)

式中，p0表示沒有目標時(i,j,l)單元處噪聲的似然函數，p1表示存在Num個目標時(i,j,l)單元目標加噪聲的似然函數。

(5)

(6)

由此得到似然比函數Lk的計算如下：

(7)

2 改進的雷達PSO-PHD-TBD算法

改進之后的PSO-PHD-TBD算法主要借助于多目標粒子群優化算法以及基于密度聚類的DBSCAN算法相結合提高對目標的位置估計精度。2.1節著重描述多目標粒子群優化算法，2.2節強調對粒子進行基于密度聚類，2.3節展示PSO-PHD-TBD的具體實施步驟，并在3.2節3)中給出了具體的算法復雜度的分析。

2.1 粒子群優化算法

粒子群優化算法(Particle Swarm Optimization, PSO)通過粒子個體通過適應度目標函數的函數值信息進行粒子優選，完成對最優條件的搜索[15]。算法中每個粒子i具有位置和搜索速度，代表最優結果的可能取值。算法首先定義適應度目標函數，通過迭代尋找最優解。迭代操作按照式(8)、式(9)更新每個粒子的運動信息，使得每個粒子移動，經過對粒子的個體最優解pi篩選，在搜索全局極值gi的同時實現對集群分布的優化[16]。

vi=wvi+c1R1·(pi-xi)+c2R2·(gi-xi)

(8)

xi+1=xi+rvi

(9)

式中，R1和R2在區間(0,1)內隨機取值，w為慣性系數，c1和c2為正實數，代表粒子向個體最優解和全局最優解移動的權重，xi和xi+1為粒子狀態，vi為搜索速度，r為更新因子。

在多目標的檢測前跟蹤過程中，由于存在多個適應度目標函數之間相互沖突，一個目標函數的最優解往往并不能滿足剩余目標函數最優解，因此每一個目標的最優解只能陷入局部最優解。而多目標粒子群優化(Multi-Object PSO, MOPSO)方法，借助NSGA-Ⅱ算法實現利用多個適應度目標函數篩選粒子，能夠將粒子引導、移動以求取多目標Pareto最優解集[17-18]。作為遺傳算法上的改進[19]，算法關鍵步驟在以下3個過程：

1) 精英選擇策略。為了避免盲目的更新迭代而采用的一種方法，將更新前后的粒子集群同時存儲進行優先選擇。

2) 快速非支配排序。將粒子群中所有粒子的解進行比較，得出粒子目標函數值之間的支配關系，對種群分層，從而指引搜索向Pareto最優解集方向進行。

3) 計算擁擠距離。計算粒子的擁擠距離是為了能夠將處于相同非支配層的粒子進行選擇性排序。對同一非支配層L中所有n個粒子的擁擠距離disi初始化。

disi=0,i=1,2,…,n

(10)

對同層的個體按照第m個目標函數值升序排列：對于排序邊緣上的粒子，令dis0=disn=D(D為一較大乘數)；對排序中間的粒子i，求取擁擠距離disi。

disi=∑mWmdisi,m

(11)

disi,m=disi,m+(disi+1,m-disi-1,m)/(fmax,m-fmin,m)

(12)

式中disi+1,m，disi-1,m為粒子i的相鄰粒子對第m個目標函數的函數值，fmax,m，fmin,m分別為層中粒子對第m個目標函數的最大值和最小值。基于所有的適應度目標函數，都需要循環上述步驟，每一個目標函數擁有與之相對應的函數權重Wm，對所有的Wmdisi,m求和得到粒子i的擁擠距離disi。

為了避免陷入局部最優解的問題，再加入變異機制增加粒子的多樣性。粒子分為三部分，一為選擇不變異，二為統一增加變異量，三為隨機添加變異量。

2.2 基于密度的數據聚類方法

粒子群優化方法使粒子都向后驗概率密度較大的區域進行移動，這樣使得基于密度進行聚類的方法得到應用。基于密度的聚類算法其實是基于距離的聚類的類型(這一點和通常使用的K-means聚類方法類似)。只要鄰近區域的粒子密度超出了閾值Thmin，就可以繼續對粒子聚類。DBSCAN是其中具有代表性的一種，該方法根據密度閾值Thmin控制聚類的增長，其主要目的是通過過濾低密度區域來發現稠密樣本點聚類。

DBSCAN方法將聚類定義為高密度相連點的最大集合，并能發現任意形狀的聚類[20]。為了實現粒子聚類，算法描述需要輸入一個如數目為N的粒子群，一個最小半徑E，密度閾值Thmin。對于尚未被聚類的粒子，算法主要統計該粒子周圍半徑E范圍內的粒子數量，如果小于Thmin，則該點為噪聲點；如果抽出的粒子其最小半徑E范圍內的粒子數目超出閾值Thmin，則該點成為一個新聚類的核心點，將所有從該點密度可達的點加入到當前聚類中，如果該點密度可達節點中有已經生成聚類的核心點，則將兩個聚類合并形成一個聚類，新聚類的核心點為原來合并的核心點的集合。

2.3 PSO-PHD-TBD算法實施步驟

借助于粒子濾波(也稱為序貫蒙特卡洛[21])的實現形式，PSO-PHD-TBD的具體實施步驟如下：

1) 初始化

(13)

2) 預測

(14)

并根據下式計算存活粒子的權重：

(15)

假設k時刻新生目標的采樣粒子數為Jk，則根據建議密度pk(·|Yk)采樣新生目標的粒子：

(16)

并根據式(17)計算新生粒子的權重：

(17)

式中,γk為新生粒子的概率假設密度函數。

3) PSO優化

在粒子的預測之后，將最新的觀測值引入優化采樣過程，根據獲得的觀測定義優化所需的適應度目標函數：

(18)

(19)

為了避免預測概率密度與似然函數的重疊區域較小的問題，在此時設置粒子群優化算法，結合量測信息，由量測驅動粒子集通過式(8)、式(9)不斷更新迭代向高似然比區域移動，提高每一個粒子的利用效率，使粒子不斷地向高似然狀態靠近[13]。具體的實施步驟如下：

①經過對父代種群粒子進行式(8)、式(9)的運算更新得到子代種群，對子代種群中的粒子進行變異操作，同時父代粒子的權重也傳遞給子代對應的粒子；

②將兩個種群合并，進行非支配排序，計算整個粒子集之間的擁擠距離；

③選擇生成新的父代種群粒子，對粒子權重重新進行歸一化處理，繼續進行下一輪的優化迭代，直到適應度函數值達到閾值Thpso。

4) 更新

對于k時刻已經過優化處理的所有粒子i=1,2,…,Lk-1+Jk，根據式(20)更新其權重：

(20)

(21)

(22)

(23)

σ2為觀測噪聲協方差，而權重表達式分母的另一項，

(24)

式中，

(25)

5) 數目估計、重采樣以及狀態提取

根據粒子權重之和估計目標數目

(26)

3 算法仿真

3.1 仿真條件設置

本節將給出多目標TBD的應用實例，仿真中假定目標運動滿足勻速模型。掃描周期T=1 s，目標運動近似為線性運動，運動方程同式(1)。假設雷達的測量區域設置為[0,2 000 m]×[0,2 000 m]，雷達位于原點位置，距離單元分辨率R=20 m，多普勒單元分辨率D=1 m/s，方位單元分辨率B=1°。為便于算法性能對比，一共產生40幀仿真數據(文中將一次掃描數據設置為一幀數據，由于掃描周期為1 s，數據信息一共歷時40 s)，經過100次蒙特卡洛實驗取均值。

仿真設置在3種不同信噪比條件下(分別為9，7和5 dB)，雷達的觀測區域存在6個目標，目標初始運動狀態如表1所示。假定目標的RCS沒有起伏，滿足Swerling 0型假設；目標存活概率PS=0.98，新生概率Pb=0.02；過程噪聲協方差為0.001；目標強度噪聲協方差為0.01；粒子數L=500；SNR=10lg(P/2σ2)，觀測噪聲協方差為1，則可依據信噪比以及觀測噪聲協方差σ2來得到目標功率。目標強度設定I=20，粒子的強度服從介于[15,25]之間的均勻分布。

表1 目標初始運動狀態

每個目標分別設置1 000個粒子，設置高斯混合模型從中采樣得出500個新生粒子，高斯混合模型一共包含3種高斯分布x(i)～N(x;m(i),P),i=1,2,3相組合，具體參數分別為m(1)=[250 m,0 m/s,1 500 m,0 m/s,I]，m(2)=[1 000 m,0 m/s,1 000 m,0 m/s,I]，m(3)=[1 500 m,0 m/s,250 m,0 m/s,I]，P=diag([20 m 20 m/s 20 m 20 m/s 1])2，3種目標的權重一致。多目標粒子群優化過程中，在一個時刻k內迭代次數gen=10，目標函數閾值Thpso=0.7，達到迭代次數后以適應度函數閾值為優化結束的最終條件。DBSCAN聚類時設置參數最小半徑E=5 m和密度閾值Thmin=10。

整個仿真場景不考慮新目標衍生。目標真實運動狀態如圖1所示，圖中三角形表示目標起始位置，正方形表示消失位置。

圖1 目標的真實航跡

圖2表示9 dB信噪比條件下，監視區域在第25幀時已轉換為笛卡爾坐標系的回波能量信息，6個目標都悉數出現在觀測區域，完全淹沒在噪聲之中。

圖2 第25幀的觀測信息

3.2 仿真結果分析

本文將所提出的PSO-PHD-TBD算法和常規PHD-TBD算法對同一仿真場景在信噪比進行性能分析。

1) 對目標的數目估計情況

經過100次蒙特卡洛的仿真實驗，得到平均每個時刻的估計數目如圖3所示。

兩種算法在隨著信噪比的降低，兩算法的估計性能都出現了相應的下降。相較之下，PSO-PHD-TBD的性能更加穩定。圖4列舉出不同信噪比條件下基于100次蒙特卡洛仿真計算而來的數目估計標準差，表2展示了40幀數據目標數目估計標準差的均值，統計值越小代表目標數目估計越準確。根據圖4和表2可以看出在數目估計精度上PSO-PHD-TBD算法相比于PHD-TBD算法也有較為明顯的提高。

(a) 9 dB的數目估計

(b) 7 dB的數目估計

(c) 5 dB的數目估計圖3 不同信噪比下算法對目標數目的估計

表2 兩種算法目標數目估計標準差均值

(a) 9 dB的數目估計標準差

(b) 7 dB的數目估計標準差

(c) 5 dB的數目估計標準差圖4 不同信噪比下算法對目標數目估計標準差

2) 對目標位置精度的估計情況

本文按照OSPA距離[24]作為誤差評判標準，設置OSPA距離參數截斷距離c=50 m，p=2。為了更好地進行對比說明粒子群優化與DBSCAN方法相結合的優勢，另外添加了結合粒子群優化但進行常規K-means聚類的PHD-TBD算法添加對照，按照相同的仿真場景，得出的結果如圖5所示。具體3種算法在不同信噪比下的OSPA誤差的均值如表3所示。

綜合來看，新提出的PSO-PHD-TBD算法主要在目標位置估計精度上較PHD-TBD算法有了明顯的進步，隨著信噪比的降低，精度上優勢更加明顯。

(a) 9 dB的目標位置估計誤差

(b) 7 dB的目標位置估計誤差

(c) 5 dB的目標位置估計誤差圖5 不同信噪比下3種算法的誤差比較

表3 不同信噪比下OSPA誤差均值 m

3) 算法實現復雜度分析

本文所提出的PSO-PHD-TBD的粒子濾波實現方式，借助了多目標粒子群優化，DBSCAN算法共同實現。DBSCAN算法的時間復雜度是O(n2)，相比于K-means聚類的時間復雜度O(n)更加復雜，另外基于NSGA-Ⅱ的MOPSO算法的時間復雜度為O(mn2)，其中m為目標函數的個數，n為粒子數目，通過理論分析PSO-PHD-TBD算法運行時間確實更長，最為理想的情況下生存粒子維持現有的分布，僅有新生粒子隨著粒子群優化而改變分布，會使得優化所用的時間最短。在Intel(R) Core(TM) i5-4590 3.30 GHz CPU，16 GB(DDR3 1600 MHz)內存，Win7 64位旗艦版計算平臺下運行不同信噪比條件下仿真的平均時長如表4所示。

表4 不同信噪比下算法運行時長 s

在信噪比降低的情況下，優化步驟會受到一定的影響，從而PSO-PHD-TBD的運算時間會隨之增加，反觀PHD-TBD算法因為不存在優化步驟，因此在時間上并沒有隨信噪比的變化而能明顯改變；因為DBSCAN的算法復雜度高于K-means聚類，因此同樣經過粒子群優化而聚類不同的算法運算效率也不同，具體而言本文提出的粒子群優化與DBSCAN算法相結合的用時更長。雖然運算成本不可避免地增加，但是很明顯的一點是本文提出的PSO-PHD-TBD的算法性能要明顯地優于PHD-TBD。在運行時間尚能容忍的情況下，利用DBSCAN聚類的PSO-PHD-TBD算法依舊是更好的選擇。

4 結束語

經過仿真對比發現，結合DBSCAN聚類的PSO-PHD-TBD算法通過粒子群優化算法結合量測優化粒子分布，再將基于密度的聚類算法運用到經過優化的粒子濾波當中可以實現較為良好的聚類效果，能夠滿足檢測目標發現其位置的要求，更是能在目標數目準確估計的情況下將估計誤差加以降低。但是DBSCAN算法對兩個初始參數鄰域半徑和鄰域最小點數還需要根據仿真過程調試，以達到最佳聚類效果。同時粒子群優化算法復雜度較高，在計算中確實犧牲了一部分計算效率，在今后的工作中則需要盡可能地解決粒子數和運算時間之間的矛盾，并希望能夠歷史信息自適應選擇粒子集群優化，避免低效率的優化。