變分貝葉斯概率數(shù)據(jù)關(guān)聯(lián)算法

惲鵬 吳盤龍 李星秀 何山

目標(biāo)跟蹤在軍事和民用領(lǐng)域均具有廣泛的應(yīng)用[1?3],而在現(xiàn)實工程中,例如,在海上目標(biāo)跟蹤[4]和空中目標(biāo)跟蹤[5]等領(lǐng)域,觀測設(shè)備在獲取與目標(biāo)相關(guān)的真實量測同時會獲取與目標(biāo)不相關(guān)的虛警數(shù)據(jù).虛警數(shù)據(jù)也被稱為雜波[6],由于無法在混有雜波信息的數(shù)據(jù)集中直接識別真實量測,基于真實量測設(shè)計的濾波器[7]無法對目標(biāo)進行跟蹤.此外密集的雜波信號會影響真實目標(biāo)的回波信號,導(dǎo)致雷達探測能力下降進而影響目標(biāo)的跟蹤精度.因此如何提高雜波環(huán)境下目標(biāo)跟蹤精度是當(dāng)前目標(biāo)跟蹤領(lǐng)域亟需解決的問題之一.

許多學(xué)者均對雜波環(huán)境下的目標(biāo)跟蹤問題展開研究,并取得了不錯的研究成果.對于單目標(biāo)跟蹤系統(tǒng),最近鄰算法[8?10]以及概率數(shù)據(jù)關(guān)聯(lián)算法(Probabilistic data association algorithm,PDA)[11?14]是較為經(jīng)典的解決方法.最近鄰算法通過構(gòu)建距離指標(biāo)函數(shù)選取與預(yù)測中心最近的數(shù)據(jù)作為真實量測,然而當(dāng)雜波距離預(yù)測中心較近時,其往往會錯誤地將雜波選取為真實量測,導(dǎo)致跟蹤精度下降.PDA算法本質(zhì)上是一種全鄰算法,其認(rèn)為關(guān)聯(lián)波門內(nèi)的每一個數(shù)據(jù)均可作為真實量測,并通過新息函數(shù)為每一個數(shù)據(jù)賦予相應(yīng)的權(quán)重.由于該方法相比于最近鄰算法可降低雜波對狀態(tài)更新的影響,因此具有更高的目標(biāo)跟蹤精度.PDA 算法的設(shè)計思想同樣被應(yīng)用于解決雜波環(huán)境下機動目標(biāo)跟蹤問題,相關(guān)的算法包括交互多模型概率數(shù)據(jù)關(guān)聯(lián)算法[15](Interactive multi-model probabilistic data association,IMM-PDA)、聯(lián)合交互多模型概率數(shù)據(jù)關(guān)聯(lián)算法[16]以及聯(lián)合交互多模型距離加權(quán)概率數(shù)據(jù)關(guān)聯(lián)算法[17]等.

從貝葉斯統(tǒng)計原理的角度來看,雜波環(huán)境下目標(biāo)真實狀態(tài)后驗概率密度函數(shù)應(yīng)為真實量測作為觀測信息情況下的狀態(tài)后驗概率密度函數(shù).由于無法直接識別真實量測,因此目標(biāo)真實狀態(tài)后驗概率密度函數(shù)無法獲取并且每一個數(shù)據(jù)作為真實量測而獲取的狀態(tài)后驗概率密度函數(shù)均可能為真實狀態(tài)后驗概率密度函數(shù).此外,當(dāng)檢測概率不為1 時,所有數(shù)據(jù)與目標(biāo)不相關(guān)時獲取的狀態(tài)后驗概率密度函數(shù)也可能成為真實狀態(tài)后驗概率密度函數(shù).代理概率密度函數(shù)表示某一變量可能的概率分布,所有數(shù)據(jù)關(guān)聯(lián)形式下獲取的狀態(tài)后驗概率密度函數(shù)可構(gòu)成一個狀態(tài)代理概率密度函數(shù)集合.基于該集合的統(tǒng)計特征以實現(xiàn)狀態(tài)近似后驗概率密度函數(shù)的獲取是一個合理的選擇.PDA 算法通過提取每一個數(shù)據(jù)?狀態(tài)后驗概率密度函數(shù)的均值與協(xié)方差并基于加權(quán)融合的方式以實現(xiàn)狀態(tài)近似后驗概率密度函數(shù)的獲取.這類近似方法雖然簡單且易于實現(xiàn),但是在數(shù)據(jù)數(shù)目較多的情形下會存在計算效率低的缺陷.本文在變分貝葉斯框架[18]下設(shè)計了一種新的概率數(shù)據(jù)關(guān)聯(lián)算法(Variational Bayesian based probabilistic data association algorithm,VB-PDA),該算法首先將關(guān)聯(lián)事件視為隨機變量并利用多項分布對兩者進行建模,隨后基于量測集、目標(biāo)狀態(tài)、關(guān)聯(lián)事件的聯(lián)合概率密度函數(shù)求取關(guān)聯(lián)事件的后驗概率密度函數(shù),最后將關(guān)聯(lián)事件的后驗概率密度函數(shù)引入變分貝葉斯框架中以獲取狀態(tài)近似后驗概率密度函數(shù).VB-PDA算法本質(zhì)上是一種在PDA 建?？蚣芟峦ㄟ^權(quán)重Kullback-Leibler (KL)平均[19]準(zhǔn)則來實現(xiàn)狀態(tài)近似后驗概率密度函數(shù)獲取的方法.相比于PDA 算法,VB-PDA 算法在提高計算效率的同時,獲取的近似后驗概率密度函數(shù)也更能體現(xiàn)狀態(tài)代理概率密度函數(shù)集的統(tǒng)計特征.相關(guān)仿真實驗結(jié)果驗證了算法的有效性.

本文內(nèi)容安排如下: 第1 節(jié)對相關(guān)問題進行了描述;第2 節(jié)對PDA 算法進行了介紹以及對VBPDA 算法進行了推導(dǎo);第3 節(jié)從權(quán)重KL 平均角度對VB-PDA 算法進行了分析;第4 節(jié)對算法計算復(fù)雜度進行了分析;第5 節(jié)通過相關(guān)仿真實驗對提出的算法進行驗證;第6 節(jié)對本文進行了總結(jié).

1 問題描述

假設(shè)線性離散系統(tǒng)狀態(tài)空間模型的狀態(tài)方程和量測方程分別為:

由于無法直接區(qū)分真實量測與雜波,狀態(tài)真實后驗概率密度函數(shù)無法獲取.如何獲取合適的狀態(tài)近似后驗概率密度函數(shù)是本文的研究目的.

2 算法介紹

不同的雜波分布會導(dǎo)致算法具備不同的數(shù)學(xué)形式.假設(shè)目標(biāo)的檢測概率為pd,門概率為pg.系統(tǒng)視場 Θ 被視為關(guān)聯(lián)門,雜波在 Θ 內(nèi)服從均勻分布,該區(qū)域體積為Vk+1.雜波數(shù)目服從泊松分布,雜波密度為λk+1.Vk+1、pd和λk+1被認(rèn)為是已知的.pg表示真實量測在視場內(nèi)的概率.其可由如下公式獲取:

式中,R 表示實數(shù)域.本文介紹的PDA 算法以及設(shè)計的VB-PDA 算法均是基于上述雜波環(huán)境給出的.

2.1 PDA 算法

2.2 VB-PDA 算法

在VB-PDA 算法中,量測集中任意一個數(shù)據(jù)為真實量測的概率相同被視為先驗信息.此外,pd、pg以及雜波相關(guān)信息同樣被視為已知的先驗信息并被使用構(gòu)建相應(yīng)的概率密度函數(shù).將關(guān)聯(lián)事件視為隨機變量,由于關(guān)聯(lián)事件的數(shù)目是有限的,因此可用多項分布進行建模.此外,關(guān)聯(lián)事件數(shù)目的確定意味著量測數(shù)目同樣確定.因此基于已有的先驗信息,關(guān)聯(lián)事件的先驗概率密度函數(shù)可被建模為:

可以發(fā)現(xiàn),基于本文建?？蚣芟峦茖?dǎo)的關(guān)聯(lián)事件權(quán)重模型與PDA 算法相同.利用VB 算法[18]對狀態(tài)近似后驗密度函數(shù)qVB(xk+1) 進行求解,其求解公式如下:

基于推導(dǎo)過程可以發(fā)現(xiàn),狀態(tài)后驗概率密度函數(shù)在貝葉斯框架下可直接被獲取,并且該概率密度函數(shù)服從高斯混合分布.然而隨著狀態(tài)更新的過程不斷進行,該后驗概率密度函數(shù)的高斯分量會呈幾何式增長,導(dǎo)致其計算復(fù)雜度過大.考慮到實時性的影響,采用VB 算法獲取狀態(tài)近似后驗概率密度函數(shù).VB-PDA 算法獲取的狀態(tài)近似后驗概率密度函數(shù)服從高斯分布,這樣的結(jié)果會為算法在后續(xù)的狀態(tài)更新過程中提供便利.此外,VB-PDA 算法獲取的狀態(tài)協(xié)方差小于狀態(tài)預(yù)測協(xié)方差以及PDA 算法獲取的狀態(tài)協(xié)方差,因此該算法可有效避免奇異現(xiàn)象的產(chǎn)生.

3 VB-PDA 算法討論

本節(jié)從代理概率密度函數(shù)集的角度說明VBPDA 算法的優(yōu)勢.權(quán)重KL 平均被認(rèn)為是最能體現(xiàn)代理概率密度函數(shù)集統(tǒng)計特征的概率密度函數(shù)[19].權(quán)重KL 平均(Ψ) 定義為:

式中,L (·) 表示下界函數(shù),VB 算法是一種基于平均場理論實現(xiàn)參數(shù)近似后驗概率密度函數(shù)獲取的方法,qVB(ΨVB) 可被分解為:

根據(jù)式(25)、式(27)和式(39),式(37)的等價形式為:

根據(jù)式(36)和式(42)可知,VB-PDA 算法獲取的參數(shù)集近似后驗概率密度函數(shù)即為PDA 算法建?？蚣芟碌臓顟B(tài)代理概率密度函數(shù)集的權(quán)重KL平均,因此該算法獲取的狀態(tài)近似后驗概率密度函數(shù)比PDA 算法獲取的狀態(tài)近似后驗概率密度函數(shù)更能體現(xiàn)代理概率密度函數(shù)集統(tǒng)計特征.

4 計算復(fù)雜度分析

算法的計算復(fù)雜度同樣是影響算法跟蹤性能的重要指標(biāo).本文將加減運算與乘除運算次數(shù)作為指標(biāo)對4種算法的計算復(fù)雜度進行分析.表1 展示了PDA算法、距離加權(quán)概率數(shù)據(jù)關(guān)聯(lián)算法[17](Distance weighted probabilistic data association,DWPDA)以及VB-PDA 算法在一步狀態(tài)更新過程中所需的加減運算與乘除運算次數(shù).可以看出,VBPDA 算法比PDA 算法以及DW-PDA 算法擁有更低的算法復(fù)雜度.當(dāng)雜波數(shù)目較大多時,VB-PDA算法計算效率優(yōu)勢更加明顯.

表1 一步狀態(tài)更新過程中所需的加減運算與乘除運算次數(shù)Table 1 The number of addition and subtraction operations and multiplication and division operations required in the process of one-step state update

5 仿真實驗

5.1 仿真場景 1

利用目標(biāo)跟蹤案例對算法的跟蹤性能進行測試.假設(shè)目標(biāo)在二維平面進行勻速直線運動,利用勻速模型(Constant velocity,CV)對該目標(biāo)進行跟蹤并且量測設(shè)備可以獲取目標(biāo)的位置信息.因此,目標(biāo)的狀態(tài)信息xk,狀態(tài)轉(zhuǎn)移矩陣Fk,量測矩陣Hk+1分別為:

式中,ak和bk表示目標(biāo)位置信息,和表示目標(biāo)速度信息,t=1 s 表示采樣時間.系統(tǒng)過程噪聲協(xié)方差和系統(tǒng)噪聲量測協(xié)方差為:

目標(biāo)真實初始狀態(tài)信息x0為[100 m,100 m,1 m/s,1 m/s],初始狀態(tài)協(xié)方差P0|0為diag{100 m2,100 m2,1 m2/s2,1 m2/s2},目標(biāo)初始狀態(tài)信息x0|0為x0與協(xié)方差為P0|0高斯白噪聲的和[20],該高斯白噪聲由Matlab 程序隨機產(chǎn)生.因此在每一次蒙特卡洛仿真實驗中,目標(biāo)初始狀態(tài)信息是不同的.

在仿真實驗中,將均方根誤差(Root mean square error,RMSE)和平均時間均方根誤差(Time average root mean square error,TRMSE)作為指標(biāo)來比較各算法的估計性能,目標(biāo)位置和速度的RMSE和TRMSE 計算公式如下:

圖2 仿真場景1 下 3種算法的速度 RMSEFig.2 The RMSEof velocity from three algorithms in scenario 1

假設(shè)雜波分布在 [ 0～400 m,0～400 m] 區(qū)域里服從均勻分布,因此Vk+1= 160 000 m2.雜波數(shù)目服從泊松分布[21].檢測概率為0.95,雜波數(shù)目統(tǒng)計期望為20.圖1～2 展示了3種算法的位置與速度RMSE,表2 展示了3種算法的位置與速度TRMSE.由表2 可以看出,VB-PDA 算法的位置估計精度在3種算法中最好,而DW-PDA 算法的速度估計精度在3種算法中最好.從TRMSE 數(shù)據(jù)可以發(fā)現(xiàn),相比于PDA 算法,VB-PDA 算法在位置和速度估計精度上分別提高了2.18%和0.11%.DW-PDA算法是在權(quán)重計算的過程中對PDA 算法的改進,而VB-PDA 算法是從狀態(tài)后驗概率密度函數(shù)獲取的過程中對PDA 算法的改進.由于兩種算法的改進方式不同,出現(xiàn)DW-PDA 算法在速度估計精度優(yōu)于VB-PDA 算法這一現(xiàn)象屬于正常情形.表3展示了4種算法一次蒙特卡洛實驗所需的計算時間,可以看出VB-PDA 算法的所需的計算時間最少,這也驗證了第4 節(jié)分析的結(jié)果.從數(shù)據(jù)中可以發(fā)現(xiàn),相比于PDA 算法,VB-PDA 算法在計算效率上提高了17.37%.

圖1 場景 1 下 3種算法的位置 RMSEFig.1 The RMSE of position from three algorithms in scenario 1

表2 場景 1 下 3種算法的 TRMSETable 2 The TRMSE of three algorithms in scenario 1

表3 場景 1 下一次蒙特卡洛仿真實驗所需的計算時間Table 3 The computational time at one Monte Carlo simulation experiment in scenario 1

5.2 仿真場景 2

考慮基于機動目標(biāo)跟蹤案例對提出算法的跟蹤性能進行測試,利用交互多模型(Interactive multimodel,IMM)算法[22]與上述3種算法結(jié)合對目標(biāo)進行跟蹤.IMM 算法中采用兩個帶不同轉(zhuǎn)彎速率勻速轉(zhuǎn)彎模型(CT)對該目標(biāo)進行跟蹤并且量測設(shè)備可以獲取目標(biāo)的位置信息.因此,目標(biāo)的狀態(tài)信息xk,子模型狀態(tài)轉(zhuǎn)移矩陣Fk,j(j=1,2),量測矩陣Hk+1,j(j=1,2) 分別為:

其中,子模型的轉(zhuǎn)彎速率分別為w1=?π/20 rad/s,w2=π/20 rad/s. 系統(tǒng)過程噪聲協(xié)方差Qk,j和量測噪聲協(xié)方差Rk+1,j設(shè)置為:

目標(biāo)真實初始狀態(tài)信息x0為 [100 m,100 m,3 m/s,3 m/s],目標(biāo)在1～50 s 做轉(zhuǎn)彎速率為π/100 rad/s的勻速轉(zhuǎn)彎運動,在51～100 s 做轉(zhuǎn)彎速率為?π/100 rad/s 的勻速轉(zhuǎn)彎運動.假設(shè)2 個子模型初始權(quán)重均為0.5,馬爾科夫轉(zhuǎn)移矩陣為子模型初始狀態(tài)協(xié)方差為diag{100 m2,100 m2,1 m2/s2,1m2/s2}.為驗證算法的有效性,對IMM-PDA 算法、IMM-DW-PDA 算法和IMM-VB-PDA 算法進行測試.

假設(shè)雜波分布在 [ 0～400 m,0～400 m] 區(qū)域里服從均勻分布,因此Vk+1= 160 000 m2.雜波數(shù)目服從泊松分布.檢測概率為0.9,雜波數(shù)目統(tǒng)計期望為10.圖3～4 給出了3種算法獲取的位置與速度RMSE.表4 給出了3種算法獲取的位置與速度TRMSE.由表4 可以看出,在此類場景下,IMMVB-PDA 算法的估計精度最佳.從TRMSE 數(shù)據(jù)可以發(fā)現(xiàn),相比于IMM-PDA 算法,IMM-VB-PDA算法在位置和速度估計精度上分別提高了5.26%和4.8%.

表4 仿真場景2 下3種算法的TRMSETable 4 The TRMSE of three algorithms in scenario 2

圖3 仿真場景 2 下 3種算法的位置 RMSEFig.3 The RMSE of position from three algorithms in scenario 2

6 結(jié)束語

圖4 仿真場景 2 下 3種算法的速度 RMSEFig.4 The RMSE of velocity from three algorithms in scenario 2

PDA 算法采用從狀態(tài)代理概率密度函數(shù)集里提取均值信息與協(xié)方差信息的方式獲取狀態(tài)近似后驗概率密度函數(shù).該方法雖然便捷,但在數(shù)據(jù)量較大情形下存在計算效率低的缺陷.針對這一問題,本文提出了一種VB-PDA 算法,該算法本質(zhì)上是在PDA 算法建?？蚣芟绿崛∠鄳?yīng)的權(quán)重KL 平均以獲取狀態(tài)近似后驗概率密度函數(shù).相比PDA 算法,VB-PDA 算法在提高計算效率的同時,獲取了近似程度更高的狀態(tài)后驗概率密度函數(shù).相關(guān)仿真實驗結(jié)果驗證了VB-PDA 算法的有效性.由于本文提出的算法是針對單目標(biāo)跟蹤問題設(shè)計的,如何將類似的思想拓展到多目標(biāo)跟蹤領(lǐng)域是未來的研究目標(biāo)之一.