方波触发勘探与开发的粒子群优化算法

邓志诚; 孙辉; 赵嘉; 王晖; 吕莉; 谢海华

doi:10.16383/j.aas.c190842

方波触发勘探与开发的粒子群优化算法

doi: 10.16383/j.aas.c190842

邓志诚^{1, 2,},
孙辉^{1, 2, 3,},
赵嘉^{1, 2, 3,},
王晖^{1, 2, 3,},
吕莉^{1, 2, 3,},
谢海华^{1, 2,}

1.
南昌工程学院信息工程学院南昌 330099
2.
江西省水信息协同感知与智能处理重点实验室南昌 330099
3.
鄱阳湖流域水工程安全与资源高效利用国家地方联合工程实验室南昌 330099

基金项目: 国家自然科学基金(52069014, 62066030, 61663029, 51669014, 61663028), 江西省教育厅科技项目(GJJ201915, GJJ180940), 江西省研究生创新专项资金项目(YC2018-S422)资助

详细信息

作者简介:
邓志诚：南昌工程学院硕士研究生. 主要研究方向为群智能算法.E-mail: deng_zc215@163.com

孙辉：南昌工程学院教授. 1988年获清华大学硕士学位, 2002年获南昌大学博士学位. 主要研究方向为群智能算法, 粗糙集和变分不等原理. 本文通信作者.E-mail: sun_hui2006@163.com

赵嘉：南昌工程学院教授. 主要研究方向为群智能算法与数据挖掘.E-mail: zhaojia925@163.com

王晖：南昌工程学院教授. 主要研究方向为群智能算法与水资源优化.E-mail: huiwang@nit.edu.cn

吕莉：南昌工程学院教授. 主要研究方向为群智能算法与目标跟踪.E-mail: lvli623@163.com

谢海华：南昌工程学院硕士研究生. 主要研究方向为群智能算法.E-mail: pxlh_xhh@163.com

计量
- 文章访问数: 540
- HTML全文浏览量: 77
- PDF下载量: 161
- 被引次数: 9
出版历程
- 收稿日期: 2019-12-11
- 录用日期: 2020-05-03
- 网络出版日期: 2022-12-23
- 刊出日期: 2022-12-23

Particle Swarm Optimization With Square Wave Triggered Exploration and Exploitation

DENG Zhi-Cheng^{1, 2
,},
SUN Hui^{1, 2, 3
,},
ZHAO Jia^{1, 2, 3
,},
WANG Hui^{1, 2, 3
,},
LV Li^{1, 2, 3
,},
XIE Hai-Hua^{1, 2
,}

1.
School of Information Engineering, Nanchang Institute of Technology, Nanchang 330099
2.
Jiangxi Province Key Laboratory of Water Information Cooperative Sensing and Intelligent Processing, Nanchang 330099
3.
National-local Joint Engineering Laboratory of Water Engineering Safety and Effective Utilization of Resources in Poyang Lake Area, Nanchang 330099

Funds: Supported by National Natural Science Foundation of China (52069014, 62066030, 61663029, 51669014, 61663028), Science and Technology Research Project of Education Department of Jiangxi Province (GJJ201915, GJJ180940), and Innovation Fund Designated for Graduate Students of Jiangxi Province (YC2018-S422)

More Information

Author Bio:
DENG Zhi-Cheng　Master stude-nt at Nanchang Institute of Technology. His main research interest is swarm intelligence algorithm

SUN Hui　Professor at Nanchang Institute of Technology. He received his master degree from Tsinghua University in 1988 and received his Ph.D. degree from Nanchang University in 2002, respectively. His research interest covers swarm intelligent algorithm, rough sets and variational inequality principle. Corresponding author of this paper

ZHAO Jia　Professor at Nanch-ang Institute of Technology. His research interest covers swarm intelligent algorithm and data mining

WANG Hui　Professor at Nanch-ang Institute of Technology. His research interest covers swarm intelligent algorithm and water resources optimization

LV Li　Professor at Nanchang Institute of Technology. Her research interest covers swarm intelligent algorithm and single target track

XIE Hai-Hua　Master student at Nanchang Institute of Technology. His main research interest is swarm intelligence algorithm

摘要

摘要: 在粒子群优化算法中, 当勘探时间持续过长, 将可能导致种群在解空间过度徘徊; 种群在开发阶段陷入局部最优后, 难以再次进行全局勘探. 针对上述问题, 提出方波触发勘探与开发的粒子群优化算法. 依据方波的周期特性, 在前半个周期内使用标准粒子群优化算法执行全局勘探, 后半个周期使用改进的更新公式执行局部开发. 经过实验验证, 在方波触发机制下, 通过为粒子提供多变步长, 可达到周期性触发勘探与开发的目的. 使用多类型测试函数, 将该算法与新改进粒子群算法、改进人工蜂群算法、改进差分算法在30、50和100维下比较, 实验结果表明, 该算法在收敛速度和精度上更具优势.
- 粒子群优化算法 /
- 方波触发机制 /
- 勘探与开发 /
- 全局优化
Abstract: In the particle swarm optimization, when the exploration time continues to be too long, it may cause the particles to be excessively lingered in the solution space. After the population has fallen into local optimum during the exploitation stage, it is difficult to conduct global exploration again. To solve these two problems, an improved particle swarm optimization called particle swarm optimization with square wave triggered exploration and exploitation is proposed. According to the periodic characteristics of the square wave, the exploration is performed using the standard particle swarm optimization equation in the first half of the period, and the exploitation is performed using the improved update equation in the second half of the period. After experimental verification, under the square wave triggered mechanism, by providing variable step sizes to particles, the purpose of triggering exploration and exploitation periodically is achieved. The square wave triggered exploration and exploitation was compared with other state-of-the-art intelligent algorithms, numerical experiment results on several types of benchmark functions of dimensions 30, 50 and 100 showed that the square wave triggered exploration and exploitation had more advantages in terms of accuracy and convergence speed.
- Particle swarm optimization algorithm /
- square wave triggered mechanism /
- exploration and exploitation /
- global optimization

HTML全文

智能电网是一种新型的电网, 它采用先进的通信网络技术和控制技术来支持更高效的能源安全传输和分配. 然而, 由于智能电网系统的复杂性和开放性, 智能电网中进行数据交换的信息网络成为易受到恶意攻击的对象^[1-2]. 例如, 2016年, 黑客攻击乌克兰国家电力部门致使国内发生了一次大规模的停电事件^[3], 造成严重经济损失. 因此, 智能电网的攻击检测研究具有重要意义.

隐蔽假数据攻击是目前恶意攻击的典型代表, 攻击者对传感器节点注入精心设计的错误数据, 接收错误数据的控制中心, 继而做出错误决策破坏系统的稳定性^[4]. 在文献[5]中, 拒绝服务攻击旨在中断电力网络通信信道的可用性. 文献[6]构建了一种对智能电网中以完整性和可用性为目标的攻击分类方法. 文献[7]解决了在电力网络中攻击检测和拓扑隔离的问题. 文献[8]提出了一种算法来识别要操作的智能电表的最优数目, 从而找到最优攻击策略, 其目的是干扰电网系统的状态估计, 影响其稳定性. 文献[9]设计了具有隐蔽特性的虚假数据攻击, 它可以使攻击前后残差基本不变, 因此, 基于卡方检测器的检测技术是无效的. 近年来, 隐蔽性攻击检测成为了研究热点之一.

针对智能电网遭受虚假数据注入攻击的检测问题, 近年来有了很多成果. 文献[10]中, 以电机的电压模型为研究对象, 使用卡尔曼滤波(Kalman filter, KF)技术得到该系统的残差序列, 通过计算攻击前后残差序列之间的Bhattacharyya距离来判断系统中是否存在攻击. 该文献的不足之处有两点: 1) Bhattacharyya距离无法辨别虚假数据注入攻击前后两个残差序列的相似性, 因为虚假数据注入攻击前后残差保持不变; 2)该电压模型是线性的, 对于实际中存在的大多数的非线性系统, 该检测算法是失效的. 文献[11]使用卡尔曼滤波技术得到系统的状态, 从系统状态角度考虑设计检测算法, 提出欧几里德检测方法, 该方法可以检测隐蔽假数据注入攻击. 该文献的不足之处是, 系统方程对于噪声是线性的, 即该噪声是加性噪声, 当系统方程对于噪声是非线性的时候, 噪声经过非线性变换, 不再服从高斯分布, 也就无法设计阈值, 无法检测攻击. 文献[12]中, 攻击者可能会缓慢改变多个传感器的测量值, 因此上述统计异常检测不会检测到个别受损的测量值, 所以提出检测思想, 这些测量值组合起来会导致状态变量远离其真值, 然而, 文中作出假设, 控制中心收集到的测量值都是服从高斯分布的, 并基于此性质设计了双边假设检验检测攻击. 该文献的不足之处是, 电网系统中的参数仅会在特定情况下服从高斯分布, 而这种情况并不常见, 且建模时没有考虑噪声, 所以该文中的检测算法局限性较大.

因此, 根据虚假数据注入攻击的特性, 考虑系统非线性和噪声统计特性未知情况, 本文提出一种基于自适应平方根无迹卡尔曼滤波器 (Unscented Kalman filter, UKF) 的智能电网隐蔽假数据攻击检测方法. 该算法可以有效地、稳定地应用于非线性系统中对系统状态作出估计, 依据系统状态设计检测算法检测攻击, 并从检测指标的角度与现有算法进行对比. 最后进行仿真实验, 实验结果证明, 所提出的基于该算法的攻击检测方法可以准确地给出状态估计值, 从而检测出隐蔽攻击.

1. 电网模型、攻击特性及问题描述

1.1 电网模型

考虑一个由3条互连总线组成的电力网络, 与之对应的拓扑图$G(V,\varepsilon )$, 如图1所示. $V = \left\{ i \right\}_1^N$表示节点集, $i \in V$与总线$i$相对应, $\varepsilon \subseteq V \times V$是图$G$的边集. 对于总线$i \in V$的相角, 建立一个连续时间二阶微分动力学方程, 称为摆动方程^[7]

图 1 3总线电网模型

Fig. 1 3-bus grid model

下载: 全尺寸图片幻灯片

$${m_i}{\ddot \delta _i}(t) + {d_i}{\dot \delta _i}(t) = {u_i}(t) + {\eta _i}(t) - \sum\limits_{j \in {N_i}} {{P_{ij}}(t)} $$

(1)

式中, ${\delta _i}$表示总线的相角, ${m_i}$为转动惯量, ${d_i}$为阻尼系数, ${u_i}$为控制输入信号, ${\eta _i}$是过程噪声, ${P_{ij}}$是从总线$i$到$j$的有功功率流

$${P_{ij}}(t) = {V_i}{V_j}{b_{ij}}\sin ({\delta _i}(t) - {\delta _j}(t))$$

式中, $V$表示总线电压, ${b_{ij}}$表示总线$i,j$之间的电纳.

引入角速度$\omega $, 将上式改写为一阶微分方程组的形式, 即

$$\begin{split} & {{\dot \delta }_i}(t) = {\omega _i}(t) \\ &{{\dot \omega }_i}(t) = - \frac{1}{{{m_i}}}({d_i}{\omega _i}(t) + {u_i} + {\varpi _i}(t) - \sum\limits_{j \in {N_i}} {{P_{ij}}(t)} ) \end{split} $$

(2)

每一条总线都安装有传感器, 将总线的数据传输回控制中心, 系统输出方程为

$${\boldsymbol{y}_{i}} = \left[ {\begin{array}{*{20}{c}} 1&0\\ 0&1 \end{array}} \right]{\boldsymbol{x}_{i}} + \boldsymbol{\rho }(t)$$

(3)

式中, ${\boldsymbol{y}_{i}}$表示总线${i}$的输出向量, $\boldsymbol{\rho }$是测量噪声.

针对上述系统, 本文采用欧拉离散化方法^[13], 即对于一个非线性系统

$$\dot{\boldsymbol x} = f(\boldsymbol{x}) + B\tau + \varpi $$

离散化后的形式为

$${\boldsymbol{x}_{k}} = {\boldsymbol{x}_{{k} - 1}} + \tau f({\boldsymbol{x}_{{k} - 1}}) + \tau B\boldsymbol{u} + \tau \varpi + {\rm{O}}({\tau ^2})$$

其中

$${\rm{O}}({\tau ^2}) = \dfrac{{{\tau ^2}}}{2}\dot f(\boldsymbol{x}) + \dfrac{{{\tau ^3}}}{3}{f^{(2)}}(\boldsymbol{x}) + \cdots + \dfrac{{{\tau ^{{{{k}}}}}}}{{k!}}{f^{({{{k}}} - 1)}}(\boldsymbol{x})$$

采样间隔值$\tau $取值充分小以保证高阶项${\rm{O}}({\tau ^2})$忽略不计, 本文中$\tau = 0.1\;{\rm{s}}$.

对式 (2) 使用上述欧拉离散化算法对该系统离散化, 得到

$$\left\{\begin{aligned} &{\boldsymbol{x}_{{k}}} = {f_{{o}}}({\boldsymbol{x}_{{k} - 1}}) + {B_o}\boldsymbol{u} + {\boldsymbol{\eta }_{k - 1}}\\\\ &{\boldsymbol{y}_{{k}}} = H{\boldsymbol{x}_{{k}}} + {\boldsymbol{\rho }_k}\\ \end{aligned} \right. $$

(4)

其中, ${f_o}( \cdot )$表示关于$x$的非线性函数, 将${\boldsymbol{x}_{{k} - 1}} + \tau {{f}}({\boldsymbol{x}_{{k} - 1}})$合并在一起, ${B_o} = \tau B$, $\eta = \tau \varpi $.

假设$\eta $和$\rho $是互不相关的非零均值高斯白噪声, 统计特性为

$$\left\{ {\begin{array}{*{20}{l}} {E({{{\boldsymbol{\eta }}}_k}) = {{\boldsymbol{q}}},\;\;{{ }}{\rm{cov}} ({{{\boldsymbol{\eta }}}_k}{{\boldsymbol{\eta }}}_j^{\rm{T}}) = {Q_k}{\delta _{kj}}} \\ {E({{{\boldsymbol{\rho }}}_k}) = {{\boldsymbol{r}}},\;\;{{ }}{\rm{cov}} ({{{\boldsymbol{\rho }}}_k}{{\boldsymbol{\rho }}}_j^{\rm{T}}) = {R_k}{\delta _{kj}}} \\ {E({\eta _k}\rho _j^{\rm{T}}) = 0} \end{array}} \right.$$

(5)

式中, Q是非负定矩阵, R是正定矩阵, ${\delta _{kj}}$是$Krone cher{\text{-}} \delta$函数, 即

$${\delta _{kj}} = \left\{ {\begin{aligned} &{0{{ }},\;k \ne j}\\ &{1{{ }},\;k = j} \end{aligned}} \right.$$

将$\eta $和$\rho $改写为${\eta _{k}} = \boldsymbol{q} + {\mu _{k}}$, ${\rho _{{k}}} = \boldsymbol{r} + {\boldsymbol{v}_{{k}}}$, 有

$$\left\{ {\begin{array}{*{20}{l}} {{\boldsymbol{x}_{k}} = {f_o}({\boldsymbol{x}_{{k} - 1}}) + {B_o}\boldsymbol{u} + \boldsymbol{q} + {\boldsymbol{\rho }_{k - 1}}}\\ {{\boldsymbol{y}_{{k}}} = H{\boldsymbol{x}_{{k}}} + \boldsymbol{r} + {\boldsymbol{v}_{k - 1}}} \end{array}} \right.$$

(6)

式中, ${\mu _k}$ 和${{\boldsymbol{v}}_k}$是互不相关的零均值高斯白噪声. ${\boldsymbol{x}_{{k}}} = {[{\delta _1} \quad {\delta _2}\quad {\delta _3}\quad {\omega _1}\quad {\omega _2}\quad {\omega _3}]^{\rm{T}}}$, ${\boldsymbol{y}_{{k}}} \in {\bf{R}}^{6}$.

1.2 攻击描述和问题提出

假设一个恶意的第三方想要破坏第1.1节中描述的系统的完整性. 本文主要考虑虚假数据注入攻击.

假设攻击者具有系统知识, 知道系统矩阵, 控制矩阵, 测量矩阵, 可以控制一系列系统中的传感器数据.

考虑假数据注入攻击模型描述为

$${\boldsymbol{z}_{a}}(k) = \left\{ {\begin{aligned} &{H{\boldsymbol{x}_{{k}}} + \xi (k),\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;k\not \subset T}\\ &{H{\boldsymbol{x}_{{a,k}}} + \xi (k) + {\boldsymbol{y}_{{a}}}(k),\;\;\;\;\;k \subset T} \end{aligned}} \right.$$

(7)

式中, $\varGamma = {{\rm{diag}}}\{{{\gamma _1},\cdots,{\gamma _n}}\}$代表传感器选择矩阵, T是攻击的时间范围, 当${\gamma _i} = 1$时, 代表第$i$个传感器遭受攻击, 否则${\gamma _i} = 0$. 且${\boldsymbol{y}_{{a}}}(k)$是攻击者精心构建的攻击向量序列, 攻击框图如图2所示.

图 2 系统遭受攻击框图

Fig. 2 block diagram of system under attack

下载: 全尺寸图片幻灯片

在隐蔽假数据注入攻击情况下, 在所考虑的攻击模型中, 攻击者构造攻击序列注入传感器, 在该序列下, 状态差值最终将发散到$\infty $, 而不会触发卡方检测器. 该攻击序列满足以下性质^[9]:

$$\begin{split} &\mathop {\lim }\limits_{k \to T} \left\| {{\boldsymbol{x}_{{a}}}(k) - \boldsymbol{x}(k)} \right\| > \alpha, \alpha > 0\\ &\mathop {\lim }\limits_{k \to T} \left\| {{\boldsymbol{r}_{{a}}}(k)} \right\| \le \varepsilon,\;\;\;\;\;\;\;\;\;\;\;\;\; 0 < \varepsilon < 1 \end{split}$$

(8)

其中, ${\boldsymbol{x}_{{a}}}({{k}})$和${\boldsymbol{r}_{{a}}}({{k}})$是遭受攻击时系统的状态变量和残差. 对攻击者克服检测机制的隐蔽性进行分析^[4]. 由输出方程可得

$$\begin{split} {\boldsymbol{r}_a} =& \left\| {(\boldsymbol{y} + \boldsymbol{a}) - H(\hat{\boldsymbol x} + \boldsymbol{c})} \right\| =\\ &\left\| {(\boldsymbol{y} - H\hat{\boldsymbol x}) + (\boldsymbol{a} - H\boldsymbol{c})} \right\|\leq \\ & \left\| {\boldsymbol{y} - H\hat{\boldsymbol x}} \right\| + \left\| {\boldsymbol{a} - H\boldsymbol{c}} \right\| = \boldsymbol{r} + {\tau _a} \end{split}$$

(9)

式中, $\boldsymbol{a}$表示攻击信号, $\hat{\boldsymbol x}$表示状态估计值, $\boldsymbol{c}$表示状态变化量. 由式(9)知, 当$\boldsymbol{a} = H\boldsymbol{c}$时, ${\boldsymbol{r}_{{a}}} = \boldsymbol{r}$, 从输出残差角度考虑, 实现了隐蔽性.

由上述分析可知, 从输出残差角度无法实现攻击检测. 所以, 需要借助滤波器从系统状态入手, 利用状态差值检测攻击信号. 本文考虑通过基于自适应平方根无迹卡尔曼滤波器的方法, 从系统状态入手来检测此种攻击.

2. 自适应平方根无迹卡尔曼滤波

本文提出将噪声估计环节^[14]加入到无迹卡尔曼滤波算法中, 实现对于噪声统计特性的在线估计, 同时改变标准UKF中状态误差协方差矩阵的迭代方式来保证滤波器的稳定性.

2.1 平方根无迹变换

本小节给出平方根无迹变换的实现方式.

步骤 1. 利用${{k}} - 1$时刻的估计状态及状态误差协方差矩阵来计算sigma采样点$\xi {}_i{{ }}$, 并给出其权值${W_{m,i}}$和${W_{c,i}}$.

$$\left\{ {\begin{array}{*{20}{l}} {{\xi _{\boldsymbol{i},k - 1}} = {{\hat{\boldsymbol x}}_{k - 1}},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;i = 1}\\ {{\xi _{\boldsymbol{i},k - 1}} = {{\hat{\boldsymbol x}}_{k - 1}} + \gamma {{({S_{k - 1}})}_i},\;\;\;i = 2, \cdots ,n + 1}\\ {{\xi _{\boldsymbol{i},k - 1}} = {{\hat{\boldsymbol x}}_{k - 1}} - \gamma {{({S_{k - 1}})}_i},\;\;\;i = n + 2, \cdots ,2n + 1} \end{array}} \right.$$

式中, $n$代表状态的维数, $\gamma = \sqrt {n + \lambda } $. 这里不再借助${\hat P_{k - 1}}$构造采样点, 而用${S_{k - 1}}$. ${S_{k - 1}}$是对状态误差协方差矩阵$P$进行QR分解得到的$P$的平方根, QR分解返回的是$P$的Cholesky分解的平方根的转置, 即${(\sqrt {(n + \lambda )P} )^{\rm{T}}}$, 故${({\boldsymbol{S}_{k - 1}})_i}$等于$(\sqrt {(n + \lambda )P} )_i^{\rm{T}}$, 这两种sigma点的选取方法在数学上是等价的.

$${W_{i,m}} = \left\{ {\begin{array}{*{20}{l}} {\frac{\lambda }{{n + \lambda }}{{ }},\;\;\;\;\;i = 1} \\\\ {\frac{1}{{2(n + \lambda )}},\;{{ }}i = 2, \cdots ,2n + 1} \end{array}} \right.\;\;\;\;\;\;\;\;$$

$$\;\;\;\;\;\;\;{W_{i,c}} = \left\{ {\begin{array}{*{20}{l}} {\frac{\lambda }{{n + \lambda }} + {{(1}} + \beta - {\alpha ^2}{{) }},i = 1} \\\\ {\frac{1}{{2(n + \lambda )}}{{ }},\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;i = 2, \cdots ,2n + 1} \end{array}} \right.$$

式中, $\lambda = {\alpha ^2}(n + \kappa ) - n$用来控制采样点之间的距离, 一般${10^{ - 4}} \leq \alpha \leq 1$, $\beta = 2$, $n + \kappa = 3$^[8].

步骤 2. 计算sigma点通过非线性函数的传播结果

$$\begin{split} &{\tau _i} = f({\xi _i}),\;\;\;i = 1, \cdots ,2n + 1\\ &\boldsymbol{z} = \sum\limits_{i = 1}^{2n + 1} {{W_{m,i}}{\boldsymbol{\tau }_i}} \\ &{P_x} = \sum\limits_{i = 1}^{2n + 1} {{W_{c,i}}({\boldsymbol{\tau }_i} - \boldsymbol{z}){{({\boldsymbol{\tau }_i} - \boldsymbol{z})}^{\rm{T}}}} \end{split}$$

(10)

式(10)是无迹变换中的计算误差协方差矩阵的公式, 平方根无迹变换将式(10)换为

$$\begin{split}{S_z} =\;& \boldsymbol{qr}(\sqrt {{\boldsymbol{W}_{c,1}}} ({\tau _1} - \boldsymbol{z})\times \cdots\times \\ &\sqrt {{\boldsymbol{W}_{c,2{n} + 1}}} ({\tau _{2n + 1}} - \boldsymbol{z}))\end{split}$$

(11)

2.2 噪声估计方法

假设噪声是非零均值的高斯噪声, 本文使用的是Sage-Husa噪声估计器, 通过极大后验估计原理, 获得次优噪声估计值, 噪声估计部分可参考文献[14].

对于系统(6), 系统方程是非线性的, 输出方程是线性的, 所以噪声估计器为

$${\hat{\boldsymbol q}_{k + 1}} = \frac{1}{k}\left[(k - 1){\hat{\boldsymbol q}_k} + {\hat{\boldsymbol x}_{k + 1}} - \sum\limits_{i = 1}^{2n + 1} {W_i^m} {f_o}({\xi _{i,k}})\right]\quad$$

$$\begin{split} {{\hat Q}_{k + 1}} =& \dfrac{1}{k}\big[({k} - 1){{\hat Q}_k} + {K_{k + 1}}{\boldsymbol{\varepsilon }_{k + 1}}{\boldsymbol{\varepsilon }_{k + 1}}K_{k + 1}^{\rm{T}} + {P_{k + 1}}-\\ & \sum\limits_{i = 1}^{2n + 1} {W_i^c} ({\boldsymbol{\tau }_{i,k + 1\left| k \right.}} - {{\hat{\boldsymbol x}}_{k + 1\left| k \right.}})\times\\ &{({\boldsymbol{\tau }_{i,k + 1\left| k \right.}} - {{\hat{\boldsymbol x}}_{k + 1\left| k \right.}})^{\rm{T}}}\big] \end{split}$$

$$\begin{split} &{{\hat{\boldsymbol r}}_{k + 1}} = \dfrac{1}{k}\left[(k - 1){{\hat{\boldsymbol r}}_k} + {\boldsymbol{y}_{k + 1}} - H{{\hat{\boldsymbol x}}_{k + 1\left| k \right.}}\right]\\ &{{\hat{\boldsymbol R}}_{k + 1}} = \dfrac{1}{k}\left[(k - 1){{\hat{\boldsymbol R}}_k} + {\boldsymbol{\varepsilon }_{k + 1}}\boldsymbol{\varepsilon }_{k + 1}^{\rm{T}} - H{P_{k + 1\left| k \right.}}H_k^{\rm{T}}\right] \end{split}$$

式中, ${\boldsymbol{\varepsilon }_k} = {\boldsymbol{y}_k} - {\hat{\boldsymbol y}_{k\left| {k - 1} \right.}}$代表残差, $K$代表滤波器增益矩阵.

2.3 自适应平方根无迹卡尔曼滤波设计

本节用到了MATLAB中的3个函数, 分别为: qr函数, cholupdate函数, diag函数.

1) $qr\;( \cdot )$ 表示QR分解, $A \in {\mathbb{{\boldsymbol{{\rm{R}}}}}^{m \times l}}$, $[Q,R] = $ $qr(A)$, 该函数生成$m \times l$上三角形矩阵R和$m \times m$酉矩阵Q, 从而, $A = Q\times R$.

2) $R{{1}} = cholupdate(R,x)$, 返回$A + xx'$的上三角Cholesky因子, $x$是具有合适长度的一个列向量, 其中, $R = chol(A)$是A的原始Cholesky分解因子.

3) ${\rm{diag}}\{ \cdot \}$函数, $D = {\rm{diag}}\{v\}$, 返回以向量v的元素为主对角线的对角矩阵. $x ={\rm{ diag}}\{A\}$, 返回A的主对角线元素的列向量. 对于式(22)中最外侧的${\rm{diag}} \{\cdot \}$, 返回一个对角矩阵, 主对角线元素是每一次迭代得到的噪声矩阵主对角线元素的平方根的实时估计值.

自适应平方根无迹卡尔曼滤波(Adaptive square-root UKF, ASRUKF)算法的具体步骤如下.

步骤 1.初始化$(k\;=\;1):$

$$\left\{ {\begin{array}{*{20}{l}} {{{\hat{\boldsymbol x}}_1} = {\rm{E}}({\boldsymbol{x}_1})}\\ {{S_1} = chol\left\{ {\rm{E}}\left(({\boldsymbol{x}_1} - {{\hat{\boldsymbol x}}_1}){{({\boldsymbol{x}_1} - {{\hat{\boldsymbol x}}_1})}^{\rm{T}}}\right)\right\} }\\ {\sqrt {{{\hat Q}_1}} = {S_1}}\\ {{{\hat{\boldsymbol q}}_1} = 0} \end{array}} \right.$$

(12)

步骤 2. 对于$k\;=\;2,3,\cdots,$进行迭代

a)时间更新

计算sigma点, 构造矩阵

$${\boldsymbol{\xi }_{k - 1}} = \left[{\hat{\boldsymbol x}_{k - 1}}{\hat{\boldsymbol x}_{k - 1}} + \gamma {({S_{k - 1}})_i}{\hat{\boldsymbol x}_{k - 1}} - \gamma {({S_{k - 1}})_i}\right]$$

(13)

$$\left\{ {\begin{array}{*{20}{l}} {{\tau _{i,k\left| {k - 1} \right.}} = f({\boldsymbol{\xi }_{i,k - 1}}) + \boldsymbol{q}}\\ {\boldsymbol{x} = \sum\limits_{i = 1}^{2n + 1} {{W_{m,i}}{\tau _{i,k\left| {k - 1} \right.}} + \hat{\boldsymbol q}} }\\ {{S_{k\left| {k - 1} \right.}} = qr\left\{ \left[ \sqrt {{W_{c,i}}} ({\tau _{2:2n + 1}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}})\sqrt {{{\hat Q}_{k - 1}}} \right]\right\} }\\ {S = cholupdate}\left\{ \left[{S_{k\left| {k - 1} \right.}},\right.\right. \\ \;\;\;\;\;\;\left.\left.\sqrt {{W_{c,i}}} ({\tau _{2:2n + 1}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}), \pm \right]\right\} \end{array}} \right.$$

(14)

b)测量更新

$${\tilde{\boldsymbol y}_k} = {\boldsymbol{y}_k} - H{\hat{\boldsymbol x}_{k\left| {k - 1} \right.}}$$

(15)

计算混合误差协方差矩阵平方根

$${P_{xy}} = S_{k\left| {k - 1} \right.}^{\rm{T}}{S_{k\left| {k - 1} \right.}}{H_k}$$

(16)

计算新息协方差矩阵

$${S_{y,k}} = qr{\bigg\{ \left[H{S_{k\left| {k - 1} \right.}}){{ }}\sqrt R \right]^{\rm{T}}}\bigg\} $$

(17)

计算滤波增益

$${K_k} = \frac{\frac{P_{xy,k}} {S_{y,k}^{\rm{T}}}}{S_{y,k}}$$

(18)

更新状态值

$${\hat{\boldsymbol x}_k} = {\hat{\boldsymbol x}_{k\left| {k - 1} \right.}} + {K_k}{\tilde{\boldsymbol y}_k}$$

(19)

更新状态误差协方差矩阵平方根

$$\left\{ {\begin{array}{*{20}{l}} {U = {K_k}S_{y,k}^{\rm{T}}} \\\\ {{S_k} = cholupdate\left\{ {S_{k\left| {k - 1} \right.}}{{ }}U{{ }} - 1\right\} } \end{array}} \right.$$

(20)

估计过程噪声均值

$${\hat{\boldsymbol q}_{k + 1}} = \dfrac{1}{k}\left[(k - 1){\hat{\boldsymbol q}_k} + {\hat{\boldsymbol x}_{k + 1}} - \sum\limits_{i = 1}^{2n + 1} {W_i^m} {f_o}({\boldsymbol{\xi }_{i,k}})\right]$$

(21)

估计过程噪声协方差矩阵平方根

$$\left\{ {\begin{array}{*{20}{l}} {\sqrt {Q1} = cholupdate\left\{ \sqrt {{{\hat Q}_{k - 1}}} ,\left| {{{\hat{\boldsymbol x}}_k} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}} \right|,\dfrac{1}{k}\right\} }\\ {\sqrt {Q2} = cholupdate\left\{ \sqrt {Q1} ,U, - \dfrac{1}{k}\right\} }\\ {\sqrt {{{\hat Q}_k}} = {\rm{diag}}\left\{ \sqrt {{\rm{diag}}\left\{\sqrt {Q2} {{\sqrt {Q2} }^{\rm{T}}}\right\}} \right\} }\\[-10pt] \end{array}} \right.$$

(22)

注1. Sage-Husa噪声估计器不能同时处理过程噪声和测量噪声都未知的情况, 否则会造成滤波发散, 故上述算法假设测量噪声的均值和方差均是已知的^[15].

注2. 式(14)的最后一个子式中, $ \pm $表示当${W_{c,1}}> 0$时, 取正; 反之, 取负, 以保证根号下不出现负数^[16].

2.4 自适应平方根无迹卡尔曼滤波稳定性

本小节给出自适应平方根无迹卡尔曼滤波算法的估计误差保持随机有界性的充分条件, 并给出证明.

对于非线性电网系统(6)

$$\begin{split} &{\boldsymbol{x}_k} = {f_o}({\boldsymbol{x}_{k - 1}}) + {B_o}\boldsymbol{u} + \boldsymbol{q} + {\boldsymbol{\mu }_{k - 1}}\\ &{\boldsymbol{y}_k} = H{\boldsymbol{x}_{k - 1}} + \boldsymbol{r} + {\boldsymbol{v}_k} \end{split}$$

(23)

定义

$${\tilde{\boldsymbol x}_k} = {\boldsymbol{x}_k} - {\hat{\boldsymbol x}_k}\;\;\;\;\;\;\;\;\;\;\;\;$$

(24)

$${\tilde{\boldsymbol x}_{k\left| {k - 1} \right.}} = {\boldsymbol{x}_k} - {\hat{\boldsymbol x}_{k\left| {k - 1} \right.}}$$

(25)

在${\hat{\boldsymbol x}_{k - 1}}$处对${{{\boldsymbol{x}}}_k}$泰勒级数展开

$$\begin{split}{\boldsymbol{x}_k} = & {\;f_o}({{\hat{\boldsymbol x}}_{k - 1}}) + {B_o}\boldsymbol{u} + \boldsymbol{q} + \nabla {f_o}({{\hat{\boldsymbol x}}_{k - 1}}){{\tilde{\boldsymbol x}}_{k - 1}}+\\ &\dfrac{1}{2}{\nabla ^2}{f_o}({{\hat{\boldsymbol x}}_{k - 1}})\tilde{\boldsymbol x}_{k - 1}^2 + \cdots + \omega \end{split}$$

(26)

其中

$${\nabla ^i}{f_{{o}}}(\hat{\boldsymbol x}){\tilde{\boldsymbol x}^i} = {\left(\sum\limits_{j = 1}^n {{{\tilde{\boldsymbol x}}_j}} \frac{\partial }{{\partial {\boldsymbol{x}_j}}}\right)^i}{\left. {{f_{{o}}}(\boldsymbol{x})} \right|_{\boldsymbol{x} = {{\hat{\boldsymbol x}}_{k - 1}}}}$$

在${\hat x_{k - 1}}$处对${\hat x_{k\left| {k - 1} \right.}}$泰勒级数展开

$$\begin{split} {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}} = &\left(1 - \frac{n}{{n + \lambda }}\right) \times (f({{\hat{\boldsymbol x}}_{k - 1}}) + B\boldsymbol{u}) + \boldsymbol{q}\; + \\ &\frac{1}{{2(n + \lambda )}} \times \sum\limits_{i = 2}^{n + 1} {(f\left[({{\hat{\boldsymbol x}}_{k - 1}}) + \gamma {{({S_{k - 1}})}_i}\right] + } \\ &B\boldsymbol{u} + \boldsymbol{q}) + \frac{1}{{2(n + \lambda )}} \times \sum\limits_{i = 2}^{n + 1} {\bigg(f\bigg[({{\hat{\boldsymbol x}}_{k - 1}})\;- } \\ &\gamma {({S_{k - 1}})_i}\bigg] + B\boldsymbol{u} + \boldsymbol{q}\bigg)=\\ \;& f({{\hat{\boldsymbol x}}_{k - 1}}) + B\boldsymbol{u} + \boldsymbol{q} + {\nabla ^2}f({{\hat{\boldsymbol x}}_{k - 1}}){P_{k - 1}}\\[-10pt] \end{split}$$

(27)

将式(26)和式(27)代入式(25), 得

$${\tilde{\boldsymbol x}_{k\left| {k - 1} \right.}} = {\beta _k}{F_k}{\tilde{\boldsymbol x}_{k - 1}} + {\boldsymbol{\mu }_k}$$

(28)

式中

$${F_k} = {\left. {\dfrac{{\partial {f_o}({\boldsymbol{x}_{k - 1}})}}{{\partial {\boldsymbol{x}_{k - 1}}}}} \right|_{\boldsymbol{x} = {{\hat{\boldsymbol x}}_{k - 1}}}}$$

(29)

其中, ${\beta _k} = {\rm{diag}}\left\{ {{\beta _{1k}}, \cdots ,{\beta _{nk}}} \right\}$为未知的矩阵, 其作用是弥补建模一阶线性化引起的模型误差.

对于式(24)中的输出方程, 定义

$${\tilde{\boldsymbol y}_{k\left| {k - 1} \right.}} = {\boldsymbol{y}_k} - {\hat{\boldsymbol y}_{k\left| {k - 1} \right.}}$$

(30)

式(30)可写成

$$\begin{split} {{\tilde{\boldsymbol y}}_{k\left| {k - 1} \right.}} =\;& H{\boldsymbol{x}_k} + \boldsymbol{r} + {\boldsymbol{v}_k} - (H{{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}} + \boldsymbol{r}) = \\ &H{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} + \boldsymbol{r} + {\boldsymbol{v}_k} \end{split}$$

(31)

由式(14)可知, ${S_{k\left| {k - 1} \right.}}$是${\hat P_{k\left| {k - 1} \right.}}$的平方根, 根据$qr( \cdot )$函数的模式$[Q,R] = qr(A)$, 此时

$$A = {\left[\sqrt {{W_{c,i}}} \left({\boldsymbol{\tau }_{2:2n + 1,k\left| {k - 1} \right.}} - {\hat{\boldsymbol x}_{k\left| {k - 1} \right.}}\right)\sqrt {{{\hat Q}_{k - 1}}} \right]^{\rm{T}}}$$

使用$qr( \cdot )$函数有$A = Q\times R$, 那么

$$\begin{split} {R^{\rm{T}}}R =\; & R{Q^{\rm{T}}}QR = {A^{\rm{T}}}A =\\ & \sum\limits_{i = 2}^{2n + 1} {\omega _i}{{\left({\boldsymbol{\tau }_{i,k\left| {k - 1} \right.}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}\right)}^{\rm{T}}}\; \times \\ &\left({\boldsymbol{\tau }_{i,k\left| {k - 1} \right.}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}\right) + {{\hat Q}_k} \end{split}$$

再根据$cholupdate$函数的模式, 此时

$$\begin{split} S_{k\left| {k - 1} \right.}^{\rm{T}}S_{k\left| {k - 1} \right.}^{} = \;\;& {R^{\rm{T}}}R + \sqrt {{W_{c,1}}} \left({\boldsymbol{\tau }_{1,k\left| {k - 1} \right.}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}\right) \times \\ &\sqrt {{W_{c,1}}} {\left({\boldsymbol{\tau }_{1,k\left| {k - 1} \right.}} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}\right)^{\rm{T}}}=\\ &P_{k\left| {k - 1} \right.}^{}\\[-10pt] \end{split}$$

(32)

与之对比, 先验状态误差协方差矩阵的实际值为

$$\begin{split} {P_{k\left| {k - 1} \right.}} = \;&{\rm{E}}\left({{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}}\tilde{\boldsymbol x}_{k\left| {k - 1} \right.}^{\rm{T}}\right)=\\ & {\rm{E}}\left(({\beta _k}{F_k}{{\hat{\boldsymbol x}}_{k - 1}} + {\boldsymbol{\mu }_k}){({\beta _k}{F_k}{{\hat{\boldsymbol x}}_{k - 1}} + {\boldsymbol{\mu }_k})^{\rm{T}}}\right)=\\ & {\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + \Delta {P_{k\left| {k - 1} \right.}} + {Q_{k - 1}} \end{split}$$

(33)

式中, $\Delta {P_{k\left| {k - 1} \right.}}$是${\rm{E}}\left(({\beta _k}{F_k}{\tilde{\boldsymbol x}_{k - 1}}){({\beta _k}{F_k}{\tilde{\boldsymbol x}_{k - 1}})^{\rm{T}}}\right)$与${\beta _k}\times {F_k} {\hat P_{k - 1}} F_k^{\rm{T}}{\beta _k}$的差值. 引入$\delta {P_{k\left| {k - 1} \right.}}$来表示误差协方差矩阵的真实值${P_{k\left| {k - 1} \right.}}$与${\hat P_{k\left| {k - 1} \right.}}$的差异, 式(32)可写为

$$\begin{split} {{\hat P}_{k\left| {k - 1} \right.}} =\;& {P_{k\left| {k - 1} \right.}} + \delta {P_{k\left| {k - 1} \right.}} = \\ &{\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + \Delta {P_{k\left| {k - 1} \right.}} + \\ &Q{}_k +\delta {P_{k\left| {k - 1} \right.}}=\\ &{\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {\Xi _{k - 1}} \\ \end{split} $$

(34)

此时, $\delta {P_{k\left| {k - 1} \right.}}$中包含了${Q_k}$与${\hat Q_k}$的差值, 因为在ASRUKF算法中噪声信息是未知的. 其中, ${\Xi _k} = \Delta {P_{k\left| {k - 1} \right.}} + {Q_k} + \delta {P_{k\left| {k - 1} \right.}}$.

预测测量协方差矩阵${\hat P_{yy}}$、互协方差矩阵${\hat P_{xy}}$以及后验状态误差协方差矩阵${\hat P_k}$可以写成

$${\hat P_{yy}} = H{\hat P_{k\left| {k - 1} \right.}}{H^{\rm{T}}} + {R_k}\;\;\;\;$$

(35)

$${\hat P_{xy}} = {\hat P_{k\left| {k - 1} \right.}}{H^{\rm{T}}}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;$$

(36)

$${\hat P_k} = {\hat P_{k\left| {k - 1} \right.}} - {\hat P_{xy}}\hat P_{yy}^{ - 1}\hat P_{xy}^{\rm{T}}$$

(37)

引理 1^[17]. 对于一个随机变量${\xi _k}$及其函数$V({\xi _k})$, 实数${\upsilon _{\min }},{\upsilon _{{\rm{\max}} }},\mu > 0,\;0 < \lambda \leq 1$对于任意的k, 若满足以下两个不等式性质:

$$\begin{split} &{\upsilon _{\min }}{\left\| {{\boldsymbol{\xi }_k}} \right\|^2} \leq V({\xi _k}) \leq {\upsilon _{\max }}{\left\| {{\boldsymbol{\xi }_k}} \right\|^2}\\ &{\rm{E}}(V({\boldsymbol{\xi }_k})\left| {{\boldsymbol{\xi }_{k - 1}}} \right.) - V\left({\boldsymbol{\xi }_{k - 1}}\right) \leq \mu - \lambda V\left({\boldsymbol{\xi }_{k - 1}}\right) \end{split}$$

(38)

则该随机变量均方有界, 即

$$\begin{split} {\rm{E}}\left({\left\| {{\xi _k}} \right\|^2}\right) \leq& \dfrac{{{\upsilon _{\max }}}}{{{\upsilon _{\min }}}}{\rm{E}}\left({\left\| {{\xi _0}} \right\|^2}\right){(1 - \lambda )^k} + \\ &\dfrac{\mu }{{{\upsilon _{\min }}}}\sum\limits_{i = 1}^{k - 1} {{{(1 - \lambda )}^i}} \\ \end{split} $$

(39)

引理 2^[17]. 假设矩阵$A \in {{\boldsymbol{{\rm{R}}}}^{n \times m}}$, $B \in {{\boldsymbol{{\rm{R}}}}^{m \times n}}$, $C \in {{\boldsymbol{{\rm{R}}}}^{n \times n}}$, 如果$A > 0,\;C > 0$, 则有

$${A^{ - 1}} > B{({B^{\rm{T}}}AB + C)^{ - 1}}{B^{\rm{T}}}$$

(40)

引理 3^[17]. 假设矩阵$A \in {{\boldsymbol{{\rm{R}}}}^{n \times n}}$, $C \in {{\boldsymbol{{\rm{R}}}}^{n \times n}}$, 如果$A > 0$, $C > 0$, 则有

$${A^{ - 1}} > {(A + C)^{ - 1}}$$

(41)

定理 1. 对于非线性系统(23), 如果有如下假设条件成立

1)对于任意的k, 存在实数${f_{\min }},{h_{\min }},{\beta _{\min }},$ ${f_{\max }}, {h_{\max }},{\beta _{\max }} \ne {{0}}$, 使得下式成立:

$$\begin{split} f_{\min }^2I \leq\; &{F_k}F_k^{\rm{T}} \leq f_{\max }^2I \\ h_{\min }^2I \leq\;& {H_k}H_k^{\rm{T}} \leq h_{\max }^2I \\ \beta _{\min }^2I \leq \;&{\beta _k}\beta _k^{\rm{T}} \leq \beta _{\max }^2I \\ \end{split} $$

(42)

2)存在实数${\Xi _{\max }},{\Xi _{\min }},{r_{\min }},{p_{\max }},{p_{\min }}$, 使得下式成立:

$$\begin{split} & {p_{\min }}I \leq {P_k} \leq {p_{{\rm{max}}}}I \\ & {\Xi _{\min }}I < {\Xi _k} \leq {\Xi _{\max }}I \\ &{r_{\min }}I \leq {R_k} \\ \end{split} $$

(43)

那么自适应平方根无迹卡尔曼滤波算法的状态估计误差将是均方有界的, 即自适应平方根无迹卡尔曼滤波器稳定.

证明. 选择函数

$${V_k}({\tilde x_k}) = \tilde x_k^{\rm{T}}\hat P_k^{ - 1}{\tilde x_k}$$

(44)

由式(43)可知

$$\dfrac{1}{{{p_{\max }}}}{\left\| {{{\tilde x}_k}} \right\|^2} \leq V({\tilde x_k}) \leq \dfrac{1}{{{p_{\min }}}}{\left\| {{{\tilde x}_k}} \right\|^2}$$

满足了引理1中第1个性质, 为满足引理1, 还需要得到式(38)中第2个性质${\rm{E}}[{V_k}({x_k})] - {V_{k - 1}}({\tilde x_{k - 1}})$的上界.

后验状态误差协方差矩阵${\hat P_k}$的另一种表达式为

$${\hat P_k} = {\hat P_{k\left| {k - 1} \right.}} - {\hat P_{xy}}\hat P_{yy}^{ - 1}\hat P_{xy}^{\rm{T}} = (I - {K_k}{H_k}){\hat P_{k\left| {k - 1} \right.}}$$

(45)

其中

$${K_k} = {\hat P_{k\left| {k - 1} \right.}}H_k^{\rm{T}}{({H_k}{\hat P_{k\left| {k - 1} \right.}}H_k^{\rm{T}} + {R_k})^{ - 1}}\;\;\;\;\;\;\;\;\;\;\;\;$$

(46)

$${\tilde x_k} = {x_k} - ({\hat x_{k\left| {k - 1} \right.}} + {\hat P_{xy}}\hat P_{yy}^{ - 1}{\tilde y_k}) = {\tilde x_{k\left| {k - 1} \right.}} - {K_k}{\tilde y_k}\;\;$$

(47)

将式(47)代入式(44), 得

$$\begin{split} {V_k}({{\tilde{\boldsymbol x}}_k}) =\;& {\left({{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} - {K_k}{{\tilde{\boldsymbol y}}_k}\right)^{\rm{T}}}\hat P_k^{ - 1}\left({{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} - {K_k}{{\tilde{\boldsymbol y}}_k}\right) =\\ & {{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}}\hat P_k^{ - 1}{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} - {\left[{H_k}{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} + {v_k}\right]^{\rm{T}}}\; \times \\ &K_k^{\rm{T}}\hat P_k^{ - 1}{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} - \tilde x_{k\left| {k - 1} \right.}^{\rm{T}}\hat P_k^{ - 1}{K_k}\big[{H_k}{{\tilde x}_{k\left| {k - 1} \right.}} + \\ &{v_k}\big] + {\left[{H_k}{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} + {v_k}\right]^{\rm{T}}}K_k^{\rm{T}} \times \\ &\hat P_k^{ - 1}{K_k}\left[{H_k}{{\tilde{\boldsymbol x}}_{k\left| {k - 1} \right.}} + {v_k}\right]\\[-10pt] \end{split}$$

(48)

式(46)的另一种等价形式为

$${K_k} = (I - {K_k}{H_k}){\hat P_{k\left| {k - 1} \right.}}H_k^{\rm{T}}R_k^{ - 1} = {\hat P_k}H_k^{\rm{T}}R_k^{ - 1}$$

(49)

式(45)用矩阵求逆公式得

$$\hat P_k^{ - 1} = \hat P_{k\left| {k - 1} \right.}^{ - 1} + H_k^{\rm{T}}R_k^{ - 1}{H_k}$$

(50)

将式(49), 式(50)和式(28)代入式(48)中, 式(48)的条件期望为

$$\begin{split} {\rm{E}}({V_k}({{\tilde x}_k})\left| {{{\tilde x}_{k - 1}}} \right.) =\;& {\rm{E}}{\bigg(} {({A_{k - 1}}{{\tilde x}_{k - 1}} + {\mu _k})^{\rm{T}}} \times\\ &\hat P_{k\left| {k - 1} \right.}^{ - 1}({A_{k - 1}}{{\tilde x}_{k - 1}} + {\mu _k})-\;\;\\ & {\left[{B_k}\left({A_{k - 1}}{{\tilde x}_{k - 1}} + {\mu _k}\right)\right]^{\rm{T}}}\\ &\left(\Sigma _k^{ - 1} - \Sigma _k^{ - 1}{B_k}{{\hat P}_k}B_k^{\rm{T}}\Sigma _k^{ - 1}\right) \times \\ &\left[{B_k}({A_{k - 1}}{{\tilde x}_{k - 1}} + {\mu _k})\right] + \\ &v_k^{\rm{T}}\Sigma _k^{ - 1}{B_k}{{\hat P}_k}B_k^{\rm{T}}\Sigma _k^{ - 1}{v_k}\left| {{{\tilde x}_{k - 1}}} \right.{\bigg)} \end{split}$$

(51)

将式(34)代入式(51)中, 并由引理3得

$$\begin{split} &{{\rm{E}}{\bigg(}{V_k}({{\tilde{\boldsymbol x}}_k})\left| {{{\tilde{\boldsymbol x}}_{k - 1}}} \right. {\bigg)} } \leq\\ &\;\;\;{ {\rm{E}}{\bigg(} {{{({\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})}^{\rm{T}}}} {{\left( {{\beta _k}{F_k}\hat P_{k{\rm{ - }}1}^{}F_k^{\rm{T}}{\beta _k}} \right)}^{ - 1}}}\; \times \\ &\;\;\;{({\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}}) + \omega _k^{\rm{T}}{{({\beta _k}{F_k}\hat P_{k{\rm{ - }}1}^{}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})}^{ - 1}} \times }\\ &\;\;\;{{\omega _k} - {{({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})}^{\rm{T}}} \times \left[ {{H_k}} \right.\left( {{\beta _k}} \right.{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}} \times }\\ &\;\;\;{{\beta _k} + \left. {{{\hat Q}_k}} \right)H_k^{\rm{T}} + {{\left. {{R_k}} \right]}^{ - 1}}({H_k}{\beta _k}{F_k} \times {{\tilde{\boldsymbol x}}_{k - 1}})\;- }\\ &\;\;\;{{{({H_k}{\omega _k})}^{\rm{T}}}\left[ {{H_k}} \right.\left({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k}\right) \times H_k^{\rm{T}} + }\\ &\;\;\;{{{\left. {{R_k}} \right]}^{ - 1}}({H_k}{\omega _k}) + v_k^{\rm{T}}R_k^{ - 1}{H_k}{{\hat P}_k}H_k^{\rm{T}}R_k^{ - 1}{v_k} {\left| {{{\tilde{\boldsymbol x}}_{k - 1}}} \right.{\bigg)} }} \end{split}$$

(52)

观察式(52)的第1项

$$\begin{split} &{\rm{E}}{\bigg(} {({\beta _k}{F_k}{{\tilde x}_{k - 1}})^{\rm{T}}}{\left( {{\beta _k}{F_k}\hat P_{k{\rm{ - }}1}^{}F_k^{\rm{T}}{\beta _k}} \right)^{ - 1}} \times \\ &\qquad({\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})\left| {{{\tilde{\boldsymbol x}}_{k - 1}}} \right.{\bigg)} = \\ & \qquad\tilde{\boldsymbol x}_{k - 1}^{\rm{T}}{{\hat P}_{k - 1}}{{\tilde{\boldsymbol x}}_{k - 1}} = {V_{k - 1}}({{\tilde{\boldsymbol x}}_{k - 1}}) \end{split}$$

(53)

两边同时减去式(53), 得

$$\begin{split} &{{\rm{E}}\bigg({V_k}({{\tilde{\boldsymbol x}}_k})\left| {{{\tilde{\boldsymbol x}}_{k - 1}}} \right.\bigg)- {V_{k - 1}}({{\tilde {\boldsymbol{x}}}_{k - 1}})}\leq\\ &\;\;\;{ {\rm{E}}\bigg( {\omega _k^{\rm{T}}} {{({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})}^{ - 1}} - }\\ &\;\;\;{{{({H_k}{\omega _k})}^{\rm{T}}} \times \left[ {{H_k}} \right.\big( {{\beta _k}} {F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} \;+ }\\ &\;\;\;{ {{{\hat Q}_k}} \big) H_k^{\rm{T}} + {{\left. {{R_k}} \right]}^{ - 1}}({H_k}{\omega _k}) + v_k^{\rm{T}}R_k^{ - 1}{H_k} \;\times }\\ &\;\;\;{{{\hat P}_k}H_k^{\rm{T}}R_k^{ - 1}{v_k} {\left| {{{\tilde{\boldsymbol{ x}}}_{k - 1}}} \right.} \bigg) - {{({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})}^{\rm{T}}} \times }\\ &\;\;\;{\left[ {{H_k}} \right.({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})H_k^{\rm{T}}}+\\ &\;\;\;{ {{\left. {{R_k}} \right]}^{ - 1}} \times \left({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}}\right)} \end{split}$$

(54)

观察式(54)的最后1项, 根据引理3, 得

$$\begin{split} &{\hat P_{k - 1}^{ - 1} > {{({H_k}{\beta _k}{F_k})}^{\rm{T}}}\left[ {({H_k}{\beta _k}{F_k}){{\hat P}_{k - 1}}{{({H_k}{\beta _k}{F_k})}^{\rm{T}}}} \right.}+\\ &\;\;{{{\left. { {H_k}{{\hat Q}_k}H_k^{\rm{T}} + {R_k}} \right]}^{ - 1}}({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})}\\[-10pt] \end{split}$$

(55)

式(55)分别左右乘${\tilde {\boldsymbol{x}}_{k - 1}}$, 有

$$\begin{split} &\tilde{\boldsymbol x}_{k - 1}^{\rm{T}}\hat P_{k - 1}^{ - 1}{{\tilde {\boldsymbol{x}}}_{k - 1}} > {({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})^{\rm{T}}}[({H_k}{\beta _k}{F_k}){{\hat P}_{k - 1}} \times\\ &\quad{({H_k}{\beta _k}{F_k})^{\rm{T}}} + {H_k}{\Xi _k}H_k^{\rm{T}} + {R_k}{]^{ - 1}}({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}}) \end{split}$$

(56)

选择

$$\begin{split} {\lambda _k} =& \; \Bigl\{ {({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})^{\rm T}}\Bigl[({H_k}{\beta _k}{F_k}){{\hat P}_{k - 1}}\; \times \\ &{({H_k}{\beta _k}{F_k})^{ {\rm{T}}}} + {H_k}{Q _k}H_k^{\rm T} + {R_k}{\Bigl]^{ - 1}}\;\times\\ & ({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})\Bigl\}/(\tilde{\boldsymbol x}_{k - 1}^{\rm{T}}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}}) \end{split}$$

(57)

由式(56)和式(57)可知, ${\lambda _k} < 1$, 得

$$\begin{split} {\lambda _k} \geq \;&{p_{\min }}{({h_{\min }}{\beta _{\min }}{f_{\min }})^2} \; \times \\ &\big[{p_{\max }}{({h_{\max }}{\beta _{\max }}{f_{\max }})^2} + \\ &{{\hat q}_{\max }}h_{\max }^2 + {r_{\max }}{\big]^{ - 1}} = {\lambda _{\min }} > 0 \end{split} $$

(58)

由式(57), 得

$$\begin{split} &- {({H_k}{\beta _k}{F_k})^{\rm{T}}}[({H_k}{\beta _k}{F_k}){{\hat P}_{k - 1}}{({H_k}{\beta _k}{F_k})^{\rm{T}}}+\\ &\quad\quad\quad{H_k}{{\hat Q}_k}H_k^{\rm{T}} + {R_k}{]^{ - 1}}({H_k}{\beta _k}{F_k}{{\tilde{\boldsymbol x}}_{k - 1}})\leq \\ &\quad\quad\quad-{\lambda _{\min }}{V_{k - 1}}({{\tilde{\boldsymbol x}}_{k - 1}}) \end{split}$$

(59)

考虑式(54)中的其他项, 得

$$\begin{split} &{{\mu _k} = {\rm{E}}\bigg({\omega _k^{\rm{T}}}\Big[{{{({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})}^{ - 1}}} }-\\ & \quad \; \quad{ {H_k}^{\rm{T}}{\left( {{H_k}} \right.}({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})H_k^{\rm{T}} + }\\ &\quad \; \quad{{\left. {{R_k}} \right)^{ - 1}} \times H_k\Big]}{\omega _k} + \\ &\quad \; \quad{\upsilon _k^{\rm{T}}R_k^{ - 1}{H_k}{{\hat P}_k}H_k^{\rm{T}}R_k^{ - 1}{v_k} {\left| {{{\tilde x}_{k - 1}}} \right.} \bigg)} \end{split}$$

(60)

考虑式(60)两边都是标量, 取迹不会改变等式, 所以根据

$${\rm{tr}}(AB) = {\rm{tr}}(BA)$$

(61)

得到

$$\begin{split} &{{\mu _k} = {\rm{tr}}\bigg\{ {\Big[ {{{({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})}^{ - 1}}} } - {H^{\rm{T}}_k} \times }\\ &\quad\quad\;{\left( {{H_k}} \right.({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})H_k^{\rm{T}} + }\\ &\quad\quad\;{{{ {{R_k}} )}^{ - 1}}{ {{H_k}} \Big]} {{Q_k}} \bigg\} + {\rm{tr}}\{[R_k^{ - 1}{H_k}{{\hat P}_k}H_k^{\rm{T}}]\left| {{{\tilde {\boldsymbol{x}}}_{k - 1}}} \right.\} } \end{split}$$

(62)

根据引理2, 得

$$\begin{split} &{({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})^{ - 1}} - \\ &\;\;\;\quad {H^{\rm{T}}_k}[{H_k}({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + \\ &\;\;\;\quad {{\hat Q}_k})H_k^{\rm{T}} + {R_k}{)^{ - 1}}{H_k}] > 0 \end{split}$$

(63)

由式(63)可知${\mu _k} > 0$, 由式(42)和式(43), 得

$$\begin{split} {\mu _k} \leq\; & {\rm{tr}}[{({\beta _k}{F_k}{{\hat P}_{k - 1}}F_k^{\rm{T}}{\beta _k} + {{\hat Q}_k})^{ - 1}}{{\hat Q}_k}]+\\ & {\rm{tr}}(R_k^{ - 1}{H_k}{{\hat P}_k}H_k^{\rm{T}})\leq\\ & \dfrac{{{q_{\max }}}}{{{{\hat q}_{\min }}}} \times L + \dfrac{{{h_{{{\max }^2}}}{p_{\max }}}}{{{r_{\min }}}} \times M = {\mu _{\max }} \end{split}$$

(64)

因此根据引理1, 得

$${\rm{E}}({V_k}({x_k})) - {V_{k - 1}}({\tilde x_{k - 1}}) \leq {\mu _{\max }} - {\lambda _{\min }}{V_{k - 1}}({\tilde x_{k - 1}})$$

(65)

保证是${\tilde x_k}$有界的, 即

$$\begin{split} {\rm{E}}({\left\| {{{\tilde{\boldsymbol x}}_k}} \right\|^2}) \leq\; & \dfrac{{{p_{\max }}}}{{{p_{\min }}}}{\rm{E}}({\left\| {{{\tilde{\boldsymbol x}}_0}} \right\|^2}){(1 - {\lambda _{\min }})^k}\;+\\ & \dfrac{{{\mu _{\max }}}}{{{p_{\min }}}}\sum\limits_{i = 1}^{k - 1} {{{(1 - {\lambda _{\min }})}^i}} \end{split}$$

(66)

因此, 自适应平方根无迹卡尔曼滤波算法状态估计误差将是均方有界的, 即自适应平方根无迹卡尔曼滤波器是稳定的. □

注3. 因为误差序列$\left\{ {{{\hat{\boldsymbol x}}_k},k \geq 1} \right\}$是随机过程, ${\rm{E}}({\left\| {{{\tilde{\boldsymbol x}}_k}} \right\|^2})$是该随机过程的均值函数, 随机有界指该随机过程有界. 即对于任意时刻$k$, 对于服从高斯分布的噪声$\omega $, 定理1确保了误差${\tilde{\boldsymbol x}_k}$的2范数平方的期望值是有界的, 即${\rm{E}}({\left\| {{{\tilde{\boldsymbol x}}_k}} \right\|^2}) \leq \alpha $上限$\alpha $的值由定理1确定.

3. 检测方法

本节分析了文献[6], [10-11]中所使用的检测方法的不足之处. 提出一种新的检测方法, 利用中心极限定理设计新的检测算法, 依据随机变量的统计特性提出双边假设检验法检测攻击, 从误检率方面对本文算法给出分析, 并与巴氏系数检测算法^[10]、欧几里德检测算法^[11]进行对比与分析.

文献[6]针对虚假数据注入攻击, 构造了CUSUM算法来检测攻击. 该算法步骤如下:

$${\chi _t} = ({y_t} - {H_k}{\hat x_{t\left| {t - 1} \right.}})\Sigma _t^{ - 1}({y_t} - {H_k}{\hat x_{t\left| {t - 1} \right.}})\;$$

$${p_t} = 1 - F({\chi _t})\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;$$

$${s_t} = \log \left(\dfrac{\alpha }{{{p_t}}}\right)\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;$$

$$\begin{split} & \Gamma = \inf \left\{ {t:{g_t} \geq h} \right\} \\ &{g_t} = \max \left\{ {0,{g_{t - 1}} + {s_t}} \right\} \\ \end{split}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; $$

文献[6]中假设噪声服从高斯分布, 可知${\chi _t}$服从${\chi ^2}$分布. $F( \cdot )$表示${\chi ^2}$分布的概率分布函数, ${p_t}$表示${\chi ^2}$分布的上(右)尾概率, $\alpha $表示系统给定的上(右)尾概率用于检验${y_t}$是否是异常值. 根据${\chi ^2}$分布的上$\alpha $分位点性质, 当${y_t}$是异常值时, ${s_t} > 0$; 当${y_t}$是正常值时, ${s_t} \leq 0$. 当${y_t}$是异常值时, ${g_t}$的值会随时间累积直到超过检测阈值$h$($h$的确定参考文献[6]), 此时发出警报. 该算法通过引入$h$, 避免了欧氏检测算法那种仅考虑误检率的阈值设计方法, 在误检率和误检周期之间做到了平衡, 既降低了噪声的影响, 又可以尽快检测出攻击.

但是, 该算法实质上是依赖与噪声的高斯分布特性, 即, 实质上还是${\chi ^2}$检测, 所以在遇到本文所研究的攻击时, 攻击前后的残差基本保持不变, 因此该算法是不可用的.

文献[10]中针对具有隐蔽性的复数域攻击, 提出一种检测思想, 滤波器在不受攻击下工作和滤波器在错误数据注入攻击下工作会产生两个服从高斯分布的残差序列, 通过计算上述两个服从高斯分布的残差序列的巴氏相似性系数$DB$来判断攻击是否发生.

当物理系统的基本动力学用线性系统和高斯统计特性的噪声进行建模时, 在正常操作期间, 残差序列${v^o}$遵循零均值高斯分布${v^o} = {\rm{N}}({r^o},{S^o})$, 攻击发生后残差序列$v_k^a = {\rm{N}}(r_k^a,S_k^a)$.

$$\begin{split} DB({v^o},v_k^a) =\;& \dfrac{1}{4}[R\{ {({r^o},r_k^a)^H}{\Lambda _k}({r^o},r_k^a)\}\;+ \\ & R\{ {({r^o},r_k^a)^H}{{\tilde \Lambda }_k}({r^o},r_k^a)\} ]\;+\\ & \dfrac{1}{2}\ln \frac{{\left| {{\Gamma _k}} \right|}}{{\sqrt {\left| {{S^o}} \right|\left| {S_k^a} \right|} }} \\ \end{split}$$

$${\Gamma _k} = \frac{1}{2}({S^o} + S_k^a) = \left[ {\begin{array}{*{20}{c}} {{\Gamma _k}}&{{{\tilde \Gamma }_k}} \\ {\tilde \Gamma _k^*}&{\Gamma _k^*} \end{array}} \right]\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;$$

该检测方法的优点同时也是缺点, 即巴氏相似性系数$DB$能检测复数域攻击, 但也仅能检测复数域攻击. 因为系统中有无复数域攻击, 差异主要体现在$DB({v^o},v_k^a)$的后两项上. 而本文的虚假数据注入攻击属于实数域范畴, 在正常操作期间, 残差序列${v^o}$遵循零均值高斯分布${v^o} = {\rm{N}}({m_1},{P_1})$, 攻击发生后残差序列$v_k^a = {\rm{N}}({m_2},{P_2})$, 实数域的巴氏距离计算式为

$$\begin{split} DB(v_k^o,v_k^a) = \;&\dfrac{1}{8}[{({m_1} - {m_2})^{\rm{T}}}{P^{ - 1}}({m_1} - {m_2}) ]\;+\\ & \dfrac{1}{2}\ln \dfrac{{\det P}}{{\sqrt {\det {P_1}\times \det {P_2}} }} \\ \end{split} $$

$$P = \frac{{{P_1} + {P_1}}}{2}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;$$

Bhattacharyya距离无法辨别虚假数据注入攻击前后两个残差序列的相似性, 因为虚假数据注入攻击前后残差保持不变. 除此之外, 该检测法的前提是线性系统, 因为线性运算虽然会改变系统中随机变量的均值和方差但并不会改变随机变量服从高斯分布的特性, 对于本文所研究的非线性系统, 这种方法是不可用的.

目前存在的检测方法大多是在加性噪声的基础上提出的, 在文献[11]中, 假设过程方程和量测方程相对于噪声是线性的, 即

$$\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\boldsymbol{x}_{k+1}} = f({\boldsymbol{x}_k},{u_k},{t_k}) + {\omega _k}$$

$$\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\boldsymbol{y}_{k+1}} = h({\boldsymbol{x}_k},{t_k}) + {v_k}$$

文献[11]中, 欧氏检测法的阈值和误检率正是基于噪声项的附加性得到的, 因为噪声并不参与到非线性变换中, 所以噪声的特性就是假设的正态分布的特性, 使用正态分布的$3\sigma $原则设计阈值. 只有在这种前提下, 欧氏距离检测算法才是正确的. 实际上, 过程和量测方程相对于噪声也可能是非线性的, 即

$$\qquad\qquad{\boldsymbol{x}_k} = f({\boldsymbol{x}_k},{u_k},{t_k},{\omega _k})$$

$$\qquad\qquad y_{k}=h(x_{k},t_{k},v_{k})$$

在这种情况下, 欧氏距离检测法的阈值设计思想是行不通的, 因为噪声不再作为附加项, 而是直接参与到非线性变换中. 此时, 噪声不再服从正态分布, 统计特性未知, 因此阈值无法设计.

综上, 本文提出一种检测思想: 在攻击未发生时, 利用本文提出的非线性滤波方法, 根据总线$i$上传感器在时间范围$[{t_0},{t_k}]$收集到的测量数据获得系统状态估计值, 进而获得后验状态误差值${\tilde{\boldsymbol x}_k}$, 此时并不知道状态误差${\tilde{\boldsymbol x}_k}$这个随机变量服从什么分布, 而且也不用对该随机变量的统计特性作出假设. 本文将时间区间$[{t_0},{t_k}]$分成 $\;l$个小区间, 对每一个小区间包含的${\tilde x_k}$使用中心极限定理^[18]构造随机变量${\boldsymbol{S}_i}\tilde x = \frac{1}{T}\sum\nolimits_{k = 1}^T {{{\tilde{\boldsymbol x}}_k}}$, 得到 $\;l$个${\boldsymbol{S}_i}\tilde x,\;i \in l$, 这 $l$个随机变量是服从标准正态分布的, 当攻击发生时, ${S_i}\tilde x$会偏离标准整体分布, 从而达到检测攻击的目的.

本文所提出的检测方法从系统内部状态误差入手, 不受限于传感器参数的概率分布未知所带来的影响, 不受限于加性噪声和非加性噪声的情况, 可以成功地完成攻击检测, 并从数学角度给出了误检率与阈值之间的关系.

3.1 攻击检测

式(24)定义了随机变量${\tilde{\boldsymbol x}_k}$, 本节定义当系统未遭受攻击时, 变量为${\tilde{\boldsymbol x}_k}$; 当系统遭受攻击时, 变量为$\tilde{\boldsymbol x}_k^a$. 下面分别计算这两个随机变量的两个统计特性, 均值和方差.

1)当系统正常运行时

$$\begin{split} {\rm{E}}({{\tilde{\boldsymbol x}}_k}) =\;& {\rm{E}}({\boldsymbol{x}_k} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}} - {K_k}({\boldsymbol{y}_k} - {H_k}{{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}))=\\ &{\rm{E}}({\boldsymbol{x}_k} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}} - {K_k}({H_k}{{\hat{\boldsymbol x}}_k} + {v_k} - {H_k}{{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}})) =\\ &{\rm{E}}(({I_k} - {K_k}{H_k})({\boldsymbol{x}_k} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}) - {K_k}{\rm{E}}({v_k}))=\\ &({I_k} - {K_k}{H_k}){\rm{E}}({\boldsymbol{x}_k} - {{\hat{\boldsymbol x}}_{k\left| {k - 1} \right.}}) - {K_k}{\rm{E}}({v_k}) \end{split}$$

(67)

当初值设置为${\hat{\boldsymbol x}_0} = {\rm{E}}({\boldsymbol{x}_0})$时, 式(67)变为

$${\rm{E}}({\tilde{\boldsymbol x}_k}) = 0$$

(68)

方差为

$$\begin{aligned} {\hat P}_k|k=\;&{\rm{E}}\left({\tilde{\boldsymbol x}}_{k}\tilde{\boldsymbol x}^{\rm{T}}_k\right)=\\ &{\rm{E}}\left[(I_{k}-K_{k}H_{k})({\boldsymbol{x}}_{k}-\hat{\boldsymbol x}_{k|k-1})-K_{k}{\rm{E}}(v_{k}\right]\times\\ &[(I_{k}-K_{k}H_{k})({\boldsymbol{x}}_{k}-{\hat{\boldsymbol x}}_{k|k-1})-K_{k}{\rm{E}}(v_{k})]^{\rm{T}}=\\ &(I_{k}-K_{k}H_{k})\hat{P}_{k|k-1}(I_{k}-K_{k}H_{k})^{\rm{T}}+K_{k}R_{k}K_{k}^{\rm{T}} \end{aligned}$$

(69)

2)当系统遭受攻击时

$$\begin{split}{\rm{E}}(\tilde{\boldsymbol x}_k^a) = \;&{\rm{E}}(\boldsymbol{x}_k^a - \hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a - {K_k}(\boldsymbol{y}_k^a + {a_k} - {H_k}\hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a))=\\ &{\rm{E}}(\boldsymbol{x}_k^a - \hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a - \\ &{K_k}({H_k}\boldsymbol{x}_k^a + {v_k} + {a_k} - {H_k}\hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a))=\\ &{\rm{ E}}(({I_k} - {K_k}{H_k})(\boldsymbol{x}_k^a - \hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a) - {K_k}{\rm{E}}({v_k})) \end{split}$$

(70)

方差为

$$\begin{split} \hat P_{k\left| k \right.}^a =\;& {\rm{E}}[(\tilde{\boldsymbol x}_k^a){(\tilde{\boldsymbol x}_k^a)^{\rm{T}}}]=\\ & {\rm{E}}\{ [({I_k} - {K_k}{H_k})(\boldsymbol{x}_k^a - \hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a) - {K_k}{\rm{E}}({v_k} + {a_k})] \;\times \\ &{[({I_k} - {K_k}{H_k})(\boldsymbol{x}_k^a - \hat{\boldsymbol x}_{k\left| {k - 1} \right.}^a) - {K_k}{\rm{E}}({v_k} + {a_k})]^{\rm{T}}} =\\ & ({I_k} - {K_k}{H_k})\hat P_{k\left| {k - 1} \right.}^a{({I_k} - {K_k}{H_k})^{\rm{T}}} + {K_k}{R_k}K_k^{\rm{T}} \end{split}$$

(71)

根据第3.1节的描述, 虽然得到了${\tilde x_k}$和$\tilde x_k^a$ 的均值和方差, 但是无法确定随机变量${\tilde x_k}$和$\tilde x_k^a$服从何种分布. 此时, 对于利用${\tilde x_k}$应用中心极限定理的性质, n个独立随机变量之和的分布函数趋于正态分布函数.

现构造随机变量$S\tilde x = \frac{1}{T}\sum\nolimits_{k = 1}^T {{{\tilde{\boldsymbol x}}_k}}$, 根据中心极限定理, 可知

$$S\tilde x \sim {\rm{N}}\left(\frac{1}{T}\sum\limits_{k = 1}^T {{\rm{E}}(\tilde{\boldsymbol x})} ,\frac{1}{T}\sum\limits_{k = 1}^T {{{\hat P}_{k\left| k \right.}}} \right)$$

同时按照上述形式构造$S{\tilde x^a} = \frac{1}{T}\sum\nolimits_{k = 1}^T {\tilde{\boldsymbol x}_k^a}$, 注意此时$S{\tilde x^a}$并不服从正态分布, 构造$S{\tilde x^a}$是为了保持变量形式上的一致, 记$S{\tilde x^a}$服从的分布为$A$, 得

$$S{\tilde x^a} \sim A\left(\frac{1}{T}\sum\limits_{k = 1}^T {{\rm{E}}(\tilde{\boldsymbol x}_k^a)} ,\frac{1}{T}\sum\limits_{k = 1}^T {\hat P_{k\left| k \right.}^a} \right)$$

为书写简便, 将$S\tilde x$和$S{\tilde x^a}$分别记为

$$S\tilde x \sim {\rm{N}}(m,\upsilon ),\;\;S{\tilde x^a} \sim A(m + r\upsilon ,q\upsilon )$$

其中, $q$和$r$分别表示系统遭受攻击后均值和方差的变化. 本文使用$m$和$\upsilon $作为系统正常行为的度量. 接下来, 建立了一对阈值${T_l}$和${T_h}$, 用于检测系统是否遭受攻击, ${T_l} = m - k\upsilon $, ${T_h} = m + k\upsilon $, 表示允许偏离原数据的最大范围. 当随机变量$S \notin [{T_l},{T_h}]$时, 则确定系统遭受攻击^[19]. 基于此, 提出一个双边假设检验

$$\left\{ {\begin{aligned} {{H_0},\;S \in [{T_l},{T_h}]} \\ {{H_1},\;S \notin [{T_l},{T_h}]} \end{aligned}} \right.\;\;$$

(72)

式中, ${H_0}$代表系统没有遭受攻击, ${H_1}$代表系统遭受攻击.

误检率${P_F}$的定义是攻击没有发生却宣布攻击发生的概率, 即$S\tilde x \sim {\rm{N}}(m,\upsilon )$, 同时, $\left| {S\tilde x - m} \right| > k\upsilon $.

$$\begin{split} {P_F}(k) &= 1 - \int\nolimits_{m - kv}^{m + kv} {\frac{{{\exp (\frac{ - {{(u - m)}^2}}{(2{v^2})})}}}{{v\sqrt {2\pi } }}} {\rm{d}}u =\\ &\;\;\;\;\;1 - \Phi (k) + \Phi ( - k) \\ \end{split} $$

4. 仿真

本节验证本文设计的自适应平方根无迹卡尔曼滤波算法估计状态及所提出检测算法检测攻击有效性, 并与欧几里得检测方法进行对比. 仿真参数和初始值为 $m = [8.9,{{ }}8.8,{{ }}8.5]$, $d = [3.1,{{ }}3.4,{{ }}3.7], $ $u = [6.3, {{ }}1.6,{{ }}8.5]^{\rm{T}}$, $B = \left[{0_{3 \times 3}};{\rm{diag}}\{1/m\}\right]$.

为了最小化由噪声引起的误报率, 欧几里得检测法的阈值设计采用$3\sigma $准则($\sigma $是噪声信号的标准差), 那么得到先验阈值$f = 3\sigma = 0.3$.

4.1 ASRUKF滤波算法的有效性

对于模型(6), 未遭受攻击时的状态估计, 从图3可以看出, ASRUKF滤波器的滤波性能很好, 滤波精度很高.

图 3 ASRUKF下的状态估计

Fig. 3 State estimation in ASRUKF

下载: 全尺寸图片幻灯片

4.2 虚假数据注入攻击下的三种检测方法

1) 欧几里得检测法^[11]. 由图 4可知, 在迭代步数$k = 30$时刻对系统注入攻击, 在$k = 80$时刻成功检测到攻击, 欧几里得检测法可以成功检测出虚假数据注入攻击.

图 4 两种检测方法针对隐蔽假数据攻击

Fig. 4 Two detection methods for covert false data attack

下载: 全尺寸图片幻灯片

2) 巴氏相似性系数检测法^[10]. 巴氏系数范围为0到1之间. 越接近0, 证明两序列越相似. 由图5可知, 巴氏相似性系数无法辨别虚假数据注入攻击前后两个残差序列的相似性, 不可以检测出虚假数据注入攻击.

图 5 巴氏相似性系数

Fig. 5 Bhattacharyya coefficient

下载: 全尺寸图片幻灯片

3) 本文所提的检测方法. 首先, 使用分位数–分位数(Quantile-quantile, Q-Q)图来确定前面根据中心极限定理构造的$S\tilde x$分布可以近似为高斯分布. 因为根据中心极限定理, 随机变量的数量越多, 它们的和也越近似服从高斯分布. 然而, 考虑实际情况, 取无限多的独立随机变量是不符合实际的, 所以在近似高斯分布的同时, 要尽量选取较少的独立随机变量. 要利用Q-Q图鉴别样本数据是否近似于高斯分布, 只需观察Q-Q图上的点是否近似地在一条直线附近.

本文的检测算法选取了$T = 30$个独立随机变量参与运算, 由图6可知, 数据的分布非常接近直线, 这证明$S\tilde x$的分布可以视为高斯分布.

图 6 $S\tilde x$的Q-Q图

Fig. 6 Quantile-quantile plot of $S\tilde x$

下载: 全尺寸图片幻灯片

由图7可以看出, 本文的检测算法对于隐蔽假数据攻击的有效性, 在迭代步数$k = 30$处注入攻击, 选取30个迭代步骤下的随机变量构造出$S\tilde x$, 所以在$k = 60$处超过了阈值, 实现了攻击检测.

图 7 本文提出的攻击检测方法

Fig. 7 Attack detection proposed in this paper

下载: 全尺寸图片幻灯片

图8给出了误检率随$k$变化的曲线.从图8可知, 一个算法的误检率过高, 是因为阈值设置得太低, 导致攻击之外的因素也会使得被检测的量超过阈值. 随着$k$值增大, 阈值变高, 误检率自然随之下降.

图 8 误检率${P_F}$

Fig. 8 False alarm rate ${P_F}$

下载: 全尺寸图片幻灯片

在图8中, 当$k = {{3}}$时, 误检率${P_F} = 0.0027$, 除此之外, 由图4可知, 欧氏检测法检测到攻击所用时间是$\Delta k = 50$; 由图8可知, 检测到攻击所用时间是$\Delta k = 30$, 检测用时更短.

综上, 对比三种检测方法可知, 本文提出的检测方法可以成功检测攻击不受限于传感器参数统计特性未知所带来的影响, 不受限于加性噪声和非加性噪声的情况, 对于实际系统的适用性更强.

5. 结束语

本文研究了智能电网系统中虚假数据注入攻击的检测问题. 针对非线性系统, 噪声统计特性未知的情况, 本文使用自适应平方根无迹卡尔曼滤波算法对系统内部状态和噪声作出估计. 针对传感器参数统计特性未知和非加性噪声的情况, 利用中心极限定理构造出符合正态分布的随机变量, 基于该随机变量提出了一种攻击检测方法, 并从评价指标的角度对算法进行分析, 该算法对于实际系统的适用性更强.

图 1 不同f 值的失败次数

Fig. 1 Number of failures for different f

下载: 全尺寸图片幻灯片

图 2 粒子在10维中的步长变化

Fig. 2 Step size changes of 10 dimensions of particles

下载: 全尺寸图片幻灯片

图 3 PSO和PSO+SW的种群多样性变化

Fig. 3 The population’s diversity of PSO and PSO + SW

下载: 全尺寸图片幻灯片

图 4 PSO和PSO + SW的收敛性

Fig. 4 The convergence performance of PSO and PSO + SW

下载: 全尺寸图片幻灯片

图 5 14个函数进化曲线图

Fig. 5 Evolution curves of 14 functions

下载: 全尺寸图片幻灯片

表 1 Penalized 2函数上使用不同的f值寻找全局最优值的迭代次数

Table 1 Numbers of iterations for finding the global optimum using different values of f on the Penalized 2 function

运行次数	0.01	0.009	0.008	0.007	0.006	0.005	0.004	0.003	0.002
1	11704	12273	11981	11625	10146	9985	11252	11198	×
2	2	2	2	2	1780	2	2	2	3056
3	13065	11979	12499	11119	2	9847	10795	9500	2
4	1405	2	3	2	10421	801	2	2	8887
5	2	13336	2	11896	1109	2	11033	×	2
6	14060	2	10895	1364	2	10132	3	9472	×
7	2758	14999	1106	1107	12761	2	1333	2	1222
8	865	2	2	2	2476	2	2	9418	426
9	2	12276	13873	11136	2	2	10078	2	2
10	14962	1407	2	2	10834	2	2	10544	×
11	1436	2109	13018	12386	1034	12742	11571	2	9060
12	2	2	2	2	2	2	2	10181	488
13	2	14577	13927	10019	2	11213	10855	1071	2
14	13346	2	2	1457	13428	2	1397	401	8615
15	1633	12661	11396	2	2	10691	2	2	2
16	2	2	2	12466	13939	2	9906	9438	8822
17	14325	13727	15113	2	2	10756	413	424	1319
18	2	1470	2	12205	12848	930	2	2	2
19	13871	2018	12345	2	2	2	11565	9610	×
20	3	505	2	11701	10586	11815	829	2	×
21	3	3	14717	2	2	978	2	11406	5039
22	2	2	2	2	12163	2	10889	601	529
23	×	14773	13978	11976	3	10870	2	2	517
24	14313	270	2	1368	2	2	10413	11753	639
25	927	2	10368	2	11707	11587	857	2	2
26	830	12428	2	12196	2	2	2	14083	9400
27	2	2	13887	2	11969	10566	12339	2	2
28	×	2	953	11848	2	2	2	9706	8173
29	14813	12621	2	2	12454	11169	9447	2	2
30	2	975	10766	12405	960	533	297	13482	14269
平均	4798	5148	6028	5277	5021	4488	4510	4907	3219

下载: 导出CSV

表 2 Levy函数上使用不同的f值寻找全局最优值的迭代次数

Table 2 Numbers of iterations for finding the global optimum using different values of f on the Levy function

运行次数	0.01	0.009	0.008	0.007	0.006	0.005	0.004	0.003	0.002
1	13874	14277	13198	10780	10421	13557	11056	9483	13787
2	2826	1183	2	1449	2	2	2	2	2
3	2	2	13689	2	11314	10167	13358	2	12695
4	×	13732	2042	11042	2	2	2	10550	2
5	×	2	2686	2	11160	11234	13494	1239	×
6	14837	13975	2	9906	2	2	1578	552	3200
7	5154	2	13005	2	11868	9903	2	2	2
8	2	14063	1244	10620	2	2	10440	9624	13760
9	12041	2	2	1362	×	10850	1376	2	2
10	2	×	12501	610	4469	603	2	2	12647
11	×	12550	2	2	2	2	10836	9370	620
12	12924	2	11695	11196	×	9632	2	2	2
13	2	13450	431	1730	4908	1492	9457	13545	13183
14	2	2	3306	2	783	2	2	2	805
15	×	12873	1861	10321	1413	9453	10773	9836	625
16	13729	3	2	2	2	2	701	2	2
17	2	2	×	12235	12042	10359	2	×	×
18	13234	14143	12546	2	2	2	8828	12012	10905
19	2780	2	2	2	10954	11268	621	375	2
20	2	14885	12626	12806	2	776	1180	2	8737
21	13953	2	2	2	2	2	2	×	2
22	4594	13092	13093	14701	11422	×	10400	9022	12969
23	2	4313	2	2	2	1381	2	2	446
24	13667	2	11992	11462	×	1442	9952	11454	2
25	2	14440	2	2	2434	2	2	3	×
26	×	754	13095	2	1446	10102	10023	2	13288
27	×	992	956	×	2	959	2	2	2
28	×	1009	2	11072	14639	2	10721	11631	8908
29	×	2	2	2	636	12115	2	2	2
30	14436	12920	12514	11453	2	2	2	10330	12355
平均	5887	5954	5259	4923	4072	4321	4494	4252	5146

下载: 导出CSV

表 3 26个基准测试函数

Table 3 Information of 26 benchmark functions

函数名	函数式	范围	最优解
Sphere	${f_1}(x) = \displaystyle\sum\limits_{i = 1}^D {{x_i}^2} $	[−100, 100]	0
Schwefel2.22	${f_2}\left( x \right) = \displaystyle\sum\limits_{i = 1}^D {\left\| { {x_i} } \right\| + \prod_{i = 1}^D {\left\| { {x_i} } \right\|} }$	[−10, 10]	0
Rosebrock	${f_3}(x) = \displaystyle\sum\limits_{i = 1}^D {\left[ {100{ {\left( { {x_{i + 1} } - x_i^2} \right)}^2} + { {\left( {1 - x_i^2} \right)}^2} } \right]}$	[−5, 10]	0
Quartic	${f_4}(x) = \displaystyle\sum\limits_{i = 1}^D {i{x_i^4} + random\left[ {0,1} \right)}$	[−1.28, 1.28]	0
Schwefel2.26	${f_5}(x) = 418.98288727243380{\cdot }D{\rm{ - } }\displaystyle\sum\limits_{i \;=\; 1}^D { {x_i}\sin \sqrt {\left\| { {x_i} } \right\|} }$$\overline { {f} } _5(x) = \displaystyle\sum\limits_{i = 1}^D { - {x_i}\sin \sqrt {\left\| { {x_i} } \right\|} },\; (D = 100)$	[−500, 500]	0 ~ 418.9826 × D
Rastrigin	${f_6}(x) = \displaystyle\sum\limits_{i = 1}^D {\left[ { {x_i^2} - 10\cos 2\pi {x_i} + 10} \right]}$	[−5.12, 5.12]	0
Ackley	${f_7}(x) = - 20{\rm{exp} }\left( { - 0.2\sqrt {\dfrac{1}{D}\displaystyle\sum\limits_{i \;=\; 1}^D { {x_i^2} } } } \right) - {\rm{exp} }\left( {\dfrac{1}{D}\displaystyle\sum\limits_{i \;=\; 1}^D {\cos 2\pi x_i} } \right) + 20 + {\rm{e}}$	[−50, 50]	0
Griewank	${f_8}(x) = \dfrac{1}{{4000}}\displaystyle\sum\limits_{i = 1}^D {{{\left( {{x_i}} \right)}^2}} - \prod _{i = 1}^D\cos \left( {\dfrac{{{x_i}}}{{\sqrt i }}} \right) + 1$	[−600, 600]	0
Penalized 1	${f_9}(x) = \dfrac{\pi }{D} \bigg\{ {\displaystyle\sum\limits_{i \;=\; 1}^{D - 1} {({y_i} - 1)} ^2}[1 + \sin (\pi {y_{i + 1} })] + {({y_D} - 1)^2} + 10{\sin ^2}(\pi {y_1}) \bigg\} +$$\displaystyle\sum\limits_{i = 1}^D {u({x_i},10,100,4),{y_i} = 1 + \dfrac{ { {x_i} + 1} }{4} }u({x_i},a,k,m) = \left\{ {\begin{aligned} &u({x_i},a,k,m),&{x_i} > a\quad\qquad \\ &0,& - a \le {x_i} \le a\;\\ &k{ {( - {x_i} - a)}^m},&{x_i} < - a \qquad\end{aligned} } \right.$	[−100, 100]	0
Penalized 2	${f_{10} }\left( x \right) = 0.1\bigg\{ { { {\sin }^2}\left( {\pi {x_1} } \right) + \displaystyle\sum\limits_{i \;=\; 1}^{D - 1} { { {\left( { {x_i} - 1} \right)}^2}[1 + { {\sin }^2}\left( {3\pi {x_{i + 1} } } \right)]+} }$${\left( { {x_D} - 1} \right)^2} {\left[ {1 + { {\sin }^2}\left( {2\pi {x_{i + 1} } } \right)} \right]} \bigg\} + \displaystyle\sum\limits_{i = 1}^D {u\left( { {x_i},5,100,4} \right)}$	[−100, 100]	0
Rotated Schwefel2.26	${f_{11} }(x) = 4.18.928\cdot D - \displaystyle\sum\limits_{i \;=\; 1}^D { {z_i} } , where \; {z_i} = \left\{ \begin{aligned} & - {y_i}\sin (\sqrt {\|{y_i}\|} ) if\|{y_i}\| \le 500 \\ & 0, \quad {\rm{否则} } \end{aligned} \right.$$ {y_i} = y_i^` + 420.96, y_i^` = M \cdot(x - 420.96), M \;is \; an \; orthogonal \; matrix $	[−500, 500]	0
Rotated Rastrigin	${f_{12} }(x) = \displaystyle\sum\limits_{i \;=\; 1}^D {[y_i^2 - 10\cos (2\pi {y_i}) + 10], where\;y = M \cdot x,\;} M \;is\; an\; orthogonal\; matrix$	[−5.12, 5.12]	0
Rotated Ackly	${f_{13} }(x) = - 20\exp \left[ - 0.2\sqrt {\dfrac{1}{D}\displaystyle\sum\limits_{i \;=\; 1}^D {y_i^2} }\right ] - \exp \left [\dfrac{1}{D}\displaystyle\sum\limits_{i \;=\; 1}^D {\cos (2\pi {y_i})} \right ] + (20 + {\rm{e} })$$where\;y = M \cdot x , M \;is\; an\; orthogonal\; matrix$	[−32, 32]	0
Rotated Griewank	${f_{14} }(x) = \dfrac{1}{ {4000} }\displaystyle\sum\limits_{i = 1}^D {y_i^2 - \prod\limits_{i = 1}^D {\cos (\dfrac{ { {y_i} } }{ {\sqrt i } })} + 1, where\; y = M \cdot x,\;} M \;is\; an \; orthogonal \; matrix$	[−600, 600]	0
Elliptic	${f_{15}}\left( x \right) = {\displaystyle\sum\limits_{i = 1}^D {\left( {{{10}^6}} \right)} ^{\dfrac{{i - 1}}{{D - 1}}}}x_{_i}^2$	[−100, 100]	0
SumSquare	${f_{16}}\left( x \right) = \displaystyle\sum\limits_{i = 1}^D {ix_{_i}^2} $	[−10, 10]	0
SumPower	${f_{17} }\left( x \right) = {\displaystyle\sum\limits_{i \;=\; 1}^D {\left\| { {x_i} } \right\|} ^{\left( {i + 1} \right)} }$	[−1, 1]	0
Schwefel2.21	${f_{18} }(x) ={\rm max}\left\{ {\left\| { {x_i} } \right\|,1 \le i \le D} \right\}$	[−100, 100]	0
Step	${f_{19} }(x) = \displaystyle\sum\limits_{i \;=\; 1}^D {\left\lfloor { {x_i} + 0.5} \right\rfloor }$	[−100, 100]	0
Exponential	${f_{20} }\left( x \right) = \exp \left( {0.5\displaystyle\sum\limits_{i\; =\; 1}^D { {x_i} } } \right)$	[−10, 10]	0
NCRastrigin	${f_{21} }\left( x \right) = \displaystyle\sum\limits_{i\; =\; 1}^D {\left[ {y_i^2 - 10\cos \left( {2\pi {y_i} } \right) + 10} \right]}$	[−5.12, 5.12]	0
Alpine	${f_{22} }\left( x \right) = \displaystyle\sum\limits_{i\; = \;1}^{D - 1} {\left\| { {x_i}\sin \left( { {x_i} } \right) + 0.1{x_i} } \right\|}$	[−10, 10]	0
Levy	${f_{23}}\left( x \right) = \displaystyle\sum\limits_{i = 1}^{D - 1} {{{\left( {{x_i} - 1} \right)}^2}\left[ {1 + {{\sin }^2}\left( {3\pi {x_{i + 1}}} \right)} \right] + {{\sin }^2}\left( {3\pi {x_1}} \right)} + \left\| {{x_D} - 1} \right\|\left[ {1 + {{\sin }^2}\left( {3\pi {x_D}} \right)} \right]$	[−10, 10]	0
Weierstrass	${f_{24} }\left( x \right) = \displaystyle\sum\limits_{i \;=\; 1}^D {\left( {\displaystyle\sum\limits_{k \;=\; 0}^{k_{max} } {\left[ { {a^k}\cos \left( {2\pi {b^k}\left( { {x_i} + 0.5} \right)} \right)} \right]} } \right)} - D\displaystyle\sum\limits_{k = 0}^{k_{max} } {\left[ { {a^k}\cos \left( {2\pi {b^k}0.5} \right)} \right]}$$a = 0.5;b = 3;k_{max} = 20$	[−1, 1]	0
Himmelballa	${f_{25}}\left( x \right) = \dfrac{1}{D}\displaystyle\sum\limits_1^D {\left( {x_{_i}^4 - 16x_i^2 + 5{x_i}} \right)} $	[−5, 5]	−78.332
Michalewice	${f_{26} }\left( x \right) = - \displaystyle\sum\limits_{i \;=\; 1}^D {\sin \left( { {x_i} } \right){ {\sin }^{20} }\left( {\dfrac{ {ix_i^2} }{\pi } } \right)}$	$\left[ {0,\pi } \right]$	−50, $D=50 $

下载: 导出CSV

表 4 30维实验结果对比

Table 4 Comparison of 30-dimensional test results

函数	指标	CLPSO	HPSO-TVAC	DMS-PSO	LFPSO	PSOLF	RPSOLF	H-PSO-SCAC	SWTPSO
f₁	MBF	8.06 × 10⁻⁹⁶	2.83 × 10⁻³³	1.53 × 10⁻¹¹³	4.69 × 10⁻³¹	0	0	0	0
	SD	3.53 × 10⁻⁹⁵	3.19 × 10⁻³³	5.14 × 10⁻¹¹³	2.50 × 10⁻³⁰	0	0	0	0
f₂	MBF	5.50 × 10⁻⁵⁷	9.03 × 10⁻²⁰	2.18 × 10²	2.64 × 10⁻¹⁷	0	0	4.09 × 10⁻²¹⁷	0
	SD	1.92 × 10⁻⁵⁶	9.58 × 10⁻²⁰	1.83 × 10²	6.92 × 10⁻¹⁷	0	0	0	0
f₃	MBF	4.18 × 10¹	2.39 × 10¹	3.49 × 10¹	2.38 × 10¹	2.68 × 10¹	1.01 × 10¹	2.36 × 10¹	3.99
	SD	3.35 × 10¹	2.65 × 10¹	2.76 × 10¹	3.17 × 10⁻¹	1.03	9.69 × 10⁻¹	1.51 × 10⁻¹	2.22 × 10¹
f₄	MBF	1.74 × 10⁻³	9.82 × 10⁻²	6.09 × 10⁻⁴	2.41 × 10⁻³	2.60 × 10⁻⁵	6.50 × 10⁻³	3.94 × 10⁻²³⁷	8.89 × 10⁻⁴
	SD	7.83 × 10⁻⁴	3.26 × 10⁻²	4.18 × 10⁻⁴	8.07 × 10⁻⁴	2.10 × 10⁻⁵	5.50 × 10⁻³	0	1.58 × 10⁻³
f₅	MBF	3.82 × 10⁻⁴	1.59 × 10³	3.21 × 10³	1.37 × 10³	2.00 × 10³	2.85 × 10³	6.64 × 10³	3.51 × 10²
	SD	1.28 × 10⁻⁷	3.26 × 10²	6.51 × 10²	6.36 × 10²	6.08 × 10²	3.79 × 10²	2.39 × 10³	1.14 × 10³
f₆	MBF	1.27 × 10¹	9.43	1.49 × 10¹	4.54	0	0	0	0
	SD	4.22	3.48	3.62	1.03 × 10¹	0	0	0	0
f₇	MBF	1.42 × 10⁻¹⁴	7.29 × 10⁻¹⁴	6.06 × 10⁻¹²	1.68 × 10⁻¹⁴	8.88 × 10⁻¹⁶	8.88 × 10⁻¹⁶	8.88 × 10⁻¹⁶	5.88 × 10⁻¹⁶
	SD	7.46 × 10⁻¹⁵	3.00 × 10⁻¹⁴	3.90 × 10⁻¹³	4.84 × 10⁻¹⁵	0	0	0	0
f₈	MBF	3.20 × 10⁻³	9.75 × 10⁻³	1.85 × 10⁻³	8.14 × 10⁻¹⁷	0	0	0	0
	SD	4.93 × 10⁻³	8.33 × 10⁻³	4.07 × 10⁻³	4.46 × 10⁻¹⁶	0	0	0	0
f₉	MBF	1.36 × 10⁻³³	2.71 × 10⁻²⁹	0	4.67 × 10⁻³¹	1.66 × 10⁻²	1.98 × 10⁻³²	2.10 × 10⁻¹	1.58 × 10⁻³²
	SD	2.82 × 10⁻³³	1.88 × 10⁻²⁹	0	9.01 × 10⁻³¹	1.27 × 10⁻²	9.23 × 10⁻³³	9.67 × 10⁻¹	5.07 × 10⁻³³
f₁₀	MBF	1.65 × 10⁻³³	2.79 × 10⁻²⁸	6.16 × 10⁻³⁵	1.51 × 10⁻²⁸	9.01 × 10⁻⁹	1.65 × 10⁻³²	1.09	1.36 × 10⁻³²
	SD	4.03 × 10⁻³³	2.18 × 10⁻²⁸	2.76 × 10⁻³⁴	8.00 × 10⁻²⁸	1.11 × 10⁻⁸	4.36 × 10⁻³³	1.42	4.84 × 10⁻³³
f₁₁	MBF	4.39 × 10³	5.32 × 10³	4.04 × 10³	5.51 × 10³	1.54 × 10³	9.98 × 10³	4.48 × 10³	6.63 × 10³
	SD	3.51 × 10²	7.00 × 10²	5.68 × 10²	5.64 × 10²	4.28 × 10²	7.58 × 10⁻¹²	1.73 × 10²	1.43 × 10³
f₁₂	MBF	8.17 × 10¹	5.29 × 10³	4.20 × 10¹	1.79	0	0	0	0
	SD	1.08 × 10¹	1.25 × 10¹	9.74	9.81	0	0	0	0
f₁₃	MBF	5.91 × 10⁻⁵	9.29	2.42 × 10⁻¹⁴	1.65 × 10⁻¹⁴	0	0	0	0
	SD	6.46 × 10⁻⁵	2.07	1.52 × 10⁻¹⁴	5.40 × 10⁻¹⁵	0	0	0	0
f₁₄	MBF	7.69 × 10⁻⁵	9.26 × 10⁻³	1.02 × 10⁻²	1.48 × 10⁻³	0	0	0	0
	SD	7.66 × 10⁻⁵	8.80 × 10⁻³	1.24 × 10⁻²	6.17 × 10⁻³	0	0	0	0

下载: 导出CSV

表 6 100维实验结果对比

Table 6 Comparison of 100-dimensional test results

函数	指标	CLPSO	HPSO-TVAC	DMS-PSO	LFPSO	PSOLF	RPSOLF	H-PSO-SCAC	SWTPSO
f₁	MBF	4.16 × 10⁻⁷⁵	5.48 × 10⁻²⁶	4.89 × 10⁻⁷⁷	2.65 × 10⁻¹⁰	0	0	0	0
	SD	1.80 × 10⁻⁷⁴	2.59 × 10⁻²⁵	1.40 × 10⁻⁷⁶	3.89 × 10⁻⁸	0	0	0	0
f₂	MBF	3.45 × 10⁻⁴⁵	3.77 × 10⁻¹⁵	5.83 × 10²	1.02	0	1.22 × 10⁻²⁵²	0	0
	SD	1.52 × 10⁻⁴⁴	1.20 × 10⁻¹⁴	3.54 × 10¹	3.21	0	0	0	0
f₃	MBF	1.46 × 10²	1.34 × 10²	1.20 × 10²	2.56 × 10¹	9.04 × 10¹	9.15 × 10¹	9.27 × 10¹	1.03 × 10¹
	SD	4.78 × 10¹	2.28 × 10²	3.62 × 10¹	5.23 × 10⁻¹	2.57 × 10¹	1.50	4.31	2.55 × 10¹
f₄	MBF	7.00 × 10⁻³	1.17 × 10⁻⁴⁸	6.40 × 10⁻³	3.21 × 10⁻³	2.35 × 10⁻⁵	8.68 × 10⁻³	3.25 × 10⁻¹²⁵	6.28 × 10⁻⁴
	SD	1.53 × 10⁻³	1.27 × 10⁻⁴⁷	2.56 × 10⁻³	3.24 × 10⁻³	1.92 × 10⁻⁵	5.47 × 10⁻³	0	1.60 × 10⁻³
f₅	MBF	5.36	2.30 × 10³	6.58 × 10⁴	3.56 × 10⁴	2.46 × 10⁴	1.29 × 10³	2.66 × 10⁴	6.16 × 10²
	SD	2.35 × 10¹	2.88 × 10³	5.21 × 10¹	3.56 × 10³	5.75 × 10³	2.05 × 10³	4.11 × 10³	2.08 × 10³
f₆	MBF	7.02	4.90 × 10¹	1.95 × 10¹	3.25 × 10¹	0	0	0	0
	SD	1.00 × 10²	4.21 × 10¹	2.59 × 10¹	5.68 × 10¹	0	0	0	0
f₇	MBF	2.74 × 10⁻¹⁴	3.21 × 10⁻¹²	9.91 × 10⁻¹	3.87 × 10⁻⁹	5.89 × 10⁻¹⁶	5.89 × 10⁻¹⁶	5.89 × 10⁻¹⁶	5.88 × 10⁻¹⁶
	SD	5.17 × 10⁻¹⁵	4.29 × 10⁻¹¹	4.43	5.89 × 10⁻⁹	0	0	0	0
f₈	MBF	3.33 × 10⁻¹⁷	2.71 × 10⁻³	7.40 × 10⁻⁴	3.89 × 10⁻²	0	0	0	0
	SD	7.29 × 10⁻¹⁷	2.32 × 10⁻²	2.28 × 10⁻³	5.78 × 10⁻²	0	0	0	0
f₉	MBF	9.33 × 10⁻³	2.39 × 10⁻²⁴	1.56 × 10⁻²	9.56 × 10⁻²	2.15 × 10⁻²	1.27 × 10⁻¹²	5.47 × 10⁻¹	1.90 × 10⁻³¹
	SD	2.87 × 10⁻²	7.19 × 10⁻²⁴	4.23 × 10⁻²	3.46 × 10⁻²	3.26 × 10⁻²	1.24 × 10⁻¹¹	7.08 × 10⁻¹	5.09 × 10⁻³⁰
f₁₀	MBF	1.10 × 10⁻³	1.74 × 10⁻²¹	5.49 × 10⁻³	8.97 × 10⁻²	6.32 × 10⁻³	3.25 × 10⁻³	4.50	2.80 × 10⁻³⁰
	SD	3.38 × 10⁻³	7.05 × 10⁻²¹	1.26 × 10⁻²	2.65 × 10⁻²	2.53 × 10⁻²	9.57 × 10⁻²	6.71	3.77 × 10⁻²⁹
f₁₁	MBF	6.32 × 10³	2.85 × 10⁴	1.12 × 10⁴	3.56 × 10⁴	4.18 × 10⁴	3.43 × 10⁴	2.07 × 10⁴	2.99 × 10⁴
	SD	2.32 × 10²	5.14 × 10³	3.65 × 10³	3.25 × 10³	1.67 × 10²	3.57 × 10³	8.63 × 10²	7.92 × 10³
f₁₂	MBF	2.56 × 10³	5.47 × 10²	9.32 × 10²	1.25 × 10²	0	0	0	0
	SD	5.32 × 10²	1.03 × 10²	1.35 × 10²	6.32 × 10¹	0	0	0	0
f₁₃	MBF	6.87 × 10⁻³	1.05 × 10⁻⁶	3.56 × 10⁻⁶	6.89 × 10⁻⁶	0	0	0	0
	SD	3.98 × 10⁻³	2.54 × 10⁻⁶	3.89 × 10⁻⁵	1.58 × 10⁻⁵	0	0	0	0
f₁₄	MBF	2.32 × 10⁻²	3.89 × 10⁻³	8.25 × 10⁻¹	9.58 × 10⁻²	0	0	0	0
	SD	1.25 × 10⁻²	5.89 × 10⁻³	9.23 × 10⁻¹	2.56 × 10⁻²	0	0	0	0

下载: 导出CSV

表 5 50维实验结果对比

Table 5 Comparison of 50-dimensional test results

函数	指标	CLPSO	HPSO-TVAC	DMS-PSO	LFPSO	PSOLF	RPSOLF	H-PSO-SCAC	SWTPSO
f₁	MBF	2.28 × 10⁻⁸⁵	5.50 × 10⁻¹¹	8.62 × 10⁻¹⁰³	9.17 × 10⁻¹⁷	0	6.57 × 10⁻²⁰¹	0	0
	SD	3.87 × 10⁻⁸⁵	1.42 × 10⁻¹⁰	2.82 × 10⁻¹⁰²	3.20 × 10⁻¹⁶	0	0	0	0
f₂	MBF	4.47 × 10⁻⁵²	1.12 × 10⁻⁷	3.09 × 10²	8.00 × 10⁻¹	0	8.01 × 10⁻¹²⁴	4.69 × 10⁻²⁴⁵	0
	SD	1.22 × 10⁻⁵¹	1.76 × 10⁻⁷	2.87 × 10¹	4.44	0	2.36 × 10⁻¹²²	0	0
f₃	MBF	6.97 × 10¹	1.36 × 10²	7.03 × 10¹	4.41 × 10¹	4.74 × 10¹	4.41 × 10¹	4.71 × 10¹	7.24
	SD	4.45 × 10¹	5.90 × 10²	4.49 × 10¹	2.82 × 10⁻¹	9.84 × 10⁻¹	4.74 × 10¹	6.75	5.33 × 10¹
f₄	MBF	3.31 × 10⁻³	3.78 × 10⁻²³	1.15 × 10⁻³	2.37 × 10⁻³	1.60 × 10⁻⁵	3.26 × 10⁻³	8.85 × 10⁻¹⁷⁶	9.08 × 10⁻⁴
	SD	9.05 × 10⁻⁴	1.39 × 10⁻²²	3.92 × 10⁻⁴	1.01 × 10⁻³	1.17 × 10⁻⁵	5.36 × 10⁻³	0	2.68 × 10⁻³
f₅	MBF	9.53	9.83 × 10²	9.76 × 10³	4.36 × 10³	5.07 × 10³	5.96 × 10²	1.20 × 10⁴	5.38 × 10²
	SD	4.76 × 10¹	1.76 × 10³	5.62 × 10²	1.27 × 10³	9.74 × 10²	1.43 × 10³	2.50 × 10³	2.01 × 10³
f₆	MBF	2.42 × 10¹	2.31 × 10¹	2.99 × 10¹	1.26 × 10¹	0	0	0	0
	SD	6.40	3.06 × 10¹	6.30	1.61 × 10¹	0	0	0	0
f₇	MBF	1.01 × 10⁻¹⁴	1.53 × 10⁻⁴	7.75 × 10⁻⁸	6.53 × 10⁻¹⁰	8.88 × 10⁻¹⁶	5.89 × 10⁻¹⁶	8.88×10⁻¹⁶	5.89 × 10⁻¹⁶
	SD	3.32 × 10⁻¹⁵	1.91 × 10⁻³	2.79 × 10⁻⁷	2.46 × 10⁻⁹	0	0	0	0
f₈	MBF	4.93 × 10⁻⁴	4.50 × 10⁻⁴	7.40 × 10⁻⁴	4.86 × 10⁻³	0	0	0	0
	SD	2.20 × 10⁻³	2.80 × 10⁻²	2.28 × 10⁻³	1.36 × 10⁻²	0	0	0	0
f₉	MBF	3.11 × 10⁻³	2.88 × 10⁻⁷	3.11 × 10⁻³	1.94 × 10⁻¹⁷	3.15 × 10⁻²	2.26 × 10⁻¹⁴	2.96 × 10⁻¹	1.15 × 10⁻³²
	SD	1.39 × 10⁻²	7.95 × 10⁻⁶	1.39 × 10⁻²	7.21 × 10⁻¹⁷	1.66 × 10⁻²	3.46 × 10⁻¹³	5.03 × 10⁻¹	1.79 × 10⁻³²
f₁₀	MBF	5.49 × 10⁻⁴	3.49 × 10⁻⁸	8.63 × 10⁻³⁴	2.39 × 10⁻⁴	4.23 × 10⁻⁸	9.67 × 10⁻¹³	2.21	1.56 × 10⁻³²
	SD	2.46 × 10⁻³	4.85 × 10⁻⁷	2.30 × 10⁻³³	1.69 × 10³	5.84 × 10⁻⁸	1.12 × 10⁻¹¹	3.38	2.39 × 10⁻³²
f₁₁	MBF	5.96 × 10³	1.07 × 10⁴	1.78 × 10³	5.72 × 10³	1.79 × 10³	1.32 × 10⁴	8.71 × 10³	1.14 × 10⁴
	SD	6.22 × 10²	2.85 × 10³	6.03 × 10²	1.64 × 10³	4.44 × 10²	2.93 × 10³	5.30 × 10²	3.27 × 10³
f₁₂	MBF	1.41 × 10³	2.61 × 10²	4.91 × 10²	2.56 × 10²	0	0	0	0
	SD	2.35 × 10²	6.73 × 10¹	1.20 × 10²	3.44 × 10¹	0	0	0	0
f₁₃	MBF	2.07 × 10⁻⁴	7.31 × 10⁻⁹	1.81 × 10⁻⁸	7.97 × 10⁻⁹	0	0	0	0
	SD	5.52 × 10⁻³	9.51 × 10⁻⁹	4.09 × 10⁻⁷	2.99 × 10⁻⁸	0	0	0	0
f₁₄	MBF	5.22 × 10⁻³	4.58 × 10⁻³	4.11 × 10⁻¹	3.74 × 10⁻³	0	0	0	0
	SD	1.19 × 10⁻²	6.37 × 10⁻³	3.72 × 10⁻¹	7.97 × 10⁻³	0	0	0	0

下载: 导出CSV

表 7 各PSO算法的Friedman 检验结果

Table 7 Friedman test results of each POS algorithm

算法	D = 30 (排名)	D = 50 (排名)	D = 100 (排名)
SWTPSO	2.86 (1)	2.46 (1)	2.61 (1)
H-PSO-SCAC	3.86 (4)	3.93 (4)	3.71 (4)
RPSOLF	3.79 (3)	3.46 (3)	3.46 (2)
PSOLF	3.57 (2)	3.39 (2)	3.68 (3)
LFPSO	5.43 7)	5.25 (5)	6.50 (8)
DMS-PSO	5.32 (6)	5.96 (7)	6.04 (7)
HPSO-TVAC	6.50 (8)	5.57 (6)	5.04 (6)
CLPSO	4.68 (5)	5.96 (7)	4.96 (5)

下载: 导出CSV

表 8 SWTPSO与9种相关算法实验结果对比

Table 8 Comparison of experimental results of SWTPSO and nine related algorithms

函数	指标	BABC	CoDE	JADE	jDEscop	SaDE	AABCLS	CMAES	ABCVSS	MPGABC	SWTPSO
f₁	MBF	1.01 × 10⁻¹⁴	1.16 × 10⁻³⁷	2.50 × 10⁻⁸⁷	1.15 × 10⁻³⁶	1.48 × 10⁻⁵⁷	4.60 × 10⁻³⁵	1.21 × 10⁻²⁸	6.68 × 10⁻²³	2.47 × 10⁻⁵¹	0
	SD	5.07 × 10⁻¹⁴	2.69 × 10⁻³⁷	8.62 × 10⁻⁸⁷	5.75 × 10⁻³⁶	6.02 × 10⁻⁵⁷	1.67 × 10⁻³⁵	2.09 × 10⁻²⁹	3.16 × 10⁻³²	9.52 × 10⁻⁵¹	0
f₂	MBF	3.14 × 10⁻²⁹	1.25 × 10⁻³⁴	4.43 × 10⁻⁷⁸	4.77 × 10⁻³⁶	8.97 × 10⁻⁵⁵	4.48 × 10⁻³³	4.88 × 10⁻²³	1.11 × 10⁻²³	4.54 × 10⁻⁴⁹	0
	SD	4.95 × 10⁻²⁹	1.50 × 10⁻³⁴	2.21 × 10⁻⁷⁷	1.85 × 10⁻³⁵	2.38 × 10⁻⁵⁴	2.16 × 10⁻³³	9.20 × 10⁻²⁴	4.17 × 10⁻²³	1.32 × 10⁻⁴⁸	0
f₃	MBF	1.23 × 10⁻¹³	1.94 × 10⁻³⁸	2.27 × 10⁻⁸⁵	1.80 × 10⁻⁴¹	5.22 × 10⁻⁵⁸	2.54 × 10⁻³⁵	4.22 × 10⁻²⁷	1.77 × 10⁻³³	6.05 × 10⁻⁵³	0
	SD	6.16 × 10⁻¹³	4.08 × 10⁻³⁸	1.14 × 10⁻⁸²	1.47 × 10⁻¹⁶	1.49 × 10⁻⁵⁷	1.24 × 10⁻³⁵	9.12 × 10⁻²⁸	7.88 × 10⁻³³	1.39 × 10⁻⁵²	0
f₄	MBF	2.96 × 10⁻⁹²	6.20 × 10⁻¹⁴³	3.01 × 10⁻¹⁰⁰	3.58 × 10⁻¹²⁴	3.38 × 10⁻⁸⁰	1.52 × 10⁻⁵⁰	6.82 × 10⁻¹³	7.21 × 10⁻⁴¹	8.51 × 10⁻¹⁰³	0
	SD	0	3.15 × 10⁻¹⁴²	1.00 × 10⁻⁹⁹	1.79 × 10⁻¹²³	1.59 × 10⁻⁷⁹	4.88 × 10⁻⁵⁰	5.39 × 10⁻¹³	3.56 × 10⁻⁴⁰	3.61 × 10⁻¹⁰²	0
f₅	MBF	2.18 × 10⁻⁶	1.21 × 10⁻²⁰	1.89 × 10⁻⁴¹	1.71 × 10⁻²²	1.71 × 10⁻³⁸	4.98 × 10⁻¹⁸	3.69 × 10⁻⁴	2.89 × 10⁻¹⁸	2.95 × 10⁻²⁸	1.13 × 10⁻³⁰⁰
	SD	1.06 × 10⁻⁵	8.94 × 10⁻²¹	8.32 × 10⁻⁴¹	8.46 × 10⁻²²	7.47 × 10⁻³⁸	1.34 × 10⁻¹⁸	1.84 × 10⁻³	6.52 × 10⁻¹⁸	4.04 × 10⁻²⁸	0
f₆	MBF	7.20	6.00 × 10⁻⁵	8.20 × 10⁻¹⁰	2.52	3.42 × 10⁻¹	1.16 × 10⁻¹	6.00 × 10⁻¹⁵	2.06	2.12 × 10¹	0
	SD	2.78	9.51 × 10⁻⁵	9.51 × 10⁻¹⁰	5.61 × 10⁻¹	4.67 × 10⁻¹	1.24 × 10⁻¹	7.45 × 10⁻¹⁶	3.77 × 10⁻¹	4.04	0
f₇	MBF	0	0	0	0	0	0	1.60 × 10⁻¹	0	0	0
	SD	0	0	0	0	0	0	4.73 × 10⁻¹	0	0	0
f₈	MBF	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	5.40 × 10⁻⁶⁶	2.67 × 10⁻¹⁰⁹	2.67 × 10⁻¹⁰⁹	2.96 × 10⁻¹⁰⁹
	SD	2.62 × 10⁻¹¹⁶	9.67 × 10⁻¹²⁵	9.60 × 10⁻¹²⁵	9.65 × 10⁻¹²⁵	3.71 × 10⁻¹²²	2.77 × 10⁻¹¹⁹	1.99 × 10⁻⁶⁵	3.08 × 10⁻¹²⁰	3.06 × 10⁻¹²²	1.14 × 10⁻¹¹¹
f₉	MBF	5.74 × 10⁻²	8.17 × 10⁻³	2.50 × 10⁻³	3.91 × 10⁻³	1.33 × 10⁻²	1.63 × 10⁻²	2.80 × 10⁻¹	6.14 × 10⁻²	5.02 × 10⁻²	6.61 × 10⁻⁴
	SD	1.11 × 10⁻²	2.79 × 10⁻³	1.54 × 10⁻³	1.47 × 10⁻³	3.33 × 10⁻³	4.68 × 10⁻³	6.70 × 10⁻²	1.38 × 10⁻²	7.98 × 10⁻³	1.37 × 10⁻²
f₁₀	MBF	6.29 × 10⁻²	3.32 × 10¹	3.19 × 10⁻¹	2.58 × 10¹	9.45	3.09	1.75 × 10⁻²⁵	1.09 × 10⁻¹	4.32 × 10⁻¹	8.19
	SD	1.24 × 10⁻¹	2.26 × 10¹	1.10	2.68 × 10¹	2.08 × 10¹	1.36 × 10¹	4.31 × 10⁻²⁶	2.52 × 10⁻¹	7.44 × 10⁻¹	5.68 × 10¹
f₁₅	MBF	0	7.34 × 10⁻¹	1.78 × 10⁻¹¹	1.03 × 10⁻¹⁴	6.77 × 10⁻¹	0	3.89 × 10²	0	0	0
	SD	0	8.82 × 10⁻¹	1.74 × 10⁻¹¹	1.31 × 10⁻¹⁴	6.87 × 10⁻¹	0	7.11 × 10¹	0	0	0
f₁₆	MBF	0	2.28 × 10¹	2.74 × 10⁻⁸	8.02 × 10⁻²	4.40 × 10⁻¹	8.00 × 10⁻²	3.78 × 10²	0	0	0
	SD	0	4.68	2.07 × 10⁻⁸	2.77 × 10⁻¹	6.51 × 10⁻¹	2.77 × 10⁻¹	4.90 × 10¹	0	0	0
f₁₇	MBF	0	2.96 × 10⁻⁴	7.88 × 10⁻⁴	1.47 × 10⁻¹⁶	6.88 × 10⁻³	4.44 × 10⁻¹⁸	1.38 × 10⁻³	3.67 × 10⁻¹⁴	0	0
	SD	0	1.48 × 10⁻³	2.82 × 10⁻³	2.60 × 10⁻¹⁶	1.22 × 10⁻²	2.22 × 10⁻¹⁷	3.34 × 10⁻³	1.84 × 10⁻¹³	0	0
f₁₈	MBF	3.78 × 10⁻¹²	4.74	4.75	6.18 × 10¹	3.64 × 10⁻¹²	1.41 × 10⁻¹¹	9.01 × 10³	4.66 × 10⁻¹²	1.05 × 10⁻¹¹	5.38 × 10²
	SD	7.28 × 10⁻¹³	2.37 × 10¹	2.32 × 10¹	1.69 × 10²	0	3.53 × 10⁻¹²	1.13 × 10³	2.47 × 10⁻¹²	3.03 × 10⁻¹²	2.01 × 10³
f₁₉	MBF	7.50 × 10⁻¹⁵	2.81 × 10⁻¹⁵	6.22 × 10⁻¹⁵	1.43 × 10¹	1.32	2.65 × 10⁻¹⁴	1.99 × 10¹	1.70 × 10⁻¹⁴	2.51 × 10⁻¹⁴	5.89 × 10⁻¹⁶
	SD	2.49 × 10⁻¹⁵	7.11 × 10⁻¹⁶	0	6.50	4.90 × 10⁻¹	3.48 × 10⁻¹⁵	2.63 × 10⁻²	5.75 × 10⁻³³	4.55 × 10⁻¹⁵	0
f₂₀	MBF	1.22 × 10⁻¹³	9.42 × 10⁻³³	2.49 × 10⁻³	1.91 × 10⁻³¹	3.24 × 10⁻²	9.42 × 10⁻³³	4.98 × 10⁻³	1.07 × 10⁻³²	9.42 × 10⁻³³	1.15 × 10⁻³²
	SD	6.09 × 10⁻¹³	1.40 × 10⁻⁴⁸	1.24 × 10⁻²	5.28 × 10⁻³¹	9.01 × 10⁻²	1.40 × 10⁻⁴⁸	1.72 × 10⁻²	5.74 × 10⁻³³	1.40 × 10⁻⁴⁸	1.79 × 10⁻³²
f₂₁	MBF	2.50 × 10⁻¹⁵	1.55 × 10⁻³³	1.01 × 10⁻⁴	5.02 × 10⁻³²	1.09 × 10⁻²	1.50 × 10⁻³³	8.01 × 10³	1.80 × 10⁻³³	1.50 × 10⁻³³	1.56 × 10⁻³²
	SD	1.25 × 10⁻¹⁴	2.47 × 10⁻³⁴	8.05 × 10⁻⁵	5.95 × 10⁻³²	3.00 × 10⁻²	0	6.25 × 10³	1.08 × 10⁻³³	0	2.39 × 10⁻³²
f₂₂	MBF	2.07 × 10⁻¹⁶	4.39 × 10⁻³	1.01 × 10⁻⁴	2.95 × 10⁻⁵	1.75 × 10⁻¹⁶	5.32 × 10⁻¹¹	8.59 × 10⁻¹	2.14 × 10⁻¹⁶	3.21 × 10⁻⁷	3.10 × 10⁻²⁹⁰
	SD	7.79 × 10⁻¹⁶	6.44 × 10⁻³	8.05 × 10⁻⁵	3.09 × 10⁻⁵	2.63 × 10⁻¹⁶	9.68 × 10⁻¹¹	8.49 × 10⁻¹	7.43 × 10⁻¹⁶	4.58 × 10⁻⁷	0
f₂₃	MBF	1.35 × 10⁻³¹	1.35 × 10⁻³¹	1.35 × 10⁻³¹	1.30 × 10⁻³⁰	7.39 × 10⁻²	1.35 × 10⁻³¹	3.77 × 10⁻¹	1.35 × 10⁻³¹	1.35 × 10⁻³¹	1.78 × 10⁻³¹
	SD	2.23 × 10⁻⁴⁷	2.47 × 10⁻³³	2.23 × 10⁻⁴⁷	3.69 × 10⁻³⁰	1.07 × 10⁻¹	2.23 × 10⁻⁴⁷	1.22	2.23 × 10⁻⁴⁷	2.23 × 10⁻⁴⁷	2.83 × 10⁻³¹
f₂₄	MBF	0	3.45	3.33 × 10⁻¹	0	4.17 × 10⁻⁴	2.56 × 10⁻⁶	9.41	0	1.56 × 10⁻²	0
	SD	0	2.94 × 10⁻¹	4.45 × 10⁻²	0	5.94 × 10⁻⁴	7.39 × 10⁻⁶	4.19	0	2.04 × 10⁻²	0
f₂₅	MBF	−7.83 × 10¹	−7.83 × 10¹	−7.83 × 10¹	−7.83 × 10¹	−7.83 × 10¹	−7.83 × 10¹	−6.43 × 10¹	−7.83 × 10¹	−7.83 × 10¹	−6.26 × 10¹
	SD	1.16 × 10⁻¹⁴	4.10 × 10⁻¹⁵	2.26 × 10⁻¹	1.12 × 10⁻¹³	1.13 × 10⁻¹	4.10 × 10⁻¹⁵	2.63	1.00 × 10⁻¹⁴	1.12 × 10⁻¹⁴	3.31
f₂₆	MBF	−5.00 × 10¹	−4.86 × 10¹	−4.98 × 10¹	−5.00 × 10¹	−4.83 × 10¹	−5.00 × 10¹	−4.10 × 10¹	−5.00 × 10¹	−5.00 × 10¹	−3.99 × 10¹
	SD	6.79 × 10⁻⁶	4.45 × 10⁻¹	3.86 × 10⁻²	9.23 × 10⁻³	2.70 × 10⁻¹	9.13 × 10⁻⁴	2.65	3.82 × 10⁻⁷	4.26 × 10⁻⁴	2.13 × 10¹

下载: 导出CSV

表 9 8种新相关算法和SWTPSO的Friedman测试结果

Table 9 Friedman test results of 8 new correlation algorithms and SWTPSO

算法	SWTPSO	MPGABC	AABCLS	JADE	ABCVSS	BABC	CoDE	jDEscop	SaDE	CMAES
平均值 (排名)	3.86 (1)	4.41 (2)	5.11 (3)	5.05 (4)	5.07 (5)	5.11 (6)	5.59 (7)	5.91 (8)	6.07 (9)	8.82 (10)

下载: 导出CSV

表 10 6种算法在15个函数上的比较结果

Table 10 Results of the fifteen functions for six algorithms

函数	指标	GDHS	ODE	IGHS	GABC	IGPSO	SWTPSO
f₁	MBF	8.07 × 10⁻³	5.67 × 10⁻⁸¹	2.14 × 10⁻¹⁰	8.58 × 10²	4.45 × 10⁻¹⁵⁶	0
	SD	8.87 × 10⁻⁴	9.55 × 10⁻⁸¹	1.22 × 10⁻¹¹	1.52 × 10³	9.29 × 10⁻¹⁵⁶	0
f₂	MBF	7.20 × 10⁻¹	3.60 × 10⁻⁶	3.96 × 10⁻¹	8.47 × 10¹	0	0
	SD	4.16 × 10⁻²	1.13 × 10⁻⁵	3.07 × 10⁻¹	3.06 × 10¹	0	0
f₃	MBF	6.17	1.81 × 10⁻⁵	2.72 × 10⁻⁶	8.49 × 10⁴	0	0
	SD	7.70 × 10⁻¹	4.69 × 10⁻⁵	9.49 × 10⁻⁷	3.51 × 10⁴	0	0
f₄	MBF	6.71 × 10⁻²	9.27 × 10⁻²	3.53	5.34 × 10¹	0	0
	SD	1.05 × 10⁻²	2.33 × 10⁻²	2.67 × 10⁻¹	1.93 × 10¹	0	0
f₅	MBF	1.07 × 10²	9.69 × 10⁻¹²	1.32 × 10²	6.07 × 10⁵	4.18 × 10⁻²	1.03 × 10¹
	SD	2.39 × 10¹	3.06 × 10⁻¹¹	4.49 × 10¹	8.23 × 10⁵	2.86 × 10⁻²	2.55 × 10¹
f₆	MBF	0	0	0	2.35 × 10³	0	0
	SD	0	0	0	2.90 × 10³	0	0
f₇	MBF	2.87 × 10⁻³	1.83 × 10⁻³	1.80 × 10⁻²	8.57 × 10⁻¹	1.40 × 10⁻⁴	6.28 × 10⁻⁴
	SD	5.74 × 10⁻⁴	2.02 × 10⁻³	2.17 × 10⁻³	6.27 × 10⁻¹	1.45 × 10⁻⁴	1.60 × 10⁻³
f₈	MBF	−4.19 × 10⁴	−4.17 × 10⁴	−3.99 × 10⁴	−4.20 × 10⁴	−4.19 × 10⁴	−3.34 × 10⁴
	SD	4.57 × 10⁻²	4.35 × 10⁻²	9.76 × 10⁻¹	4.57 × 10⁻²	1.33 × 10⁻²	1.28 × 10⁴
f₉	MBF	1.16 × 10⁻²	4.52 × 10⁻⁸	4.41 × 10⁻⁸	4.64 × 10¹	0	0
	SD	1.14 × 10⁻³	1.41 × 10⁻⁷	3.47 × 10⁻⁹	1.63 × 10¹	0	0
f₁₀	MBF	1.44 × 10⁻²	1.72 × 10⁻⁴	5.91 × 10⁻⁶	7.03	3.55 × 10⁻¹⁵	5.88 × 10⁻¹⁶
	SD	6.71 × 10⁻⁴	5.44 × 10⁻⁴	1.63 × 10⁻⁷	2.47	0	0
f₁₁	MBF	4.50 × 10⁻³	7.40 × 10⁻⁴	2.28 × 10⁻³	2.14 × 10¹	0	0
	SD	5.18 × 10⁻⁴	2.28 × 10⁻³	4.27 × 10⁻³	3.00 × 10¹	0	0
f₁₂	MBF	6.83 × 10⁻⁶	3.49 × 10⁻²⁵	5.26 × 10⁻¹³	2.93 × 10⁻¹	1.50 × 10⁻⁷	1.90 × 10⁻³¹
	SD	6.46 × 10⁻⁷	8.66 × 10⁻²⁵	4.25 × 10⁻¹⁴	9.27 × 10⁻¹	3.70 × 10⁻⁸	5.09 × 10⁻³⁰
f₁₃	MBF	3.33 × 10⁻⁴	1.07 × 10⁻¹⁵	1.86 × 10⁻⁶	2.04 × 10⁵	1.75 × 10⁻⁵	2.80 × 10⁻³⁰
	SD	4.66 × 10⁻⁵	5.60 × 10⁻⁶	5.89 × 10⁻⁶	4.44 × 10⁵	5.60 × 10⁻⁶	3.77 × 10⁻²⁹
f₁₄	MBF	1.94 × 10¹	1.23 × 10¹	4.06 × 10¹	3.17 × 10¹	0	0
	SD	1.20	7.70	7.70	5.15	0	0
f₁₅	MBF	1.33 × 10⁻¹	4.29 × 10²	3.15 × 10¹	2.18 × 10⁵	0	0
	SD	1.61 × 10⁻²	8.19 × 10²	6.83 × 10¹	2.42 × 10⁴	0	0

下载: 导出CSV

表 11 6种算法的Friedman检验结果

Table 11 Friedman test results of 6 algorithms

算法	平均值 (排名)
SWTPSO	1.83 (1)
IGPSO	1.97 (2)
ODE	3.13 (3)
IGHS	3.93 (4)
GDHS	4.20 (5)
GABC	5.93 (6)

下载: 导出CSV

表 12 CEC2015函数集

Table 12 CEC2015 test suite

函数序号	函数类型	描述	最优值
1	单峰函数	Bent Cigar 旋转函数	100
2	单峰函数	Discus 旋转函数	200
3	简单多峰函数	Schwefel 偏移旋转函数	300
4		Schwefel偏移旋转函数	400
5		Katsuura 偏移旋转函数	500
6		HappyCat 偏移旋转函数	600
7		HGBat 偏移旋转函数	700
8		Griewank + Rosenbrocl 扩展偏移旋转函数	800
9		Schffer＇s F6 扩展偏移旋转函数	900
10	混合函数	混合函数 1 (N = 3)	1000
11		混合函数 2 (N = 4)	1100
12		混合函数 3 (N = 5)	1200
13	组合函数	组合函数 1 (N = 5)	1300
14		组合函数 2 (N = 3)	1400
15		组合函数 3 (N = 5)	1500

下载: 导出CSV

表 13 CEC2015实验结果对比

Table 13 Test results comparison of CEC2015 test suite

函数符号		PSO	FA	FFPSO	HPSOFF	HFPSO	SWTPSO
1	MBF	3.9049 × 10⁹	2.8899 × 10¹⁰	9.2383 × 10¹⁰	4.7539 × 10⁹	1.1795×10⁹	3.0501 × 10⁸
2	MBF	9.9760 × 10⁴	1.3418 × 10⁵	6.9430 × 10⁶	9.7376 × 10⁴	8.5653 × 10⁴	5.6540 × 10⁴
3	MBF	3.3113 × 10²	3.3850 × 10²	3.4771 × 10²	3.3059 × 10²	3.2638 × 10²	3.2631 × 10²
4	MBF	7.7928 × 10³	8.0330 × 10³	9.6696 × 10³	6.8199 × 10²	5.1202×10³	5.8445×10³
5	MBF	5.0418 × 10²	5.0425 × 10²	5.0586 × 10²	5.0422 × 10²	5.0410 × 10²	5.0403 × 10²
6	MBF	6.0096 × 10²	6.0410 × 10²	6.0755 × 10²	6.0090 × 10²	6.0076 × 10²	6.0053 × 10²
7	MBF	7.0405 × 10²	7.6997 × 10²	8.9401 × 10²	7.0750 × 10²	7.0074 × 10²	7.0046 × 10²
8	MBF	4.0013 × 10³	2.3491 × 10⁶	1.6032 × 10⁸	1.7191 × 10⁴	2.6354 × 10³	1.5070 × 10³
9	MBF	9.1358 × 10²	9.1381 × 10²	9.1427 × 10²	9.1365 × 10²	9.1337 × 10²	9.1328 × 10²
10	MBF	7.5602 × 10⁶	2.9938 × 10⁷	3.9391 × 10⁸	1.1337 × 10⁷	5.4690 × 10⁶	1.1662 × 10⁷
11	MBF	1.1614 × 10³	1.2889 × 10³	2.1094 × 10³	1.1551 × 10³	1.1336 × 10³	1.1267 × 10³
12	MBF	2.2417 × 10³	2.9655 × 10³	4.9876 × 10⁵	2.0617 × 10³	1.7752 × 10³	1.7819 × 10³
13	MBF	1.7719 × 10³	1.9596 × 10³	3.6329 × 10³	1.7390 × 10³	1.6866×10³	1.6552 × 10³
14	MBF	1.6644 × 10³	1.7522 × 10³	2.1683 × 10³	1.6711 × 10³	1.6469 × 10³	1.6510 × 10³
15	MBF	2.5522 × 10³	2.8256 × 10³	3.9013 × 10³	2.6529 × 10³	2.4467 × 10³	2.4639 × 10³

下载: 导出CSV

表 14 Friedman检验结果

Table 14 Results of Friedman test

算法	秩均值 (排名)
SWTPSO	1.53 (1)
HFPSO	1.73 (2)
HPSOFF	3.33 (3)
PSO	3.40 (4)
FA	5.00 (5)
FFPSO	6.00 (6)

下载: 导出CSV

参考文献(51)

[1]	Kennedy J, Eberhart R C. Particle swarm optimization. In: Proceedings of the 1995 International Conference on Neural Networks. Perth, Australia: 1995. 1942−1948
[2]	贾承丰, 韩华, 吕亚楠, 张路.基于word2vec和粒子群的链路预测算法[J/OL].自动化学报, 已录用 Jia Cheng-Feng, Han Hua, Lv Ya-Nan, Zhang Lu. Link prediction algorithm based on word2vec and particle swarm. Acta Automatica Sinica, to be published
[3]	Chou J S, Pham A D Nature-inspired metaheuristic optimization in least squares support vector regression for obtaining bridge scour information. Information Sciences, 2017, 399: 64-80 doi: 10.1016/j.ins.2017.02.051
[4]	Tran B, Xue B, Zhang M Variable-Length Particle Swarm Optimization for Feature Selection on High-Dimensional Classification. IEEE Transactions on Evolutionary Computation, 2018, 23(3): 473-487
[5]	Mousavi-Avval S H, Rafiee S, Sharifi M, Hosseinpour S, Notarnicola B, Tassielli G, et al. Application of multi-objective genetic algorithms for optimization of energy, economics and environmental life cycle assessment in oilseed production. Journal of Cleaner Production, 2017, 140: 804-815 doi: 10.1016/j.jclepro.2016.03.075
[6]	Chen S M, Huang Z C. Multiattribute decision making based on interval-valued intuitionistic fuzzy values and particle swarm optimization techniques. Information Sciences, 2017, 397: 206-218
[7]	吴庆, 赵涛, 佃松宜, 郭锐, 李胜川, 方红帏, 等.基于FPSO的电力巡检机器人的广义二型模糊逻辑控制. 自动化学报, 已录用 Wu Qing, Zhao Tao, Dian Song-Yi, Guo Rui, Li Sheng-Chuang, Fang Hong-Wei, et al. General type-2 fuzzy logic control for a power-line inspection robot based on FPSO. Acta Automatica Sinica, to be published
[8]	Shi Y, Eberhart R. A modified particle swarm optimizer. In: Proceedings of the 1998 IEEE International Conference on Evolutionary Computation. Anchorage, USA: IEEE, 1998. 69−73
[9]	Zhan Z H, Zhang J, Li Y, Chung H S H. Adaptive particle swarm optimization. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics), 2009, 39(6): 1362- 1381 doi: 10.1109/TSMCB.2009.2015956
[10]	Zhang L, Tang Y, Hua C, Guan X. A new particle swarm optimization algorithm with adaptive inertia weight based on Bayesian techniques. Applied Soft Computing, 2015, 28: 138-149 doi: 10.1016/j.asoc.2014.11.018
[11]	Ratnaweera A, Halgamuge S K, Watson H C. Self-organizing hierarchical particle swarm optimizer with time-varying acceleration coefficients. IEEE Transactions on Evolutionary Computation, 2004, 8(3): 240-255 doi: 10.1109/TEVC.2004.826071
[12]	Chen K, Zhou F, Yin L, Wang S, Wang Y, Wan F. A hybrid particle swarm optimizer with sine cosine acceleration coefficients. Information Sciences, 2018, 422: 218-241 doi: 10.1016/j.ins.2017.09.015
[13]	Li S F, Cheng C Y. Particle Swarm Optimization with Fitness Adjustment Parameters. Computers & Industrial Engineering, 2017, 113:831-841
[14]	Isiet M, Gadala M. Self-adapting control parameters in particle swarm optimization. Applied Soft Computing, 2019, 83: 105653 doi: 10.1016/j.asoc.2019.105653
[15]	Kennedy J. Bare bones particle swarms. In: Proceedings of the 2003 IEEE Swarm Intelligence Symposium. Indianapolis, USA: 2003. 80−87
[16]	Wang H, Wu Z, Rahnamayan S, Liu Y, Ventresca M. Enhancing particle swarm optimization using generalized opposition-based learning. Information Sciences, 2011, 181(20): 4699-4714 doi: 10.1016/j.ins.2011.03.016
[17]	Yan B, Zhao Z, Zhou Y, Yuan W, Li J, Wu J, et al. A particle swarm optimization algorithm with random learning mechanism and Levy flight for optimization of atomic clusters. Computer Physics Communications, 2017, 219:79-86 doi: 10.1016/j.cpc.2017.05.009
[18]	Xia X, Xing Y, Wei B, Zhang Y, Li X, Deng X, et al. A fitness-based multi-role particle swarm optimization. Swarm and Evolutionary Computation, 2019, 44: 349-364 doi: 10.1016/j.swevo.2018.04.006
[19]	Liu H R, Cui J C, Lu Z D, Liu D Y, Deng Y J. A hierarchical simple particle swarm optimization with mean dimensional information. Applied Soft Computing, 2019, 76: 712-725 doi: 10.1016/j.asoc.2019.01.004
[20]	Kennedy J, Mendes R. Neighborhood topologies in fully informed and best-of-neighborhood particle swarms. IEEE Transactions on Systems, Man, and Cybernetics, Part C (Applications and Reviews), 2006, 36(4): 515-519 doi: 10.1109/TSMCC.2006.875410
[21]	Mendes R, Kennedy J, Neves J. The fully informed particle swarm: simpler, maybe better. IEEE Transactions on Evolutionary Computation, 2004, 8(3): 204-210 doi: 10.1109/TEVC.2004.826074
[22]	Lim W H, Isa N A M, Particle swarm optimization with adaptive time-varying topology connectivity. Applied Soft Computing, 2014, 24: 623-642 doi: 10.1016/j.asoc.2014.08.013
[23]	Lin A, Sun W, Yu H, Wu G, Tang H. Global genetic learning particle swarm optimization with diversity enhancement by ring topology. Swarm and Evolutionary Computation, 2019, 44: 571-583 doi: 10.1016/j.swevo.2018.07.002
[24]	Zhang K, Huang Q, Zhang Y. Enhancing comprehensive learning particle swarm optimization with local optima topology. Information Sciences, 2019, 471: 1-18 doi: 10.1016/j.ins.2018.08.049
[25]	Alatas B, Akin E, Ozer A B. Chaos embedded particle swarm optimization algorithms. Chaos, Solitons & Fractals, 2009, 40(4): 1715-1734
[26]	Kao Y T, Zahara E. A hybrid genetic algorithm and particle swarm optimization for multimodal functions. Applied Soft Computing, 2008, 8(2): 849-857 doi: 10.1016/j.asoc.2007.07.002
[27]	Aydilek İ B. A hybrid firefly and particle swarm optimization algorithm for computationally expensive numerical problems. Applied Soft Computing, 2018, 66: 232-249 doi: 10.1016/j.asoc.2018.02.025
[28]	Lynn N, Suganthan P N. Ensemble particle swarm optimizer. Applied Soft Computing, 2017, 55: 533-548 doi: 10.1016/j.asoc.2017.02.007
[29]	Wang F, Zhang H, Li K, et al. A hybrid particle swarm optimization algorithm using adaptive learning strategy. Information Sciences, 2018, 436: 162-177
[30]	Cui L, Li G, Lin Q, Du Z, Gao W, Chen J, et al. A novel artificial bee colony algorithm with depth-first search framework and elite-guided search equation. Information Sciences, 2016, 367: 1012-1044
[31]	Liang J J, Qin A K, Suganthan P N, Baskar S. Comprehensive learning particle swarm optimizer for global optimization of multimodal functions. IEEE Transactions on Evolutionary Computation, 2006, 10(3): 281-295 doi: 10.1109/TEVC.2005.857610
[32]	Haklı H, Uğuz H. A novel particle swarm optimization algorithm with Levy flight. Applied Soft Computing, 2014, 23: 333-345 doi: 10.1016/j.asoc.2014.06.034
[33]	Liang J J, Suganthan P N. Dynamic multi-swarm particle swarm optimizer. In: Proceedings of the 2005 IEEE Swarm Intelligence Symposium. Pasadena, USA: 2005. 124−129
[34]	Jensi R, Jiji G W. An enhanced particle swarm optimization with levy flight for global optimization. Applied Soft Computing, 2016, 43: 248-261 doi: 10.1016/j.asoc.2016.02.018
[35]	Gao W, Chan F T S, Huang L, Liu S. Bare bones artificial bee colony algorithm with parameter adaptation and fitness-based neighborhood. Information Sciences, 2015, 316: 180-200 doi: 10.1016/j.ins.2015.04.006
[36]	Jadon S S, Bansal J C, Tiwari R, Sharma H. Accelerating artificial bee colony algorithm with adaptive local search. Memetic Computing, 2015, 7(3): 215-230 doi: 10.1007/s12293-015-0158-x
[37]	Kiran M S, Hakli H, Gunduz M, Uguz H. Artificial bee colony algorithm with variable search strategy for continuous optimization. Information Sciences, 2015, 300: 140-157 doi: 10.1016/j.ins.2014.12.043
[38]	Cui L, Zhang K, Li G, Fu X, Wen Z, Lu N, et al. Modified Gbest-guided artificial bee colony algorithm with new probability model. Soft Computing, 2018, 22(7): 2217-2243 doi: 10.1007/s00500-017-2485-y
[39]	Wang Y, Cai Z, Zhang Q. Differential evolution with composite trial vector generation strategies and control parameters. IEEE Transactions on Evolutionary Computation, 2011, 15(1): 55-66 doi: 10.1109/TEVC.2010.2087271
[40]	Zhang J, Sanderson A C. JADE: adaptive differential evolution with optional external archive. IEEE Transactions on Evolutionary Computation, 2009, 13(5): 945-958 doi: 10.1109/TEVC.2009.2014613
[41]	Brest J, Maucec, M S. Self-adaptive differential evolution algorithm using population size reduction and three strategies. Soft Computing, 2011, 15(11): 2157-2174 doi: 10.1007/s00500-010-0644-5
[42]	Qin A K, Huang V L, Suganthan P N. Differential evolution algorithm with strategy adaptation for global numerical optimization. IEEE Transactions on Evolutionary Computation, 2009, 13(2): 398-417 doi: 10.1109/TEVC.2008.927706
[43]	Hansen N, Ostermeier A. Completely derandomized self-adaptation in evolution strategies. Evolutionary Computation, 2001, 9(2): 159-195 doi: 10.1162/106365601750190398
[44]	Khalili M, Kharrat R, Salahshoor K, Sefat M H. Global dynamic harmony search algorithm: GDHS. Applied Mathematics and Computation, 2014, 228: 195-219 doi: 10.1016/j.amc.2013.11.058
[45]	El-Abd M. An improved global-best harmony search algorithm. Applied Mathematics and Computation, 2013, 222: 94-106 doi: 10.1016/j.amc.2013.07.020
[46]	Rahnamayan S, Tizhoosh H R, Salama M M A. Opposition-based differential evolution. IEEE Transactions on Evolutionary computation, 2008, 12(1): 64-79 doi: 10.1109/TEVC.2007.894200
[47]	Zhu G, Kwong S. Gbest-guided artificial bee colony algorithm for numerical function optimization. Applied Mathematics and Computation, 2010, 217(7): 3166-3173 doi: 10.1016/j.amc.2010.08.049
[48]	Ouyang H, Gao L, Li S, Kong X Y. Improved global-best-guided particle swarm optimization with learning operation for global optimization problems. Applied Soft Computing, 2017, 52: 987-1008 doi: 10.1016/j.asoc.2016.09.030
[49]	Liang J J, Qu B Y, Suganthan P N, Chen Q. Problem Definitions and Evaluation criteria for the CEC 2015 Competition on Learning-based Real-parameter Single Objective Optimization. Technical Report 201411A, Computational Intelligence Laboratory, Zhengzhou University, China, Nanyang Technological University, Singapore, 2014. 29: 625−640
[50]	Kora P, Krishna K S R. Hybrid firefly and particle swarm optimization algorithm for the detection of bundle branch block. International Journal of the Cardiovascular Academy, 2016, 2(1): 44-48 doi: 10.1016/j.ijcac.2015.12.001
[51]	Arunachalam S, AgnesBhomila T, Babu M R. Hybrid particle swarm optimization algorithm and firefly algorithm based combined economic and emission dispatch including valve point effect. In: Proceedings of the 2014 International Conference on Swarm, Evolutionary, and Memetic Computing. Odisha, India: 2014. 647−660

施引文献

期刊类型引用(3)

1.	周玉，裴泽宣，王培崇，陈博. 引入相量算子和流向算子的天鹰优化算法. 浙江大学学报(工学版). 2024(02): 304-316 . 百度学术
2.	李大海，曾能智，王振东. 基于相对距离和历史成功率机制的增强麻雀搜索算法. 计算机应用研究. 2024(06): 1640-1648 . 百度学术
3.	赵嘉，陈文平，肖人彬，王晖. 面向多峰优化问题的自主学习萤火虫算法. 控制与决策. 2022(08): 1971-1980 . 百度学术

其他类型引用(6)

资源附件(0)

访问统计

图(5) / 表(14)

计量

文章访问数: 540
HTML全文浏览量: 77
PDF下载量: 161
被引次数: 9

1. 电网模型、攻击特性及问题描述
1.1 电网模型
1.2 攻击描述和问题提出
2. 自适应平方根无迹卡尔曼滤波
2.1 平方根无迹变换
2.2 噪声估计方法
2.3 自适应平方根无迹卡尔曼滤波设计
2.4 自适应平方根无迹卡尔曼滤波稳定性
3. 检测方法
3.1 攻击检测
4. 仿真
4.1 ASRUKF滤波算法的有效性
4.2 虚假数据注入攻击下的三种检测方法
5. 结束语

1. 电网模型、攻击特性及问题描述
1.1 电网模型
1.2 攻击描述和问题提出
2. 自适应平方根无迹卡尔曼滤波
2.1 平方根无迹变换
2.2 噪声估计方法
2.3 自适应平方根无迹卡尔曼滤波设计
2.4 自适应平方根无迹卡尔曼滤波稳定性
3. 检测方法
3.1 攻击检测
4. 仿真
4.1 ASRUKF滤波算法的有效性
4.2 虚假数据注入攻击下的三种检测方法
5. 结束语

参考文献(51)

施引文献

资源附件(0)

访问统计

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

留言板

方波触发勘探与开发的粒子群优化算法

doi: 10.16383/j.aas.c190842

计量

Particle Swarm Optimization With Square Wave Triggered Exploration and Exploitation

1. 电网模型、攻击特性及问题描述

1.1 电网模型

1.2 攻击描述和问题提出

2. 自适应平方根无迹卡尔曼滤波

2.1 平方根无迹变换

2.2 噪声估计方法

2.3 自适应平方根无迹卡尔曼滤波设计

2.4 自适应平方根无迹卡尔曼滤波稳定性

3. 检测方法

3.1 攻击检测

4. 仿真

4.1 ASRUKF滤波算法的有效性

4.2 虚假数据注入攻击下的三种检测方法

5. 结束语

期刊类型引用(3)

其他类型引用(6)

计量

目录

1. 电网模型、攻击特性及问题描述

1.1 电网模型

1.2 攻击描述和问题提出

2. 自适应平方根无迹卡尔曼滤波

2.1 平方根无迹变换

2.2 噪声估计方法

2.3 自适应平方根无迹卡尔曼滤波设计

2.4 自适应平方根无迹卡尔曼滤波稳定性

3. 检测方法

3.1 攻击检测

4. 仿真

4.1 ASRUKF滤波算法的有效性

4.2 虚假数据注入攻击下的三种检测方法

5. 结束语

留言板

方波触发勘探与开发的粒子群优化算法

doi: 10.16383/j.aas.c190842

计量

出版历程

Particle Swarm Optimization With Square Wave Triggered Exploration and Exploitation

1. 电网模型、攻击特性及问题描述

1.1 电网模型

1.2 攻击描述和问题提出

2. 自适应平方根无迹卡尔曼滤波

2.1 平方根无迹变换

2.2 噪声估计方法

2.3 自适应平方根无迹卡尔曼滤波设计

2.4 自适应平方根无迹卡尔曼滤波稳定性

3. 检测方法

3.1 攻击检测

4. 仿真

4.1 ASRUKF滤波算法的有效性

4.2 虚假数据注入攻击下的三种检测方法

5. 结束语

期刊类型引用(3)

其他类型引用(6)

计量

出版历程

目录

1. 电网模型、攻击特性及问题描述

1.1 电网模型

1.2 攻击描述和问题提出

2. 自适应平方根无迹卡尔曼滤波

2.1 平方根无迹变换

2.2 噪声估计方法

2.3 自适应平方根无迹卡尔曼滤波设计

2.4 自适应平方根无迹卡尔曼滤波稳定性

3. 检测方法

3.1 攻击检测

4. 仿真

4.1 ASRUKF滤波算法的有效性

4.2 虚假数据注入攻击下的三种检测方法

5. 结束语