CN110677140A - 一种含未知输入和非高斯量测噪声的随机系统滤波器 - Google Patents

Abstract

本发明涉及一种含未知输入和非高斯量测噪声的随机系统滤波器。首先，针对含未知输入/干扰的离散时间线性随机系统，基于状态方程完成状态的预测估计；其次，针对实际系统的量测噪声通常含有较多野值，噪声的概率分布往往具有尖峰肥尾等强非高斯特性，在求解干扰估计与状态滤波时基于Huber函数代替传统的均方误差准则构造指标函数；再次，通过最小化指标函数，基于不动点迭代法解算出干扰估计与状态的滤波估计；最后，将状态预测、干扰估计、状态滤波按时刻递推，给出滤波器的设计流程。本发明可推广应用于组合导航、目标跟踪、信号处理等领域，解决实际系统的高精度状态估计问题。

Description

一种含未知输入和非高斯量测噪声的随机系统滤波器

技术领域

本发明涉及一种含未知输入/干扰和非高斯量测噪声的随机系统滤波器，针对系统所受到的未知输入或未知动态干扰，采用状态预测、干扰估计、状态滤波的递推滤波器结构，同时考虑到量测噪声的非高斯特性，基于Huber函数构造干扰与状态估计的性能指标函数，增强了滤波器对于量测野值的鲁棒性，本发明可应用于飞行器、船舶、车辆等运动体的组合导航系统，也可推广应用于目标跟踪、信号处理等相关领域。

背景技术

飞行器、船舶、车辆等运动体依靠导航系统实时获取自身的运动信息和姿态信息，而滤波算法又是导航系统信息处理的核心所在。传统的以卡尔曼滤波为典型代表的滤波算法在实际工程应用中已经取得了很大的成功，发挥了积极的作用。然而，运动体的任务需求日益提高，所面临的环境也日益复杂多样，这些都对导航系统的实时性、可靠性和精确性提出了进一步的挑战。为了提高导航系统的性能，以便能够快速、可靠、精准地获取和估计运动信息和姿态信息，一方面需要从硬件角度提高传感器性能或者研制新的智能感知机制，另一方面则要从软件角度提升滤波算法的适应性和鲁棒性。特别是在硬件资源有限的条件下，后者更是具有重要意义。

实际系统总是或多或少受到未知输入或未知干扰的影响。以惯性导航系统为例，陀螺仪的漂移、加速度计的常值偏置、外部环境的振动、安装误差等都可以视为干扰。因此，研究存在未知输入/干扰下的滤波方法即抗干扰滤波，具有积极的理论意义和工程应用价值。作为一种典型的抗干扰滤波方法，未知输入滤波无需对干扰信号的特性做出假设，能够在干扰信号完全未知的情形下给出状态的最优估计，一般为最小方差无偏估计。未知输入滤波最初由美国学者Kitanidis于上世纪八十年代提出，经过三十余年的发展，已建立起了一套较为成熟的理论体系。同时，这些理论成果也成功应用在了导航定位、环境监测、故障检测等相关领域。例如，专利1(火星动力下降段非线性三步滤波方法，申请号为201310341821.7)利用未知输入滤波技术对火星探测器动力下降段的位置和速度进行精确估计；专利2(一种基于智能手机移动检测的城市快速路交通状态估计方法，申请号为201610056587.7)将未知输入滤波应用于城市快速路交通状态估计问题当中，实现交通密度和边界流量的同步估计。

然而，现有的研究成果都是基于系统的过程噪声和量测噪声均服从高斯分布的假设，在最小方差意义下给出状态的最优估计。事实上，非高斯噪声广泛存在于实际系统当中。当随机变量的概率分布具有偏斜非对称、尖峰肥尾等强非高斯特性时，一般不能够再近似为高斯分布，见专利3(一种基于负熵的非高斯线性随机系统滤波方法，申请号为201410124598.5)。特别是对于量测噪声，其概率分布经常出现肥尾现象，这意味着量测噪声取得极端异常值的可能性很大。或者说，系统的量测通道经常受到野值的影响。在这种情况下，现有滤波器的估计性能会变差，甚至有可能导致滤波发散。因此，有必要研究一种针对非高斯随机系统的未知输入滤波器。

发明内容

本发明的技术解决问题是：针对现有的未知输入滤波方法局限于高斯系统的不足，提供一种适合于非高斯系统的新型随机抗干扰滤波器，解决复杂系统的高精度状态估计问题。

本发明的技术解决方案为：一种含未知输入和非高斯量测噪声的随机系统新型抗干扰滤波器，其特征在于，包括以下步骤：首先，基于状态方程得到状态的预测估计与预测误差的协方差矩阵；其次，基于Huber函数构造干扰与状态联合估计的性能指标函数；进而，通过最小化性能指标函数，借助矩阵相关运算和不动点迭代法给出干扰估计与状态的滤波估计，以及相应的估计误差和协方差矩阵；最后，总结归纳出包括状态预测、干扰估计、状态滤波在内的递推滤波器设计流程；具体步骤如下：

第一步，状态预测

(1)一类含有未知输入/干扰的离散时间线性随机系统，其状态空间模型表示如下：

其中，为系统状态变量，为量测输出变量，

为系统的未知输入，也可认为是动态特性未知的干扰，

和

分别表示过程噪声和量测噪声，为互不相关的零均值白噪声且协方差分别为正定矩阵Q_k和R_k，过程噪声ω_k认为服从高斯分布，而量测噪声υ_k为非高斯随机变量，尤其是具有尖峰肥尾的强非高斯特性，A_k，G_k，C_k和H_k为已知的时变矩阵，假设(A_k,C_k)可观且H_k列满秩，系统初始状态x₀服从高斯分布

且与噪声ω_k和υ_k均不相关。

(2)基于状态方程，根据系统运行的前一时刻的干扰估计和状态滤波估计，可得到当前时刻的状态预测估计、预测误差、预测误差的协方差矩阵，分别如下三式所示：

需要指出的是，在k＝0时刻，状态预测不用以上公式，直接由初始条件给定，即

第二步，基于Huber函数的性能指标函数构造

考虑到量测噪声的尖峰肥尾非高斯特性，这里采用Huber函数构造相应的性能指标函数。当量测噪声含有野值时，Huber函数等价于残差的l₁范数，因此对量测野值具有较强的鲁棒性。定义

其第i个分量表示为ξ_k[i]。给出如下的性能指标函数：

其中，ρ(·)为Huber函数，定义如下：

其中τ是函数自变量，γ为一个可调节的参数。

第三步，干扰估计与状态滤波

(1)基于第二步构造的性能指标函数，将其最小化则可得到干扰的估计与状态的滤波估计，即估计问题转化为优化问题为了求解这一优化问题，将指标函数分别对x_k和d_k求导，可以得到以下方程组：

进一步，定义函数ψ(·)为ψ(τ)＝ρ′(τ)/τ，同时定义矩阵Ψ_k＝diag{ψ(ξ_k[1]),…,ψ(ξ_k[m])}，

因此，可以得到如下所示的方程组：

根据分块矩阵求逆公式、矩阵逆引理以及代数运算，可以求解上述方程组。为了表述简便，进一步定义矩阵

则方程组的解表示如下：

(2)注意到在上述结果中，M_k，K_k均与Ψ_k有关，而Ψ_k是x_k和d_k的函数，故上式可看作是关于x_k和d_k的不动点方程。基于不动点迭代法，我们可以得到下述的迭代求解过程：

其中，t表示迭代时标，

和

由前述定义结合给出，迭代初值可以取为

迭代终止条件或者可以通过直接设定迭代次数给出，或者可以给出相应判断准则及阈值。迭代终止后，记

其中t_f表示最后一次迭代的时标。这样，干扰的估计与状态的滤波估计表示为：

同时，我们可以得到相应的估计误差与协方差矩阵如下两式所示：

第四步，滤波器按时刻递推

基于前三个步骤的结果，在系统运行的每一时刻首先进行状态预测获取

接着进行干扰估计与状态滤波获取

进而递推至系统的下一时刻，重复这些步骤，即可得到滤波器的完整工作流程。

本发明与现有技术相比的优点在于：针对同时含有未知输入/干扰和非高斯量测噪声的复杂随机系统，提出了一种新型随机抗干扰滤波方法，通过采用基于Huber函数的未知输入滤波算法，相对于现有的基于最小方差准则的未知输入滤波方法而言，提高了滤波器对量测野值的鲁棒性，增强了滤波器对受扰非高斯随机系统的适应性，使得滤波器能够在复杂情况下实现高精度的状态估计。

附图说明

图1为本发明一种含未知输入和非高斯量测噪声的随机系统滤波器的设计流程图；

图2为本发明一种含未知输入和非高斯量测噪声的随机系统滤波器的不动点迭代法过程示意图。

具体实施方式

如图1所示，本发明具体实现步骤如下(以某一做近似匀加速直线运动的运动体为例来说明方法的具体实现)：

1、状态预测

(1)当运动体做近似匀加速直线运动时，运动学方程如下所示：

量测方程为：

y_k＝s_k+d_k+υ_k

其中s_k，v_k和a_k分别表示运动体在k时刻的位置、速度和加速度，量测信息y_k由位置传感器提供，T为采样周期，d_k表示未知输入/干扰，

(i＝1,2,3)为相互独立的高斯白噪声，υ_k表示零均值方差为R_k的噪声，其具有学生t分布、Laplace分布或者其他尖峰肥尾分布的特性。将运动学方程和量测方程写成紧凑形式，则有：

其中，相关矩阵和向量定义如下：

C_k＝[1 0 0]，H_k＝1

显然，(A_k,C_k)满足可观性条件，H_k满足列满秩要求，ω_k为零均值高斯向量方差阵为

初始状态

且与噪声均不相关。

(2)基于状态方程，根据运动体前一时刻的干扰估计和状态滤波，可得到当前时刻的状态预测、预测误差、预测误差的协方差矩阵分别如下三式所示：

需要指出的是，在k＝0时刻，状态预测不用以上公式，直接由初始条件给定，即

2、基于Huber函数的性能指标函数构造

考虑到量测噪声的尖峰肥尾非高斯特性，使用Huber函数构造相应的性能指标函数。当量测噪声含有野值时，Huber函数等价于残差的l₁范数，因此对量测野值具有较强的鲁棒性。定义

其第i个分量表示为ξ_k[i]。我们给出如下的性能指标函数：

其中，ρ(·)为Huber函数，定义如下：

其中τ是函数自变量，γ为一个可调节的参数。

3、干扰估计与状态滤波

(1)基于上一步构造的性能指标函数，将其最小化则可得到干扰估计与状态滤波，即需要求解优化问题：将指标函数J_k对x_k和d_k分别求导，可以得到以下方程组：

定义函数ψ(τ)＝ρ′(τ)/τ，同时定义矩阵Ψ_k＝diag{ψ(ξ_k[1]),…,ψ(ξ_k[m])}，则上述方程组具有如下形式：

根据分块矩阵求逆公式、矩阵逆引理以及相关代数运算，我们可以求解上述方程组。为了表述简便，进一步定义矩阵

则方程组的解可表示为：

(2)在上述结果中，注意到M_k，K_k均与Ψ_k有关，而Ψ_k是x_k和d_k的函数，故上式可看作是关于x_k和d_k的不动点方程。基于不动点迭代法，我们可以得到下述迭代过程：

其中，t为迭代时标，

和由前述定义结合

给出，迭代初值可以取为

迭代终止条件可以通过直接设定迭代次数给出，或者也可以给出相应判断准则及阈值。迭代终止后，我们记

其中t_f是最后一次迭代的时标，则干扰的估计和状态的滤波估计表示为：

同时，可以计算得到干扰估计误差、状态滤波误差，以及相应的协方差矩阵，如下所示：

4、滤波器按时刻递推

基于前述步骤中的结果，在系统运行的每一时刻首先进行状态预测获取接着进行干扰估计与状态滤波获取

进而递推至下一时刻，重复这些步骤，可得到滤波器的完整工作流程，即可给出复杂环境下运动体的运动信息(包括位置、速度和加速度)的实时高精度估计。

本发明说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (5)

1.一种含未知输入和非高斯量测噪声的随机系统滤波器，其特征在于，包括以下步骤：

第一步，基于系统状态方程进行状态预测，得到状态的预测估计与相应的协方差矩阵；

第二步，基于Huber函数构造干扰与状态联合估计的性能指标函数；

第三步，通过最小化第二步中的性能指标函数，基于矩阵运算和不动点迭代法解算出干扰估计与状态的滤波估计，以及相应的协方差矩阵；

第四步，将前述步骤中的状态预测、干扰估计、状态滤波按时刻递推，最终得到滤波器。

2.根据权利要求1所述的含未知输入和非高斯量测噪声的随机系统滤波器，其特征在于：所述第一步，状态预测过程如下：

(1)一类含有未知输入/干扰的离散时间线性随机系统，其状态空间模型表示如下：

其中，

为系统状态变量，

为量测输出变量，

为系统的未知输入，即动态特性未知的干扰，

和

分别表示过程噪声和量测噪声，为互不相关的零均值白噪声且方差为正定阵Q_k和R_k，过程噪声ω_k认为是高斯分布，而量测噪声υ_k具有强非高斯特性，特别是尖峰肥尾现象，A_k，G_k，C_k和H_k为已知的时变矩阵，(A_k,C_k)可观且H_k列满秩，系统初始状态x₀服从高斯分布

且与噪声均不相关；

(2)基于状态方程，根据前一时刻得到的干扰和状态的估计，得到当前时刻状态的预测估计、预测误差与协方差矩阵分别如下三个公式所示：

在初始k＝0时刻，状态预测直接由系统的初始条件给定，即

3.根据权利要求1所述的含未知输入和非高斯量测噪声的随机系统滤波器，其特征在于：所述第二步，基于Huber函数的性能指标函数构造过程如下：

基于量测噪声的非高斯特性，使用Huber函数构造相应的性能指标函数，Huber函数等价于残差的l₁范数，对于噪声野值具有抑制能力，定义

其第i个分量表示为ξ_k[i]，构造出如下的性能指标函数：

其中，ρ(·)为Huber函数，定义如下：

其中τ是函数自变量，γ为一个可调节的参数。

4.根据权利要求1所述的含未知输入和非高斯量测噪声的随机系统滤波器，其特征在于：所述第三步，干扰估计与状态滤波过程如下：

(1)基于第二步构造的性能指标函数，将其最小化可以得到干扰的估计与状态的滤波估计，即需要求解如下方程：

定义函数ψ(τ)＝ρ′(τ)/τ，同时定义矩阵Ψ_k＝diag{ψ(ξ_k[1]),…,ψ(ξ_k[m])},

经过相关矩阵运算，方程的解表示为：

(2)M_k、K_k均与Ψ_k有关，Ψ_k是x_k和d_k的函数，上式是关于x_k和d_k的不动点方程，基于不动点迭代法，得到下述迭代过程：

其中，t为迭代时标，

和

由前述定义结合给出，迭代初值可以取为

迭代终止条件可直接设定迭代次数或者给出判断准则及阈值；迭代终止后，记

其中t_f是最后一次迭代的时标，则干扰估计与状态滤波表示为：

同时，相应的估计误差与协方差矩阵分别如下两式所示：

5.根据权利要求1所述的含未知输入和非高斯量测噪声的随机系统滤波器，其特征在于：所述第四步，滤波器的按时刻递推过程如下：

在每一时刻首先进行状态预测获取状态的预测估计

及相应的协方差

接着进行干扰估计与状态滤波获取

及相应的协方差

进而递推至下一时刻，重复这些步骤，最终即得到了滤波器的完整工作流程。

Images