基于高斯粒子滤波算法的改进及应用

孙翌晨 李军

摘 要：针对粒子滤波存在的样本贫化现象，提出了一种优化重选样本粒子的粒子滤波算法。这种方法在引入最新量测后将状态后验概率密度逼近为一个高斯分布，在粒子贫化问题逐渐凸显后，通过该分布重新采集粒子后再进行运算，有效缓解了传统方法在粒子贫化后出现的滤波精度下降的问题。仿真结果表明，新的粒子滤波算法有更高的滤波精度和运行效率。

关键词：粒子滤波；后验概率；粒子贫化；重新选取

近年来，粒子滤波在目标跟踪领域得到了越来越广泛的应用。常见的非线性滤波方法，如扩展卡尔曼滤波（Extended Kalman Filter，EKF），无迹卡尔曼滤波（Unscented Kalman Filter，UKF）都是针对非线性系统的线性卡尔曼滤波方法的变形与改进，因此使用条件也受到卡尔曼滤波算法的条件限制[1]。而粒子滤波算法通过蒙特卡罗仿真手段产生大量粒子，随着采样粒子数不断增大，其散布情况将逐渐逼近状态的后验概率密度。粒子滤波在解决非高斯分布系统问题上具有明显的优势，可以说它是目前非高斯非线性系统状态估计的“最优”滤波器[2]。但是，随着时间的递推，会出现粒子的退化问题。通常，有两种方法可以减轻粒子退化问题：一是增加重采样环节；二是选择合适的重要密度函数进行更有效的采样[3-5]。常规的重采样方法随着迭代次数的增加，会出现粒子贫化问题，为此，人们提出了许多不同的方法来解决这个问题，如高斯粒子滤波算法（Gaussian particle filter），重采样粒子移动算法（Resample-Move Alogrithm）[6]，增加马尔可夫链蒙特卡罗（Markov Chain Monte Carlo）移动步骤[7-9]，对粒子进行正则（Regula—rization）重采样[10]。

笔者将标准粒子滤波算法和高斯粒子滤波相结合，引入一个重新选择粒子的过程，即粒子优化重选粒子滤波（Optimized Repicking Particle Filter）。当粒子贫化问题出现时，通过由当前滤波值和方差生成的高斯密度函数来逼近后验概率分布，重新采样后再進行运算。由于引入粒子优化重选的过程，缓解了粒子的退化与贫化问题，所以滤波精度要高于标准算法。最后的仿真结果也证明了这一点。

1 高斯粒子滤波算法

GPF通过由高斯密度函数来逼近后验概率分布，其基本思想是利用描述k时刻目标状态xk的后验概率分布 ，是对应权值为的粒子集，其中是0到k时刻的状态集。权值被归一化为 。由此，k时刻的后验概率密度可表示为

根据文献[11]，选择=为IDF可使权值的方差最小化，但通常情况下很难求得的表达式，因此一种简单常用的替代方案是选择先验作为重要密度函数。因此，可将重要性权值写为：

1.1 量测更新

当接受到第k个观测值zk之后，可以利用样本及其权值来计算滤波值和方差pk，将状态后验概率密度逼近为一个高斯分布，可表示为

其中

1.2 时间更新

由于已经将近似为高斯函数，故可将状态预测概率密度近似为：

其蒙特卡罗逼近为

式中，是从 采样得到的粒子，状态预测概率密度函数的均值和方差可表示为

2 样本优化重选粒子滤波

与一般的粒子滤波相比，GPF不需要对粒子进行重采样，降低了计算量和复杂度，但是滤波精度较之一般的粒子滤波稍差，因而考虑将两者结合。

由此，引入样本优化重选粒子滤波算法ORPF。在粒子贫化问题出现后，引入高斯粒子滤波的思想将后验概率密度逼近为一个高斯分布。由通过蒙特卡罗仿真得到粒子集各个粒子的权值为1/N。进而，在的基础上通过状态方程产生新的粒子集 进行K+1时刻的计算。这样，后一时刻的粒子根据前一时刻的新采样的粒子集更新而来，这样就可以保留粒子的多样性。所以，当在出现粒子贫化问题时，可以采用这种方法来解决这个问题。

新算法的描述如下：

（1）由p（x0）得到N个采样点 。

（2）计算权值 并对其归一化，得到

（3）输出K时刻的状态估计为：

方差为

（4）判断是否需要重采样：计算的值，如果Neff

（5）进行重采样，将原来的带权样本映射为等权样本

（6）判断是否要重新采集粒子，设重采样后被复制次最多的粒子繁殖出了a个子代，当a>0.2N时则认为粒子多样性丧失，此时进入步骤（7），否则，跳入步骤（8）。

（7）根据通过蒙特卡罗采样得到粒子集

（8）通过状态方程产生新粒子集

（9）重复步骤（2）～（8）。

3 实验仿真结果及分析

3.1 仿真一

选用单变量非静态增长模型（UNGM），仿真对象的过程模型和量测模型如下。

过程模型：

量测模型：

式中，ω（t）和v（t）为零均值高斯噪声。仿真取噪声方差Q为10，量测噪声方差R为1。分别用增加马尔可夫链蒙特卡罗粒子滤波（MCMCPF）和样本优化重选粒子滤波算法（ORPF）来对这样的非线性系统进行状态估计和跟踪。时间步长50次。

其中均方根误差

具体数据如表1所示：

状态估计曲线如图1、2所示：

如表1所示，过程噪声方差Q=10、量测噪声方差R=1时。在N=100，MCMCPF的误差均值为3.6946，而ORPF为3.5602；MCMCPF的RMSE为5.2191，而ORPF为4.0066；在状态估计时间上，MCMCPF为0.009502s，而ORPF为0.010143s。在N=500时，MCMCPF的误差均值为2.6765，而ORPF为2.4154；MCMCPF的RMSE为3.6747，而ORPF为3.2293；在时间估计上，MCMCPF为0.085130s，而ORPF为0.028229s。由此，说明ORPF的估计精度比MCMCPF略高，在粒子数较少时，两种算法的估计时间相差不大，但随着粒子数的增加，MCMCPF的估计时间急剧上升，远超ORPF。

如圖1、2、3、4所示，ORPF算法的估计误差是略小于MCMCPF的，滤波精度是略优于MCMCPF的。

综上，可以得出，在高度非线性模型下，ORPF在滤波精度上略优于MCMCPF，在实时性上要优于MCMCPF。

3.2 仿真二

选用被动定位系统中二维纯方位目标跟踪模型，仿真对象的过程模型和量测模型如下：

过程模型：

量测模型：

其中为系统k时刻的状态值，即目标在x和y轴上的位置和速度，为k-1时刻x，y方向的系统噪声，为k时刻的观测噪声。为k时刻的观测角度。

如图5、6所示，ORPF在X轴方向和Y轴方向上的估计误差都小于PF的估计误差，说明ORPF的估计精度高于PF。图7展示了2种算法的状态估计曲线，显然ORPF的跟踪精度要好。

4 结束语

本文提出了一种基于GPF的加入了粒子优化重选步骤的新的粒子滤波。传统的重采样方法中权值较大的粒子会被多次复制，造成粒子的贫化。在贫化问题出现时，利用当前时刻的估计值和方差构建高斯密度函数来逼近后验概率分布，从中重新采样并更新粒子。该方法能够有效的缓解粒子贫化带来的问题。仿真结果表明，提出的ORPF在估计性能上要好于传统的PF。

[参考文献]

[1]占荣辉，张军.非线性滤波理论与目标跟踪应用[M].北京：国防工业出版社，2013：108-167.

[2]朱志宇.粒子滤波算法及其应用[M].北京：科学出版社，2010.

[3]胡士强，敬忠良.粒子滤波算法综述[J].2005，20（4）：361-365.

[4]杨小军，潘泉，王睿，等.粒子滤波进展与展望[J].控制理论与应用， 2006，23（2）：261-267.

[5]张俊根.粒子滤波及其在目标跟踪中的应用研究[D].西安电子科技大学，2011.

[6]Gilks W R，Berzuini C.Following a moving target—Monte Carlo inference for dynamic Bayesian models[J].Journal of the Royal Statistical Society：Series B （Statistical Methodology）， 2001，63（1）：127-146.

[7]Fearnhead P.Markov chain Monte Carlo，sufficient statistics， and particle filters[J].Journal of Computational and Graphical Statistics，2002，11（4）：848-862.

[8]Carlin B P，Polson N G，Stoffer D S.A Monte Carlo approach to nonnormal and nonlinear state-space modeling[J].Journal of the American Statistical Association，1992，87（418）：493-500.

[9]Berzuini C，Best N G，Gilks W R，et al.Dynamic conditional independence models and Markov chain Monte Carlo methods[J]. Journal of the American Statistical Association，1997， 92（440）： 1403-1412.

[10]Sarkar P.Sequential Monte Carlo Methods in Practice[J]. Technometrics，2003，45（1）：106-106.

[11]Sanjeev Arulampalam M，Maskell S，Gordon N，et al.A tutorial on particle filters for online nonlinear/non-Gaussian Bayesian tracking[J].Signal Processing，IEEE Transactions on， 2002，50（2）：174-188.