CN110649911B - 一种基于α散度的分布式非线性卡尔曼滤波方法 - Google Patents

Abstract

本发明本发明属于信号处理领域，具体提出一种基于α散度的分布式非线性卡尔曼滤波方法，能有效减小系统自身噪声和外部噪声的影响，得到对真实值高精度的估计，可用于解决信号处理中的滤波问题和参数估计问题，尤其是目标跟踪问题。本发明中，首先，通过最小化每个节点的中间近似后验分布函数与真实后验分布函数之间的α散度以获取中间状态估计；然后，通过最小化每个节点的最终后验分布近似函数与其邻居节点的中间后验分布近似函数之间前向KL散度的凸组合来计算每个节点的最终状态估计结果。本发明不同于现有的基于最小方差准则的分布式粒子滤波和基于最小均方误差准则的分布式扩展卡尔曼滤波，本有效提高了非线性分布式卡尔曼滤波的稳态性能。

Description

一种基于α散度的分布式非线性卡尔曼滤波方法

技术领域

本发明属于信号处理领域，涉及信号处理领域的目标跟踪问题，尤其涉及到分布式无线 传感器网络上目标跟踪问题，具体为一种基于α散度的分布式非线性卡尔曼滤波方法。

背景技术

卡尔曼滤波算法不仅可以对平稳的一维随机过程进行估计，也可以对非平稳的多维随机 过程进行估计，而且卡尔曼滤波算法是递推的，具有存储量小,收敛速度快，实时处理速度快 等优点，因此卡尔曼滤波在复杂系统中有更广泛的应用，如导航、目标跟踪、定位等。此外，卡尔曼滤波算法也被用于动态系统的预测，如星体运动轨迹，商品交换价格的变化趋势。

目前，对于线性卡尔曼滤波的研究已经取得了巨大的进展，然而对非线性卡尔曼滤波的 研究仍有所欠缺。在卡尔曼滤波中常用后验分布表示其状态估计，由于非线性卡尔曼滤波的 高斯后验分布不具有封闭形式，无法直接求解其后验分布，这就为非线性卡尔曼的研究带来了巨大困难。

随着无线传感器技术的快速发展，近年来，分布式网络成为焦点，大量低成本的无限传 感器被广泛应用于分布式网络中。在分布式网络中，每个节点仅同其邻点进行信息交流，如 此可以大大提高网络的可扩展性和灵活性，能大量节省通信资源，在环境监测，救灾管理，参数估计等领域有着广泛应用。然而目前，分布式的扩散融合技术在分布式线性卡尔曼滤波 中的应用较为成熟，而对于非线性动力系统中的研究尚不够充分。现有的分布式非线性卡尔 曼滤波算法有分布式扩展卡尔曼滤波和分布式粒子滤波等。扩展卡尔曼滤波算法对非线性动 力学系统做了线性近似处理，再利用线性卡尔曼滤波系统的方法进行处理，引入了线性化误差，具有收敛速度慢等局限性。

发明内容

本发明的目的在于提出一种基于α散度的分布式非线性卡尔曼滤波方法，采用一个指数 族的分布函数q(x)，来近似真实的后验分布，并采用α散度衡量两个分布之间的差距，α散 度越小表示两个分布越接近；本发明通过直接优化α散度的方法来计算近似分布函数，从而得到对目标真实状态的估计结果。

为实现上述目的，本发明采用的技术方案如下：

一种基于α散度的分布式非线性卡尔曼滤波方法，包括以下步骤：

步骤1.针对节点r，根据其t-1的后验分布函数与状态转移方程计算t 时刻的先验分布函数/>

步骤2：设定t时刻的采样的建议分布函数以及初始化中间近似 分布函数/>

步骤3：节点r根据建议分布进行采样，使得每个粒子都独立同分布于建议分布：

其中，S为输入采样粒子数，表示节点r在t时刻采样的第s个粒子的状态向量；

步骤4：根据下式计算节点k的粒子权重，并进行归一化处理：

其中，表示节点r在t时刻第s个采样粒子的权重，W_r,t表示节点r在t时刻所有采样粒子 的权重之和，y_u,t表示节点r的邻居节点u在t时刻的观测值；/>表示邻居节点u对/>观测值的预测，由观测方程计算得到；/>表示节点r对/>的预测，由状态转移方 程计算得到；

步骤5：根据粒子权重计算近似中间后验分布函数其服从高斯分布

步骤6：在分布式网络中，将每个节点得到的近似中间后验分布函数进行扩散， 计算得到节点r的最终状态估计q_r(x_t)：

其中，a_u,r表示在扩散时邻点u对当前节点r的权重，满足

并计算q_r(x_t)的均值μ_r,t和协方差矩阵

μ_r,t即为节点r的状态估计向量。

本发明的有益效果在于：

本发明提出的一种基于α散度的分布式非线性卡尔曼滤波方法，具有如下优点：

1.本发明提出的方法可用于非线性动力系统，相比经典的线性卡尔曼滤波算法，应用范 围更广；

2.本发明提出的是一种分布式扩散策略，每个节点处只需同其邻点交流信息，允许各个 节点同时进行处理，不需要将各个节点的信息送至融合中心进行处理，需要的通信能量少， 运算效率高；

3.本发明提出的分布式算法比相应的集中式算法的稳健性强。对于集中式而言，当处理 中心出现问题时，整个系统将无法正常工作，分布式算法能有效避免融合中心出现问题时导 致整个系统崩溃的风险；

4.本发明在计算每个节点的中间状态估计时，基于蒙特卡洛技术直接优化的无偏算法直 接优化近似分布与真实后验分布之间的α散度，不需要对非线性函数做线性近似；

5.本发明使用的α散度可以看作是KL散度的推广，其中涉及的参数α满足α∈(-∞,∞)。 当α＝1时，α散度可退化为后向KL散度；当α＝0时，α散度可退化为前向KL散度；当α＝0.5， α散度可退化为海林格距离(Hellinger distance)。本发明可以根据具体的应用场景选择 合适的α的取值，通过调整α的大小控制散度范围，可能获得更优的估计结果(不同α值近 似分布结果如图1所示)。

附图说明

图1为本发明基于α散度的分布式非线性卡尔曼滤波方法中每个节点的流程示意图。

图2为不同α值的近似分布结果。

图3为实施例中采用的分布式网络拓扑结构(以σ_Q＝0.05网络中有30个节点为例)。

图4为实施例中本发明方法与分布式扩展卡尔曼滤波，分布式粒子滤波和集中式α散度 方法在某一次蒙特卡罗实验中一个节点的跟踪结果图。

图5为实施例中过程噪声取σ_Q＝0.05时，本发明提出的基于α散度的分布式非线性卡尔 曼滤波与分布式扩展卡尔曼滤波，分布式粒子滤波和集中式α散度方法的位置均方误差(MS E)的对比图。

图6为实施例中过程噪声取σ_Q＝0.05时，本发明提出的基于α散度的分布式非线性卡尔 曼滤波与分布式扩展卡尔曼滤波，分布式粒子滤波和集中式α散度方法的速度MSE对比图。

表1为实施例中不同的过程噪声参数下的基于α散度的分布式非线性卡尔曼滤波，分布 式扩展卡尔曼滤波，分布式粒子滤波和集中式α散度方法的平均误差结果。

具体实施方式

下面结合附图和实施例对本发明作进一步说明。

本发明中，考虑有个节点的分布式网络结构，每个节点r都有一个传感器可对目标进 行跟踪，与节点r直接相连的节点称之为他的邻点(每个节点都与自身相连)表示为u，将r 的邻点网络记为N_r；本发明提出的基于α散度的分布式非线性卡尔曼滤波方法主要基于扩散 策略，可分为两步：在自适应阶段得到每个节点的中间状态估计，再通过组合阶段将中间状 态估计在每个节点的邻域内进行扩散；每个节点根据其所有邻点观测信息累积量得到的的真实后验分布为/>其中，x_t表示t时刻的状态向量；/>表示节点k的所有邻居 节点在1～t时刻观测数据的集合：

每个节点的中间状态估计用中间近似后验分布函数来描述，每个节点的最终状态 估计用近似后验分布函数q_r(x_t)来描述，q_r(x_t)的均值即为对目标的状态估计值。

用α散度来衡量每个节点的真实后验分布与近似后验分布函数之间的差距，定义每个节 点的α散度表达式为：

其中，α为预设参数，如图2所示，α的取值范围为{-∞,+∞}/{0}。假设各邻节的观测信息 彼此独立，即之间彼此独立。

优化每个节点的α散度，得到在α散度意义下与真实后验分布最接近的中间近似后验分 布函数令/>使α散度最小即对/>与做矩匹配，得到：/>E表示期望；由于此等式两边均有/>不能直接求解，因此本发明通过一次迭代的方法得到/>

本实施例提供一种基于α散度的分布式非线性卡尔曼滤波方法，其流程如图1所示，具 体包括以下步骤：

步骤1.针对节点r，根据其t-1的后验分布函数与状态转移方程计算t 时刻的先验分布函数/>

步骤2：设定t时刻的采样的建议分布函数以及初始化中间近似 分布函数/>

步骤3：节点r根据建议分布进行采样，使得每个粒子都独立同分布于建议分布：

其中，S为输入采样粒子数，表示节点r在t时刻采样的第s个粒子的状态向量；

步骤4：根据下式计算节点k的粒子权重，并进行归一化处理：

其中，表示节点r在t时刻第s个采样粒子的权重，W_r,t表示节点r在t时刻所有采样粒子 的权重之和，y_u,t表示节点r的邻居节点u在t时刻的观测值；/>表示邻居节点u对/>观测值的预测，由观测方程计算得到；/>表示节点r对/>的预测，由状态转移方 程计算得到；

步骤5：根据上式求得的粒子权重计算中间近似后验分布其服从高斯分布/>

步骤6：在分布式网络中，将每个节点得到的中间后验分布近似函数进行扩散， 最小化节点r的最终近似后验分布与其所有邻点的中间近似后验分布函数的凸组合之间的前 向KL散度，定义其KL散度为：

其中，a_u,r表示在扩散时邻点u对当前节点r的贡献，满足

令上述前向KL散度最小，得到：

因为上式中均服从高斯分布，因此q_r(x_t)也服从高斯分布，其均值和协方差矩阵计算 如下：

得到每个节点的后验分布近似函数q_r(x_t)，其均值即为该节点的状态估计向量；

步骤7：进行时间迭代，得到下一时刻的状态估计结果。

仿真测试

将本发明提出的方法用于分布式网络的目标跟踪中，与分布式扩展卡尔曼滤波，分布式 粒子滤波和集中式α散度方法进行比较。整个分布式网络的拓扑结构共30个节点，拓扑结构 如图3所示，每个节点处一个传感器分别对目标进行跟踪，各节点的观测噪声功率一样。目 标运动的状态空间模型如下：

x_t＝F_tx_t-1+w_t,w_t～N(0,Q_t)

y_t＝h(x_t)+v_t,v_t～N(0,R_t)

其中，F_t表示状态转移矩阵，h(x_t)表示观测函数，w_t是状态噪声，v_t是观测噪声，w_t与v_t是均值为零，协方差矩阵分别为Q_t和R_t的高斯白噪声：

其中s表示传感器的位置信息。本发明中采取恒定的测量速率，令Δt＝1；观测噪声的协方差 矩阵为σ_R＝20；状态向量是四维的：/>包括两个方向的位置信 息，第二维和第四维表示两个方向的速度信息。采样粒子数为1000，迭代次数500次，蒙特卡洛实验100次，α散度参数的取值为α∈{0.5,1}，过程噪声的取值为 σ_Q∈{0.001,0.005,0.01,0.05,0.1,0.3}。

将本发明提出的基于α散度的分布式卡尔曼滤波方法，α分别取0.5和1(图中分别标注 为“分布式AKF,α＝0.5”和“分布式AKF,α＝1”)与“Nonlinear Kalman FilteringWith Divergence Minimization”所提出的集中式基于α散度卡尔曼滤波方法(图中标注为“集中式 AKF,α＝0.5”)，“Diffusive Particle Filtering For DistributedMultisensor Estimation”提出的分布式粒子滤波方法(图中标注为“分布式PF”)和“Distributed nonlinear Kalman filtering with applications to wirelesslocalization”提出的分布式扩展卡尔曼滤波(图中标注为“分 布式EKF”)仿真结果如图4、5、6和表1所示。

图4表明本发明提出的基于α散度分布式非线性卡尔曼滤波方法同其他几种方法都能有 效地对目标的位置进行跟踪。图5表明本发明提出的分布式方法与其他方法相比，本发明方 法的位置状态估计MSE曲线的收敛速度与其他方法相差不大，稳态性能略差于集中式的方 法，但比其他两种分布式方法的稳态误差都小，并且在α＝0.5时比α＝1时(退化为KL散度) 稳态误差更小。

图6表明在过程噪声取值为0.05时，本发明提出的方法与其他三种方法相比，他们的速 度状态估计MSE曲线的收敛速度大致都相同，但本发明提出的方法的稳态误差比其他方法都小，性能更好，甚至比集中式的方法都好，并且在α＝0.5时效果更好。表1表明在过程噪声 较小时，取值在一定范围内，本发明提出的方法相比其他两种分布式方法的稳态误差更好， 但会略小于集中式的方法，α＝0.5比α＝1的平均误差更小。

表1

以上所述，仅为本发明的具体实施方式，本说明书中所公开的任一特征，除非特别叙述， 均可被其他等效或具有类似目的的替代特征加以替换；所公开的所有特征、或所有方法或过 程中的步骤，除了互相排斥的特征和/或步骤以外，均可以任何方式组合。

Claims (1)

1.一种基于α散度的分布式非线性卡尔曼滤波方法，该方法用于分布式网络的目标跟踪中，设置分布式网络有个节点，每个节点r都有一个传感器对目标进行跟踪；该方法包括以下步骤：

步骤1.针对节点r，根据其t-1时刻的后验分布函数与状态转移方程计算t时刻的先验分布函数/>

步骤2：设定t时刻的采样的建议分布函数以及初始化中间近似分布函数/>

步骤3：节点r根据建议分布进行采样，使得每个粒子都独立同分布于建议分布：

其中，S为输入采样粒子数，表示节点r在t时刻采样的第s个粒子的状态向量；

步骤4：根据下式计算节点r的粒子权重，并进行归一化处理：

其中，表示节点r在t时刻第s个采样粒子的权重，W_r,t表示节点r在t时刻所有采样粒子的权重之和，y_u,t表示节点r的邻居节点u在t时刻的观测值；/>表示邻居节点u对观测值的预测，由观测方程计算得到；/>表示节点r对/>的预测，由状态转移方程计算得到；

步骤5：根据粒子权重计算近似中间后验分布函数其服从高斯分布

步骤6：在分布式网络中，将每个节点得到的近似中间后验分布函数进行扩散，计算得到节点r的最终状态估计q_r(x_t)：

其中，a_u,r表示在扩散时邻点u对当前节点r的权重，满足

并计算q_r(x_t)的均值μ_r,t和协方差矩阵

μ_r,t即为节点r的状态估计向量。