多机动目标跟踪的IMM-GMPHD滤波算法
时间：2022-10-23 01:26:32
文章编号：1003-6199(2011)04-0089-06
摘 要：针对现有多机动目标跟踪算法精度低、计算量大、约束条件苛刻等问题，本文将高斯混合概率假设密度（Gaussian Mixture PHD，GMPHD）滤波器和交互式多模型（Interacting Multiple Model，IMM）相结合，提出交互式多模型GM-PHD（Interacting Multiple Model GMPHD，IMMGMPHD）滤波算法。算法不仅避免了多目标跟踪中的数据关联问题，而且在漏检、目标密集、目标机动、航迹交叉、目标数目未知的杂波环境下能够稳定、精确地估计目标数目和状态。100次蒙特卡洛（Monte Carlo，MC）仿真结果表明，IMMGMPHD滤波器能在不增加额外计算负担的基础上，体现出较高的精确度和较强的鲁棒性。
关键词：多机动目标跟踪;高斯混合概率假设密度;交互式多模型;蒙特卡洛
中图分类号： TN953;TP302 文献标识码：A
IMMGMPHD Filter for Multiple Maneuvering Targets Tracking
YAN Kang, YAN Yude
(Nanjing University of Science and Technology, Nanjing 210094, China)
Abstract:Considering the traditional data association algorithm of multiple maneuvering targets tracking being of hard constraint condition, lower estimated accuracy, and higher computational complexity. In this paper an interacting multiplemodel (IMM) implementation of the GMPHD filter is proposed to estimate the states of multiple targets which may have different maneuverability in measurements system. The new algorithm avoids the difficult problem of data association. It is able to deal with clutter, miss detection, false alarm, dense, and cross targets tracking effectively. 100 Monte Carlo (MC) simulation results show that the proposed MMGMPHD filter significantly outperforms in estimating the number and states of the multiple maneuvering targets. The proposed algorithm has higher tracking accuracy and more steady tracking performance.
Key words:multiple maneuvering targets tracking;gaussian mixture;probability hypothesis density;interacting multiple model;monte carlo
1 引 言
多目标跟踪技术在军事和民用领域都有着广泛的应用，如空防、精确制导，以及汽车防碰撞、机器人视觉等。在多目标跟踪场景中，目标的出现和消亡，目标间的交叉、编队、迂回等各种协同和非协同机动，还有大量欺骗性干扰和复杂电磁环境和强杂波背景的复合干扰都对目标跟踪带来大量的不确定性。
随着多目标跟踪理论研究的日益深入，形成了一些成熟的多目标跟踪方法，概括来讲可分两大类［1］：①极大似然类滤波算法②贝叶斯类滤波算法。这两类方法都需要数据互联理论和滤波技术共同完成，目标数目的不确定性让传统的状态空间方法难以处理，且随着目标量测的增长，数据互联的计算复杂度也呈指数级增长。近年来，基于随机有限集（Random Finite Set，RFS）［2］理论的多目标跟踪算法受到了国内外学者的广泛关注。Mahler提出概率假设密度（Probability Hypothesis Density，PHD）［3］的概念。Vo等给出了PHD粒子滤波的实现（Sequential Monte Carlo PHD，SMC-PHD）［4］滤波器。虽然PHD滤波器可以联合估计目标个数和目标状态，但目标状态的提取需要专门的算法获得。Vo等在线性高斯混合（Gaussian Mixture，GM）条件下，给出了有闭式解的递推方程即高斯混合概率假设密度（Gaussian Mixture PHD，GMPHD）［5］滤波器，不需要特别的聚类算法，而且简化了运算。以上提到的基于随机集合理论的跟踪算法都是针对整个多目标状态随机集合的运算，并不区分个体目标自身的特性，当目标发生机动时，并不能识别每个目标的机动情况。而交互式多模型（Interacting Multiple Model，IMM）［6］算法是跟踪单个机动目标的有效方法，它具有明显的并行结构，能有效地实现多模型的自动切换，在计算精度和复杂度上可获得比较好的平衡。
针对上述问题，本文将IMM算法和GMPHD滤波器相结合，提出了IMM-GMPHD滤波算法。算法能快速跟踪数目变化的多机动目标，完成对目标数目和状态的准确估计。本文最后给出仿真结果，并进行了分析。
2 多目标跟踪模型
2.1 随机有限集多目标模型
在单目标跟踪系统中，目标的状态和量测值是两个向量，且向量维数不随时间改变。在多目标跟踪系统中，状态和量测值是各单个目标的状态和量测值的集合，集合的维数随时间改变。记在时刻k，监控区域内的目标个数为Nk，单个目标的状态为xk，则多目标的状态可以表示为一个RFS Xk=xnknxkn=1，其中xnk表示第n(n=1,…,Nk)个目标的状态向量；量测RFS表示为Zk=znknzkn=1，znk表示第n(n=1,…,Nk)个目标的量测。
每个目标都服从线性高斯的运动模型和量测模型
fk|k－1(x|ζ)=N(x;Fk－1ζ,Qk－1) （1）
gk(z|x)=N(z;Hkx,Rk) （2）
其中，N:;m,p表示其密度均值为m，协方差为P；Fk－1为状态转移矩阵，Qk－1为模型噪声协方差矩阵；Hk为量测矩阵，Rk为量测噪声协方差矩阵。
2.2 交互式多模型原理
交互式多模型算法是目前机动目标跟踪中比较优越的混合式估计算法，在合理给定运动模型集合的前提下，IMM可以跟踪目标的任意机动。IMM的贝叶斯递归过程可表示为：
1）初始化模型混合概率
μrk－1|z1:k－1(rk－1|rk,z1:k－1) （3）
2）初始化状态概率密度
pxk－1|rk－1,z1:k－1(xk－1|rk,z1:k－1)（4）
3）模型概率预测
rk－1|rk,z1:k－1μrk|z1:k－1（5）
4）模型r下的状态概率密度预测
xk－1|rk,z1:k－1pxk|rk,z1:k－1 （6）
5）模型r的概率更新
μrk|z1:k－1μrk|z1:k （7）
6）模型r下的状态概率密度更新
pxk|rk,z1:k－1pxk|rk,z1:k（8）
7）估计融合
pxk|z1:k=∑rkp(xk|rk,z1:k)μ(rk|z1:k) （9）
其中rk∈1,…,Nr表示k时刻模型的指示参数，Nr表示总的模型个数；xk和zk分别表示k时刻目标的状态向量和观测向量；式（9）描述了IMM算法的实质，即并行处理每个模型下的目标跟踪并计算该模型当前概率，最后按模型的概率权重对各模型下的估计结果进行融合。
3 IMMGMPHD滤波
为了获得IMMGMPHD滤波器，我们首先对PHD滤波器作简要介绍。
3.1 PHD滤波器
因为随机有限集的条件期望［7］是不存在的，那么RFS的一阶矩也就不存在，而随机点过程中的一阶矩是可以定义的，它表示空间中随机点数的期望值，一个自然的想法是借用它作为RFS一阶矩的定义。有限集X可等价地表示为广义函数∑x∈Xδx，其中δx是中心在x的狄拉克δ函数。因而随机有限集Ξ可表示为随机密度∑x∈Ξδx。则随机有限集Ξ的一阶矩密度即PHD可表示为：
υ(x)=E∑y∈Ξδy(x)
=∫∑y∈Ξδy(x)fΞ(x)δX（10）
υ(x)类似于常增益Kalman滤波器，通过递推预测和更新PHD获取目标个数的分布信息，进而获得目标状态的分布信息［8］。
3.2 IMMGMPHD滤波递推
假设任意时刻每个模型下的GMPHD包含相同数量的高斯分量，且每个模型的每个高斯分量具有一一对应关系，分别计算每个模型下GMPHD以及每个高斯分量下模型的条件概率，最终将PHD融合，这样即可建立具有封闭解的IMMGMPHD算法的递归模型。
3.2.1 模型条件初始化
1）模型混合概率初始化
计算多目标状态随机集合的混合转移概率
μ(l)(i,j)k－1|k－1Jk－1l=1：μ(l)(i,j)k－1|k－1=1jπijμ(l)(i)k－1（11）
其中μ(l)(i,j)k－1|k－1表示对于任意高斯分量l，模型i到模型j的混合转移概率，式中j满足：
j=∑ri=1πijμ(l)(i)k－1 （12）
πij为模型i到模型j的先验转移概率。
2）混合条件PHD的初始化
初始化的混合条件PHD可表示为：
(j)k－1|k－1(x)=∑Jk－1l=1(l)(j)kN(x;(l)(j)k－1,(l)(j)k－1) （13）
其中任意高斯分量l在k－1时刻的均值估计满足：
(l)(j)k－1=∑ri=1m(l)(j)k－1|k－1(l)(i,j)k－1|k－1 （14）
高斯分量l的权值和协方差矩阵满足：
(l)(j)k－1=∑ri=1ω(l)(i)k－1(l)(i,j)k－1 （15）
(l)(j)k－1|k－1=∑ri=1(l)(i)k－1|k－1+(m(l)(i)k－1|k－1－
(l)(j)k－1)(m(l)(i)k－1|k－1－(l)(j)k－1)T×μ(l)(i,j)k－1 (16)
3.2.2 条件PHD的预测和更新
1）条件PHD预测
υ(j)k|k－1(x)=∑Jk|k－1l=1ω(l)(j)k|k－1Nx;m(l)(j)k|k－1,P(l)(j)k|k－1
=υ(j)s,k|k－1(x)+υ(j)β,k|k－1(x)+υ(j)γ,k|k－1(x)（17）
其中υ(j)s,k|k－1(x)，υ(j)β,k|k－1(x)，υ(j)γ,k|k－1(x)分别代表旧目标、衍生目标和新生目标PHD强度。
旧目标PHD强度
υ(j)s,k|k－1(x)=
ps,k∑Jk－1l=1ω(l)(j)k－1Nx;m(l)(j)s,k|k－1,P(l)(j)s,k|k－1 （18）
m(l)(j)s,k|k－1=F(j)k－1(l)(j)k－1 （19）
P(l)(j)s,k|k－1=F(j)k－1(l)(j)k－1F(j)k－1T+Q(j)k－1 （20）
衍生目标PHD强度
υ(j)β,k|k－1=
∑Jk－1l=1∑Jβ,km=1(l)(j)k－1ω(m)(j)β,kNx;m(l,m)(j)β,k|k－1,P(l,m)(j)β,k|k－1 （21）
m(l,m)(j)β,k|k－1=F(m)(j)β,k－1(l)(j)k－1（22）
P(l,m)(j)β,k|k－1=F(m)(j)β,k－1(l)(j)k－1(F(m)(j)β,k－1)T+Q(m)(j)β,k－1（23）
其中，Jβ,k，ω(m)(j)β,k为衍生目标状态随机集的先验参数，这些参数决定了衍生目标集合的PHD。
2）条件PHD更新
计算目标的后验条件PHD为：
υ(j)k(x)=
(1－p(j)D,k)υ(j)k|k－1(x)+∑z∈Zkυ(j)D,k(x|z) （24）
其中，p(j)D,k表示在模型j下的目标检测概率； υ(j)D,k(x|z)表示模型j下，基于观测z的后验条件PHD，满足
υ(j)D,k(x|z)=∑Jk|k－1(1+nzk)l=Jk|k－1nzk+1ω(l)(j)D,k|kN(x;m(l)(j)D,k,P(l)(j)D,k) （25）
m(l)(j)D,k=m(l)(j)k|k－1+Kkz－Hkm(l)(j)k|k－1（26）
P(l)(j)D,k=(I－KkHk)P(l)(j)k|k－1 （27）
K(j)k=P(l)(j)k|k－1HT(HkP(l)(j)k|k－1(Hk)T+Rk)－1 （28）
对应高斯分量权值为
ω(l)(j)D,k=Λ(l)(j)kκk(z)+∑Jk|k－1(1+nzk)l=Jk|k－1nzk+1Λ(l)(j)k （29）
Λ(l)(j)k=
(1－P(j)D,k)ω(l)(j)k|k－1l=1,2,…,Jk|k－1P(j)D,kω(l)(j)k|k－1ρ(l)(j)k(z|m(l)(j)D,k,P(l)(k)D,k)l=(Jk|k－1+1),…,Jk （30）
ρ(l)(j)k(•)表示似然函数［9］；κk(z)表示杂波随机集的PHD。
3.2.3 模型概率预测与更新
1）模型概率预测
计算与υ(j)k|k－1中每个高斯分量相对应的预测模型概率，它的集合表示如下
μ(j)k|k－1=μ(1)(j)k|k－1,…,μ(Jk－1)(j)k|k－1,
μ(Jk－1+1)(j)k|k－1,…μ(Jk－1×(Jβ,k+1))(j)k|k－1,
μ(Jk－1×(Jβ,k+1)+1)(j)k|k－1,…,μ(Jk－1×(Jβ,k+1)+Jγ,k)(j)k|k－1 （31）
其中
μ(l)(j)k|k－1=(l)(j)k－1 （32）
μ(Jk－1×β+l)(j)k=(l)(j)k－1 （33）
μ(Jk－1×Jβ,k+1+g)k|k－1=μ(g)(j)γ（34）
μ(g)(j)γ表示模型j下新生PHD中的每个高斯分量的初始概率。
2）模型概率更新
模型概率集合μ(l)(j)kJkl=1中每个元素分量的更新可通过下式计算
μ(l)(j)k=1c(l)Λl(j)k(l)(j) （35）
c(l+m×Jk|k－1)=∑rj=1(l+m×Jk|k－1)j （36）
(l+m×Jk|k－1)j=Λ(l+m×Jk|k－1)k,j(l)j（37）
(l)j=∑ri=1πi,jμ(l)(i)k|k－1 （38）
3.2.4 PHD融合估计
k时刻总体PHD估计，即
υk(x)=∑jυ(j)(x)μ(j)k
=∑Jklω(l)kNx;m(l)k,P(l)k（39）
其中
m(l)k=∑rj=1m(l)(j)kμ(l)(j)k （40）
P(l)k=
∑rj=1P(l)(j)k+(m(l)k－
m(l)(j)k)(m(l)k－m(l)(j)k)Tμ(l)(j)k（41）
ω(l)k=∑rj=1ω(l)(j)kμ(l)(j)k （42）
3.2.5 剪枝与合并
在递归过程中，后验PHD的高斯分量会随时间无限增加。为了控制每时刻增加的高斯分量数目，可预置一个门限τ，将小于权值τ的高斯分量滤除掉，这种操作通常称为剪枝。此外，可预置一个合并阈值U，当后验PHD中的一些高斯分量由于彼此接近程度到达该阈值时，对这些高斯分量进行合并处理。
3.2.6 目标数目、状态估计
目标的数目估计是通过对融合后PHD的各高斯分量的权值求和取整得到，而目标的状态估计是根据估计的目标数目，依次提取混合权重最大的高斯分量作为状态估计值集合。这样便可完成当前时刻的IMMGMPHD滤波估计。
4 仿真分析
4.1 仿真场景设置
考虑杂波环境下的二维多目标运动场景，比较IMMGMPHD滤波算法和CSGMPHD滤波算法的优劣。设定监控时间长度为T=50s，采样间隔Δt=1s，方位角误差σφ=0.01rad，距离误差σρ=30m，检测概率pD=0.99，生存概率pS=0.9，剪枝阈值τ=1E－5，合并阈值U=4。
设杂波为泊松点过程，强度为κk(zk)=λcfc(zk)令λc=50，fc(zk)=ρ(zk)，其中ρ(•)表示服从均匀分布的概率密度。
在监控时间段内有四个不同机动方式的目标相继出现在监控区域内，如图1所示。
在图1中，Target 1 出生于第1秒，消亡于第45秒，全程作匀加速运动。 Target 2 出生于第10秒，在10～17秒作匀加速运动，在第18秒分散成Target 21和Target 22两个目标；Target 21在18～32秒作匀速运动在第32秒被Target 3击中消亡；Target 22在18～30秒作转弯机动，31～40秒作匀速运动，在第40秒消亡。Target 3出生于第20秒，在20～23秒作匀加速运动，在25～32秒作匀速运动，在第32秒击中Target 21一并消亡。
IMM模型设定如下：
模型1为CV模型，模型2为CA模型。
模型3为“当前”模型，大气扰动α=0.5，最大加速度amax =50 m/s2。
模型间的马尔可夫转移矩阵
πi,j=0.980.010.010.0250.950.0250.030.010.96
新生目标随机集合的PHD为：
υ(j)k=0.5×［Nx;m(1)γ,Pγ+Nx;m(2)γ,Pγ+Nx;m(3)γ,Pγ］
m(1)γ=6000－10005050Tmm(2)γ=－6000－10005050Tm
m(3)γ=6000－10007070TmPγ=diag(600600600600)m2
4.2 仿真结果与分析
图2和图3分别给出了IMMGMPHD和CSGMPHD滤波算法在相同的观测情况下，对空间中相继出现的多个机动目标的位置估计。
可以看出，不管目标是否发生机动，两种算法均能估计出多目标的位置，进一步比较可以看出，由于运动模型不匹配CSGMPHD滤波算法出现了明显的目标丢失现象和较大的估计误差；受杂波影响，CSGMPHD比IMMGMPHD产生更多的虚假目标估计，而IMMGMPHD表现出跟踪的连贯性和完整性。
图4和图5分别展示了100次MC仿真实验后，两种算法对目标数目估计随时间变化曲线。两种算法对目标数目的估计均较接近于真实值，但是由于杂波和目标机动的影响，CSGMPHD算法会产生显著的偏差；相反，IMMGMPHD算法受目标机动的影响较小，估计曲线波动相对较小，更接近于真实值。
从表1中可以看出，在本文设定的仿真条件下，IMMGMPHD滤波器的整体性能明显优于CSGMPHD滤波器，额外的时间消耗也是在合理范围之内。
5 结 论
针对现有多机动目标跟踪算法精度低、计算量大、约束条件苛刻等问题，本文将GMPHD滤波器和交互式多模型相结合，提出了IMMGMPHD滤波算法。新算法充分利用了IMM算法在机动目标跟踪方面的优势，结合基于随机集理论的GMPHD算法，不仅避免了多目标跟踪的数据关联问题，而且在复杂的背景环境下能够稳定、精确地估计目标数目和状态。IMMGMPHD滤波算法性能稳定，跟踪精度高，适应性强，有一定理论意义和工程应用价值。
参考文献
［1］ BAR SHALOM Y, LI X. Multitarget_multisensor tracking: Principles and Techniques ［M］.Storrs, CT: YBS Publishing, 1995.
［2］ MAHLER R P S. Multitarget bayes filtering via first-order multitarget moments ［J］. IEEE Transactions on Aerospace and Electronic Systems, 2003, 39(4): 11521178.
［3］ SHOENFELD PETER S. A particle filter algorithm for the multi-target probability hypothesis density ［J］. Signal Processing, Sensor Fusion, and Target Recognition XIII, 2004: 315-325.
［4］ VO B N, MA W K. The Gaussian mixture probability hypothesis densities filter ［J］. IEEE Transactions on Signal Processing, 2006, 54(11): 4091-4104.
［5］ PANTA K, VO BA NGU, SINGH S. Novel data association schemes for the probability hypothesis density filters ［J］. IEEE Transactions on Aerospace and Electronic Systems, 2007, 43(2):556-570.
［6］ SIMO SARKKA Vehtari1, Jouko Lampinen. Rao-Blackwellized particle filter for multiple target tracking［C］. 7th International Conference on Information Fusion, 2007 8(1):2-15.
［7］ LIU H, DOU L H, PAN F, DONG L X. Research on data association in three passive sensors network ［C］. 2007 IEEE International Conference on Control and Automation, ICCA, 2007: 3235-3238.
［8］ ARULAMPALAM M S,MASKELL S,GORDON N,CLAPP T.A tutorial on particle filters for online nonlinear/nongaussian bayesian tracking［J］.IEEE Trans. Signal Processing, 2002, 50:174-188.
［9］ B.Vo, S.SINGH, A. Doucet. Sequential monte carlo methods for multitarget filtering with random finite sets ［J］.IEEE Transactions on Aerospace and Electronic Systems, 2005.8, 41(4):2855-2865.省略)；闫玉德（1954―），男，江苏南京人，副教授，硕士生导师，研究方向：机动目标跟踪，数据融合技术，嵌入式系统设计。
