交互多模型的Rao—Blackwellized粒子滤波算法在多目标跟踪中的应用
张万里+何金刚
摘 要：在雷达信号处理领域，运动目标跟踪问题一直是研究的重点，而非线性滤波问题则是难点所在。本文以Rao-Blackwellized粒子滤波（RBPF）算法为基础，针对其应用于雷达目标跟踪时产生的算法收敛情况较差问题进行研究。为提高算法的收敛性能及滤波精度，对RBPF算法中线性部分采用交互多模型的方法进行状态预测和校正，利用数字仿真对更改后的算法进行验证，结果证明改进后的RBPF算法在各个方面均有较大程度改善。
关键词：雷达目标跟踪；制导信息估计；粒子滤波；RBPF算法；交互多模型
中图分类号：TP273 文献标识码：A 文章编号：1673-5048（2014）04-0003-05
0 引 言
在雷达信号处理领域，运动目标跟踪问题一直是长期以来的难点所在，实践证明在线性高斯系统中利用最小均方根误差准则进行目标状态估计的Kalman滤波方法[1]是最优的估计方法，但针对非线性非高斯系统，尽管采用基于局部线性化近似的扩展Kalman滤波（ExtendKF，EKF）[2]以及基于确定性采样的Unscented卡尔曼滤波（UnscentedKF，UKF）[3]方法可以解决一定形式的弱非线性、弱高斯条件下的目标跟踪问题，但由于其对系统模型的限制过强，在实际应用中大多无法满足。20世纪90年代出现了以粒子滤波（Particlefilter，PF）[4]为代表的非线性滤波方法，即利用蒙特卡罗采样得到的随机样本（也称为粒子）的加权和来近似状态的整个后验概率密度，其本质是采用蒙特卡罗仿真来获得高维积分的近似数值解，并解决各种估计问题。
粒子滤波面临的两个最大问题：一是粒子退化问题，即经过若干次迭代后，重要性权值可能集中到少数粒子上，这些粒子已经不能有效表达后验概率密度函数，为解决此问题，Gordon[5]等人提出了重采样方法，其思想是减少权值较小的粒子数，增加权值较大的粒子数。另一个问题是采用粒子数目过多导致计算的复杂度增加，当前的解决方法主要是从系统模型出发，利用模型自身的特性来提高滤波器性能。RaoBlackwellized方法[6]即将线性状态从系统中分离出来，利用Kalman滤波器对线性状态进行估计，利用粒子滤波（PF）对剩余的非线性状态进行估计，然后利用贝叶斯定理求取状态的后验概率。由于RBPF降低了粒子滤波状态的维数，与使用相同粒子数的传统PF算法相比，可以获得更优的性能。
当前的RaoBlackwellized粒子滤波（RBPF）中采用单一系统模型作为其中近似线性状态的估计，在跟踪机动目标时与真实飞行轨迹存在偏差。为解决此类问题，本文采用交互多模型的方法来获得状态的重要性分布，从而实现粒子状态的估计，并以仅有角度信息的雷达双目标跟踪问题为例，对改进的算法进行验证。
1 RaoBlackwellized粒子滤波算法
1.1 粒子滤波算法原理
解决目标跟踪问题的最优方法是贝叶斯滤波方法，它通过两个步骤来实现：状态预测和状态更新。贝叶斯滤波的实质是通过获得目标的后验概率密度，根据某些准则（如最大后验估计）近似地计算出目标状态值。定义系统模型如下：
其中：xk为目标在k时刻的状态，如目标的位置、速度、加速度等信息；yk为k时刻的测量值，如目标的位置、弹目距离、目标与传感器的相对角度等；p（xk|xk-1）为目标的动态模型，表征目标状态的动态变化情况；p（yk|xk）为系统的测量模型，表征目标在干扰情况下的测量变化情况。最优滤波的目的就是为了在已知观测信息的前提下获得目标的后验概率。
利用ChapmanKolmogoroff公式可得目标的后验概率密度为
上式从理论意义上提供了最优滤波问题的解决方法，但在非线性系统求解过程中无穷维积分的运算极为困难，无法得到其精确最优解。1.2 RaoBlackwellized粒子滤波算法流程
在RaoBlackwellized粒子滤波算法中，引入任意潜在变量λ，系统的动态模型和测量模型分别变为p（xk|xk-1，λk-1）和p（yk|xk，λk），已知重要性分布为π（λk|λ（i）1∶k-1，y1∶k），对当前粒子群{ω（i）k，λ（i）k，m（i）k，P（i）k：i=1，2，…，N}进行处理，其中m为均值，P为协方差，ω为粒子权重，N为粒子数。在k时刻，Rao-Blackwellized粒子滤波算法的流程如下：
（1）对粒子均值m和协方差P作卡尔曼滤波预测
k，P（i）k） （8）
利用RBPF算法可将多目标跟踪问题分为两个部分：多目标数据关联中后验概率分布的估计和基于数据关联单个目标跟踪的估计。可以分别通过序列重要性采样及Kalman滤波进行最小均方误差估计来解决，将跟踪过程简化为目标判别，即判别当前得到的测量值是目标还是杂波，并在此基础上对目标进行跟踪。
通过设定数据关联指标ck，当ck=j时表示当前测量值对应第j个目标，当ck=0时表示当前测量值经判别为杂波。
为使用RBPF滤波算法，必须首先确定一个重要性分布用以计算不同时刻k各个粒子的权值，即确定分布π（ck|y1∶k，c（i）1∶k-1），利用贝叶斯公式可以方便求取概率密度p（ck|y1∶k，c（i）1∶k-1），因此RBPF算法默认将p（ck|y1∶k，c（i）1∶k-1）作为最优的重要性分布π（ck|y1∶k，c（i）1∶k-1）来计算。 2 采用内部多模型的RBPF算法改进
当前RBPF中对数据关联后线性部分的处理一般采用EKF或UKF算法进行，但EKF或UKF算法多是针对于单一模型进行的，当系统状态在不同模型间发生变化时（如目标在平飞和机动间切换），利用RBPF算法中的线性部分采样单一的模型进行预测和校正往往会带来较大误差，为此对原始的RBPF算法进行改进，采样交互多模型的方法来进行辨识，判断目标当前是处于机动状态还是非机动状态，进而采取不同的模型来代入运算，具体过程如下：
（1）建立采用连续Wiener过程加速度模型（CWPA），系统状态为
Xk=（xk yk x·k y·k x¨k y¨k）T（11）
动力学模型为
（15）
其中：v（t）～N（0，σ2w）。此系统模型虽为非线性，但可采用EKF或UKF求解。
（3）对于每个粒子，分别采用两种模型对RBPF算法的线性部分进行预测及校正，求取其均值和方差，然后利用其每步的信息及误差方差的统计特性求取其对于每个模型的概率密度，即权重。
μik=f（vik，Sik） （16）
其中：μik为粒子权重；vik为测量信息；Sik为误差方差的统计特性。
（4）按照权重比例对计算的粒子均值及方差进行加权，作为下一次计算的初值，以此类推，可求取不同时刻的粒子状态。
mk=∑2
下面针对此问题分别采用线性Kalman滤波及RBPF粒子滤波来仿真，采用粒子数目为10。
图1所示为采用Kalman滤波的估计结果，可以看出，采用Kalman滤波得到的目标运动轨迹输出完全不能够跟踪上目标的真实运动轨迹，这是因为目标的观测模型中不只存在高斯噪声，而且在整个视场内存在均匀散布的杂波测量值，这样导致Kalman滤波算法很快失效。由图2～3可以看出，采用粒子滤波可以将真实目标轨迹与噪声
3.2 双目标跟踪的数字仿真实现
下面对仅有方位角测量信息的雷达双目标跟踪问题进行仿真验证，图5为雷达测量的示意图，此时目标的动态方程与上例中相同，但测量模型不同，此时测量量为角度值，使用两个固定位置的传感器对于两个目标进行测量，测量方程如下：
θik=arctan（yj，k-siy
xj，k-six）+rik（22）
其中：xj，k，yj，k为目标j的位置；six，siy为第i个传感器的位置；rik～N（0，σ2）为测量噪声，此时测量方程为非线性形式，因此需采用EKF或UKF配合使用RBPF算法。
图6为雷达的角度测量值随时间的变化情况，
从图中可以看出，针对两个传感器及两个目标可测量得到4组测量值，同时在视场范围内存在一定数量的杂波测量值。
针对当前双目标雷达跟踪的问题，采用单一Wiener过程加速度模型及机动转弯模型均能完成目标跟踪的过程，即滤波输出能够跟踪目标的运动轨迹，但与真实轨迹均有一定的偏离，按照本文交互多模型算法对RBPF的线性部分进行加权处理，不同时间段内选取不同的模型进行计算，同时按照加权结果生成新的粒子，得到的轨迹及误差如图7所示。图8给出采用交互多模型与采用单一模型的误差对比，从图中可见，采用交互多模型算法可以有效降低单一模型带来的估计误差。
4 结 论
本文首先给出Rao-Blackwellized粒子滤波算法的基本原理及算法流程，然后利用数字仿真验证了其在单目标跟踪问题中的应用效果以及相比Kalman滤波算法的优越性，针对其应用在多目标跟踪问题中存在的滤波局部不收敛的现象，采用交互多模型的方法对粒子的估值进行预测与更新，数字仿真验证可以看出，相比原始的Rao-Black wellized粒子滤波算法，更改后滤波算法的收敛性更好且跟踪精度更高，能够较好地完成杂波干扰下的双目标跟踪任务。
参考文献：
[1]KalmanRE.ANewApproachtoLinearFilteringandPredictionProblems[J].JournalofFluidsEngineering，1960，82（1）：35-45.
[2]SunaharaY.AnApproximateMethodofStateEstimation forNonlinearDynamicalSystems[J].JournalofFluids Engineering，1969，92（2）：385-393.
[3]JulierS，UhlmanJ，DurrantWhyteHF.ANewMethod fortheNonlinearTransformationofMeansandCovariances inFiltersandEstimators[J].IEEETransactionsonAutomaticControl，2000，45（3）：477-482.
[4]CappeO，GodsillSJ，MoulinesE.AnOverviewofExisting MethodsandRecentAdvancesinSequentialMonteCarlo[J].ProceedingsoftheIEEE，2007，95（5）：899-924.
[5]GordonN，SalmondD，SmithAFM.NovelApproachto Nonlinear/NonGaussianBayesianStateEstimation[J]. IEEProceedingsF（RadarandSignalProcessing），1993，140（2）：107-113.
[6]SarkkaS，VehtariA，LampinenJ.RaoBlackwellizedParticleFilterforMultipleTargetTracking[J].Information Fusion，2007，8（1）：2-15.
（1）建立采用连续Wiener过程加速度模型（CWPA），系统状态为
Xk=（xk yk x·k y·k x¨k y¨k）T（11）
动力学模型为
（15）
其中：v（t）～N（0，σ2w）。此系统模型虽为非线性，但可采用EKF或UKF求解。
（3）对于每个粒子，分别采用两种模型对RBPF算法的线性部分进行预测及校正，求取其均值和方差，然后利用其每步的信息及误差方差的统计特性求取其对于每个模型的概率密度，即权重。
μik=f（vik，Sik） （16）
其中：μik为粒子权重；vik为测量信息；Sik为误差方差的统计特性。
（4）按照权重比例对计算的粒子均值及方差进行加权，作为下一次计算的初值，以此类推，可求取不同时刻的粒子状态。
mk=∑2
下面针对此问题分别采用线性Kalman滤波及RBPF粒子滤波来仿真，采用粒子数目为10。
图1所示为采用Kalman滤波的估计结果，可以看出，采用Kalman滤波得到的目标运动轨迹输出完全不能够跟踪上目标的真实运动轨迹，这是因为目标的观测模型中不只存在高斯噪声，而且在整个视场内存在均匀散布的杂波测量值，这样导致Kalman滤波算法很快失效。由图2～3可以看出，采用粒子滤波可以将真实目标轨迹与噪声
3.2 双目标跟踪的数字仿真实现
下面对仅有方位角测量信息的雷达双目标跟踪问题进行仿真验证，图5为雷达测量的示意图，此时目标的动态方程与上例中相同，但测量模型不同，此时测量量为角度值，使用两个固定位置的传感器对于两个目标进行测量，测量方程如下：
θik=arctan（yj，k-siy
xj，k-six）+rik（22）
其中：xj，k，yj，k为目标j的位置；six，siy为第i个传感器的位置；rik～N（0，σ2）为测量噪声，此时测量方程为非线性形式，因此需采用EKF或UKF配合使用RBPF算法。
图6为雷达的角度测量值随时间的变化情况，
从图中可以看出，针对两个传感器及两个目标可测量得到4组测量值，同时在视场范围内存在一定数量的杂波测量值。
针对当前双目标雷达跟踪的问题，采用单一Wiener过程加速度模型及机动转弯模型均能完成目标跟踪的过程，即滤波输出能够跟踪目标的运动轨迹，但与真实轨迹均有一定的偏离，按照本文交互多模型算法对RBPF的线性部分进行加权处理，不同时间段内选取不同的模型进行计算，同时按照加权结果生成新的粒子，得到的轨迹及误差如图7所示。图8给出采用交互多模型与采用单一模型的误差对比，从图中可见，采用交互多模型算法可以有效降低单一模型带来的估计误差。
4 结 论
本文首先给出Rao-Blackwellized粒子滤波算法的基本原理及算法流程，然后利用数字仿真验证了其在单目标跟踪问题中的应用效果以及相比Kalman滤波算法的优越性，针对其应用在多目标跟踪问题中存在的滤波局部不收敛的现象，采用交互多模型的方法对粒子的估值进行预测与更新，数字仿真验证可以看出，相比原始的Rao-Black wellized粒子滤波算法，更改后滤波算法的收敛性更好且跟踪精度更高，能够较好地完成杂波干扰下的双目标跟踪任务。
参考文献：
[1]KalmanRE.ANewApproachtoLinearFilteringandPredictionProblems[J].JournalofFluidsEngineering，1960，82（1）：35-45.
[2]SunaharaY.AnApproximateMethodofStateEstimation forNonlinearDynamicalSystems[J].JournalofFluids Engineering，1969，92（2）：385-393.
[3]JulierS，UhlmanJ，DurrantWhyteHF.ANewMethod fortheNonlinearTransformationofMeansandCovariances inFiltersandEstimators[J].IEEETransactionsonAutomaticControl，2000，45（3）：477-482.
[4]CappeO，GodsillSJ，MoulinesE.AnOverviewofExisting MethodsandRecentAdvancesinSequentialMonteCarlo[J].ProceedingsoftheIEEE，2007，95（5）：899-924.
[5]GordonN，SalmondD，SmithAFM.NovelApproachto Nonlinear/NonGaussianBayesianStateEstimation[J]. IEEProceedingsF（RadarandSignalProcessing），1993，140（2）：107-113.
[6]SarkkaS，VehtariA，LampinenJ.RaoBlackwellizedParticleFilterforMultipleTargetTracking[J].Information Fusion，2007，8（1）：2-15.
（1）建立采用连续Wiener过程加速度模型（CWPA），系统状态为
Xk=（xk yk x·k y·k x¨k y¨k）T（11）
动力学模型为
（15）
其中：v（t）～N（0，σ2w）。此系统模型虽为非线性，但可采用EKF或UKF求解。
（3）对于每个粒子，分别采用两种模型对RBPF算法的线性部分进行预测及校正，求取其均值和方差，然后利用其每步的信息及误差方差的统计特性求取其对于每个模型的概率密度，即权重。
μik=f（vik，Sik） （16）
其中：μik为粒子权重；vik为测量信息；Sik为误差方差的统计特性。
（4）按照权重比例对计算的粒子均值及方差进行加权，作为下一次计算的初值，以此类推，可求取不同时刻的粒子状态。
mk=∑2
下面针对此问题分别采用线性Kalman滤波及RBPF粒子滤波来仿真，采用粒子数目为10。
图1所示为采用Kalman滤波的估计结果，可以看出，采用Kalman滤波得到的目标运动轨迹输出完全不能够跟踪上目标的真实运动轨迹，这是因为目标的观测模型中不只存在高斯噪声，而且在整个视场内存在均匀散布的杂波测量值，这样导致Kalman滤波算法很快失效。由图2～3可以看出，采用粒子滤波可以将真实目标轨迹与噪声
3.2 双目标跟踪的数字仿真实现
下面对仅有方位角测量信息的雷达双目标跟踪问题进行仿真验证，图5为雷达测量的示意图，此时目标的动态方程与上例中相同，但测量模型不同，此时测量量为角度值，使用两个固定位置的传感器对于两个目标进行测量，测量方程如下：
θik=arctan（yj，k-siy
xj，k-six）+rik（22）
其中：xj，k，yj，k为目标j的位置；six，siy为第i个传感器的位置；rik～N（0，σ2）为测量噪声，此时测量方程为非线性形式，因此需采用EKF或UKF配合使用RBPF算法。
图6为雷达的角度测量值随时间的变化情况，
从图中可以看出，针对两个传感器及两个目标可测量得到4组测量值，同时在视场范围内存在一定数量的杂波测量值。
针对当前双目标雷达跟踪的问题，采用单一Wiener过程加速度模型及机动转弯模型均能完成目标跟踪的过程，即滤波输出能够跟踪目标的运动轨迹，但与真实轨迹均有一定的偏离，按照本文交互多模型算法对RBPF的线性部分进行加权处理，不同时间段内选取不同的模型进行计算，同时按照加权结果生成新的粒子，得到的轨迹及误差如图7所示。图8给出采用交互多模型与采用单一模型的误差对比，从图中可见，采用交互多模型算法可以有效降低单一模型带来的估计误差。
4 结 论
本文首先给出Rao-Blackwellized粒子滤波算法的基本原理及算法流程，然后利用数字仿真验证了其在单目标跟踪问题中的应用效果以及相比Kalman滤波算法的优越性，针对其应用在多目标跟踪问题中存在的滤波局部不收敛的现象，采用交互多模型的方法对粒子的估值进行预测与更新，数字仿真验证可以看出，相比原始的Rao-Black wellized粒子滤波算法，更改后滤波算法的收敛性更好且跟踪精度更高，能够较好地完成杂波干扰下的双目标跟踪任务。
参考文献：
[1]KalmanRE.ANewApproachtoLinearFilteringandPredictionProblems[J].JournalofFluidsEngineering，1960，82（1）：35-45.
[2]SunaharaY.AnApproximateMethodofStateEstimation forNonlinearDynamicalSystems[J].JournalofFluids Engineering，1969，92（2）：385-393.
[3]JulierS，UhlmanJ，DurrantWhyteHF.ANewMethod fortheNonlinearTransformationofMeansandCovariances inFiltersandEstimators[J].IEEETransactionsonAutomaticControl，2000，45（3）：477-482.
[4]CappeO，GodsillSJ，MoulinesE.AnOverviewofExisting MethodsandRecentAdvancesinSequentialMonteCarlo[J].ProceedingsoftheIEEE，2007，95（5）：899-924.
[5]GordonN，SalmondD，SmithAFM.NovelApproachto Nonlinear/NonGaussianBayesianStateEstimation[J]. IEEProceedingsF（RadarandSignalProcessing），1993，140（2）：107-113.
[6]SarkkaS，VehtariA，LampinenJ.RaoBlackwellizedParticleFilterforMultipleTargetTracking[J].Information Fusion，2007，8（1）：2-15.
