闪烁噪声下的自适应无迹卡尔曼算法

陈 政，刘 康

（昆明理工大学信息工程与自动化学院，云南 昆明 650500）

0 引言

卡尔曼最初提出的滤波基本理论只适用于线性系统，并且要求量测也必须是线性的[1]。卡尔曼算法可以看成是一种贝叶斯方法，贝叶斯方法通过对模型引入先验分布，并利用似然函数求解后验分布的方法来对模型中的不确定性进行建模[2]。在线性估计算法的发展历程中，卡尔曼滤波扮演了重要的角色，但该算法实现的前提是系统模型已知，与实际工业中的状态估计问题不相符[3]；同时，真实环境中卫星雷达受到的非高斯闪烁噪声进一步加重了后验概率密度的非高斯特性[4]，于是，传统的卡尔曼滤波已经不再适用[5][6]。由此带来了改变后的使用到非线性系统观测中的滤波理论，如无迹卡尔曼滤波器（UKF），可用于导航、目标状态估计、飞行器轨道确定等领域并获得了良好效果，但噪声统计特性未知时，UKF滤波精度下降[7]。

全球卫星导航系统在许多领域都发挥出了重要的作用，并且仍然在不断地改进与创新[8]。近年来，随着我国空间技术的发展，在轨运行的卫星越来越多[9]。卫星侦察范围广，可对目标进行长时间、大范围的连续侦察和监视，获取情报的时效性强，是现代侦察中不可或缺的重要手段[10]。人造卫星在空中飞行，当自身携带的能源耗尽时便被回收。基于此，设计出来一个复杂计算度低、精确度高的非线性动态系统的滤波算法尤为重要，本文结合EM算法和无迹卡尔曼滤波算法，提出一种新的方法，仿真结果表明，本文提出的方法滤波精度较高，能够有效的跟踪目标。

1 问题描述

1.1 闪烁噪声

随着计算机技术的不断发展，利用计算机来实现目标追踪功能成为目前计算机领域中最热门的课题之一[11]。目标跟踪是计算机视觉领域的一个重要研究方向，被广泛应用于智能监控、目标识别、交通监视等方面[12]。

闪烁噪声分布跟高斯分布的主要区别在于尾部较长，其中心区域则类似高斯分布[2]。本文采用混合高斯模型来模拟闪烁噪声，其概率密度如下：

其中，N(x; μi,Pi)表示均值为μi，方差为 Pi的高斯分布，闪烁强度 ε ∈ [ 0,1]。

1.2 运动模型

目标跟踪的只要目的就是估计移动目标的状态轨迹[13]。考虑一般的目标跟踪问题，在笛卡尔坐标系中，目标做匀速运动的离散时间状态系统方程和观测方程分别为：

其中，kx表示目标在K时刻的状态向量，yk表示K时刻状态的观测值，F为状态转移矩阵，L为过程噪声驱动矩阵，在一般的雷达目标跟踪中，()·h为有界的非线性函数，过程噪声1-kw 是均值为零方差为Q的高斯白噪声，实际情况中，Q经常是未知的，观测噪声kv是均值为零方差为R的高斯白早噪声，实际情况中，kv经常是闪烁噪声。

2 算法描述

2.1 无迹卡尔曼滤波算法

考虑如下含有加性高斯噪声的离散非线性系统[14-15]：

其中， k ≥ 0 为离散时间变量， x ∈ Rn为状态向量，z ∈ Rm为输出向量。非线性函数 fn∈Rm，h ∈ Rn。过程噪声 wk-1和量测噪声 vk分别为n维和m维的高斯白噪声，并具有以下统计特性：

在很多情况下，系统噪声都是非零均值。因此，给出噪声均值非零情况下的 UKF 算法。由于系统过程噪声和量测噪声为非零均值的高斯噪声且其均值分别为q和r。则系统噪声可写为如下形式：

显然，kµ和ηk为互不相关的零均值高斯白噪声，且其方差分别为Q和R。

此时, 离散非线性系统 (4) 可以被等价表述为以下形式：

则UKF 算法的计算步骤如下：

步骤1. 初始化系统变量

步骤2. 时间更新

根据确定的采样策略，得到服从均值 xk-1、方差 pk-1的Sigma 点 { χi,k-1},i = 0,...,2n，为：

Sigma点的对应权值为：

其中，α描述了 Sigma 点的散布程度,通常取一小的正值（如 0.01），k通常取0；β用来描述 x的分布信息（Gauss情况下β的最优值为 2）；()i表示矩阵 ( n + λ )Pk-1的平方根矩阵的第k-1列； w(m)(i = 0,1,2,...,2n)为求一阶统计特性时

i

的权系数 w(c)(i = 0,1,2,...,2n)为求二阶统计特性时

i的权系数.计算 Sigma 点 {χi,k−1}, i = 0,··· ,2n 沿非线性函数 f ( .)的传播结果， 并考虑非零均值噪声的作用。

从而得到随机向量沿非线性函数 (.)f 传播的后验均值和协方差：

步骤3. 量测更新

其对应权值与式 (9) 的权值相同。

计算 Sigma 点｛χi,k|k-1｝, i = 0,··· ,2n 沿非线性函数 h( .)的传播结果，并考虑非零均值噪声的作用。

得到随机向量沿非线性函数 (.)h 传播的后验均值、自协方差和互协方差：

步骤4. 滤波更新

根据量测数据kz得到最小方差估计结果为：

步骤4 完成后，则返回步骤2，开始下一时刻的滤波过程。

2.2 EM算法在目标跟踪中的应用

首先给出一般形式的线性状态空间模型的参数估计问题：

其中F为状态转移矩阵，H为观测矩阵， wt，vt为过程噪声和观测噪声，它们分别服从以下分布：

对运动目标模型的参数估计即在知道tx和ty的情况下（将tx和ty合起来看做完全数据）对目标的运动参数F、H、Q、R进行估计[16]。设参数的F、H、Q、R的对数似然函数为：

多维高斯分布的概率密度函数：

且 xt+1～ N(F xt,Q)， yt～ N(H xt,R),即完全数据的似然也服从高斯分布，将这两个式子带入式（26）中，得：

将上式展开，我们就可以得到对数似然函数的最终表达式：

由EM算法的基本原理可知，需要对（28）式进行最大化处理，即：

化简后得：

EM 算法是一种迭代算法，因此我们可以先假定F,H,Q,R的值，得到tx和ty的值后进行平滑滤波（E步），得到新的tx和ty的值后再利用式（30-33）估计出参数F,H,Q,R新的值（M步），然后再回到E步，这样不断的重复此过程直到参数收敛[17]。

3 仿真实验

3.1 仿真场景

为了验证本文提出的基于EM算法的改进无迹卡尔曼滤波算法（下面简称EMUKF）的滤波效果，在过程噪声统计特性未知和带有闪烁噪声的情况下，分别采用UKF算法和EMUKF算法进行仿真。二维情况下，设雷达观测站位于坐标原点，目标做近似匀速直线运动，运动方程描述如下：

向量 x = [ x x˙ y y˙]′，其中x和x˙分别表示X轴上的位置分量和速度分量，y和y˙分别表示Y轴上的位置分量和速度分量，目标的初始状态为 x0=[50 1 50 - 1 ]′， zk是观测向量， rk和θk分别是目标到观测站的距离和角度。仿真中过程噪声是均值为零的高斯白噪声，观测噪声是闪烁噪声，过程噪声的统计特性都是未知。仿真中过程噪声统计特性如下：

仿真中先验过程噪声 Qk=diag ｛ 0. 01 , 0.01｝。仿真时间为100 s，采样周期T=1s。

3.2 仿真结果

针对过程噪声统计特性未知，过程噪声是一般高斯白噪声和观测噪声是闪烁噪声这两种情况分别做 100次蒙特卡洛仿真。仿真实验的均方根误差（Root mean square error, RMSE）定义为：

其中N为蒙特卡洛仿真次数。

图1为仿真中真实轨迹和两种算法滤波轨迹的对比，图2、3为仿真下两种算法的均方根误差曲线比较。表1是仿真下两种滤波算法均方根误差的均值和方差的统计数据。

图1 三条不同轨迹Fig.1 Three different tracks

图2 位置均方根误差Fig.2 Root mean square error about position

通过上述仿真结果可以看出闪烁噪声和过程噪声统计特性未知会造成滤波效果差甚至发散的情况。本文提出的算法在闪烁噪声下比标准的无迹卡尔曼滤波具有更良好的滤波效果，其均方根误差也明显小于标准的无迹卡尔曼滤波算法。

图3 速度均方根误差Fig.3 Root mean square error about velocity

表1 算法性能比较Tab.1 Algorithm performance comparison

4 结论

针对雷达闪烁噪声情况下，过程噪声统计特性未知的非线性目标跟踪问题，采用无迹卡尔曼滤波的同时运用EM算法对过程噪声进行参数估计，获得较准确的参数后明显的提高了滤波的精度和稳定性。特别是在闪烁噪声的情况下，本文提出的自适应无迹卡尔曼滤波算法要明显优于标准的无迹卡尔曼滤波，表现出很强的抗噪声性能。

[1] 石章松, 刘忠等编著, 目标跟踪与数据融合理论及方法[M]. 北京: 国防工业出版社, 2010. 7.

[2] 侯禹腾. 数据增广求解贝叶斯Logistic 回归模型的方法研究[J]. 软件, 2014, 35(7): 109-115.

[3] 王写, 刘妹琴, 张森林.卡尔曼滤波器、H∞滤波器和鲁棒混合Kalman/H∞滤波器的比较[J].新型工业化, 2012, 2(7):36-45.

[4] I. Bilik; J. Tabrikian. Target tracking in glint noise environment using nonlinear non-Gaussian Kalman filter. 2006 IEEE Conference on Radar. 2006: 6 pp.

[5] Sorenson H W. Kalman filtering: Theory and Application.New York: IEEE, 1985.

[6] F. Daum. Nonlinear filters: beyond the Kalman filter. IEEE Aerospace and Electronic Systems Magazine.2005, 20(8):57 –69.

[7] 秦永元, 张洪钺, 汪淑华等编著, 卡尔曼滤波与组合导航原理(第3版)[M]. 西安:西北工业大学出版社, 2015. 6.

[8] 陈杨毅, 陈凌宇.基于最优几何的改进型多星座卫星选星算法[J]. 软件, 2014,35(1): 26-31.

[9] 杨冬, 孙剑伟. 基于令牌桶算法的卫星数据地面传输流量控制方法研究[J]. 软件, 2016, 37(3): 99-103.

[10] 张志恒, 尹路明, 王茂磊. 采用组合优化策略的电子侦察卫星规划方法[J]. 软件, 2014,35(4): 143-149

[11] 邱亚钦, 刘渭滨, 户磊, 等.动态背景下的运动目标检测方法[J].新型工业化, 2011,1(9): 22-29.

[12] 王江, 付强, 全权, 等.基于Kalman 滤波和直方图匹配的双目视觉跟踪[J].新型工业化, 2013,3(2): 23-33.

[13] [美] Anton J。Haug著, 王欣, 于晓译 Bayesian Estimation and Tracking A Practical Guide(贝叶斯估计与跟踪实用指南)[M]. 北京: 国防工业出版社, 2014.

[14] 武飞, 柳炳利. Kalman 滤波技术在地球化学数据处理中的应用[J]. 软件, 2013,34(9): 70-74.

[15] 韩崇昭, 朱洪艳, 段战胜等著, 多源信息融合(第二版)[M].北京: 清华大学出版社, 2010.9.

[16] Ming Lei; Chongzhao Han; Panzhi Liu. Expectation Maximization (EM) algorithm-based nonlinear target tracking with adaptive state transition matrix and noise covariance. 2007 10th International Conference on Information Fusion. 2007:1–7.

[17] A.P. Dempster; N.M.Laird; D.B.Rubin. Maximum Likelihood from Incomplete Data via the EM Algorithm. Journal of the Royal Statistical Society. SeriesB (Methodological), 1977,39(1): 1-38.