CN104268567A - 一种观测数据聚类划分的扩展目标跟踪方法 - Google Patents

Abstract

一种观测数据聚类划分的扩展目标跟踪方法。其包括：1)获得目标数的估计值；2)观测数据椭球门限处理；3)观测数据K-means聚类划分；4)概率假设密度滤波扩展目标跟踪。本发明提供的观测数据聚类划分的扩展目标跟踪方法首先利用极大似然估计方法获得目标数的估计值，利用椭球门限对观测数据进行处理以降低杂波观测对聚类划分的影响后，利用K-means对观测数据集合进行聚类划分，进而得到观测数据聚类划分的扩展目标跟踪方法。本发明方法能够快速获得正确的划分观测集合，计算复杂度较距离门限划分方法有较大降低，同时在多目标跟踪效果方面优于观测数据距离门限划分的概率假设密度滤波扩展目标跟踪算法。

Description

一种观测数据聚类划分的扩展目标跟踪方法

技术领域

本发明属于多扩展目标跟踪技术领域，特别是涉及一种观测数据聚类划分的扩展目标跟踪方法。

背景技术

多目标跟踪是多传感器数据融合领域的一个重要理论和实践问题。最近几十年来，国内外众多专家学者对它进行了深入的研究，取得了丰硕的成果。目前常用的多目标跟踪方法主要分为两类。一类是传统的基于数据关联的多目标跟踪方法，另一类是非关联的多目标跟踪算法，主要是基于随机有限集理论的概率假设密度滤波方法。无论是传统的基于数据关联的多目标跟踪方法还是基于随机有限集理论的方法，绝大多数方法都需要满足这样一个假设：每个观测都来自于一个可以被看作质点的目标。然而，由于传感器分辨率、目标存在多个反射点以及传感器与目标相对位置关系等原因，这一假设并不一定能成立，从而出现了由于传感器分辨率不断提高或者目标距离传感器较近，单个目标的不同反射点同时产生多个观测形成的扩展目标问题。在随机有限集框架下，2003年Mahler通过传递目标状态的后验概率密度的一阶统计矩的方法得到了概率假设密度(Probability hypothesis density，PHD)滤波器。2006年，Vo提出了适合线性高斯条件的高斯混合PHD(Gaussian Mixture PHD,GM-PHD)滤波器。2009年，Mahler提出了概率假设密度滤波算法解决扩展目标的跟踪问题，并给出了相关算法理论框架，即ET-PHD(Extended Target PHD)算法。2012年，Granstrom提出了基于上述框架对于扩展目标的线性条件下的高斯混合模型的实现方法，即ET-GM-PHD算法。Mahler提出的ET-PHD扩展目标跟踪算法理论框架中，需要将观测数据(包括目标及杂波)划分为多种可能的集合，然后利用观测数据划分的所有可能组合进行目标状态更新计算，这种方法对于目标数较多和杂波较多的情况几乎是无法计算的。Granstrom针对这个问题提出了一种利用距离门限实现观测数据集合划分的方法，但距离门限范围的设置依赖于实验数据，导致算法性能与观测数据集合划分距离范围参数选取情况密切相关，跟踪估计性能很不稳定。

发明内容

为了解决上述问题，本发明的目的在于提供一种能够快速获得正确的划分观测集合，计算复杂度较距离门限划分方法低，同时在多目标跟踪效果方面优于距离门限划分方法的概率假设密度滤波扩展目标跟踪方法。

为了达到上述目的，本发明提供的观测数据聚类划分的扩展目标跟踪方法包括按顺序进行的下列步骤：

1)获得目标数的估计值的S1阶段；

2)观测数据椭球门限处理的S2阶段；

3)观测数据K-means聚类划分的S3阶段；

4)概率假设密度滤波扩展目标跟踪的S4阶段；

在步骤1)中，所述的获得目标数的估计值的方法是利用极大似然估计方法从观测数据中获得目标数的估计值，为观测数据K-means聚类划分做准备。

在步骤2)中，所述的观测数据椭球门限处理的方法是利用椭球门限对进行更新计算的观测数据进行门限处理以减少观测数据中的杂波成分，为下一步观测数据K-means聚类划分做准备。

在步骤3)中，所述的观测数据K-means聚类划分的方法是将S1阶段中获得的目标数的估计值作为聚类数，利用K-means对S2阶段中椭球门限处理后的观测数据集合进行聚类划分，为下一步概率假设密度滤波扩展目标跟踪做准备。

在步骤4)中，所述的概率假设密度滤波扩展目标跟踪的方法是基于S3阶段中K-means聚类划分后的观测数据对扩展目标进行跟踪滤波。本发明提供的一种观测数据聚类划分的扩展目标跟踪方法是利用极大似然估计的思想，获取目标数的估计值，并采用K-means聚类算法对观测数据集合划分，进而得到概率假设密度滤波的扩展目标跟踪结果。本发明方法能够获得正确的观测集合划分，计算复杂度较距离门限划分方法有较大降低，同时在多目标跟踪效果方面优于观测数据距离门限划分的概率假设密度滤波扩展目标跟踪算法。

附图说明

图1为本发明提供的观测数据聚类划分的扩展目标跟踪方法流程图。

图2为目标的运动轨迹及观测值。

图3为目标数估计结果随时间变化曲线图。

图4为扩展目标状态估计随时间变化曲线图。

图5为OSPA距离随时间变化曲线图。

具体实施方式

下面结合附图和具体实施例对本发明提供的观测数据聚类划分的扩展目标跟踪方法进行详细说明。

图1为本发明提供的观测数据聚类划分的扩展目标跟踪方法流程图。其中的全部操作都是在以计算机为核心的计算机系统中完成的，操作的主体均为计算机系统。

如图1所示，本发明提供的观测数据聚类划分的扩展目标跟踪方法包括按顺序进行的下列步骤：

1)获得目标数的估计值的S1阶段：本阶段利用极大似然估计方法从观测数据中获得目标数的估计值，为观测数据K-means聚类划分做准备，然后进入下一步S2阶段；

在此阶段中，每个周期获得的观测数据中包括目标的扩展观测数以及杂波观测数，每个目标的扩展观测数服从γ为均值的Poisson分布，杂波观测数服从λ为均值的Poisson分布。由于多个独立的Poisson分布和仍然服从Poisson分布，因此每一个周期获得的观测数服从均值为γN_X,k+λ的Poisson分布，即：

$p\left(M_k\right|(\mathrm{\γ}{N}_{x,k}+\mathrm{\λ}))=\frac{{(\mathrm{\γ}{N}_{x,k}+\mathrm{\λ})}^{M_k}}{M_k!}{e}^{-(\mathrm{\γ}{N}_{x,k}+\mathrm{\λ})}---\left(1\right)$

式中，N_x,k表示k时刻的目标数，M_k表示k时刻的观测数。

每个目标的扩展观测数均值γ已知，杂波观测数均值λ已知，对于k时刻获得的M_k个观测数据可通过极大似然估计方法获得目标数N_X,k的估计值。

${\hat{N}}_{X,k}=\arg \max p\left(M_k\right|(\mathrm{\γ}{N}_{x,k}+\mathrm{\λ}))---\left(2\right)$

2)观测数据椭球门限处理的S2阶段：本阶段利用椭球门限对进行更新计算的观测数据进行门限处理，以减少观测数据中的杂波成分，为下一步观测数据K-means聚类划分做准备，然后进入下一步S3阶段；

在此阶段中，为了降低杂波观测数对观测数据聚类划分的影响，首先利用椭球门限对进行更新计算的观测数据进行门限处理，以减少观测数据中的杂波成分。在k时刻预测步计算后，第i个观测值与第j个状态预测值(高斯项均值)的残差向量可以表示为：

${\mathrm{\ϵ}}^{\left(\mathrm{ji}\right)}=z_k^{\left(i\right)}-H_k{m}_{k|k-1}^{\left(j\right)}---\left(3\right)$

残差协方差矩阵可以表示为：

$S_k^{\left(j\right)}=H_k{P}_k^{\left(j\right)}{H}_k^T+R_k---\left(4\right)$

其中，H_k为测量矩阵，为第j个预测协方差矩阵，R_k为观测噪声矩阵。

利用椭球门限判别对观测数据进行门限处理，以减少观测数据中的杂波成分：

${\mathrm{\ϵ}}^{\left(\left(\mathrm{ji}\right)\right)T}{\left(S_k^{\left(j\right)}\right)}^{-1}{\mathrm{\ϵ}}^{\left(\mathrm{ji}\right)}\≤\mathrm{Tg}---\left(5\right)$

当观测空间的维数为2时，椭球门限阈值可用下式表示：

Tg＝-2ln(1-P_g)    (6)

其中，P_g为目标观测数据落入椭球门限的概率。

对观测数据进行椭球门限处理后，距离目标观测较近的杂波被剔除掉，从而降低了杂波观测数对观测数据聚类划分的影响。

3)观测数据K-means聚类划分的S3阶段：本阶段将S1阶段获得的目标数的估计值作为聚类数，利用K-means对S2阶段中椭球门限处理后的观测数据集合进行聚类划分，然后进入下一步S4阶段；

在此阶段中，所述的观测数据K-means聚类划分的具体方法如下：首先将目标数的估计值设定为每个周期的划分聚类数，然后使用K-means聚类算法通过相似度准则将观测数据集合划分为个观测聚类子集。

4)概率假设密度滤波扩展目标跟踪的S4阶段：基于S3阶段中观测数据K-means聚类划分结果对扩展目标进行跟踪滤波；

在此阶段中，所述的概率假设密度滤波扩展目标跟踪的具体方法如下：使用线性高斯条件下的扩展目标(ET-GM-PHD)滤波器对S3阶段中划分的聚类观测数据进行跟踪滤波。

针对扩展目标的ET-GM-PHD算法的递推公式可以分为预测步和更新步两个步骤。其中预测步与传统点目标GM-PHD方法相同即：

预测步：

${\mathrm{\υ}}_{k|k-1}\left(x\right)=\underset{j=1}{\overset{J_{k|k-1}}{\mathrm{\Σ}}}{\mathrm{\ω}}_{k|k-1}^{\left(j\right)}N(x,m_{k|k-1}^{\left(j\right)},P_{k|k-1}^{\left(j\right)})---\left(7\right)$

其中，J_k|k-1为目标强度函数高斯项个数；是目标强度函数高斯项权值；和是目标强度函数高斯项的均值和方差；N(x；m,P)表示高斯分布。

更新步：

${\mathrm{\υ}}_{k|k}\left(x\right)={\mathrm{\υ}}_{k|k}^{\mathrm{ND}}\left(x\right)+\underset{p\∠Z_k}{\mathrm{\Σ}}\underset{W\∈p}{\mathrm{\Σ}}{\mathrm{\υ}}_{k|k}^D(x,W)---\left(8\right)$

表示传感器未探测到的目标的高斯形式强度函数。

$\begin{array}{c}{\mathrm{\υ}}_{k|k}^{\mathrm{ND}}\left(x\right)=\underset{j=1}{\overset{J_{k|k-1}}{\mathrm{\Σ}}}{\mathrm{\ω}}_{k|k}^{\left(j\right)}N(x,m_{k|k}^{\left(j\right)},P_{k|k}^{\left(j\right)})\\ {\mathrm{\ω}}_{k|k}^{\left(j\right)}=(1-(1-e^{{\mathrm{\γ}}^{\left(i\right)}})){P_D}^{\left(i\right)}{\mathrm{\ω}}_{k|k-1}^{\left(j\right)}\\ m_{k|k}^{\left(j\right)}=m_{k|k-1}^{\left(j\right)}\\ P_{k|k}^{\left(j\right)}=P_{k|k-1}^{\left(j\right)}\end{array}---\left(9\right)$

其中γ(i)，分别为的缩写，是线性高斯条件下的近似。一般情况下假设扩展观测数均值为常数即γ⁽ⁱ⁾＝γ，目标检测概率也为常数，即 $P_D^{\left(i\right)}=P_D.$

表示传感器探测到的目标的高斯形式强度函数。

$\begin{array}{c}{\mathrm{\υ}}_{k|k}^D(x,W)=\underset{j=1}{\overset{J_{k|k-1}}{\mathrm{\Σ}}}{\mathrm{\ω}}_{k|k}^{\left(j\right)}N(x,m_{k|k}^{\left(j\right)},P_{k|k}^{\left(j\right)})\\ {\mathrm{\Γ}}^{\left(i\right)}=e^{{-\mathrm{\γ}}^{\left(i\right)}}{\left({\mathrm{\γ}}^{\left(i\right)}\right)}^{\left|W\right|}\\ {\mathrm{\φ}}_W^{\left(i\right)}=N(z_W;H_W{m}_{k|k}^{\left(j\right)},H_W{P}_{k|k}^{\left(j\right)}{H}_W^T+R_W)\\ {\mathrm{\Φ}}_W^{\left(i\right)}={\mathrm{\φ}}_W^{\left(i\right)}\underset{z_k\∈W}{\mathrm{\Π}}\frac{1}{{\mathrm{\λ}}_k{c}_k\left(z_k\right)}\end{array}---\left(10\right)$

$\begin{array}{c}{m}_{k|k}^{\left(j\right)}=m_{k|k-1}^{\left(j\right)}+K_k^{\left(j\right)}(z_W-H_W{m}_{k|k-1}^{\left(j\right)})\\ P_{k|k}^{\left(j\right)}=(I-K_k^{\left(j\right)}{H}_W)P_{k|k-1}^{\left(j\right)}\\ K_k^{\left(j\right)}=P_{k|k-1}^{\left(j\right)}{H}_W^T{(H_W{P}_{k|k-1}^{\left(j\right)}{H}_W^T+R_W)}^{-1}\end{array}---\left(11\right)$

其中，观测z_k属于观测子集W，集合的势为|W|，即估计的目标数N_X,k；每个周期观测数据中杂波观测数服从Poisson分布，其均值为λ_k，一般情况下假设杂波观测数均值为常数，即λ_k＝λ；c_k(z_k)是观测空间中杂波的分布函数，假设其服从均匀分布，c_k(z_k)＝1/V，V是观测空间。

其中，H_W，R_W分别表示扩维后的观测矩阵及观测噪声协方差矩阵；blkdiag(·)表示块对角矩阵计算。

$\begin{array}{c}{d}_W={\mathrm{\δ}}_{\left|W\right|,1}+\underset{l=1}{\overset{J_{k|k-1}}{\mathrm{\Σ}}}{\mathrm{\Γ}}^{\left(l\right)}{P}_D^{\left(l\right)}{\mathrm{\Φ}}_W^l{\mathrm{\ω}}_{k|k-1}^{\left(l\right)}\\ {\mathrm{\ω}}_p=\frac{\underset{W\∈p}{\mathrm{\Π}}{d}_W}{\underset{p^{\′}\∠Z_k}{\mathrm{\Σ}}\underset{W^{\′}\∈p^{\′}}{\mathrm{\Π}}{d}_{W^{\′}}}\\ {\mathrm{\ω}}_{k|k}^{\left(j\right)}={\mathrm{\ω}}_p\frac{{\mathrm{\Γ}}^{\left(i\right)}{P}_D^{\left(i\right)}}{d_W}{\mathrm{\Φ}}_W^{\left(i\right)}{\mathrm{\ω}}_{k|k-1}^{\left(j\right)}\end{array}---\left(13\right)$

其中，δ_|W|,1为克罗内克(Kronecker delta)函数；p'表示观测集合划分的每一种可能；W'表示每一个观测集合划分可能中的观测子集；ω_p表示当前观测划分可能的权重。

本发明提供的观测数据聚类划分的扩展目标跟踪方法的效果可以通过以下仿真结果进一步说明。

仿真数据描述：实验设置一个[-100,100]×[-100,100]的监视区域，4个目标做匀速直线运动，采样周期T＝1s，仿真40步，目标状态向量为[x,x',y,y']，其中(x,y)表示目标的位置，(x',y')表示目标的速度，观测空间维数n_z＝2，目标运动方程及观测方程满足线性高斯模型，检测概率P_D＝0.9，目标存活概率P_s＝0.99，杂波强度λ＝10，服从Poisson分布，杂波位置均匀分布在[-100,100]×[-100,100]的观测空间。每个目标扩展观测数服从γ＝5的Poisson分布。新生目标的概率假设密度函数为：

${\mathrm{\γ}}_k\left(x\right)=0.1N(x;m_{\mathrm{\γ}}^{\left(1\right)},P_{\mathrm{\γ}})+0.1N(x;m_{\mathrm{\γ}}^{\left(2\right)},P_{\mathrm{\γ}})+0.1N(x;m_{\mathrm{\γ}}^{\left(3\right)},P_{\mathrm{\γ}})+0.1N(x;m_{\mathrm{\γ}}^{\left(4\right)},P_{\mathrm{\γ}})$

其中, $m_{\mathrm{\γ}}^{\left(1\right)}={[-\mathrm{10,0,20,0}]}^T,m_{\mathrm{\γ}}^{\left(2\right)}={[\mathrm{20,0},-\mathrm{20,0}]}^T,m_{\mathrm{\γ}}^{\left(3\right)}={[-\mathrm{40,0,20,0}]}^T,$ P_γ＝[5,1,5,1]^T。处理观测数据的椭球门限参数P_g＝0.9。高斯项修剪阈值T＝10^-5，合并门限U＝4，高斯项最大个数J_max＝200，误差调节因子c＝20。进行100次Monte Carlo实验。

图2为4个目标的运动轨迹图，‘—o—’表示目标的观测值，‘——’表示目标真实的状态值。图3为目标状态估计随时间变化曲线图，‘—o—’表示目标的估计值，‘——’表示目标真实的状态值。

图4为目标数估计随时间变化曲线图，图5为OSPA距离随时间变化曲线图。如图所示，‘—*—’为采用本发明方法估计结果，‘—△—’为采用距离门限划分方法的估计结果，‘——’为目标数真实值。本发明提供的观测数据聚类划分的扩展目标跟踪方法的目标数估计结果及OSPA距离都要优于距离门限划分的扩展目标跟踪方法，说明本发明方法获得了较距离门限划分方法更准确的扩展目标观测集合划分。在相同运行环境下两种方法的平均处理时间如表1所示，由于本文方法利用K-means方法进行观测数据集合划分，其计算复杂度为O(M_k·n·k)，与观测数M_k成正比，k为迭代次数；而距离门限划分方法进行观测数据集合划分时其计算复杂度为与观测数M_k的4次方成正比。可以看出本发明方法远快于距离门限划分方法。

表1

Claims (5)

1.一种观测数据聚类划分的扩展目标跟踪方法，其特征在于：所述的观测数据聚类划分的扩展目标跟踪方法包括按顺序进行的下列步骤： 

1)获得目标数的估计值的S1阶段； 

2)观测数据椭球门限处理的S2阶段； 

3)观测数据K-means聚类划分的S3阶段； 

4)概率假设密度滤波扩展目标跟踪的S4阶段。 

2.根据权利要求1所述的一种观测数据聚类划分的扩展目标跟踪方法，其特征在于：在步骤1)中，所述的获得目标数的估计值的方法是利用极大似然估计方法从观测数据中获得目标数的估计值，为观测数据K-means聚类划分做准备。 

3.根据权利要求1所述的观测数据聚类划分的扩展目标跟踪方法，其特征在于：在步骤2)中，所述的观测数据椭球门限处理的方法是利用椭球门限对进行更新计算的观测数据进行门限处理，以减少观测数据中的杂波成分，为下一步观测数据K-means聚类划分做准备。 

4.根据权利要求1所述的观测数据聚类划分的扩展目标跟踪方法，其特征在于：在步骤3)中，所述的观测数据K-means聚类划分的方法是将S1阶段中获得的目标数的估计值作为聚类数，利用K-means对S2阶段中椭球门限处理后的观测数据集合进行聚类划分，为下一步概率假设密度滤波扩展目标跟踪做准备。 

5.根据权利要求1所述的观测数据聚类划分的扩展目标跟踪方法，其特征在于：在步骤4)中，所述的概率假设密度滤波扩展目标跟踪的方法是基于S3阶段中K-means聚类划分后的观测数据对扩展目标进行跟踪滤波。