-
摘要: 在统计流形空间中, 从信息几何角度考虑非线性状态后验分布近似的实质是后验分布与相应参数化变分分布之间的Kullback-Leibler (KL)散度最小化问题, 同时也可以转化为变分置信下界的最大化问题. 为了提升非线性系统状态估计的精度, 在高斯系统假设条件下结合变分贝叶斯(Variational Bayes, VB)推断和Fisher信息矩阵推导出置信下界的自然梯度, 并通过分析其信息几何意义, 阐述在统计流形空间中置信下界沿其方向不断迭代增大, 实现变分分布与后验分布的“紧密”近似; 在此基础上, 以状态估计及其误差协方差作为变分超参数, 结合最优估计理论给出一种基于自然梯度的非线性变分贝叶斯滤波算法; 最后, 通过天基光学传感器量测条件下近地轨道卫星跟踪定轨和纯角度被动传感器量测条件下运动目标跟踪仿真实验验证, 与对比算法相比, 所提算法具有更高的精度.
-
关键词:
- 非线性滤波 /
- 信息几何 /
- 变分贝叶斯推断 /
- 自然梯度 /
- Fisher信息矩阵
Abstract: In statistical manifold space, the essence of nonlinear state posterior distribution approximation from the perspective of information geometry is minimizing Kullback-Leibler (KL) divergence between posterior distribution and the corresponding approximated distribution; Meanwhile, it is equivalent to maximizing evidence low bound. Aiming at the problem of improving the estimation accuracy of nonlinear system state, the natural gradient of evidence lower bound is derived under Gaussian system assumption by combining with Fisher information matrix and variational Bayesian (VB) inference, which produces a faster movement direction to the posterior distribution, and realizing a close approximation between variational distribution and the posterior. On this basis, a variational Bayesian Kalman filtering algorithm using natural gradient is proposed for updating the variational hyperparameters of state estimation and the associated error covariance. Simulations in low earth orbit target tracking system with space-based optical sensors and bearing-only target tracking system are presented verifying that the proposed algorithm has higher accuracy than the comparison algorithms. -
在过去的十几年内, 随着传感器网络[1]、智能电网[2]、机器人集群[3]等多智能体系统的广泛应用, 各种分布式优化算法[4−5]应运而生. 其中, 最优一致性算法[6−7]因其能够同时满足一致性和最优性的特点而备受关注. 具体而言, 在最优一致性问题的场景下, 所有智能体都拥有一个局部代价函数, 通过使用自身信息以及与邻居的交互信息来设计一类分布式算法, 使每个智能体的状态都能够收敛到全局代价函数的最优解.
在早期的研究中, 针对最优一致性问题, 无论是离散时间算法[8−9]还是连续时间算法[10−11], 大多基于线性系统. 例如, 文献[8−9]针对一阶系统分别提出基于消失步长的次梯度下降算法和基于固定步长的EXTRA算法. 文献[10]针对有向图下的一阶系统提出一种基于分布式比例积分(Proportional-integral, PI)的最优一致性算法; 文献[11]在文献[10]的基础上进一步优化, 消除了对辅助变量通信的需求. 然而, 大多数实际系统(如无人机集群)的动力学模型都是非线性的. 因此, 近年来, 越来越多的学者致力于非线性多智能体系统[12−18]的最优一致性研究. 例如, 文献[12]通过内模法实现了受扰非线性多智能体系统的最优一致性. 进一步地, 文献[13]提出一种自适应内模法, 解决了一类最小相位多智能体系统的最优一致性问题. 文献[14]提出一类小增益方法来处理非线性多智能体系统. 在嵌入式的框架下, 文献[15−17]分别研究了具有静态不确定性、动态不确定性以及复杂动力学模型的多智能体系统的最优一致性. 此外, 文献[18]将内模法和嵌入式框架相结合, 解决了二阶非线性多智能体系统的最优一致性问题. 需要注意的是, 上述研究均依赖额外的参考信息才能实现非线性多智能体系统的最优一致性. 例如, 在文献[15−18]中, 由于采用了嵌入式框架, 必须依赖由一阶参考系统生成的最优信号, 才能达到最优一致性. 然而, 在某些只能使用相对状态信息的应用场景中, 获取额外的参考信息并不现实, 如仅配备相对位置传感器的多无人机系统. 因此, 对于非线性多智能体系统而言, 设计一类仅使用相对状态信息的最优一致性算法仍是一项具有挑战性的任务.
最近, 通过采用PI调节的方式, 文献[19]解决了无向图下带有状态时延和扰动的纯反馈多智能体系统的最优一致性问题. 文献[19]引入一类分布式变量, 将最优一致性问题转化为PI调节问题, 并在仅使用相对状态信息的情况下, 实现了近似最优一致性. 值得指出的是, 文献[19]所提出的分布式变量是基于无向图的, 但在实际应用中, 由于通信带宽的限制, 网络链接通常为单链路, 这表明研究适用于权重非平衡有向图的分布式变量具有重要的实际意义. 此外, 死区输入是实际控制系统中最常见的非线性源之一, 如单连杆机械臂、液压伺服阀、压电转换器等系统. 死区非线性的存在常与闭环系统的性能下降、严重振荡甚至更糟的不稳定现象直接相关. 因此, 设计合适的死区输入非线性补偿, 已成为解决上述控制系统最优一致性问题的当务之急.
在上述讨论的启发下, 本文研究了有向图下一类带有死区输入非线性和有界扰动的严格反馈多智能体系统的最优一致性问题. 首先, 我们提出一种新的分布式PI变量, 将最优一致性问题转化为PI调节问题. 然后, 针对权重非平衡有向图下的严格反馈多智能体系统, 结合所提出的分布式PI变量和预设性能控制, 设计一类基于PI调节的近似最优一致性算法. 本文的主要创新点总结如下:
1)提出一种新的分布式PI变量, 将最优一致性问题转化为PI调节问题, 使得经典的控制技术能够通过调节PI变量的方式来处理更加复杂的非线性多智能体系统. 相较于文献[19]中基于无向图的分布式变量, 本文提出的分布式PI变量形式更加简洁, 且适用于更一般的权重非平衡有向图. 此外, 与采用嵌入式框架的文献[18]相比, 本文无需依赖额外的参考信息, 仅使用相对状态信息即可实现最优一致性.
2)研究了权重非平衡有向图下带有死区输入非线性和有界扰动的严格反馈多智能体系统的最优一致性问题, 并基于PI调节提出一类新的近似最优一致性算法. 与文献[19]中的理想控制输入相比, 本文对控制输入施加了死区约束, 所考虑的情况更加复杂.
本文的剩余部分安排如下: 第1节介绍基础知识和研究问题. 第2节引入新的分布式PI变量并构建基于PI调节的重要引理. 第3节设计基于PI调节的严格反馈多智能体系统的近似最优一致性算法. 第4节进行仿真实验. 第5节对全文工作进行总结与展望.
符号说明: 令$ {\bf{R}} $为实数集; $ {\bf{R}}_{\geq 0} $为非负实数集; $ {\bf{R}}^N $和$ {\bf{R}}^n $分别为$ N $维和$ n $维列向量空间; $ {\bf{R}}^{N \times n} $为$ N $行$ n $列的矩阵空间; $ 0_N $, $ 1_N \in {\bf{R}}^N $分别表示元素全为$ 0 $和元素全为1的列向量; $ I_n $表示维度为$ n $的单位矩阵; $ \lambda_i(A) $表示矩阵$ A \in {\bf{R}}^{N \times N} $的第$ i $个特征值, $ i=1,\; 2,\; \cdots ,\; N $; 由元素$ d_i $构成的对角阵和列向量分别记为$ \mathrm{diag}\{d_1,\; d_2,\; \cdots,\; d_N\} $, $ \mathrm{col}(d_1,\; d_2,\; \cdots, d_N) $; $ {\cal{L}}_{\infty} $表示一致有界集; $ \|\cdot\| $表示欧几里德范数.
1. 基础知识和问题描述
1.1 图论
定义$ {\cal{G}} = ({\cal{V}},\; {\cal{E}},\; {\cal{A}}) $为包含$ N $个智能体的有向图[20−22], 其中, $ {\cal{V}}=1,\; 2,\; \cdots,\; N $表示所有智能体的集合; $ {\cal{E}} \subseteq {\cal{V}} \times {\cal{V}} $表示智能体之间的有向边集. 边$ (j,\; i) $表示不同智能体$ j,\; i $之间存在有向链路, 且智能体$ i $能收到智能体$ j $发送的信息. 换言之, 智能体$ j $是智能体$ i $的入邻居. $ {\cal{A}} = [a_{ij}] \in {\bf{R}}^{N \times N} $表示邻接矩阵, 当且仅当$ (j,\; i)\in{\cal{E}} $时, $ a_{ij} > 0 $; 其他情况下, $ a_{ij} $均为$ 0 $. 智能体$ i $的入度记为$ d_i = \mathop{\sum}_{j = 1}^{N}a_{ij} $. $ {\cal{L}} = {\cal{D}} - {\cal{A}} $为有向图$ {\cal{G}} $的Laplacian矩阵, 其中, $ {\cal{D}} = \operatorname{diag}\{d_1,\; d_2,\; \cdots,\; d_N\} $. 此外, 如果一个有向图中任意两个不同智能体之间都存在一条有向路径, 则该图是强连通的.
引理 1. 假设有向图$ {\cal{G}} $是强连通的[23−24], 则具有如下性质:
1)存在一个属于Laplacian矩阵零特征值的左特征向量$ r=(r_1,\; r_2,\; \cdots,\; r_N)^{\mathrm{T}} $, $ r_i > 0 $, 使得等式$ r^{\operatorname{T}}{\cal{L}} = 0_N^{\operatorname{T}} $和 $ \mathrm{\sum_{ }^{ }}_{i=1}^Nr_i=1 $成立;
2)存在一个半正定矩阵$ \tilde{{\cal{L}}} = ({R{\cal{L}} + {\cal{L}}^{\operatorname{T}}R})/{2} $, 其中, $ R = \mathrm{diag}\{r_1,\; r_2,\; \cdots,\; r_N\} $, 且$ \tilde{{\cal{L}}} $的特征值可以按序排列为: $ 0 = \lambda_1(\tilde{{\cal{L}}}) < \lambda_2(\tilde{{\cal{L}}}) \leq \lambda_3(\tilde{{\cal{L}}}) \leq \cdots \leq \lambda_N(\tilde{{\cal{L}}}) $;
3)存在一个对角线元素均大于$ 0 $的半正定矩阵 $ \mathrm{exp}({-{\cal{L}}t}) $, $ t > 0 $, 且当$ t \to \infty $时, $ \mathrm{exp}({-{\cal{L}}t}) \to 1_Nr^{\operatorname{T}} $.
1.2 近似最优一致性问题
近似最优一致性问题旨在为包含$ N $个智能体的多智能体系统设计一个恰当的控制算法, 使得:
$$ \begin{align} \mathop{\lim}\limits_{t \to \infty}\| s_i(t)-s^{*} \|=\kappa\epsilon ,\; \quad \forall i\in{\cal{V}} \end{align} $$ (1) 其中, $ s_i(t): {\bf{R}}_{\geq 0} \to {\bf{R}}^n $表示每个智能体的状态变量;$ s^* \in {\bf{R}}^n $表示全局代价函数$ f(s(t)) = \mathop{\sum}_{i = 1}^{N}f_i(s_i(t)) \in {\bf{R}} $的最小值点, 即最优状态, $ f_i(s_i(t)): {\bf{R}}^n \to {\bf{R}} $表示每个智能体的局部代价函数; $\epsilon>0 $为设计参数, 其取值可以为任意小; $\kappa $表示大于零的常数且独立于$\epsilon $.
2. 基于PI调节的最优一致性
2.1 分布式PI变量
本节提出一种新的分布式PI变量$ q_i(t): {\bf{R}}_{\geq 0} \to {\bf{R}}^n $, $ i \in {\cal{V}} $, 将最优一致性问题转化为调节问题, 其定义如下:
$$ \begin{align} q_{i}(t) = \;&s_{i}(t) + a\int_{0}^{t}\hat{r}_{i}^{-1}(\tau)(\nabla{f}_i(s_i(\tau)) \otimes {1}_N){\mathrm{d}}\tau\; +\\ &b\int_{0}^{t}\sum_{j = 1}^{N}a_{ij}(s_{i}(\tau)-s_{j}(\tau)){\mathrm{d}}\tau + \int_{0}^{t}v_i(\tau){\mathrm{d}}\tau \end{align} $$ (2) 并且
$$ \begin{align} \dot{v}_i(t) = ab\mathop{\sum}_{j = 1}^{N}a_{ij}(s_{i}(t)-s_{j}(t)) \end{align} $$ (3) $$ \begin{align} \dot{z}_i(t) = -\mathop{\sum}_{j = 1}^{N}a_{ij}(z_i(t) - z_j(t)) \end{align} $$ (4) 其中, $ a,\; b > 0 $为调节参数; $ z_i(t) \in {\bf{R}}^N $; $ \hat{r}_i(t)=z_i^i(t) $表示$ z_i(t) $的第 $ i $ 个元素. 此外, $ s_i(0) \in {\bf{R}}^n $; $ v_i(0) = 0_n $; $ \hat{r}_i(0) = z_i^i(0) = 1 $, $ z_i(t) $的其他元素初值均为$ 0 $.
注 1. 由文献[23]的定理1可知, 在权重非平衡有向图下, 式(2)右侧的微分形式能够渐近收敛至最优状态$ s^{*} $. 此外, 从式(2)的结构可以看出, 所提出的分布式变量是由状态变量与状态变量的积分项相加而成, 故命名为PI变量. 因此, 当时间趋于无穷大时, 若PI变量收敛至零, 则系统的所有状态变量将收敛至最优状态$ s^{*} $, 即实现最优一致性.
注 2. 分布式PI变量$q_i(t) $的设计灵感来源于PI一致性误差变量[25]. 然而, 与文献[25]研究的一致性不同, 本文提出的分布式PI变量解决的是最优一致性问题. 此外, 相较于文献[19]中的分布式变量, 式(2)不包含二重积分, 形式更加简洁, 并且通过构造Laplacian矩阵的零左特征向量估计器, 使得所提出的分布式PI变量能够在不依赖图的全局信息的情况下, 实现权重非平衡有向图下的最优一致性.
2.2 收敛性分析
本节将对多智能体系统的收敛性进行分析证明. 在此背景下, 多智能体系统的网络拓扑和代价函数的假设如下:
假设 1. 有向图$ {\cal{G}} $是强连通的.
假设 2. 每个智能体的局部代价函数$ f_i: {\bf{R}}^n \to {\bf{R}} $, $ i \in {\cal{V}} $是连续可微且强凸的, 即存在$ c_i > 0 $, 使得$ (X \,-\, Y)^{\operatorname{T}}[\nabla f_i(X) - \nabla f_i(Y)] \geq c_i\|X\, -\, Y\|^2 $, 其中, $ X,\; Y \in {\bf{R}}^n $. 此外, 每个函数$ f_i $都具有Lipschitz梯度, 即存在$ C_i > 0 $, 使得$ \|\nabla f_i(X) - \nabla f_i(Y)\| \leq C_i\|X - Y\| $.
为了进一步明确分布式PI变量$ q_i(t) $与多智能体系统最优一致性之间的关系, 本文提出以下定理:
定理 1. 考虑一个包含$ N $个智能体的多智能体系统, 在满足假设1和假设2的前提下, 当时间趋于无穷大时, 若所有的PI变量满足$\|q_i(t)\|\leq \epsilon,\;i\in \cal{V} $, 则该系统能够实现如式(1)所述的近似最优一致性. 此外, 式(2) ~ (4)中的所有信号有界. 换言之, $ s_i(t), \; v_i(t),\; \dot{v}_i(t),\; z_i(t),\; \dot{z}_i(t) \in {\cal{L}}_{\infty},\; i \in {\cal{V}} $. 进一步地, 如果PI变量趋向于零, 即当$ t \to \infty $时, 若$ q_i(t) \to 0_n $, 则$ \mathop{\lim}_{t \to \infty}s_i(t) = s^*,\; i \in {\cal{V}} $.
证明. 将式(2) ~ (4)中的变量统一改写为向量和矩阵的形式, 具体如下:
$$\left\{ \begin{aligned} &q(t) = \operatorname{col}(q_1(t),\; q_2(t),\; \cdots,\; q_N(t))\\ &\notag s(t) = \operatorname{col}(s_1(t),\; s_2(t),\; \cdots,\; s_N(t))\\ &\notag v(t) =\operatorname{col}(v_1(t),\; v_2(t),\; \cdots,\; v_N(t))\\ &\notag z(t) =\operatorname{col}(z_1(t),\; z_2(t),\; \cdots,\; z_N(t))\\ &\notag \hat{R}(t) = \operatorname{diag}\{\hat{r}_1(t),\; \hat{r}_2(t),\; \cdots,\; \hat{r}_N(t)\} \end{aligned}\right. $$ 进一步地, 定义向量$ \nabla f(s(t))=\mathrm{col}(\nabla f_1(s_1(t))\ \otimes 1_N,\; \nabla f_2(s_2(t))\otimes1_N,\; \cdots,\; \nabla f_N(s_N(t))\otimes1_N) $, 则PI变量可改写为:
$$ \begin{split} q(t)=\; & s(t)+a\int_0^t(\hat{R}^{-1}(\tau)\otimes I_n)\nabla f(s(\tau)){\mathrm{d}}\tau\; + \\ &b\int_0^t({\cal{L}}\otimes I_n)s(\tau){\mathrm{d}}\tau+\int_0^tv(\tau){\mathrm{d}}\tau \end{split} $$ (5) 并且
$$ \begin{align} \dot{v}(t) = ab({\cal{L}}\otimes{I_n})s(t) \end{align} $$ (6) $$ \begin{align} \dot{z}(t) = -({\cal{L}}\otimes{I_n})z(t) \end{align} $$ (7) 首先, 为了简化PI变量$ q(t) $, 定义辅助变量$ \eta(t) $为:
$$ \begin{split} \eta(t) =\; &-a\int_{0}^{t}(\hat{R}^{-1}(\tau)\otimes{I_n})\nabla{f}(s(\tau)){\mathrm{d}}\tau\; - \\ &b\int_{0}^{t}({\cal{L}}\otimes{I_n})s(\tau){\mathrm{d}}\tau - \int_{0}^{t}v(\tau){\mathrm{d}}\tau \end{split} $$ (8) 结合式(5), 可知$ q(t) = s(t) - \eta(t) $. 为了便于表示, 下文将时间变量$ t $省略. 根据式(6)和(8), 可得:
$$ \begin{align} \begin{pmatrix} \dot{\tilde{\eta}}\\ \dot{\tilde{v}} \end{pmatrix} = &\begin{pmatrix} -a(\hat{R}^{-1}\otimes{I_n})\nabla{f}(s) - b({\cal{L}}\otimes{I_n})s - v \\ ab({\cal{L}}\otimes{I_n})s \end{pmatrix} \end{align} $$ (9) 然后, 定义最优变量$ \eta_e = {1}_N \otimes s^* $, $ v_e = -a(\hat{R}^{-1}\, \otimes I_n)\nabla f(\eta_e) $. 进一步地, 令误差变量$ \tilde{\eta} = \eta - \eta_e $, $ \tilde{v} = v - v_e $, 联立$ q = s - \eta $可知, $ s = q + \tilde{\eta} + \eta_e $. 将$ v = \tilde{v} + v_e $和$ s = q + \tilde{\eta} + \eta_e $代入式(9), 并结合$ ({\cal{L}}\; \otimes I_n)\eta_e=({\cal{L}}\otimes I_n)(1_N\otimes s^*)=0 $, 则可求得误差变量$ \tilde{\eta} $和$ \tilde{v} $的动力学方程为:
$$ \begin{split} \begin{pmatrix} \dot{\tilde{\eta}}\\ \dot{\tilde{v}} \end{pmatrix} =\;& \begin{pmatrix}-a(\hat{R}^{-1}\otimes{I}_n)\nabla{f}(q + \tilde{\eta} + \eta_e) \\ ab({\cal{L}}\otimes{I}_n)\tilde{\eta} \end{pmatrix} + \\&\begin{pmatrix}a(\hat{R}^{-1}\otimes{I}_n)\nabla{f}(\eta_e) - b({\cal{L}}\otimes{I}_n)\tilde{\eta} \\ ab({\cal{L}}\otimes{I}_n)q \end{pmatrix} + \\ &\begin{pmatrix} - \tilde{v} - b({\cal{L}}\otimes{I}_n)q \\ {0}_{Nn} \end{pmatrix}\\[-1pt] \end{split} $$ (10) 最后, 为了便于统一表示, 定义$ \xi = \operatorname{col}(\tilde{\eta},\; \tilde{v}) $, 并利用$ \nabla{f}(\eta) = \nabla{f}(\tilde{\eta} + \eta_e) $. 于是, 式(10)可改写为:
$$ \begin{split} \underbrace{\begin{pmatrix} \dot{\tilde{\eta}} \\ \dot{\tilde{v}}\end{pmatrix}}_{\dot{\xi}} =\;& \underbrace{\begin{pmatrix} -a(R^{-1}\otimes I_n)g - b({\cal{L}}\otimes I_n)\tilde{\eta} - \tilde{v} \\ ab({\cal{L}}\otimes I_n)\tilde{\eta} \end{pmatrix}}_{{\bf{f}}(t,\; \xi)} + \\ &\;\;\underbrace{\begin{pmatrix} a((R^{-1}-\hat{R}^{-1})\otimes I_n)k \\ ab({\cal{L}}\otimes{I}_n)q \end{pmatrix}}_{{\bf{g}}(t,\; \xi)} + \\ &\;\;\underbrace{\begin{pmatrix} -a(R^{-1}\otimes I_n)h - b({\cal{L}}\otimes I_n)q \\{0}_{Nn} \end{pmatrix}}_{{\bf{g}}(t,\; \xi)} \\[-1pt]\end{split} $$ (11) 其中, $ g \,= \,\nabla f(\eta)-\nabla f(\eta_e) $; $ k \,=\, \nabla{f}(q + \tilde{\eta} \,+ \,\eta_e)\;- \nabla{f}(\eta_e) $; $ h = \nabla{f}(q + \tilde{\eta} + \eta_e) - \nabla{f}(\tilde{\eta} + \eta_e) $; $ R $的定义见引理1. 因此, 本定理的证明将分为如下两步: 标称系统的收敛性证明和受扰系统的收敛性证明.
步骤 1. 标称系统$ \dot{\xi} = {\bf{f}}(t,\; \xi) $的收敛性.
由式(11)可知, 标称系统$ \dot{\xi} = {\bf{f}}(t,\; \xi) $可以表示为:
$$ \begin{align} \begin{pmatrix} \dot{\tilde{\eta}}\\ \dot{\tilde{v}} \end{pmatrix} = &\begin{pmatrix} -a(R^{-1}\otimes I_n)g - b({\cal{L}}\otimes I_n)\tilde{\eta} - \tilde{v} \\ ab({\cal{L}}\otimes I_n)\tilde{\eta} \end{pmatrix} \end{align} $$ (12) 设Lyapunov函数为:
$$ \begin{align} V(\tilde{\eta},\; \tilde{v}) = \frac{1}{2}\begin{pmatrix}\tilde{\eta}^{\operatorname{T}} \tilde{v}^{\operatorname{T}}\end{pmatrix}P \begin{pmatrix} \tilde{\eta}\\ \tilde{v} \end{pmatrix} \end{align} $$ (13) 其中,
$$ \begin{align} P = \begin{pmatrix} a(\phi+1) \quad &1 \\ 1 \quad &\displaystyle\frac{1}{a} \end{pmatrix}\otimes(R\otimes I_n) \end{align} $$ (14) $ \phi > 0 $是设计参数, 具体范围将在下文给出. 求$ V $关于时间的导数, 可得:
$$ \begin{split}\dot{V}=\; & -a^2(\phi+1)\tilde{\eta}^{\mathrm{T}}g-ab\phi\tilde{\eta}^{\mathrm{T}}(R{\cal{L}}\otimes I_n)\tilde{\eta}\ - \\ & a(\phi+1)\tilde{\eta}^{\mathrm{T}}(R\otimes I_n)\tilde{v}-a\tilde{v}^{\mathrm{T}}g-\tilde{v}^{\mathrm{T}}(R\otimes I_n)\tilde{v}\end{split} $$ (15) 由于$ -\; ab\phi\tilde{\eta}^{\mathrm{T}}(R{\cal{L}}\otimes I_n)\tilde{\eta}=-ab\phi\tilde{\eta}^{\mathrm{T}}(\tilde{\cal{L}}\otimes I_n)\tilde{\eta} $, 其中, $ \tilde{{\cal{L}}} $的定义见引理1, 所以有:
$$ \begin{split}\dot{V}=\; & -a^2(\phi+1)\tilde{\eta}^{\mathrm{T}}g-ab\phi\tilde{\eta}^{\mathrm{T}}(\tilde{\cal{L}}\otimes I_n)\tilde{\eta}\ - \\ & a(\phi+1)\tilde{\eta}^{\mathrm{T}}(R\otimes I_n)\tilde{v}-a\tilde{v}^{\mathrm{T}}g-\tilde{v}^{\mathrm{T}}(R\otimes I_n)\tilde{v}\end{split} $$ (16) 下面对$ \dot{V} $进行放缩处理. 首先, 根据假设2中局部代价函数的强凸性, 可得:
$$ \begin{align} -a^2(\phi + 1)\tilde{\eta}^{\operatorname{T}}g \leq -a^2(\phi + 1)c\|\tilde{\eta}\|^2 \end{align} $$ (17) 其中, $ c = \operatorname{min}\{c_1,\, c_2,\, \cdots,\, c_N\} > 0 $ 是最小的强凸常数.
其次, 使用Rayleigh quotient不等式, 可得:
$$ \begin{align} -ab\phi\tilde{\eta}^{\operatorname{T}}(\tilde{{\cal{L}}} \otimes I_n)\tilde{\eta} \leq -ab\phi\lambda_2(\tilde{{\cal{L}}})\|\tilde{\eta}\|^2 \end{align} $$ (18) $$ \begin{align} -\tilde{v}^{\operatorname{T}}(R \otimes I_n)\tilde{v} \leq -r_{\mathrm{min}}\|\tilde{v}\|^2 \end{align} $$ (19) 其中, $ \lambda_2(\tilde{{\cal{L}}}) $是矩阵$ \tilde{{\cal{L}}} $的第二小特征值; $ r_{\mathrm{min}} $是矩阵$ R $的最小特征值. 基于上述放缩结果, 使用Young不等式继续处理式(16)的剩余项, 则有:
$$ \begin{split} -a(\phi + 1)&\tilde{\eta}^{\operatorname{T}}(R \otimes I_n)\tilde{v} \leq \frac{r_{\mathrm{min}}}{3}\|\tilde{v}\|^2 \;+ \\ &\frac{3a^2(\phi+1)^2r_{\mathrm{max}}^2}{4r_{\mathrm{min}}}\|\tilde{\eta}\|^2 \end{split} $$ (20) $$ \begin{align} -a\tilde{v}^{\operatorname{T}}g \leq \frac{3a^2C^2}{4r_{\mathrm{min}}}\|\tilde{\eta}\|^2 + \frac{r_{\mathrm{min}}}{3}\|\tilde{v}\|^2 \end{align} $$ (21) 其中, $ r_{\mathrm{max}} $是矩阵$ R $的最大特征值; $ C = \operatorname{max}\{C_1,\; C_2,\; \cdots,\; C_N\} > 0 $是最大的Lipschitz常数.
然后, 将式(17) ~ (21)代入式(16)得:
$$ \begin{split} \dot{V}& \leq - \left(a^2(\phi + 1)c - \frac{3a^2C^2}{4r_{\mathrm{min}}}\right)\|\tilde{\eta}\|^2 \;- \\ &\left(ab\phi\lambda_2(\tilde{{\cal{L}}}) - \frac{3a^2(\phi+1)^2r_{\mathrm{max}}^2}{4r_{\mathrm{min}}}\right)\|\tilde{\eta}\|^2 - \frac{r_{\mathrm{min}}}{3}\|\tilde{v}\|^2 \end{split} $$ (22) 令
$$ \begin{align} \begin{cases} a^2(\phi + 1)c - \displaystyle\frac{3a^2C^2}{4r_{\mathrm{min}}} > 0 \\ ab\phi\lambda_2(\tilde{{\cal{L}}}) - \displaystyle\frac{3a^2(\phi + 1)^2r_{\mathrm{max}}^2}{4r_{\mathrm{min}}} > 0 \end{cases} \end{align} $$ (23) 解得:
$$ \begin{align} \begin{cases} \phi + 1 > \displaystyle\frac{3C^2}{4cr_{\mathrm{min}}} \\ \phi > \displaystyle\frac{3a(\phi + 1)^2r_{\mathrm{max}}^2}{4b\lambda_2(\tilde{{\cal{L}}})r_{\mathrm{min}}} \end{cases} \end{align} $$ (24) 显然, 当$\phi,\;a,\;b $满足上述条件时, $ \dot{V}\le0 $. 进一步地, $ \dot{V} $可表示为:
$$ \begin{align} \dot{V} \leq -\alpha\left(\|\tilde{\eta}\|^2 + \|\tilde{v}\|^2\right) \end{align} $$ (25) 其中, $ \alpha = \operatorname{min}\{\frac{r_{\mathrm{min}}}{3},\, a^2(\phi + 1)c - \frac{3a^2C^2}{4r_{\mathrm{min}}} + ab\phi\lambda_2(\tilde{{\cal{L}}}) \;- \frac{3a^2(\phi \;+\; 1)^2r_{\mathrm{max}}^2}{4r_{\mathrm{min}}}\} $.
最后, 由式(13)和(14)可得, $ V \;\leq\; \beta(\|\tilde{\eta}\|^2 \;+ \|\tilde{v}\|^2) $, 其中, $ \beta $是矩阵$ P $的最大特征值. 综上所述, 可推得$ \dot{V} \leq -({\alpha}/{\beta})\times V $. 因此, 标称系统$ \dot{\xi} = {\bf{f}}(t,\; \xi) $指数收敛于原点, 即当$ t \to \infty $时, $ \tilde{\eta} \to {0}_{Nn},\; \tilde{v} \to {0}_{Nn} $, 且收敛速率不低于$ \alpha/\beta $.
步骤 2. 受扰系统$ \dot{\xi} = {\bf{f}}(t,\; \xi) + {\bf{g}}(t,\; \xi) $的收敛性.
下面对扰动项$ {\bf{g}}(t,\; \xi) $进行放缩处理, 可得:
$$ \begin{split} {\bf{g}}(t,\; \xi) \leq\; & {\begin{pmatrix} aC\|(R^{-1}-\hat{R}^{-1})\otimes I_n\|(\|\tilde{\eta}\| + \|q\|) \\ ab\|{\cal{L}}\otimes I_n\|\|q\| \end{pmatrix}} \;+ \\ &{\begin{pmatrix} aC\|R^{-1}\otimes I_n\|\|q\| + b\|{\cal{L}}\otimes I_n\|\|q\| \\ {0}_{Nn} \end{pmatrix}} \end{split} $$ (26) 由引理1可知, $ \mathop{\lim}_{t \to \infty}\mathrm{exp}({-{\cal{L}}t}) = {1}_Nr^{\operatorname{T}} $, 所以有$ \mathop{\lim}_{t \to \infty}z(t) = \mathop{\lim}_{t \to \infty}\mathrm{exp}({-\;({\cal{L}}\;\otimes\; I_{N})\;t})z(0) = ({1}_{N}r^{\operatorname{T}} \otimes I_{N})z(0) = {1}_{N}\otimes r $. 因此, 当$ t \to \infty $时, $ \hat{R}^{-1} \to R^{-1} $.
此外, 当$ t \to \infty $时, 若$ q(t) $有界, 则$ {\bf{g}}(t,\; \xi) $有界. 由文献[26]的引理$ 9.4 $、推论$ 9.1 $以及引理$ 9.5 $的证明和结论可知, 变量$ \xi(t) $也有界, 即$ (\tilde{\eta}(t),\; \tilde{v}(t)) $有界. 相似地, 当$ t \to \infty $时, 若$ q(t) \to {0}_{Nn} $, 则$ {\bf{g}}(t,\; \xi) \to {0}_{2Nn} $, 那么受扰系统指数收敛于原点, 即当$ t \to \infty $时, $ (\tilde{\eta}(t),\; \tilde{v}(t)) \to {0}_{2Nn} $.
综上所述, 当时间$ t $趋于无穷大时, 如果$ q(t) $有界, 即$\mathop{\lim}_{t \to \infty}||q(t)|| \leq \epsilon $, 那么$ \tilde{\eta}(t) $和$ \tilde{v}(t) $也有界, 相应地, $ \mathop{\lim}_{t \to \infty}||\tilde{\eta}(t)|| \leq \kappa_1\epsilon,\ \kappa_1> 0$. 再根据$ \tilde{\eta}(t) $和$ \tilde{v}(t) $的定义, 可知$ \eta(t) $和$ v(t) $有界. 进一步地, 由$ q(t) = s(t) - \eta(t) $可知, $ s(t) $有界, 故可推得$ \mathop{\lim}_{t \to \infty}||s(t)\; - {1}_N \otimes s^*||$ $=\; \mathop{\lim}_{t \to \infty}||s(t) \,-\, \eta_e||$ $=\; \mathop{\lim}_{t \to \infty}||s(t)\; - \eta(t) + \tilde{\eta}(t)|| $ $= \mathop{\lim}_{t \to \infty}||q(t) + \tilde{\eta}(t)|| $ $\leq \mathop{\lim}_{t \to \infty}||q(t)|| + \mathop{\lim}_{t \to \infty}||\tilde{\eta}(t)||\leq \kappa\epsilon $, 其中, $\kappa = 1 + \kappa_1 $. 再结合式(6), 可推得$ \dot{v}(t) $有界. 此外, 从上述证明可知, 当$ t $趋于无穷大时, $ \hat{R}^{-1} $趋于$ R^{-1} $, 说明$ z(t) $有界, 从而根据式(7), 可得$ \dot{z}(t) $有界. 因此, 式(2) ~ (4)中的所有信号有界. 再者, 当时间$ t $趋于无穷大时, 如果$ q(t) $趋于$ {0}_{Nn} $, 那么$ (\tilde{\eta}(t),\; \tilde{v}(t)) $趋于$ {0}_{2Nn} $. 再通过$ \tilde{\eta}(t) $和$ \tilde{v}(t) $的定义, 可知$ \mathop{\lim}_{t \to \infty}\eta\;(t) = \eta_e = {1}_N \;\otimes s^* $, $ \mathop{\lim}_{t \to \infty}v\;(t) = v_e $, 又因为$ q\;(t) = s\;(t) - \eta\;(t) $, 所以$ \mathop{\lim}_{t \to \infty}q(t) \;= \;\mathop{\lim}_{t \;\to\; \infty}(s(t)\; -\; \eta(t)) \;=\; 0 $, 即$ \mathop{\lim}_{t \to \infty}s(t) =\mathop{\lim}_{t \to \infty} \eta(t) = {1}_N \otimes s^* $.
□ 注 3. 由定理1的证明可知, 所提出的分布式PI变量$q_i(t) $独立于图的全局信息, 且所需的交互信息量更少, 仅需相对状态信息即可解决最优一致性问题(1). 相比之下, 文献[18]中所提方法需要依赖最优信号发生器生成的参考信息. 因此, 本文提出的分布式PI变量在减少信息量方面具有显著的优势.
3. 基于PI调节的严格反馈多智能体系统最优一致性
为了验证所提出分布式PI变量的适用性, 本节针对N个智能体组成的多智能体系统, 提出一种基于PI调节的近似最优一致性算法, 其中, 每个智能体的动力学模型具有严格反馈的形式, 且带有死区输入和有界扰动.
3.1 严格反馈多智能体系统
考虑具有如下动力学的多智能体系统[27]:
$$\left\{ \begin{aligned} &\dot{x}_{il} = f_{il}(t,\; \bar{x}_{il}) + g_{il}(t,\; \bar{x}_{il})x_{i,\; l+1} \\ &\dot{x}_{im} = f_{im}(t,\; \bar{x}_{im}) + g_{im}(t,\; \bar{x}_{im})D_{i}\big(u_{i} + d_{i}(t,\; \bar{x}_{im})\big) \end{aligned}\right. $$ (27) 其中, $ t\ge0;\; i\in{\cal{V}};\; l=1,\; 2,\; \cdots,\; m-1 $; $ m \geq 2 $表示阶数; $ x_{il}: {\bf{R}}_{\geq 0} \to {\bf{R}} $表示第$ i $个智能体的第 $ l $ 阶状态;$ \bar{x}_{il} = [x_{i1},\; x_{i2},\; \cdots,\; x_{il}]^{\operatorname{T}} \;\in\; {\bf{R}}^l $ 表示状态向量; $ u_i \in {\bf{R}} $表示第 $ i $个智能体的控制输入; 函数$ f_{il},\; g_{il}: {\bf{R}}_{\geq 0} \times {\bf{R}}^l \to {\bf{R}} $在$ t $上是分段连续的, 且对$ \bar{x} _{il} $是Lipschitz的, $ l=1,\; 2,\; \cdots,\; m $; 有界扰动$ d_i: {\bf{R}}_{\geq 0}\; \times {\bf{R}}^m \to {\bf{R}} $在$ t $上是分段连续的, 且对$ \bar{x}_{im} $是Lipschitz的. 此外, 死区函数$ D_i: {\bf{R}} \to {\bf{R}} $定义如下:
$$ \begin{align} D_{i}(u_i) = \begin{cases} h_{i}^{1}(u_i),\;&u_{i}\leq-b_{i}^1 \\ 0,\;&-b_{i}^1 \leq u_{i} \leq b_{i}^2 \\ h_{i}^{2}(u_i),\;&u_{i}\geq b_{i}^2 \end{cases} \end{align} $$ (28) 其中, $ b_{i}^1,\; b_{i}^2 \geq 0 $表示死区断点, $ |b_{i}^2 - b_{i}^1| $表示死区带大小. $ h_i^1: (-\infty,\; -b_i^1] \to {\bf{R}} $, $ h_i^2: [b_i^2,\; +\infty) \to {\bf{R}} $表示死区空间的非线性函数, $ h_i^1 $和$ h_i^2 $都是局部Lipschitz的. 再者, $ h_i^1(-b_i^1) = 0,\; h_i^2(b_i^2) = 0 $, 且满足以下性质:
$$ \begin{align} &\mathop{\lim}_{u_i \to -\infty}h_i^1(u_i) = -\infty \end{align} $$ (29) $$ \begin{align} &\mathop{\lim}_{u_i \to +\infty}h_i^2(u_i) = +\infty \end{align} $$ (30) 因此, 死区函数$ D_i $也是局部Lipschitz的.
接下来, 本文将为系统(27)设计一个分布式控制器, 以实现如下控制目标:
$$ \mathop{\lim}\limits_{t \to \infty} ||x_{i1}(t) - x^{*}|| \leq \kappa\epsilon $$ (31) 其中, $i \in \mathcal{V}$; $ x^* \;\in\; {\bf{R}} $表示全局代价函数$ f(x(t)) = \mathop{\sum}_{i = 1}^{N}f_i(x_{i1}(t)) \in {\bf{R}} $的最小值点; $\kappa,\ \epsilon$的定义参考式(1).
注 4. 上述关于$ D_i $的定义代表一类相当一般的死区函数, 即$ D_i $既是非对称的, 也是非线性的. 其中, 函数$ h_i^1 $和$ h_i^2 $无需是线性的、严格增的, 甚至不需要是可微的. 此外, 函数$ h_i^1 $和$ h_i^2 $未知其具体表达式, 且断点$ b_{i}^1 $和$ b_{i}^2 $不必相同, 亦未必已知.
3.2 分布式控制器设计
为了设计一个能够调节PI变量的分布式控制器, 本文基于预设性能控制[28−29], 提出如下控制方案:
首先, 定义辅助函数$ A_f: (-1,\; 1) \to {\bf{R}} $, 其表达式为: $ A_f(\triangle) = \mathrm{ln}\left(({1 + \triangle})/({1 - \triangle})\right) $.
其次, 构造严格正的性能函数$ \rho_{i1}(t),\; t \geq 0 $, 使其满足以下两点特征:
1) $ \rho_{i1}(0) > |q_i(0)| $;
2) 包含动态和静态性能指标, 如收敛速度、稳态误差等, 即$ \rho_{i1}(t)=(\rho_{i1}(0)-\rho_{i1}(\infty))\mathrm{e}^{-\beta_{i1}t}+\rho_{i1}(\infty) $, 其中, $ \beta_{i1} $表示最小指数收敛速度; $ \rho_{i1}(\infty) =\epsilon_{i1}>0$表示最大稳态误差.
再次, 设置中间控制信号$ \varphi_{i1} $为:
$$ \begin{align} \varphi_{i1}(t,\; x_{i1}) = -\mathrm{sgn}(g_{i1})k_{i1}A_f\left(\frac{q_i(t)}{\rho_{i1}(t)}\right) \end{align} $$ (32) 其中, $ \mathrm{sgn} $表示符号函数, 用来获取函数$ g_{i1} $的正负; $ k_{i1} > 0 $表示控制增益.
然后, 设置中间控制信号$ \varphi_{il},\; l=2,\; \cdots,\; m\ - 1 $为:
$$ \begin{split}\varphi_{il}(t,\; \overline{x}_{il})=\; & -\mathrm{sgn}(g_{il})k_{il}\ \times \\ & A_f\left(\frac{x_{il}-\varphi_{i,\; l-1}(t,\; \overline{x}_{i,\; l-1})}{\rho_{il}(t)}\right)\end{split} $$ (33) 其中, $ k_{il} > 0 $; 性能函数$ \rho_{il} $的构造除初值外与$ \rho_{i1} $ 相同, 其初值满足 $ \rho_{il}\; (0)\; > \; |x_{il}(0)-\varphi_{i,\; l-1}(0, \overline{x}_{i,\; l-1}(0))| $.
最后, 控制输入$ u_i $设置为:
$$ \begin{split}u_i(t,\; \overline{x}_{im})=\; & -\mathrm{sgn}(g_{im})k_{im}\ \times \\ & A_f\left(\frac{x_{im}-\varphi_{i,\; m-1}(t,\; \overline{x}_{i,\; m-1})}{\rho_{im}(t)}\right)\end{split} $$ (34) 其中, $ k_{im} > 0 $; 性能函数$ \rho_{im} $的构造除初值外与$ \rho_{il} $相同, 其初值满足$ \rho_{im}(0) \;>\; |x_{im}(0)\;-\;\varphi_{i,\; m-1}(0, \overline{x}_{i,\; m-1}(0))| $.
对于智能体的状态$ x_{il} $、函数$ g_{il}\ 和\ f_{il} $以及扰动$ d_i,\; t \in [0,\; \infty),\; i \in {\cal{V}},\; l = 1,\; 2,\; \cdots,\; m $, 本文作出如下假设:
假设 3. 智能体$ i $的各阶状态$ x_{i1},\; x_{i2},\; \cdots ,\; x_{im} $均可通过测量获得.
假设 4. 函数$ g_{il} $严格为正或者严格为负, 且符号已知.
假设 5. 存在连续的非负函数$ \bar{d}_{i}: {\bf{R}}^m \to {\bf{R}}_{\geq 0} $和$ \bar{f}_{il}: {\bf{R}}^l \to {\bf{R}}_{\geq 0} $使得:
$$ \begin{align} &|d_{i}(t,\; \bar{x}_{im})| \leq \bar{d}_{i}(\bar{x}_{im}) \end{align} $$ (35) $$ \begin{align} &|f_{il}(t,\; \bar{x}_{il})| \leq \bar{f}_{il}(\bar{x}_{il}) \end{align} $$ (36) 此外, 存在连续且严格正的函数$ \underline{g}_{il},\; \overline{g}_{il}: {\bf{R}}^l \to {\bf{R}}_{\geq 0} $, 使得:
$$ \begin{align} \underline{g}_{il}(\bar{x}_{il}) \leq |g_{il}(t,\; \bar{x}_{il})| \leq \overline{g}_{il}(\bar{x}_{il}) \end{align} $$ (37) 下面, 为了实现严格反馈多智能体系统(27)的最优一致性, 提出本文的第二个定理:
定理 2. 假设严格反馈多智能体系统(27)的通信拓扑满足假设1, 局部代价函数$ f_i,\; i \in {\cal{V}} $满足假设2, 并且假设3 ~ 5均成立, 则使用所设计的分布式控制器(32) ~ (34), 可以使得系统(27)按照预先设定的性能指标实现控制目标(31).
证明. 见附录A.
注 5. 所提出的分布式控制器(32) ~ (34)不需要依赖任何系统非线性的先验知识, 也没有采用近似结构去估计系统中的非线性因素, 如神经网络、模糊系统等, 从而降低了算法设计的复杂度. 此外, 由定理2的证明可知, 分布式控制器(32) ~ (34)的主要限制在于归一化误差变量式$\zeta_{il}$的初值大小, 即要求$ |\zeta_{il}(0)| < 1 $, $ l=1,\; 2,\; \cdots,\; m $. 因此, 需选择适当的控制增益$ k_{il} $和性能函数$ \rho_{il} $, 以确保$ \zeta_{il} $的初值在(−1, 1)之内, 即控制目标的初值在$(- \rho_{il} $, $ \rho_{il} $)之内, 从而可以保证系统按照预设的性能指标收敛.
注 6. 由定理2可知, 基于PI变量$ q_i(t) $所设计的分布式控制器(32) ~ (34)能够实现有向非平衡图下带有死区输入非线性和有界扰动的严格反馈多智能体系统的近似最优一致性. 此外, 鉴于文献[19]的定理2, 通过使用类似的假设和处理方式, 所提出的基于PI调节的分布式控制器也能适用于含有状态时延的情况.
4. 仿真实验
本节通过仿真案例验证了所提出的基于PI调节的严格反馈多智能体系统近似最优一致性算法的有效性. 具体而言, 考虑一个由$ N = 5 $个智能体(单连杆机械臂)组成的多智能体系统, 其中, 第$ i, i \in {\cal{V}} $个智能体的动力学模型[30]描述如下:
$$ \begin{align} \begin{cases} E_i\ddot{p}_i + H_i\dot{p}_i+G_i\mathrm{sin}(p_i) = J_i \\ M_i\dot{J}_i+A_{i}J_i = u_i-K_i\dot{p}_i \end{cases} \end{align} $$ (38) 其中, $ p_i,\; \dot{p}_i,\; \ddot{p}_i $分别表示连杆位置, 速度, 加速度; $ E_i $表示机械惯性; $ H_i $表示接头处的粘性摩擦系数; $ G_i $表示与负载质量和重力系数相关的正常数; $ J_i $表示电气子系统产生的扭矩; $ M_i $表示电枢电感; $ A_i $表示电枢电阻; $ K_i $表示反电动势系数; $ u_i $表示控制输入, 其物理含义为机电转矩. 令$ x_{i1} = p_i $, $ x_{i2} = \dot{p}_i $, $ x_{i3} = J_i $, 则系统(38)可改写为:
$$ \begin{align} \begin{cases} \dot{x}_{i1} = x_{i2} \\ \dot{x}_{i2} = -\displaystyle\frac{G_i}{E_i}\mathrm{sin}(x_{i1}) - \displaystyle\frac{H_i}{E_i}x_{i2} + \displaystyle\frac{1}{E_i}x_{i3} \\ \dot{x}_{i3} = -\displaystyle\frac{K_i}{M_i}x_{i2} - \displaystyle\frac{A_i}{M_i}x_{i3} + \displaystyle\frac{1}{M_i}D_i(u_i + d_i) \end{cases} \end{align} $$ (39) 其中, $ D_i $和$ d_i $分别是施加到 $ u_i $上的死区函数和有界扰动, 具体定义见式(27)和(28).
首先, 给定上述系统的仿真参数. 令$ x_{11}(0)= 1.5 $, $ x_{21}(0) = 2 $, $ x_{31}(0) = -3.7 $, $ x_{41}(0) = 5.5 $, $ x_{51}(0) = -1 $. $ x_{i2} $, $ x_{i3} $的初值均为$ 0 $, 这表明在初始时刻, 所有机械臂呈静止状态且扭矩为零. $ E_i = 1 \ \mathrm{kg}\cdot\mathrm{m}^2 $, $ H_i = 1\ \mathrm{Nm}\cdot\mathrm{s}/ \mathrm{rad} $, $ G_i = 10 $, $ M_i = 0.05\ \mathrm{H} $, $ A_i=0.5\ \Omega $, $ K_i = 10 \ \mathrm{Nm/A} $. 此外, 有界扰动$ d_i = 0.2\mathrm{cos}(t)+2x_{i2}+ x_{i3} $, 死区函数$ D_i(u_i) $的表达式如下所示:
$$ D_i(u_i)=\left\{\begin{aligned} & 0.8u_i+0.6\mathrm{sin}(0.9\pi)+2.346\ 4, \\ & \qquad\qquad\qquad\qquad\; \; \; \; \;u_i > 6.4 \\ & 0.12u_i^2+0.54u_i+0.6\mathrm{sin}(0.9\pi)-0.904\ 8 ,\\ & \qquad\qquad\qquad\qquad\; \; \;\; \; 1.3 < u_i\le6.4 \\ & 0.6\mathrm{sin}((u_i-0.4)\pi),\; \; \; \, 0.4 < u_i\le1.3 \\ & 0,\qquad\qquad\qquad\qquad\; \, -0.2\le u_i\le0.4 \\ & \mathrm{sin}(2(u_i+0.2)\pi),\; \; \; \; \; \; -0.8\le u < -0.2 \\ & -0.28u_i^2+0.6u_i+0.659\ 2-\mathrm{sin}(1.2\pi), \\ & \qquad\qquad\qquad\qquad\, \; \; \;\; -3.6\le u_i < -0.8 \\ & 2u_i+2.070\ 4-\mathrm{sin}(1.2\pi), \\ & \qquad\qquad\qquad\qquad\; \; \; \; \;u_i < -3.6\; \end{aligned}\right. $$ (40) 其次, 分布式控制器(32) ~ (34)的参数分别设置为: $ k_{i1} = 60 $, $ k_{i2} = 600 $, $ k_{i3} = 1\ 600 $, 性能函数$ \rho_{i1}=\rho_{i2}=\rho_{i3}=69.95\mathrm{e}^{-0.1t}+0.05 $, 分布式PI变量$q_i(t) $的调节参数设置为$ a=1,\ b=2 $.
然后, 每个智能体的代价函数设置为$ f_1=\mathrm{ln}(2 +\mathrm{e}^{x_{11}})+0.5x_{11}^2 $, $f_2 = (x_{21} + 2)^2$, $ f_3=\mathrm{ln}(1+\mathrm{e}^{x_{31}})\; + x_{31}^2 $, $ f_4 = 0.5x_{41}^2 $, $f_5 = (x_{51}+2.5)^2$. 显然, 所有代价函数均满足假设2.
最后, 给出各智能体之间的通信拓扑图, 即图1所示的权重非平衡有向图.
从图2(b)可以看出所有智能体的状态变量都能够收敛到最优状态$x^*=-1.17 $的邻域. 进一步地, 为了突出所提出基于PI变量调节的近似最优一致性算法的优势, 在相同的动力学模型、假设条件、调节参数以及控制方式下, 将文献[18] 采用嵌入式框架的最优一致性算法与本文算法进行对比, 实验结果如图2所示. 从图2(a)和图2(b)可以看出, 本文所提出的算法响应速度更快、收敛时间更短.
本文算法的其他仿真结果如图3 ~ 图7所示. 由图3可知, 状态变量$ x_{i2} $能够收敛至零, 并结合状态变量$ x_{i1} $的收敛情况, 可以发现所提出的近似最优一致性算法能够满足控制目标(31). 图4显示, 状态变量$ x_{i3} $最终收敛至$ -9.22 $, 这表明系统(38)的收敛使所有智能体维持相同的扭矩. 从图5可以看出, 全局代价函数的梯度$ \sum_{i = 1}^{5}\nabla f_i(x_{i1}(t)) $收敛至零的邻域, 验证了所提出的基于PI调节的控制算法能够实现严格反馈多智能体系统(27)的近似最优一致性问题.
从图6可以发现, PI变量$ q_i(t) $按照指数衰减的趋势收敛至零的邻域, 表明了所设计的控制器(32) ~ (34)以及选择的性能函数能够将PI变量控制在预设性能范围内. 进一步地, 结合图2 ~ 图5可以发现, 状态变量$ x_{i1},\; x_{i2},\; x_{i3} $以及全局代价函数的梯度$ \sum_{i = 1}^{5} \nabla f_i(x_{i1}(t)) $的收敛趋势和收敛时间均与PI变量$ q_i(t) $类似, 这表明PI变量不仅决定系统能否收敛, 而且还影响其收敛性能. 由图7可知, 变量$ \hat{r}_i(t) $收敛至Laplacian矩阵零特征值对应的左特征向量, 且收敛后的$ \hat{r}_i(t) $满足引理1中的性质, 证明了本文设计的Laplacian矩阵零左特征向量估计器(4)的有效性.
5. 结束语
本文通过PI调节的方式, 解决了有向图下严格反馈多智能体系统的最优一致性问题. 首先, 提出一种新的分布式PI变量, 该变量独立于图的全局信息和系统的动力学模型, 并将最优一致性问题转化为PI调节问题, 使得经典的控制技术能够通过PI调节来实现复杂系统的最优一致性. 然后, 本文将所提出的分布式PI变量与预设性能控制相结合, 设计一类新的分布式控制算法, 有效实现了权重非平衡有向图下带有死区输入非线性和有界扰动的严格反馈多智能体系统的近似最优一致性. 最后, 对所提出的算法进行仿真验证, 证明了该算法的有效性. 未来的研究将考虑更一般的切换拓扑和时变代价函数.
附录 A. 定理 2 的证明
证明.
首先, 定义归一化误差变量:
$$ \begin{align} \zeta_{i1} = \frac{q_i(t)}{\rho_{i1}(t)} \end{align} $$ (A1) $$ \begin{align} \zeta_{il} = \frac{x_{il} - \varphi_{i,\; l-1}(t,\; \bar{x}_{i,\; l-1})}{\rho_{il}(t)} \end{align} $$ (A2) 其中, $ t \geq 0;\; i \in {\cal{V}};\; l = 2,\; \cdots,\; m $. 于是, 中间控制信号$ \varphi_{il} $和控制输入$ u_i $可重新表述为:
$$ \begin{align} &\varphi_{il}(\zeta_{il}) = -\mathrm{sgn}(g_{il})k_{il}A_f(\zeta_{il}) \end{align} $$ (A3) $$ \begin{align} &u_{i}(\zeta_{im}) = -\mathrm{sgn}(g_{im})k_{im}A_f(\zeta_{im}) \end{align} $$ (A4) 其中, $ i \in {\cal{V}};\; l = 1,\; 2,\; \cdots,\; m-1 $. 结合式(2)、(8)、(A1)以及(A2), 状态变量可改写为:
$$ \begin{align} &x_{i1} = \zeta_{i1}\rho_{i1}(t) + \eta_{i} \end{align} $$ (A5) $$ \begin{align} &x_{il} = \zeta_{il}\rho_{il}(t) + \varphi_{i,\; l-1}(\zeta_{i,\; l-1}) \end{align} $$ (A6) 然后, 联立式(2)、(8)、(27)、(A1)、(A5)以及(A6)对$ \zeta_{i1} $关于时间$ t $求导, 可得:
$$ \begin{split} \dot{\zeta}_{i1} =\;& \phi_{i1}(t,\; \zeta_{i1},\; \zeta_{i2}) = \frac{1}{\rho_{i1}(t)}\;\times\\&\Big[f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i) + g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t) \;+\\ & \eta_i)(\zeta_{i2}\rho_{i2}(t) + \varphi_{i1}(\zeta_{i1})) - \dot{\eta}_i - \zeta_{i1}\dot{\rho}_{i1}(t)\Big]\\ \end{split} $$ (A7) 同理可得:
$$ \begin{split} \dot{\zeta}_{i2} =\;&\phi_{i2}(t,\; \zeta_{i1},\; \zeta_{i2},\; \zeta_{i3}) =\frac{1}{\rho_{i2}(t)}\Big[f_{i2}(t,\; \zeta_{i1}\rho_{i1}(t) \;+\\ & \eta_i,\; \zeta_{i2}\rho_{i2}(t) + \varphi_{i1}(\zeta_{i1}))+ g_{i2}(t,\; \zeta_{i1}\rho_{i1}(t)\;+\\ &\eta_i,\; \zeta_{i2}\rho_{i2}(t) + \varphi_{i1}(\zeta_{i1}))(\zeta_{i3}\rho_{i3}(t)+\varphi_{i2}(\zeta_{i2})) \;- \\ &\frac{\mathrm{d}\varphi_{i1}}{\mathrm{d}\zeta_{i1}} \times \phi_{i1}(t,\; \zeta_{i1},\; \zeta_{i2})- \zeta_{i2}\dot{\rho}_{i2}(t)\Big] \end{split} $$ (A8) $$ \begin{split} & \dot{\zeta}_{il}= \phi_{il}(t,\; \zeta_{i1},\; \zeta_{i2},\; \cdots,\; \zeta_{i,\;l+1}) = \frac{1}{\rho_{il}(t)}\;\times\\ &\Big[f_{il}(t,\, \zeta_{i1}\rho_{i1}(t)+\eta_i,\, \cdots,\, \zeta_{il}\rho_{il}(t) + \varphi_{i,\, l-1}(\zeta_{i,\, l-1})) \;+\\ & g_{il}(t,\; \zeta_{i1}\rho_{i1}(t) +\eta_i,\; \cdots,\; \zeta_{il}\rho_{il}(t)+\varphi_{i,\; l-1}(\zeta_{i,\; l-1})) \;\times\\& (\zeta_{i,\; l+1}\rho_{i,\; l+1}(t)+\varphi_{il}(\zeta_{il})) - \frac{\mathrm{d}\varphi_{i,\; l-1}}{\mathrm{d}\zeta_{i,\; l-1}} \;\times \\&\phi_{i,\; l-1}(t,\; \zeta_{i1},\; \cdots,\; \zeta_{il}) - \zeta_{il}\dot{\rho}_{il}(t)\Big] \end{split} $$ (A9) 其中, $ l=3,\; \cdots,\; m-1 $.
$$ \begin{split}& \dot{\zeta}_{im} = \phi_{im}(t,\; \zeta_{i1},\; \zeta_{i2},\; \cdots,\; \zeta_{im}) =\frac{1}{\rho_{im}(t)}\Big[f_{im}(t,\;\\ & \;\;\;\zeta_{i1}\rho_{i1}(t)+\eta_i,\; \cdots,\; \zeta_{im}\rho_{im}(t) + \varphi_{i,\;m-1} (\zeta_{i,\; m-1})) \;+\\ & \;\;\; g_{im}(t,\; \zeta_{i1}\rho_{i1}(t) +\eta_i,\; \cdots,\; \zeta_{il}\rho_{il}(t) \;+\\ & \;\;\;\varphi_{i,\; m-1}(\zeta_{i,\; m-1}))\times D_i(u_{i}(\zeta_{im}) + d_i(\zeta_{i1}\rho_{i1}(t)\;+\\ & \;\;\;\eta_i,\; \cdots,\; \zeta_{im}\rho_{im}(t) + \varphi_{i,\; m-1}(\zeta_{i,\; m-1}))) \;-\\ & \;\;\; \frac{\mathrm{d}\varphi_{i,\; m-1}}{\mathrm{d}\zeta_{i,\; m-1}}\phi_{i,\; m-1}(t,\; \zeta_{i1},\; \cdots,\; \zeta_{im}) - \zeta_{im}\dot{\rho}_{im}(t)\Big]\\ \end{split} $$ (A10) 最后, 令$ \zeta_i=\mathrm{col}(\zeta_{i1},\; \cdots,\; \zeta_{im}),\; i\in\cal{V} $. 进一步地, 令$ \zeta=\mathrm{col}(\zeta_1,\; \cdots,\; \zeta_N) $. 于是, 联立式(A7) ~ (A10), 建立如下的初值问题:
$$ \begin{align} \dot{\zeta} = \phi\left(t,\; \zeta\right) = \begin{bmatrix}\phi_{11}\left(t,\; \zeta_{11},\; \zeta_{12}\right)\\ \vdots\\ \phi_{Nm}\left(t,\; \zeta_{N1},\; \cdots,\; \zeta_{Nm}\right)\end{bmatrix} \end{align} $$ (A11) 此外, $ \zeta(0) \in {{\Omega}}_{\zeta} $, 其中, 开集$ {{\Omega}}_{\zeta} $为:
$$ \begin{align} {{\Omega}}_{\zeta} = \underbrace{(-1,\; 1)\times\cdots\times(-1,\; 1)}_{N \times m} \end{align} $$ (A12) 显然, 由式(32) ~ (34)中性能函数$ \rho_{il},\; l= 1, 2,\; \cdots,\; m $的初值设置可知, $ \zeta(0) \in {{\Omega}}_{\zeta} $恒成立.
因此, 定理2的证明可分为如下两个阶段: 初值问题(A11)解的有界性证明和$ \tau_{\mathrm{max}} $大小的确定.
阶段 1. 初值问题(A11)解的有界性.
根据式(A7) ~ (A10)可知, $ \phi_{il} $的性质与$ \rho_{il} $, $ f_{il} $, $ g_{il} $, $ u_i $, $ d_i $, $ \varphi_{il} $, 以及$ \dot{\eta}_i $有关, $ i \in {\cal{V}},\; l = 1,\; 2,\; \cdots, m $. 由性能函数的特征可知, $ \rho_{il} $在$ t $上是连续可微且有界的. 动力学模型(27)表明$ f_{il},\; g_{il} $以及$ d_i $在$ t $上是连续可微的, 且对$ \zeta_{il} $是Lipschitz的. 由式(A3)和(A4)可知, $ \varphi_{il} $和$ u_i $是光滑的, 即$ \varphi_{il} $和$ u_i $对$ \zeta_{il} $是Lipschitz的. 此外, 式(10)表明$ \dot{\eta}_{i} $对$ \zeta_{il} $也是Lipschitz的. 综上所述, $ \phi_{il} $在$ t $上是连续可微的, 并对$ \zeta_{il} $是Lipschitz的. 因此, 由文献[28]中的定理1可知, 对于初值问题(A11), 在$ t \in [0,\; \tau_{\mathrm{max}}) $上, $ \zeta(t) \in {{\Omega}}_{\zeta} $, 其中$ \tau_{\mathrm{max}} > 0 $.
阶段 2. $ \tau_{\mathrm{max}} $大小的确定.
由$ \zeta(t) \in {{\Omega}}_{\zeta} $和式(A12)可知, $ \forall t \in [0,\; \tau_{\mathrm{max}}) $,$ \zeta_{il}(t) \in (-1,\; 1),\; i \in {\cal{V}},\; l = 1,\; 2,\; \cdots,\; m $, 即$ \zeta_{il}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的. 此外, 令$ \varepsilon_ {il}(t) = A_f(\zeta_{il}) $, 即$ \forall t \in [0,\; \tau_{\mathrm{max}}) $, $ \zeta_{il}(t) \in (-1,\; 1) $有:
$$ \begin{align} \varepsilon_{il}(t) = \mathrm{ln}\left(\frac{1+\zeta_{il}(t)}{1-\zeta_{il}(t)}\right) \end{align} $$ (A13) 其中, $ i \in {\cal{V}};\; l = 1,\; 2,\; \cdots,\; m $.
下面, 对$ \varepsilon_{il}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上的有界性进行分析.
当$ l = 1 $时, 令$ V_{i1} = \frac{1}{2}\varepsilon_{i1}^2(t) $, 并对时间$ t $求导, 联立式(A7)和(A13)可得:
$$ \begin{split} \dot{V}_{i1} =\; &\frac{2\varepsilon_{i1}(t)}{(1-\zeta_{i1}^2)\rho_{i1}(t)}\Big[f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i) \;+ \\ &g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i)(\zeta_{i2}\rho_{i2}(t) \;+\\ &\varphi_{i1}(\zeta_{i1})) - \dot{\eta}_i - \zeta_{i1}\dot{\rho}_{i1}(t)\Big] \end{split} $$ (A14) 令$ F_{i1}(t) = f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t) + \eta_i) - \dot{\eta}_i - \zeta_{i1}\dot{\rho}_{i1}(t) $, 由上述分析可知, $ \rho_{i1}(t) $和$ \zeta_{i1} $是有界的, 再根据式(A1)和定理1可推得$ \eta_i $有界. 结合假设5中函数$ \bar{f}_{i1} $的连续性、极值定理以及式(36), 可进一步推导出:
$$ \begin{align} |f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t) +\; \eta_i)| \leq \bar{f}_{i1}(\zeta_{i1}\rho_{i1}(t) +\; \eta_i) \leq f_{i1}^U \end{align} $$ (A15) 其中, $ f_{i1}^U > 0 $. 所以$ f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\ \eta_i) $在$ t \in [0, \tau_{\mathrm{max}}) $上是有界的. 相似地, 由式(10)可推得$ \dot{\eta}_i $是有界的. 此外, 由性能函数的表达式可知, $ \dot{\rho}_{i1}(t) $也是有界的. 因此, $ F_{i1}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即存在$ \bar{F}_{i1} > 0 $, 使得$ |F_{i1}(t)| \leq \bar{F}_{i1},\; \forall t \in [0,\; \tau_{\mathrm{max}}) $. 所以有:
$$ \begin{split} & \dot{V}_{i1} \leq\frac{2|\varepsilon_{i1}(t)|}{(1-\zeta_{i1}^2)\rho_{i1}(t)}\Big[\bar{F}_{i1} + \mathrm{sgn}(\varepsilon_{i1}) \;\times \\ &\quad g_{i1}(t,\;\zeta_{i1}\rho_{i1}(t)+\eta_i) \times(\zeta_{i2}\rho_{i2}(t)+\varphi_{i1}(\zeta_{i1}))\Big]\\ \end{split} $$ (A16) 进一步地, 由式(A3)和(A13)可得:
$$ \begin{split} \dot{V}_{i1} \leq\;& \frac{2|\varepsilon_{i1}(t)|}{(1-\zeta_{i1}^2)\rho_{i1}(t)}\Big[\bar{F}_{i1} + |g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t) + \eta_i)|\;\times\\ &\big(\mathrm{sgn}(\varepsilon_{i1})\mathrm{sgn}(g_{i1})\zeta_{i2}\rho_{i2}(t) - k_{i1}|\varepsilon_{i1}(t)|\big)\Big] \end{split} $$ (A17) 由$ \zeta_{il} \in (-1,\; 1),\; l = 1,\; 2,\; \cdots,\; m $和性能函数的表达式可知, $ \zeta_{i2} $和$ \rho_{i2}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即存在$ \bar{B}_{i2} \;>\; 0 $, 使得$ |\; \mathrm{sgn}(\varepsilon_{i1})\mathrm{sgn}(g_{i1})\zeta_{i2}\rho_{i2}(t)|\le \overline{B}_{i2},\ \forall t\in[0,\; \tau_{\mathrm{max}}) $, 则有:
$$ \begin{split} \dot{V}_{i1} \leq\;& \frac{2|\varepsilon_{i1}(t)||g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i)|}{(1-\zeta_{i1}^2)\rho_{i1}(t)} \;\times \\ &\bigg[\frac{\bar{F}_{i1}}{|g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i)|} \;+\; \bar{B}_{i2} -\; k_{i1}|\varepsilon_{i1}(t)|\bigg]\\ \end{split}$$ (A18) 由假设5可知, 函数$ \underline{g}_{il},\; l = 1,\; 2,\; \cdots,\; m $是连续且严格正的. 再结合$ \zeta_{i1} $, $ \rho_{i1}(t) $的有界性和极值定理, 可推得函数$ \underline{g}_{il} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上也是有界的, 即存在$ g_{i1}^L > 0 $, 使得$ 0 < g_{i1}^L \leq \underline{g}_{i1}(\zeta_{i1}\rho_{i1}(t)+ \eta_i) \leq |g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i)|,\; \forall t \in [0,\; \tau_{\mathrm{max}}) $. 所以有:
$$ \begin{split} \dot{V}_{i1} \leq\; &\frac{2|\varepsilon_{i1}(t)||g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i)|}{(1-\zeta_{i1}^2)\rho_{i1}(t)} \;\times \\ &\bigg[\frac{\bar{F}_{i1}}{g_{i1}^L} + \bar{B}_{i2} - k_{i1}|\varepsilon_{i1}(t)|\bigg] \end{split} $$ (A19) 因为$ \zeta_{i1} \in (-1,\; 1) $, 所以$ 1-\zeta_{i1}^2(t) > 0 $. 此外, 由于$ \rho_{i1}(t) > 0 $, 则在$ t \in [0,\; \tau_{\mathrm{max}}) $上, 当$ {\bar{F}_{i1}}/{g_{i1}^L}+\bar{B}_{i2}\; - k_{i1}|\varepsilon_{i1}(t)| \,<\, 0 $, 即 $ |\varepsilon_{i1}(t)| > {\bar{F}_{i1}}/({g_{i1}^{L}k_{i1}}) + {\bar{B}_{i2}}/{k_{i1}}, \forall t \in [0, \tau_{\mathrm{max}}) $时, $ \dot{V}_{i1} < 0. $
因此, 经上述分析可知, $ \varepsilon_{i1} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即$ \forall t \in [0,\; \tau_{\mathrm{max}}) $,
$$ |\varepsilon_{i1}(t)| \leq \bar{\varepsilon}_{i1} = \mathrm{max}\left\{|\varepsilon_{i1}(0)|,\; \frac{\bar{F}_{i1}}{g_{i1}^{L}k_{i1}} + \frac{\bar{B}_{i2}}{k_{i1}}\right\} $$ (A20) 进一步地, 根据(A13)可求得$ \varepsilon_{i1} $的反函数, 并结合$ \varepsilon_{i1} $的有界性可推得:
$$ \begin{align} -1 < \underline{\zeta}_{i1} \leq \zeta_{i1}(t) \leq \overline{\zeta}_{i1} < 1,\; \quad \forall t \in [0,\; \tau_{\mathrm{max}}) \end{align} $$ (A21) 其中, $ \underline{\zeta}_{i1}=(\mathrm{e}^{-\overline{\varepsilon}_{i1}}-1)/(\mathrm{e}^{-\overline{\varepsilon}_{i1}}+1) $, $ \overline{\zeta}_{i1}=(\mathrm{e}^{\overline{\varepsilon}_{i1}}-1)/ (\mathrm{e}^{\overline{\varepsilon}_{i1}}+1) $.
此外, 由式(A3)和(A13)可知, $ \varphi_{i1} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即$ |\varphi_{i1}| \leq k_{i1}\bar{\varepsilon}_{i1},\; \forall t \in [0,\; \tau_{\mathrm{max}}) $. 对$ \varphi_{i1} $求关于时间$ t $的导数, 则有:
$$ \begin{split} \dot{\varphi}_{i1} =\;& \frac{\mathrm{d}\varphi_{i1}}{\mathrm{d}\zeta_{i1}}\phi_{i1}(t,\; \zeta_{i1},\; \zeta_{i2}) = \frac{-2\mathrm{sgn}(g_{i1})k_{i1}}{(1-\zeta_{i1}^2)\rho_{i1}(t)}\;\times\\ &\Big[f_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)+\eta_i) + g_{i1}(t,\; \zeta_{i1}\rho_{i1}(t)\;+\\ &\eta_i)(\zeta_{i2}\rho_{i2}(t)+\varphi_{i1}(\zeta_{i1})) - \dot{\eta}_i - \zeta_{i1}\dot{\rho}_{i1}(t)\Big]\\ \end{split} $$ (A22) 采用与之前相同的分析方式可推得$ \dot{\varphi}_{i1} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上也是有界的.
当$ l = 2,\; 3,\; \cdots,\; m-1 $时, 令$ V_{il} = \frac{1}{2}\varepsilon_{il}^2(t),\; \forall t \in [0,\; \tau_{\mathrm{max}}) $, 采用与$ V_{i1} $相同的分析可推得, $ \varepsilon_{il}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上也是有界的, 即$ \forall t \in [0,\; \tau_{\mathrm{max}}) $,
$$ |\varepsilon_{il}(t)| \leq \bar{\varepsilon}_{il} = \mathrm{max}\left\{|\varepsilon_{il}(0)|,\; \frac{\bar{F}_{il}}{g_{il}^{L}k_{il}} + \frac{\bar{B}_{i,\; l+1}}{k_{il}}\right\} $$ (A23) 同样地, 可得到:
$$ \begin{align} -1 < \underline{\zeta}_{il} \leq \zeta_{il}(t) \leq \overline{\zeta}_{il} < 1,\; \quad \forall t \in [0,\; \tau_{\mathrm{max}}) \end{align} $$ (A24) 其中, $ \underline{\zeta}_{il}=(\mathrm{e}^{-\overline{\varepsilon}_{il}}-1)/(\mathrm{e}^{-\overline{\varepsilon}_{il}}+1) $, $ \overline{\zeta}_{il}=(\mathrm{e}^{\overline{\varepsilon}_{il}}-1)/ (\mathrm{e}^{\overline{\varepsilon}_{il}}+1) $, 以及$ \varphi_{il},\; \dot{\varphi}_{il} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的.
当$ l = m $时, 令$ V_{im} = \frac{1}{2}\varepsilon_{im}^2(t),\; \forall t \in [0,\; \tau_{\mathrm{max}}) $, 并对时间$ t $求导, 可得:
$$ \begin{split} \dot{V}_{im} =\;& \frac{2\varepsilon_{im}(t)}{(1-\zeta_{im}^2)\rho_{im}(t)}\Big[f_{im}(t,\; \bar{\mu}_{im}) + g_{im}(t,\; \bar{\mu}_{im}) \;\times\\ & D_{i}\big(u_i(\zeta_{im}) + d_{i}(t,\; \bar{\mu}_{im})\big) - \frac{\mathrm{d}\varphi_{i,\; m-1}}{\mathrm{d}\zeta_{i,\; m-1}} \;\times\\ & \phi_{i,\; m-1}(t,\; \zeta_{i1},\; \cdots,\; \zeta_{im}) - \zeta_{im}\dot{\rho}_{im}(t)\Big] \end{split} $$ (A25) 其中, $ \bar{\mu}_{im}\;\;=\;\;\zeta_{i1}\rho_{i1}(t)\;\;+\;\;\eta_i,\; \cdots,\; \zeta_{im}\rho_{im}(t)\;+ \varphi_{i,\; m-1} (\zeta_{i,\; m-1}) $. 同样地, 令
$$ \begin{split} F_{im}(t) =\;&f_{im}(t,\; \bar{\mu}_{im}) - \frac{\mathrm{d}\varphi_{i,\; m-1}}{\mathrm{d}\zeta_{i,\; m-1}}\; \times \\&\phi_{i,\; m-1}(t,\; \zeta_{i1},\; \cdots,\; \zeta_{im}) - \zeta_{im}\dot{\rho}_{im}(t)\\ \end{split} $$ (A26) 由$ l = m-1 $的分析结果可知, $ \dot{\varphi}_{i,\; m-1} $有界, 即$ \frac{\mathrm{d}\varphi_{i,\; m-1}}{\mathrm{d}\zeta_{i,\; m-1}}\times\phi_{i,\; m-1}\; (\; t,\; \zeta_{i1},\; \cdots\ ,\; \zeta_{im}\; ) $有界, 剩余部分采用与之前相同的证明, 可推得$ |F_{im}(t)| \leq \bar{F}_{im}, \forall t \in [0,\; \tau_{\mathrm{max}}) $, 则有:
$$ \begin{split} \dot{V}_{im} \leq\;& \frac{2|\varepsilon_{im}(t)||g_{im}(t,\; \bar{\mu}_{im})|}{(1-\zeta_{im}^2)\rho_{im}(t)}\bigg[\frac{\bar{F}_{im}}{|g_{im}(t,\; \bar{\mu}_{im})|} \;+ \\&\mathrm{sgn}(\varepsilon_{im})\mathrm{sgn}(g_{im})D_i\big(u_i(\zeta_{im}) + d_i(t,\; \bar{\mu}_{im})\big)\bigg]\\ \end{split} $$ (A27) 相似地, 存在$ g_{im}^L > 0 $, 使得$ 0 < g_{im}^L \leq \underline{g}_{im}(\bar{\mu}_{im}) \leq |g_{im}(t,\; \bar{\mu}_{im})|,\; \forall t \in [0,\; \tau_{\mathrm{max}}) $, 再结合式(A4)和(A13)可推得:
$$ \begin{split} &\dot{V}_{im} \leq \frac{2|\varepsilon_{im}(t)||g_{im}(t,\; \bar{\mu}_{im})|}{(1-\zeta_{im}^2)\rho_{im}(t)}\bigg[\frac{\bar{F}_{im}}{g_{im}^L} + \mathrm{sgn}(\varepsilon_{im}) \; \times \\&\;\;\;\;\;\;\mathrm{sgn}(g_{im})D_i\big(d_i(t,\; \bar{\mu}_{im}) - \mathrm{sgn}(g_{im}) \times k_{im}\varepsilon_{im}(t)\big)\bigg]\\ \end{split} $$ (A28) 为了证明$ \varepsilon_{im}(t) $在$ t \in [0,\; \tau_{\mathrm{max}}) $上的有界性, 需要先找出$ \dot{V}_{im}<0 $的条件, 即找出合适的$ \varepsilon_{im}(t) $, 使得:
$$ \begin{split}\frac{\bar{F}_{im}}{g_{im}^L} \;+\; & \mathrm{sgn}(\varepsilon_{im})\mathrm{sgn}(g_{im})D_i\big(d_i(t,\; \bar{\mu}_{im}) \;- \\ &\mathrm{sgn}(g_{im})k_{im}\varepsilon_{im}(t)\big) < 0 \end{split} $$ (A29) 接下来, 根据$ \mathrm{sgn}(g_{im}) $取值, 分两种情况讨论.
情况 1. $ \mathrm{sgn}(g_{im}) = 1 $.
由于$ \mathrm{sgn}(g_{im}) = 1 $, 所以式(A29)左侧可改写为:
$$ \begin{align} \frac{\bar{F}_{im}}{g_{im}^L} + \mathrm{sgn}(\varepsilon_{im})D_i\big(d_i(t,\; \bar{\mu}_{im}) - k_{im}\varepsilon_{im}(t)\big) \end{align} $$ (A30) 由式(29)可知, 存在$ c_{i1} \in {\bf{R}} $, $ \forall u_i < c_{i1} $, 使得:
$$ \begin{align} \mathop{\lim}_{u_i\to-\infty}D_{i}(u_i) & = -\infty \Rightarrow D_i(u_i) < -\frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A31) 因此, 令$ D_i(d_i(t,\; \bar{\mu}_{im}) - k_{im}\varepsilon_{im}(t)) < -{\bar{F}_{im}}/{g_{im}^L} $, 即要求存在合适的$ \varepsilon_{im}(t) $, 使得$ d_i(t,\; \overline{\mu}_{im})-k_{im}\ \times \varepsilon_{im}(t) < c_{i1} $. 此外, 由假设5可知, 函数$ \bar{d}_{i} $是非负且连续的. 通过$ l = 2,\; \cdots,\; m-1 $的分析可知, $ \varphi_{il} $是有界的. 结合$ \eta_i,\; \zeta_{il},\; \rho_{il},\; l = 1,\; \cdots,\; m $的有界性以及极值定理, 可推得函数$ \bar{d}_{i} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即存在$ d_{i}^U > 0 $, 使得$ |d_{i}(t,\; \bar{\mu}_{im})| \leq \bar{d}_{i}(\bar{\mu}_{im}) \leq d_{i}^U ,\; \forall t \in [0,\; \tau_{\mathrm{max}}) $. 所以, 解得$ \varepsilon_{im}(t) > ({d_{i}^U - c_{i1}})/ {k_{im}} $. 进一步地, 上式还可表达为:
$$ \begin{align} \varepsilon_{im}(t) > \frac{|d_{i}^U - c_{i1}|}{k_{im}},\; \quad \forall t \in [0,\; \tau_{\mathrm{max}}) \end{align} $$ (A32) 同理, 由式(30)可知, 存在$ c_{i2} \in {\bf{R}} $, $ \forall u_i > c_{i2} $, 使得:
$$ \begin{align} \mathop{\lim}_{u_i\to+\infty}D_{i}(u_i) & = +\infty \Rightarrow D_i(u_i) > \frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A33) 令$ D_i(d_i(t,\; \bar{\mu}_{im}) - k_{im}\varepsilon_{im}(t)) > {\bar{F}_{im}}/{g_{im}^L} $, 并能推得$ \varepsilon_{im}(t) < ({-d_{i}^U - c_{i2}})/{k_{im}} $. 进一步地, 上式还可表达为:
$$ \begin{align} \varepsilon_{im}(t) < -\frac{|-d_{i}^U - c_{i2}|}{k_{im}},\; \quad \forall t \in [0,\; \tau_{\mathrm{max}}) \end{align} $$ (A34) 联合式(A30) ~ (A34)可知, 在$ t \in [0,\; \tau_{\mathrm{max}}) $上有
$$ \begin{align} \mathrm{sgn}(\varepsilon_{im})D_i(d_i(t,\; \bar{\mu}_{im}) - k_{im}\varepsilon_{im}(t)) < -\frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A35) 其中, $ \varepsilon_{im}(t) $的取值范围为$ |\varepsilon_{im}(t)| > {\sigma_i}/{k_{im}} $, $ \sigma_i = \mathrm{max} \left\{|d_{i}^U - c_{i1}|,\; |-d_{i}^U - c_{i2}|\right\} $.
情况 2. $ \mathrm{sgn}(g_{im}) = -1 $.
由于$ \mathrm{sgn}(g_{im}) = -1 $, 所以式(A29)可改写为:
$$ \begin{align} \frac{\bar{F}_{im}}{g_{im}^L} - \mathrm{sgn}(\varepsilon_{im})D_i\big(d_i(t,\; \bar{\mu}_{im}) + k_{im}\varepsilon_{im}(t)\big) \end{align} $$ (A36) 采用与情况1相同的分析可知, 在$ t \in [0,\; \tau_{\mathrm{max}}) $上, 当$ \varepsilon_{im}(t) < -{|c_{i1}-d_{i}^U|}/{k_{im}} $时, 有:
$$ \begin{align} D_i\big(d_i(t,\; \bar{\mu}_{im}) + k_{im}\varepsilon_{im}(t)\big) < -\frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A37) 当$ \varepsilon_{im}(t) > {|-c_{i2} - d_{i}^U|}/{k_{im}} $时, 有:
$$ \begin{align} D_i\big(d_i(t,\; \bar{\mu}_{im}) + k_{im}\varepsilon_{im}(t)\big) > \frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A38) 联合式(A36) ~ (A38)可得:
$$ \begin{align} \mathrm{sgn}(\varepsilon_{im})D_i(d_i(t,\; \bar{\mu}_{im}) + k_{im}\varepsilon_{im}(t)) > \frac{\bar{F}_{im}}{g_{im}^L} \end{align} $$ (A39) 其中, $ \varepsilon_{im}(t) $的取值范围为$ |\varepsilon_{im}(t)| > {\sigma_i}/{k_{im}} $, $ {\sigma}_i = \mathrm{max} \left\{|c_{i1} - d_{i}^U|,\; |- c_{i2} - d_{i}^U |\right\} $.
综上所述, 当$ |\varepsilon_{im}(t)| > {\sigma_i}/{k_{im}} $时, 式(A29)成立, 进而$ \dot{V}_{im} < 0 $. 因此, $ \varepsilon_{im} $在$ t \in [0,\; \tau_{\mathrm{max}}) $上是有界的, 即$ \forall t \in [0,\; \tau_{\mathrm{max}}) $,
$$ \begin{align} |\varepsilon_{im}(t)| \leq \bar{\varepsilon}_{im} = \mathrm{max}\left\{|\varepsilon_{im}(0)|,\; \frac{\sigma_i}{k_{im}}\right\} \end{align} $$ (A40) 进一步地, 根据(A13)可求得$ \varepsilon_{im} $的反函数, 并结合$ \varepsilon_{im} $的有界性推得:
$$ \begin{align} -1 < \underline{\zeta}_{im} \leq \zeta_{im}(t) \leq \overline{\zeta}_{im} < 1,\;\;\; \forall t \in [0,\; \tau_{\mathrm{max}}) \end{align} $$ (A41) 其中, $ \underline{\zeta}_{i1}=(\mathrm{e}^{-\overline{\varepsilon}_{i1}}-1)/(\mathrm{e}^{-\overline{\varepsilon}_{i1}}+1) $, $ \overline{\zeta}_{i1}=(\mathrm{e}^{\overline{\varepsilon}_{i1}}-1)/ (\mathrm{e}^{\overline{\varepsilon}_{i1}}+1) $.
下面证明$ \tau_{\mathrm{max}} = +\infty $. 经上述分析可知, 在$ t \in [0,\; \tau_{\mathrm{max}}) $上, $ \zeta_{il}(t) \in {{\Omega}}_{\zeta'},\; l = 1,\; \cdots,\; m $, $ {{\Omega}}_{\zeta'} \subset {{\Omega}}_{\zeta} $, 其中,
$$ \begin{align} {{\Omega}}_{\zeta'} = \underbrace{(-\underline{\zeta}_{i1},\; \overline{\zeta}_{i1})\times\cdots\times(-\underline{\zeta}_{im},\; \overline{\zeta}_{im})}_{N \times m} \end{align} $$ (A42) 即不存在$ t'\in [0,\; \tau_{\mathrm{max}}) $使得$ \zeta_{il}(t') \notin {{\Omega}}_{\zeta'} $. 显然, 这与文献[28]中提议1的结论相矛盾, 故$ \tau_{\mathrm{max}} = +\infty $. 因此, $ \forall t \in [0,\; +\infty) $,
$$ \begin{align} \zeta_{il}(t) \in {{\Omega}}_{\zeta'} \subset {{\Omega}}_{\zeta},\;\quad i \in {\cal{V}},\; l = 1,\; \cdots,\; m \end{align} $$ (A43) 所以, 在$ t \in [0,\; +\infty) $上, 有:
$$ \begin{align} \begin{cases} -\rho_{i1}(t) \leq q_i(t) \leq \rho_{i1}(t) \\ -\rho_{il}(t) \leq x_{il} - \varphi_{i,\; l-1}(t,\; \bar{x}_{i,\; l-1}) \leq \rho_{il}(t) \end{cases} \end{align} $$ (A44) 其中, $ l=2,\; 3,\; \cdots,\; m $. 由性能函数$ \rho_{i1}(t) $的表达式可知, $\mathop{\lim}_{t\to+\infty}\rho_{i1}(t) = \epsilon_{i1} > 0 $, 结合$ -\rho_{i1}(t) \leq q_i(t) \leq \rho_{i1}(t) $可知, $\mathop{\lim}_{t\to+\infty}||q_i(t)|| \leq \epsilon_{i1}$. 进一步可推得$ \mathop{\lim}_{t\to+\infty}||q(t)|| \;\leq\; \epsilon$, 其中, $\epsilon \;=\; \mathop{\max}\{\epsilon_{11},\; \epsilon_{21},\; \cdots, \epsilon_{N1}\} $. 令 $s_i(t) \;\;= \;\;x_{i1}(t) $, 根据定理1可知, $\mathop{\lim}_{t\to+\infty} ||x_{i1}(t) - x^*|| \leq \kappa\epsilon$. 由此可见, 系统(27)能按照预先设定的性能指标实现控制目标(31), 定理2成立.
□ -
表 1 文中变量和符号含义
Table 1 The meaning of variables and symbols
变量 符号含义 $ {\boldsymbol x}_k $ $ k $时刻目标状态真实值 $ {\boldsymbol x}_{k|k} $ $ k $时刻目标状态估计值 $ {\boldsymbol P}_{k|k} $ $ k $时刻目标状态估计误差协方差 $ {\boldsymbol z}_k $ 传感器在$ k $时刻的量测值 $ {\boldsymbol\omega}_k $ $ k $时刻的系统噪声 $ {\boldsymbol\upsilon}_k $ $ k $时刻的量测噪声 $ {\boldsymbol Q}_k $ $ k $时刻系统噪声方差 $ {\boldsymbol R}_k $ $ k $时刻量测噪声方差 $ d_x $ 目标状态向量的维数 $ d_z $ 量测向量的维数 $ {\boldsymbol F}_{k|k-1} $ $ k-1 $时刻到$ k $时刻的状态转移矩阵 $ {\boldsymbol H}_{k} $ $ k $时刻量测矩阵 $ {\boldsymbol\psi}_{k} $ 变分分布参数 $ p\left({\boldsymbol x}_k|{\boldsymbol z}_{k}\right) $ $ k $时刻目标状态后验分布 $ q\left({\boldsymbol x}_k|{\boldsymbol\psi}_{k}\right) $ 以$ {\boldsymbol\psi}_{k} $为参数的变分分布 $ \mathcal{L}\left({\boldsymbol\psi}_{k}\right) $ 以$ {\boldsymbol\psi}_{k} $为变分分布参数的置信下界 ${\cal{D}}\left(q\left({\boldsymbol x}_k|{\boldsymbol\psi}_{k}\right)|| p\left({\boldsymbol x}_k|{\boldsymbol z}_{k}\right)\right)$ 变分分布$ q\left({\boldsymbol x}_k|{\boldsymbol\psi}_{k}\right) $与状态后验分布$ p\left({\boldsymbol x}_k|{\boldsymbol z}_{k}\right) $的KL散度 $ {\boldsymbol J}_{{\boldsymbol\psi}_k} $ 以$ {{\boldsymbol\psi}_k} $为参数的 Fisher 信息矩阵 $ \mathcal{M} $ 流形空间 $ \mathcal{S} $ 流形空间中的概率分布集合 $ \mathcal{F} $ 流形空间中的平滑映射函数 ${{\boldsymbol v} }_{OP}$ 流形空间中的$ O $点处指向$ P $点的切向量 $|{{\boldsymbol v} }_{OP}|$ 切向量${{\boldsymbol v} }_{OP}$的模 表 2 目标的轨道根数
Table 2 The orbital elements of target
半长轴 (km) 离心率 倾角 (deg) 近地点角 (deg) 升交点赤经 (deg) 7500 0.1 15 30 12 表 3 算法平均估计误差均值的对比
Table 3 Comparison of the mean estimation errors of the algorithm
算法 EKF UKF IEKF VBKF-NG x 轴位置 6.8731 6.8493 6.8290 6.6025 x 轴速度 − 0.0382 0.0418 − 0.0368 − 0.0241 y 轴位置 2.8793 2.8770 2.8688 2.7563 y 轴速度 − 0.0257 − 0.0243 − 0.0222 − 0.0103 z 轴位置 4.7665 4.6759 4.5546 4.3286 z 轴速度 0.1272 0.1097 0.1076 0.1062 表 4 算法平均运行时间$( \times 10^{-4}\;{\rm{s}}) $的对比
Table 4 The comparison of the average run time $( \times 10^{-4}\;{\rm{s}}) $ of algorithms
算法 EKF UKF IEKF VBKF-NG 时间 0.2580 0.9542 1.0989 1.1205 表 5 算法平均RMSE的对比
Table 5 Comparison of average RMSE between algorithms
算法 径向距(km) 径向速度(km/min) EKF 1.5052 0.0813 UKF 1.1833 0.0670 VBKF-COD 0.9181 0.0432 VBKF-NG 0.2990 0.0161 -
[1] 潘泉, 胡玉梅, 兰华, 孙帅, 王增福, 杨峰. 信息融合理论研究进展: 基于变分贝叶斯的联合优化. 自动化学报, 2019, 45(7): 1207−1233Pan Quan, Hu Yu-Mei, Lan Hua, Sun Shuai, Wang Zeng-Fu, Yang Feng. Information fusion progress: Joint optimization based on variational Bayesian theory. Acta Automatica Sinica, 2019, 45(7): 1207−1233 [2] Ito K, Xiong K Q. Gaussian filters for nonlinear filtering problems. IEEE Transactions on Automatic Control, 2000, 45(5): 910−927 doi: 10.1109/9.855552 [3] Hu Y M, Wang X Z, Pan Q, Hu Z T, Moran B. Variational Bayesian Kalman filter using natural gradient. Chinese Journal of Aeronautics, 2022, 35(5): 1−10 doi: 10.1016/j.cja.2021.08.033 [4] Gu D B. A game theory approach to target tracking in sensor networks. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics), 2011, 41(1): 2−13 doi: 10.1109/TSMCB.2010.2040733 [5] 潘泉, 胡玉梅, 马季容. 基于变分贝叶斯联合优化的情报监视与侦察. 指挥信息系统与技术, 2020, 11(2): 1−8Pan Quan, Hu Yu-Mei, Ma Ji-Rong. Intelligence, surveillance and reconnaissance based on variational Bayesian joint optimization. Command Information System and Technology, 2020, 11(2): 1−8 [6] Hu X Y, Yang L Q, Xiong W. A novel wireless sensor network frame for urban transportation. IEEE Internet of Things Journal, 2015, 2(6): 586−595 doi: 10.1109/JIOT.2015.2475639 [7] Spantini A, Baptista R, Marzouk Y. Coupling techniques for nonlinear ensemble filtering. SIAM Review, 2022, 64(4): 921−953 doi: 10.1137/20M1312204 [8] Silva B, Fisher R M, Kumar A, Hancke G P. Experimental link quality characterization of wireless sensor networks for underground monitoring. IEEE Transactions on Industrial Informatics, 2015, 11(5): 1099−1110 doi: 10.1109/TII.2015.2471263 [9] Vu T, Rahmani A. Distributed consensus-based Kalman filter estimation and control of formation flying spacecraft: Simulation and validation. In: Proceedings of the AIAA Guidance, Navigation, and Control Conference. Kissimmee, USA: AIAA, 2015. 7−12 [10] Wang B H, Chen W S, Zhang B, Shi P, Zhang H Y. A nonlinear observer-based approach to robust cooperative tracking for heterogeneous spacecraft attitude control and formation applications. IEEE Transactions on Automatic Control, 2023, 68(1): 400−407 doi: 10.1109/TAC.2022.3143082 [11] Tronarp F, García-Fernández Á F, Särkkä S. Iterative filtering and smoothing in nonlinear and non-Gaussian systems using conditional moments. IEEE Signal Processing Letters, 2018, 25(3): 408−412 doi: 10.1109/LSP.2018.2794767 [12] Humpherys J, West J. Kalman filtering with Newton's method. IEEE Control Systems Magazine, 2010, 30(6): 101−106 doi: 10.1109/MCS.2010.938485 [13] Alessandri A, Gaggero M. Moving-horizon estimation for discrete-time linear and nonlinear systems using the gradient and newton methods. In: Proceedings of the IEEE 55th Conference on Decision and Control. Las Vegas, USA: IEEE, 2016. 2906−2911 [14] Akyildiz Ö D, Chouzenoux É, Elvira V, Míguez J. A probabilistic incremental proximal gradient method. IEEE Signal Processing Letters, 2019, 26(8): 1257−1261 doi: 10.1109/LSP.2019.2926926 [15] Gultekin S, Paisley J. Nonlinear Kalman filtering with divergence minimization. IEEE Transactions on Signal Processing, 2017, 65(23): 6319−6331 doi: 10.1109/TSP.2017.2752729 [16] Hoffman M D, Blei D M, Wang C, Paisley J. Stochastic variational inference. The Journal of Machine Learning Research, 2013, 14(1): 1303−1347 [17] Salimans T, Kingma D P, Welling M. Markov chain Monte Carlo and variational inference: Bridging the gap. In: Proceedings of the 32nd International Conference on International Conference on Machine Learning. Lille, France: JMLR.org, 2015. 1218−1226 [18] Acerbi L. Variational Bayesian Monte Carlo. In: Proceedings of the 32nd International Conference on Neural Information Processing Systems. Montréal, Canada: Curran Associates Inc., 2018. 8223−8233 [19] Acerbi L. Variational Bayesian Monte Carlo with noisy likelihoods. In: Proceedings of the 34th International Conference on Neural Information Processing Systems. Vancouver, Canada: Curran Associates Inc., 2020. Article No. 688 [20] Wang P Y, Blunsom P. Collapsed variational Bayesian inference for hidden Markov models. In: Proceedings of the 16th International Conference on Artificial Intelligence and Statistics. Scottsdale, USA: AISTATS, 2013. 599−607 [21] Amari S I. Natural gradient works efficiently in learning. Neural Computation, 1998, 10(2): 251−276 doi: 10.1162/089976698300017746 [22] 李宇波. 基于信息几何理论的滤波方法研究 [博士学位论文], 国防科技大学, 中国, 2017.Li Yu-Bo. Study of Filtering Methods via Information Geometry [Ph.D. dissertation], National University of Defense Technology, China, 2017. [23] Zhang G D, Sun S Y, Duvenaud D, Grosse R B. Noisy natural gradient as variational inference. In: Proceedings of the 35th International Conference on Machine Learning. Stockholm, Sweden: ICML, 2018. 5847−5856 [24] Ollivier Y. Online natural gradient as a Kalman filter. Electronic Journal of Statistics, 2018, 12(2): 2930−2961 [25] Cheng Y Q, Wang X Z, Morelande M, Moran B. Information geometry of target tracking sensor networks. Information Fusion, 2013, 14(3): 311−326 doi: 10.1016/j.inffus.2012.02.005 [26] Schmitt L, Fichter W. Globally valid posterior Cramér-Rao bound for three-dimensional bearings-only filtering. IEEE Transactions on Aerospace and Electronic Systems, 2019, 55(4): 2036−2044 doi: 10.1109/TAES.2018.2881352 [27] 胡玉梅, 潘泉, 胡振涛, 郭振. 基于自然梯度的噪声自适应变分贝叶斯滤波算法. 自动化学报, 2023, 49(10): 2094−2108Hu Yu-Mei, Pan Quan, Hu Zhen-Tao, Guo Zhen. A novel noise adaptive variational Bayesian filter using natural gradient. Acta Automatica Sinica, 2023, 49(10): 2094−2108 [28] Duan T, Avati A, Ding D Y, Thai K K, Basu S, Ng A, et al. NGBoost: Natural gradient boosting for probabilistic prediction. In: Proceedings of the 37th International Conference on Machine Learning. JMLR.org. Vienna, Austria: 2020. Article No. 252 [29] Tseng P. An analysis of the EM algorithm and entropy-like proximal point methods. Mathematics of Operations Research, 2004, 29(1): 27−44 doi: 10.1287/moor.1030.0073 [30] Khan M E, Baqué P, Fleuret F, Fua P. Kullback-Leibler proximal variational inference. In: Proceedings of the 28th International Conference on Neural Information Processing Systems. Montreal, Canada: MIT Press, 2015. 3402−3410 [31] Smidl V, Quinn A. Variational Bayesian filtering. IEEE Transactions on Signal Processing, 2008, 56(10): 5020−5030 doi: 10.1109/TSP.2008.928969 [32] Beal M J, Ghahramani Z. The Variational Kalman Smoother, Technical Report GCNU TR, 2001, 3, Computational Neuroscience Unit, University College London, UK, 2001. [33] Hu Y M, Wang X Z, Lan H, Wang Z F, Moran B, Pan Q. An iterative nonlinear filter using variational Bayesian optimization. Sensors, 2018, 18(12): Article No. 4222 [34] Ble D M, Kucukelbir A, McAuliffe J D. Variational inference: A review for statisticians. Journal of the American Statistical Association, 2017, 112(518): 859−877 doi: 10.1080/01621459.2017.1285773 [35] Sain R, Mittal V, Gupta V. A comprehensive review on recent advances in variational Bayesian inference. In: Proceedings of the International Conference on Advances in Computer Engineering and Applications. Ghaziabad, India: IEEE, 2015. 488−492 [36] Beal M J. Variational Algorithms for Approximate Bayesian Inference [Ph.D. dissertation], Cambridge University, UK, 2003. [37] Amari S I. Information Geometry and Its Applications. Tokyo: Springer, 2016. [38] Amari S, Douglas S C. Why natural gradient? In: Proceedings of International Conference on Acoustics, Speech and Signal Processing. Seattle, USA: IEEE, 1998. 1213−1216 [39] Rezende D J, Mohamed S, Wierstra D. Stochastic backpropagation and approximate inference in deep generative models. In: Proceedings of the 31st International Conference on International Conference on Machine Learning. Beijing, China: JMLR.org, 2014. II-1278−II-1286 [40] Lan H, Liang Y, Zhang W, Yang F, Pan Q. Iterated minimum upper bound filter for tracking orbit maneuvering targets. In: Proceedings of the 16th International Conference on Information Fusion. Istanbul, Turkey: IEEE, 2013. 1051−1057 -