李文娟, 呂 婧, 顧 紅, 蘇衛民
(南京理工大學電子工程與光電技術學院, 江蘇 南京 210094)
多目標跟蹤是為了在復雜背景下獲得準確的目標數目和狀態。近年來,基于隨機有限集(random finite set,RFS)的多目標跟蹤算法[1-3]發展迅猛,主要有序貫蒙特卡羅(sequential Monte Carlo, SMC)[4-6]和高斯混合(Gaussian mixture, GM)[7-9]兩種實現方式。概率假設密度(probability hypothesis density,PHD)[10-13]傳播一階矩,是最簡單的RFS類跟蹤算法。
占據多個分辨單元的目標稱為擴展目標[14]。在多擴展目標跟蹤問題中,對擴展目標的量測數目和擴展建立合理模型,能夠提高跟蹤性能[15-17]。本文聚焦擴展目標的量測數目。文獻[14]提出擴展目標的量測數目服從泊松分布。隨后,文獻[18]提出一種泊松擴展目標概率假設密度(Poisson extended target PHD, PET-PHD)濾波算法。然而,這種算法在跟蹤鄰近擴展目標時不能精確估計目標數目[19-20]。針對這一問題,文獻[20]提出擴展目標的量測數目服從二項分布的假設。基于以下兩點原因,文獻[20]認為二項分布比泊松分布更為合理。一方面,一個擴展目標占據多個分辨單元,一般假設分辨單元互相獨立且檢測概率相同[14],擴展目標的量測產生過程顯然是一個多伯努利實驗。另一方面,泊松分布是二項分布在擴展目標占據的分辨單元數趨于無窮大,檢測概率趨于無窮小的極限情況下近似得到的一種分布[21]。然而,在實際應用中,擴展目標占據的分辨單元數有限,不可能趨于無窮大。因此,擴展目標的量測數目服從泊松分布與實際情況不符,服從二項分布更為合理。文獻[20]的仿真結果表明,服從基于二項分布的ET-PHD(binominal distribution based ET-PHD,BET-PHD)濾波算法能夠獲得比PET-PHD更好的跟蹤性能。由于二項分布是一種特殊的多伯努利分布[1],文獻[20]提出的濾波算法也稱為多伯努利ET-PHD(multi-Bernoulli ET-PHD, MB-ET-PHD)。文獻[22]首次將二項分布量測數目模型和多伯努利濾波器結合,用于跟蹤雜波環境下單個擴展目標。
在BET-PHD濾波算法中,擴展目標的檢測概率和量測數目最大值是非常重要的先驗參數。然而,在實際應用中這兩個參數事先未知。參數嚴重不匹配會導致BET-PHD算法跟蹤性能急劇下降。由于文獻[22]已經給出量測數目最大值的估計方法,本文聚焦未知目標檢測概率的估計方法:通過二項分布的共軛先驗貝塔分布和似然函數二項分布構造貝葉斯迭代公式,同時將檢測概率的估計方法引入BET-PHD的高斯混合濾波器,提出一種貝塔高斯ET-PHD(beta Gaussian ET-PHD,BG-ET-PHD)濾波器。通過引入貝塔分布估計檢測概率已經不是首次。文獻[23]提出一種聚焦點目標跟蹤問題、自適應估計檢測概率的BG-PHD濾波器。在PET-PHD濾波器中,作為先驗信息的目標量測數目在實際應用中也是事先未知的。相應的自適應估計目標量測數目的濾波器稱為伽馬高斯ET-PHD(gamma Gaussian ET-PHD, GG-ET-PHD)[24]。仿真結果表明,BG-ET-PHD濾波器估計參數準確,能夠獲得比GG-ET-PHD更好的跟蹤性能。
給定k時刻的任一擴展目標xk,假設該目標占據γk(xk)個分辨單元,分辨單元之間互相獨立且檢測概率pD,k(xk)相同[14]。量測的產生過程可以看作一個分辨單元上的γk(xk)次獨立重復檢測實驗,并且每次實驗只有兩種結果:檢測出量測和沒有檢測出量測。顯然,這是統計學的多伯努利實驗[21]。因此,擴展目標的量測數目應該服從二項分布,其最大值為γk(xk),平均值為λk(xk)=γk(xk)pD,k(xk)。
假設擴展目標的量測數目服從二項分布,相應的濾波算法稱為BET-PHD[20]。其GM實現方式要求系統線性高斯,狀態方程和測量方程分別為
fk|k-1(xk|k-1|xk-1)=N(xk|k-1;Fk-1xk-1,Qk-1)
(1)
gk(zk|xk|k-1)=N(zk;Hkxk|k-1,Rk)
(2)
式中,Fk-1和Hk分別是狀態轉移矩陣和觀測矩陣;Qk-1和Rk分別是過程噪聲協方差和量測噪聲協方差。
這樣的濾波器稱為BET-GMPHD,下面給出k時刻BET-GMPHD的預測強度和更新強度。

Dk|k-1(xk|k-1|Zk-1)=Db,k(xk)+

(3)

如果k時刻BET-GMPHD的預測強度是




(4)

(5)
(6)
(7)
(8)
(9)

(10)

(11)
(12)

在這一節中,首先利用貝塔分布是二項分布的共軛先驗,得到估計檢測概率的貝葉斯迭代公式。然后,對包含檢測概率和運動參數的增廣狀態建立貝塔高斯分布模型,得到新的濾波器BG-ET-PHD。
二項分布的共軛先驗貝塔分布[25]是一個關于檢測概率pD的二值分布。

(13)



p(pD,k|k-1|Wk-1)=β(pD,k|k-1;pk|k-1,qk|k-1)
(14)
相應地,k時刻似然函數p(Wk|pD,k|k-1)是二項分布Bin(γk|k-1,|Wk|,pD,k|k-1),表達式為
p(Wk|pD,k|k-1)=p(|Wk||pD,k|k-1)=
C(γk|k-1,|Wk|)(pD,k|k-1)|Wk|(1-pD,k|k-1)γk|k-1-|Wk|
(15)

通過貝葉斯準則[25],檢測概率pD,k的后驗分布為
p(pD,k|Wk)∝p(pD,k|k-1|Wk-1)p(Wk|pD,k|k-1)=
β(pD,k;pk,qk)·Lβ
(16)
式中
pk=|Wk|+pk|k-1,qk=γk|k-1+qk|k-1-|Wk|

(17)



≈
β(pD,k|k-1;pk|k-1,qk|k-1)
(18)

μβ,k|k-1
(19)
(20)
假設1參考關于gamma預測分布假設[24,26],有理由相信假設1同樣合理。一般情況下,相比k-1時刻的后驗分布,變量k時刻的先驗分布均值不變,協方差增大[24]。為了確保pk|k-1>0和qk|k-1>0,預測協方差選擇為
,μβ,k-1(1-μβ,k-1))
一般,假設pD,k(xk)與擴展目標xk無關[27],即pD,k(xk)=pD,k。將pD,k作為增廣變量,令ξk=(pD,k,xk)表示擴展目標xk的增廣狀態,k時刻ξk服從貝塔高斯分布,即
p(ξk|Wk)=p(pD,k|Wk)p(xk|Wk)=
β(pD,k;pk,qk)N(xk;mk,Pk)=BG(ξk;ζk)
(21)
式中,p(pD,k|Wk)服從貝塔分布β(pD,k;pk,qk),在線性系統下p(xk|Wk)服從高斯分布N(xk;mk,Pk),ζk={pk,qk,mk,Pk}。

(1) 預測等式
定理1假設在k-1時刻BG-ET-PHD濾波器的更新強度是一個貝塔高斯的混合形式

(22)
那么k時刻BG-ET-PHD的預測強度為

(23)
式中


(24)
(25)
并且

證明對于k時刻的存活目標,相應的目標PHD強度[27]為

(26)
假設pS,k(ξk|k-1)與狀態ξk|k-1無關[27],即pS,k(ξk|k-1)=pS,k,p(ξk|k-1|ξk-1)為狀態轉移矩陣,可分解為
p(ξk|k-1|ξk-1)=p(pD,k|k-1|pD,k-1)p(xk|k-1|xk-1)
(27)
其中,p(xk|k-1|xk-1)=N(xk|k-1;Fk-1xk-1,Qk)。
將Dk-1(ξk-1|Wk-1)代入式(26),得到

(28)
式(28)分別使用了文獻[1]的結論



和假設1。
證畢
(2) 更新等式
定理2如果BG-ET-PHD濾波器的預測強度是一個貝塔高斯混合,表示為

(29)



(30)

(31)
(32)
(33)
(34)
(35)
(36)
(37)
(38)
(39)

(40)
(41)
(42)


(43)
(44)

證明將ξk取代xk,通過文獻[20]的推導過程,可以推導得到
Dk(ξk|Zk)=L(Zk|ξk|k-1)Dk|k-1(ξk|k-1|Zk-1)
(45)
式中

(46)
wρ見式(11);dW見式(12)。
(47)

將Dk|k-1(ξk|k-1|Zk-1)的貝塔高斯混合形式代入式(45)的無量測更新部分
(1-pD)γ(x)Dk|k-1(ξk|k-1|Zk-1)=

(48)

對于式(45)的量測更新部分,有


(49)
式中,φW還可以表示為
γk,|W|,pD,k)W
(50)
假設2

(51)



(52)

將Dk|k-1(ξk|k-1|Zk-1)的貝塔高斯混合形式代入φWDk|k-1(ξk|k-1|Zk-1),得到

(53)
式(53)主要采用了以下兩個分布乘積公式。首先,根據式(16)可以得到


(54)

然后,利用兩個高斯分布的乘積公式[1]得到




(55)

將式(53)代入式(49),得到式(30)的量測更新部分。
證畢




(56)
式中



圖1 場景1的目標真實軌跡Fig.1 True target trajectories for scenario 1


圖2 檢測概率的估計值Fig.2 Estimated value of detection probability

圖3 最大量測數目的估計值Fig.3 Estimated value of maximum measurement number


圖4 場景2的目標真實軌跡Fig.4 True target trajectories for scenario 2

圖不同的BET-PHD的目標數目估計Fig.5 Target number estimation of BET-PHD under different
采用最多同時存在10個目標的場景3評估兩種GG-ET-PHD和BG-ET-PHD濾波器性能,目標真實軌跡如圖6所示。其中,在第40~80 s,目標6和目標7相差10 m同向平行運動。

圖6 場景3的目標真實軌跡Fig.6 True target trajectories for scenario 3
圖7給出了兩種濾波器100次蒙特卡羅平均得到的目標數目估計和最佳子模式分配(optimal sub-pattern assignment, OSPA)距離誤差。可以看到,在平行運動的時間段(第40~80 s),相比Poisson模型的GG-ET-PHD濾波器,基于二項分布的BG-ET-PHD濾波器能夠獲得較為準確的目標數目估計和較小的OSPA誤差。這意味著基于二項分布模型的濾波器能更好地處理鄰近擴展目標的跟蹤問題。同時,還可以看出,BG-ET-PHD濾波器參數估計值與真實值之間的較小誤差沒有影響濾波器的跟蹤性能。此外,從圖7(b)可以看出,GG-ET-PHD和BG-ET-PHD濾波器在目標新生時刻獲得的OSPA誤差距離較大。這是因為這兩個濾波器延遲一個時間步長估計新生目標的數目和狀態,導致新生時刻新生目標的OSPA誤差較大。

圖7 兩個ET-PHD濾波器的跟蹤性能對比Fig.7 Tracking performance of two ET-PHD filters
針對基于二項分布的BET-PHD濾波算法中檢測概率和量測數目最大值在實際應用中未知的問題,本文提出基于貝塔高斯混合分量的BG-ET-PHD濾波器。該濾波器能夠自適應估計各個擴展目標的檢測概率和量測數目最大值。仿真結果表明BG-ET-PHD濾波器能夠準確估計未知參數,能獲得比基于泊松模型的GG-ET-PHD濾波器更好的跟蹤性能。