一种改进的粒子滤波方法

一种改进的粒子滤波方法

0基于线性化的一般方法非线性滤波问题是通过在线获得受干扰噪声的观测量，并对非线性系统的不观测值进行加权。它在自动控制、信号处理、目标跟踪、人工智能以及导航制导等诸多领域具有广泛的应用。最著名的非线性滤波方法是扩展卡尔曼滤波(extenedKalmanfilter,EKF),其基本思想是使用泰勒展开对非线性系统进行线性化,但这种线性化误差较大,而且许多实际问题中很难得到非线性函数的Jacobian矩阵。近年来出现了一些无需计算Jacobian矩阵的非线性滤波方法,如无迹卡尔曼滤波器(unscentedKalmanfilter,UKF)、均差滤波器(divideddifferencefilter,DDF)、中心差分滤波器(centraldifferencefilter,CDF)等,这些方法的共同问题是在非线性、非高斯性较强时收敛性急剧下降甚至发散。粒子滤波器(particlefilter,PF)使用序贯蒙特卡罗方法,通过寻找一组在状态空间中传播的随机采样对后验概率密度进行近似,它可以完整反映状态的后验分布并容易得到如均值、模和方差等统计特征,适用于任何分布的非线性系统。重要性采样(importancesampling,IS)是粒子滤波的重要部分,因此要求建议分布(proposaldistribution,PD)能够准确地逼近后验概率分布。通常设计这样的建议分布函数比较困难,比较常用的方法是对系统状态的转移概率进行采样,但如果最新的观测信息位于先验概率分布的尾部或者似然函数相比先验概率是峰化的,就会导致粒子选择的盲目性,降低估计的精度。为克服这一问题,出现了一些基于线性化的建议分布改进方法。例如,文献提出了一种利用EKF高斯逼近产生建议分布函数的卡尔曼粒子滤波(extendedKalmanparticlefilter,EKPF),使滤波性能有所改善。但由于EKF使用泰勒一阶展开并需要计算Jacobian矩阵,因此算法的精度和效率依然受到限制。文献提出了一种基于UKF的无迹粒子滤波(unscentedparticlefilter,UPF),UKF的估计精度高于EKF,可以获得更加接近真实后验分布的建议分布,因此无迹粒子滤波的性能有了很大改善。但UKF在选取Sigma点集时,采样点数与状态维数相关,计算量随着维数增大上升较快。集合卡尔曼滤波(ensembleKalmanfilter,EnKF)是一种蒙特卡罗方法的卡尔曼滤波,它通过集合采样避免了Jacobian矩阵的计算,同时具有处理高维非线性系统的能力。文献利用EnKF对粒子滤波的采样权值进行改进,用于处理非高斯数据同化。文献使用两个采样集合,通过EnKF和PF并行处理,利用PF计算EnKF采样点的权值,并将这种方法应用于洪水预测模型。本文利用集合卡尔曼滤波产生粒子滤波每一时刻的建议分布函数,这样的建议分布函数能够更接近真实的后验概率密度,且具有较小的计算量。仿真实验验证了这种方法的优越性。1标准粒子滤波算法的基本原理对于如下非线性离散系统xk=f(xk-1,vk-1)(1)zk=h(xk,wk)(2)式中,xk∈Rnx为系统在k时刻的状态;zk∈Rnz为系统在k时刻的观测向量;vk∈Rn为系统过程噪声;wk∈Rm为k时刻的观测噪声;映射f:Rnx×Rnv→Rnx和h:Rnx×Rnw→Rnz都是有界非线性函数,分别代表系统的状态和观测模型。滤波的目的就是要获得系统后验分布p(xk|z1:k),继而得到系统状态的统计特性,如均值、最大后验概率和置信区间等。粒子滤波又称自举滤波(bootstrapfiltering,BF)或重要性采样重采样(samplingimportanceresampling,SIR)滤波,其本质是采用序贯蒙特卡罗(sequentialMonteCarlo,SMC)方法。基本思想是:从后验概率密度p(xk|z1:k)中独立的抽取N个采样点,通过加权求和,近似表示后验概率分布。在很多情况下,后验分布可能是多变量、高维、多峰、非解析的,因此很难直接从中采样,为此引入了重要性采样。其基本思路是使用易于采样的建议分布函数q(xk|z1:k)代替后验密度函数进行采样。重要性采样的一个重要缺陷是“蜕化”问题,这时就需要用到“重采样”技术。其基本思想是抑制或剔除小权值粒子,对大权值粒子依其权值进行复制。标准粒子滤波算法可总结如下(N为粒子个数):(1)初始化k=0,Fori=1:N,从先验概率p(x0)中随机采样x(i)0。(2)Fork=1,2…(a)重要性采样Fori=1:N,抽取采样点ˆx(i)k~q(xk|x(i)1:k-1,z1:k),计算每个粒子的权值ω(i)k=ω(i)k-1p(zk|ˆx(i)k)p(ˆx(i)k|x(i)k-1)q(ˆx(i)k|x(i)0:k-1,z1:k)归一化权值˜ω(i)k=ω(i)k/Ν∑j=1ω(j)k(b)重采样得到x(i)k并将权值置为1Ν。(c)输出p(x0:k|z1:k)≈1ΝΝ∑i=1δ(x(i)0:k)(dx0:k)2k的分类提高粒子滤波算法性能最有效的方法是选择一个适当的建议分布函数,也就是使建议分布函数和似然函数具有更大的重叠区域。针对建议分布选择问题,最优选择标准是最小化重要性权值的方差,并有如下重要命题:p(x0:k|x(i)0:k-1,z1:k)是最小化的基于x(i)0:k-1和z1:k的重要性权方差Varq(xk|x(i)0:k-1,z1:k)(ω(i)k)的最优建议分布。但上述分布至少存在两个缺陷,首先是必须从可能是非标准的分布中采样粒子比较困难;其次是涉及积分计算,而该积分通常是非解析的。目前,通常选取p(x0:k|x(i)0:k-1)作为建议分布函数,但由于该分布未融入最新观测信息,会导致较高的权值方差,影响滤波精度。文献分别提出了利用EKF和UKF产生建议分布,使算法精度得到提高。本文使用EnKF来产生建议分布函数,由于EnKF能够得到状态的最大后验概率估计,将最新观测信息融入建议分布,使得产生的采样样本更接近真实的样本,又避免了对非线性系统进行线性化,同时有效控制了计算量,提高了算法性能。2.1改进的卡尔曼增益估计算法EnKF的本质是一种基于蒙特卡罗方法的卡尔曼滤波。EnKF因其较高的计算效率和处理高维非线性系统的能力,被广泛应用于大气预报领域。其基本思想是,初始化一组系统的状态采样作为背景集合,利用观测信息通过卡尔曼滤波对背景数据集中每个个体进行更新,得到分析集合。分析集合用来估计状态的真实均值和方差。通过系统模型传递采样集合,可以得到下一时刻的背景数据集。这样使用集合估计真实统计值,提高了估计精度,同时使计算量明显降低。定义集合Xbk={xbk,i,i=1,2,…,n}为k时刻状态的背景集合,它由k-1时刻的分析集合Xak-1传递而来,n为集合样本数。采样均值和方差可由下式计算ˆxbk=1nn∑i=1xbk,i(3)ˆΡbk=1n-1n∑i=1(xbk,i-ˆxbk)(xbk,i-ˆxbk)Τ(4)实际应用中,不必计算ˆΡbk,使用如下公式代替ˆΡkxh=1n-1n∑i=1(xbk,i-ˆxbk)(h(xbk,i)-h(ˆxbk))Τ(5)ˆΡkhh=1n-1n∑i=1(h(xbk,i)-h(ˆxbk))(h(xbk,i)-h(ˆxbk))Τ(6)卡尔曼增益的计算为Κk=ˆΡkxh(ˆΡkhh+Rk)-1(7)式中,Rk表示k时刻的观测误差协方差矩阵。利用最新观测信息,可以对背景集合进行更新,得到k时刻的分析集合Xak={xak,i,i=1,2,…,n}如下xak,i=xbk,i+Κk(yk,i-h(xbk,i)),i=1,2,⋯,n(8)式中,yk,i是以观测值yk为均值,Rk为方差的高斯分布的采样。相应的,分析集合的均值和方差为ˆxak=1nn∑i=1xak,i(9)ˆΡak=1n-1n∑i=1(xak,i-ˆxak)(xak,i-ˆxak)Τ(10)2.2在关注高估计精度的同时进行实验结果由前述算法流程可以看出,不同于EKF,EnKF无须对非线性系统进行线性化,避免了Jacobian矩阵的计算,使用采样法近似非线性分布,提高了计算精度,并且可以处理非可导的非线性系统。文献指出,对于受高斯噪声扰动的非线性系统,EnKF对均值和方差的估计至少可以达到二阶水平。所以相比传统的EKF,EnKF产生的分布与后验分布的重叠区域更大,可以对状态分布进行更好的高斯估计。另外,虽然UKF对于非线性问题具有很高的估计精度,但其Sigma点的采样数目与系统的维数相关,当维数增大时,算法的计算量上升比较快,而EnKF的采样点数目是启发式的,可以灵活设定,使算法在保证较高估计精度的同时有效控制计算量。这使得EnKF成为计算粒子滤波建议分布函数的更好选择。使用EnKF产生粒子滤波的建议分布函数,主要是通过传播后验分布的高斯估计以及将每个时刻的最近观测与之结合。换句话说,就是使用EnKF对如下后验概率密度进行递归估计p(xk|z1:k)≈pΝ(xk|z1:k)=Ν(ˆxk,ˆΡk)(11)在粒子滤波的框架中,对每个粒子使用一个独立的EnKF产生和传递高斯建议分布q(x(i)k|x(i)1:k-1,z1:k)=Ν(ˆx(i)k,ˆΡ(i)k)(12)在k-1时刻,使用EnKF和最近观测信息对每个粒子的建议分布的均值和方差进行更新,在k时刻从这个分布中采样得到新的粒子。2.3集合卡尔曼粒子滤波器本文提出的这种使用EnKF产生粒子滤波建议分布的新滤波器叫做集合卡尔曼粒子滤波(ensembleKalmanparticlefilter,EnKPF)。算法步骤如下:3pf、enkf和opf关于状态估计的比较为比较EnKPF和现有滤波算法如PF、PF-EKF和UPF等对非线性系统的估计性能,这里对算法进行仿真实验。仿真硬件环境是Intel(R)Core2DuoCPUT7250@2.00GHz,2GRAM,WindowsXP操作系统,软件利用MatlabR2007b编写。考虑如下合成尺度估计模型xk=1+sin(απk)+0.5xk-1+vk-1yk={0.2x2k+wk,k≤300.5xk-2+wk,k＞30(13)式中,α=6e-2;过程噪声vk服从Gamma分布Gamma(4,3);观测噪声wk服从正态分布N(0,1e-4)。给定观测信息yk,使用不同滤波器对系统状态进行估计。粒子数为100,使用残差采样法,将采样均值作为每一时刻算法输出值。令EnKF的集合采样个数为5,运行60个时刻,每次随机起始并进行100次独立重复实验。图1为一次独立实验条件下得到的4种滤波算法的状态估计情况,可以看出PF-EKF、UPF和EnKPF都是用最大后验估计产生粒子滤波的建议分布函数,通过融入最新观测信息将先验分布向似然分布移动,估计性能优于标准PF算法。但由于EKF使用一阶泰勒展开,忽略高阶项,这种线性化过程产生了较大误差,所以PF-EKF的估计性能低于UPF和EnKPF。同时,UKF和EnKF都是通过近似非线性概率密度分布代替非线性函数的方法,在这一单次独立实验当中EnKPF的估计性能与UPF较为相似。图2为4种粒子滤波算法的均方根估计误差,标准PF和PF-EKF会出现不稳定的情况,且估计误差较大,而UPF和EnKPF由于建议分布函数更加接近真实后验分布,其状态估计的均方根误差明显降低,并且EnKPF的估计效果最好。图3为各粒子滤波算法状态值和观测值概率分布情况,可以观察到使用标准PF、PF-EKF和EnKPF对观测值概率分布和状态后验概率分布的估计情况,由于EnKPF在产生建议函数时既融入了最新观测信息,又避免了对非线性系统的线性化处理,直接通过真实非线性系统进行传递,因此它的估计效果最好。表1列出了各滤波器状态估计的均方根误差的均值和方差以及平均运行时间,“—”表示未进行考察对比。由表1可知,对于非线性系统状态估计问题,EnKF的估计性能优于EKF,因而利用EnKF产生的建议分布函数更加接近真实后验概率分布,所以EnKPF的估计性能优于标准PF和PF-EKF算法。相比UPF,在当前的集合采样点数下,EnKPF取得了更好的估计效果。在运行时间上,由于UPF和EnKPF都使用了采样点集的策略,因此高于传统PF及PF-EKF,而UPF的运行时间较之EnKPF更长,这主要是由于在利用UKF产生建议分布函数时,因为噪声项的存在,需要对状态进行扩维处理,导致Sigma粒子个数比较多。另外,EnKPF不同于UPF采样点数固定,其采样点数灵活可变,可以通过增加采样点数进一步提高估计精度。为比较不同粒子数情况下各粒子滤波算法的性能,分别选取粒子数为10、30、50进行实验。图4给出了PF、PF-EKF、UPF和EnKPF在不同粒子数时状态估计的均方根误差均值,可以看出,在粒子数较少时EnKPF仍然具有较高的估计性能。为考察不同集合采样点数情况下EnKPF的性能,分别选取采样点数为10、20、30、40、50和100进行仿真。表2给出了不同采样点数情况下状态估计的均方根误差的均值和平均运行时间。由表2可知,随着采样点数增加,EnKPF对于系统状态估计的均方根误差呈递减趋势,但是采样点数为50时的估计性能与采样点数为100时的估计性能相差不大,而算法运行时间却随着采样点数变化而显著增加。由此可知,当采样点数达到一定数目后,增加采样点数对进一步提高算法估计精度效果甚微,反而会降低算法的计算效率。4标准粒子