基于ET-PHD濾波器和變分貝葉斯近似的擴展目標(biāo)跟蹤算法

何祥宇，李 靜，楊數(shù)強，夏玉杰

（1.洛陽師范學(xué)院物理與電子信息學(xué)院，河南洛陽 471934；2.洛陽師范學(xué)院信息技術(shù)學(xué)院，河南洛陽 471934）

（?通信作者電子郵箱xyuhe@foxmail.com）

0 引言

目標(biāo)跟蹤技術(shù)是一個備受國內(nèi)外學(xué)者關(guān)注的研究課題，應(yīng)用范圍十分廣泛［1-3］。多目標(biāo)跟蹤（Multi-Target Tracking，MTT）技術(shù)研究是目標(biāo)跟蹤技術(shù)領(lǐng)域內(nèi)的一個重要研究分支，多目標(biāo)跟蹤的目的是利用受雜波和噪聲影響的傳感器測量值估計目標(biāo)的個數(shù)及狀態(tài)。隨機有限集（Random Finite Set，RFS）理論被認(rèn)為是處理多目標(biāo)跟蹤問題的一種簡潔高效的方法，通過把多目標(biāo)狀態(tài)集和測量集定義為隨機有限集，該理論建立了嚴(yán)謹(jǐn)?shù)淖顑?yōu)貝葉斯多目標(biāo)濾波框架。

在基于RFS 理論的多目標(biāo)跟蹤方法中，概率假設(shè)密度（Probability Hypothesis Density，PHD）濾波器［4］是一種很有潛力的多目標(biāo)貝葉斯濾波算法，它不用復(fù)雜的數(shù)據(jù)關(guān)聯(lián)過程就可以處理復(fù)雜環(huán)境中的多目標(biāo)跟蹤問題。標(biāo)準(zhǔn)的PHD 濾波器是基于點目標(biāo)假設(shè)的多目標(biāo)跟蹤算法，點目標(biāo)假設(shè)即：在每一個時間步，每一個目標(biāo)假設(shè)只能產(chǎn)生至多一個量測值。與點目標(biāo)不同的是，一個擴展目標(biāo)在每一個時間步上可產(chǎn)生多個量測值［5］。因此，應(yīng)用標(biāo)準(zhǔn)的PHD 濾波器跟蹤多個擴展目標(biāo)會產(chǎn)生不可靠的跟蹤結(jié)果。

為在RFS理論框架內(nèi)處理多擴展目標(biāo)跟蹤問題，文獻(xiàn)［6］中提出了標(biāo)準(zhǔn)PHD 濾波器的一個擴展版本，稱之為擴展目標(biāo)PHD（Extended Target PHD，ET-PHD）濾波器。目前，ET-PHD濾波器實現(xiàn)方法主要有兩種，分別是擴展目標(biāo)序貫蒙特卡洛PHD（Extended Target Sequential Monte Carlo PHD，ET-SMCPHD）濾波器［7］和擴展目標(biāo)高斯混合PHD（Extended Target Gaussian Mixture PHD，ET-GM-PHD）濾波器［8］。ET-SMC-PHD濾波器使用隨機樣本粒子的加權(quán)和形式近似后驗強度函數(shù)，利用聚類及其優(yōu)化算法實現(xiàn)目標(biāo)狀態(tài)的抽取與估計，是ETPHD 濾波器一種通用實現(xiàn)。ET-GM-PHD 濾波器是ET-PHD濾波器在線性高斯假設(shè)條件下的封閉解，利用高斯混合函數(shù)描述目標(biāo)的后驗強度函數(shù)。此外，近期提出的擴展目標(biāo)箱粒子PHD（Extended Target Box-Particle PHD，ET-Box-PHD）濾波器［9］和橢圓擴展目標(biāo)箱粒子EET-BP-PHD（Ellipse Extended Target Box Particle PHD，）濾波器［10］是兩種新穎的擴展目標(biāo)跟蹤算法，克服了ET-SMC-PHD 濾波過程中使用大量樣本粒子近似后驗強度函數(shù)時存在計算量大的不足，提高了算法的運行效率。

標(biāo)準(zhǔn)的ET-PHD 濾波器需要測量噪聲的先驗信息，但這種測量噪聲的先驗信息在許多實際應(yīng)用環(huán)境中都是未知的。測量噪聲協(xié)方差的顯著不匹配會導(dǎo)致ET-PHD 濾波器的目標(biāo)跟蹤性能退化，從而限制了ET-PHD 濾波器的實際應(yīng)用范圍。處理未知測量噪聲問題的方法主要有交互多模型（Interactive Multiple Models，IMM）方法［11］、粒子方法［12］和變分貝葉斯（Variational Bayesian，VB）近似方法［13］等。由于具有較低的計算復(fù)雜度，VB 近似方法已經(jīng)成功地被用于基于RFS理論的濾波器中處理未知測量噪聲協(xié)方差環(huán)境中的多點目標(biāo)跟蹤問題［14-17］，并且都取得了令人滿意的跟蹤效果。

為處理未知測量噪聲協(xié)方差情況下的多擴展目標(biāo)跟蹤問題，本文提出了一種基于ET-PHD 濾波器和VB 近似的擴展目標(biāo)跟蹤算法。本文根據(jù)VB 近似理論及標(biāo)準(zhǔn)的ET-PHD 濾波器，提出了ET-PHD 濾波器的一種擴展版本，并利用相對熵最小化方法給出了其在線性高斯假設(shè)條件下的解析實現(xiàn)。本文算法是一種遞歸的多擴展目標(biāo)跟蹤算法，利用預(yù)測和更新過程實現(xiàn)后驗強度函數(shù)和測量噪聲協(xié)方差的計算與估計，進(jìn)而獲取相應(yīng)的目標(biāo)跟蹤結(jié)果。

1 擴展的ET-PHD濾波器

定義fk|k-1(xk|xk-1)和gk(zk|xk)分別為目標(biāo)狀態(tài)轉(zhuǎn)移函數(shù)和測量似然，在線性高斯假設(shè)條件下，ET-PHD 濾波器的系統(tǒng)模型可分別表示為下面形式：

其中：N(x；m，P)表示均值為m、協(xié)方差為P的高斯函數(shù)；xk和zk分別表示目標(biāo)狀態(tài)和測量值；Fk-1和Hk分別表示目標(biāo)狀態(tài)轉(zhuǎn)移矩陣和測量矩陣；Qk-1和Rk分別表示過程噪聲協(xié)方差矩陣和測量噪聲協(xié)方差矩陣。

根據(jù)ET-PHD 濾波器的系統(tǒng)模型，為描述基于VB 近似的ET-PHD 濾波器的公式體系，限定目標(biāo)狀態(tài)xk和測量噪聲協(xié)方差Rk是互相獨立的。若定義(xk，Rk)為k時間步擴展目標(biāo)的增廣狀態(tài)變量，則增廣狀態(tài)變量(xk，Rk)的聯(lián)合轉(zhuǎn)移函數(shù)

其中，fk|k-1(xk|xk-1)和fk|k-1(Rk|Rk-1)分別表示單目標(biāo)狀態(tài)和測量噪聲協(xié)方差的轉(zhuǎn)移函數(shù)。

此外，增廣狀態(tài)的新生目標(biāo)RFS限定為泊松隨機集，定義Db，k(xk，Rk)為新生目標(biāo)RFS的強度函數(shù)，擴展的ET-PHD 濾波器的遞歸公式描述如下。

假設(shè)Dk-1(xk-1，Rk-1)為k-1 時間步的目標(biāo)聯(lián)合后驗強度函數(shù)，則k時間步預(yù)測的強度函數(shù)Dk|k-1(xk，Rk)可寫為：

若在k時間步的測量集合為Zk，則聯(lián)合后驗強度函數(shù)

其中：pD，k和κk(zk)分別表示傳感器的檢測概率和泊松雜波RFS 的強度；γ(xk)表示擴展目標(biāo)xk產(chǎn)生的測量個數(shù)期望值；g(zk|xk，Rk)為測量zk的似然函數(shù)；p∠Zk表示劃分；|W|表示子集合W中元素個數(shù)；ωp和dW分別為劃分p和子集W的非負(fù)系數(shù)。ωp和dW計算式分別為：

其中：當(dāng)|W|=1時δ|W|，1=1，否則δ|W|，1=0。

在更新式（5）～（8）中，似然g(zk|xk，Rk)是目標(biāo)狀態(tài)xk和測量噪聲協(xié)方差變量Rk的耦合函數(shù)，聯(lián)合強度函數(shù)Dk(xk，Rk)是一種非解析形式。

2 基于VB近似的解析實現(xiàn)

擴展的ET-PHD 濾波器解析實現(xiàn)的一次循環(huán)過程包含預(yù)測和更新兩個過程，解析實現(xiàn)的基礎(chǔ)除線性高斯模型之外，還需限定目標(biāo)的后驗強度函數(shù)和新目標(biāo)強度函數(shù)為高斯與逆伽馬（Inverse-Gamma，IG）混合分布的形式。由此，定義k時間步的新目標(biāo)強度Db，k(xk，Rk)為下面的函數(shù)形式：

根據(jù)上述的假設(shè)與定義，擴展的ET-PHD 濾波器解析實現(xiàn)的預(yù)測過程和更新過程詳述如下。

2.1 預(yù)測過程

若k-1 時間步的強度函數(shù)Dk-1(xk-1，Rk-1)可表示為下面的函數(shù)形式：

則k時間步預(yù)測的強度函數(shù)Dk|k-1(xk，Rk)可寫作：

其中Db，k(xk，Rk)的表達(dá)式與式（9）相同，且有：

式中，ρ(i，j)表示取值范圍為(0，1]的衰減因子。

根據(jù)式（11）～（17），式（11）描述的預(yù)測強度函數(shù)Dk|k-1(xk，Rk)可以表示為下面的函數(shù)形式：

利用式（18）給出的預(yù)測的聯(lián)合強度函數(shù)Dk|k-1(xk，Rk)的表達(dá)式，預(yù)測的測量噪聲協(xié)方差的計算公式為：

其中diag(?)表示構(gòu)建對角矩陣的運算。

2.2 更新過程

在更新過程，k時間步的測量集利用距離劃分方法［8］進(jìn)行劃分，進(jìn)而實現(xiàn)后驗強度函數(shù)Dk(xk，Rk)的計算，將式（6）代入式（5）可得：

若k時間步的預(yù)測的聯(lián)合強度函數(shù)Dk|k-1(xk，Rk)的函數(shù)形式與式（18）相同，將式（18）代入式（22）可得：

其中：為方便描述強度函數(shù)(xk，Rk)的函數(shù)表達(dá)式，定義如下符號：

其中⊕表示垂直向量連接。

利用文獻(xiàn)［13］給出的相對熵最小化結(jié)論，聯(lián)合強度函數(shù)，Rk)可表示為：

為獲取最優(yōu)的計算結(jié)果，聯(lián)合強度函數(shù)(xk，Rk)表達(dá)式中的各參數(shù)采用定點迭代算法進(jìn)行計算，定點迭代算法描述如下：

3）利用下面公式進(jìn)行若干次的迭代計算：

其中blkdiag(?)表示構(gòu)建分塊對角矩陣。

定點迭代算法的迭代次數(shù)可根據(jù)實際需要進(jìn)行設(shè)置，此外式（32）中的其他參數(shù)計算如下：

計算得到的強度函數(shù)Dk(xk，Rk)用于獲取擴展目標(biāo)的狀態(tài)和個數(shù)，并通過修剪和合并限制Dk(xk，Rk)中高斯與逆伽馬混合分量的個數(shù)，進(jìn)而提高運算效率，實現(xiàn)未知噪聲環(huán)境中多個擴展目標(biāo)的實時有效的跟蹤。

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

利用兩個仿真實驗驗證本文算法的有效性，實驗場景設(shè)置為包含雜波且目標(biāo)個數(shù)隨機變化的一個二維平面。實驗環(huán)境為：英特爾八代酷睿i5-8250 處理器、8 GB 運行內(nèi)存、Windows 10操作系統(tǒng)和R2012a版本的Matlab仿真軟件。

k時間步的目標(biāo)狀態(tài)向量xk由目標(biāo)位置和速度構(gòu)成，測量值zk定義為受噪聲影響的目標(biāo)位置信息。仿真實驗中采用式（1）和式（2）中描述的目標(biāo)跟蹤模型，模型中的各個參數(shù)定義如下：

抽樣間隔定義為Δ=1 s，目標(biāo)的存活和檢測概率分別設(shè)置為pS，k=0.99 和pD，k=0.99，每個目標(biāo)產(chǎn)生的測量個數(shù)的泊松率設(shè)置為γ=10。標(biāo)準(zhǔn)差σq已知，值設(shè)置為σq=1.95 m；跟蹤過程中標(biāo)準(zhǔn)差σr的值是未知的。仿真實驗中用到的新目標(biāo)強度函數(shù)定義為：

強度函數(shù)中高斯分量的修剪閾值為Tw=10-3，合并閾值為U=4，高斯分量個數(shù)最大值設(shè)置為Jmax=100。逆伽馬分布的參數(shù)初始值設(shè)置為α0=β0=1，退化因子的值定義為ρ=0.9。

不同算法的目標(biāo)位置估計精度利用最優(yōu)子模式分配（Optimal Sub-Pattern Assignment，OSPA）距離［18］進(jìn)行評估，計算公式為：

其中：X和Y分別表示目標(biāo)的真實狀態(tài)集和狀態(tài)估計集；p和c分別設(shè)置為p=2和c=3 000 m。

仿真場景中設(shè)定有三個目標(biāo)運動，圖1 為三個目標(biāo)的真實運動軌跡。圖1所示的仿真場景中，目標(biāo)1的運動時間設(shè)置為50個持續(xù)時間步，目標(biāo)2的運動時間在5到45時間步之間，目標(biāo)3在時間步10出現(xiàn)。

圖1 真實目標(biāo)軌跡Fig.1 True target tracks

3.1 固定參數(shù)實驗

本實驗中雜波設(shè)置為速率為r=10泊松隨機有限集，ETGM-PHD 濾波器中的標(biāo)準(zhǔn)差設(shè)置為σr=0.65。圖2 給出的是一次實驗運行中50 個時間步的真實目標(biāo)軌跡和包含雜波量測的傳感器測量值。

圖2 真實目標(biāo)軌跡與測量值Fig.2 True target tracks and measurements

圖3為不同算法在一次實驗運行中目標(biāo)的PHD權(quán)值。從圖3 中可以看出，在17 和22 等時間步上，ET-GM-PHD 濾波器估計的不同目標(biāo)的PHD 權(quán)值存在不同程度的減小，但本文算法的PHD 權(quán)值估計卻十分可靠。原因是ET-GM-PHD 濾波器使用了不匹配的測量噪聲參數(shù)進(jìn)行濾波計算，導(dǎo)致了濾波性能的下降；但本文算法在跟蹤過程能夠獲得合適的測量噪聲參數(shù)，進(jìn)而取得了可靠的目標(biāo)PHD權(quán)值。

圖3 本文算法和ET-GM-PHD濾波器中不同目標(biāo)的PHD權(quán)值Fig.3 PHD weights of different targets of proposed algorithm and ET-GM-PHD filter

圖4給出的是一次實驗運行中不同算法在50個時間步上的目標(biāo)位置估計。從圖4 中可以看出，在17、22 和37 時間步上，ET-GM-PHD濾波器未能獲取相應(yīng)的目標(biāo)位置，原因是ETGM-PHD濾波器中目標(biāo)PHD權(quán)值的大幅下降造成了這種目標(biāo)位置估計丟失現(xiàn)象的產(chǎn)生，但本文算法在跟蹤過程中沒有出現(xiàn)目標(biāo)丟失的現(xiàn)象。

圖4 真實目標(biāo)軌跡與位置估計Fig.4 True target tracks and position estimations

在實驗中，不同算法的平均性能從平均的目標(biāo)個數(shù)估計和平均的OSPA 距離這兩個方面進(jìn)行分析，相應(yīng)的仿真結(jié)果分別由圖5 和圖6 給出。從圖5 中可以看出，在11～50 時間步上，ET-GM-PHD 濾波器平均的目標(biāo)個數(shù)估計與目標(biāo)個數(shù)的真實值有著明顯的偏差，該偏差產(chǎn)生的原因是ET-GM-PHD濾波器中的不合適的目標(biāo)PHD 權(quán)值估計。通過對比發(fā)現(xiàn)，本文算法平均的目標(biāo)個數(shù)估計非常接近真實的目標(biāo)個數(shù)，偏差很小。圖6 顯示，本文算法的平均OSPA 距離遠(yuǎn)小于ET-GM-PHD 濾波器的平均OSPA 距離，這種情況在11～50 時間步之間尤為明顯。

圖5 不同算法平均的目標(biāo)個數(shù)估計Fig.5 Average target number estimations of different algorithms

圖6 不同算法的平均OSPA距離Fig.6 Average OSPA distances of different algorithms

3.2 變化參數(shù)實驗

本仿真實驗中，不同算法的仿真結(jié)果為不同參數(shù)取值情況下進(jìn)行的1 000 次仿真實驗中獲取的平均值，ET-GM-PHD濾波器中標(biāo)準(zhǔn)差σr取值分別為0.65、0.75、0.95 和1.65 四種情況。首先，目標(biāo)量測的泊松率固定為γ=10，雜波速率在5到50 之間以等變化間隔5 進(jìn)行取值，不同算法在不同的雜波速率取值情況下時間上平均的OSPA 距離的對比結(jié)果如圖7所示。然后，雜波速率固定為r=10，目標(biāo)量測的泊松率在3到12 之間變化，不同目標(biāo)量測泊松率時的時間平均的OSPA距離由圖8給出。

圖7 顯示，不同算法的跟蹤精度隨著雜波速率的增大而下降，但本文算法時間平均的OSPA 距離在不同的雜波速率上總是一直小于ET-GM-PHD 濾波器相應(yīng)的OSPA 距離。圖8顯示，當(dāng)目標(biāo)量測的泊松率在10 附近取值時，所有算法都可以獲取最佳的跟蹤精度；但對于所有的目標(biāo)量測泊松率，本文算法都可以獲得精度更好的目標(biāo)位置估計。根據(jù)圖7 和圖8給出的仿真對比結(jié)果可以得到，在未知測量噪聲協(xié)方差環(huán)境中，本文算法能夠獲得可信賴的多擴展目標(biāo)跟蹤結(jié)果。

圖7 不同雜波速率下的時間平均的OSPA距離Fig.7 Time averaged OSPA distances with different clutter rates

圖8 不同目標(biāo)量測泊松率下的時間平均的OSPA距離Fig.8 Time averaged OSPA distances with different Poisson rates of target measurement

4 結(jié)語

為處理未知測量噪聲協(xié)方差情況下的多擴展目標(biāo)跟蹤問題，利用VB近似理論和高斯與逆伽馬混合分布表示的聯(lián)合強度函數(shù)，本文提出了一種擴展的ET-PHD 濾波器及其線性高斯解析實現(xiàn)。仿真結(jié)果表明：本文算法能夠以可靠的性能跟蹤未知測量噪聲協(xié)方差環(huán)境中個數(shù)未知且隨時間變化的多個擴展目標(biāo)。在本文研究的基礎(chǔ)之上，下一步將研究ET-SMCPHD 濾波器的算法優(yōu)化方法，以提高ET-SMC-PHD 濾波器的算法運行速度。