多目标跟踪的混合高斯概率假设密度滤波算法
吴盘龙,陈风
南京理工大学自动化学院,南京 (210094)
摘 要:为解决目标数未知或随时间变化时的多目标跟踪问题,将多目标状态和观测信
息表示为随机集的形式,建立了多目标跟踪的混合高斯概率假设密度(PHD)滤波方法。当
目标初始的先验概率密度满足高斯分布的形式时,通过将状态噪声、观测噪声、目标的繁衍、
新目标的产生、目标的存活概率和检测概率表示成混合高斯的形式,之后每个时刻的后验概
率密度均能表示成混合高斯的形式。线性混合高斯 PHD滤波方法将Kalman滤波引入到 PHD
滤波中,利用混合高斯成分预测和更新随机集的 PHD,并估计出目标的状态。实验结果表
明,在杂波环境下混合高斯 PHD 滤波方法可以有效地跟踪目标状态。
关键词:多目标跟踪;随机集;混合高斯;概率假设密度
中图分类号:TN911
多目标跟踪技术在军事和民用方面有着广泛应用,如弹道导弹防御,空中和海上监视,
空中交通管制等。在杂波环境中的多目标进行跟踪比较困难,当目标或杂波出现或消失时,
目标数及目标产生的观测值可能随时间改变。对多目标跟踪的一种有效方法是对各单个目标
分别滤波,这要求对单个目标及其观测值正确互联。但如果存在两目标轨迹相似或小范围内
目标数多时,再用数据关联进行多目标跟踪又要牵涉到失跟问题的处理[1]。
利用随机集方法可以解决目标数未知或随时间变化的情况。在多目标跟踪问题中,随机
集实际上就是元素及元素的个数都是随机变量的集合。当目标的数目未知或不断变化时,目
标数是一个离散随机变量,状态空间的维数也会随目标数的不同取值而变化。于是,多目标
的状态模型和观测模型可以表示为随机有限集形式。Mahler[2,3]提出的有限集统计理论给出
了随机集的统计特性,考虑到了不同维状态空间之间的比较,完善了多目标系统规范的贝叶
斯方法的有关内容,将单传感器单目标中的近似方法推广到多传感器多目标系统的研究中,
提出了概率假设密度(Probability Hypothesis Density, PHD),即目标状态后验密度的一阶矩,
实现对目标状态和目标数的估计。
本文就 Mahler 提出的这种统计多目标跟踪算法进行研究,讨论了这种有限集统计跟踪
算法的混合高斯实现问题,并将其应用到线性高斯模型下的多目标跟踪。实验结果表明,这
种方法可以解决计算过大的问题,且能很好的估计目标的状态及目标数。本文对水下三维情
形的被动 TMA 问题进行探讨,利用被动声纳的方位角、俯仰角和频率的量测信息为基础,
采用 UKF 算法对水下被动跟踪目标的状态进行估计,并与 EKF 算法进行了仿真比较。
1. PHD 多目标跟踪算法
在单目标跟踪系统中,目标的状态和观测值是两个向量,且向量的维数不随时间改变。
在多目标跟踪系统中,状态和观测值是各单个目标的状态和观测值的集合,它们的维数随时
间改变。在多目标运动的随机有限集模型中,状态集和观测集可以分别用状态空间的随机有
- 1 -
中国科技论文在线
限集 { },1 , ( ), ,k k k M kX x x= L 和观测空间的随机有限集 { },1 , ( ), ,k k k N kY y y= L 表示。其中,
( )M k 及 分别表示 时刻目标数及观测值数,由于部分观测值可能源于杂波所以
。
( )N k k
( ) ( )N k M k≥
随机有限集的概率假设密度与随机变量的期望类似,但集合之间无法定义加法运算,所
以随机集的期望没有意义[4]。有限集 X 在物理上可等价地表示为广义函数 xx X δ∈∑ ,其中
xδ 是中心在 x的狄拉克δ 函数。因而随机有限集Ξ可表示为随机密度 xx δ∈Ξ∑ 。随机有限
集 的一阶矩密度(或PHD)定义为Ξ [5]。
( ) ( ) ( ) ( )y yy yD x E x x f X Xδ δ Ξ∈Ξ ∈Ξ⎡ ⎤= =⎣ ⎦∑ ∑∫ δ (1)
在多目标跟踪问题中, 是在0( )D xΞ 0x 点期望的目标密度。 ( )S D x dx∫ 是出现在S内的
期望的目标数。而 的极值(图形的峰点)给出了目标的状态估计。 ( )D x
PHD滤波分为预测和更新两个步骤,PHD预测方程是:
| 1 , | 1 | 1 1( | ) ( ) ( ( ) ( | ) ( | )) ( )k k k s k k k k k kD x Y b x p x f x x D dζ β ζ ζ− − −= + +∫ ζ− (2)
( )kb x 表示 时刻新出现目标的RFS的强度函数,k | 1( | )k k xβ ζ− 表示由 时刻状态为1k −
ζ 的目标衍生的RFS的强度函数, , ( )s kp ζ 表示 1k − 时刻状态为ζ 的目标在 时刻仍存活的
概率,
k
| 1( | )k kf x ζ− 表示单个目标的转移概率密度。
假设目标相互独立运动,则PHD更新方程可以用下式近似计算:
( )
1
, | 1
, | 1
, | 1
( ) ( | ) ( )
( ) 1 ( ) ( )
( ) ( ) ( | ) ( )k
D k k k k
k D k k k
y Y k D k k k k
P x f y x D x
D x p x D x
k y P f y D dζ ζ ζ ζ+
−
−
∈ −
= − + +∑ ∫ (3)
其中, 表示 时刻噪声密度函数,( )kk y k ( ) ( )k k kk y c yλ= , 1(D kP x )+ 是检测概率, 1( | )kP y x +
是单个目标的似然函数, kλ 是单位面积内每帧扫描时的平均杂波数,服从Poisson分布,
是跟踪区域内杂波的概率密度。 ( )kc y
对更新的 PHD 积分,取最接近的整数值,就可以得到 时刻的目标数的期望值,即
。目标的位置可以从更新 PHD 的最大 个峰值点所在位置得到。
k
( )k kN D x d= ∫ x kN
2. 混合高斯 PHD 多目标跟踪
跟踪假设
对于线型高斯多目标模型,采用PHD滤波,需对目标的产生、衍生、和探测做以下假设[6]:
(1)每个目标都服从线性高斯的运动模型和测量模型:
( )| 1 1 1( | ) ; ,k k k kf x N x F Qζ ζ− −= − (4)
( ) ( )| ; ,k kg y x N y H x R=
P
k (5)
其中, 表示其密度均值为 ,协方差为 ;(:; , )N m p m 1kF − 为状态转移矩阵, 为
系统噪声协方差矩阵; 为观测矩阵,
1kQ −
kH kR 为测量噪声协方差矩阵。
(2)目标的存活概率和探测概率相互独立:
, ,( )S k S kP x P= , , ( ) ,D kP x PD k= (6)
- 2 -
中国科技论文在线
(3)新生目标随机集和衍生目标集概率密度都为高斯混合形态:
,
( ) ( ) ( )
, ,
1
( ) ( ; , )
b kJ
i i
k b k b k
i
b x w N x m P
=
= ∑ ,ib k
− −
(7)
,
( ) ( ) ( ) ( )
| 1 , , 1 , 1 , 1
1
( | ) ( ; , )
kJ
j j j j
k k k k k k
j
x w N x F m Q
β
β β β ββ ζ ζ− −
=
= +∑ (8)
其中: ,b kJ , ; , ,( ),ib kw ( ),ib km ( ),ib kP ,1, , b ki J= L 等参数决定了新生目标随机集的概率密
度, , , 为新生目标随机集概率密度的第 i个高斯成分的权值,均值和方差,( ),ib kw ( ),ib km ( ),ib kP
,b kJ 是 时刻新生目标高斯成分的个数。同理:k ,kJβ , ( ),jkwβ ; ( ), 1jkFβ − , , ,( ), 1jkdβ − ( ), 1jkQβ −
,1, , kj Jβ= L 等参数决定了由目标ζ 衍生的目标随机集的概率密度,衍生的目标一般都在
ζ 的附近。
混合高斯PHD滤波算法
混合高斯PHD滤波算法由以下六个步骤组成[7]:
第一步:初始化。用 个高斯成分的加权和对混合高斯PHD算法进行初始化 0J
0
( ) ( ) ( )
0|0 0 0 0
1
( ; , )
J
i i
i
iD w N x m P
=
=∑ (9)
其中 是均值为m,方差为 的高斯分布。则权值的和是期望的初始目标数,
。
( ; , )N x m P P
0
( )
0 0
1
ˆ
J
i
i
w T
=
=∑
第二步:预测
1| 1 , 1( ) ( ) ( )k k k S kD x b x D x+ + += + (10)
( ) ( ) ( )
1 , 1 , 1
1
( ) ( ; , )
bJ
i i
k b k b k
i
b x w N x m P+ + +
=
=∑ , 1ib k+
+ +
(11)
( ) ( ) ( )
, 1| , 1| , 1|
1
( ) ( ; , )
tJ
i i i
S k k s t s k k s k k
i
D x p w N x m P+
=
= ∑ (12)
其中, 是新产生目标的PHD, 是幸存目标的PHD;1( )kb x+ , 1( )S kD x+ ( ) ( ), 1|i is k k k km F m+ = ,
( ) ( )
, 1|
i i T
s k k k k k kP Q F P+ = + F
}
。
第三步:更新。当 的测量信息1k +
11 1,1 1,|
{ , ,
tk k k Y
Y y y ++ + += L 有效时,后验密度由下式
计算:
1|
1
( ) ( ) ( )
1| 1 1| 1 1| 1 1| 1
1
( ) (1 ) ( ) ( ) ( ; ( ), )
k k
k
J
i i
k k D k k k k k k k
y Y i
D x p D x w y N x m y P
+
+
+ + + + + + + +
∈ =
= − + ∑ ∑ i (13)
1|
( ) ( )
1| 1( )
1 ( ) ( )
1| 11
( )
( )
( ) ( )k k
i i
D k k ki
k J j j
D k t kj
p w q y
w y
c y p w q yλ +
+ +
+
+ +=
= + ∑ (14)
( ) ( ) ( )
1 1 1| 1 1 1| 1( ) ( ; , )
i i
k k k t k k kq y N y H m R H P H+ + + + + += + i Tk k+
i
i
)
(15)
( ) ( ) ( ) ( )
1| 1 1| 1 1 1|( ) ( )
i i i
k k k k k k k km y m K y H m+ + + + + += + − (16)
( ) ( ) ( )
1| 1 1 1 1|
i i
k k k k k kP I K H P+ + + + +⎡ ⎤= −⎣ ⎦ (17)
( ) ( ) ( ) 1
1 1| 1 1 1| 1 1(
i i T i T
k k k k k k k k kK P H H P H R
−
+ + + + + + += + (18)
- 3 -
中国科技论文在线
第四步:剪枝。将小于权值τ 的高斯成分滤除掉。
1
1
1
( )
1 ( ) ( ) ( )1
1| 1 1 1 1( )
111
( ; , )
t
t
t
P
P
J l J
k i il
k k k k kJ j
i Nkj N
w
D w N
w
+ +
+
+=
+ + + + +
= ++= +
= ∑ ∑∑
ix m P
1k+
(19)
其中 是大于阈值的权值。 ( )1, 1 ,ik pw i N J+ = + L
第五步:合并。当一些高斯成分间的距离非常接近时(小于阈值U ),可以将这些高
斯成分进行合并
设置 , 且 0l = ( )1 1{ 1, , | }ik kI i J w τ+ += = L >
}
。重复以下步骤:
1l l= +
( )
1arg max
i
ti I
j w +∈=
( ) ( ) ( ) 1 ( ) ( )
1 1 1 1 1{ | ( ) ( ) ( )
i j T i i j
k k k k kL i I m m P m m U
−
+ + + + += ∈ − − ≤
( ) ( )
1 1
l i
k k
i L
w w+ +
∈
=∑%
( ) ( ) ( )
1 1( )
1
1l i
k kl
i Ik
m w
w+ +∈+
= ∑% % 1ikm +
( ) ( ) ( ) ( ) ( ) ( ) ( )
1 1 1 1 1 1( )
1
1 ( ( )( )l i i l i lk k k k k kl
i Lk
P w P m m m m
w+ + + + + +∈+
= + − −∑% % %% 1 )i Tk+
\I I L=
直至 I φ= 。
如果 ,用 个具有最大权值的高斯成分替代maxl J> maxJ ( ) ( ) ( )1 1 1{ , , }i i i lk k k iw m P 1+ + + =%% % 。否则
即为合并后的高斯成分。 ( ) ( ) ( )1 1 1{ , , }i i i lk k k iw m P+ + +%% % 1=
第六步:状态估计。提取权值大于 1τ 的高斯成分,则
( ) ( )
1 1 1
ˆ { :i ik k kX m w 1}τ+ + += >
= × S
(20)
3. 仿真实验
设 PHD 目标跟踪范围为 , 内有三个目标,每个目标的状
态包括位置和速度。状态变量
[01000] [01000]S m m
[ ], , , Tk k k k kX x y x y= & & 。目标运动符合线性高斯模型,目标状
态方程为:
2
2
1
1 0 0 / 2 0
0 1 0 0 / 2
0 0 1 0 0
0 0 0 1 0
k k
T T
T T
X X
T
T
kω+
⎡ ⎤⎡ ⎤ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥= + ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦
(21)
目标测量方程为:
1 0 0 0
0 1 0 0k k k
y X v⎡ ⎤= ⎢ ⎥⎣ ⎦ +
⎤⎦
(22)
其中,T 为采样时间。 ,为目标加速度引起的系统过程噪声。过程噪声
为零均值高斯白噪声,其方差矩阵为 。
, ,,
T
k x k y kω ω ω⎡= ⎣
2 2([ , ])
x y
Q diag ω ωσ σ= , ,, Tk x k y kv v v⎡ ⎤= ⎣ ⎦ 为测量噪声,
其方差矩阵为R, 。 2 2([ , ])
x yv v
R diag σ σ=
- 4 -
中国科技论文在线
仿真条件和相关参数为: , , ,检测概率
,生存概率 ,裁剪阈值
1T s= 2 2
x yω ωσ σ= = 2 2 yv vσ σ= =
= = 1 5eτ = − ,合并阈值 4U = ,状态估计阈值
1 τ = 。 , ,65 10k mλ − −= × 2 6 210kc m= max 20J = 。新目标的出现符合Poisson过程,
, ,
。
( ) ( , , )k b bb x N x m P= [500 ,400 ,10 / ,10 / ]m m m m s m s=
2(10 [10,10,1,1])P diag −= ×
b
b
0 5 10 15 20 25 30 35 40 45 50
0
1
x/
km
t/s
0 5 10 15 20 25 30 35 40 45 50
0
1
y/
km
t/s
图1 目标轨迹的测量值和真实值 图2 目标轨迹估计值和真实值
仿真结果如图1至图3所示。图1为目标轨迹的测
量值,图2是混合高斯PHD估计的目标轨迹,图3是
目标数目的估计结果。从实验结果可以看出,对目
标数变化的多目标情况,混合高斯PHD滤波对目标
状态的点估计效果很好,且基本能正确估计目标数。
0 5 10 15 20 25 30 35 40 45 50
0
1
2
3
t/s
Ta
rg
et
n
um
be
r
True target number
Estimated number
4 结论
多目标跟踪问题中,当目标数未知或随时间变
化时,PHD滤波可对杂波环境下的目标状态和目标
数同时进行估计。本文研究了用混合高斯PHD滤波
实现跟踪随机集的方法。仿真实验表明,在杂波环
境下,混合高斯PHD滤波可以跟踪线型高斯模型下
的目标状态和目标数,跟踪效果较好。
图3 目标数目估计
参考文献
[1] Shin J, Guibas L J, Zhao F. A Distributed Algorithm for Managing Multi-Target Identities in
Wireless Ad-hoc Sensor Networks[C]// Proceeding of 2nd Workshop on Information
Processing in Sensor Networks, Palo Alto, CA, 2003.
[2] Goodman I, Mahler R, Nguyen H. Mathematics of Data Fusion[M]. Kluwer Academic
Publishers,1997.
[3] Mahler R. Multi-target Bayes Filtering via First-order Multi-target Moments[J]. IEEE Transactions on
- 5 -
中国科技论文在线
Aerospace and Electronic Systems, 2003, 39: 1152-1178.
[4] Daley D, Vere-Jones D. An Introduction to the Theory of Point Processes[M]. Springer-Verlag, 1988.
[5] Vo B, Singh S, Doucet A. Sequential Monte Carlo Implementation of the PHD Filter for
Multi-target[C]// Proc. International Conference on Information Fusion, Cairns, Australia,
2003:792-799.
[6] Vo B, Ma W K. A Closed Form Solution for the Probability Hypothesis Density Filter[J]. Proc
FUSION 2005, 2005,2:1-8.
[7] Vo B, Ma W K. The Gaussian Mixture Probability Hypothesis Density Filter[J]. IEEE
Transactions on Signal Processing ,2006,54(11):4091-4014.
Gaussian Mixture Probability Hypothesis Density Filter for Multiple
Target Tracking
WU Pan-long, CHEN Feng
School of Automation, Nanjing University of Science and Technology, Nanjing (210094)
Abstract
When the number of targets is unknown or varied with time, the target state and measurements
can be represented as random sets. The Gaussian mixture probability hypothesis density(PHD)
filter is implemented to track the multi-targets. The analytical analysis of the method show that the
posterior intensity at any subsequent time step remains a Gaussian mixture under the assumption
that the state noise, the measurement noise, target spawn intensity, new birth intensity, target
survival probability, and detection probability are all Gaussian mixture. The Kalman filter is
embedded in the Gaussian mixture PHD filter. This method uses Gaussian components to predict
and update the PHD of random sets, and estimates targets states. Experiments show that the
Gaussian mixture PHD filter can be used to track multi-target in clutter effectively.
Keywords: Multi-target tracking; Random set; Gaussian Mixture; Probability hypothesis density
- 6 -
中国科技论文在线
1. PHD多目标跟踪算法
2. 混合高斯PHD多目标跟踪
3. 仿真实验
4 结论