2.845

2023影响因子

(CJCR)

  • 中文核心
  • EI
  • 中国科技核心
  • Scopus
  • CSCD
  • 英国科学文摘

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

带时间相关乘性噪声多传感器系统的分布式融合估计

马静 杨晓梅 孙书利

马静, 杨晓梅, 孙书利. 带时间相关乘性噪声多传感器系统的分布式融合估计. 自动化学报, 2023, 49(8): 1745−1757 doi: 10.16383/j.aas.c210147
引用本文: 马静, 杨晓梅, 孙书利. 带时间相关乘性噪声多传感器系统的分布式融合估计. 自动化学报, 2023, 49(8): 1745−1757 doi: 10.16383/j.aas.c210147
Ma Jing, Yang Xiao-Mei, Sun Shu-Li. Distributed fusion estimation for multi-sensor systems with time-correlated multiplicative noises. Acta Automatica Sinica, 2023, 49(8): 1745−1757 doi: 10.16383/j.aas.c210147
Citation: Ma Jing, Yang Xiao-Mei, Sun Shu-Li. Distributed fusion estimation for multi-sensor systems with time-correlated multiplicative noises. Acta Automatica Sinica, 2023, 49(8): 1745−1757 doi: 10.16383/j.aas.c210147

带时间相关乘性噪声多传感器系统的分布式融合估计

doi: 10.16383/j.aas.c210147
基金项目: 国家自然科学基金(61573132, 61873087), 黑龙江省信息融合估计与检测重点实验室资助
详细信息
    作者简介:

    马静:黑龙江大学数学科学学院教授. 主要研究方向为信息融合状态估计. E-mail: majing@hlju.edu.cn

    杨晓梅:黑龙江大学数学科学学院硕士研究生. 主要研究方向为信息融合状态估计. E-mail: 15615088651@163.com

    孙书利:黑龙江大学电子工程学院教授. 主要研究方向为状态估计, 多传感器信息融合和网络化系统滤波. 本文通信作者. E-mail: sunsl@hlju.edu.cn

Distributed Fusion Estimation for Multi-sensor Systems With Time-correlated Multiplicative Noises

Funds: Supported by National Nature Science Foundation of China (61573132, 61873087) and Information Fusion Estimation and Detection Key Laboratory in Heilongjiang Province
More Information
    Author Bio:

    MA Jing Professor at the School of Mathematics Science, Heilongjiang University. Her main research interest is information fusion state estimation

    YANG Xiao-Mei Master student at the School of Mathematics Science, Heilongjiang University. Her main research interest is information fusion state estimation

    SUN Shu-Li Professor at the School of Electronic Engineering, Heilongjiang University. His research interest covers state estimation, multi-sensor information fusion, and networked systems filtering. Corresponding author of this paper

  • 摘要: 研究带时间相关乘性噪声多传感器系统的分布式融合估计问题, 其中时间相关的乘性噪声满足一阶高斯−马尔科夫过程. 通过引入虚拟状态和虚拟过程噪声, 构建了虚拟状态的递推方程. 首先, 基于新息分析方法, 分别对系统状态和虚拟状态设计局部一步预报器. 然后, 基于一步预报器设计状态的局部线性滤波器、多步预报器和平滑器. 推导了任意两个局部状态估计误差之间的互协方差矩阵. 接着, 基于线性最小方差意义下的矩阵加权、对角矩阵加权和标量加权融合算法, 给出相应的分布式融合状态估值器. 最后, 分析算法的稳定性. 仿真研究验证了该算法的有效性.
  • 卡尔曼滤波广泛应用于各种领域, 如惯性导航、制导、全球定位系统、目标跟踪、通信、信号处理、控制等[1-2]. 阿波罗登月计划和C-5A飞机导航系统的设计是卡尔曼滤波算法早期应用中最成功的实例. 经典的卡尔曼滤波算法是在精确已知系统模型的前提下给出的, 但在许多重要的物理过程中, 通过引入乘性噪声来刻画模型参数的不确定性, 使建立的数学模型更符合实际系统. 乘性噪声广泛存在于网络化系统[3-4]和通信系统[5-6]中. 在通信系统[5-6]中, 信号源的传输信道是一个衰减系数可由增强现实描述的时变衰减信道. 这种信道的时变衰减性可看作系统的乘性扰动[7-9]. 乘性噪声的存在使系统呈现了非线性. 当乘性噪声为白噪声时, 由于白噪声具有相同时刻相关且不同时刻独立的性质, 因此这种噪声的扰动不具有持久性, 只会影响到当前时刻的观测[3-4]. 但当乘性噪声为有色噪声时, 扰动会持续一段时间或一直存在, 因此滤波器的设计具有一定的挑战性[7-15].

    对带有时间相关的乘性噪声系统, 文献[7]假设加性噪声为独立的白噪声, 通过构建一个递推的新息过程, 提出了线性最小均方意义下的递推滤波器. 文献[8]将文献[7]的结果推广到非零均值乘性噪声的情形. 文献[9]进一步将文献[7]的结果推广到带有相关加性噪声的系统中, 提出了相应的递推最优滤波算法, 但文献[9]需要计算乘性噪声系数矩阵的逆. 文献[10]考虑了非高斯乘性噪声系统的滤波问题. 进一步, 当系统同时受到多个乘性噪声扰动时, 文献[11]通过引入一些新的递推项, 提出了最小均方最优递推滤波算法. 基于文献[11]提出的方法, 文献[12]考虑了网络传输中的随机数据包丢失和时间相关的乘性噪声, 提出了线性最优滤波器、预报器和平滑器. 文献[13]处理的是带时间相关通道噪声系统的最优滤波问题. 但文献[1113]考虑的加性噪声仍然为相互独立的白噪声. 文献[1415]对带时间相关的乘性噪声系统, 提出了Tobit卡尔曼滤波算法. 然而, 以上研究成果都是基于单传感器系统进行的研究, 鲜有对带时间相关乘性噪声的多传感器系统的信息融合估计问题的研究.

    目前, 多传感器信息融合状态估计有集中式融合、分布式融合和序贯融合三种基本方法. 集中式融合虽然可以得到全局最优融合估计, 但计算负担大、容错性差. 序贯融合更便于处理异步采样系统[16]. 分布式融合具有并行结构、易于故障的检测与隔离的特点, 得到了大家的广泛关注. 目前比较流行的分布式融合算法有联邦卡尔曼滤波[17]、极大似然最优融合准则[18]、线性最小方差意义下的加权融合算法[19]、协方差交叉融合算法[20-21]和文献[22-23]提出的一种类卡尔曼结构的递推分布式融合滤波算法. 文献[24]对乘性噪声为白噪声的多传感器系统, 提出了基于线性最小方差意义下的矩阵加权融合算法; 考虑网络传输中的随机滞后和丢包, 提出了分布式融合滤波器. 文献[25]考虑了有限步相关的加性噪声, 提出了分布式融合估计. 基于三种融合方法和改进的协方差交叉融合方法, 文献[26]对带有不确定噪声方差系统, 提出了鲁棒加权融合卡尔曼估值器. 但以上结果是在乘性噪声为白噪声或在有限步相关加性噪声下提出的, 而没有考虑时间相关乘性噪声. 文献[27]考虑了时间相关的乘性噪声, 在线性最小方差意义下, 提出了分布式标量加权融合滤波器, 但没有考虑预报和平滑问题.

    基于以上分析, 本文研究带时间相关乘性噪声和相关加性噪声多传感器系统的分布式融合估计问题. 主要贡献如下: 1)通过引入虚拟状态, 提出一种新的增广方法, 与文献[11]注释6基于Kronecker乘积提出的增广方法相比, 具有更小的计算负担; 2)为进一步减轻计算负担, 采用非增广方法, 在线性最小方差意义下, 提出了系统状态的局部线性最优滤波器、预报器和平滑器, 推广了现有文献中的结果; 3)推导了任意两个局部估计误差之间的互协方差矩阵, 来确定分布式融合权重; 4)对提出的分布式融合估值器进行了稳态分析, 给出稳态估计存在的一个充分条件. 当稳态估计存在时, 方差矩阵和互协方差矩阵均可任选初值离线迭代计算, 即融合权重可以离线计算. 因此分布式稳态融合估计可明显减轻融合中心的在线计算负担, 便于实时应用.

    符号描述: $L({{\boldsymbol{y}}_i}(1),{{\boldsymbol{y}}_i}(2), \cdots ,{{\boldsymbol{y}}_i}(t))$ 表示由量测序列$({{\boldsymbol{y}}_i}(1),{{\boldsymbol{y}}_i}(2), \cdots ,{{\boldsymbol{y}}_i}(t))$所生成的线性空间, 其中下标$i$表示第$i$个传感器. ${{\bf{R}}^n}$表示$n$维列向量空间, 上标T表示转置, ${\rm{tr}}$表示矩阵的迹, ${\rm{E}}[ \cdot ]$表示数学期望, $\rho ( \cdot )$ 表示矩阵“$\cdot$”的谱半径, ${\rm{diag}}\{{ \cdot}\}$表示对角块为“$\cdot$”的块对角矩阵. ${\boldsymbol{x}} \bot {\boldsymbol{y}}$表示随机变量${\boldsymbol{x}}$与${\boldsymbol{y}}$正交, 即${\rm{E}}[{\boldsymbol{x}}{{\boldsymbol{y}}^{\rm{T}}}] = 0$, ${\delta _{tk}}$表示Kronecker delta函数. ${\hat {\boldsymbol{x}}_i}( \cdot | \circ )$表示基于线性空间$L({{\boldsymbol{y}}_i}(1), {{\boldsymbol{y}}_i}(2), \cdots , {{\boldsymbol{y}}_i}( \circ ))$对随机变量${\boldsymbol{x}}( \cdot )$的估计, 其中“$\cdot$” 和 “$\circ$” 表示时刻. ${\tilde {\boldsymbol{x}}_i}( \cdot| \circ) = {{\boldsymbol{x}}_i}( \cdot ) - {\hat {\boldsymbol{x}}_i}( \cdot | \circ )$ 表示估计误差, ${\boldsymbol{P}}_{ij}^{xy}( \cdot | \circ ) = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}( \cdot | \circ )\tilde {\boldsymbol{y}}_j^{\rm{T}}( \cdot | \circ )]$表示随机变量${\tilde {\boldsymbol{x}}_i}( \cdot | \circ )$和${\tilde {\boldsymbol{y}}_j}( \cdot | \circ )$之间的互协方差矩阵, 其中${\boldsymbol{P}}_{ij}^{xx}( \cdot | \circ ) = {\boldsymbol{P}}_{ij}^x( \cdot | \circ ),$ ${\boldsymbol{P}}_{ii}^{xx}( \cdot | \circ ) = {\boldsymbol{P}}_i^x( \cdot | \circ ),$ ${\boldsymbol{P}}_{ii}^{xy}( \cdot | \circ ) = {\boldsymbol{P}}_i^{xy}( \cdot | \circ) ,$ $1 \leq i,j \leq L.$

    考虑带时间相关乘性噪声的多传感器线性离散随机系统:

    $${\boldsymbol{x}}(t + 1) = {\boldsymbol{Ax}}(t) + {{{\boldsymbol{\Gamma}} w}}(t)$$ (1)
    $${{\boldsymbol{y}}_i}(t) = ({{\boldsymbol{C}}_i} + {\tilde {\boldsymbol{C}}_i}{r_i}(t)){\boldsymbol{x}}(t) + {{\boldsymbol{v}}_i}(t)$$ (2)
    $${r_i}(t + 1) = {h_i}{r_i}(t) + {p_i}(t),\;i = 1, \cdots ,L$$ (3)

    式中, ${\boldsymbol{x}}(t) \in {{\bf{R}}^n}$是系统状态, ${{\boldsymbol{y}}_i}(t) \in {{\bf{R}}^{{m_i}}}$是第$i$个传感器的量测. ${\boldsymbol{w}}(t)$、${{\boldsymbol{v}}_i}(t)$和 ${p_i}(t) \in {\bf{R}}$, $i = 1, \cdots , L$, 为零均值白噪声, ${r_i}(t) \in {\bf{R}}$为满足一阶高斯−马尔科夫过程的乘性噪声. ${\boldsymbol{A}}$、${\boldsymbol{\Gamma }}$、${{\boldsymbol{C}}_{{i}}}$、${\tilde {\boldsymbol{C}}_{{i}}}$、${h_i}$为已知的适当维数的常值矩阵.

    假设1. 加性噪声${\boldsymbol{w}}(t)$、${{\boldsymbol{v}}_i}(t)$和${p_i}(t)$具有如下统计信息:

    $$ \begin{split} &{\rm{E}}[{\boldsymbol{w}}(t)] = 0,\;{\rm{E}}[{{\boldsymbol{v}}_i}(t)] = 0,\;{\rm{E}}[{p_i}(t)] = 0,\;\; i=1,\cdots,L\\ &{\rm{E}}\left[ {\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{w}}(t)} \\ {{{\boldsymbol{v}}_i}(t)} \\ {{p_i}(t)} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{w}}^{\rm{T}}}(k)}&{{\boldsymbol{v}}_j^{\rm{T}}(k)}&{{p_j}(k)} \end{array}} \right]} \right] =\\ &\qquad\left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{Q}}^w}}&{{{\boldsymbol{R}}_j}}&0 \\ {{\boldsymbol{R}}_i^{\rm{T}}}&{{\boldsymbol{Q}}_{ij}^v}&0 \\ 0&0&{Q_{ij}^p{\delta _{ij}}} \end{array}} \right]{\delta _{tk}},\;\;\;i,j = 1, \cdots ,L \end{split} $$

    式中, 当$i = j$时, 有${\boldsymbol{Q}}_{ii}^v = {\boldsymbol{Q}}_i^v$, $Q_{ii}^p = Q_i^p$.

    假设2. ${\boldsymbol{x}}(t)$的初值${\boldsymbol{x}}(0)$的均值为 ${\rm{E}}[{\boldsymbol{x}}(0)] = {{\boldsymbol{\mu }}_0}$, 协方差矩阵为${\rm{E}}[({\boldsymbol{x}}(0) - {{\boldsymbol{\mu }}_0}) {({\boldsymbol{x}}(0) - {{\boldsymbol{\mu }}_0})^{\rm{T}}}] = {P_0}$, 且与加性噪声 ${\boldsymbol{w}}(t)$、${{\boldsymbol{v}}_i}(t)$ 和 ${p_i}(t)$, $t = 0,1, 2, \cdots$不相关.

    假设3. ${r_i}(t)$的初值${r_i}({{0}})$的均值为${\rm{E}}[{r_i}(0)] = 0$, 方差为${\rm{E}}[r_i^2(0)] = Q_i^{{r_0}}$, 且与 ${\boldsymbol{w}}(t)$, ${{\boldsymbol{v}}_i}(t)$, ${p_i}(t)$, $t = 0,1,2, \cdots$和${\boldsymbol{x}}(0)$不相关. 且${r_i}({{0}})$与${r_j}({{0}})$, $i \ne j$ 也不相关.

    本文目的是基于量测序列$({{\boldsymbol{y}}_i}(1),{{\boldsymbol{y}}_i}(2), \cdots , {{\boldsymbol{y}}_i}(t))$, $i = 1, \cdots ,L$, 求线性最小方差意义下的分布式融合估值器${\hat {\boldsymbol{x}}_o}(t|t + N)$, 包括滤波器$(N = 0)$、预报器$(N < 0)$和平滑器$(N > 0)$.

    注 1. 在实际应用中, 当系统和传感器在相同的背景噪声下或受相同噪声源影响时, 会导致相关的加性噪声; 当广义系统转化为正常系统[28]且时滞系统转化为无时滞系统时[29], 也会产生相关的加性噪声. 时间相关乘性噪声可由传输信道的不理想引起, 例如文献[5-6]将时变衰减信道的信道系数描述为一个增强现实过程; 信号在采样、选通和幅值调制的过程中, 可能产生时间相关的乘性噪声[7-8]; 利用脉冲雷达对移动车辆的位置和方位角进行观测时, 通常会受到时间相关乘性噪声的影响[11]. 文献[10]将时间相关的乘性噪声理解为传感器系数调整不够准确或作为传感器的一种损耗.

    已知获得状态估值器的关键在于获取观测的一步预报误差, 即新息${{\boldsymbol{\epsilon }}_i}(t) = {{\boldsymbol{y}}_i}(t) - {\hat {\boldsymbol{y}}_i}(t|t - 1)$. 由观测方程(2)和乘性噪声模型(3)可以看出, 观测的一步预报器${\hat {\boldsymbol{y}}_i}(t|t - 1)$中包含了${r_i}(t){\boldsymbol{x}}(t)$的一步预报器. 因此获得${r_i}(t){\boldsymbol{x}}(t)$的一步预报器是获得状态估值器的关键. 为此, 引入虚拟状态向量${\boldsymbol{x}}_i^r(t) = {r_i}(t){\boldsymbol{x}}(t)$. 那么系统(1) ~ (3)可重新写为如下新系统:

    $${\boldsymbol{x}}(t + 1) = {{{\boldsymbol{A}}}}{\boldsymbol{x}}(t) + {\boldsymbol{\Gamma w}}(t)$$ (4)
    $${\boldsymbol{x}}_i^r(t + 1) = {{\boldsymbol{\Phi }}_i}{\boldsymbol{x}}_i^r(t) + {{\boldsymbol{\xi }}_i}(t)$$ (5)
    $${{\boldsymbol{y}}_i}(t) = {{\boldsymbol{C}}_i}{\boldsymbol{x}}(t) + {\tilde {\boldsymbol{C}}_i}{\boldsymbol{x}}_i^r(t) + {{\boldsymbol{v}}_i}(t),\;i = 1, \cdots ,L$$ (6)

    式中, ${{\boldsymbol{\Phi }}_i} = {h_i}{\boldsymbol{A}}$, 虚拟过程噪声${{\boldsymbol{\xi }}_i}(t)$计算如下:

    $${{\boldsymbol{\xi }}_i}(t) = {h_i}{\boldsymbol{\Gamma }}{r_i}(t){\boldsymbol{w}}(t) + {\boldsymbol{A}}{p_i}(t){\boldsymbol{x}}(t) + {\boldsymbol{\Gamma }}{p_i}(t){\boldsymbol{w}}(t)$$ (7)

    为便于后续引理和定理的推导, 先给出如下正交关系.

    注 2. 由式(1)和式(3)可得$L({\boldsymbol{x}}(1), \cdots ,{\boldsymbol{x}}(t)) \subset L({\boldsymbol{x}}(0),{\boldsymbol{w}}(0), \cdots ,{\boldsymbol{w}}(t - 1))$ 和 $L({r_i}(1), \cdots ,{r_i}(t)) \subset L({r_i}(0),{p_i}(0),\cdots ,{p_i}(t - 1))$. 因此, 利用假设1 ~ 3有${\boldsymbol{w}}(t + k) \bot {\boldsymbol{x}}(t)、$${{\boldsymbol{v}}_i}(t + k) \bot {\boldsymbol{x}}(t)、$${p_i}(t \;+\; k) \bot {r_i}(t),$ $k \geq 0$、${p_i}(t) \bot {\boldsymbol{x}}(k)$、${r_i}(t) \bot {\boldsymbol{x}}(k)$、${r_i}(t) \bot {\boldsymbol{w}}(k)$ 和${r_i}(t) \bot {{\boldsymbol{v}}_i}(k)$, $\forall t,k$.

    注3. 对式(3)迭代得${r_i}(t) = {({h_i})^t}{r_i}(0)\; +$$\sum\nolimits_{k = 0}^{t - 1}\cdot {{{({h_i})}^{t - k - 1}}} {p_i}(k)$, 由假设1和假设3可知, ${\rm{E}}[{r_i}(t)] = {({h_i})^t}{\rm{E[}}{r_i}{\rm{(0)]}} \;+\;\sum\nolimits_{k = 0}^{t - 1} {{{({h_i})}^{t - k - 1}}} {\rm{E[}}{p_i}(k)] \;=\; 0$. 因此, 有${\rm{E}}[{\boldsymbol{x}}_i^r(t){\boldsymbol{v}}_j^{\rm{T}}(t + k)] = {\rm{E}}[{r_i}(t)]$${\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{v}}_j^{\rm{T}}(t + k)] = 0$, ${\rm{E}}[{\boldsymbol{x}}_i^r(t){\boldsymbol{w}}_{}^{\rm{T}}(t + k)] \;=\; {\rm{E}}[{r_i}(t)]$${\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{w}}_{}^{\rm{T}}(t + k)] \;= \; 0$, $\forall k$; ${\rm{E}}[{\boldsymbol{x}}_i^r(t)p_j^{\rm{T}}(t + k)] = $${\rm{E}}[{r_i}(t)p_j^{\rm{T}}(t + k)]{\rm{E}}[{\boldsymbol{x}}(t)] = 0$, $k \geq 0$, 即${\boldsymbol{x}}_i^r(t) \bot {\boldsymbol{v}}_j^{}(t + k)$, ${\boldsymbol{x}}_i^r(t) \bot {\boldsymbol{w}}(t + k)$, $\forall k$和${\boldsymbol{x}}_i^r(t) \bot p_j^{}(t + k)$, $k \geq 0$. 由式(4)和式(6)有$L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t)) \subset L({{\boldsymbol{v}}_i}(1), \cdots ,{{\boldsymbol{v}}_i}(t),{\boldsymbol{x}}(0),{\boldsymbol{w}}(0), \cdots ,$ ${\boldsymbol{w}}(t - 1), {\boldsymbol{x}}_i^r(1), \cdots ,{\boldsymbol{x}}_i^r(t))$. 从而有 ${\boldsymbol{w}}\left(t \;+\; k\right)\; \bot $$L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t)),$ ${{\boldsymbol{v}}_i}(t + k) \bot L({{\boldsymbol{y}}_i}(1), \cdots ,{{\boldsymbol{y}}_i}(t))$, $k \geq 1$, ${p_i}(t) \bot L({{\boldsymbol{y}}_i}(1), \cdots ,{{\boldsymbol{y}}_i}(t))$.

    引理 1. 虚拟过程噪声${{\boldsymbol{\xi }}_i}(t)$具有如下统计信息:

    $$\left\{\begin{aligned} &{\rm{E}}[{{\boldsymbol{\xi }}_i}(t)] = 0, \;\;\\ &{\rm{E}}[{{\boldsymbol{\xi }}_i}(t){{\boldsymbol{w}}^{\rm{T}}}(k){{]}} = 0,\\ &{\rm{E}}[{{\boldsymbol{\xi }}_i}(t){\boldsymbol{v}}_j^{\rm{T}}(k){{]}} = 0, \end{aligned}\right.{\;\;\forall t,k,\;\;\;i,j = 1, \cdots ,L}$$ (8)
    $$\begin{split} {\boldsymbol{Q}}_{ij}^{\boldsymbol{\xi }}(t) = \;&{\rm{E[}}{{\boldsymbol{\xi }}_i}{\rm{(}}t{\rm{)}}{\boldsymbol{\xi }}_j^{\rm{T}}{\rm{(}}k{\rm{)]}} = {\rm{(}}h_i^2{\boldsymbol{Q}}_i^r(t){\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}} \;+ \\ &{{\boldsymbol{Q}}_i^p{\boldsymbol{A\Xi }}(t){{\boldsymbol{A}}^{\rm{T}}} + {\boldsymbol{Q}}_i^p{\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}){\delta _{ij}}{\delta _{tk}}} \end{split}$$ (9)

    式中, ${\boldsymbol{Q}}_{ii}^{\boldsymbol{\xi}} (t) = {\boldsymbol{Q}}_i^{\boldsymbol{\xi}} (t)$. ${\boldsymbol{\Xi }}(t) = {\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{x}}^{\rm{T}}}(t)]$ 和$Q_i^r(t) =$ ${\rm{E}}[r_i^2(t)]$分别为${\boldsymbol{x}}(t)$和${r_i}(t)$的二阶矩阵, 可递推计算为:

    $${\boldsymbol{\Xi }}(t + 1) = {\boldsymbol{A\Xi }}(t){{\boldsymbol{A}}^{\rm{T}}} + {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}$$ (10)
    $$Q_i^r(t + 1) = h_i^2Q_i^r(t) + Q_i^p$$ (11)

    式中, 初值为${\boldsymbol{\Xi }}(0) = {{\boldsymbol{P}}_0} + {{\boldsymbol{\mu }}_0}{\boldsymbol{\mu }}_0^{\rm{T}}$, $Q_i^r(0) = Q_i^{{r_0}}$, $i = 1, \cdots ,L$.

    证明. 见附录A.

    注4. 由引理1可知, ${{\boldsymbol{\xi }}_i}(t)$为独立于${\boldsymbol{w}}(k)$和${{\boldsymbol{v}}_j}(k), \forall k$的零均值白噪声. 由式(4) ~ (6)有$L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t)) \subset L({{\boldsymbol{v}}_i}(1), \cdots ,{{\boldsymbol{v}}_i}(t),{\boldsymbol{x}}_i^r(0),{{\boldsymbol{\xi }}_i}(0),$ $\cdots ,{{\boldsymbol{\xi }}_i}(t - 1), {\boldsymbol{x}}(0),{\boldsymbol{w}}(0), \cdots ,{\boldsymbol{w}}(t - 1)) .$ 从而有 ${{\boldsymbol{\xi }}_i}(t + k) \bot L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t)) ,$ $k \geq 0,$ ${{\boldsymbol{\xi }}_j}(t + k) \bot $$L({{\boldsymbol{y}}_i}(1), \cdots ,{{\boldsymbol{y}}_i}(t)),$ $\forall k,$ $j \ne i$ .

    注 5. 对系统(4) ~ (6), 可以引入增广状态向量${\bar {\boldsymbol{x}}_i}(t) = {[ {{{\boldsymbol{x}}^{\rm{T}}}(t)}\;\;{({\boldsymbol{x}}_i^r} (t))^{\rm{T}}}{]^{\rm{T}}}$和增广过程噪声向量${\bar {\boldsymbol{w}}_i}(t) = {[ {{{\boldsymbol{w}}^{\rm{T}}}(t)}\;\;{{\boldsymbol{\xi }}_i^{\rm{T}}(t)} ]^{\rm{T}}}$, 则有如下增广系统:

    $${\bar {\boldsymbol{x}}_i}(t + 1) = {\bar {\boldsymbol{\Phi }}_i}{\bar {\boldsymbol{x}}_i}(t) + \bar {\boldsymbol{\Gamma }}{\bar {\boldsymbol{w}}_i}(t)$$ (12)
    $${{\boldsymbol{y}}_i}(t) = {\bar {\boldsymbol{C}}_i}{\bar {\boldsymbol{x}}_i}(t) + {{\boldsymbol{v}}_i}(t),\;\;i = 1, \cdots ,L$$ (13)

    其中

    $$\left\{\begin{aligned} &{\bar {\boldsymbol{\Phi }}_i} = \left[ {\begin{array}{*{20}{c}} {\boldsymbol{A}}&0 \\ 0&{{{\boldsymbol{\Phi }}_i}} \end{array}} \right] \;\;\\ &\bar {\boldsymbol{\Gamma }} = \left[ {\begin{array}{*{20}{c}} {\boldsymbol{\Gamma }}&0 \\ 0&{{{\boldsymbol{I}}_n}} \end{array}} \right] \\ &{\bar {\boldsymbol{C}}_i} = [\begin{array}{*{20}{c}} {{{\boldsymbol{C}}_i}}&{{{\tilde {\boldsymbol{C}}}_i}} \end{array}] \end{aligned}\right.$$ (14)

    根据引理1可知, ${\bar {\boldsymbol{w}}_i}(t)$和${{\boldsymbol{v}}_i}(t)$为零均值的相关白噪声, 且满足如下噪声统计信息:

    $$\left\{\begin{aligned} &{\rm{E}}[{\bar {\boldsymbol{w}}_i}(t)] = 0 \;\;\\ &{\rm{E}}[{{\boldsymbol{v}}_i}(t)] = 0\\ &{{{\boldsymbol{Q}}_{ij}}(t) = {\rm{E}}[{\bar {\boldsymbol{w}}_i}(t)\bar {\boldsymbol{w}}_j^{\rm{T}}(k)] = }\\ &\qquad{\rm{E}}\left\{ {\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{w}}(t)} \\ {{{\boldsymbol{\xi }}_i}(t)} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{w}}^{\rm{T}}}(k)}&{{\boldsymbol{\xi }}_j^{\rm{T}}(k)} \end{array}} \right]} \right\} =\\ &\qquad\left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{Q}}^w}}&0 \\ 0&{{\boldsymbol{Q}}_{ij}^\xi (t)} \end{array}} \right]{\delta _{tk}}\\ & {{\boldsymbol{S}}_j}{\rm{ = E}}[{\bar {\boldsymbol{w}}_i}(t){\boldsymbol{v}}_j^{\rm{T}}(k)] = \\ &\qquad{\rm{E}}\left\{ {\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{w}}(t)} \\ {{{\boldsymbol{\xi }}_i}(t)} \end{array}} \right]{\boldsymbol{v}}_j^{\rm{T}}(k)} \right\} = \left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{R}}_j}} \\ 0 \end{array}} \right]{\delta _{tk}} \end{aligned}\right.$$ (15)

    式中, $ {{\boldsymbol{Q}}_{ii}}(t) = {{\boldsymbol{Q}}_i}(t)$.

    这样, 原系统(1) ~ (3)的分布式融合估计问题已经被转化为带相关白噪声增广系统(12)和(13)的分布式融合估计问题. 而增广系统(12)和(13)为标准卡尔曼滤波算法所对应的系统, 其分布式融合估计算法可参见文献[19]. 本文忽略计算负担小的加/减法的运算, 通过计算算法乘/除法的次数, 给出增广方法的计算量级. 因为${\bar {\boldsymbol{x}}_i}(t) \in {{\bf{R}}^{2n}}$, 故增广方法的计算量级为${\rm{O}}{(2n)^3} = {\rm{O}}(8{n^3}),$与文献[11]注释6基于Kronecker乘积提出的增广方法${\rm{O}}({n^2} + 3n + 2)^3$相比, 具有更小的计算负担.

    为了进一步减轻计算负担, 下面将采用直接对式(4) ~ (6)的系统状态和虚拟状态设计一步预报器的思想, 提出如下非增广估计方法.

    引理2. 对于系统(4) ~ (6), 在假设1 ~ 3条件下, 新息${{\boldsymbol{\epsilon }}_i}(t)$和新息方差${\boldsymbol{Q}}_i^\epsilon (t)$计算如下:

    $${{\boldsymbol{\epsilon }}_i}(t) = {{\boldsymbol{y}}_i}(t) - {{\boldsymbol{C}}_i}{\hat {\boldsymbol{x}}_i}(t|t - 1) - {\tilde {\boldsymbol{C}}_i}\hat {\boldsymbol{x}}_i^r(t|t - 1)$$ (16)
    $$\begin{split} {\boldsymbol{Q}}_i^\epsilon (t) =\;&{{\boldsymbol{C}}_i}{\boldsymbol{P}}_i^x(t|t - 1){\boldsymbol{C}}_i^{\rm{T}}{\rm{ + }}{\tilde {\boldsymbol{C}}_i}{\boldsymbol{P}}_i^{{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_i^{\rm{T}}\; + \\ &{{\boldsymbol{C}}_i}{\boldsymbol{P}}_i^{x{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_i^{\rm{T}} \;+ \\ &{\tilde {\boldsymbol{C}}_i}{({\boldsymbol{P}}_i^{x{x^r}}(t|t - 1))^{\rm{T}}}{\boldsymbol{C}}_i^{\rm{T}} + {\boldsymbol{Q}}_i^v\\[-10pt] \end{split}$$ (17)

    证明. 见附录B.

    由引理2可知, 新息${{\boldsymbol{\epsilon }}_i}(t)$依靠虚拟状态${\boldsymbol{x}}_i^r(t)$的一步预报器$\hat {\boldsymbol{x}}_i^r(t|t - 1)$, ${\boldsymbol{Q}}_i^\epsilon (t)$依靠系统状态和虚拟状态的一步预报误差的方差矩阵${\boldsymbol{P}}_i^x(t|t - 1)$和${\boldsymbol{P}}_i^{{x^r}}(t|t - 1)$以及协方差矩阵${\boldsymbol{P}}_i^{x{x^r}}(t|t - 1)$. 下面给出系统状态${\boldsymbol{x}}(t)$和虚拟状态${\boldsymbol{x}}_i^r(t)$的一步预报器和相应的方差矩阵以及协方差矩阵.

    定理 1. 在假设1 ~ 3条件下, 系统状态${\boldsymbol{x}}(t)$和虚拟状态${\boldsymbol{x}}_i^r(t)$的局部一步预报器为:

    $${\hat {\boldsymbol{x}}_i}(t + 1|t) = {\boldsymbol{A}}{\hat {\boldsymbol{x}}_i}(t|t - 1) + {{\boldsymbol{L}}_i}(t){{\boldsymbol{\epsilon }}_i}(t)$$ (18)
    $$\hat {\boldsymbol{x}}_i^r(t + 1|t) = {{\boldsymbol{\Phi }}_i}\hat {\boldsymbol{x}}_i^r(t|t - 1){\rm{ + }}{\boldsymbol{L}}_i^{{x^r}}(t){{\boldsymbol{\epsilon }}_i}(t)$$ (19)

    式中, 增益矩阵为:

    $$\begin{split} {{\boldsymbol{L}}_i}(t) =\;& [{\boldsymbol{A}}({\boldsymbol{P}}_i^x(t|t - 1){\boldsymbol{C}}_i^{\rm{T}}\; + \\ &{\boldsymbol{P}}_i^{x{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_i^{\rm{T}}) + {\boldsymbol{\Gamma }}{{\bf{R}}_i}]{({\boldsymbol{Q}}_i^\epsilon (t))^{ - 1}}\end{split}\;\;$$ (20)
    $$\begin{split} {\boldsymbol{L}}_i^{{x^r}}(t) =\;& [{{\boldsymbol{\Phi }}_i}{({\boldsymbol{P}}_i^{x{x^r}}(t|t - 1))^{\rm{T}}}{\boldsymbol{C}}_i^{\rm{T}}\; + \\ &{\boldsymbol{P}}_i^{{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_i^{\rm{T}}]{({\boldsymbol{Q}}_i^\epsilon (t))^{ - 1}} \end{split}\qquad\qquad$$ (21)

    一步预报误差方差矩阵为:

    $$\begin{split} {\boldsymbol{P}}_i^x(t + 1|t) =\;& {\boldsymbol{AP}}_i^x(t|t - 1){{\boldsymbol{A}}^{\rm{T}}}\; - \\ &{{\boldsymbol{L}}_i}(t){\boldsymbol{Q}}_i^\epsilon (t){\boldsymbol{L}}_i^{\rm{T}}(t) + {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}} \end{split}\;\;\;$$ (22)
    $$\begin{split} {\boldsymbol{P}}_i^{{x^r}}(t + 1|t) =\;& {{\boldsymbol{\Phi }}_i}{\boldsymbol{P}}_i^{{x^r}}(t|t - 1){\boldsymbol{\Phi }}_i^{\rm{T}} \;- \\ &{\boldsymbol{L}}_i^{{x^r}}(t){\boldsymbol{Q}}_i^\epsilon (t){({\boldsymbol{L}}_i^{{x^r}}(t))^{\rm{T}}} + {\boldsymbol{Q}}_i^\xi (t) \end{split}$$ (23)

    ${\tilde {\boldsymbol{x}}_i}(t + 1|t)$和$\tilde {\boldsymbol{x}}_i^r(t + 1|t)$的协方差矩阵为:

    $$\begin{split} {\boldsymbol{P}}_i^{x{x^r}}(t + 1|t) =\;& {\boldsymbol{AP}}_i^{x{x^r}}(t|t - 1){\boldsymbol{\Phi }}_i^{\rm{T}}\; -\\ &{{\boldsymbol{L}}_i}(t){\boldsymbol{Q}}_i^\epsilon (t){({\boldsymbol{L}}_i^{{x^r}}(t))^{\rm{T}}} \end{split}$$ (24)

    式中, 初值为${\hat {\boldsymbol{x}}_i}(0| - 1) = {{\boldsymbol{\mu }}_0}$, $\hat {\boldsymbol{x}}_i^r(0| - 1) = 0$, ${{\boldsymbol{P}}}_{i}^{x}(0|\;- 1)= {{\boldsymbol{P}}}_{0}, {{\boldsymbol{P}}}_{i}^{{x}^{r}}({0}|-1)={{\boldsymbol{P}}}_{0}{Q}_{i}^{{r}_{0}}$和${\boldsymbol{P}}_i^{x{x^r}}(0| - 1) = 0$.

    证明. 见附录C.

    注6. 由定理1可以看出, 与以往算法不同, 新息计算中同时涉及了系统状态和虚拟状态的一步预报器, 这给后续方差矩阵和互协方差矩阵, 尤其是多步平滑器的方差矩阵的推导带来了难度(见定理6和证明).

    下面将基于一步预报器${\hat {\boldsymbol{x}}_i}(t + {\rm{1}}|t)$设计多步预报器${\hat {\boldsymbol{x}}_i}(t|t + N)\;(N < - {\rm{1}})$、滤波器$(N = 0)$和多步平滑器$(N > 0)$.

    定理2. 在假设1 ~ 3条件下, 系统状态${\boldsymbol{x}}(t)$的局部, $N(N < - 1)$步预报器为:

    $${\hat {\boldsymbol{x}}_i}(t{\rm{|}}t + N) = {{\boldsymbol{A}}^{ - N - 1}}{\hat {\boldsymbol{x}}_i}(t + N + 1|t + N)$$ (25)

    $N$步预报误差方差矩阵为:

    $$\begin{split} {\boldsymbol{P}}_i^x(t|t{\rm{ + }}N) =\;& {{\boldsymbol{A}}^{ - N - 1}}{\boldsymbol{P}}_i^x(t + N + 1|t + N)\;\times\\ &{({{\boldsymbol{A}}^{ - N - 1}})^{\rm{T}}} +\sum\limits_{k = 1}^{ - N - 1} {{\boldsymbol{A}}^{ - N - k - 1}}\;\times\\ &{\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}{{({{\boldsymbol{A}}^{ - N - k - 1}})}^{\rm{T}}} \end{split}$$ (26)

    式中, ${\hat {\boldsymbol{x}}_i}(t + N + 1|t + N)$ 和 ${\boldsymbol{P}}_i^x(t + N + 1|t + N)$由定理1给出.

    证明. 见附录D.

    定理 3. 在假设1 ~ 3的条件下, 系统状态${\boldsymbol{x}}(t)$的局部滤波器和平滑器${\hat {\boldsymbol{x}}_i}(t|t + N)$, $N \geq 0$为:

    $$\begin{split} {\hat {\boldsymbol{x}}_i}(t|t + N) =\;& {\hat {\boldsymbol{x}}_i}(t|t + N - 1) \;+ \\ &{{\boldsymbol{K}}_i}(t|t + N){{\boldsymbol{\epsilon }}_i}(t + N) \end{split}$$ (27)

    式中, 增益矩阵为:

    $$\begin{split} &{{\boldsymbol{K}}_i}(t|t + N) = [{\boldsymbol{\Omega }}_i^x(t + N|t + N - 1){\boldsymbol{C}}_i^{\rm{T}} \;+ \\ &\;\;\;\;\;{\boldsymbol{\Omega }}_i^{x{x^r}}(t + N|t + N - 1)\tilde {\boldsymbol{C}}_i^{\rm{T}}]{({\boldsymbol{Q}}_i^\epsilon (t + N))^{ - 1}} \end{split}$$ (28)
    $$\begin{split} &{\boldsymbol{\Omega }}_i^x(t + N + 1|t + N) = {\boldsymbol{\Omega }}_i^x(t + N|t + N - 1){{\boldsymbol{A}}^{\rm{T}}}\; -\\ &\;\;\;\;\;\;\;{{\boldsymbol{K}}_i}(t|t + N){\boldsymbol{Q}}_i^\epsilon (t + N){\boldsymbol{L}}_i^{\rm{T}}(t + N)\\[-10pt] \end{split}$$ (29)
    $$\begin{split} &{\boldsymbol{\Omega }}_i^{x{x^r}} (t + N + 1|t + N) = {\boldsymbol{\Omega }}_i^{x{x^r}}(t + N|t + N - 1){\boldsymbol{\Phi }}_i^{\rm{T}} \;- \\ &\;\;\;\;\;\;\;\;{{\boldsymbol{K}}_i}(t|t + N){\boldsymbol{Q}}_i^\epsilon (t + N){({\boldsymbol{L}}_i^{{x^r}}(t + N))^{\rm{T}}}\\[-10pt] \end{split}$$ (30)

    误差方差矩阵为:

    $$\begin{split} {\boldsymbol{P}}_i^x(t|t + N) =\;& {\boldsymbol{P}}_i^x(t|t + N - 1) \;-\\ &{{\boldsymbol{K}}_i}(t|t + N){\boldsymbol{Q}}_i^\epsilon (t + N){\boldsymbol{K}}_i^{\rm{T}}(t|t + N) \end{split}$$ (31)

    式中, ${\boldsymbol{\Omega }}_i^x(t + N|t + N - 1) = {\rm{E}}[{\boldsymbol{x}}(t)\tilde {\boldsymbol{x}}_i^{\rm{T}}(t+N|t+ N - 1)]$ 的初值为 ${\boldsymbol{\Omega }}_i^x(t|t - 1) = {\boldsymbol{P}}_i^x(t|t \;-\; 1),$ ${\boldsymbol{\Omega }}_i^{x{x^r}}(t\; + N|t + N - 1) = {\rm{E}}[{\boldsymbol{x}}(t){(\tilde {\boldsymbol{x}}_i^r(t+N|t+N - 1))^{\rm{T}}}]$的初值为${\boldsymbol{\Omega }}_i^{x{x^r}} (t|t - 1) = {\boldsymbol{P}}_i^{x{x^r}}(t|t - 1)$由定理1给出.

    证明. 见附录E.

    注7. 文献[11, 27]滤波算法需要计算${r_i}(t){\boldsymbol{x}}(t)$的滤波器和一步预报器, 而本文算法只需要计算虚拟状态${\boldsymbol{x}}_i^r(t) = {r_i}(t){\boldsymbol{x}}(t)$的一步预报器, 降低了算法的复杂性. 且文献[11, 27]仅设计了滤波器, 没有考虑多步预报、平滑和相应的融合估计问题.

    为了获得融合权重, 下面推导任意两个局部估计误差之间的互协方差矩阵.

    定理 4. 在假设1 ~ 3的条件下, 第$i$和第$j$个传感器子系统之间的一步预报误差互协方差矩阵计算如下:

    $$\begin{split} {\boldsymbol{P}}_{ij}^x(t + 1|t) =\;& {\boldsymbol{AP}}_{ij}^x(t|t - 1){{\boldsymbol{A}}^{\rm{T}}} - {\boldsymbol{AP}}_{ij}^{x\epsilon }(t){\boldsymbol{L}}_j^{\rm{T}}(t) \;- \\ &{{\boldsymbol{L}}_i}(t){({\boldsymbol{P}}_{ji}^{x\epsilon }(t))^{\rm{T}}}{{\boldsymbol{A}}^{\rm{T}}} + {{\boldsymbol{L}}_i}(t){\boldsymbol{Q}}_{ij}^\epsilon (t){\boldsymbol{L}}_j^{\rm{T}}(t) \;- \\ &{{\boldsymbol{L}}_i}(t){\boldsymbol{R}}_i^{\rm{T}}{{\boldsymbol{\Gamma }}^{\rm{T}}} - {\boldsymbol{\Gamma }}{{\boldsymbol{R}}_j}{\boldsymbol{L}}_j^{\rm{T}}(t) + {\boldsymbol{\Gamma Q}}_{}^w{{\boldsymbol{\Gamma }}^{\rm{T}}} \end{split}$$ (32)
    $${\boldsymbol{P}}_{ij}^{x\epsilon }(t) = {\boldsymbol{P}}_{ij}^x(t|t - 1){\boldsymbol{C}}_j^{\rm{T}} + {\boldsymbol{P}}_{ij}^{x{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_j^{\rm{T}}$$ (33)

    新息互协方差矩阵计算如下:

    $$\begin{split} {\boldsymbol{Q}}_{ij}^\epsilon (t) =\;& {{\boldsymbol{C}}_i}{\boldsymbol{P}}_{ij}^x(t|t - 1){\boldsymbol{C}}_j^{\rm{T}} + {{\boldsymbol{C}}_i}{\boldsymbol{P}}_{ij}^{x{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_j^{\rm{T}} \;+ \\ &{\tilde {\boldsymbol{C}}_i}{({\boldsymbol{P}}_{ji}^{x{x^r}}(t|t - 1))^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}} \;+ \\ &{\tilde {\boldsymbol{C}}_i}{\boldsymbol{P}}_{ij}^{{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_j^{\rm{T}} + {\boldsymbol{Q}}_{ij}^v\\[-10pt] \end{split}$$ (34)

    $\tilde {\boldsymbol{x}}_i^r(t + 1|t)$与$\tilde {\boldsymbol{x}}_j^r(t + 1|t)$的互协方差矩阵为:

    $$\begin{split} {\boldsymbol{P}}_{ij}^{{x^r}}(t + 1|t) =\;& {{\boldsymbol{\Phi }}_i}{\boldsymbol{P}}_{ij}^{{x^r}}(t|t - 1){\boldsymbol{\Phi }}_j^{\rm{T}} \;- \\ &{{\boldsymbol{\Phi }}_i}{\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t){({\boldsymbol{L}}_j^{{x^r}}(t))^{\rm{T}}} \;-\\ &{\boldsymbol{L}}_i^{{x^r}}(t){({\boldsymbol{P}}_{ji}^{{x^r}\epsilon }(t))^{\rm{T}}}{\boldsymbol{\Phi }}_j^{\rm{T}}\; + \\ &{\boldsymbol{L}}_i^{{x^r}}(t){\boldsymbol{Q}}_{ij}^\epsilon (t){({\boldsymbol{L}}_j^{{x^r}}(t))^{\rm{T}}} \end{split}$$ (35)
    $${\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t) = {({\boldsymbol{P}}_{ji}^{x{x^r}}(t|t - 1))^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}} + {\boldsymbol{P}}_{ij}^{{x^r}}(t|t - 1)\tilde {\boldsymbol{C}}_j^{\rm{T}}$$ (36)

    ${\tilde {\boldsymbol{x}}_i}(t + 1|t)$与$\tilde {\boldsymbol{x}}_j^r(t + 1|t)$的互协方差矩阵为:

    $$\begin{split} &{\boldsymbol{P}}_{ij}^{x{x^r}}(t + 1|t) = {\boldsymbol{AP}}_{ij}^{x{x^r}}(t|t - 1){\boldsymbol{\Phi }}_j^{\rm{T}} \;- \\ &\;\;\;\;\;\;\;\;\;{\boldsymbol{AP}}_{ij}^{x\epsilon }(t){({\boldsymbol{L}}_j^{{x^r}}(t))^{\rm{T}}} - {{\boldsymbol{L}}_i}(t){({\boldsymbol{P}}_{ji}^{{x^r}\epsilon }(t))^{\rm{T}}}{\boldsymbol{\Phi }}_j^{\rm{T}}\; +\\ &\;\;\;\;\;\;\;\;\;{{\boldsymbol{L}}_i}(t){\boldsymbol{Q}}_{ij}^\epsilon (t){({\boldsymbol{L}}_j^{{x^r}}(t))^{\rm{T}}} - {\boldsymbol{\Gamma }}{{\boldsymbol{R}}_j}{({\boldsymbol{L}}_j^{{x^r}}(t))^{\rm{T}}}\\[-10pt] \end{split}$$ (37)

    初值为 ${\boldsymbol{P}}_{ij}^x(0|\; -\; 1) \;=\; {{\boldsymbol{P}}_0}$, ${\boldsymbol{P}}_{ij}^{{x^r}}(0| \;-\; 1)\;=\;0$ 和 ${\boldsymbol{P}}_{ij}^{x{x^r}}(0| - 1) = 0$.

    证明. 见附录F.

    定理5. 第$i$和第$j$个传感器子系统的$N$$( N < - 1 )$步预报误差之间的互协方差矩阵计算如下:

    $$\begin{split} &{\boldsymbol{P}}_{ij}^x(t|t + N) = {{\boldsymbol{A}}^{ - N - 1}}{\boldsymbol{P}}_{ij}^x(t + N + 1|t +N)\;\times\\ &\qquad {({{\boldsymbol{A}}^{ - N - 1}})^{\rm{T}}} \;+ \\ &\qquad\sum\limits_{k = 1}^{ - N - 1} {{{\boldsymbol{A}}^{ - N - k - 1}}} {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}{({{\boldsymbol{A}}^{ - N - k - 1}})^{\rm{T}}} \end{split}$$ (38)

    式中, ${\boldsymbol{P}}_{^{ij}}^x(t + N + 1|t + N)$由定理4给出.

    证明. 利用附录式(D2), 类似式(26)的推导过程, 可直接获得式(38).

    定理6. 第$i$和第$j$个传感器子系统的滤波$(N= 0)$误差之间和平滑$ (N>0)$误差之间的互协方差矩阵计算如下:

    $$\begin{split} {\boldsymbol{P}}_{ij}^x(t|t + N) = \;&{\boldsymbol{P}}_{ij}^x(t|t - 1)\; - \\ &\sum\limits_{{k_j} = 0}^N {{\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j})} {\boldsymbol{K}}_j^{\rm{T}}(t|t + {k_j})\; - \\ &\sum\limits_{{k_i} = {{0}}}^N {{{\boldsymbol{K}}_i}(t|t + {k_i})} {({\boldsymbol{P}}_{ji}^{x\epsilon }(t + {k_i}))^{\rm{T}}}\; + \\ &\sum\limits_{{k_i},{k_j} = 0}^N {{{\boldsymbol{K}}_i}(t{\rm{|}}t + {k_i})}\; \times \\ &{\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_j}){\boldsymbol{K}}_j^{\rm{T}}(t{\rm{|}}t + {k_j})\\[-10pt] \end{split}$$ (39)
    $$\begin{split} {\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j}) =\;& {\boldsymbol{P}}_{ij}^x(t|t - 1){({{\boldsymbol{A}}^{{k_j}}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}} \;+ \\ &{\boldsymbol{P}}_{ij}^{x{x^r}}(t|t - 1){({\boldsymbol{\Phi }}_j^{{k_j}})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}}\; - \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j} - {m_j})} \; \times \\ &[{\boldsymbol{L}}_j^{\rm{T}}(t + {k_j} - {m_j}){({{\boldsymbol{A}}^{{m_j} - 1}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}}\; + \\ &{({\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}))^{\rm{T}}}\; \times \\ &{({\boldsymbol{\Phi }}_j^{{m_j} - 1})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}}],\;{k_j} > 0\\[-10pt] \end{split}$$ (40)
    $$\begin{split} {\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j}) =\;& {({\boldsymbol{P}}_{ji}^{x{x^r}}(t|t - 1))^{\rm{T}}}{({{\boldsymbol{A}}^{{k_j}}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}}\; + \\ &{\boldsymbol{P}}_{ij}^{{x^r}}(t|t - 1){({\boldsymbol{\Phi }}_j^{{k_j}})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}} \;- \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j} - {m_j})}\;\times\\ &[({\boldsymbol{L}}_j^{\rm{T}}(t + {k_j} - {m_j})\; \times {({{\boldsymbol{A}}^{{m_j} - 1}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}}\; +\\ & {({\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}))^{\rm{T}}}\; \times \\ &{({\boldsymbol{\Phi }}_j^{{m_j} - 1})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}}],{k_j} > 0\\[-10pt] \end{split}$$ (41)

    新息协方差${\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_j})$为:

    $$\begin{split} {\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},&t + {k_j}) = {({\boldsymbol{P}}_{ji}^{x\epsilon }(t + {k_i}))^{\rm{T}}}{({{\boldsymbol{A}}^{{k_j}}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}} \;+ \\ &{({\boldsymbol{P}}_{ji}^{{x^r}\epsilon }(t + {k_i})^{\rm{T}}}{({\boldsymbol{\Phi }}_j^{{k_j}})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}} \;- \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_j} - {m_j})} \; \times \\ &[{\boldsymbol{L}}_j^{\rm{T}}(t + {k_j} - {m_j}){({{\boldsymbol{A}}^{{m_j} - 1}})^{\rm{T}}}{\boldsymbol{C}}_j^{\rm{T}} \;+ \\ &{\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}){^{\rm{T}}}{({\boldsymbol{\Phi }}_j^{{m_j} - 1})^{\rm{T}}}\tilde {\boldsymbol{C}}_j^{\rm{T}}] \;+ \\ &\sum\limits_{{m_j} = {k_j} - {k_i}}^{{k_j}} {{\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i},t + {k_j} - {m_j}} {\rm{)}}\; \times \\ &{{\rm{(}}{{\boldsymbol{C}}_j}{\boldsymbol{A}}_{}^{{m_j} - 1}{\boldsymbol{\Gamma }})^{\rm{T}}},\;\;\;{k_i} < {k_j}\\[-10pt] \end{split}$$ (42)

    其中

    $$\begin{split} {\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i},&t + {k_j} - {m_j}{\rm{)}} \;= \\ & -\sum\limits_{{m_i} = 1}^{{k_i}} {[{{\boldsymbol{C}}_i}{{\boldsymbol{A}}^{{m_i} - 1}}{{\boldsymbol{L}}_i}(t + {k_i} - {m_i})} \;+ \\ &{\tilde {\boldsymbol{C}}_i}{\boldsymbol{\Phi }}_i^{{m_i} - 1}{\boldsymbol{L}}_i^{{x^r}}(t + {k_i} - {m_i})]\; \times \\ &{\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i} - {m_i},t + {k_j} - {m_j}{\rm{)\; + }}\\ &\sum\limits_{{m_i} = 1}^{{k_i}} {{{\boldsymbol{C}}_i}{\boldsymbol{A}}_{}^{{m_i} - 1}} {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{\delta _{{k_i} - {m_i},{k_j} - {m_j}}} \;+ \\ &{\boldsymbol{R}}_i^{\rm{T}}{\delta _{{k_i},{k_j} - {m_j}}}\\[-10pt] \end{split}$$ (43)

    当${k_i} = {k_j}$时, 有${\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_i}) = {\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i})$; 当${k_i} > {k_j}$时, 有 ${\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_j}) = ({\boldsymbol{Q}}_{ji}^\epsilon (t + {k_j}, t + {k_i}))^{\rm{T}}$, ${k_i},{k_j} = {{0}}, \cdots ,N$, 其中${\boldsymbol{P}}_{ij}^{x\epsilon }(t)$和${\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t)$分别由式(33)和式(36)计算.

    证明. 见附录G.

    注8. 由增广状态定义${\bar {\boldsymbol{x}}_i}(t) = {[ {{{\boldsymbol{x}}^{\rm{T}}}(t)}\;{({\boldsymbol{x}}_i^r} (t))^{\rm{T}}}{]^{\rm{T}}}$可以看出, ${\bar {\boldsymbol{x}}_i}(t)$包含了虚拟状态${\boldsymbol{x}}_i^r(t)$, 因此在增广方法中, 求系统状态${\boldsymbol{x}}(t)$的滤波器、预报器和平滑器${\hat {\boldsymbol{x}}_i}(t|t + N)$, 都需求解虚拟状态${\boldsymbol{x}}_i^r(t)$的相应估值器$\hat {\boldsymbol{x}}_i^r(t|t + N)$. 而非增广方法(定理1 ~ 6)只需要求解虚拟状态的一步预报器$\hat {\boldsymbol{x}}_i^r(t + 1|t)$来计算新息, 不需要计算其滤波器、多步预报器和多步平滑器, 以及相应的方差矩阵和互协方差矩阵. 所以, 非增广方法可避免很多对虚拟状态的不必要计算, 简化了算法. 因此, 非增广方法更具有优越性.

    基于局部估值器${\hat {\boldsymbol{x}}_i}(t|t + N)$、局部估计误差方差矩阵${\boldsymbol{P}}_i^x(t|t + N)$、任意两个局部估计误差之间的互协方差矩阵${\boldsymbol{P}}_{ij}^x(t|t + N)$和线性最小方差意义下的三种加权最优融合算法[19]给出分布式融合估值器:

    $${\hat {\boldsymbol{x}}_o}(t|t + N) = \sum\limits_{i = 1}^L {{\boldsymbol{A}}_i^N(t)} {\hat {\boldsymbol{x}}_i}(t|t + N)$$ (44)

    式中, ${{\boldsymbol{A}}^N}(t) = [{\boldsymbol{A}}_1^N(t), \cdots ,{\boldsymbol{A}}_L^N(t)]$, ${\boldsymbol{A}}_i^N(t)$为最优融合权重, 可以为矩阵权重、对角矩阵权重和标量权重.

    当${\boldsymbol{A}}_i^N(t)$为矩阵权重时:

    $$\left\{\begin{aligned} &{{\boldsymbol{A}}^N}(t) = {({{\boldsymbol{F}}^{\rm{T}}}{{\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}})^{ - 1}}{{\boldsymbol{F}}^{\rm{T}}}{{\boldsymbol{P}}^{ - 1}}(t|t + N)\\ &{{\boldsymbol{P}}_o}(t|t + N) = {({{\boldsymbol{F}}^{\rm{T}}}{{\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}})^{ - 1}} \end{aligned}\right.$$ (45)

    式中, ${\boldsymbol{F}} = [{{\boldsymbol{I}}_n}, \cdots ,{{\boldsymbol{I}}_n}]_{n \times nL}^{\rm{T}}$, ${\boldsymbol{P}}(t|t + N)$ $= [{\boldsymbol{P}}_{ij}^x(t|t \;+ N)]_{nL \times nL}$, $i,j = 1, \cdots ,L$, 且有 ${{\boldsymbol{P}}_o}(t|t + N) \; \leq {\boldsymbol{P}}_i^x (t|t + N)$.

    当${\boldsymbol{A}}_i^N(t)$为对角矩阵权重时:

    $$\left\{\begin{aligned} &{{\boldsymbol{A}}^N}(t) = {({{\boldsymbol{F}}^{\rm{T}}}{\bar {\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}})^{ - 1}}{{\boldsymbol{F}}^{\rm{T}}}{\bar {\boldsymbol{P}}^{ - 1}}(t|t + N)\\ &{{\boldsymbol{P}}_o}(t|t + N) = {({{\boldsymbol{F}}^{\rm{T}}}{\bar {\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}})^{ - 1}}\; \times \\ &\qquad{{\boldsymbol{F}}^{\rm{T}}}{\bar {\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{P}}(t|t + N){\bar {\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}}\; \times \\ &\qquad{({{\boldsymbol{F}}^{\rm{T}}}{\bar {\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{F}})^{ - 1}} \end{aligned}\right.$$ (46)

    式中, $\bar {\boldsymbol{P}}(t|t + N) = {[\bar {\boldsymbol{P}}_{ij}^x(t|t + N)]_{nL \times nL}}$, 且有${\bar {\boldsymbol{P}}_o}(t|t \;+ N) \leq \bar {\boldsymbol{P}}_i^x(t|t + N)$, 其中${\bar {\boldsymbol{P}}_o}(t|t + N)$和$\bar {\boldsymbol{P}}_{ij}^x(t|t + N)$分别为${{\boldsymbol{P}}_o}(t|t + N)$和${\boldsymbol{P}}_{ij}^x(t|t + N)$的对角元素构成的对角矩阵.

    当${\boldsymbol{A}}_i^N(t)$为标量权重时:

    $$\left\{\begin{aligned} &{{{\boldsymbol{A}}^N}(t) = \frac{{{{\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{e}}}}{{{{\boldsymbol{e}}^{\rm{T}}}{{\boldsymbol{P}}^{ - 1}}(t|t + N){\boldsymbol{e}}}}}\\ &{{{\boldsymbol{P}}_o}(t|t + N) = \sum\limits_{i,j = 1}^L {{\boldsymbol{A}}_i^N(t)} {\boldsymbol{A}}_j^N(t){\boldsymbol{P}}_{ij}^x(t|t + N)} \end{aligned}\right.$$ (47)

    式中, ${\boldsymbol{e}} = [1, \cdots ,1]_{L \times 1}^{\rm{T}}$, ${\boldsymbol{P}}(t|t + N) \;=\; [{\rm{tr}}{\boldsymbol{P}}_{ij}^x(t|t + N)]_{L \times L}$, 且有 ${\rm{tr}}{{\boldsymbol{P}}_o}(t|t + N) \leq {\rm{tr}}{\boldsymbol{P}}_i^x(t|t + N)$.

    结合式(14) ~ (19)可得, 增广的一步预报器${{\hat{ \bar {\boldsymbol{x}}}}_i}(t \;+ 1|t) = [ {\hat {\boldsymbol{x}}_i^{\rm{T}}(t + 1|t)}\;{{{(\hat {\boldsymbol{x}}_i^r(t + 1|t))}^{\rm{T}}}} ]^{\rm{T}}$为:

    $${{\hat{ \bar {\boldsymbol{x}}}}_i}(t + 1|t) = ({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}(t){\bar {\boldsymbol{C}}_i}){{\hat{ \bar {\boldsymbol{x}}}}_i}(t|t - 1) + {\bar {\boldsymbol{L}}_i}(t){{\boldsymbol{y}}_i}(t)$$ (48)

    结合式(14)和附录式(B1) ~ (C4)可得, 误差方程${{\hat{ \bar {\boldsymbol{x}}}}_i}(t + 1|t) = {[ {\tilde {\boldsymbol{x}}_i^{\rm{T}}(t + 1|t)}\;\;{{{(\tilde {\boldsymbol{x}}_i^r(t + 1|t))}^{\rm{T}}}} ]^{\rm{T}}}$为:

    $$\begin{split} { {\tilde { \bar {\boldsymbol{x}}}}_i}(t + 1|t) = \;&({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}(t){\bar {\boldsymbol{C}}_i}){ {\tilde { \bar {\boldsymbol{x}}}}_i}(t|t - 1)\; +\\ &\bar {\boldsymbol{\Gamma }}{\bar {\boldsymbol{w}}_i}(t) - {\bar {\boldsymbol{L}}_i}(t){{\boldsymbol{v}}_i}(t) \end{split}$$ (49)

    式中, ${\bar {\boldsymbol{L}}_i}(t) = {[ {{\boldsymbol{L}}_i^{\rm{T}}(t)}\;\;{({\boldsymbol{L}}_i^{{x^r}}(t)} )^{\rm{T}}}{]^{\rm{T}}}$. 由式(49)可得方差矩阵:

    $$\begin{split} {{\boldsymbol{P}}_i}(t + 1|t) =\;& ({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}(t){\bar {\boldsymbol{C}}_i}){{\boldsymbol{P}}_i}(t|t - 1)\; \times \\ &{({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}(t){\bar {\boldsymbol{C}}_i})^{\rm{T}}} + \bar {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}_i}(t){\bar {\boldsymbol{\Gamma }}^{\rm{T}}}\; -\\ &\bar {\boldsymbol{\Gamma }}{{\boldsymbol{S}}_i}\bar {\boldsymbol{L}}_i^{\rm{T}}(t) -\bar {\boldsymbol{L}}_i^{}(t){\boldsymbol{S}}_i^{\rm{T}}{\bar {\boldsymbol{\Gamma }}^{\rm{T}}}\; + \\ &\bar {\boldsymbol{L}}_i^{}(t){\boldsymbol{Q}}_i^v\bar {\boldsymbol{L}}_i^{\rm{T}}(t) \end{split}$$ (50)

    定理7. 如果$\rho ({\boldsymbol{A}}) < 1$和$\left| {{h_i}} \right| < 1$, $i = 1, \cdots ,L$, 则在任意初始条件${\boldsymbol{\Xi }}(0) \geq 0$和$Q_i^r(0) \geq 0$下, 式(10)和式(11)的解${\boldsymbol{\Xi }}(t)$和$Q_i^r(t)$将收敛到唯一的半正定矩阵${\boldsymbol{\Xi }}$和$Q_i^r$且满足:

    $${\boldsymbol{\Xi }} = {\boldsymbol{A\Xi }}{{\boldsymbol{A}}^{\rm{T}}} + {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}$$ (51)
    $$Q_i^r = h_i^2Q_i^r + Q_i^p$$ (52)

    即${\boldsymbol{\Xi }} = \mathop {\lim }\nolimits_{t \to \infty } {\boldsymbol{\Xi }}(t)$, $Q_i^r = \mathop {\lim }\nolimits_{t \to \infty } Q_i^r(t)$. 由式(9)和式(15), 有:

    $$\left\{\begin{aligned} &{{\boldsymbol{Q}}_i} = \mathop {\lim }\limits_{t \to \infty } {{\boldsymbol{Q}}_i}(t) = {\rm{diag}}\{{{\boldsymbol{Q}}^w},{\boldsymbol{Q}}_i^\xi \}\\ &{\boldsymbol{Q}}_i^\xi = \mathop {\lim }\limits_{t \to \infty } {\boldsymbol{Q}}_i^\xi (t) = h_i^2Q_i^r{\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}}\; +\\ &\qquad\;\; Q_i^p{\boldsymbol{A\Xi }}{{\boldsymbol{A}}^{\rm{T}}} + Q_i^p{\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}} \end{aligned}\right.$$ (53)

    定理8. 如果$\rho ({\boldsymbol{A}}) < 1$和$\left| {{h_i}} \right| < 1$, $i = 1, \cdots ,L$, 且$({\bar {\boldsymbol{\Phi }}_i} - \bar {\boldsymbol{\Gamma }}{{\boldsymbol{S}}_i}{({\boldsymbol{Q}}_i^v)^{ - 1}}{\bar {\boldsymbol{C}}_i},\bar {\boldsymbol{\Gamma }}{\bar {\boldsymbol{Q}}_a})$ 是可稳对, 其中${\bar {\boldsymbol{Q}}_a}\bar {\boldsymbol{Q}}_a^{\rm{T}} = {{\boldsymbol{Q}}_i} - {{\boldsymbol{S}}_i}{({\boldsymbol{Q}}_i^v)^{ - 1}}{\boldsymbol{S}}_i^{\rm{T}}$, 则在任意初始条件${{\boldsymbol{P}}_i}(0| - 1) \geq 0$时, 方程(50)的解${{\boldsymbol{P}}_i}(t + 1|t)$将收敛到如下方程:

    $$\begin{split} {{\boldsymbol{P}}_i} =\;& ({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}{\bar {\boldsymbol{C}}_i}){{\boldsymbol{P}}_i}{({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}{\bar {\boldsymbol{C}}_i})^{\rm{T}}}\; +\\ &\bar {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}_i}{\bar {\boldsymbol{\Gamma }}^{\rm{T}}} - \bar {\boldsymbol{\Gamma }}{{\boldsymbol{S}}_i}\bar {\boldsymbol{L}}_i^{\rm{T}} - \bar {\boldsymbol{L}}_i^{}{\boldsymbol{S}}_i^{\rm{T}}{\bar {\boldsymbol{\Gamma }}^{\rm{T}}} + \bar {\boldsymbol{L}}_i^{}{\boldsymbol{Q}}_i^v\bar {\boldsymbol{L}}_i^{\rm{T}} \end{split}$$ (54)

    即${{\boldsymbol{P}}_i} = \mathop {\lim }\nolimits_{t \to \infty } {{\boldsymbol{P}}_i}(t + 1|t)$ 和${\bar {\boldsymbol{L}}_i} = \mathop {\lim }\nolimits_{t \to \infty } {\bar {\boldsymbol{L}}_i}(t)$, $i = $ $ 1, \cdots ,L.$

    且局部稳态一步预报器:

    $${{\hat{ \bar {\boldsymbol{x}}}}_i}(t + 1|t) = ({\bar {\boldsymbol{\Phi }}_i} - {\bar {\boldsymbol{L}}_i}{\bar {\boldsymbol{C}}_i}){{\hat{ \bar {\boldsymbol{x}}}}_i}(t|t - 1) + {\bar {\boldsymbol{L}}_i}{{\boldsymbol{y}}_i}(t)$$ (55)

    是渐近稳定的.

    证明. 见附录H.

    注9. 由${{\hat{ \bar {\boldsymbol{x}}}}_i}(t + 1|t) = {[ {\hat {\boldsymbol{x}}_i^{\rm{T}}(t + 1|t)}\;{{{(\hat {\boldsymbol{x}}_i^r(t + 1|t))}^{\rm{T}}}} ]^{\rm{T}}}$得 ${\hat {\boldsymbol{x}}_i}(t + 1|t) = [ {{{\boldsymbol{I}}_n}}\;\;0 ]{{\hat{ \bar {\boldsymbol{x}}}}_i}(t + 1|t) ,$ ${\boldsymbol{\varSigma }}_i^{} = \mathop {\lim }\nolimits_{t \to \infty } {\boldsymbol{P}}_i^x(t \;+ 1|t)$${\rm{ = }}\;[ {{{\boldsymbol{I}}_n}}\;\;0 ]{{\boldsymbol{P}}_i}{[ {{{\boldsymbol{I}}_n}}\;\;0 ]^{\rm{T}}}$和${{\boldsymbol{L}}_i} = [ {{{\boldsymbol{I}}_n}}\;\;0 ]{\bar {\boldsymbol{L}}_i}$. 表明本文提出的局部稳态一步预报器是渐近稳定的. 由局部稳态一步预报器的渐近稳定性可知, 局部稳态多步预报器、滤波器和平滑器也是渐近稳定的. 由文献[24]可知, 如果所有的局部稳态估值器是渐近稳定的, 那么任两个传感器子系统之间存在稳态互协方差矩阵且稳态互协方差矩阵可任选初值离线迭代计算. 从而由第2.4节分布式融合估计算法可知, 分布式融合稳态估值器也是渐近稳定的.

    注10. 集中式融合算法的观测和相应的增广状态为${{\boldsymbol{y}}_c}(t) = {[ {{\boldsymbol{y}}_1^{\rm{T}}(t)}\; \cdots \;{{\boldsymbol{y}}_L^{\rm{T}}(t)} ]^{\rm{T}}} \in {{\bf{R}}^{\sum\nolimits_{i = 1}^L {{m_i}} }}$, $\bar {\boldsymbol{x}}(t + 1) = {[ {{{\boldsymbol{x}}^{\rm{T}}}(t + 1)}\;\;{{{({\boldsymbol{x}}_1^r(t + 1))}^{\rm{T}}}}\; \cdots \;{{{({\boldsymbol{x}}_L^r(t + 1))}^{\rm{T}}}} ]^{\rm{T}}} \in$${{\bf{R}}^{n(L + 1)}}$, 相应的观测方程为${{\boldsymbol{y}}_c}(t) = {H_c}\bar {\boldsymbol{x}}(t) + {{\boldsymbol{v}}_c}(t)$, 其中${H_c}$为适当的矩阵, ${{\boldsymbol{v}}_c}(t)$为增广的观测噪声. 通常情况下, ${m_i} \leq n$即$n(L + 1) > \sum\nolimits_{i = 1}^L {{m_i}} $. 因此, 集中式融合算法的计算量级为${{\rm{O}}_c}({n^3}{(L + 1)^3})$. 根据式(44)和式(45), 有分布式矩阵加权融合算法的计算量级为${{\rm{O}}_d}({n^3}{L^3})$, 因此在融合中心有${{\rm{O}}_d}({n^3}{L^3}) \; < {{\rm{O}}_c} ({n^3}{(L + 1)^3})$. 在稳态估计存在时, 集中式和分布式的增益和方差矩阵都可离线计算, 在线只需估值的计算. 易知, 稳态集中式融合的在线计算负担为${{\rm{O}}_c}({n^2}{(L + 1)^2})$, 而稳态分布式融合的在线计算负担为${{\rm{O}}_d}({n^2}{L^2})$, 因此可明显减小融合中心的在线计算负担, 便于实时应用.

    仿真实验采用文献[11, 14]的雷达跟踪系统以及文献[24]的不间断电源系统, 基于Matlab仿真软件, 验证本文算法的有效性和可行性. 同时与现有文献算法进行比较研究.

    例1. 考虑带三个传感器的雷达跟踪系统[11, 14]:

    $${\boldsymbol{x}}(t + 1) = \left[ {\begin{array}{*{20}{c}} 1&T&0&0&0&0 \\ 0&1&1&0&0&0 \\ 0&0&{{\rho _1}}&0&0&0 \\ 0&0&0&1&T&0 \\ 0&0&0&0&1&1 \\ 0&0&0&0&0&{{\rho _2}} \end{array}} \right]{\boldsymbol{x}}(t) + {\boldsymbol{w}}(t)$$ (56)
    $${{\boldsymbol{y}}_i}(t) = ({{\boldsymbol{C}}_i} + {\tilde {\boldsymbol{C}}_i}{r_i}(t)){\boldsymbol{x}}(t) + {{\boldsymbol{v}}_i}(t),\;\;i = 1,2,3$$ (57)
    $${r_i}(t + 1) = {h_i}{r_i}(t) + {p_i}(t),\;\;i = 1,2,3$$ (58)

    式中, ${\boldsymbol{x}}(t) = {[ {s(t)}\;{\dot s(t)}\;{{U_1}(t)}\;{\theta (t)}\;{\dot \theta (t)}\;{{U_2}(t)} ]^{\rm{T}}}$表示目标的状态向量, $s(t)$表示目标的位置, $\dot s(t)$表示目标的速度, ${U_1}(t)$和${U_2}(t)$表示机动相关过程噪声; $\theta (t)$表示目标的方位角, $\dot \theta (t)$表示方位角速度; $T$表示采样周期. 过程噪声${\boldsymbol{w}}(t) = {[ 0\;0\;{{w_1}(t)}\;0\;0\;{{w_2}(t)} ]^{\rm{T}}}$和观测噪声${{\boldsymbol{v}}_i}(t)$, $i = 1,2,3$, 满足关系${{\boldsymbol{v}}_i}(t) = {{\boldsymbol{d}}_i}{\boldsymbol{w}}(t) \;+ {{\boldsymbol{\eta }}_i}(t)$, ${{\boldsymbol{d}}_i}$表示相关系数, ${p_i}(t)$和${{\boldsymbol{\eta }}_i}(t)$是零均值方差为$Q_i^p$和${\boldsymbol{Q}}_i^\eta $的独立白噪声. 仿真中, 取: ${\rho _1} = 0.4$, ${\rho _2} = 0.5$, $T = 1$s, ${{\boldsymbol{C}}_1} = \left[ {\begin{aligned} 1\;0\;0\;0\;0\;0 \\ 0\;0\;0\;1\;0\;0 \end{aligned}} \right]$

    $Q_1^p = {(0.05)^2}(1 - {{\rm{e}}^{- \frac {2T}{5}}})$, ${{\boldsymbol{C}}_2} = \left[ {\begin{aligned} 1\;0\;0\;0\;0\;0 \\ 0\;0\;0\;0\;0\;0 \end{aligned}} \right]$

    $Q_2^p = {(0.05)^2}(1 - {{\rm{e}}^{ - \frac {2T}{5}}})$, ${{\boldsymbol{C}}_3} = \left[ {\begin{aligned} 0\;0\;0\;0\;0\;0 \\ 0\;0\;0\;1\;0\;0 \end{aligned}} \right]$

    $Q_3^p = {(0.04)^2}(1 - {{\rm{e}}^{ - \frac {2T}{5}}})$, ${\tilde {\boldsymbol{C}}_1} = \left[ {\begin{aligned} 1\;0\;0\;0\;0\;0 \\ 0\;0\;0\;{\rm{1}}\;0\;0 \end{aligned}} \right]$

    ${\boldsymbol{Q}}_1^\eta = {\rm{diag}}\{{{(1000)^2},{(0.01)^2}}\}$, ${\tilde {\boldsymbol{C}}_2} = \left[ {\begin{aligned} 1\;0\;0\;0\;0\;0 \\ 0\;0\;0\;0\;0\;0 \end{aligned}} \right]$

    ${\boldsymbol{Q}}_2^\eta = {\rm{diag}\{{(900)^2},{(0.02)^2}}\}$, ${\tilde {\boldsymbol{C}}_{\rm{3}}} = \left[ {\begin{aligned} {{0}}\;0\;0\;0\;0\;0 \\ 0\;0\;0\;{\rm{1}}\;0\;0 \end{aligned}} \right]$

    ${\boldsymbol{Q}}_3^\eta = {\rm{diag}}\{{(800)^2},{(0.01)^2}\}$, ${d_1} = 0$

    ${{\boldsymbol{d}}_{\rm{2}}} = \left[ {\begin{aligned} 0\;0\;{0.1}\;0\;0\;0 \\ 0\;0\;0\;{{0}}\;0\;{0.1} \end{aligned}} \right]$, ${{\boldsymbol{d}}_3} = \left[ {\begin{aligned} 0\;0\;0\;0\;0\;{0.1} \\ 0\;0\;{0.1}\;{{0}}\;0\;0 \end{aligned}} \right]$

    ${h_1} = {h_2} = {{\rm{e}}^{\frac{T}{5}}}$, ${h_3} = 0.9{{\rm{e}}^{\frac {T}{5}}}$${{\boldsymbol{Q}}^w} = {\rm{diag}}\left\{0,0, {\left(\dfrac{103}{3}\right)^2}, 0, 0,1.1 \times {\rm{1}}{{{0}}^{{\rm{ - 8}}}}\right\}$ ${\boldsymbol{x}}(0) = [ {100}\;\;{400}\;\;0\;\;1\;\;1.07 \times {\rm{1}}{{{0}}^{{\rm{ - 3}}}}\;\; 0 ]^{\rm{T}}$

    目的是求位置$s(t)$的分布式矩阵加权融合滤波器${\hat s_o}(t|t)$.

    应用定理3$(N = 0)$, 可得位置$s(t)$的局部滤波器${\hat s_i}(t|t)$. 根据定理6$(N = 0)$, 可得任意两个局部滤波误差之间的互协方差矩阵${\boldsymbol{P}}_{ij}^s(t|t)$. 最后, 应用矩阵加权算法, 可得矩阵加权融合器滤波器${\hat s_o}(t|t)$. 图1给出了分布式矩阵加权融合滤波器的位置跟踪曲线. 图2为本文局部滤波算法、文献[11]滤波算法(第1个传感器子系统过程噪声与观测噪声独立${d_1} = 0 )$和本文矩阵加权融合算法、状态增广方法基于1000次蒙特卡洛实验的滤波均方根误差(Root mean squared error, RMSE)比较曲线:

    $$ {{{\rm{RMSE}}}}_{}^l(t) = \sqrt {{\frac{\left[ {\sum\limits_{l = 1}^{1000} {{{({s^l}(t) - \hat s_{}^l(t|t))}^2}} } \right]} {1000}}}$$

    式中, $l$表示第$l$次实验. 由图2可知, 3种滤波算法具有相同的估计精度, 且均低于矩阵加权融合滤波的估计精度.

    图 1  分布式矩阵加权融合滤波的位置跟踪性能
    Fig. 1  The position tracking performance of distributed fusion filter weighted by matrices
    图 2  RMSEs比较图
    Fig. 2  Comparisons of RMSEs

    例2. 考虑如下1 kVA不间断电源系统[24]:

    $$\begin{split} {\boldsymbol{x}}(t + 1) = \;&\left[ {\begin{array}{*{20}{c}} {{{0}}{\rm{.9226}}}&{{\rm{ - 0}}{\rm{.6330}}}&0 \\ 1&0&0 \\ 0&1&0 \end{array}} \right]{\boldsymbol{x}}(t)\;+\\ & \left[ {\begin{array}{*{20}{c}} {0.5} \\ 0 \\ {0.2} \end{array}} \right]w(t) \end{split}$$

    观测方程同式(57), 其中, ${{\boldsymbol{C}}_1} = {{\boldsymbol{C}}_2} = {{\boldsymbol{C}}_3} = [{23.738} {20.287}\;\;0 ],$ ${Q^w} = 5,$ ${\tilde {\boldsymbol{C}}_1} = [ 8\;\;{12}\;\;3 ],$ ${\tilde {\boldsymbol{C}}_2} = [ {10}\;\;8\;\;5 ],$ ${\tilde {\boldsymbol{C}}_3} = [ 8\;\;6\;\;{10} ],$ ${h_1} = 0.6,$ ${h_2} = - 0.7,$ ${h_3} = 0.8,$ $Q_1^p = 2,$ $Q_2^p = 3,$ $Q_3^p = 0.5,$ ${d_1} = - 0.5,$ ${d_2} = 0.4,$ ${d_3} = - 1,$ $Q_1^\eta = 5,$ $Q_2^\eta = 15,$ $Q_3^\eta = 10,$ ${P_0} = 0.1{I_3}$.

    图3为局部滤波(Local filtering, LF)、矩阵加权(Matrix weighted, MF)、对角矩阵加权(Diagonal matrix weighting, DF)、标量加权(Scalar weighted, SF) (文献[27]的融合算法)和集中式融合滤波(Centralized filtering, CF)误差方差的比较图. 图3(a) ~ 图3(c)表示相应的状态分量. 可以看出, 稳态融合滤波是存在的, 估计精度由高到低依次为: CF、MF、DF、SF、LF. 图4为矩阵加权融合预报、滤波和平滑的误差方差比较. 由图4可以看出, 两步平滑的误差方差最小, 估计精度最高; 两步预报的误差方差最大, 估计精度最低.

    图 3  LF、SF、DF、MF、CF的误差方差
    Fig. 3  Variances of LFs, SF, DF, MF, and CF
    图 4  矩阵加权融合预报器、滤波器和平滑器的估计误差方差
    Fig. 4  Variances of fusion predictor, filter and smoother weighted by matrices

    本文针对具有时间相关乘性噪声的线性离散随机多传感器系统, 首先, 通过引入虚拟状态和虚拟过程噪声, 构建了虚拟状态所满足的虚拟状态方程. 在线性最小方差意义下设计了相应的局部估值器, 包括滤波器、预报器和平滑器. 推导了任两个局部估计误差之间的互协方差矩阵的计算公式, 进而提出基于线性最小方差意义下的分布式加权融合算法, 给出了相应的分布式融合估值器. 最后, 分析了算法的稳定性, 给出了稳态滤波器存在的一个充分条件. 在后续研究工作中, 将结合文献[30]提出的模型参数辨识方法, 研究系统含有未知模型参数, 尤其是乘性噪声参数未知时的自校正估计问题.

    证明. 由式(7)和注2可得:${} $

    $$\begin{split} {\rm{E}}[{{\boldsymbol{\xi }}_i}(t)] =\;& {h_i}{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}(t)]{\rm{E}}[{\boldsymbol{w}}(t)]\; + \\ &{\rm{E}}[{p_i}(t)]{\rm{E}}[{\boldsymbol{x}}(t + 1)] \end{split}\tag{A1}$$
    $$\begin{split} {\rm{E}}[{{\boldsymbol{\xi }}_i}(t){{\boldsymbol{w}}^{\rm{T}}}(k)] =\;& {h_i}{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}(t)]{\rm{E}}[{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k)]\; + \\ &{\boldsymbol{A}}{\rm{E}}[{p_i}(t)]{\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(k)] \;+ \\ &{\boldsymbol{\Gamma }}{\rm{E}}[{p_i}(t)]{\rm{E}}[{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k)] \end{split}\tag{A2}$$
    $$\begin{split} {\rm{E}}[{{\boldsymbol{\xi }}_i}(t){\boldsymbol{v}}_j^{\rm{T}}(k)] =\;& {h_i}{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}(t)]{\rm{E}}[{\boldsymbol{w}}(t){\boldsymbol{v}}_j^{\rm{T}}(k)] \;+ \\ &{\boldsymbol{A}}{\rm{E}}[{p_i}(t)]{\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{v}}_j^{\rm{T}}(k)]\; + \\ &{\boldsymbol{\Gamma }}{\rm{E}}[{p_i}(t)]{\rm{E}}[{\boldsymbol{w}}(t){\boldsymbol{v}}_j^{\rm{T}}(k)] \end{split}\tag{A3}$$

    利用${\rm{E}}[{r_i}(t)] = 0$和${\rm{E}}[{p_i}(t)] = 0$, 可得式(8)成立.

    将式(7)代入${\rm{E}}[{{\boldsymbol{\xi }}_i}(t){\boldsymbol{\xi }}_i^{\rm{T}}(k)]$, 并利用假设1, ${r_i}(t) \bot {\boldsymbol{w}}(k)$和${p_i}(t) \bot {\boldsymbol{x}}(k)$, 可得:

    $$\begin{split} {\rm{E}}[{{\boldsymbol{\xi }}_i}(t){\boldsymbol{\xi }}_i^{\rm{T}}(k)] =\;& h_i^2{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{r_i}{\rm{(}}k{\rm{)]{\rm{E}}[}}{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){\rm{]}}{{\boldsymbol{\Gamma }}^{\rm{T}}} \;+ \\ &{\boldsymbol{A}}{\rm{E}}[{p_i}{\rm{(}}t{\rm{)}}{p_i}{\rm{(}}k{\rm{)]{\rm{E}}[}}{\boldsymbol{x}}(t){{\boldsymbol{x}}^{\rm{T}}}(k){\rm{]}}{{\boldsymbol{A}}^{\rm{T}}} \;+ \\ &{\boldsymbol{\Gamma }}{\rm{E}}[{p_i}{\rm{(}}t{\rm{)}}{p_i}{\rm{(}}k{\rm{)]{\rm{E}}[}}{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){\rm{]}}{{\boldsymbol{\Gamma }}^{\rm{T}}}\; + \\ &\{ {h_i}{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{\boldsymbol{w}}(t){{\boldsymbol{x}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)]}}{{\boldsymbol{A}}^{\rm{T}}}\; + \\ &{h_i}{\boldsymbol{\Gamma }}{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)]}}{{\boldsymbol{\Gamma }}^{\rm{T}}} \;+ \\ &{\boldsymbol{A}}{\rm{E}}[{p_i}{\rm{(}}t{\rm{)}}{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)]}}{{\boldsymbol{\Gamma }}^{\rm{T}}}\} + {\{ *\} ^{\rm{T}}} \end{split}\tag{A4}$$

    式中, ${\{ *\} ^{\rm{T}}}$表示与邻近 $\{ \cdot \}$中的项互为转置. 利用${r_i}(t) \bot {p_i}(k)$和${\boldsymbol{x}}(t) \bot {\boldsymbol{w}}(k)$, $k \geq t$, 可得:

    $$\left\{\begin{aligned} &{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)] = }}\\ &\qquad{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{p_i}{\rm{(}}k{\rm{)}}]{\rm{E}}[{\boldsymbol{w}}(t){{\boldsymbol{w}}^{\rm{T}}}(k)] = {{0}}\\ &{\rm{E}}[{p_i}{\rm{(}}t{\rm{)}}{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)] = }}\\ &\qquad{\rm{E}}[{p_i}{\rm{(}}t{\rm{)}}{p_i}{\rm{(}}k{\rm{)}}]{\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(k)] = 0\\ &{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{\boldsymbol{w}}(t){{\boldsymbol{x}}^{\rm{T}}}(k){p_i}{\rm{(}}k{\rm{)] = }}\\ &\qquad{\rm{E}}[{r_i}{\rm{(}}t{\rm{)}}{p_i}{\rm{(}}k{\rm{)}}]{\rm{E}}[{\boldsymbol{w}}(t){{\boldsymbol{x}}^{\rm{T}}}(k)] = {{0}} \end{aligned}\right.\tag{A5}$$

    类似地, 当$i \ne j$时, 利用${r_i}(t) \bot {r_j}(t)$和${p_i}(t) \bot {p_j}(t)$, 可得${\rm{E}}[{{\boldsymbol{\xi }}_i}(t){\boldsymbol{\xi }}_j^{\rm{T}}(k)] = 0$; 并将式(A5)代入式(A4), 可得式(9).

    利用式(1), 可得${\boldsymbol{x}}(t)$的二阶矩阵为:

    $$\begin{split} {\boldsymbol{\Xi }}(t + 1) =\;& {\rm{E}}[{\boldsymbol{x}}(t + 1){{\boldsymbol{x}}^{\rm{T}}}(t + 1)] = \\ &{\boldsymbol{A\Xi }}(t){{\boldsymbol{A}}^{\rm{T}}} + {\boldsymbol{A}}{\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(t)]{{\boldsymbol{\Gamma }}^{\rm{T}}} \;+ \\ &{\boldsymbol{\Gamma }}{\rm{E}}[{\boldsymbol{w}}(t){{\boldsymbol{x}}^{\rm{T}}}(t)]{{\boldsymbol{A}}^{\rm{T}}} + {\boldsymbol{\Gamma }}{{\boldsymbol{Q}}^w}{{\boldsymbol{\Gamma }}^{\rm{T}}} \end{split}\tag{A6}$$

    利用${\boldsymbol{x}}(t) \bot {\boldsymbol{w}}(t)$, 可得式(10). 类似地, 由式(3)和${r_i}(t) \bot {p_i}(t)$, 可得式(11).

    证明. 由射影理论[1]可知, 新息表达式为${{\boldsymbol{\epsilon }}_i}(t) = {{\boldsymbol{y}}_i}(t) \; - {\hat {\boldsymbol{y}}_i} (t|t - 1)$, 对观测方程式(6)等号两边同时在线性空间$L({{\boldsymbol{y}}_i}(1), \cdots ,{{\boldsymbol{y}}_i}(t))$上取投影, 利用${{\boldsymbol{v}}_i}(t) \bot L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t - 1)),$可得${\hat {\boldsymbol{y}}_i}(t|t - 1) = {{\boldsymbol{C}}_i}{\hat {\boldsymbol{x}}_i}(t|t - 1) + {\tilde {\boldsymbol{C}}_i}\hat {\boldsymbol{x}}_i^r(t|t - 1)$. 式(16)得证. 将式(6)代入式(16), 新息${{\boldsymbol{\epsilon }}_i}(t)$可重写为:

    $${{\boldsymbol{\epsilon }}_i}(t) = {{\boldsymbol{C}}_i}{\tilde {\boldsymbol{x}}_i}(t|t - 1) + {\tilde {\boldsymbol{C}}_i}\tilde {\boldsymbol{x}}_i^r(t|t - 1) + {{\boldsymbol{v}}_i}(t)\tag{B1}$$

    将式(B1)代入${\boldsymbol{Q}}_i^\epsilon (t) = {\rm{E}}[{{\boldsymbol{\epsilon }}_i}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)]$, 注意到${\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1)\cdot {(\tilde {\boldsymbol{x}}_i^r(t|t - 1))^{\rm{T}}}{\rm{]}} = {\boldsymbol{P}}_i^{x{x^r}}(t|t - 1)$、${\tilde {\boldsymbol{x}}_i}(t|t - 1) \bot {{\boldsymbol{v}}_i}(t)$和 $\tilde {\boldsymbol{x}}_i^r(t|t \;- 1) \bot {{\boldsymbol{v}}_i}(t)$, 可得式(17).

    证明. 由递推射影公式[1], 有:

    $${\hat {\boldsymbol{x}}_i}(t + 1|t) = {\hat {\boldsymbol{x}}_i}(t + 1|t - 1) + {{\boldsymbol{L}}_i}(t){{\boldsymbol{\epsilon }}_i}(t)\tag{C1}$$

    式中, 预报增益${{\boldsymbol{L}}_i}(t) = {\rm{E}}[{\boldsymbol{x}}(t + 1){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)]{({\boldsymbol{Q}}_i^\epsilon (t))^{ - 1}}$. 对式(4)等号两边在$L({{\boldsymbol{y}}_i}(1), \cdots ,{{\boldsymbol{y}}_i}(t - 1))$取投影, 注意到${\boldsymbol{w}}(t) \bot L({{\boldsymbol{y}}_i}(1), \cdots , {{\boldsymbol{y}}_i}(t - 1))$, 有${\hat {\boldsymbol{x}}_i}(t + 1|t - 1) = {\boldsymbol{A}}{\hat {\boldsymbol{x}}_i}(t|t \;- 1)$, 将其代入式(C1), 可得式(18). 下面求解增益矩阵${{\boldsymbol{L}}_i}(t)$.

    由状态方程(4), 可得:

    $$\begin{split} {\rm{E}}[{\boldsymbol{x}}(t + 1){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)] =\;& {\boldsymbol{A}}{\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)]\; + \\ &{\boldsymbol{\Gamma }}{\rm{E}}[{\boldsymbol{w}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)] \end{split}\tag{C2}$$

    式中, ${\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)] = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)]$, ${\rm{E}}[{\boldsymbol{w}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t)] = {\rm{E}}[{\boldsymbol{w}}(t)\;\times {\boldsymbol{y}}_i^{\rm{T}}(t)] = {\rm{E}}[{\boldsymbol{w}}(t){\boldsymbol{v}}_i^{\rm{T}}(t)] = {{\boldsymbol{R}}_i}$. 上文推导中利用了${\boldsymbol{x}}(t) = {\tilde {\boldsymbol{x}}_i}(t|t \;- 1) + {\hat {\boldsymbol{x}}_i}(t|t - 1)$、${\hat {\boldsymbol{x}}_i}(t|t - 1) \bot {{\boldsymbol{\epsilon }}_i}(t)$和 ${\boldsymbol{w}}(t) \bot {\hat {\boldsymbol{y}}_i}(t|t - 1)$. 由式(B1)和式(C2), 可直接得到式(20). 类似上面的推导过程, 并注意到${{\boldsymbol{\xi }}_i}(t) \bot {{\boldsymbol{\epsilon }}_i}(t)$, 则式(19)和式(21)可得证.

    用式(4)减式(18), 用式(5)减式(19), 得一步预报误差方程为:

    $${\tilde {\boldsymbol{x}}_i}(t + 1|t) = {\boldsymbol{A}}{\tilde {\boldsymbol{x}}_i}(t|t - 1) + {\boldsymbol{\Gamma w}}(t) - {{\boldsymbol{L}}_i}(t){{\boldsymbol{\epsilon }}_i}(t)\tag{C3}$$
    $$\tilde {\boldsymbol{x}}_i^r(t + 1|t) = {{\boldsymbol{\Phi }}_i}\tilde {\boldsymbol{x}}_i^r(t|t - 1) + {{\boldsymbol{\xi }}_i}(t) - {\boldsymbol{L}}_i^r(t){{\boldsymbol{\epsilon }}_i}(t)\tag{C4}$$

    移项得:

    $${\tilde {\boldsymbol{x}}_i}(t + 1|t) + {{\boldsymbol{L}}_i}(t){{\boldsymbol{\epsilon }}_i}(t) = {\boldsymbol{A}}{\tilde {\boldsymbol{x}}_i}(t|t - 1) + {\boldsymbol{\Gamma w}}(t)\tag{C5}$$
    $$\tilde {\boldsymbol{x}}_i^r(t + 1|t) + {\boldsymbol{L}}_i^r(t){{\boldsymbol{\epsilon }}_i}(t) = {{\boldsymbol{\Phi }}_i}\tilde {\boldsymbol{x}}_i^r(t|t - 1) + {{\boldsymbol{\xi }}_i}(t)\tag{C6}$$

    对式(C5)等号两边取方差, 利用${\tilde {\boldsymbol{x}}_i}(t + 1|t) \bot {{\boldsymbol{\epsilon }}_i}(t)$, 整理可得式(22). 类似地, 利用$\tilde {\boldsymbol{x}}_i^r(t + 1|t) \bot {{\boldsymbol{\epsilon }}_i}(t)$、${{\boldsymbol{\xi }}_i}(t) \bot \tilde {\boldsymbol{x}}_i^{}(t|t \;- 1)$、${{\boldsymbol{\xi }}_i}(t) \bot \tilde {\boldsymbol{x}}_i^r(t|t - 1)$和${{\boldsymbol{\xi }}_i}(t) \bot {\boldsymbol{w}}(t)$, 式(23)和式(24)得证.

    证明. 对式(1)迭代, 有:

    $${\boldsymbol{x}}(t) = {{\boldsymbol{A}}^{ - N - 1}}{\boldsymbol{x}}(t + N + 1) + \sum\limits_{k = 1}^{ - N - 1} {{{\boldsymbol{A}}^{ - N - k - 1}}} {\boldsymbol{\Gamma w}}(t - k)\tag{D1}$$

    对式(D1)等号两边同时在线性空间$L({{\boldsymbol{y}}_i}(1), {{\boldsymbol{y}}_i}(2), \cdots , {{\boldsymbol{y}}_i}(t + N))$上取射影, 由注1有${\boldsymbol{w}}(t - k) \bot L({{\boldsymbol{y}}_i}(1),{{\boldsymbol{y}}_i}(2), \cdots , {{\boldsymbol{y}}_i}(t + N))$, $k = 1, \cdots , - N - 1$, 得式(25)成立. 由式(D1)减式(25), 可得:

    $$\begin{split} {\tilde {\boldsymbol{x}}_i}(t{\rm{|}}t + N) =\;& {{\boldsymbol{A}}^{ - N - 1}}{\tilde {\boldsymbol{x}}_i}(t + N + 1|t + N) \;+ \\ &\sum\limits_{k = 1}^{ - N - 1} {{{\boldsymbol{A}}^{ - N - k - 1}}} {\boldsymbol{\Gamma w}}(t - k) \end{split}\tag{D2}$$

    将式(D2)代入$N$步预报误差方差矩阵${\boldsymbol{P}}_i^x(t|t{\rm{ + }}N) = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t + N)\tilde {\boldsymbol{x}}_i^{\rm{T}}(t|t + N)]$中, 利用假设1和${\boldsymbol{w}}(t - k) \bot {\tilde {\boldsymbol{x}}_i}(t\; + N + 1|t + N)$, $k = 1, \cdots , - N - 1$, 得式(26).

    证明. 由递推射影公式[1]可直接得式(27). 将式(B1)代入$N$步平滑增益矩阵${{\boldsymbol{K}}_i}(t|t + N){\rm{ = {\rm{E}}}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t + N)]({\boldsymbol{Q}}_i^\epsilon (t + N))^{ - 1}$,利用${\boldsymbol{x}}(t) \bot {{\boldsymbol{v}}_i}(t + N)$, 可得:

    $$\begin{split} {{\boldsymbol{K}}_i}(t|t + N) =\;& ({\rm{E}}[{\boldsymbol{x}}(t)\tilde {\boldsymbol{x}}_i^{\rm{T}}(t{\rm{ + }}N|t{\rm{ + }}N - 1)]{\boldsymbol{C}}_i^{\rm{T}} \;+ \\ &{\rm{E}}[{\boldsymbol{x}}(t){(\tilde {\boldsymbol{x}}_i^r(t{\rm{ + }}N|t{\rm{ + }}N - 1))^{\rm{T}}}]\tilde {\boldsymbol{C}}_i^{\rm{T}})\; \times \\ &{({\boldsymbol{Q}}_i^\epsilon (t + N))^{ - 1}} \end{split}\tag{E1}$$

    定义${\boldsymbol{\Omega }}_i^x(t + N|t + N - 1) = {\rm{E}}[{\boldsymbol{x}}(t)\tilde {\boldsymbol{x}}_i^{\rm{T}}(t{\rm{ + }}N|t{\rm{ + }}N - 1)]$和${\boldsymbol{\Omega }}_i^{x{x^r}}(t + N|t + N - 1) = {\rm{E}}[{\boldsymbol{x}}(t){(\tilde {\boldsymbol{x}}_i^r(t{\rm{ + }}N|t{\rm{ + }}N - 1))^{\rm{T}}}]$, 并将其代入式(E1), 可得式(28)成立.

    将式(C3)中的$t$替换为$t + N$, 并代入${\rm{E}}[{\boldsymbol{x}}(t)\tilde {\boldsymbol{x}}_i^{\rm{T}}(t{\rm{ + }}N + 1|t{\rm{ + }}N)]$, 可得:

    $$\begin{split} {\boldsymbol{\Omega }}_i^x(t + N + 1|t + N) =\;& {\boldsymbol{\Omega }}_i^x(t + N|t + N - 1){{\boldsymbol{A}}^{\rm{T}}}\; -\\ &{\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t{\rm{ + }}N)]{\boldsymbol{L}}_i^{\rm{T}}(t{\rm{ + }}N) \;+ \\ &{\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(t{\rm{ + }}N)]{{\boldsymbol{\Gamma }}^{\rm{T}}} \end{split}\tag{E2}$$

    利用平滑增益定义, 可得${\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t + N)] =$${{\boldsymbol{K}}_i}(t|t\; + N){\boldsymbol{Q}}_i^\epsilon (t + N)$, 由注1可知${\rm{E}}[{\boldsymbol{x}}(t){{\boldsymbol{w}}^{\rm{T}}}(t{\rm{ + }}N)] = 0$, 代入式(E2), 得式(29). 由式(C4)和${\rm{E}}[{\boldsymbol{x}}(t){({{\boldsymbol{\xi }}_i}(t + N))^{\rm{T}}}] = 0$, 可得:

    $$\begin{split} &{\boldsymbol{\Omega }}_i^{x{x^r}}(t + N + 1|t + N) = {\boldsymbol{\Omega }}_i^{x{x^r}}(t + N|t + N - 1){\boldsymbol{\Phi }}_i^{\rm{T}} \;- \\ &\;\;\;\;\;\;{\rm{E}}[{\boldsymbol{x}}(t){\boldsymbol{\epsilon }}_i^{\rm{T}}(t + N)]{({\boldsymbol{L}}_i^{{x^r}}(t + N))^{\rm{T}}} \end{split}\tag{E3}$$

    由式(E3), 可得式(30).

    用$ {\boldsymbol{x}}(t)$减去式(27), 可得:

    $${\tilde {\boldsymbol{x}}_i}(t|t + N) = {\tilde {\boldsymbol{x}}_i}(t|t + N - 1) - {{\boldsymbol{K}}_i}(t|t + N){{\boldsymbol{\epsilon }}_i}(t + N)\tag{E4}$$

    将式(E4)中的${{\boldsymbol{K}}_i}(t|t + N){{\boldsymbol{\epsilon }}_i}(t + N)$移项, 利用${\tilde {\boldsymbol{x}}_i}(t|t \;+ N) \bot {{\boldsymbol{\epsilon }}_i}(t + N)$, 可得式(31).

    证明. 利用式(C3)和式(C4), 并注意到${\boldsymbol{P}}_{ij}^{x\epsilon }(t) = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t\;- 1){\boldsymbol{\epsilon }}_j^{\rm{T}}(t)]$和${\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t) =$ ${\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1){\boldsymbol{\epsilon }}_j^{\rm{T}}(t)]$, 可直接得式(32)、式 (35)和式(37). 利用式(B1)和假设1, 可直接得式(33)、式 (34)和式(36).

    证明. 将式(E4)迭代, 有:

    $${\tilde {\boldsymbol{x}}_i}(t|t + N) = {\tilde {\boldsymbol{x}}_i}(t|t - 1) - \sum\limits_{{k_i} = 0}^N {{{\boldsymbol{K}}_i}(t|t + {k_i}){{\boldsymbol{\epsilon }}_i}(t + {k_i})} \tag{G1}$$

    将式(G1)代入第$i$和第$j$个传感器子系统的$N \; (N > 0)$步平滑误差协方差矩阵${\boldsymbol{P}}_{ij}^x(t|t + N) = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t + N) \tilde {\boldsymbol{x}}_j^{\rm{T}}(t|t\; + N)]$, 可得:

    $$\begin{split} &{\boldsymbol{P}}_{ij}^x(t|t + N) = {\boldsymbol{P}}_{ij}^x(t|t - 1)\; + \\ &\;\;\;\;\;\;\;\sum\limits_{{k_i},{k_j} = 0}^N {{{\boldsymbol{K}}_i}(t|t + {k_i})} {\rm{E}}[{{\boldsymbol{\epsilon }}_i}(t + {k_i}) \;\times \\ &\;\;\;\;\;\;\;{\boldsymbol{\epsilon }}_j^{\rm{T}}(t + {k_j})]{\boldsymbol{K}}_j^{\rm{T}}(t|t + {k_j}) \;- \\ &\;\;\;\;\;\;\;\sum\limits_{{k_j} = 0}^N {{\rm{E}}[{{\tilde {\boldsymbol{x}}}_i}(t|t - 1){\boldsymbol{\epsilon }}_j^{\rm{T}}(t + {k_j})]} {\boldsymbol{K}}_j^{\rm{T}}(t|t + {k_j}) \;- \\ &\;\;\;\;\;\;\;{\kern 1pt} \sum\limits_{{k_i} = 0}^N {{{\boldsymbol{K}}_i}(t|t + {k_i})} {\rm{E}}[{{\boldsymbol{\epsilon }}_i}(t + {k_i})\tilde {\boldsymbol{x}}_j^{\rm{T}}(t|t - {\rm{1}})] \end{split}\tag{G2}$$

    定义${\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j}) = {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1){\boldsymbol{\epsilon }}_j^{\rm{T}}(t + {k_j})]$、${\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j}) = {\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1){\boldsymbol{\epsilon }}_j^{\rm{T}}(t + {k_j})]$并将其代入式(G2), 得式(39). 当${k_j} = 0$时, ${\boldsymbol{P}}_{ij}^{x\epsilon }(t)$由式(33)计算, 下面推导${k_j} > 0$的情况.

    将式(B1)中的$t$换成$t + {k_j}$, 可得$t + {k_j}$时刻的新息为:

    $$\begin{split} {{\boldsymbol{\epsilon }}_j}(t+{k_j}) =\;& {{\boldsymbol{C}}_j}{\tilde {\boldsymbol{x}}_j}(t+{k_j}|t+{k_j} - 1)\; + \\ &{\tilde {\boldsymbol{C}}_j}\tilde {\boldsymbol{x}}_j^r(t+{k_j}|t+{k_j} - 1) + {{\boldsymbol{v}}_j}(t\;{ + }\;{k_j}) \end{split}\tag{G3}$$

    由式(G3)、${\tilde {\boldsymbol{x}}_i}(t|t - 1) \bot {{\boldsymbol{v}}_j}(t{\rm{ + }}{k_j})$ 和 $\tilde {\boldsymbol{x}}_i^r(t|t - 1)\bot {{\boldsymbol{v}}_j}(t\;+ {k_j})$, 可得:

    $$\begin{split} {\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j}) =\;& {\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1) \;\times \\ &\tilde {\boldsymbol{x}}_j^{\rm{T}}(t + {k_j}{\rm{|}}t + {k_j} - 1)]{\boldsymbol{C}}_j^{\rm{T}}\; + \\ &{\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1){(\tilde {\boldsymbol{x}}_j^r(t + {k_j}{\rm{|}}t + {k_j} - 1))^{\rm{T}}}]\tilde {\boldsymbol{C}}_j^{\rm{T}} \end{split}\tag{G4}$$
    $$\begin{split} {\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j}) =\;& {\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1) \;\times \\ &\tilde {\boldsymbol{x}}_j^{\rm{T}}(t + {k_j}{\rm{|}}t + {k_j} - 1)]{\boldsymbol{C}}_j^{\rm{T}}\; + \\ &{\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1){(\tilde {\boldsymbol{x}}_j^r(t + {k_j}{\rm{|}}t \;+ {k_j} - 1))^{\rm{T}}}]\tilde {\boldsymbol{C}}_j^{\rm{T}} \end{split}\tag{G5}$$

    利用式(C3)和式(C4)将一步预报误差${\tilde {\boldsymbol{x}}_j}(t + {k_j}| t\; + {k_j} - 1)$和$\tilde {\boldsymbol{x}}_j^r(t + {k_j}|t + {k_j} - 1)$分别递推到$t$时刻, 有:

    $$\begin{split} &{\tilde {\boldsymbol{x}}_j}(t + {k_j}|t + {k_j} - 1) = {{\boldsymbol{A}}^{{k_j}}}{\tilde {\boldsymbol{x}}_j}(t|t - 1)\; - \\ &\;\;\;\;\;\;\sum\limits_{{m_j} = 1}^{{k_j}} {{{\boldsymbol{A}}^{{m_j} - 1}}{{\boldsymbol{L}}_j}(t + {k_j} - {m_j})} \,{{\boldsymbol{\epsilon }}_j}(t + {k_j} - {m_j})\; +\\ &\;\;\;\;\;\;\sum\limits_{{m_j} = 1}^{{k_j}} {{{\boldsymbol{A}}^{{m_j} - 1}}} {\boldsymbol{\Gamma w}}(t + {k_j} - {m_j}) \end{split}\tag{G6}$$
    $$\begin{split} &\tilde {\boldsymbol{x}}_j^r(t + {k_j}|t + {k_j} - 1) = {\boldsymbol{\Phi }}_j^{{k_j}}\tilde {\boldsymbol{x}}_j^r(t|t - 1)\; -\\ &\;\;\;\;\;\;\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{\Phi }}_j^{{m_j} - 1}} {\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}){{\boldsymbol{\epsilon }}_j}(t + {k_j} - {m_j})\; + \\ &\;\;\;\;\;\;\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{\Phi }}_j^{{m_j} - 1}} {{\boldsymbol{\xi }}_j}(t + {k_j} - {m_j}) \end{split}\tag{G7}$$

    由 ${\boldsymbol{x}}(t) \bot {\boldsymbol{w}}(t + {k_j} - {m_j})$ 和$L({{\boldsymbol{y}}_i}(1),\cdots,{{\boldsymbol{y}}_i}(t - 1)) \bot {\boldsymbol{w}}(t \;+ {k_j} - {m_j})$, 可得${\tilde {\boldsymbol{x}}_i}(t|t - 1) \bot {\boldsymbol{w}}(t + {k_j} - {m_j})$、${k_j} - {m_j} \geq 0$. 类似地, 有${\tilde {\boldsymbol{x}}_i}(t|t - 1) \bot {{\boldsymbol{\xi }}_j}(t + {k_j} - {m_j})$, ${k_j} - {m_j} \geq 0$. 由式(G6)和式(G7), 可得:

    $$\begin{split} &{\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1)\tilde {\boldsymbol{x}}_j^{\rm{T}}(t{\rm{ + }}{k_j}|t{\rm{ + }}{k_j} - 1)] =\\ &\;\;\;\;\;\;\;\;{\boldsymbol{P}}_{ij}^x(t|t - 1){({{\boldsymbol{A}}^{{k_j}}})^{\rm{T}}} - \sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j} - {m_j})} \;\times \\ &\;\;\;\;\;\;\;\;{\boldsymbol{L}}_j^{\rm{T}}(t + {k_j} - {m_j}){({{\boldsymbol{A}}^{{m_j} - 1}})^{\rm{T}}} \end{split}\tag{G8}$$
    $$\begin{split} &{\rm{E}}[{\tilde {\boldsymbol{x}}_i}(t|t - 1){(\tilde {\boldsymbol{x}}_j^r(t{\rm{ + }}{k_j}|t{\rm{ + }}{k_j} - 1))^{\rm{T}}}] = \\ &\;\;\;\;\;\;\;\;{\boldsymbol{P}}_{ij}^{x{x^r}}(t|t - 1){({\boldsymbol{\Phi }}_j^{{k_j}})^{\rm{T}}} - \sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{x\epsilon }(t + {k_j} - {m_j})}\; \times \\ &\;\;\;\;\;\;\;\;{({\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}))^{\rm{T}}}{({\boldsymbol{\Phi }}_j^{{m_j} - 1})^{\rm{T}}} \end{split}\tag{G9}$$

    将式(G8)和式(G9)代入式(G4), 整理得式(40). 类似地, 有:

    $$\begin{split} &{\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1){({\tilde {\boldsymbol{x}}_j}(t + {k_j}{\rm{|}}t + {k_j} - 1))^{\rm{T}}}] = \\ &\;\;\;\;\;\;\;\;\;{({\boldsymbol{P}}_{ji}^{x{x^r}}(t|t - 1))^{\rm{T}}}{({{\boldsymbol{A}}^{{k_j}}})^{\rm{T}}} \;- \\ &\;\;\;\;\;\;\;\;\;\sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j} - {m_j})} \; \times \\ &\;\;\;\;\;\;\;\;\;{\boldsymbol{L}}_j^{\rm{T}}(t + {k_j} - {m_j}){({{\boldsymbol{A}}^{{m_j} - 1}})^{\rm{T}}} \end{split}\tag{G10}$$
    $$\begin{split} &{\rm{E}}[\tilde {\boldsymbol{x}}_i^r(t|t - 1){(\tilde {\boldsymbol{x}}_j^r(t + {k_j}{\rm{|}}t + {k_j} - 1))^{\rm{T}}}] = \\ &\;\;\;\;\;\;\;{\boldsymbol{P}}_{ij}^{{x^r}}(t|t - 1){({\boldsymbol{\Phi }}_j^{{k_j}})^{\rm{T}}} - \sum\limits_{{m_j} = 1}^{{k_j}} {{\boldsymbol{P}}_{ij}^{{x^r}\epsilon }(t + {k_j} - {m_j})} \; \times \\ &\;\;\;\;\;\;\;{({\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j}))^{\rm{T}}}{{\rm{(}}{\boldsymbol{\Phi }}_j^{{m_j} - 1})^{\rm{T}}}\\[-5pt] \end{split}\tag{G11}$$

    由式(G10)和式(G11), 可得式(41).

    最后, 求解新息互协方差矩阵${\boldsymbol{Q}}_{ij}^\epsilon (t + {k_i},t + {k_j}) = {\rm{E}}[{{\boldsymbol{\epsilon }}_i}(t + {k_i}){\boldsymbol{\epsilon }}_j^{\rm{T}}(t + {k_j})],{k_i} < {k_j}$. 将式(G6)和式(G7)代入式(G3), 整理得新息:

    $$\begin{split} {{\boldsymbol{\epsilon }}_j}(t + {k_j}) =\;& {{\boldsymbol{C}}_j}{{\boldsymbol{A}}^{{k_j}}}{\tilde {\boldsymbol{x}}_j}(t|t - 1) + {\tilde {\boldsymbol{C}}_j}{\boldsymbol{\Phi }}_j^{{k_j}}\tilde {\boldsymbol{x}}_j^r(t|t - 1) \;- \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {[{{\boldsymbol{C}}_j}{{\boldsymbol{A}}^{{m_j} - 1}}{{\boldsymbol{L}}_j}(t + {k_j} - {m_j})} \; + \\ &{\tilde {\boldsymbol{C}}_j}{\boldsymbol{\Phi }}_j^{{m_j} - 1}{\boldsymbol{L}}_j^{{x^r}}(t + {k_j} - {m_j})]{{\boldsymbol{\epsilon }}_j}(t + {k_j} - {m_j})\; + \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {{{\boldsymbol{C}}_j}{\boldsymbol{A}}_{}^{{m_j} - 1}} {\boldsymbol{\Gamma w}}(t + {k_j} - {m_j}) \;+ \\ &\sum\limits_{{m_j} = 1}^{{k_j}} {{{\tilde {\boldsymbol{C}}}_j}{\boldsymbol{\Phi }}_j^{{m_j} - 1}} {{\boldsymbol{\xi }}_j}(t + {k_j} - {m_j}) + {{\boldsymbol{v}}_j}(t + {k_j}) \end{split}\tag{G12}$$

    由注4得${{\boldsymbol{\epsilon }}_i}(t + {k_i}) \bot {{\boldsymbol{\xi }}_j}(t + {k_j} - {m_j})$. 由注3得${{\boldsymbol{\epsilon }}_i}(t + {k_i}) \bot {{\boldsymbol{v}}_j}(t + {k_j})$, ${k_i} < {k_j}$. 当${k_i} < {k_j} - {m_j}$ 时, ${{\boldsymbol{\epsilon }}_i}(t + {k_i}) \bot {\boldsymbol{w}}(t + {k_j} \;- {m_j})$; 当${k_i} \geq {k_j} - {m_j}$时, ${{\boldsymbol{\epsilon }}_i}(t + {k_i})$与${\boldsymbol{w}}(t + {k_j} - {m_j})$相关. 令${\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i},t + {k_j} - {m_j}{\rm{)}} = {\rm{E}}[{{\boldsymbol{\epsilon }}_i}(t + {k_i}){{\boldsymbol{w}}^{\rm{T}}}(t + {k_j} - {m_j})]$. 利用此定义和(G12), 可得不同时刻的新息互协方差矩阵(42).

    下面推导${\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i},t + {k_j} - {m_j}{\rm{)}}$. 利用${\boldsymbol{w}}(t + {k_j} \;- {m_j}) \bot {\tilde {\boldsymbol{x}}_i}(t|t - 1)$、${\boldsymbol{w}}(t + {k_j} - {m_j}) \bot \tilde {\boldsymbol{x}}_i^r(t|t - 1)$ 和 ${\boldsymbol{w}}(t + {k_j} \;- {m_j}) \bot {{\boldsymbol{\xi }}_i}(t + {k_i} - {m_i})$, 有:

    $$\begin{split} &{\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i},t + {k_j} - {m_j}{\rm{) = }}\\ &\;\;\;\;\;\;\;\;\;\; - \sum\limits_{{m_i} = 1}^{{k_i}} {[{{\boldsymbol{C}}_i}{{\boldsymbol{A}}^{{m_i} - 1}}{{\boldsymbol{L}}_i}(t + {k_i} - {m_i})}\; + \\ &\;\;\;\;\;\;\;\;\;\;{\tilde {\boldsymbol{C}}_i}{\boldsymbol{\Phi }}_i^{{m_i} - 1}{\boldsymbol{L}}_i^{{x^r}}(t + {k_i} - {m_i})]\; \times \\ &\;\;\;\;\;\;\;\;\;\;{\boldsymbol{\Omega }}_i^{\epsilon w}{\rm{(}}t + {k_i} - {m_i},t + {k_j} - {m_j}{\rm{) \;+ }}\\ &\;\;\;\;\;\;\;\;\;\;\sum\limits_{{m_i} = 1}^{{k_i}} {{{\boldsymbol{C}}_i}{\boldsymbol{A}}_{}^{{m_i} - 1}} {\boldsymbol{\Gamma }}\; \times \\ &\;\;\;\;\;\;\;\;\;\;{\rm{E}}[{\boldsymbol{w}}(t + {k_i} - {m_i}){{\boldsymbol{w}}^{\rm{T}}}(t + {k_j} - {m_j})]\; + \\ &\;\;\;\;\;\;\;\;\;\;{\rm{E}}[{{\boldsymbol{v}}_i}(t + {k_i}){{\boldsymbol{w}}^{\rm{T}}}(t + {k_j} - {m_j})] \end{split}\tag{G13}$$

    由式(G13)和假设1, 可得式(43).

    证明. 因为$\rho ({\boldsymbol{A}}) < 1$和 $\left| {{h_i}} \right| < 1$, 由${\bar {\boldsymbol{\Phi }}_i} ={\rm{diag}}\{{\boldsymbol{A}},{h_i}\}$, 则有$\rho ({\bar {\boldsymbol{\Phi }}_i}) < 1$. 因此$({\bar {\boldsymbol{\Phi }}_i},{\bar {\boldsymbol{C}}_i})$ 为可检对. 由于$({\bar {\boldsymbol{\Phi }}_i} - \bar {\boldsymbol{\Gamma }}{{\boldsymbol{S}}_i}{({\boldsymbol{Q}}_i^v)^{ - 1}}{\bar {\boldsymbol{C}}_i}, \bar {\boldsymbol{\Gamma }}{\bar {\boldsymbol{Q}}_a})$是可稳对, 则由经典卡尔曼滤波理论[1]可知稳态一步预报器式(55)存在, 且是渐近稳定的.

  • 图  1  分布式矩阵加权融合滤波的位置跟踪性能

    Fig.  1  The position tracking performance of distributed fusion filter weighted by matrices

    图  2  RMSEs比较图

    Fig.  2  Comparisons of RMSEs

    图  3  LF、SF、DF、MF、CF的误差方差

    Fig.  3  Variances of LFs, SF, DF, MF, and CF

    图  4  矩阵加权融合预报器、滤波器和平滑器的估计误差方差

    Fig.  4  Variances of fusion predictor, filter and smoother weighted by matrices

  • [1] 邓自立, 王欣, 高媛. 建模与估计(第二版). 北京: 科学出版社, 2016.

    Deng Zi-Li, Wang Xin, Gao Yuan. Modeling and Estimation (Second Edition). Beijing: The Science Press, 2016.
    [2] 徐宁寿. 随机信号估计与系统控制. 北京: 北京工业大学出版社, 2001.

    Xu Ning-Shou. Random Signal Estimation and System Control. Beijing: Beijing University of Technology Press, 2001.
    [3] 祁波, 孙书利. 带未知通信干扰和丢包补偿的多传感器网络化不确定系统的分布式融合滤波. 自动化学报, 2018, 44(6): 1107--1114.

    Qi Bo, Sun Shu-Li. Distributed fusion filtering for multi-sensor networked uncertain systems with unknown communication disturbances and compensations of packet dropouts. Acta Automatica Sinica, 2018, 44(6): 1107--1114.
    [4] 李娜, 马静, 孙书利. 带多丢包和滞后随机不确定系统的最优线性估计. 自动化学报, 2015, 41(3): 611--619.

    Li Na, Ma Jing, Sun Shu-Li. Optimal linear estimation for stochastic uncertain systems with multiple packet dropouts and delays. Acta Automatica Sinica, 2015, 41(3): 611--619.
    [5] Kim K J, Iltis R A. Joint detection and channel estimation algorithms for QS-CDMA signals over time-varying channels. IEEE Transaction on Communication, 2002, 50(5): 845--855. doi: 10.1109/TCOMM.2002.1006565
    [6] Iltis R A. Joint estimation of PN code delay and multipath using the extended Kalman filter. IEEE Transactions on Automatic Control, 1990, 38(10): 1677--1685.
    [7] Chow B S, Birkemeier W P. A new structure of recursive estimator. IEEE Transaction on Automatic Control, 1989, 34(5): 568--572.
    [8] Chow B S, Birkemeier W P. A new recursive filter for systems with multiplicative noise. IEEE Transactions on Information Theory, 1990, 36(6): 1430--1435. doi: 10.1109/18.59939
    [9] Zhang L, Zhang X. An optimal filtering algorithm for systems with multiplicative/additive noises. IEEE Signal Processing Letters, 2007, 14(7): 469--472. doi: 10.1109/LSP.2006.891331
    [10] Yu X, Jin G, Li J. Target tracking algorithm for system with Gaussian/non-Gaussian multiplicative noise. IEEE Transactions on Vehicular Technology, 2020, 69(1): 90--100. doi: 10.1109/TVT.2019.2952368
    [11] Liu W. Optimal filtering for discrete-time linear systems with time-correlated multiplicative measurement noises. IEEE Transactions on Automatic Control, 2016, 61(7): 1972--1978. doi: 10.1109/TAC.2015.2480238
    [12] Wang S, Wang Z, Dong H, F E Alsaadi. Recursive state estimation for linear systems with lossy measurements under time-correlated multiplicative noises. Journal of the Franklin Institute, 2020, 57: 1887--1908.
    [13] Liu W, Shi P. Optimal linear filtering for networked control systems with time-correlated fading channels, Automatica, 2019, 101: 345--353. doi: 10.1016/j.automatica.2018.11.042
    [14] Li W, Jia Y, Du J. Tobit Kalman filter with time-correlated multiplicative measurement noise. IET Control Theory and Applications, 2017, 11(1): 122--128. doi: 10.1049/iet-cta.2016.0624
    [15] Geng H, Wang Z, Liang Y, Cheng Y, Alsaadi F E. Tobit Kalman filter with time-correlated multiplicative sensor noises under redundant channel transmission. IEEE Sensors Journal, 2017, 17(24): 8367--8377. doi: 10.1109/JSEN.2017.2766077
    [16] Lin H, Sun S. Globally optimal sequential and distributed fusion state estimation for multi-sensor systems with cross-correlated noises. Automatica, 2019, 101: 128--137. doi: 10.1016/j.automatica.2018.11.043
    [17] Carlson N A. Federated square root filter for decentralized parallel processes. IEEE Transactons on Aerospace and Electronic Systems, 1990, 26(3): 517--525. doi: 10.1109/7.106130
    [18] Kim K H. Development of track to track fusion algorithm. In: Proceedings of the American Control Conference. Maryland, USA: IEEE, 1994. 1037−1041
    [19] Sun S, Deng Z. Multi-sensor optimal information fusion Kalman filter. Automatica, 2004, 40(6): 1017--1023. doi: 10.1016/j.automatica.2004.01.014
    [20] 从金亮, 李银伢, 戚国庆, 盛安冬. 快速协方差交叉融合算法及应用. 自动化学报, 2020, 46(7): 1433--1444.

    Cong Jin-Liang, Li Yin-Ya, Qi Guo-Qing, Sheng An-Dong. A fast covariance intersection fusion algorithm and its application. Acta Automatica Sinica, 2020, 46(7): 1433--1444.
    [21] 王雪梅, 刘文强, 邓自立. 不确定系统改进的鲁棒协方差交叉融合稳态Kalman预报器, 自动化学报, 2016, 42(8): 1198--1206.

    Wang Xue-Mei, Liu Wen-Qiang, Deng Zi-Li. Modified robust covariance intersection fusion steady-state Kalman predictor for uncertain systems. Acta Automatica Sinica, 2016, 42(8): 1198--1206.
    [22] Sun S. Distributed optimal linear fusion estimators. Information fusion, 2020, 63: 56--73. doi: 10.1016/j.inffus.2020.05.006
    [23] Sun S. Distributed optimal linear fusion predictors and filters for systems with random parameter matrices and correlated noise. IEEE Transactons on Signal Processing, 2020, 68: 1064--1074. doi: 10.1109/TSP.2020.2967180
    [24] Ma J, Sun S. Distributed fusion filter for networked stochastic uncertain systems with transmission delays and packet dropouts. Signal Processing, 2017, 130: 268--278. doi: 10.1016/j.sigpro.2016.07.004
    [25] Tian T, Sun S, Lin H. Distributed fusion filter for multi-sensor systems with finite-step correlated noises. Information Fusion, 2019, 46: 128--140. doi: 10.1016/j.inffus.2018.05.002
    [26] Wang X, Liu W, Deng Z. Robust weighted fusion Kalman estimators for multi-model multisensor systems with uncertain-variance multiplicative and linearly correlated additive white noises. Signal Processing, 2017, 137: 339--355. doi: 10.1016/j.sigpro.2017.02.015
    [27] Yang X, Song Z, Ma J. Distributed fusion filter for systems with time-correlated multi-plicative noise. In: Proceedings of the 39th Chinese Control Conference. Shenyang, China: IEEE, 2020. 2741−2745
    [28] 孙书利, 马静. 随机奇异系统的分布式降阶最优融合Kalman滤波器. 自动化学报, 2006, 32(20): 286--290.

    Sun Shu-Li, Ma Jing, Distributed reduced-order optimal fusion Kalman filters for stochastic singular system. Acta Automatica Sinica, 2006, 32(20): 286--290.
    [29] Sun Shu-Li, Wang Guang-Hui, Modeling and estimation for networked systems with multiple random transmission delays and packet losses. Systems & Control Letters, 2014, 73, 6--16.
    [30] 段广全, 孙书利. 带未知模型参数和衰减观测率系统自校正分布式融合估计. 自动化学报, 2021, 47(2): 423--431.

    DUAN Guang-Quan, SUN Shu-Li. Self-tuning distributed fusion estimation for systems with unknown model parameters and fading measurement rates. Acta Automatica Sinica, 2021, 47(2): 423--431.
  • 期刊类型引用(3)

    1. Yongpeng CUI,Xiaojun SUN. Multi-Sensor Fusion Adaptive Estimation for Nonlinear Under-observed System with Multiplicative Noise. Chinese Journal of Electronics. 2024(01): 282-292 . 必应学术
    2. 崔永鹏,孙小君,张扬. 带乘性噪声的欠观测系统无迹增量Kalman融合估计. 黑龙江大学工程学报(中英俄文). 2024(02): 66-74 . 百度学术
    3. 王锦. 面向分布式多传感器的FOA大数据融合算法研究. 北部湾大学学报. 2024(04): 60-67 . 百度学术

    其他类型引用(10)

  • 加载中
图(4)
计量
  • 文章访问数:  1500
  • HTML全文浏览量:  487
  • PDF下载量:  215
  • 被引次数: 13
出版历程
  • 收稿日期:  2021-02-18
  • 网络出版日期:  2021-08-16
  • 刊出日期:  2023-08-21

目录

/

返回文章
返回