一种基于熵惩罚的EM未知杂波估计的PHD多目标跟踪方法与流程

本发明属于目标跟踪技术领域，特别是未知杂波估计的PHD跟踪方法。

背景技术：

多目标跟踪(Multi-Target Tracking，MTT)是目标跟踪的研究热点，其在国防、军事、环境监测、安全监控、智能交通等领域具有广泛的应用价值。多目标跟踪方法主要有基于数据关联和随机有限集两类，前者的主要思想是通过数据关联建立观测信息与目标身份信息的匹配，从而将多目标跟踪问题转化为对单一目标状态的估计。但由于数据关联需要大量的计算，所以这类方法处理起来极为复杂。针对这一问题，Mahler将随机有限集理论(Random Finite Set，RFS)引入到多目标跟踪问题的求解中，提出了基于概率假设密度(Probability Hypothesis Density，PHD)滤波方法。PHD滤波器是RFS理论框架下多目标状态集合的后验概率密度的一阶统计矩，其实现主要分为适用于非线性非高斯情形下的序贯蒙特卡洛PHD(Sequential Monte Carlo PHD，SMC-PHD)滤波器和适用于线性高斯情形下的高斯混合PHD(Gaussian Mixture PHD，GM-PHD)滤波器。这类方法避免了数据关联，从而降低了计算难度，成为了近年来国内外研究热点。

在传统的PHD滤波器中，通常假定杂波的数目服从泊松分布，且在整个区域分布均匀。此种假定在杂波干扰较少的跟踪环境中不会有太大问题，但在复杂环境下，由于杂波个数和分布未知，这种假设使方法的鲁棒性和跟踪精度会严重降低。针对此类问题，目前的研究文献中有很多解决方法，但都存在着各种不足，林晓东等人于2011年提出了改进的基于概率假设密度滤波的多目标跟踪方法，该方法在目标状态中加入标签，除去权重较大标签代表的杂波干扰，此方法能更好地避免杂波干扰，但对新的目标不易发现。Lian F等人于2010年提出了未知杂波强度估计的PHD滤波方法，该方法运用有限混合模型和期望最大化(Expectation Maximization，EM)方法，较全面地解决了多目标跟踪中未知杂波估计的问题，但是杂波模型对初始值要求较敏感，不合理的初始值将影响估计精度。另外，传统EM方法计算量大，收敛速度慢，在复杂杂波背景下进行跟踪难以满足实时性的要求。

技术实现要素：

本发明所解决的技术问题在于提供一种基于熵惩罚EM的未知杂波估计的PHD(Entropy-penalized Expectation Maximization PHD，EPEM-PHD)滤波方法。

实现本发明目的的技术解决方案为：首先，提出了将带有熵惩罚因子的EM方法应用到PHD滤波器中，以加速PHD方法执行效率，提高了方法的实时性。其次，通过动态系数自适应调节使得加速混合模型分量消亡，在减少方法迭代次数的同时，又避免陷入局部最优解。此外，熵惩罚EM方法还具有对初始参数不敏感的优点。通过SMC-PHD滤波器的实现，保证了该方法在非线性情况下也有很好的应用。

本发明与现有技术相比具有以下优点：

1.针对多目标跟踪方法对强实时性的需求，首次提出了将带有熵惩罚因子的EM方法应用到PHD滤波器中，以加速PHD方法执行速度，使得改进后的PHD方法在多目标跟踪中具有实时性强的优势。

2.在熵惩罚EM方法中赋予熵惩罚因子合适的动态调节系数，通过动态系数自适应调节使得加速混合模型分量消亡，在减少方法迭代次数的同时，又避免陷入局部最优解。

3.提出了基于序贯Monte Carlo改进EM的未知杂波估计PHD方法，在完成杂波强度估计的同时，保证了在非线性情况下也能有很好的应用。

下面结合附图对本发明作进一步详细描述。

附图说明

图1为本发明方法的流程图。

图2为无线传感器网络协同跟踪仿真场景图中X坐标真实值与估计值对比。

图3为无线传感器网络协同跟踪仿真场景图中Y坐标真实值与估计值对比。

图4为不同方法跟踪误差对比图。

具体实施方式

依据附图，对本发明的技术方案作具体说明。

所述基于熵惩罚EM的未知杂波估计的PHD滤波方法，包括以下步骤：

步骤1进行基于熵惩罚EM方法的未知杂波估计。

采用随机有限集为杂波强度建模，杂波强度可表示为

k(z)＝λc(z) (1)

其中，λ为杂波数目，c(z)为杂波分布。

(1)杂波数目的估计

记H＝{h₁，h₂…h_k-1}，h_i表示i时刻的量测个数，i＝1，2，…，k-1。设杂波的运动为泊松点过程，即

其中，λ即为平均每帧的杂波个数。其无偏最小方差估计为

即为杂波数目的估计值。

(2)杂波分布的估计

采用有限混合模型描述杂波分布

其中，M为分量数目，θ_m＝{μ_m，φ_m}为分量均值μ_m与分量方差φ_m构成的分量m参数集合，θ＝{θ₁，…，θ_M，π₁，…，π_M}为混合参数集合，混合权重π_m表示混合模型中各分量的比例系数，并且

在有限混合模型的初始步骤中增加一个均匀分布分量，用于拟合来源于背景噪声的稀疏杂波分布，既可实现初始分量对目标区域较为全面的表示，同时又能达到缩减初始分量数目的目的。新的有限混合模型表示如下

由量测集合Z的n个元素相互独立，得到随机有限集合参数θ的极大似然函数为

采用EM方法求参数θ的极大似然估计。令EM方法的缺失参数集合为E＝{e₁，e₂，…，e_n}，其中e_j为M维向量，e_mj＝0或1，分别表示第j个测量属于或不属于混合模型中的第m分量，且满足

由完整数据集Y＝{z₁，z₂，…，z_n，E}，可得最大似然法求解的相关log似然函数为

由混合模型的各分量均为高斯分布，可得混合模型为

其中，d表示量测维数，δ(z，μ；φ)表示z与μ的马氏距离，e_mi为缺失参数，μ_m、φ_m、π_m为混合参数。

步骤2进行基于熵惩罚的改进EM的参数估计

未施加惩罚因子的EM方法通常步骤为

E-Step：求解参数e_mi，有

M-Step：将似然函数式(7)最大化，以获得混合参数π_m、μ_m和φ_m的估计值，有

混合权重

分量均值

分量方差

对参数e_mi与混合权重π_m施加惩罚因子，得到新的EM方法。具体为

(1)对混合权重π_m施加熵惩罚因子，得到新的M-Step。

混合权重π_m是迭代过程中分量删减的依据，π_m值越大，则该分量在下一次迭代中越可能被留存下来。而当π_m越小，直至小于预先设定的阈值ξ时，该分量就会被删减。取表示所有分量的信息熵。其中，-lnπ_m表示第m分量中一个数据点所散发的信息量，当时，则认为没有关于π_m的信息量。这里通过最大化W，即最小化信息熵实现更好更快地获得π_m的信息。在似然函数式(7)中添加惩罚因子以便对EM方法的迭加过程进行调节，可得

其中，α为自适应调节系数，且满足α∈(0，1)，有

采用拉格朗日乘子法，将式(5)代入式(15)，并对π_m求一阶偏导，有

如果则混合模型的第m分量在第t+1次迭代时，其混合权重将比没有施加惩罚因子时要小。按照迭代规则，越小，则减小的趋势越大。当减小到阈值ξ时，认为第m分量是无效的，将此分量进行删减。

分量数目M每减少一次，为保证约束条件和的成立，对和调整为

由式(17)、式(13)与式(14)形成新的M-Step。

(2)对缺失参数e_mi施加熵惩罚因子，得到新的E-Step。

以-ln e_mi表示量测z_i属于第m分量的信息量，则为z_i属于各分量的信息熵。当e_mi＝1/M，即z_i属于任何一个分量的概率相同时，则认为信息量不存在。对e_mi施加惩罚因子，有

其中，为关于的加权平均值。β_i为自适应调节系数，且满足β_i∈(0，1)，有

当时，将小于没有施加惩罚因子时的值，加快了量测信息向相关分量聚拢的速度，从而加速了分量数M的删减。

由式(20)形成施加熵惩罚后的新的E-Step。

利用改进后的E-Step和M-Step交替迭代直到收敛，即得到有限混合模型缺失参数e_mi和混合参数μ_m、φ_m和π_m的估计结果。

步骤3，得到基于熵惩罚EM的未知杂波估计的PHD滤波方法

根据式(10)和有限混合模型参数估计结果，可得杂波分布c(z)。

PHD算法中，SMC-PHD主要步骤为：

(1)粒子预测

对k-1时刻的粒子进行重要性采样，符号“～”表示随机抽样。则粒子权值为

其中，φ(x_k，x_k-1)＝e_k|k-1(x_k-1)f_k|k-1(x_k|x_k-1)+β(x_k|x_k-1)，为k时刻第i个目标的状态，J_k为出生目标粒子数目。Z_k＝{z_k，1，z_k，2，...，z_k，j，...，z_k，n(k)}，z_k，j表示k时刻第j个量测值，f_k|k-1(x_k|x_k-1)为状态转移密度，e_k|k-1(x_k)为目标生存概率，β(x_k|x_k-1)为源自状态x_k-1的衍生目标强度，γ_k(x_k)为出生目标强度；

(2)粒子更新

对于每个量测z∈Z_k，计算目标量测强度C_k(z)，并且更新预测步的粒子权值

(3)重采样

第1步：对粒子权值求和

第2步：对重采样得

第3步：与相乘得PHD近似即是k时刻目标数估计值

将杂波分布c(z)代入式(1)可得估计值k(z)，将k(z)代入式(23)和式(24)，即可得到EPEM-PHD滤波器。

下面对本发明的方法进行仿真验证：

仿真实验采用序贯Monte Carlo概率假设密度滤波器(SMC-PHD)与所提出的EPEM-PHD比较跟踪性能。仿真实验在AMD FX-7500 2.10GHz处理器、1024MB内存的PC上，使用MATLAB R2012b平台实现。模拟一个[-120，120]×[-100，100]m的x-y二维平面区域，先后出现4个作无规则运动的目标。目标的运动模型为

x_k＝F·x_k-1+G·ω_k (22)

其中，x_k-1＝[x_x，v_x，x_y，v_y]^T为目标状态，由x，y方向的目标位置及速度[x_x，v_x]^T和[x_y，v_y]^T组成；F为状态转移矩阵，G为噪声输入矩阵，ω_k为过程噪声，取高斯白噪声ω_k～N(0，Q)，F与G分别表示为

目标量测方程Z_k为

其中，{V_1，k}和{V_2，k}为相互独立的零均值高斯白噪声。系统采样周期T＝1s，仿真时间步长为40。出生目标强度为γ_k＝0.2N(·；x，Q)，其正态分布的均值和方差为

目标存活概率P_s(x)取0.95，目标检测概率P_D(x)设为1。令平均每帧杂波个数的真实值为λ＝50；杂波在设定区域内随机产生。每个目标赋予1000个粒子，随机生成一个随机数作为每帧杂波个数。

实验中所选择的四个目标之中，目标1在第1秒出现，终止于第7秒；目标2在第8秒出现，终止于第25秒；目标3在第12秒出现，终止于第37秒；目标4在第26秒出现，终止于第40秒。实验结果如图2至图4所示。

图2和图3显示了单次Monte Carlo仿真中的X坐标和Y坐标的跟踪效果，从图中可以看出，改进方法能够正确估计多目标状态，且优于PHD方法。跟踪第三个目标在第25秒至第33秒时PHD方法偏离了目标的真实位置，但EPEM-PHD却很好地进行目标跟踪，性能提高的原因在于进行了杂波强度的在线估计。采用Wasserstein距离评价多目标跟踪性能，其表达形式为：

其中，为多目标状态的估计值和真实值，C^ij为正实数且p取2。

图4显示了通过50次Monte Carlo仿真后，本发明方法和传统PHD多目标跟踪方法的Wasserstein距离。从图中可以看出，本发明方法的Wasserstein距离较小，说明EPEM-PHD方法在多目标跟踪精度上的准确性高于PHD，并且能够在杂波未知的环境下，使得估计结果更接近于真实的分布，具有更好的稳定性和鲁棒性，同时也证明了本发明方法的有效性。

同时，系统通过50次Monte Carlo仿真后由CPU平均运算时间得到，PHD方法的平均耗时为7.24秒，而EPEM-PHD平均耗时为7.37秒，耗时增幅为1.80％，是因为对未知杂波进行了估计。增幅并不显著，所以本发明方法保持了较好的实时性。

综上所述，本发明公布的一种基于熵惩罚的EM未知杂波估计的PHD多目标跟踪方法在杂波强度未知的环境下，与传统的PHD滤波方法相比，保持了较好的实时性，具有跟踪精度高且稳定的优势。本发明的基于熵惩罚的EM未知杂波估计的PHD多目标跟踪方法在处理杂波未知时的多目标跟踪方面具有积极的意义。