高斯过程回归模型多扩展目标多伯努利滤波器

陈 辉 ,李国财 ,韩崇昭 ,杜金瑞

(1.兰州理工大学电气工程与信息工程学院,甘肃兰州 730050;2.西安交通大学电子与信息工程学院综合自动化研究所,陕西西安 710049)

1 引言

传统多目标跟踪(multiple target tracking,MTT)方法[1–4]大多基于点目标模型,即假定目标是一个点,量测只源于目标的质心,假定每个目标在给定的时间最多生成一个量测值.近些年来,随着传感器分辨率的提高和对邻近目标的跟踪需求,目标不能再宽泛的假设为一个点,而是有一定空间范围的扩展,它表现在量测不只源于目标的质心,而是目标在每个采样周期内产生多个量测,学者们习惯上称这样的目标为扩展目标(extended target,ET),其对应的跟踪问题为扩展目标跟踪(extended target tracking,ETT)[5–6].ETT 问题的提出具有非常重大的现实意义,最近几年引起了国内外学者的高度关注.此外,在实际的跟踪问题中通常要对聚集目标(例如蜂群无人机)进行跟踪,即群目标跟踪(group target tracking,GTT)问题,如果把群目标看作一个整体,则这个广义目标往往具有相似的运动属性.这时候,就无需对群中的子目标进行单独跟踪,而跟踪的重点则集中于群目标整体运动特性和空间扩展.所以,群目标跟踪(GTT)与扩展目标跟踪(ETT)在跟踪方法上具有类似性.

众所周知,多目标跟踪(MTT)中关联关系的确定是解决最优估计问题的首要难题,对于多扩展目标(multi extended target,MET)而言,由于每个目标在采样周期内对应多个量测,即传统多“点目标”跟踪的假设不再成立,这将进一步复杂化目标和量测之间的关联关系,这使得多扩展目标跟踪(multi extended target tracking,METT)问题的解决变得异常困难.近些年来,为了避免去处理复杂的数据关联问题,Mahler基于有限集统计理论(finite set statistics,FISST)[2,7]的多目标跟踪方法受到了学者们的极大关注.一些近似的最优多目标滤波器陆续被提出,主要包括多目标矩递推滤波器[3]和多伯努利滤波器[8–11].由于多伯努利滤波器是利用多个伯努利过程直接近似递推了多目标后验概率密度,这使得贝叶斯多目标跟踪问题的求解更为直观.

METT另外一个核心议题就是如何利用目标的多个量测获取目标形状(轮廓)的最优估计,实际上最近几年对扩展目标形状估计的研究已经取得了一些引人瞩目的成果[12–16].在诸多具有形状信息的扩展目标跟踪建模方法中,随机矩阵模型(random matrix model,RMM)[17–18]和随机超曲面模型(random hypersurface model,RHM)[19–22]是其中最具代表性的两种算法.RMM用椭圆表征目标的形状,并使用高斯逆Wishart分布来表示含形状信息的扩展目标状态密度[23].然而,仅利用椭圆表示目标的基本轮廓特征,对目标形状的描述具有相对的局限性.不同于RMM,RHM使用随机超曲面表示扩展目标的形状[24–26],将径向函数作为扩展目标的形状函数,进而利用径向函数傅里叶级数展开的有限阶系数所构成的参数向量去表示星凸形扩展目标的不规则形状.但是,仅利用傅里叶级数展开只能在一定程度上去描述复杂的局部轮廓特征.

高斯过程回归(Gaussian process regression,GPR)易于处理的后验计算和极具吸引力的分析性质已经广泛的应用于信号处理、机器学习、概率回归和分类统计的研究当中[27–30].近几年,GPR作为一种非参数贝叶斯学习方法,已被用于扩展目标形状建模[31–33],并取得不错的跟踪估计效果.与已有星凸形多扩展目标跟踪方法不同的是,本文方法不使用傅立叶级数建模目标形状轮廓,而是在RHM方法中,借助GPR学习目标形状的径向函数,从而完成对星凸形扩展目标不规则形状的精确建模.

本文针对不规则形状星凸形多扩展(群)目标的跟踪问题,首先在利用Mahler基于FISST将星凸形多扩展目标的状态集和量测集分别建模为伯努利随机有限集(random finite set,RFS)和泊松随机有限集的基础上,利用GPR对星凸形扩展(群)目标的量测源建模.然后将GPR建模的量测空间分布模型和泊松观测模型[35–36]与势均衡MET 多伯努利(cardinality balanced multi extended target multi-Bernoulli,ETCBMeMBer)[37]滤波算法相结合,并进一步研究其滤波算法的高斯混合(Gaussian mixture,GM)实现方式,并结合具有高滤波精度以及具备解决含高维状态向量非线性问题能力的容积卡尔曼滤波(cubature Kalman filter,CKF)[34]对GPR模型的非线性量测方程进行求解.此外,本文利用拟Jaccard(quasi-Jaccard)[38]距离与最优次模式分配(optimal subpattern assignment,OSPA)距离[39]分别对星凸形多扩展(群)目标的形状估计与运动参数估计进行评测.最后,通过构造星凸形多扩展(群)目标的跟踪实验,验证了本文算法的有效性.

2 问题描述

2.1 高斯过程

2.2 高斯过程回归

首先,给出如式(8)所示的量测模型:

2.3 高斯过程回归建模

星凸形(star-convex shape)满足:如果集合E(x)⊂RN中至少存在一个点m与该集合内任意点x连线上的所有点构成的集合仍是该集合的子集,则集合E(x)所构成的形状为星凸形.为了便于计算,用径向函数r=f(θ)来描述星凸形,如图1所示.

图1 星凸形的径向函数描述Fig.1 Radius function representation of star-convex shape

径向函数是目标轮廓点与质心连线的半径关于对应角度的函数.为了同时估计星凸形扩展目标的状态,首先给出其联合估计模型.

1) 目标运动模型.

其中:T为采样周期,α为遗忘因子,决定动态速度和历史量测信息占的权重.即当α=0时,所有量测值同样重要,随着α的减小,历史量测信息的权重将增大.

2) 目标量测模型.

本文主要研究对象为星凸形扩展目标.基于该目标的径向函数描述,根据量测源的位置主要表明两种量测模型:目标轮廓模型和目标表面模型.

a) 目标轮廓模型:量测源于目标轮廓,则量测方程可以表示为

根据式(24)并结合式(33)利用GPR对星凸形扩展目标的径向函数描述进行建模,可得如式(36)标准量测方程.其中hk,l(xk)为

需要注意的是,式(39)已经包含了量测噪声协方差,所以式(22)中的量测噪声协方差可以忽略不计.

以上的量测模型只适合于量测源于目标轮廓边缘的情况.当量测源于目标的表面时(目标内部),通过利用尺度因子缩放轮廓模型可得目标表面量测所满足的模型.其中,尺度因子sk,l服从如下均匀分布:

为了满足高斯滤波器的推理条件,则采用一个高斯随机变量去近似尺度因子[27],即

以上假设是进行高斯滤波器递推的一种必要性假设.文献[27]特别证明了这种假设的合理性,当量测均匀分布在物体表面,而尺度因子设置为高斯分布函数时,基于随机超曲面模型的滤波器也能获得较好的估计效果.进而目标表面量测方程可整理为

由尺度因子的分布可知,当尺度因子sk,l=1时,该方程与式(36)相同,进而无论量测源于目标何处,目标的量测方程都如式(45)所示.

3 基于GPR的多扩展目标多伯努利滤波算法

基于以上利用GPR建立的扩展目标运动模型和量测模型,重点研究利用ETCBMeMBer滤波算法对星凸形扩展目标(ET)的运动参数和形状参数进行联合递推估计,并且详细推导基于GPR的ETCBMeMBer滤波算法的具体GM实现过程.但由于GPR量测方程的非线性,这使得问题的求解变得更加复杂.在解决非线性滤波的问题中,本文利用CKF算法对GPR非线性量测方程的高斯加权积分进行近似.CKF采用三阶球面–径向容积准则对非线性函数的高斯加权积分进行近似,具有比无损卡尔曼滤波(unscented Kalman filter,UKF)更高的滤波精度和数值稳定性,而且适用于求解包含高维状态向量的非线性问题.下面将提出并详细推导高斯过程回归高斯混合多扩展目标多伯努利(GPR–GM–ETCBMeMBer)滤波算法.

3.1 GPR–GM–ETCBMeMBer滤波算法

1) 预测步.

GPR–GM–ETCBMeMBer滤波算法的预测步与点目标CBMeMBer的预测步相同,即如果k −1时刻的多扩展目标后验概率密度是一个多伯努利形式:

其中前一部分为存活目标多伯努利参数,后一部分表示k时刻新生目标多伯努利参数,且概率密度函数可由高斯混合加权和近似表示为

2) 更新步.

对于扩展目标而言,每个目标对应多个量测,其量测更新步是研究的重点.首先本文采用距离划分算法对每个周期传感器得到的所有量测进行划分,进而将所有扩展目标的量测划分为量测集(簇)的形式,然后利用泊松似然模型进行更新.

由预测步可知,k时刻的多目标后验概率密度由漏检部分和量测更新部分组成,即

其中:℘为量测集Zk的划分方式,W为每种划分方式中的量测单元数.另外,漏检部分如下:

3.2 GPR–GM–ETCBMeMBer算法实现伪码

为了更具体地说明本文算法的递推滤波过程,以下给出该算法的高斯混合实现伪码.

4 多扩展(群)目标的运动参数估计与形状参数估计性能评价

多扩展(群)目标的估计状态向量由运动状态和扩展状态(形状)两部分组成.因此,本文引入OSPA距离和拟Jaccard(quasi-Jaccard)距离分别对扩展目标质心位置估计和形状估计的性能进行客观全面的评估.

1) OSPA距离评价质心位置估计:首先由OSPA距离的定义可知:对于两个任意的集合和,若集合X的元素个数小于集合Y 的元素个数,即mn,则这两个集合之间的OSPA距离为

其中:d(c)(x,y)=min(c,‖x,y‖),Πk表示集合{1,2,···,k}的全排列,π(i)表示该全排列的第i种排列.距离阶次p ∈(1,∞);截断系数c>0,用于调节集合元素数误差的影响;当m>n时,

2) 拟Jaccard距离评价形状估计:为评估滤波器对多扩展目标形状的估计性能,文献[38]针对多扩展目标依据Jaccard距离的思想提出拟Jaccard距离方法对多扩展目标形状的估计效果进行评价,具体计算过程如下:

首先,对于多扩展目标而言,利用拟Jaccard距离的关键是真实目标与估计目标之间关联矩阵的建立.即利用估计目标质心位置和真实目标质心位置找到估计目标形状和真实目标形状的匹配关系,即

最后,根据式(74)分别计算k时刻单个扩展目标的拟Jaccard距离,再对各个扩展目标的拟Jaccard距离求均值作为k时刻整体形状估计的评价指标.

5 仿真实验

通过构造具有不规则形状,且做匀速运动(CV)的多扩展目标和多群目标跟踪的仿真场景,同时选取基于RHM的高斯混合多扩展目标多伯努利(RHM–GM–ETCBMeMBer)滤波算法与本文所提出的GPR–GM–ETCBMeMBer滤波算法相比较,验证所提算法对多扩展(群)目标跟踪的有效性.仿真场景设置监测时间长度为50 s,采样周期为T=1 s,观测区域为[−500,500]m×[−500,500]m,杂波个数服从均值为λ=5的泊松分布,即杂波密度为5×10−6/m2.每个扩展目标的运动状态演化服从相互独立的线性高斯运动模型,即k时刻第i个扩展目标的状态演化可描述为如下模型:

其中:过程噪声标准差σq=0.01 m/s2;量测噪声协方差为Rk=diag{0.12,0.12};分别由式(31)和式(32)计算得到;遗忘因子的选择参考文献[27]的做法,即遗忘因子α=1×10−5.形状基本点个数(形状参数维数)NB=50,即相应的角度为

高斯过程的超参数分别为:σr=0.2,σf=0.5,l=2π.每个扩展目标量测个数服从均值为γ=20的泊松分布.尺度因子sk,l服从为了客观验证算法的性能,根据统计指标的选定,将构造100 次蒙特卡罗(Monte Carlo,MC)仿真实验评价滤波器性能.

扩展目标存活概率pS=0.99,传感器的检测概率为pD=0.99.新生目标服从伯努利分布,概率密度

表1 多扩展目标初始参数Table 1 Initial parameters of MET

图2给出了多扩展(群)目标随时间运动的真实轨迹.在本文的多伯努利滤波算法中,每条航迹的高斯分量的截断阈值为TP=10−3,其合并阈值为U=4,最大项数为Jmax=100.为了避免假设航迹数的不断增加,设其存在概率阈值为Tr=10−2,最大航迹数Tmax=10.根据已有扩展目标跟踪文献场景中OSPA距离参数的配置,选择c=1.5,p=1,轮廓点个数为50.

图2 目标真实运动轨迹Fig.2 Actual target trajectories

5.1 星凸形多扩展目标跟踪

依据上述仿真场景,采用本文提出的GPR–GM–ETCBMeMBer滤波算法和RHM–GM–ETCBMeMBer滤波算法分别对如图3所示的多个扩展目标进行跟踪估计,同时为了验证算法对星凸不规则形状的估计能力,实验构造出几种颇具代表性的星凸型扩展目标,进而构造MC仿真实验对比两种滤波器对多扩展目标的跟踪估计性能.其中目标1和目标2形状如图3(a),目标3形状如图3(b)所示.

图3 扩展目标形状及量测分布Fig.3 Shape and measurements distribution of ET

图4为单次仿真实验下两种滤波器对多个具有不规则形状的扩展目标整体跟踪效果图.图中可得,两种滤波器都能够精确地跟踪估计多个扩展目标.为了能清晰地反映出两种滤波器对目标形状的估计效果,图5给出了在目标新生时刻(图5(a),5(c)和5(e))和目标消亡时刻(图5(b),5(d)和5(f))对每个星凸形扩展形状估计的局部放大效果图.由图5可知,由于在新生时刻形状先验未知和量测信息累积较少,对MET形状估计效果比较差,但随着量测信息的累积,两种滤波器在目标消亡时刻都能够较为准确的跟踪估计多个星凸形扩展目标的形状且稳定收敛于目标的真实形状,并且从局部估计放大图可以看出本文算法对MET形状的估计要优于RHM–GM–ETCBMeMBer算法.

图4 两种滤波器对多扩展目标的跟踪效果图Fig.4 The tracking result to MET of the two filters

图6–8分别给出了在100次MC仿真实验下,两种滤波器对多扩展目标的数目(势)、质心位置和形状的估计结果对比.图6显示出两种滤波器的势估计都是无偏的,这与ETCBMeMBer滤波器对目标势的无偏估计一致.图7是两种滤波器对MET质心位置估计的OSPA距离,可以看出两种OSPA距离都非常小,进而表明两种算法都能对多扩展目标进行比较准确的定位,但本文算法对扩展目标的定位要更加精确.随着量测信息的不断累积,对多扩展目标的质心估计将更加准确,因而OSPA距离值将随着时间逐渐减小并趋于稳定.

图5 两种滤波器对多扩展目标形状估计的局部放大图Fig.5 The partial enlarged effect of the two filters for shape estimation of MET

图6 两种滤波器对多扩展目标的势估计Fig.6 Cardinality estimation statistics of the two filters to MET

图7 多扩展目标质心位置估计的OSPAFig.7 OSPA statistics of the MET centroid position estimation

图8显示了两种滤波器对多个星凸形扩展目标形状估计的拟Jaccard距离统计值.统计值表明两种滤波器都能够对MET的形状进行比较精确的跟踪估计.但相较于RHM–GM–ETCBMeMBer滤波器,本文算法对MET形状的估计更加精确.该统计值分别在t=15 s,t=25 s和t=41 s时突然增大,之后随时间推移又逐渐减小至趋于稳定,这可以得到很好的解释,即t=15 s,t=25 s时刻出现了新生目标,导致拟Jaccard距离值增大.t=41 s出现了消亡目标,形状已经估计准确的目标在此刻消失导致整体的拟Jaccard距离值变大.这一系列现象的出现进一步验证了本文算法对星凸形多扩展目标形状估计性能的适应性和有效性.

图8 多扩展目标形状估计的拟Jaccard距离Fig.8 Quasi-Jaccard distance statistics of the MET shape estimation

5.2 星凸形多群目标跟踪

为了验证本文所提算法对多个星凸形群目标的跟踪估计性能,本节构造对形如图9的多个星凸形群目标进行跟踪估计的实验.群目标航迹依然如图2所示,每个群目标由20个具有相同运动特性的子目标组成.其中目标1和目标3形如图9(a),目标2形如图9(b).为了对比RHM–GM–ETCBMeMBer 与本文所提出的GPR–GM–ETCBMeMBer 滤波算法对多群目标的跟踪估计性能,本文构造MC仿真实验,结果如图10–14所示.

图9 群目标形状及量测分布Fig.9 Shape and measurements of the group target

图10–11分别为单次MC仿真实验下两种滤波器对多个群目标的整体跟踪效果图和每个群目标在新生时刻(图11(a),11(c)和11(e))以及消亡时刻(图11(b),11(d)和11(f))轮廓估计效果的局部放大图.

图10 两种滤波器对多群目标的跟踪效果图Fig.10 The tracking result to multiple group target of the two filters

图11 两种滤波器对多群目标形状估计的局部放大图Fig.11 The partial enlarged effect of the two filters for shape estimation of multiple group target

由图10可以看出,两种滤波器都能同时实现对多个具有复杂轮廓的群目标的实时跟踪.且由图11可以看出随着量测信息的累积,在消亡时刻目标的轮廓估计愈加准确且收敛于目标的真实轮廓,同时可以明显的看出本文所提算法对群目标的轮廓估计较为准确,进而验证了该算法对于不规则轮廓群目标跟踪的有效性.

图12–14是在100次独立MC仿真实验下,两种滤波器对多个群目标的目标数、质心位置、轮廓的估计结果对比.由图12可知两种滤波器都能够准确地估计出群目标数的变化.图13可以看出,两种滤波器都能较为准确的对群目标进行跟踪定位.注意到,t=15 s时刻OSPA距离突然增大,这是因为在该时刻出现了新生目标,目标新生时刻,质心位置估计不准确导致OSPA距离值变大.t=41 s已经定位准确的扩展目标消失,导致OSPA距离值出现略微增大的现象.

图12 两种滤波器对多群目标的势估计Fig.12 Cardinality estimation statistics of the two filters to multiple group target

图13 多群目标质心位置估计的OSPAFig.13 OSPA statistics of the multiple group target centroid position estimation

由图14给出的两种滤波算法对多个不规则形状群目标轮廓估计的拟Jaccard距离统计值可以看出:两种滤波器都能够对多个群目标轮廓形状进行有效的估计.且由该值的统计规律可知,本文提出的GPR–GM–ETCBMeMBer算法对群目标轮廓形状的估计要优于RHM–GM–ETCBMeMBer算法.而且,多目标形状的拟Jaccard距离值与群目标数的变化一致,即拟Jaccard距离统计值在群目标新生和消亡时刻都将会发生明显的变化.

图14 多群目标形状估计的拟Jaccard距离Fig.14 Quasi-Jaccard distance statistics of the multiple group target shape estimation

综合以上两个仿真,本文所提的滤波算法能够同时对杂波、漏检、噪声共存环境下不规则形状星凸形多扩展(群)目标的目标数、运动状态和形状进行较为准确的估计.与RHM–GM–ETCBMeMBer算法相比,本文所提算法较RHM–GM–ETCBMeMBer算法在形状估计等诸多性能上都有比较明显的提升.

实际上,GPR是计算开销较大的算法,在本仿真实验中,确实也存在着GPR学习径向函数计算开销较大的问题.另外,本文是一个GPR高维情形应用的实例,有重要的理论价值和意义.所以,希望在今后的研究工作中,将进一步提高该算法的计算效率,开展工程算法的应用性研究.

6 结论

本文在采用高斯过程回归对星凸形扩展目标进行建模的基础上,利用多伯努利随机有限集和泊松随机有限集分别描述多扩展目标的状态和观测,结合多扩展目标多伯努利滤波器,提出了GPR–GM–ETCBMe-MBer滤波算法,并详细推导并给出其高斯混合具体实现过程.最后,在杂波、漏检、噪声共存的跟踪环境下,分别对多个具有不规则形状的扩展目标和群目标进行跟踪估计,验证了该算法的有效性.仿真结果表明,相对于已有的基于傅里叶级数近似的RHM–GM–ETCBMeMBer滤波器,本文算法能够更为有效地对多扩展(群)目标的运动状态、形状轮廓进行精确地估计.对于未来的研究,有必要将文中算法扩展到三维空间,这更符合实际的跟踪需求.另外,本文算法的计算效率有待提高,结合GPR方法进行形状特征分析与高效识别的工程算法应用性研究将刻不容缓.