-
摘要: 研究了信息物理系统中假数据注入(False data injection, FDI)攻击信号的检测问题. 在分布式融合框架下, 首先将FDI攻击信号建模为信息物理系统模型中的未知输入, 从而使得攻击信号的检测问题转化为对FDI攻击信号的实时估计问题. 其次, 在每个传感器端设计基于自适应卡尔曼滤波的FDI攻击信号的局部估计器; 在融合中心端引入补偿因子, 设计分布式信息融合准则以导出攻击信号的融合估计器. 特别地, 当FDI攻击信号是时变情况时, 融合过程中补偿因子的引入可以大大提高对攻击信号的估计精度. 最后, 通过两个仿真算例验证所提算法的有效性.
-
关键词:
- 自适应Kalman滤波 /
- 假数据注入攻击 /
- 攻击信号估计 /
- 信息融合
Abstract: This paper is concerned with the estimation problem of false data injection (FDI) attacks in cyber-physical systems, where the false data is injected into the actuator. Under the distributed fusion framework, the attack signal is modeled as an unknown parameter in the addressed system models, and each local estimate of the attack signal at the sensor is obtained based on adaptive Kalman filtering method. Subsequently, by introducing the compensation factor, an optimal weighing fusion criterion is designed to improve estimation precision of attack signals. Particularly, when the FDI attack signal is time-varying, the introduction of compensation factor during fusing process can indeed improve the estimation precision largely. Finally, two illustrative examples are used to show the effectiveness of the proposed methods. -
信息物理系统(Cyber-physical systems, CPSs)是一个综合了计算、网络和物理环境的多维智能化复杂系统, 它借助有线或无线通信网络将各个关键设施整合在一起, 使得人机和物理进程的交互更加便捷. 随着网络通信、嵌入式系统、计算机控制及相关硬件技术的不断发展, 信息物理系统已引起了工业界的广泛关注[1-3]. 然而在此框架下, 原有系统的封闭性被打破, 使其面临着来自网络攻击的安全威胁. 例如: 2010年, 伊朗的首座核电站−布什核电站被震网病毒攻击, 使得伊朗的第一座核设施推迟发电, 严重损害了伊朗的工业设施[4]; 2011年美国伊利诺伊州一处水利控制系统遭到网络攻击, 险些造成大面积的供水中断. 当信息物理系统遭受网络攻击时, 准确、及时地监测到攻击信号对于监测中心采取高效的防御策略是至关重要的[5].
从控制角度来看, 信息物理系统是传统数字控制系统融合通信技术的下一代网络化控制系统[6]. 信息物理系统中被控对象的量测信号和控制信号均通过网络进行传输, 因此无论是传感器到控制器端, 还是控制器到执行器端, 均有受到恶意网络攻击的可能, 从而影响系统的稳定性, 甚至引发系统崩溃, 造成严重的生产事故与经济损失[7-8]. 典型的网络攻击有三种, 分别是拒绝服务攻击(Denial of service, DoS)、欺骗攻击和重放攻击. 针对拒绝服务攻击, 文献[9]将DoS攻击建模为一类能量约束问题, 针对攻击者与防御者设计递归分布式卡尔曼估计器进行双边优化; 文献[10]假定DoS攻击有界的情况下, 通过构建嵌套切换模型得到了CPS在基于包控制方法下的稳定性条件. 针对重放攻击, 文献[11]从防御者的角度提出一种带补偿策略的数学模型来描述重放攻击和带宽约束, 并在线性最小方差意义下设计了递归分布式卡尔曼融合估计器. 欺骗攻击又称为假数据注入(False data injection, FDI)攻击, 它通过向系统注入错误的控制信号或测量信号影响信号数据的准确性. 注意到FDI攻击可以通过欺骗攻击检测机制来影响信息物理系统, 从而造成攻击检测器的漏报或虚警[12]. 文献[13]中研究了电力系统中针对FDI攻击的状态估计问题, 揭示了现有错误测量检测算法中存在的脆弱性. 该研究表明,即使攻击者的资源受限, 依旧可以改变状态估计的结果; 而在假定攻击者可以获取系统参数及所有数据流的基础上, 此类攻击则可以对系统产生一定影响的同时却不被检测到. 目前信息物理系统中的FDI攻击主要采用异常信号检测方法进行检测[14], 如基于二元假设的贝叶斯检测[15-18]、基于卡尔曼滤波器的
$ \chi^2 $ 检测[19-21]和加权最小平方检测[22-24]. 后两种方法均是利用观测值构造新息残差, 然后通过与一个给定的阈值比较来判决是否受到攻击. 攻击检测的阈值对检测精度起着重要作用, 然而阈值通常在检测前根据已有经验进行选择, 这无疑会降低检测精度. 而且对于CPSs而言, 特别是电力、医疗系统, 仅仅检测或识别到攻击是不够的, 因为这些基础设施系统无法快速对其关闭、重启、恢复, 以对抗攻击. 为此, 不同于文献[15-24]中基于阈值的FDI攻击信号检测方法, 本文将对CPSs中控制器与执行器间的通信网络遭受的FDI攻击信号进行实时估计, 不仅可以检测攻击信号是否存在, 而且还可以掌握攻击信号的基本特征, 从而使防御方采取简单有效的补救措施以降低系统性能的损失程度.为了对FDI攻击信号进行实时估计, 本文将攻击信号建模为状态动态方程中的一个未知参数, 然后在每个传感器端利用自适应卡尔曼滤波设计FDI攻击信号的局部估计. 由于基于递归最小二乘的自适应Kalman滤波方法适应于被估计的参数是时不变或者变化缓慢的, 而恶意的攻击信号往往是时变的. 因此, 本文充分利用多传感器融合所提供的冗余信息, 通过引入补偿因子和设计分布式融合准则, 来提高时变攻击信号的估计精度. 注意到补偿因子是通过影响局部互协方差的信息来提高估计性能, 而局部互协方差矩阵存在于分布式融合结构中. 最后, 通过两个仿真验证所设计的融合算法对FDI攻击信号估计的有效性. 特别地, 当攻击信号是时变时, 仿真结果表明所引入的补偿因子可以明显提高攻击信号的融合估计性能.
符号说明.
${\rm{E}}$ 表示数学期望,$ {\rm diag}\{N_1,\cdots,N_n\} $ 表示由$ N_1,\cdots,N_n $ 组成的块对角矩阵,$ I $ 表示单位矩阵,$ \perp $ 表示正交.1. 问题描述与分析
考虑如下离散线性时变状态空间模型
$$ {\boldsymbol x}(k) = A(k){\boldsymbol x}(k-1)+B(k){\boldsymbol u}(k-1)+{\boldsymbol w}(k)$$ (1) 其中,
$ {\boldsymbol x}(k)\in{\bf { R}}^n $ 为系统状态,$ {\boldsymbol u}(k)\in{\bf { R}}^l $ 为控制输入;$ A(k) $ 与$ B(k) $ 分别为系统状态转移矩阵和控制输入矩阵. 当控制输入信号$ {\boldsymbol u}(k) $ 受到FDI攻击时, 被攻击后的控制输入信号描述为$$ {\boldsymbol u}_a(k) = {\boldsymbol u}(k)+{\boldsymbol \theta}(k) $$ (2) 其中,
$ {\boldsymbol \theta}(k)\in{\bf { R}}^l $ 表示FDI攻击信号. 因此, 控制信号遭受攻击后的系统(1)建模为$$ {\boldsymbol x}(k) = A(k){\boldsymbol x}(k-1)+B(k)({\boldsymbol u}(k-1)+{\boldsymbol \theta}(k-1))+{\boldsymbol w}(k) $$ (3) 然后, 利用
$ L $ 个传感器对系统(3)的状态进行实时监测, 且每个传感器的量测$ {\boldsymbol y}_i(k) $ 建模为$$ {\boldsymbol y}_i(k) = C_i(k){\boldsymbol x}(k)+{\boldsymbol v}_i(k),\;\; i = 1,2,\cdots,L $$ (4) 其中,
$ C_i(k) $ 为量测矩阵.$ {\boldsymbol w}(k)\in{\bf { R}}^n $ ,$ {\boldsymbol v}_i(k)\in{\bf { R}}^{m_i} $ 是相互独立的零均值高斯白噪声, 且满足$$\begin{split} &{ {\rm{E}}}\left\{\left[ \begin{aligned} {\boldsymbol w}(k)\\ {\boldsymbol v}_i(k)\\ {\boldsymbol v}_j(k) \end{aligned}\right]\left[ \begin{aligned} {\boldsymbol w}(k)\\ {\boldsymbol v}_i(k)\\ {\boldsymbol v}_j(k) \end{aligned}\right]^{\rm T}\right\} = \\ &\qquad{\rm diag}\{Q(k), R_i(k),R_j(k)\},\;\;i\ne j \end{split}$$ (5) 本文要解决的问题是: 针对控制信号遭受FDI攻击篡改的信息物理系统(3), 根据各个传感器的量测信息, 计算得到FDI攻击信号的局部估计
$ \hat{{\boldsymbol \theta}}_i(k);$ 然后, 基于局部估计$\hat{{\boldsymbol \theta}}_i(k),\;i = 1,2,\cdots,L,$ 设计满足$\sum\nolimits^L_{i = 1}W^\theta_i(k) = I$ 最优加权融合矩阵$ W^\theta_1(k),\cdots,W^\theta_L(k) ,$ 使得攻击信号的分布式融合估计器$$ \hat{{\boldsymbol \theta}}_0(k) = \sum\limits_{i = 1}^L W^\theta_i(k)\hat{{\boldsymbol \theta}}_i(k) $$ (6) 在线性最小方差意义下是最优的.
注1. 本文的主要思想是将攻击信号
$ {\boldsymbol \theta}(k) $ 建模为系统(3)中的未知输入, 然后在信息融合框架下利用多组传感器的冗余量测信息, 实现对$ {\boldsymbol \theta}(k) $ 的实时估计. 特别地, 分布式信息融合不仅可以保证所设计估计器的鲁棒性, 而且也能够提高对攻击信号的估计精度. 注意到攻击信号的融合估计器与原闭环系统(1)中的控制信号设计是相互独立的, 即: 本文的控制信号是原有闭环控制系统根据未被攻击情况下的期望性能而提前设计给出的, 不依赖于融合中心的状态估计结果. 此外, 由于文中所给出的局部估计算法是以递推形式给出, 所以实现局部的状态估计不需要局部可观性条件. 但是, 为了保证局部估计器可以提供较好的估计性能, 则需要$ (A(k),C_i(k)) $ 是可观的.2. 攻击信号的估计 − 分布式融合策略
在分布式融合结构下, 首先基于每个传感器量测信息
$ \{{\boldsymbol y}_i(1),\cdots,{\boldsymbol y}_i(k)\} $ , 结合Kalman滤波和带有遗忘因子的递归最小二乘估计方法导出攻击信号的局部估计$ \hat{{\boldsymbol \theta}}_i(k) $ , 主要结果由定理1给出.定理 1. 给定初始值
$P_i(0\vert 0) = P_0$ ,$\Upsilon_i(0) = 0$ ,$ S_i(0) =$ $ \omega I_p $ ,$ \hat{{\boldsymbol \theta}}_i(0) = {\boldsymbol \theta}_0, $ $ \hat{{\boldsymbol x}}_i(0\vert 0) = {\boldsymbol x}_0 $ 和遗忘因子$ 0 < {\lambda _i} \le 1 $ , 则局部估计$ \hat{{\boldsymbol \theta}}_i(k) $ 由以下递推公式计算$$ \left\{ \begin{aligned} &\widetilde{{\boldsymbol y}}_i(k) = {\boldsymbol y}_i(k)-C_i(k)[A(k)\hat{{\boldsymbol x}}_i(k-1\vert k-1)\;+ \\ &\qquad \quad\;\; B(k)({\boldsymbol u}(k-1)+\hat{{\boldsymbol \theta}}_i(k-1))] \\ &\hat{{\boldsymbol \theta}}_i(k) = \hat{{\boldsymbol \theta}}_i(k-1)+\Gamma_i(k)\widetilde{{\boldsymbol y}}_i(k) \\ &\hat{{\boldsymbol x}}_i(k\vert k) = A(k)\hat{{\boldsymbol x}}_i(k-1\vert k-1)\;+ \\ &\qquad \;\;\;\;\quad B(k)({\boldsymbol u}(k-1)+\hat{{\boldsymbol \theta}}_i(k-1))\;+\\ &\qquad \;\;\;\;\quad K_i(k)\widetilde{{\boldsymbol y}}_i(k)+ \Upsilon_i(k)[\hat{{\boldsymbol \theta}}_i(k)\;-\\ &\qquad \;\;\;\;\quad\hat{{\boldsymbol \theta}}_i(k-1)] \end{aligned}\right. $$ (7) $$ \left\{ \begin{aligned} &P_i(k\vert k-1) = A(k)P_i(k-1\vert k-1)A^{\rm T}(k)+Q(k) \\&\Sigma_i(k) = C_i(k)P_i(k\vert k-1)C_i^{\rm T}(k)+R_i(k) \\&K_i(k) = P_i(k\vert k-1)C_i^{\rm T}(k)\Sigma_i^{-1}(k) \\&P_i(k\vert k) = [I_n-K_i(k)C_i(k)]P_i(k\vert k-1) \end{aligned}\right. $$ (8) $$ \left\{ \begin{aligned} &\Upsilon_i(k) = [I_n-K_i(k)C_i(k)]A(k)\Upsilon_i(k-1)\;+ \\&\qquad \;\;\;\quad[I_n-K_i(k)C_i(k)]B(k) \\&\Omega_i(k) = C_i(k)A(k)\Upsilon_i(k-1)+C_i(k)B(k) \\&\Lambda_i(k) = [\lambda_i\Sigma_i(k)+\Omega_i(k)S_i(k-1)\Omega_i^{\rm T}(k)]^{-1} \\&\Gamma_i(k) = S_i(k-1)\Omega_i^{\rm T}(k)\Lambda_i(k) \\&S_i(k) = \frac{1}{\lambda_i}S_i(k-1)-\frac{1}{\lambda_i}\Gamma_i(k)\Omega_i(k)S_i(k-1) \end{aligned}\right.\;\;\; $$ (9) 证明. FDI攻击信号
$ {\boldsymbol \theta}(k) $ 的局部估计器的推导类似于文献[25]与文献[26]的推导过程, 故在此省略. □注 2. 局部FDI攻击信号估计
$ \hat{{\boldsymbol \theta}}_i(k) $ 和局部状态估计$ \hat{{\boldsymbol x}}_i(k\vert k) $ 由递推式(7)计算, 其中, 未知输入$ {\boldsymbol \theta}(k-1) $ 的真实值由$ \hat{{\boldsymbol \theta}}_i(k-1) $ 给出. 当$ {\boldsymbol \theta}(k-1) = \hat{{\boldsymbol \theta}}_i(k-1) $ 时, 则$$ \widetilde{{\boldsymbol y}}_i(k) = C_i(k)A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)+C_i(k){\boldsymbol w}(k)+{\boldsymbol v}_i(k) $$ (10) 其中,
$ \widetilde{{\boldsymbol x}}_i(k-1\vert k-1) = {\boldsymbol x}(k-1)-\hat{{\boldsymbol x}}_i(k-1\vert k-1) .$ 显然新息(10)与标准Kalman滤波的新息相同, 此时将$ B(k)\hat{{\boldsymbol \theta}}_i(k-1) $ 作为局部状态估计器的输入信号, 则可直接由Kalman滤波得到状态的估计值; 当${\boldsymbol \theta}(k-1)\ne $ $ \hat{{\boldsymbol \theta}}_i(k-1)$ 时, 将无法避免攻击信号的估计误差, 此时则用$ \Upsilon_i(k)[\hat{{\boldsymbol \theta}}_i(k)-\hat{{\boldsymbol \theta}}_i(k-1)] $ 补偿这一误差. 另一方面, 遗忘因子$ \lambda_i\;(0 < {\lambda _i} \le 1) $ 是用来控制过去量测信息被遗忘的速度, 以体现新近数据的作用, 即: 小的遗忘因子$ \lambda_i $ 代表过去的量测信息将被快速遗忘. 根据递归最小二乘法的核心思想, 当$ {\boldsymbol \theta}(k) $ 是缓慢变化时, 式(7)中的局部估计器$ \hat{{\boldsymbol \theta}}_i(k) $ 的性能随着$ \lambda_i $ 的减小而提高; 而当$ {\boldsymbol \theta}(k) $ 是时不变参数时, 则式(7)中的局部估计器$ \hat{{\boldsymbol \theta}}_i(k) $ 的性能随着$ \lambda_i $ 的减小而降低[26]. 特别地, 当$ {\boldsymbol \theta}(k) $ 是变化速率快时, 不能通过调节遗忘因子保证局部估计器的收敛性.定义
$$ \left\{ \begin{aligned} & \widetilde{{\boldsymbol \theta}}_i(k) = {\boldsymbol \theta}(k)-\hat{{\boldsymbol \theta}}_i(k) \\ & P^\theta_{ij}(k) = {\rm E}(\widetilde{{\boldsymbol \theta}}_i(k)\tilde{{\boldsymbol \theta}}^{\rm T}_j(k),\;\;\forall i,j \end{aligned} \right. $$ (11) 当在
$ k $ 时刻$ L $ 个局部估计$ \hat{{\boldsymbol \theta}}_i(k) $ 由定理1给出时, 根据文献[27]的结论可以导出满足式(6)最优融合估计$\hat{{\boldsymbol \theta}}_0(k) = $ $ \sum^L_{i = 1}W^\theta_i(k)\hat{{\boldsymbol \theta}}_i(k) $ 的权重$W^\theta_i(k),\;i = 1,2,\cdots,L,$ 由下式计算得到$$ W^\theta(k) = \Sigma_\theta^{-1}(k)e_\theta\left(e_\theta^{\rm T}\Sigma_\theta^{-1}(k)e_\theta\right)^{-1} $$ (12) 其中,
$W^\theta(k) = [W^\theta_1(k),W^\theta_2(k),\cdots,W^\theta_L(k)]^{\rm T}$ 与$e_\theta = $ $ [I_p,\cdots, I_p]^{\rm T}$ 均为$ pL\times p $ 的矩阵;$ \Sigma_\theta(k) = (P^\theta_{ij}(k)) $ 为$ pL\times pL $ 的对称正定矩阵. 在此情形下, 最优分布式融合估计器$ \hat{{\boldsymbol \theta}}_0(k) $ 的协方差矩阵为$$ P^\theta_0(k) = (e_\theta^{\rm T}\Sigma_\theta^{-1}(k)e_\theta)^{-1} $$ (13) 且满足
${\rm tr}(P^\theta_{ii}(k))\ge{\rm tr}(P^\theta_0(k))$ .由式(11)可知, 最优加权矩阵
$W^\theta_i(k)\;(i = 1,2,\cdots, $ $ L)$ 的计算需要求解每个协方差矩阵$P^\theta_{ij}(k)\;(\forall i,j)$ . 因此, 定理2将给出$ P^\theta_{ij}(k) $ 的递推形式. 在给出主要结果之前, 需要定义如下变量$$ \left\{ \begin{aligned} &\widetilde{{\boldsymbol x}}_i(k\vert k) = {\boldsymbol x}(k)-\hat{{\boldsymbol x}}_i(k\vert k) \\&P^x_{ij}(k) = { {\rm{E}}}(\widetilde{{\boldsymbol x}}_i(k\vert k)\tilde{{\boldsymbol x}}^{\rm T}_j(k\vert k)) \\&\Psi_{ij}(k) = { {\rm{E}}}(\widetilde{{\boldsymbol x}}_i(k\vert k)\tilde{{\boldsymbol \theta}}^{\rm T}_j(k)) \end{aligned}\right. $$ (14) 定理 2. 给定补偿因子
$ \eta>0 $ , FDI攻击信号的协方差矩阵$ P^\theta_{ij}(k) $ 由以下的递推公式计算$$ \begin{split} P^\theta_{ij}(k) =\;& [I_p-\Gamma_i(k)C_i(k)B(k)]P^\theta_{ij}(k-1)\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;+\\ &\Gamma_i(k)C_i(k)A(k)P^x_{ij}(k-1)\;\times \\ &A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)-[I_p-\Gamma_i(k)C_i(k)B(k)]\;\times \\ &\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;-\\ &\Gamma_i(k)C_i(k)A(k)\Psi_{ij}(k-1)\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;+\\ &\Gamma_i(k)C_i(k)Q(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &\Gamma_i(k)\delta_{ij}R_i\Gamma^{\rm T}_j(k)+\eta I\\[-10pt] \end{split} $$ (15) 其中,
$$ \begin{split} P^x_{ij}(k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)P^x_{ij}(k-1)A^{\rm T}(k)+B(k)P^\theta_{ij}(k-1)\;\times\\ &B^{\rm T}(k)+ A(k)\Psi_{ij}(k-1)B^{\rm T}(k)\;+\\ &B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)\;+ \\ &Q(k)][I_n-(K_j(k)+\Upsilon_j(k)\Gamma_j(k))C_j(k)]^{\rm T}\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]\delta_{ij}R_i[K_j(k)\;+\\ &\Upsilon_j(k)\Gamma_j(k)]^{\rm T} \\[-10pt]\end{split} $$ (16) $$ \begin{split} {\Psi}_{ij}(k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\&[A(k)\Psi_{ij}(k-1)+B(k)P^\theta_{ij}(k-1)]\;\times \\&[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;- \\&[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\&[A(k)P^x_{ij}(k-1)A^{\rm T}(k)\;+ \\ &B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)+Q(k)]C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\&[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]\delta_{ij}R_i\Gamma^{\rm T}_j(k)\\[-10pt] \end{split} $$ (17) 其中,
$\delta_{ii} = 1,\delta_{ij} = 0,\;i\ne j$ .证明. 定义
$ {\boldsymbol \phi}(k) = {\boldsymbol \theta}(k)-{\boldsymbol \theta}(k-1) $ , 则由式(7)可导出攻击信号$ {\boldsymbol \theta}(k) $ 的局部估计误差为$$ \begin{split} \widetilde{{\boldsymbol \theta}}_i(k) = \;&{\boldsymbol \phi}(k)+{\boldsymbol \theta}(k-1)-[\hat{{\boldsymbol \theta}}_i(k-1)+\Gamma_i(k)\widetilde{{\boldsymbol y}}_i(k)] = \\ &\widetilde{{\boldsymbol \theta}}_i(k-1)-\Gamma_i(k)C_i(k)\;\times \\ &[A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)+B(k)\widetilde{{\boldsymbol \theta}}_i(k-1)]\;- \\ &\Gamma_i(k)C_i(k){\boldsymbol w}(k)-\Gamma_i(k){\boldsymbol v}_i(k)+{\boldsymbol \phi}(k) = \\ &[I_p-\Gamma_i(k)C_i(k)B(k)]\widetilde{{\boldsymbol \theta}}_i(k-1)\;- \\ &\Gamma_i(k)C_i(k)A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)\;- \\ &\Gamma_i(k)C_i(k){\boldsymbol w}(k)-\Gamma_i(k){\boldsymbol v}_i(k)+{\boldsymbol \phi}(k)\\[-10pt] \end{split} $$ (18) 由式(3)和式(7)可得状态估计误差为
$$ \begin{split} \widetilde{{\boldsymbol x}}_i(k\vert k) =\;& [I_n-K_i(k)C_i(k)][A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)\;+ \\ &B(k)\widetilde{{\boldsymbol \theta}}_i(k-1)]-\Upsilon_i(k)[\hat{{\boldsymbol \theta}}_i(k)-\hat{{\boldsymbol \theta}}_i(k-1)]\;+ \\ &[I_n-K_i(k)C_i(k)]{\boldsymbol w}(k)-K_i(k){\boldsymbol v}_i(k) = \\ &[I_n-K_i(k)C_i(k)][A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)\;+ \\ &B(k)\widetilde{{\boldsymbol \theta}}_i(k-1)]+\Upsilon_i(k)[\widetilde{{\boldsymbol \theta}}_i(k)\;-\\ &\widetilde{{\boldsymbol \theta}}_i(k-1)-{\boldsymbol \phi}(k)]\;+ \\ &[I_n-K_i(k)C_i(k)]{\boldsymbol w}(k)-K_i(k){\boldsymbol v}_i(k) \\[-10pt]\end{split} $$ (19) 将式(18)代入式(19), 导出
$$ \begin{split} \widetilde{{\boldsymbol x}}_i(k\vert k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)\widetilde{{\boldsymbol x}}_i(k-1\vert k-1)+B(k)\widetilde{{\boldsymbol \theta}}_i(k-1)]\;+ \\ &[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]{\boldsymbol w}(k)\;- \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\boldsymbol v}_i(k) \\[-10pt]\end{split}$$ (20) 注意到
$ {\boldsymbol \phi}(k) $ 是未知变量, 在此假设它与其他统计变量无关, 且它的协方差为$\eta I = { {\rm{E}}}({\boldsymbol \phi}(k){\boldsymbol \phi}^{\rm T}(k))$ . 在此情况下, 调整参数$ \eta $ 可以改变所设计融合估计器对攻击信号$ {\boldsymbol \theta}(k) $ 的估计性能, 即通过$ \eta $ 实现对估计误差$ \widetilde{{\boldsymbol \theta}}_i(k) $ 的协方差进行补偿, 因而$ \eta $ 称为补偿因子. 根据${\boldsymbol w}(k)\perp\widetilde{{\boldsymbol x}}_j(k- $ $ 1\vert k-1),$ ${\boldsymbol v}_i(k)\;\perp$ $\widetilde{{\boldsymbol x}}_j\;(k-1\vert k-1),$ ${\boldsymbol \phi}_i(k)\;\perp\;\widetilde{{\boldsymbol x}}_j(k- $ $ 1\vert k- 1),$ $ {\boldsymbol w}(k)\perp\widetilde{{\boldsymbol \theta}}_i(k-1) ,$ $ {\boldsymbol v}_i(k)\perp\widetilde{{\boldsymbol \theta}}_i(k-1) ,$ ${\boldsymbol \phi}_i(k)\perp $ $ \widetilde{{\boldsymbol \theta}}_i(k-1),$ 由式(18)可得$P ^\theta_{ij}(k)$ 由下式计算$$ \begin{split} P^\theta_{ij}(k) =\;& [I_p-\Gamma_i(k)C_i(k)B(k)]P^\theta_{ij}(k-1)\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;+ \\ &\Gamma_i(k)C_i(k)A(k)P^x_{ij}(k-1)A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &\Gamma_i(k)C_i(k){\rm E}({\boldsymbol w}(k){\boldsymbol w}^{\rm T}(k))C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &\Gamma_i(k){\rm E}({\boldsymbol v}_i(k){\boldsymbol v}_j^{\rm T}(k))\Gamma^{\rm T}_j(k)+{\rm E}({\boldsymbol \phi}(k){\boldsymbol \phi}^{\rm T}(k))\;- \\ &[I_p-\Gamma_i(k)C_i(k)B(k)]\Psi^{\rm T}_{ji}(k-1)\;\times \\ &A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)-\Gamma_i(k)C_i(k)A(k)\;\times \\ &\Psi_{ij}(k-1)[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;+ \\ &\Gamma_i(k)C_i(k){\rm E}({\boldsymbol w}(k){\boldsymbol v}^{\rm T}_j(k))\Gamma^{\rm T}_j(k)\;+ \\ &\Gamma_i(k){\rm E}({\boldsymbol v}_i(k){\boldsymbol w}^{\rm T}(k))C^{\rm T}_j(k)\Gamma^{\rm T}_j(k) \\[-10pt]\end{split} $$ (21) 根据式(5), 式(21)简化为
$$ \begin{split} P^\theta_{ij}(k) = \;&[I_p-\Gamma_i(k)C_i(k)B(k)]P^\theta_{ij}(k-1)\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}+\Gamma_i(k)C_i(k)A(k)\;\times \\ &P^x_{ij}(k-1)A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;- \\ &[I_p-\Gamma_i(k)C_i(k)B(k)]\Psi^{\rm T}_{ji}(k-1)\;\times \\ &A^{\rm T}(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)-\Gamma_i(k)C_i(k)A(k)\;\times \\ &\Psi_{ij}(k-1)[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;+ \\ &\Gamma_i(k)C_i(k)Q(k)C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &\Gamma_i(k){\rm E}({\boldsymbol v}_i(k){\boldsymbol v}_j^{\rm T}(k))\Gamma^{\rm T}_j(k)+\eta I\\[-10pt] \end{split} $$ (22) 从而导出式(15)成立. 另一方面, 由式(20)可得
$ P^x_{ij}(k) $ 由下式计算$$ \begin{split} P^x_{ij}(k) = \;&[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)P^x_{ij}(k-1)A^{\rm T}(k)\;+\\ &B(k)P^\theta_{ij}(k-1)B^{\rm T}(k)\;+ \\ &A(k)\Psi_{ij}(k-1)B^{\rm T}(k)\;+\\ &B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)+ \\&{\rm E}({\boldsymbol w}(k){\boldsymbol w}^{\rm T}(k))]\;\times \\ &[I_n-(K_j(k)+\Upsilon_j(k)\Gamma_j(k))C_j(k)]^{\rm T}\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol v}_j^{\rm T}(k))\;\times \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]^{\rm T}\;+ \\ &[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &{\rm E}({\boldsymbol w}(k){\boldsymbol v}_j^{\rm T}(k))[K_j(k)+\Upsilon_j(k)\Gamma_j(k)]^{\rm T}\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol w}^{\rm T}(k))\;\times \\ &[I_n-(K_j(k)+\Upsilon_j(k)\Gamma_j(k))C_j(k)]^{\rm T} \\[-10pt]\end{split} $$ (23) 根据式(5), 式(23)转换为
$$ \begin{split} P^x_{ij}(k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)P^x_{ij}(k-1)A^{\rm T}(k)\;+\\ &B(k)P^\theta_{ij}(k-1)B^{\rm T}(k)\;+ \\ &A(k)\Psi_{ij}(k-1)B^{\rm T}(k)\;+ \\ &B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)+Q(k)]\;\times \\ &[I_n-(K_j(k)+\Upsilon_j(k)\Gamma_j(k))C_j(k)]^{\rm T}\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol v}_j^{\rm T}(k))\;\times \\ &[K_j(k)+\Upsilon_j(k)\Gamma_j(k)]^{\rm T} \end{split} $$ (24) 故式(16)成立. 进一步, 根据式(18)和式(20)导出
$ \Psi_{ij}(k) $ 由下式计算$$ \begin{split} \Psi_{ij}(k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)\Psi_{ij}(k-1)+B(k)P^\theta_{ij}(k-1)]\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;- \\ &[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)P^x_{ij}(k-1)A^{\rm T}(k)\;+\\ &B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)]\;\times \\ &C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)+[I_n-(K_i(k)\;+\\ &\Upsilon_i(k)\Gamma_i(k))C_i(k)]{\rm E}({\boldsymbol w}(k)\;\times \\ &{\boldsymbol w}^{\rm T}(k))C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)+[I_n-(K_i(k)+ \Upsilon_i(k)\;\times \\ &\Gamma_i(k))C_i(k)]{\rm E}({\boldsymbol w}(k){\boldsymbol v}^{\rm T}_j(k))\Gamma^{\rm T}_j(k)+ [K_i(k)\;+\\ &\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol w}^{\rm T}(k))C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol v}^{\rm T}_j(k))\Gamma^{\rm T}_j(k) \end{split} $$ (25) 由式(5)和式(25)导出
$$ \begin{split} \Psi_{ij}(k) =\;& [I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)]\;\times \\ &[A(k)\Psi_{ij}(k-1)+B(k)P^\theta_{ij}(k-1)]\;\times \\ &[I_p-\Gamma_j(k)C_j(k)B(k)]^{\rm T}\;- \\ &[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))C_i(k)][A(k)\;\times \\ &P^x_{ij}(k-1)A^{\rm T}(k)+B(k)\Psi^{\rm T}_{ji}(k-1)A^{\rm T}(k)]\;\times \\ &C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)+[I_n-(K_i(k)+\Upsilon_i(k)\Gamma_i(k))\;\times \\ &C_i(k)]{\rm E}({\boldsymbol w}(k){\boldsymbol w}^{\rm T}(k))C^{\rm T}_j(k)\Gamma^{\rm T}_j(k)\;+ \\ &[K_i(k)+\Upsilon_i(k)\Gamma_i(k)]{\rm E}({\boldsymbol v}_i(k){\boldsymbol v}^{\rm T}_j(k))\Gamma^{\rm T}_j(k) \end{split} $$ (26) 故式(17)成立. □
注 3. 当补偿因子
$ \eta = 0 $ 时, 意味着$ \widetilde{{\boldsymbol \theta}}_i(k) $ 的协方差信息无法得到补偿. 在此情况下, 如果攻击信号$ {\boldsymbol \theta}(k) $ 恒等于一个常值, 则所设计的融合估计器将提供满意的估计性能; 但是如果攻击信号$ {\boldsymbol \theta}(k) $ 是时变的, 则$ \eta $ 取$ 0 $ 时相应的融合估计性能将会因为无法得到协方差的补偿信息而变差. 由于$ {\boldsymbol \theta}(k) $ 是未知的, 所以补偿因子$ \eta $ 作为一个可调参数可以提高对攻击信号的估计精度. 但是如何设计补偿因子的选取准则以保证满意的融合性能非常具有挑战性, 将作为今后在攻击信号实时估计方面的重要研究方向之一.注 4. 当获得FDI攻击信号的估计值时, 攻击对系统性能的破环仍然存在. 为此, 可以采取在系统(1)中直接减去攻击信号融合估计值的方法来尽可能地降低攻击所带来的损失, 即带有补偿策略的系统(1)由如下动态方程描述
$$ \begin{split} {\boldsymbol x}(k) =\;& A(k){\boldsymbol x}(k-1)+B(k){\boldsymbol u_a}(k-1)\;+ \\&{\boldsymbol w}(k)-B(k)\hat\theta_0(k-1) = \\&A(k){\boldsymbol x}(k-1)+B(k){\boldsymbol u}(k-1)\;+ \\&{\boldsymbol w}(k)+B(k)(\theta(k-1)-\hat\theta_0(k-1)) \end{split} $$ (27) 其中,
$ \hat\theta_0(k-1) $ 是攻击信号的融合估计值. 由式(27)可知, 融合估计器精度越高, 系统所受到的破坏就越小.根据定理1和定理2, 给定遗忘因子
$\lambda_1,\cdots,\lambda_L$ 和补偿因子$ \eta $ , 则攻击信号的融合估计$ \hat{{\boldsymbol \theta}}_0(k) $ 实现的算法如下:算法1. 最优融合估计算法
步骤1. 根据式(8), 计算卡尔曼增益
$ K_i(k) $ , 然后根据式(9), 计算参数估计增益$ \Gamma_i(k) $ ;步骤2. 将步骤1计算的结果代入式(7), 得到攻击信号的局部估计
$ \hat{{\boldsymbol \theta}}_i(k) $ 和局部状态估计$ \hat{{\boldsymbol x}}_i(k\vert k) $ ;步骤3. 根据式(15)计算协方差矩阵
$ P^\theta_{ii}(k) $ 和互协方差矩阵$ P^\theta_{ij}(k) $ ;步骤4. 将步骤3计算的结果代入式(12)得到融合权重
$ W^\theta_i(k) $ ;步骤5. 将步骤2和步骤4的结果代入式(6), 计算得到融合估计
$ \hat{{\boldsymbol \theta}}_0(k) $ ;步骤6. 重复步骤1 ~ 5, 得到
$k+1 $ 时刻的最优融合估计$ \hat{{\boldsymbol \theta}}_0(k+1) $ .3. 仿真算例
算例 1. 考虑一个IEEE四路配电线路的电网系统. 本文采用互连的分布式发电器 (Distributed energy generator, DEG)模型, 4个DEG被建模为电压源. 4个DEG通过对应的公共耦合点(Point of common coupling, PCC)与主电网连接, PCC的电压为
$ {\boldsymbol v}_s = (v_1,v_2,v_3, $ $ v_4)^{\rm T} ,$ 其中$ v_i $ 为第$ i $ 个PCC电压. 为了保持DEG的正常工作, 这些PCC电压需要保持在一定的参考值. 在每个DEG和其余的电力网络之间存在耦合电感. 然后, 节点电压方程可以用下述线性状态空间模型表示为[3]$$ \dot{{\boldsymbol x}}(t) = \left[\begin{array}{cccc} \;175.9&\;\,176.8&\;511.0&\;1\;036.0 \\ \,-350.0\;\,&\;\,0&\;\,0&\;\,\;\;\,0 \\ -544.2&-474.8&-408.8&\;\;-828.8 \\ -119.7&-554.6&-968.8&-1\;077.5 \end{array}\right]{\boldsymbol x}(t) $$ (28) 其中,
$ {\boldsymbol x}(t) = {\boldsymbol v}_s(t)-{\boldsymbol v}_{\rm ref}(t) $ 为PCC状态电压偏差,$ {\boldsymbol v}_{\rm ref}(t) $ 为PCC参考电压; 针对系统(28), 当采样周期取$T_0 = $ $ 0.001\;{\rm{s}}$ 时, 则离散化系统的状态转移矩阵为$$ A_d = \left[\begin{array}{cccc} \;\;\,1.0779&-0.0691&\;\;\,0.1412&\;\;\,0.6022 \\ -0.3665&\;\;\,1.0012&-0.0417&-0.1247 \\ -0.4245&-0.2761&\;\;\,0.8447&-0.6105 \\ \;\;\,0.1515&-0.2261&-0.5207&\;\;\,0.5588 \end{array}\right] $$ 考虑到系统中存在的噪声, 得到如下离散系统方程:
$$ {\boldsymbol x}(k) = A_d{\boldsymbol x}(k-1)+{\boldsymbol w}(k) $$ (29) 注意到在没有控制信号输入的情况下, 系统(29)无法稳定. 因此, 将设计控制器以确保系统稳定. 考虑到控制输入信号受到FDI攻击, 利用两组传感器对系统进行量测, 得到如下动态方程
$$ \begin{split} {\boldsymbol x}(k) =\;& A_d{\boldsymbol x}(k-1)+B_d({\boldsymbol u}(k-1)+\\ &{\boldsymbol \theta}(k-1))+{\boldsymbol w}(k) \end{split}$$ (30) $$ {\boldsymbol y}_i(k) = C_i{\boldsymbol x}(k)+{\boldsymbol v}_i(k),\;\;i = 1,2\;\;\;\;\;$$ (31) 其中,
$B_d \;=\; [2.6890,\;-4.3035,\;-8.3410,\;-7.0725]^{\rm T},$ ${\boldsymbol u}(k) = -[5.7495,3.4187,-7.3869,8.3114]{\boldsymbol x}(k)$ 是通过极点配置设计的稳定控制器. 系统过程噪声$ {\boldsymbol w}(k) $ 的方差取$ Q = {\rm diag}\{0.1,0.2,0.2,0.1\} $ . 量测矩阵为$$ C_1 = \left[\begin{array}{ccccccc} 0&0&0&0 \\ 0&1&0&0 \\ 0&0&1&0 \\ 0&1&1&0 \end{array}\right] $$ $$ C_2 = \left[ \begin{array}{ccccccccc} 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 1 \end{array}\right] $$ 量测噪声的协方差矩阵分别取
$R_1 = {\rm diag}\{0.5,0.6, 0.3, $ $ 0.2\}$ ,$ R_2 = {\rm diag}\{0.8,0.3,0.5,0.9\} $ . 在仿真中考虑以下两种情况:情况 1. 攻击信号保持常值不变, 取
$$ {\boldsymbol \theta}(k)\equiv1 $$ (32) 令补偿因子
$ \eta = 0 $ 及遗忘因子$ \lambda = 0.95 $ . 执行算法1, 仿真结果如图1和图2所示.
图 1 算例1: 情况1中攻击信号和融合估计的轨迹Fig. 1 Example 1: The trajectories of attack signal and its fusion estimation under Case 1
图 2 算例1: 情况1中攻击信号的局部估计器与融合估计器之间的性能比较Fig. 2 Example 1: The performance comparison between local estimators and fusion estimators under Case 1由图1可知, 当遭受的攻击信号保持常值不变时, 所提出的融合算法可以很好地估计到攻击信号. 图2给出了300次Monte Carlo实验的局部估计器和融合估计器的均方误差(Mean square errors, MSEs)曲线. 由此图可以看到融合估计器的性能优于局部估计器的性能, 从而验证了基于多传感器信息融合的攻击信号估计方法优于基于单一传感器的估计方法.
情况 2. 攻击信号
$ {\boldsymbol \theta}(k) $ 是时变的, 这里假设$ {\boldsymbol \theta}(k) $ 按如下正弦函数变化$$ {\boldsymbol \theta}(k) = {\rm sin}(0.3k) $$ (33) 令补偿因子
$ \eta = 0.2 $ 和遗忘因子$ \lambda = 0.7 $ . 执行算法1, 仿真结果如图3 ~ 6所示.
图 3 算例1: 情况2中攻击信号和融合估计轨迹Fig. 3 Example 1: The trajectories of attack signal and its fusion estimation under Case 2
图 4 算例1: 情况2中攻击信号的局部估计器与融合估计器之间的性能比较Fig. 4 Example 1: The performance comparison between local estimators and fusion estimators under Case 2
图 5 算例1: 情况2中不同补偿因子下攻击信号融合估计性能的比较Fig. 5 Example 1: The comparison of fusion estimation performance of attack signal under different compensation factors under Case 2
图 6 算例1: 情况2中不同补偿因子下攻击信号融合估计的轨迹Fig. 6 Example 1: The trajectories of fusion estimation of attack signal under different compensation factors under Case 2由图3可知, 当攻击信号是时变的情况时, 所设计的分布式融合器仍可以准确地估计攻击信号. 图4给出了100次Monte Carlo实验的局部估计器和融合估计器的MSE曲线. 由此图可知, 融合估计性能明显优于局部估计性能, 进一步说明了融合策略可以提高攻击信号的估计精度. 另一方面, 当取不同的补偿因子
$ \eta $ 时, 相应的融合估计性能如图5和图6所示, 即$ \eta = 0.2 $ 和$ \eta = 0.5 $ 的融合性能明显优于$ \eta = 0 $ 时的融合估计器, 说明了当FDI攻击信号是时变的情形时, 补偿因子$ \eta $ 对融合估计性能发挥着重要的作用, 即: 不考虑补偿策略的融合方法, 相应的估计性能会变得很差.算例 2. 考虑一个由两个传感器监控的机动目标, 其运动可以由它的位置和速度矢量来描述. 定义状态向量
${\boldsymbol{x}}(k) = $ $ {\rm col}\{X_s(k), \dot X_s(k), Y_s(k), \dot Y_s(k)\}$ , 其中,$(X_s(k), Y_s(k))$ 分别为X与Y轴的位置坐标, 同时$ (\dot X_s(k), \dot Y_s(k)) $ 是对应的速度. 引入控制信号改变机动目标的速度以控制移动目标的运动, 考虑控制信号遭受FDI攻击, 机动目标的运动轨迹可由如下时变状态空间模型描述[28]:$$ \boldsymbol{x}(k)= \left[ \begin{array}{ccccccc} 1&f(k)&0&0\\ 0&1&0&0\\ 0&0&1&f(k)\\ 0&0&0&1 \end{array}\right]\boldsymbol{x}(k-1)+\boldsymbol{w}(k)\;+ $$ $$ \left[ \begin{array}{ccc} 0\\ 10f(k)\\ 0 \\ 10f(k) \end{array}\right](\boldsymbol{u}(k-1)+\boldsymbol{\theta}(k-1))\\[-20pt] $$ (34) $$ {\boldsymbol y}_i(k) = C_i{\boldsymbol x}(k)+{\boldsymbol v}_i(k),\;\;i = 1,2 \;\;\qquad\qquad\quad $$ (35) 其中,
$ f(k) $ 表示采样周期取$ f(k) = 0.9+0.1{\rm sin}(k) $ ; 控制信号$ {\boldsymbol u}(k) $ 取$ {\boldsymbol u}(k) = -[0,0.01,0,0.01]x(k-1) ,$ 表示改变机动目标的速度;$ {\boldsymbol \theta}(k) $ 表示未知的FDI攻击信号, 且取$$ {\boldsymbol \theta}(k) = {\boldsymbol \theta}(k-1)+c(k) $$ (36) 其中,
$ c(k) $ 为方差为1的零均值高斯随机信号. 系统过程噪声$ {\boldsymbol w}(k) $ 和量测噪声$ {\boldsymbol v}_1(k) $ ,$ {\boldsymbol v}_2(k) $ 的方差分别取$Q = $ $ {\rm diag}\{0.1,\;0.2,\;0.3,\;0.2\},$ $R_1 = {\rm diag}\{0.3,\;0.2\},$ $R_2 = $ ${\rm diag}\{0.2,\;0.5\}.$ 量测矩阵取[28]$$ C_1 = \left[\begin{array}{ccccccc} 0.5&1&0&0 \\ 0&1&0.9&0.6 \end{array}\right] $$ $$ C_2 = \left[\begin{array}{ccccccc} 0.9&0.8&0&0 \\ 0&1&0.5&0.1 \end{array}\right] $$ 最后, 令遗忘因子
$ \lambda = 0.6 $ , 补偿因子$ \eta_1 = 1, \eta_2 = 0 $ . 执行算法1, 仿真结果如图7 ~ 9所示.
图 7 算例2: 不同补偿因子下攻击信号融合估计的轨迹Fig. 7 Example 2: The trajectories of fusion estimation of attack signal under different compensation factors
图 8 算例2: 不同补偿因子下攻击信号融合估计性能的比较Fig. 8 Example 2: The comparison of fusion estimation performance of attack signal under different compensation factors
图 9 算例2: 攻击信号的局部估计器与融合估计器之间的性能比较Fig. 9 Example 2: The performance comparison of attack signal between local estimators and fusion estimators由图7可知, 系统参数时变的情况时, 所设计的估计算法仍可以很好地估计攻击信号. 图8基于Monte Carlo实验给出了不同补偿因子下融合估计器的MSE曲线. 图7和图8表明, 补偿因子
$ \eta = 1 $ 的情况下, 算法的性能优于$ \eta = 0 $ 的情况(无补偿的情况). 图9基于Monte Carlo实验给出了局部估计器和融合估计器的MSE曲线. 图9表明融合性能明显优于局部估计器. 仿真结果验证了算法的有效性.4. 结束语
本文研究了信息物理系统中遭受的FDI攻击信号的分布式融合估计问题. 首先, 根据FDI攻击策略对系统的影响, 将攻击信号建模为系统状态方程中的一个未知输入, 然后基于自适应卡尔曼滤波, 利用量测信息对系统状态和攻击参数进行联合估计, 从而得到攻击参数的局部变化轨迹. 在此基础上, 在分布式融合结构下, 引入补偿因子并设计了最优融合估计方法, 提高了攻击信号估计的精度. 最后, 通过例子验证了所提算法的有效性. 今后的研究方向包括如下几点: 首先, 选择合适的补偿因子对本文算法估计精度的提高有着显著的影响, 如何选择补偿因子将是未来研究的重要方向之一; 此外, 许多实际系统无法简单地利用线性模型描述, 因此将该算法推广至非线性情况的攻击信号估计具有重要的意义; 最后, 研究控制器与攻击信号融合估计器及攻击信号实时补偿策略的协同设计, 以降低攻击所带来的损失是另一个重要的研究方向.
-
图 1 算例1: 情况1中攻击信号和融合估计的轨迹
Fig. 1 Example 1: The trajectories of attack signal and its fusion estimation under Case 1
图 2 算例1: 情况1中攻击信号的局部估计器与融合估计器之间的性能比较
Fig. 2 Example 1: The performance comparison between local estimators and fusion estimators under Case 1
图 3 算例1: 情况2中攻击信号和融合估计轨迹
Fig. 3 Example 1: The trajectories of attack signal and its fusion estimation under Case 2
图 4 算例1: 情况2中攻击信号的局部估计器与融合估计器之间的性能比较
Fig. 4 Example 1: The performance comparison between local estimators and fusion estimators under Case 2
图 5 算例1: 情况2中不同补偿因子下攻击信号融合估计性能的比较
Fig. 5 Example 1: The comparison of fusion estimation performance of attack signal under different compensation factors under Case 2
图 6 算例1: 情况2中不同补偿因子下攻击信号融合估计的轨迹
Fig. 6 Example 1: The trajectories of fusion estimation of attack signal under different compensation factors under Case 2
图 7 算例2: 不同补偿因子下攻击信号融合估计的轨迹
Fig. 7 Example 2: The trajectories of fusion estimation of attack signal under different compensation factors
图 8 算例2: 不同补偿因子下攻击信号融合估计性能的比较
Fig. 8 Example 2: The comparison of fusion estimation performance of attack signal under different compensation factors
-
[1] Johansson K H, Pappas G J, Tabuada P, Tomlin C J. Guest editorial special issue on control of cyber-physical systems. IEEE Transactions on Automatic Control, 2014, 59(12): 3120-3121 doi: 10.1109/TAC.2014.2363896 [2] Fink J, Ribeiro A, Kumar V. Robust control for mobility and wireless communication in cyber-physical systems with application to robot teams. Proceedings on the IEEE, 2012, 100(1): 164-178 doi: 10.1109/JPROC.2011.2161427 [3] Li H, Lai L, Poor H V. Multicast routing for decentralized control of cyber physical systems with an application in smart grid. IEEE Journal on Selected Areas in Communications, 2012, 30(6): 1097-1107 doi: 10.1109/JSAC.2012.120708 [4] Langner R. Stuxnet: dissecting a Cyberwarfare Weapon. IEEE Security & Privacy, 2011, 9(3): 49-51 [5] Hu L, Wang Z, Han Q, Liu X. State estimation under false data injection attacks: Security analysis and system protection. Automatica, 2018, 87: 176-183 doi: 10.1016/j.automatica.2017.09.028 [6] 艾纳安德尔. 信息物理系统和工业4.0. 智能制造, 2015, 9: 10-12 doi: 10.3969/j.issn.1671-8186.2015.07.002Einar Andel. Cyber-physical systems and industry 4.0. Intelligent Manufacturing, 2015, 9: 10-12 doi: 10.3969/j.issn.1671-8186.2015.07.002 [7] Pasqualetti F, Dorfler F, Bullo F. Attack detection and identification in cyber-physical systems. IEEE Transactions on Automatic Control, 2013, 58(11): 2715-2729 doi: 10.1109/TAC.2013.2266831 [8] Pang Z, Liu G, Zhou D, Hou F, Sun D. Two-channel false data injection attacks against output tracking control of networked systems. IEEE Transactions on Industrial Electronics, 2016, 63(5): 3242-3251 doi: 10.1109/TIE.2016.2535119 [9] Chen B, Ho D W C, Zhang W, Yu L. Distributed dimensionality reduction fusion estimation for cyber-physical systems under DoS attacks. IEEE Transactions on Systems, Man, and Cybernetics: Systems, 2019, 49(2): 455-468 doi: 10.1109/TSMC.2017.2697450 [10] Lai S Y, Chen B, Li T X, Yu L. Packet-based state feedback control under DoS attacks in cyber-physical systems. IEEE Transactions on Circuits and Systems II: Express Briefs, 2019. 66(8): 1421−1425 [11] Chen B, Ho D W C, Hu G, Yu L. Secure fusion estimation for bandwidth constrained cyber-physical systems under replay attacks. IEEE Transactions on Cybernetics, 2018, 48(6): 1862-1876 doi: 10.1109/TCYB.2017.2716115 [12] Ding D, Han Q, Xiang Y, Ge X, Zhang X. A survey on security control and attack detection for industrial cyber-physical systems. Neurocomputing, 2018, 275: 1674-1683 doi: 10.1016/j.neucom.2017.10.009 [13] Liu Y, Ning P, Reiter M K. False data injection attacks against state estimation in electric power grids. ACM Transactions on Information and System Security, 2011, 14(1): 1-33 [14] Ye N, Zhang Y, Borror C M. Robustness of the Markov-chain model for cyber-attack detection. IEEE Transactions on Reliability, 2004, 53(1): 116-123 doi: 10.1109/TR.2004.823851 [15] Chorppath A K, Alpcan T, Boche H. Bayesian mechanisms and detection methods for wireless network with malicious users. IEEE Transactions on Mobile Computing, 2016, 15(10): 2452-2465 doi: 10.1109/TMC.2015.2505724 [16] Kailkhura B, Han Y S, Brahma S, Varshney P K. Distributed Bayesian detection in the presence of byzantine data. IEEE Transactions on Signal Processing, 2015, 63(9): 5250-5263 [17] Kailkhura B, Han Y S, Brahma S, Varshney P K. Asymptotic analysis of distributed Bayesian detection with byzantine data. IEEE Signal Processing Letters, 2015, 22(5): 608-612 doi: 10.1109/LSP.2014.2365196 [18] Rawat A S, Anand P, Chen H, Varshney P K. Collaborative spectrum sensing in the presence of byzantine attacks in cognitive radio networks. IEEE Transactions on Signal Processing, 2011, 59(2): 774-786 doi: 10.1109/TSP.2010.2091277 [19] Manandhar K, Cao X, Hu F, Liu Y. Detection of faults and attacks including false data injection attack in smart grid using kalman filter. IEEE Transactions on Control of Network Systems, 2014, 1(4): 370-379 doi: 10.1109/TCNS.2014.2357531 [20] Mo Y, Chabukswar R, Sinopoli B. Detecting integrity attacks on SCADA systems. IEEE Transactions on Control Systems Technology, 2014, 22(4): 1396-1407 doi: 10.1109/TCST.2013.2280899 [21] Rawat D B, Bajracharya C. Detection of false data injection attacks in smart grid communication systems. IEEE Signal Processing Letters, 2015, 22(10): 1652-1656 doi: 10.1109/LSP.2015.2421935 [22] Liu L, Esmalifalak M, Ding Q, Emesih V A, Han Z. Detecting false data injection attacks on power grid by sparse optimization. IEEE Transactions on Smart Grid, 2014, 5(2): 612-621 doi: 10.1109/TSG.2013.2284438 [23] Deng R, Xiao G, Lu R. Defending against false data injection attacks on power system state estimation. IEEE Transactions on Industrial Informatics, 2017, 13(1): 198-207 doi: 10.1109/TII.2015.2470218 [24] Huang Y, Tang J, Cheng Y, Li H, Campbell K A, Han Z. Real-time detection of false data injection in smart grid networks: An adaptive CUSUM method and analysis. IEEE Systems Journal, 2016, 10(2): 532-543 doi: 10.1109/JSYST.2014.2323266 [25] Zhang Q, Basseville M. Statistical detection and isolation of additive faults in linear time-varying systems. Automatica, 2014, 50(10): 2527-2538 doi: 10.1016/j.automatica.2014.09.004 [26] Zhang Q. Adaptive Kalman filter for actuator fault diagnosis. Automatica, 2018, 93: 333-342 doi: 10.1016/j.automatica.2018.03.075 [27] 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 [28] Chen B, Ho D W C, Zhang W, Yu L. Networked Fusion Estimation With Bounded Noises. IEEE Transactions on Automatic Control, 2017, 62(10): 5415-5421 doi: 10.1109/TAC.2017.2696746 -
下载:
下载:
计量
- 文章访问数: 1123
- HTML全文浏览量: 181
- PDF下载量: 277
- 被引次数: 0
