一种基于拒绝采样的粒子滤波算法的非线性动态系统信号处理方法技术方案

一种基于拒绝采样的粒子滤波算法的非线性动态系统信号处理方法，采用粒子描述动态系统的状态空间，设非线性动态系统的状态空间模型为：ｘｋ＝ｆ（ｘｋ－１）＋ｖｋ－１ｚｋ＝ｈ（ｘｋ）＋ｎｋ，其中，ｘｋ和ｚｋ分别表示系统在ｋ时刻的状态和观测值，ｆ（ｘｋ－１）和ｈ（ｘｋ）分别表示系统的状态转移方程和观测方程，ｖｋ－１和ｎｋ分别表示系统噪声和观测噪声；对于每一个新采样的粒子，首先计算其被接受的概率，然后判断其是否被接受。本发明专利技术能有效减小粒子权重退化、粒子多样性好，有效处理非线性和非高斯的信号。

【技术实现步骤摘要】

本专利技术涉及基于粒子滤波算法的信号处理方法，要求保护的技术方案属于信号处 理、人工智能和计算机视觉领域。
技术介绍
动态系统的状态估计问题涉及很多领域，尤其是信号处理、人工智能和计算机视 觉领域。传统的卡尔曼滤波只适用于线性高斯系统，而扩展卡尔曼滤波也只能应对系统的 弱非线性。因此，适用于非线性、非高斯系统的粒子滤波备受关注。粒子滤波是一种基于蒙特卡洛模拟和递推贝叶斯估计的滤波方法。它采用粒子 描述状态空间，使用一组带权重的粒子近似表示系统的后验概率密度，并通过模型方程和 观测信息实现递推的估计过程。常见的粒子滤波算法包括SIR粒子滤波、辅助粒子滤波 (APF)、正则粒子滤波(RPF)、高斯粒子滤波(GPF)和无迹粒子滤波(UPF)。粒子滤波算法中的两个关键技术是建议分布的选取和重采样算法。虽然好的建议 分布可以延缓粒子权重的退化问题，但最优的建议分布在一般情况下是无法获取的。同时， 虽然重采样算法解决了粒子权重的退化问题，但常见的重采样算法只是对粒子的简单复制 和剔除，从而又导致了粒子匮乏问题(粒子多样性的丧失)。
技术实现思路
为了克服现有基于粒子滤波算法的动态系统信号处理方法中的粒子权重退化和 粒子匮乏、无法有效处理非线性和非高斯信号的不足，本专利技术提出一种能有效减小粒子权 重退化、粒子多样性好，有效处理非线性和非高斯的信号的基于拒绝采样的粒子滤波算法 的非线性动态系统信号处理方法。为了解决上述技术问题提出的技术方案为—种基于拒绝采样的粒子滤波算法的非线性动态系统信号处理方法，采用粒子描 述动态系统的状态空间，设非线性动态系统的状态空间模型为xk = f (xj +vk_1zk = h (xk) +nk其中，xk和zk分别表示系统在k时刻的状态和观测值，f (x^)和h (xk)分别表示 系统的状态转移方程和观测方程，Vh和nk分别表示系统噪声和观测噪声；对于每一个新采样的粒子，首先计算其被接受的概率，然后判断其是否被接受，递 推过程包括以下步骤第一步，根据k-1时刻的N个粒子xj，i = 1,2,…，N，N为自然数，估计观测状 态zk的发生概率P(zk)p{zk) = \p(zk ^^{X^ I^^Jt/On)3 其中，< 为给定xj时xk的某个特征值，取为均值EUklxJ)或一个采样值i p(xk|xk_11)。第二步，初始化k时刻接受的粒子数R为0 ；第三步，若R < N，则重复下述过程对于每一个粒子&_/，经采样得到xk、计算 xkiJ的接受概率为pi ! 、其中系数？1与粒子的采样次数d成反比，且一的计算方法 如下 第四步，从接受的R个粒子中选取N个代表系统状态的后验概率密度。第五步，输出系统状态的估计值J^1: 本专利技术的技术构思为该算法的一个递推过程包括以下基本步骤1)、根据k_l时刻的N个粒子估计观测状态zk的发生概率。2)、初始化k时刻接受的粒子数R为0。3)、若R < N，则重复下述过程对于每一个粒子&_/，经采样得到xkij，计算xkij的 接受概率为pi r、其中系数pi与粒子的采样次数！^成反比。4)、从接受的R个粒子中选取N个代表系统状态的后验概率密度。5)、输出。本专利技术具有以下优点1、采用拒绝采样算法使接受的粒子具有相同的权重，省去了重采样算法，避免了 传统粒子滤波算法中的粒子权重退化和粒子匮乏问题。2、通过选择不同的采样接受策略可以减少算法的运行时间。附图说明图1为基于拒绝采样的粒子滤波算法流程图。图2为基于拒绝采样的粒子滤波算法与其他粒子滤波算法在不同粒子数下的均 方根误差均值对比曲线示意图。图3为基于拒绝采样的粒子滤波算法与其他粒子滤波算法在不同粒子数下的运 行时间对比曲线示意图。图4为基于拒绝采样的粒子滤波算法在不同mi-pi下的均方根误差均值对比曲线示意图。图5为基于拒绝采样的粒子滤波算法在不同一^1下的运行时间对比曲线示意图。 具体实施例方式下面结合附图和实施例对本专利技术作进一步说明。参照图1 图5，一种基于拒绝采样的粒子滤波算法的非线性动态系统信号处理 方法，采用粒子描述动态系统的状态空间，设非线性动态系统的状态空间模型为xk = f (Xh) +VHzk = h(xk)+nk其中，xk和zk分别表示系统在k时刻的状态和观测值，f (x^)和h (xk)分别表示 系统的状态转移方程和观测方程，Vh和nk分别表示系统噪声和观测噪声；对于每一个新采样的粒子，首先计算其被接受的概率，然后判断其是否被接受，递 推过程包括以下步骤第一步，根据k_l时刻的N个粒子‘/(i = 1,2,…，N)估计观测状态zk的发生 概率P(zk) ' 其中，< 为给定xj时xk的某个特征值，取为均值EUklxJ)或一个采样值i p(xk|xk_11)。第二步，初始化k时刻接受的粒子数R为0 ；第三步，若R < N，则重复下述过程对于每一个粒子，经采样得到xkij，计算 xkiJ的接受概率为pi ! 、其中系数？1与粒子的采样次数d成反比，且一的计算方法 如下“ K) k) 第四步，从接受的R个粒子中选取N个代表系统状态的后验概率密度。第五步，输出系统状态的估计值J^1: 本实施例通过一个非线性动态系统的状态估计对本专利技术和其它几种粒子滤波算 法进行比较。系统的状态空间模型如下 其中，系统噪声方差和观测噪声方差分别取为10和1。p(zk)的计算采用如下方 法 粒子的采样接受使用如下分治策略对于p(zk| y ki) <p(zk)的粒子，m1 = mi = 1，pi = Pl = 5 ；对于其他粒子，mi = m2 = 5，pi = p2 = 1。设定观测时间为100，运行次数 为100，在粒子数N分别取50、100、150、200、250、300时，本专利技术提出的算法与其他粒子滤波 算法所产生的均方根误差(RMSE)均值和所运行的时间分别如图2和图3所示。图2中，本专利技术算法的RMSE均值在粒子数目小于150时明显优于其它算法，在粒 子数目大于150时和UPF算法接近。图3中，本专利技术算法在不同粒子数目下的运行时间略 多于APF算法，约为UPF算法的55%。图4、图5分别为本专利技术算法在上述分治策略中mi、p2取1，Pl、m2取不同值时的 RMSE均值和运行时间对比曲线。在不同的Pl、m2取值下，RMSE均值基本不变，而运行时间 可以部分减少。权利要求，采用粒子描述动态系统的状态空间，其特征在于设非线性动态系统的状态空间模型为xk＝f(xk-1)+vk-1zk＝h(xk)+nk其中，xk和zk分别表示系统在k时刻的状态和观测值，f(xk-1)和h(xk)分别表示系统的状态转移方程和观测方程，vk-1和nk分别表示系统噪声和观测噪声；对于每一个新采样的粒子，首先计算其被接受的概率，然后判断其是否被接受，递推过程包括以下步骤第一步，根据k-1时刻的N个粒子xk-1i，i＝1，2，...，N，估计观测状态zk的发生概率p(zk)p(zk)＝∫p(zk|xk-1)p(xk-1|zk-1)d(xk-1) <mrow><mo>&Proportional;</mo><mundero本文档来自技高网...

【技术保护点】
一种基于拒绝采样的粒子滤波算法的非线性动态系统信号处理方法，采用粒子描述动态系统的状态空间，其特征在于：设非线性动态系统的状态空间模型为：ｘ↓［ｋ］＝ｆ（ｘ↓［ｋ－１］）＋ｖ↓［ｋ－１］ｚ↓［ｋ］＝ｈ（ｘ↓［ｋ］）＋ｎ↓［ｋ］其中，ｘ↓［ｋ］和ｚ↓［ｋ］分别表示系统在ｋ时刻的状态和观测值，ｆ（ｘ↓［ｋ－１］）和ｈ（ｘ↓［ｋ］）分别表示系统的状态转移方程和观测方程，ｖ↓［ｋ－１］和ｎ↓［ｋ］分别表示系统噪声和观测噪声；对于每一个新采样的粒子，首先计算其被接受的概率，然后判断其是否被接受，递推过程包括以下步骤：第一步，根据ｋ－１时刻的Ｎ个粒子ｘ↓［ｋ－１］↑［ｉ］，ｉ＝１，２，．．．，Ｎ，估计观测状态ｚ↓［ｋ］的发生概率ｐ（ｚ↓［ｋ］）：ｐ（ｚ↓［ｋ］）＝∫ｐ（ｚ↓［ｋ］｜ｘ↓［ｋ－１］）ｐ（ｘ↓［ｋ－１］｜ｚ↓［ｋ－１］）ｄ（ｘ↓［ｋ－１］）∝＊ｐ（ｚ↓［ｋ］｜ｘ↓［ｋ－１］↑［ｉ］）＝＊∫ｐ（ｚ↓［ｋ］｜ｘ↓［ｋ］）ｐ（ｘ↓［ｋ］｜ｘ↓［ｋ－１］↑［ｉ］）ｄ（ｘ↓［ｋ］）≈ｐ（ｚ↓［ｋ］｜μ↓［ｋ］↑［ｉ］）其中，μ↓［ｋ］↑［ｉ］为给定ｘ↓［ｋ－１］↑［ｉ］时ｘ↓［ｋ］的某个特征值，取为均值Ｅ（ｘ↓［ｋ］｜ｘ↓［ｋ－１］↑［ｉ］）或一个采样值ｘ↓［ｋ］↑［ｉ］：ｐ（ｘ↓［ｋ］｜ｘ↓［ｋ－１］↑［ｉ］）。第二步，初始化ｋ时刻接受的粒子数Ｒ为０；第三步，若Ｒ＜Ｎ，则重复下述过程：对于每一个粒子ｘ↓［ｋ－１］↑［ｉ］，经采样得到ｘ↓［ｋ］↑［ｉｊ］，计算ｘ↓［ｋ］↑［ｉｊ］的接受概率为ｐ↑［ｉ］.ｒ↑［ｉｊ］，其中系数ｐ↑［ｉ］与粒子ｘ↓［ｋ－１］↑［ｉ］的采样次数ｍ↑［ｉ］成反比，且ｒ↑［ｉｊ］的计算方法如下：ｒ↑［ｉｊ］＝ｐ（ｘ↓［ｋ］↑［ｉｊ］｜ｚ↓［ｋ］）／ｐ（ｘ↓［ｋ］↑［ｉｊ］）＝ｐ（ｚ↓［ｋ］｜ｘ↓［ｋ］↑［ｉｊ］）／ｐ（ｚ↓［ｋ］）第四步，从接受的Ｒ个粒子中选取Ｎ个代表系统状态的后验概率密度。第五步，输出系统状态的估计值ｘ↓［ｋ］↑［μ］：ｘ↓［ｋ］↑［μ］≈１／Ｎｘ↓［ｋ］↑［ｉ］。...

【技术特征摘要】

【专利技术属性】
技术研发人员：潘赟，郑宁，严晓浪，王一木，孙纲德，
申请(专利权)人：浙江大学，
类型：发明
国别省市：86[中国|杭州]