面向监视系统的鲁棒群目标跟踪方法及跟踪系统与流程

本发明涉及鲁棒群目标跟踪方法技术领域，具体地，涉及一种面向监视系统的鲁棒群目标跟踪方法，尤其是一种基于随机矩阵模型的鲁棒群目标跟踪方法。

背景技术：

群目标的定义为：满足一定目标间距的条件下，在一定长的时间内保持空间关系(如位置等)相对固定(或具有相似的状态参数)，并且至少两个以上可分辨或不可分辨目标混合形成的目标簇。群目标跟踪需要同时解决群的整体运动和形状的联合估计问题。基于随机矩阵模型的群目标跟踪方法，利用对称正定矩阵来表征群的形状。即在二维平面，群的形状为具有一定大小和方位的椭圆。

对现有的文献检索发现，尽管基于随机矩阵模型的群目标跟踪方法优于以往的群目标跟踪方法，但真实的群目标跟踪场景中往往面临着传感器噪声和环境杂波的干扰，已有的基于随机矩阵模型的群目标跟踪结果会受到严重影响。当观测噪声不服从高斯分布而服从某种重尾分布，传感器往往返回一些野值引起状态估计方差变大，进而造成群目标估计性能的下降。其次，需要考虑环境中杂波对群目标跟踪的影响，对环境杂波进行建模可提高估计性能。实际上，在应用场景中需要一种鲁棒的群目标跟踪方法来解决服从重尾分布的观测噪声、环境杂波带来的问题，还能在对群目标进行跟踪的同时对杂波密度进行估计。

技术实现要素：

针对现有技术中的缺陷，本发明的目的是提供一种面向监视系统的鲁棒群目标跟踪方法及跟踪系统，解决服从重尾分布的观测噪声、环境杂波带来的问题，同时在对群目标进行跟踪的时候对杂波密度进行估计，实现了在传感器观测噪声具有重尾分布特性而产生异常量测，监视环境中存在杂波干扰情况下的联合估计群目标形状、运动状态和环境杂波密度。

根据本发明提供的一种面向监视系统的鲁棒群目标跟踪方法，包括以下步骤：

s1、混合观测模型的建立：建立监视系统量测的混合观测模型，所述混合观测模型由源自群目标的观测模型和源自环境杂波的观测模型组成，根据混合观测模型构建监视系统量测的似然函数；

s2、联合概率密度的计算：根据监视系统量测的似然函数，计算监视系统的各状态变量和量测的完全数据联合概率密度；

s3、近似后验概率密度的计算：基于共轭分布假设，最小化各状态变量的近似后验概率密度的kullback-leibler散度，然后基于步骤s2得到的各状态变量和量测的完全数据联合概率密度计算各状态变量的近似后验概率密度；

s4、量测中的杂波密度参数的计算：保持各状态变量的近似后验概率密度不变，最大化监视系统量测的似然函数下界计算量测中的杂波密度参数；

s5、各状态变量的估计值计算：根据各状态变量的近似后验概率密度和量测中的杂波密度参数，计算各状态变量的数学期望从而得到各状态变量的估计值，重复执行步骤s2－s4指定次数，认为跟踪结束；

s6、群目标的跟踪估计：根据最终得到的各状态变量的估计值和杂波密度参数，实现群目标运动状态、形状和环境杂波密度的估计。

本发明还提供一种面向监视系统的鲁棒群目标跟踪系统，包括，

混合观测模型的构建模块：用于构建监视系统量测的混合观测模型，所述混合观测模型由源自群目标的观测模型和源自环境杂波的观测模型组成，根据混合观测模型构建监视系统量测的似然函数；

联合概率密度的计算模块：根据监视系统量测的似然函数，计算监视系统的各状态变量和量测的完全数据联合概率密度；

近似后验概率密度的计算模块：最小化各状态变量的近似后验概率密度的kullback-leibler散度，然后基于联合概率密度计算模块得到的各状态变量和量测的完全数据联合概率密度计算各状态变量的近似后验概率密度；

杂波密度参数的计算模块：保持各状态变量的近似后验概率密度不变，最大化监视系统量测的似然函数下界计算量测中的杂波密度参数；

各状态变量估计值的计算模块：根据各状态变量的近似后验概率密度和量测中的杂波密度参数，计算各状态变量的数学期望从而得到各状态变量的估计值，重复执行计算指定次数，认为跟踪结束；

群目标的跟踪估计模块：根据最终得到的各状态变量的估计值和杂波密度参数，实现群目标运动状态、形状和环境杂波密度的估计。

在步骤s1或混合观测模型的构建模块中，源自群目标的观测模型的量测服从重尾学生t分布，其均值为零，方差为群目标形状xt，自由度为v。时刻t的第j个量测ztj服从以下学生分布：

其中，xt为t时刻群目标运动状态，，xt为t时刻群目标的形状状态，ztj为时刻t的第j个量测，参数v≥1控制着学生t分布重尾的程度，c是系统的观测矩阵，t表示为矩阵转置，d是ztj的维度，可以为2(二维平面)或者3(三维空间)。量测ztj的等价高斯分布表达为：

p(ztj|xt)＝n(cxt,xt/wtj)

其中变量代表了一组nt个的权重变量。每个权重变量都一一对应时刻t的一个量测，权重变量wtj的取值表示量测对后验估计的贡献；

源自环境杂波的观测模型的量测服从均匀分布假设杂波量测在空间中均匀地分布在监视区域θ内，即杂波量测服从以下分布：

其中u(·)为均匀密度分布。监视区域的杂波量测数量假设服从已知均值的泊松分布。对于服从非均匀分布的环境背景，可以采用有限混合模型来模拟环境噪声的多峰分布，典型的混合模型如高斯混合模型。

综上，量测服从以下混合生成模型：

其中c是监视系统的观测矩阵，参数v≥1控制着学生t分布重尾的程度，在本应用中，推荐选取参数v＝4，θ代表监视区域，当量测源自杂波时，杂波量测在空间中均匀地分布在监视区域θ内，标识符是一组nt量测的二值变量，ctj为量测类别变量，取值为1代表量测由群目标产生并服从学生分布，取值为0代表量测为环境杂波产生并服从均匀分布，ztj为时刻t的第j个量测。标识符ct的概率密度如下：

其中参数θ为量测中杂波密度，基于上述定义，量测似然函数为一组源自群目标且服从重尾学生t分布的量测与源自环境杂波且服从均匀分布的量测求乘积得到，可以表达为：

其中代表二项分布，是时刻t源自群目标和环境杂波的nt个观测，变量代表了一组nt个的权重变量，每个权重变量wtj都一一对应时刻t的一个量测，权重变量wtj的取值表示第j个量测对后验估计的贡献。

在步骤s2或联合概率密度的计算模块中，监视系统的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)由监视系统量测的似然函数p(zt|xt,xt,wt,ct)、xt和xt的联合概率密度p(xt,xt|z^t-1)、权重变量wt的概率密度p(wt)，一组量测的二值变量ct的概率密度p(ct)，相乘得到。监视系统的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)的具体计算过程如下：

n(xt；xt|t-1,pt|t-1)iw(xt；vt|t-1,vt|t-1)

其中xt为t时刻群目标运动状态，服从高斯分布n(·)，xt|t-1是群目标运动状态t-1时刻的预测状态，pt|t-1是群目标运动状态t-1时刻的预测状态的协方差矩阵。gamma(·)是伽玛分布，b(·)是二项分布，xt为t时刻群目标的形状状态，群目标的形状预测状态xt|t-1服从逆维希特分布iw，iw分布自由度为vt|t-1，协方差矩阵为vt|t-1，z^t-1为从开始时刻到时刻t-1的累积量测，参数v≥1控制着学生t分布重尾的程度。

在步骤s3或近似后验概率密度的计算模块中，最小化各状态变量的近似后验概率密度的kullback-leibler散度，然后基于步骤s2得到的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)通过迭代解法计算各状态变量的近似后验概率密度q(xt)，q(xt)，q(wt)和q(ct)。具体计算方法为，对系统各状态变量的联合概率密度p(xt,xt,wt,ct|z^t)进行因子分解，并根据各状态变量分布的共轭假设定义每个状态变量的先验分布和后验分布都属于同一种分布，可以得到：

其中q(xt)，q(xt)，q(wt)和q(ct)分别是状态变量xt，xt，wt和ct的近似后验概率密度。根据变分贝叶斯方法，各状态变量的近似后验概率密度可以通过最小化近似概率密度q(xt)，q(xt)，q(wt)和q(ct)与系统各状态变量的联合概率密度p(xt,xt,wt,ct|z^t)之间的kullback-leibler散度(kld)求得，即：

q(xt),q(xt),q(wt),q(ct)＝argminkld(q(xt)q(xt)q(wt)q(ct)||p(xt,xt,wt,ct|z^t))，argmin是最小化符号，kld是kullback-leibler散度(kullback-leiblerdivergence)的缩写，argminkld表示最小化kl散度。

对上式进行优化求解可以得到各状态变量的近似后验概率密度迭代解法，对其中一个状态变量概率密度进行优化，保持其他状态变量概率密度不变可以得到监视系统的各状态变量的近似后验概率密度的以下解析公式：

其中e[·]代表随机变量的期望。const是相对于状态变量xt，xt，wt和ct的常量。假设变量wt和ct独立于系统状态xt和xt，上述公式中的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)可以分解为

计算得到各状态变量和量测的完全数据联合概率密度的对数logp(xt,xt,wt,ct,zt|z^t-1)的表达如下：

第(i+1)次迭代后qx(·)、qx(·)、qw(·)和qc(·)表达式分别为和在下文中，上标(i)和(i+1)表示第i次和第i+1次的表达式。通过计算各状态变量和量测的完全数据联合概率密度对数logp(xt,xt,wt,ct,zt|z^t-1)的期望，可给出以上各状态概率密度的表达如下，

基于量测ztj的wtj的近似后验概率密度对数表达为：

其中和分别代表变量xt、和ctj的期望，即和是时刻t第i次迭代的目标协方差矩阵，t表示为矩阵转置，tr代表矩阵的迹，为t时刻群目标的形状状态xt的逆矩阵。

对上式两端取指数并进行归一化，概率密度服从分布。

状态xt的概率密度对数表达可以计算为：

其中

其中pt|t-1为时刻t第i次迭代的目标预测协方差矩阵，是pt|t-1的逆矩阵。

状态xt的概率密度对数表达可以计算为：

是时刻t第i次迭代的目标协方差矩阵，t表示为矩阵转置，tr代表矩阵的迹。

以上公式最后一个等式是逆维希特分布的对数形式，将等式中的系数相加可以得到：

形状变量xt服从分布，并且iw分布自由度为协方差矩阵为

其中

基于量测ztj的变量ctj的概率密度对数表达可以计算为：

其中，q是中间变量，

在以上公式的推导中，可以看出和不存在解析解，故采用一阶泰勒展开进行近似估计，进而得到和

根据以上公式最后一个等式右边可以得到ctj服从伯努利分布，q(ctj)和q(ct)的表达式如下：

基于上述推导即可对进行迭代求解。

在步骤s4或杂波密度参数的计算模块中，保持各状态变量的近似后验概率密度不变，通过对监视系统量测的似然函数下界求杂波密度参数θ的导数并令其等于0，求解得到杂波密度参数θ的值：对数似然函数下界可表示为：

求解可得：

其中const是相对于参数θ的常数，是从开始时刻到时刻t的累积观测

在步骤s5或各状态变量估计值的计算模块中，分别计算群目标运动状态xt、群目标的形状状态xt、一组量测的权重变量wt、一组量测的二值变量ct的数学期望从而得到各状态变量的估计值，

由于wt的第j个变量wtj服从伽马分布，取其均值作为状态估计，即由于xt服从高斯分布，取其均值作为状态估计，即由于xt分布服从逆维希特分布，取其均值作为状态估计，则根据逆维希特分布性质由于ct的第j个变量服从二项分布，取其均值作为状态估计，即

与现有技术相比，本发明具有如下的有益效果：

1、本发明提出了一种新的基于随机矩阵模型的鲁棒群目标跟踪算法。采用学生t分布不仅可以实现对群目标形状的建模，而且通过定义每个量测的权重系数可以刻画该量测对于群目标估计的贡献。通过构建量测的混合生成模型，在群目标系统状态模型中集成环境杂波密度估计，利用环境杂波密度估计结果提高了群目标状态估计的性能。在变分贝叶斯方法的框架下，采用期望最大化的数学方法，通过估计和最大化两步骤不断迭代的方式计算出系统状态各变量的概率密度及估计值，除此之外，还可得到杂波密度参数。所提出的算法可以采用对量测迭代处理的实现方式，使得计算复杂度符合实际工程需求。

2、本发明的目的在于提出了一种在传感器观测噪声具有重尾分布特性而产生异常量测，监视环境中存在杂波干扰情况下的联合估计群目标形状、运动状态和环境杂波密度的方法。本方法针对群目标量测建立混合生成模型，基于学生t分布和变分贝叶斯方法，得到群目标跟踪系统完全数据对数似然函数期望的计算方法，在最大似然准则下，得到系统各状态变量的概率密度及杂波密度参数，可同时实现群目标运动状态、形状和环境杂波密度的估计。

3、本发明利用学生t分布和变分贝叶斯方法，提供了一种基于随机矩阵模型的鲁棒群目标跟踪算法，该算法可以采用对量测迭代处理的实现方式，并且对群目标运动状态、形状和环境杂波密度实现同时估计。此外，本算法框架清晰利于实现，从而为复杂环境下的群目标跟踪系统提供了重要的技术支持。

附图说明

通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：

图1为本发明中理想情况下群目标运动轨迹；

图2为本发明中理想情况下群目标位置估计比较；

图3为本发明中理想情况下群目标形状估计比较；

图4为本发明中观测噪声服从重尾分布情况下群目标运动轨迹；

图5为本发明中观测噪声服从重尾分布情况下群目标位置估计比较；

图6为本发明中观测噪声服从重尾分布情况下群目标形状估计比较；

图7为本发明中杂波及重尾分布情况下群目标运动轨迹；

图8为本发明中杂波及重尾分布情况下群目标位置估计比较；

图9为本发明中杂波及重尾分布情况下群目标形状估计比较；

图10为本发明中杂波及重尾分布情况下杂波密度估计。

具体实施方式

下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变化和改进。这些都属于本发明的保护范围。

从图1可以看出群目标的方位随着运动轨迹的变化而随之变化。

图2和图3分别给出了在没有考虑观测误差的重尾分布和环境杂波影响的理想情况下三种算法的位置估计rms和形状估计rms比较结果。koch方法在此理想的仿真场景中具有最小的估计误差，对于位置估计我们的方法非常接近于feldmann方法，对于形状估计我们的方法优于feldmann方法。

图4给出了重尾分布条件下的群目标运动轨迹图,“+”代表量测。从图5可以看出重尾分布对于量测的影响，少量的量测有着非常大的误差，这对跟踪算法提出了很高的要求。传统的跟踪算法由于没有考虑此类异常的势必会造成估计性能的下降。

图5和图6分别给出了在重尾分布情况下三种算法的位置估计rms和形状估计rms比较结果。我们的方法由于使用了学生t分布来拟合量测异常点从而获得了最好的位置估计和形状估计，koch方法和feldmann方法由于对于量测中的异常点考虑不足故性能都有大的下降，koch方法由于只针对理想情况故估计误差最大。

图7给出了杂波条件下的群目标运动轨迹图，量测用+号表示，而杂波用点表示。从图8可以看出监视区域内的异常点更多分布更为广泛。

图8和图9分别给出了在杂波情况下三种算法的位置估计rms和形状估计rms比较结果。我们的方法由于使用了学生t分布来拟合量测异常点从而获得了最好的位置估计和形状估计，尤其是针对形状估计koch方法和feldmann方法几乎发散。这就说明了在真实场景中非常常见的杂波对于koch方法和feldmann的影响巨大。

图10给出了杂波密度估计，我们的方法对于杂波密度的估计值非常接近与真实值。

实施例

本实施例中，本发明的面向监视系统的鲁棒群目标跟踪方法，介绍如下：

s1、混合观测模型的建立：建立监视系统量测的混合观测模型，所述混合观测模型由源自群目标的观测模型和源自环境杂波的观测模型组成，根据混合观测模型构建监视系统量测的似然函数p(zt|xt,xt,wt,ct)；

s2、联合概率密度的计算：根据监视系统量测的似然函数，计算监视系统的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)；

s3、近似后验概率密度的计算：基于共轭分布假设，最小化各状态变量的近似后验概率密度q(xt)，q(xt)，q(wt)和q(ct)的kullback-leibler散度，然后基于步骤s2得到的各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)通过迭代解法计算各状态变量的近似后验概率密度q(xt)，q(xt)，q(wt)和q(ct)；

s4、量测中的杂波密度参数的计算：保持各状态变量的近似后验概率密度不变，最大化监视系统量测的似然函数下界计算量测中的杂波密度参数θ，

s5、各状态变量的估计值计算：根据各状态变量的近似后验概率密度q(xt)，q(xt)，q(wt)和q(ct)和量测中的杂波密度参数θ，计算各状态变量的数学期望从而得到各状态变量的估计值，重复执行步骤s2－s4指定次数，认为跟踪结束；

s6、群目标的跟踪估计：根据最终得到的各状态变量的估计值和杂波密度参数，实现群目标运动状态、形状和环境杂波密度的估计。

接下来对本发明进行详细的描述。

本实施例提供的基于随机矩阵模型的鲁棒群目标联合运动状态和形状跟踪方法，经过利用模拟真实环境的传感器和背景环境产生的测试数据进行测试，实施步骤如下：

步骤1，在理想情况下，观测噪声可假设服从高斯分布。在实际应用中，由于环境恶劣、干扰较多、仪器精度等因素，观测噪声往往并不服从理想的高斯分布而服从一些重尾分布，造成传感器往往返回一些野值，引起状态估计方差变大，进而造成群目标状态估计性能的下降。其次，考虑到环境中存在杂波，在跟踪过程中需要考虑杂波对群目标状态估计的影响。故采用学生t分布来对观测噪声的重尾分布和量测野值进行建模，采用均匀分布或有限混合模型对环境杂波进行建模。引入一组权重变量wt来表征各个观测对状态估计的贡献。同样地，引入一组二进制变量ct来表征该量测属于目标或者杂波。

步骤2，计算系统各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)。根据贝叶斯定理，根据变量wt和ct独立于系统状态xt和xt，各状态变量和量测的完全数据联合概率密度p(xt,xt,wt,ct,zt|z^t-1)可分解为

n(xt；xt|t-1,pt|t-1)iw(xt；vt|t-1,vt|t-1)

步骤3，根据变分贝叶斯方法，各状态变量的近似后验概率密度可以通过最小化近似概率密度q(xt)，q(xt)，q(wt)和q(ct)与系统的各状态变量联概率密度p(xt,xt,wt,ct|z^t)之间的kullback-leibler散度求得以下解析公式：

其中e[·]代表随机变量的期望。const是相对于状态变量xt，xt，wt和ct的常量。上述公式中的系统各状态变量和量测的完全数据联合概率密度对数logp(xt,xt,wt,ct,zt|z^t-1)的完整表达如下：

第(i+1)次迭代后qx(·)、qx(·)、qw(·)和qc(·)表达式分别为和在下文中，上标(i)和(i+1)表示第i次和第i+1次的表达式。其中和分别代表变量wtj、xt、和ctj的期望，即和

权重变量wtj服从分布，其中，γtj为临时变量，为了简化表达式，

运动状态xt服从分布，其中

形状变量xt服从分布，其中

量测类别变量ctj服从伯努利分布，即其中

基于上述推导得到的概率密度的迭代求解过程。

步骤4，在保持各状态变量后验概率密度不变的情况下，通过最大化对数似然函数下界对未知参数θ进行估计。通过对似然函数下界求参数θ的导数并令其等于0，求解可得参数θ的值。

在计算过程中，当ctj出现极小值时会影响程序的稳定性并造成算法发散。因此，在程序中定义ctj的门限值，当超过预设门限值则量测来源于目标。在每一轮的迭代过程中，大于门限值的ctj被置于1；而小于门限值的ctj则被置于0。未知参数θ的计算公式为：

步骤5，重复执行步骤2－步骤4指定次数，认为算法已经收敛，计算数学期望得到各变量估计值，即

本发明提供的基于随机矩阵模型的鲁棒群目标跟踪算法，可在传感器观测噪声具有重尾分布特性而产生异常量测，监视环境中存在杂波干扰情况下实现对群目标运动状态、形状和环境杂波密度的联合估计。针对模拟真实环境的传感器和背景环境产生的测试数据，基于学生t分布和变分贝叶斯方法，得到群目标跟踪系统完全数据对数似然函数的计算方法。采用期望最大化的数学方法，通过估计和最大化两步骤不断迭代的方式计算出系统状态各变量的概率密度及估计值。所提出的算法非常灵活，可以采用对量测迭代处理的实现方式，使得计算复杂度符合实际工程需求。除此之外，通过极大似然估计可得到杂波密度参数。本算法框架清晰利于实现，可以达到实时计算的性能要求。从而可在真实场景群目标跟踪监视系统中得到广泛应用，并为良好的信息融合提供重要的技术支持。

本领域技术人员知道，除了以纯计算机可读程序代码方式实现本发明提供的系统及其各个装置、模块、单元以外，完全可以通过将方法步骤进行逻辑编程来使得本发明提供的系统及其各个装置、模块、单元以逻辑门、开关、专用集成电路、可编程逻辑控制器以及嵌入式微控制器等的形式来实现相同功能。所以，本发明提供的系统及其各项装置、模块、单元可以被认为是一种硬件部件，而对其内包括的用于实现各种功能的装置、模块、单元也可以视为硬件部件内的结构；也可以将用于实现各种功能的装置、模块、单元视为既可以是实现方法的软件模块又可以是硬件部件内的结构。

以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变化或修改，这并不影响本发明的实质内容。在不冲突的情况下，本申请的实施例和实施例中的特征可以任意相互组合。