粒子滤波的基本原理笔记

粒子滤波的基本原理

粒子滤波算法广泛应用在视觉跟踪领域、通信与信号处理领域、机器人、图像处理、金融经济、以及目标定位、导航、跟踪领域，其本质是利用当前和过去的观测量来估计未知量的当前值。在粒子滤波算法中使用了大量随机样本，采用蒙特卡洛仿真来完成递推贝叶斯滤波过程，其核心是使用一组具有相应权值的随机样本（粒子）来表示状态的后验分布。该方法的基本思路是选取一个重要性概率密度函数并从中进行随机抽样，得到一些带有相应权值的随机样本后，在状态观测的基础上调节权值的大小和粒子的位置，再使用这些样本来逼近状态后验分布，最后通过这组样本的加权求和作为状态的估计值。粒子滤波不受系统模型的线性和高斯假设约束，采用样本形式而不是函数形式对状态概率密度进行描述，使其不需要对状态变量的概率分布作过多的约束，适用于任意非线性非高斯动态系统，是目前最适合于非线 性、非高斯系统状态的滤波方法【Arulampalam M S, Maskell S, Gordon N, et al. A tutorial on particle filters for online nonlinear/non-Gaussian Bayesian tracking[J]. IEEE Transactions on Signal Processing, 2002, 50(2):174-188.】

1 动态系统的状态空间模型

状态空间模型包括系统状态方程和观测方程，其通用的表示方法分别为【梁军. 粒子滤波算法及其应用研究[D]. 哈尔滨工业大学, 2009.】【黄小平, 王岩, 廖鹏程. 粒子滤波原理及应用——MATLAB 仿真[M].电子工业出版社. 2017】

()1,k k k f -=X X W （1） (),k k k h =Z X V （2）

其中()f ∙和()h ∙为已知函数，

k W 和k V 是概率密度已知的随机变量，k X 代表k 时刻的状态量，k Z 代表k 时刻的观测量，k W 和k V 是相互独立的。

关于系统的状态方程和观测方程，通常也可用()1k k p -X X 表示状态转移模型；()k k p Z X 表示观测似然模型；()0p X 表示初始状态的先验分布；()0:1:k k p X Z 表示系统的后验密度；

()1:k k p X Z 表示边沿后验密度，或称为后验滤波密度。卡尔曼滤波以及粒子滤波算法的本质即是利用观测序列1:k Z 对当前状态进行优化，从而得到k 时刻的后验滤波密度，进而得到k 时刻的状态值。【Merwe R V D, Doucet A, Freitas N D, et al. The unscented particle filter[C]// International Conference on Neural Information Processing Systems. MIT Press, 2000:563-569.】

2 贝叶斯滤波原理

贝叶斯估计是利用先验知识和实际观测值构造系统的后验概率来估计系统的状态量。若系统的先验概率密度为()0:k p X ，观测值为1:k Z ，则系统的后验概率密度可表示为

()()

()

()()1:0:0:1:0:1:0:0:0:k k k k k k k k k

p p p p p d =⎰Z X X Z X Z X X X （3）

从式（3）可以看出，后验概率实际上是利用观测值对先验概率进行修正。由于进行了修正，后验概率密度较先验概率密度更加接近被估计值的真实概率密度。滤波问题是为了计算后验滤波概率密度()1:k k p X Z ，它实际上是()0:1:k k p X Z 的边沿概率密度即

()()1:0:1:01......k k k k k p p d d -=⎰⎰X Z X Z X X （4）

当每一个新的观测量k Z 到来时，即需要完成式（4）得到后验滤波概率密度。

为清晰的给出后验滤波概率密度的逐级递推形式。首先，通过预测步骤得到不包含k 时刻观测值的k 时刻系统的先验滤波概率密度

()()()1:1111:11k k k k k k k p p p d -----=⎰X Z X X X Z X （5）

其中()1k k p -X X 由系统状态方程决定，()11:1k k p --X Z 是1k -时刻的后验滤波概率密度。然后，根据贝叶斯公式，利用k 时刻的先验滤波概率密度()1:1k k p -X Z 推导k 时刻的后验概率密度的表示形式

()()()

()

1:11:1:1k k k k k k k k p p p p --=

Z X X Z X Z Z Z （6）

其中()k k p Z X 由观测方程决定。

()()()1:11:1k k k k k k k p p p d --=⎰Z Z Z X X Z X （7）

一般为归一化常数。因此贝叶斯滤波的递推框图如下所示

()0p X

∙∙∙∙∙∙∙∙∙∙∙

图1 贝叶斯递推滤波框图

贝叶斯递推过程中的主要的几个公式

()()()1:1111:11k k k k k k k p p p d -----=⎰X Z X X X Z X （利用上一时刻的后验滤波密度，结合状态转移概率对这一时刻的先验滤波概率密度进行求解）

()()()1:11:1k k k k k k k p p p d --=⎰Z Z Z X X Z X （利用先验滤波概率密度以及观测似然密度得到观测概率）

()()()

()

1:11:1:1k k k k k k k k p p p p --=

Z X X Z X Z Z Z （利用这一时刻的先验滤波概率密度，结合观测似然概率

密度对这一时刻的后验滤波概率密度进行求解）

以上三个式子是贝叶斯滤波递推求解后验滤波概率密度的核心公式，称为贝叶斯递推滤波算法。这里给出贝叶斯滤波原理的一个简单例子。

一个简单的例子：假设动态系统状态空间模型中状态方程和观测方程分别为

1k k k -=+X X W （8） k k k =+Z X V （9）

其中()~0,1k N W ，()~0,0.5k N V ，并给出初始状态00=X 。由此，可以得到：

0X 是以概率1取值为0，即

()()00p δ=X X （10）

状态转移概率

()()2112k k k k p --⎛⎫-=- ⎪ ⎪⎝⎭

X X X X （11）

观测似然函数密度

()()220.5k k k k p ⎛⎫

-=- ⎪ ⎪⨯⎝⎭

Z X Z X （12）

该系统等效于某信号经过一个()0,1N 的高斯信道进行传输，然后再利用一个测量误差为

()0,0.5N 的测量设备对该信号进行测量。当从0时刻向1时刻转移时

()()()(

)()101000

1000212p p p d p d δ==⎛⎫=- ⎪

⎝⎭

⎰⎰X Z X X X X X X X X X （13）