简单粒子滤波与扩展卡尔曼滤波性能对比MATLAB仿真系统
项目介绍
本项目是一个基于MATLAB的仿真环境,专门用于定量对比和视觉化展示两种非线性滤波算法:扩展卡尔曼滤波(EKF)与粒子滤波(PF)。在处理高度非线性系统时,传统的线性滤波往往失效。本项目通过模拟一个包含强非线性状态转移方程和测量方程的时变系统,直观地展现了粒子滤波如何通过蒙特卡洛采样克服线性化误差,以及在相同参数条件下,相比于依赖雅可比矩阵的EKF,PF在估计精度上的优势。
功能特性
- 强非线性模型建模:包含复杂的三角函数、分式结构以及二次项观测方程,能够有效测试滤波器在极端非线性条件下的稳定性。
- 双算法同步仿真:在同一时间基准下运行EKF与支持重采样的SIR粒子滤波。
- 动态粒子分布监控:提供粒子群随时间演化的分布图,反映滤波器在状态空间中的覆盖程度。
- 定量误差评价机制:自动计算并输出两种算法的均方根误差(RMSE),通过数据对比性能。
- 多维度图形化输出:包含状态轨迹对比图、瞬时绝对误差波动图以及粒子密集度演变图。
使用方法
- 启动MATLAB软件。
- 将包含本项目逻辑的所有函数文件置于当前工作路径。
- 运行主函数。程序将自动执行60个步长的仿真迭代。
- 观察命令行窗口输出的RMSE数值。
- 分析弹出的图形窗口,对比蓝色的EKF轨迹、红色的PF轨迹与黑色真实状态的吻合程度。
系统要求
- 软件环境:MATLAB R2016b 及以上版本。
- 硬件环境:具备基础绘图能力的通用计算机。无需特殊图形加速支持。
功能实现逻辑说明
主程序严格遵循以下逻辑步骤进行非线性估计仿真:
- 参数初始化
设置仿真总时长为60个时间步,粒子总数预设为100个。定义过程噪声方差Q=1,测量噪声方差R=1。初始化系统起始状态为0.1。
- 状态与观测数据生成
程序递归生成真实状态轨迹。状态方程结合了线性项、分式非线性项以及基于时间的余弦扰动,并叠加高斯过程噪声。观测方程采用状态平方的分数形式,并叠加高斯测量噪声。
- 扩展卡尔曼滤波(EKF)实现
在每一个时间步,程序首先根据非线性函数进行状态预测。随后,计算状态转移方程关于前一状态的偏导数(Jacobian矩阵 F),以及观测方程关于预测状态的偏导数(Jacobian矩阵 H)。通过计算卡尔曼增益 K,利用观测残差对预测值进行修正,并同步更新误差协方差矩阵。
- 粒子滤波(PF)实现
采用顺序重要性重采样(SIR)架构:
- 采样阶段:从状态转移分布中抽取100个新粒子,模拟噪声对状态转移的影响。
- 权重计算:根据当前观测值,计算每个粒子对应的观测似然概率(基于高斯核函数),并对所有权重进行归一化处理。
- 重采样阶段:调用系统重采样算法,根据权重比例对粒子进行复制或剔除,解决粒子耗尽问题,重采样后重置所有粒子权重。
- 状态估计:取所有粒子位置的算术平均值作为该时刻的最优估计结果。
- 性能度量与可视化
计算整个时间序列内估计值与真实值的均方根误差。绘图阶段将结果分为三个区域:上方区域展示三条轨迹的对比;左下方区域对比两种算法的绝对误差波形;右下方区域以散点图形式展示粒子群分布随时间的迁移过程。
关键函数与算法细节
- 状态转移雅可比矩阵计算
在EKF逻辑中,程序通过手动导出的公式计算 F = 0.5 + 25 * (1 - x^2) / (1 + x^2)^2。这是EKF处理非线性的核心,通过局部线性化来近似非稳态变化。
- 系统重采样辅助算法
为了保证粒子的多样性并避免权重退化,该算法在[0, 1/N]区间产生一个随机起点,并以1/N为步长在累积分布函数(CDF)上进行等间距采样。这种方式相比随机采样具有更低的方差,确保了高权重粒子被准确保留。
- 观测似然估计
PF在更新权重时,利用了高斯概率密度函数。通过计算观测值与粒子预测观测值差值的平方,衡量粒子与真实测量值的匹配程度,实现从粒子群中筛选“优良品种”的功能。