CN104467742A - 基于高斯混合模型的传感器网络分布式一致性粒子滤波器 - Google Patents

Abstract

本发明提供了一种非线性非高斯系统中多目标情况下的分布式一致性状态估计方法IC-GMM-DPF。该方法基于期望最大化的高斯混合模型GMM，并结合联合概率数据互联策略，设计了GMM中参数一致性机制。本发明方法联合解决了全网一致性估计、多目标数据互联和网络盲节点问题，在错误量测与杂波较多的复杂场景下仍然具有较强的状态估计精度及鲁棒性。

Description

基于高斯混合模型的传感器网络分布式一致性粒子滤波器

技术领域

本发明涉及传感器网络的信息融合系统，提供一种非线性非高斯系统中多目标环境下的分布式一致性状态估计方法。

背景技术

在实际应用环境中，分布式传感器网络的监测区域通常会包含多个目标。当目标较为分散时，可以采用单目标状态估计方法分别进行处理。然而，被估对象也可能以密集多目标甚至编队目标的形式出现，比如敌军探测与跟踪系统对敌方装甲车队的监控、无人机机群对空中编队目标的跟踪等。

现有大多数一致性状态估计方法很少涉及多目标场景，或者假设具有完全理想的数据互联结果。针对线性高斯系统中多目标环境下的状态估计问题，有关文献推导了基于联合概率数据互联的卡尔曼一致性滤波器，提出了多目标分布式一致性状态估计方法JPDA-KCF，但是未考虑网络中的盲节点问题；有关学者提出了多目标信息一致性方法(Multi-targetInformation Consensus,MTIC)，联合处理了数据互联及网络盲节点等问题，但是MTIC假设网络中的传感器节点和目标数量已知，只适用于先验预知的固定网络，方法可扩展性较弱，且没有考虑非线性非高斯系统中的一致性状态估计问题。

发明内容

1.要解决的技术问题

本发明针对多目标环境中的一致性状态估计问题，基于期望最大化及高斯混合滤波思想，并结合概率数据互联方法JPDA，提出了基于高斯混合模型的传感器网络分布式一致性粒子滤波器IC-GMM-DPF，联合解决了多目标数据互联、网络盲节点及一致性估计问题。

2.技术方案

本发明所述的基于信息一致性及高斯混合模型的分布式粒子滤波器，如图2所示，包括以下技术措施：1.基于提示分布产生当前时刻所需的滤波粒子；2.在每个节点执行分布式EM方法计算当前时刻全局分布的参数集；3.针对每个目标的GMM分布产生新的滤波粒子集；4.节点集内的每个传感器都执行一致性信息迭代与融合，获取到每个目标对应的新的滤波粒子集。

3.有益效果

本发明相比背景技术具有如下的优点：

(1)通过增加网络节点可以促使状态估计性能的提高，IC-GMM-DPF对网络节点失效具有较强的鲁棒性；

(2)IC-GMM-DPF将网络的目标状态估计任务以分布式的形式分配到每个传感器节点，降低了计算复杂度，单个传感器节点仅需少量滤波粒子即可获取预期的估计性能；

附图说明

图1：基于EM的GMM参数估计；

图2：IC-GMM-DPF方法流程图；

图3：仿真场景；

图4：一致性迭代次数变化时的估计性能比较；

图5：杂波变化时的估计性能比较；

图6：网络节点个数变化时的性能比较；

图7：目标个数变化时的性能比较。

具体实施方式

以下结合说明书附图对本发明作进一步详细描述。参照说明书附图，本发明目标状态更新的单次循环方式分以下几个步骤：

1、问题描述

考虑某一网络具有N_S个同构传感器节点，其节点之间的通信网络可有无向图来表示。集合包含了图中的所有顶点，代表网络中的所有传感器节点；集合包含图中所有的边，代表不同节点之间的可行通信链接。所有与节点S_i具有直接通信链接的节点称为S_i的邻居节点，并构成集合假设网络中存在N_T个可观测目标，构成的集合可表示为

$x_{k+1}^j=\mathrm{\Φ}{x}_k^j+w_k^j---\left(1\right)$

其中，为状态转移矩阵，过程噪声服从高斯分布

k时刻，传感器S_i获取到Γ_i,k个关于目标的量测，记为传感器并不预知量测与目标的对应关系。假设量测来源于目标T^j，由以下观测方程得到

$z_{i,k}^{\mathrm{\τ}}=H_{i,k}^j{x}_{i,k}^j+v_{i,k}^j---\left(2\right)$

其中，为传感器S_i对于目标T^j的观测矩阵，量测噪声设定为服从零均值高斯分布的随机变量，且方差为

每个传感器节点保持着对感兴趣目标的先验或预测状态估计状态估计方差为需要说明的是，本章仍然以状态方差的逆(称为信息矩阵)应用于滤波的整个过程中，表示形式为 $J_{i,k}^j={\left(P_{i,k}^j\right)}^{-1}.$

2、IC-GMM-DPF方法推导

JPDA的主要思想就是通过顺序贝叶斯估计迭代方程采用递归的方式更新每个目标的边缘滤波分布：

$p^j\left(x_k^j\right|z_{k-1}^j)=\∫p^j\left(x_k^j\right|x_{k-1}^j)p^j\left(x_{k-1}^j\right|z_{k-1}^j)d{x}_{k-1}^j---\left(3\right)$

其中表示k时刻目标j的状态，表示直到k-1时刻所获取的关于目标j的量测数据，表示目标j的先验概率密度，表示关于目标j的似然函数。目标j的联合似然函数表示为

$p^j\left(x_k^j\right|z_k^j)=\underset{i=1}{\overset{N_S}{\mathrm{\Π}}}[{\mathrm{\β}}_i^{\mathrm{j\τ}}+\underset{\mathrm{\τ}=1}{\overset{{\mathrm{\Γ}}_i}{\mathrm{\Σ}}}{\mathrm{\β}}_i^{\mathrm{j\τ}}{P}_{i,T}\left(z_{i,k}^j\right|x_k^j)]---\left(4\right)$

其中表示传感器i中目标j与量测τ互联的概率，τ＝1,...,Γ_i，表示没有量测与目标j关联(此时刻没有观测到目标j)。是传感器i关于目标j的似然函数。假定各传感器关于目标j的似然函数相互独立，基于式(4)所定义的目标似然函数，有

$p^j\left(x_k^j\right|z_k^j)\∝p^j\left(z_k\right|x_k^j)p^j\left(x_k^j\right|z_{1:k-1})---\left(5\right)$

然后，计算边缘互联的后验概率

${\mathrm{\β}}_i^{\mathrm{j\τ}}=p({\tilde{r}}_{i,k}^j=\mathrm{\τ}|z_{1:k})=\underset{\{{\widetilde{\mathrm{\λ}}}_{i,k}\∈{\widetilde{\mathrm{\Λ}}}_{i,k};{\tilde{r}}_{i,k}^j=\mathrm{\τ}\}}{\mathrm{\Σ}}p\left({\widetilde{\mathrm{\λ}}}_{i,k}\right|z_{1:k})---\left(6\right)$

其中，是传感器i的联合目标量测互联假设序列。联合互联假设的后验概率可以表示为：

其中，表达式在文献中已经给出；表示为第i个传感器量测空间的体积，定义式为R_i,S是节点i的传感半径， $I_i=\{\mathrm{\τ}\∈\{1,...,{\mathrm{\Γ}}_i\}:r_i^j\≠0\}.$ 表达式(7)中的表示传感器i基于目标j的信息计算得到关于量测τ的预测似然函数，其标准形式为

$p^j\left(z_{i,k}^j\right|z_{1:k-1})=\∫p_{i,T}\left(z_{i,k}^j\right|x_k^j)p^j\left(x_k^j\right|z_{1:k-1})d{x}_k^j---\left(8\right)$

在Monte Carlo滤波框架中，式(7)中的预测似然函数可以基于提议分布的Monte Carlo样本进行近似估计。假定可以基于前一时刻边缘滤波分布的估计分布以分布式的方式近似得到关于目标j的样本序列那么，在当前时间步长中，可以基于优化设计的提议分布产生关于目标状态的滤波粒子，即

$x_k^{j,\left(m\right)}~q^j\left(x_k^j\right|x_{k-1}^{j,\left(m\right)},z_k),m=1,...,M_p---\left(9\right)$

为了降低状态估计的计算复杂度，将先验分布作为目标j的提议分布，可以节省方法的运行时间。因此，提议分布可以表示为

$q^j\left(x_k^j\right|x_{k-1}^{j,\left(m\right)},z_k)=p^j\left(x_k^{j,\left(m\right)}\right|x_{k-1}^{j,\left(m\right)})---\left(10\right)$

于是，基于Monte Carlo样本，式(7)中预测似然函数的估计为

$p^j\left(z_{i,k}^j\right|z_{1:k-1})\≈\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_k^{j,\left(m\right)}{p}_{i,T}\left(z_{i,k}^j\right|x_k^{j,\left(m\right)})---\left(11\right)$

根据式(10)中的定义，重要性权值等于预测权重，即重要性权值的更新方程为

${\mathrm{\α}}_k^{j,\left(m\right)}\∝{\mathrm{\ω}}_{k-1}^{j,\left(m\right)}\frac{p^j\left(x_k^{j,\left(m\right)}\right|x_{k-1}^{j,\left(m\right)})}{q^j\left(x_k^j\right|x_{k-1}^{j,\left(m\right)},z_k)}---\left(12\right)$

且有

$\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_k^{j,\left(m\right)}=1---\left(13\right)$

此时，可以将预测似然函数的近似估计代入式(7)，得到联合互联后验概率的近似估计。更进一步地，基于式(6)可以计算得到关于边缘目标与量测的互联后验概率的近似估计。而上述近似估计可以用于式(5)中，从而估算出目标似然函数。最后，设置新的重要性权值为

$\left\{\begin{array}{c}{\mathrm{\ω}}_k^{j,\left(m\right)}\∝{\mathrm{\ω}}_{k-1}^{j,\left(m\right)}\left(z_k\right|x_k^{j,\left(m\right)})\\ \underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\ω}}_k^{j,\left(m\right)}=1\end{array}\right.---\left(14\right)$

由上所述，基于新的重要性权值以及当前时刻的边缘滤波分布样本序列的计算方式转化为分布式解法。

3、高斯混合模型中的参数一致性

以下为了标记简化，隐去离散系统中的时间指数k。用于高斯混合的EM方法能够以分布式的方式在每个节点计算每个目标的局部充分统计特征，而基于前一时刻本地及邻居关于目标的全局统计特征可以在每个节点中计算得到下一时刻关于目标的全局综合统计特征。在此，假定传感器网络由N_S个独立节点组成，跟踪J个目标，其中传感器i具有Γ_i个量测，并且在每个滤波时刻产生M_p个滤波粒子用于估计每个目标。假定关于目标状态后验分布可由M_g个高斯项混合得到，如下式所示：

$p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\Ψ}}^j)=\underset{g=1}{\overset{M_g}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g}^{j}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^j,{\mathrm{\Σ}}_g^j)---\left(15\right)$

其中表示关于目标j的所需估计的分布式参数集合。此外，高斯混合模型中第g个高斯项的估计状态服从均值为方差为的高斯分布：

$p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^j,{\mathrm{\Σ}}_g^j)=\frac{1}{\sqrt{2\mathrm{\π}}{\left|{\mathrm{\Σ}}_g^j\right|}^{1/2}}{e}^{-\frac{1}{2}{(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_g^j)}^T{\left({\mathrm{\Σ}}_g^j\right)}^{-1}(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_g^j)}---\left(16\right)$

传感器i关于目标j的局部统计量为

${\mathrm{\α}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{n=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}---\left(17\right)$

${\mathrm{\β}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}{x}_i^{j,\left(m\right)}---\left(18\right)$

${\mathrm{\γ}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_{i,g}^{j,q}){(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_{i,g}^{j,q})}^T---\left(19\right)$

其中式(17)、(18)、(19)中的指数q表示分布式GMM方法的迭代次数，是传感器i中第g个高斯项所估计的局部均值。后验概率由以下方程进行解算：

${\mathrm{\α}}_{i,g,m}^{j,q}=\frac{{\mathrm{\α}}_{i,g}^{j,q-1}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^{j,q-1},{\mathrm{\Σ}}_g^{j,q-1})}{\underset{g=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g}^{j,q-1}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^{j,q-1},{\mathrm{\Σ}}_g^{j,q-1})}---\left(20\right)$

各节点局部统计量的计算称为期望步骤(E-step)。在此基础上，基于节点自身局部统计量及邻居的全局统计量，可以计算得到该节点的全局统计量。不妨以表示j的局部统计量，表示j的全局统计量，于是，在一致性处理阶段，采用离散一致性滤波器对每个节点的全局统计量进行近似估计。

其中ζ_i是表示传感器i的信息更新速率，其取值方式类似于前几章所述的一致性速率因子，在此不在赘述。在最大化步骤(M-step)中，每个传感器节点采用以下方程计算得到全局分布所需的参数。

$\left\{\begin{array}{c}{\mathrm{\μ}}_g^{j,q+1}=\frac{{\mathrm{\χ}}_{i,g}^{j,q+1}\left(2\right)}{M_p{\mathrm{\χ}}_{i,g}^{j,q+1}\left(1\right)}\\ {\mathrm{\Σ}}_g^{j,q+1}=\frac{{\mathrm{\χ}}_{i,g}^{j,q+1}\left(3\right)}{M_p{\mathrm{\χ}}_{i,g}^{j,q+1}\left(1\right)}\end{array}\right.---\left(22\right)$

基于EM方法的GMM参数估计过程如图1所示。

4、IC-GMM-DPF方法描述

图2显示了IC-GMM-DPF的方法流程。在IC-GMM-DPF中，首先从每个传感器中的目标初始边缘分布开始，基于式(9)所示的提议分布产生当前时刻所需的滤波粒子；然后，在每个节点执行分布式EM方法(迭代q次)并计算当前时刻全局分布的参数集；接下来，针对每个目标的GMM分布产生新的滤波粒子集；之后，一致性节点集内的每个传感器都执行带有JPDA门控阈值的一致性信息迭代与融合，进而获取到每个目标对应的新的滤波粒子集(包括粒子对应的权重)。由于滤波粒子集的获取过程主要基于对应目标的边缘后验分布并通过了一致性信息迭代及处理，因而此过程是分布式的。

5、仿真比较与分析

本章设置非线性非高斯系统进行性能比较，参与仿真比较的方法包括：

(1)基于信息一致性的联合概率数据互联滤波方法IC-GMM-DPF；

(2)基于联合概率数据互联的卡尔曼一致性滤波方法JPDA-KCF；

(3)基于真实数据互联的集中式卡尔曼滤波方法GT-CF-KF；

5.1仿真设置

采用典型的二维观测传感器网络空中管制系统作为仿真场景，其中某目标在平面区域内作机动转弯运动，并定义目标状态为(x_k,y_k)为k时刻目标位置，为k时刻目标移动速度；为未知的角速度；目标的非线性状态转移模型为

x_k＝Φx_k-1+w_k-1

且有

$\mathrm{\Φ}=\left(\begin{array}{ccccc}1& \sin \mathrm{\Ω}{T}_s& 0& \frac{\cos \mathrm{\Ω}{T}_s-1}{\mathrm{\Ω}}& 0\\ 0& \cos \mathrm{\Ω}{T}_s& 0& -\sin \mathrm{\Ω}{T}_s& 0\\ 0& \frac{1-\cos \mathrm{\Ω}{T}_s}{\mathrm{\Ω}}& 1& \frac{\sin \mathrm{\Ω}{T}_s}{\mathrm{\Ω}}& 0\\ 0& \sin \mathrm{\Ω}{T}_s& 0& \cos \mathrm{\Ω}{T}_s& 0\\ 0& 0& 0& 0& 1\end{array}\right)$

其中过程噪声Q＝diag(10，10，1，1,1)；目标数量N_T＝2，单个时间步长为T_s＝1s，初始状态为

$P_i^{1-}\left(1\right)=P_i^{2-}\left(1\right)=\mathrm{diag}\left(\mathrm{100,25,100,25,1.7}\right)$

网络中各传感器的非线性量测模型为

$z_{i,k}=\left(\begin{array}{c}\sqrt{{(x_k-x_{S_i})}^2+{(y_k-y_{S_i})}^2}\\ \frac{(x_k-x_{S_i}){\stackrel{\·}{x}}_k+(y_k-y_{S_i}){\stackrel{\·}{y}}_k}{\sqrt{{(x_k-x_{S_i})}^2+{(y_k-y_{S_i})}^2}}\end{array}\right)+v_{i,k}$

采用两项零均值高斯混合项来表示v_i,k，并且用一项方差较大的高斯项表示异常值，例如

图3给出了仿真场景示例。假定网络中每个传感器可以获取到关于两个目标的初始状态及初始方差，其他参数的设置及方法比较说明与第一种情况相同。

5.2性能分析与比较

图4至图7分别给出了各方法估计性能随一致性迭代次数、杂波数量、网络节点数量、网络目标个数等参数的变化情况。从整体上看，本发明方法IC-GMM-DPF仍然在同等条件下具有相对较高的估计精度，在恶劣环境下的估计性能明显优于同类方法。同时可以注意到，IC-GMM-DPF依然能够实现近似最优化估计，其性能保持收敛于集中式估计GT-CF-KF。上述现象说明，IC-GMM-DPF通过EM及GMM方法对非高斯噪声的近似，且以粒子滤波方法解决了非线性模型的近似问题，因此在局部滤波估计上仍然能够保持较高的估计精度，同时也降低了数据互联给最终估计带来的负面影响(状态估计与数据互联相互牵制，精确的状态估计结果将提高数据互联的成功率)。总体上看，IC-GMM-DPF在非线性非高斯系统中具有较强的适应性。

Claims (2)

1.基于高斯混合模型的传感器网络分布式一致性粒子滤波方法，其特征在于包括以下步骤：

(1)基于提示分布产生当前时刻所需的滤波粒子；

(2)在每个节点执行分布式EM方法计算当前时刻全局分布的参数集，实现高斯混合模型中的参数一致性；

(3)针对每个目标的GMM分布产生新的滤波粒子集；

(4)节点集内的每个传感器都执行一致性信息迭代与融合，获取到每个目标对应的新的滤波粒子集。

2.根据权利要求1所述的方法，其特征在于步骤(2)采用以下技术措施：

假定传感器网络由N_S个独立节点组成，跟踪J个目标，其中传感器i具有Γ_i个量测，并且在每个滤波时刻产生M_p个滤波粒子用于估计每个目标，假定关于目标状态后验分布可由M_g个高斯项混合得到，如下式所示：

$p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\Ψ}}^j)=\underset{g=1}{\overset{M_g}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g}^{j}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^j,{\mathrm{\Σ}}_g^j)---\left(1\right)$

其中表示关于目标j的所需估计的分布式参数集合，高斯混合模型中第g个高斯项的估计状态服从均值为方差为的高斯分布：

$p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^j,{\mathrm{\Σ}}_g^j)=\frac{1}{\sqrt{2\mathrm{\π}}|{\mathrm{\Σ}}_g^j{|}^{1/2}}{e}^{-\frac{1}{2}{(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_g^j)}^T{\left({\mathrm{\Σ}}_g^j\right)}^{-1}(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_g^j)}---\left(2\right)$

传感器i关于目标j的局部统计量为

${\mathrm{\α}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{n=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}---\left(3\right)$

${\mathrm{\β}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}{x}_i^{j,\left(m\right)}---\left(4\right)$

${\mathrm{\γ}}_{i,g}^{j,q}=\frac{1}{M_p}\underset{m=1}{\overset{M_p}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g,m}^{j,q}(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_{i,g}^{j,q}){(x_i^{j,\left(m\right)}-{\mathrm{\μ}}_{i,g}^{j,q})}^T---\left(5\right)$

其中式(3)、(4)、(5)中的指数q表示分布式GMM方法的迭代次数，是传感器i中第g个高斯项所估计的局部均值，后验概率由以下方程进行解算：

${\mathrm{\α}}_{i,g,m}^{j,q}=\frac{{\mathrm{\α}}_{i,q}^{j,q-1}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^{j,q-1},{\mathrm{\Σ}}_g^{j,q-1})}{\underset{g=1}{\overset{M_g}{\mathrm{\Σ}}}{\mathrm{\α}}_{i,g}^{j,q-1}p\left(x_i^{j,\left(m\right)}\right|{\mathrm{\μ}}_g^{j,q-1},{\mathrm{\Σ}}_g^{j,q-1})}---\left(6\right)$

各节点局部统计量的计算称为期望步骤(E-step)，在此基础上，基于节点自身局部统计量及邻居的全局统计量，可以计算得到该节点的全局统计量，不妨以表示j的局部统计量，表示j的全局统计量，于是，在一致性处理阶段，采用离散一致性滤波器对每个节点的全局统计量进行近似估计，

其中ζ_i是表示传感器i的信息更新速率。