本发明属于动态电力系统态势感知技术领域，公开了一种融合递推状态估计的动态电力系统异常检测方法，其包括以下步骤：首先根据电力系统动态特性建立简化的系统节点电压的非线性递推模型，合理表征电力系统负荷动态变化对系统节点电压的影响；然后基于递推状态估计滤波算法实现基于非线性递推模型的系统节点电压动态估计，在此基础上，进一步构造系统节点电压残差随机矩阵；最后构造基于特征谱均值‑方差统计分析的系统残差动态性能指标，以有效反映动态电力系统异常对系统节点电压残差矩阵特征值分布的影响，进而根据自适应统计阈值判定，最终实现对动态电力系统的有效状态评估和异常检测。
一种融合递推状态估计的动态电力系统异常检测方法
技术领域
[0001] 本发明属于动态电力系统态势感知技术领域，尤其涉及一种融合递推状态估计的动态电力系统异常检测方法。
背景技术
[0002] 作为一种新兴技术，动态电力系统态势感知可以有效获取、理解和预测导致系统状况变化的关键因素，为复杂电力系统智能运维和决策提供依据，在电力调度、输配电管理
等领域具有广阔的应用前景。它可以促进电网自动化系统的集成，及时发现系统的弱点和
威胁，并显着提高电力系统的智能化水平，从而保证电力系统安全稳定运行。此外，动态电
力系统态势感知还可以与有效的系统决策和支持模块相结合，从而提高电力系统的可靠
性。因此，动态电力系统态势感知是实现未来智能电网不可或缺的一部分。
[0003] 近年来，随着电力系统规模的不断增大和大数据、人工智能等先进数据分析技术的发展，对动态电力系统态势感知的研究受到广泛关注。电力系统态势感知的时间尺度可
以从瞬态过程的几毫秒、长时状态评估的几分钟到电网规划的几个月甚至几年变化。
[0004] 在现代电力系统中，以实时向量测量单元为基础的监测数据集，其关联的数据管理至关重要。基于PMU数据流的数据驱动方法被认为非常适合于电力系统态势感知应用。其
中，基于随机矩阵理论分析框架，提出了一种数据驱动分析架构来感知高维复杂电力网格
的稳态运行状态，构造了一种基于圆环定理的平均谱半径统计指标进行系统运行异常事件
检测。在此基础上，大量研究集中在基于RMT的电力系统稳态状态评估中，提出了多种基于
线性特征值统计量的改进指标体系来评估电力系统稳定性，并分析了关联统计量的统计相
关性。然而，线性特征值统计量指标的异常检测性能依赖于节点功率变化灵敏度。因此，必
须根据经验设置其异常检测阈值，并且一些初期异常可能无法有效察觉。此外，在动态电力
系统的长时状态监测中，由于电力负荷的正常变化，特别是当智能电网存在分布式可再生
能源时，通常无法保证系统节点电压恒定不变。因此，上述电力系统态势感知性能指标是不
合适的，无法实现对动态电力系统正常和异常运行状态的有效区分。值得注意的是，有效的
系统状态观测器能够很好地估计动态电力系统的运行特征参数，有效地消除原始PMU数据
流中存在的复杂时空关联，从而获得系统正常运行状态下具有显著随机分布特征的系统残
差矩阵。
发明内容
[0005] 本发明的目的在于提出一种融合递推状态估计的动态电力系统异常检测方法，以有效地实现动态电力系统长时状态评估和异常检测。
[0006] 本发明为了实现上述目的，采用如下技术方案：
[0007] 一种融合递推状态估计的动态电力系统异常检测方法，包括如下步骤：
[0008] s1.对复杂电力系统进行多节点电气等效，进而对由电气等效后的多节点构成的复杂电力系统，根据电力系统长时动态特性建立简化的系统节点电压的非线性递推模型；
[0009] s2.基于递推状态估计滤波算法实现基于非线性递推模型的系统节点电压动态递推估计，得到系统节点电压原始残差矩阵；
[0010] 在此基础上，进一步构造随机扰动注入的系统节点电压残差随机矩阵；
[0011] s3.构造基于特征谱均值‑方差统计量的系统残差动态性能指标，进而根据构造的自适应统计判定阈值，实现对动态电力系统的有效状态评估和异常检测。
[0012] 优选地，步骤s1中，经过n节点电气等效的复杂电力系统可表示为：
[0013] xk+1＝f(xk)+qk (1)
[0014] zk+1＝h(xk+1)+rk+1 (2)
[0015] 其中，xk为(2n‑1)×1维电压幅值状态向量，包括n个系统节点电压幅值与相角；
[0016] zk为(2n+2m)×1维系统量测向量，包括n个系统节点电压幅值与相角，以及m条线路有功功率潮流与无功功率潮流；
[0017] f(x)为等效非线性过程方程，用以表征电压幅值的状态向量；
[0018] h(x)为等效非线性过程方程，用以表征系统量测向量的函数关系；
[0019] qk～N(0,Qk)为n×1维系统状态扰动向量；
[0020] rk+1～N(0,Rk+1)为m×1维系统量测扰动向量；
[0021] Qk和Rk+1分别为系统状态扰动向量和系统量测扰动向量对应的协方差矩阵。
[0022] 优选地，非线性递推模型的表达式为：
[0023] xk+1|k＝Sk+bk (3)
[0024] Sk＝αHxk|k+(1‑αH)xk|k‑1 (4)
[0025] bk＝βH(Sk‑Sk‑1)+(1‑βH)bk‑1 (5)
[0026] 其中，xk+1|k表示k时刻对k+1时刻系统状态的预测值；
[0027] xk|k表示k时刻当前系统状态的量测值；
[0028] xk|k‑1表示k‑1时刻对k时刻系统状态的预测值；
[0029] Sk、Sk‑1、bk以及bk‑1均为中间过程变量；αH和βH分别为经验拟合系数。
[0030] 优选地，步骤s2中，递推状态估计滤波算法包括无迹卡尔曼滤波算法或粒子滤波算法。
[0031] 优选地，步骤s2中，系统节点电压原始残差矩阵为系统节点电压幅值与系统节点电压幅值递推估计值的差值矩阵 表示为：
[0032]
[0033] 其中，
[0034] xi,k为系统节点电压幅值矩阵的第i行、第k列元素；
[0035] 为系统节点电压幅值递推估计矩阵的第i行、第k列元素；
[0036] 为系统节点电压原始残差矩阵的第i行、第k列元素；
[0037] 差值矩阵 的列数为w，表示系统节点电压原始残差矩阵的采样窗口大小。
[0038] 优选地，步骤s2中，随机扰动注入的系统节点电压残差随机矩阵Xk为：
[0039]
[0040] 其中，e为满足高斯分布的0均值随机注入扰动，其方差为
[0041] 优选地，步骤s3中基于特征谱均值‑方差统计量的系统残差动态性能指标SAk的表达式为：
[0042]
[0043] 其中，λj,k为k时刻随机扰动注入的系统节点电压残差随机矩阵的协方差矩阵特征值；
[0044] 和 分别为系统正常状态下节点j对应特征值序列λj的均值和标准差，该系统正常状态节点特征值序列由满足高斯分布的0均值随机注入扰动e的均值和标准差近似替
代。
[0045] 优选地，步骤s3中，自适应统计判定阈值δ为：
[0046]
[0047] 其中，
[0048] 其中，f为施加于单一节点电压的异常分量值；
[0049] λn为实非相关Wishart矩阵w‑1EET的最小特征值；
[0050] E为(2n‑1)×w维矩阵，其中每一列分别对应随机注入扰动向量e，其方差为
[0051] 表示异常分量值f的平方的数学期望；‑1 T
[0052] 表示实非相关Wishart矩阵w EE的最小特征值的数学期望；‑1 T
[0053] 为实非相关Wishart矩阵w EE的第i特征值的方差。
[0054] 优选地，自适应统计判定阈值δ为单一异常幅值下异常分量值f检测的充分条件，当系统中多个节点叠加异常分量大于异常分量值f时，自适应统计判定阈值δ均能够可靠检
测。
[0055] 本发明具有如下优点：
[0056] 如上所述，本发明融合递推状态估计理论，在对动态电力系统运行时空关联状态进行有效估计基础上，基于随机矩阵特征值统计分析方法，在设计有效的态势感知性能指
标的同时，给出了明确的异常检测控制限，具有清晰严谨的数学物理意义，能更加有效的量
化动态电力系统潜在异常特性，提高动态电力系统状态评估和异常检测的有效性和合理
性。
附图说明
[0057] 图1为本发明实施例中融合递推状态估计的动态电力系统异常检测方法实施流程图。
[0058] 图2为本发明实施例中系统节点电压原始数据与残差生成数据的正态性检验对比图。
[0059] 图3为本发明实施例中融合递推状态估计的动态电力系统异常检测方法实施效果图。
具体实施方式
[0060] 下面结合附图以及具体实施方式对本发明作进一步详细说明：
[0061] 如图1所示，一种融合递推状态估计的动态电力系统异常检测方法，包括如下步骤：
[0062] I.对复杂电力系统进行多节点电气等效，进而对由电气等效后的多节点构成的复杂电力系统，根据电力系统长时动态特性建立简化的系统节点电压的非线性递推模型。
[0063] 其中，经过n节点电气等效的复杂电力系统可表示为：
[0064] xk+1＝f(xk)+qk (1)
[0065] zk+1＝h(xk+1)+rk+1 (2)
[0066] 其中，xk为(2n‑1)×1维电压幅值状态向量，包括n个系统节点电压幅值与相角，由于系统中有且仅有一个平衡节点，故该平衡节点的相角恒为0。
[0067] zk为(2n+2m)×1维系统量测向量，包括n个系统节点电压幅值与相角，以及m条线路有功功率潮流与无功功率潮流。
[0068] f(x)为等效非线性过程方程，用以表征系统节点电压幅值的状态向量。
[0069] h(x)为等效非线性过程方程，用以表征系统量测向量的函数关系。
[0070] qk～N(0,Qk)为n×1维系统状态扰动向量。
[0071] rk+1～N(0,Rk+1)为m×1维系统量测扰动向量。
[0072] Qk和Rk+1分别为系统状态扰动向量和系统量测扰动向量对应的协方差矩阵。
[0073] 其中，系统节点电压的非线性递推模型表示为：
[0074] xk+1|k＝Sk+bk (3)
[0075] Sk＝αHxk|k+(1‑αH)xk|k‑1 (4)
[0076] bk＝βH(Sk‑Sk‑1)+(1‑βH)bk‑1 (5)
[0077] 其中，xk+1|k表示k时刻对k+1时刻系统状态的预测值；
[0078] xk|k表示k时刻当前系统状态的量测值；
[0079] xk|k‑1表示k‑1时刻对k时刻系统状态的预测值；
[0080] Sk、Sk‑1、bk以及bk‑1均为中间过程变量；αH和βH分别为经验拟合系数。
[0081] 以上非线性递推模型可表征电力系统负荷动态变化对系统节点电压时间序列的影响。
[0082] II.基于递推状态估计滤波算法实现基于非线性递推模型的系统节点电压动态递推估计，得到系统节点电压原始残差矩阵。
[0083] 在此基础上，进一步构造随机扰动注入的系统节点电压残差随机矩阵。
[0084] 其中，I步骤中的非线性递推模型可由递推状态估计滤波算法进行有效的状态估计。递推状态估计滤波算法包括无迹卡尔曼滤波或粒子滤波等非线性状态估计算法。
[0085] 下面以自适应无迹卡尔曼滤波状态估计算法为例，包括三个步骤：
[0086] II.1：系统状态与协方差初始化。
[0087] 对于如式(3)到(5)中的两参数指数平滑递推模型，需要首先确定初始两步长对应状态，因此自适应无迹卡尔曼滤波状态估计算法的实际迭代过程从第三步开始，前两步对
应状态采用实际系统真值，其协方差矩阵初始化为系统量测扰动向量对应的协方差矩阵
Rk+1。
[0088] II.2：状态与协方差预测。
[0089] 令 和Pk分别为k时刻对应的状态向量和协方差矩阵的估计值；
[0090] 为状态向量 的估计均值。
[0091] 根据对称采样原来，其Sigma点集{χi,k}定义为：
[0092]
[0093] 其中，n为系统节点数量，即系统状态数量；
[0094] λ＝α2(n+κ)–n，κ为用于捕获特征分布的高阶矩信息，α为经验相关系数，通常取值为[0.0001,1]； 是 的第i列。
[0095] 由此，进一步可得均值和协方差向量的权值{Wm}和{Wc}分别为：
[0096]
[0097]
[0098] Wim＝Wic＝(n+λ)/2,i＝1,...,2n (9)
[0099] 其中，β＝2对于随机扰动注入下的高斯分布是最优的。
[0100] 进一步地，状态向量均值 及其关联协方差Pk的一步预测估计为：
[0101] χi,k+1|k＝f(χi,k)+qk (10)
[0102]
[0103]
[0104] 其中，χi,k+1|k表示k时刻对k+1时刻系统第i个节点状态预测的Sigma点集。
[0105] f(χi,k)表示k时刻系统第i个节点状态量测值构造Sigma点集的状态函数。
[0106] qk表示为n×1维系统状态扰动向量。
[0107] 表示k时刻对k+1时刻系统状态预测的Sigma点集的加权平均估计值。
[0108] Pk+1|k表示k时刻对k+1时刻系统状态预测的Sigma点集的加权协方差矩阵。
[0109] 其中， 为系统噪声协方差矩阵的估计值，可由次优Sage‑Husa滤波器近似估计得到，通过引入信息修正方差矩阵 表示k时刻系统输出量测值的估计均值。
[0110] 系统噪声的无偏递推协方差估计值表示为：
[0111]
[0112] 其中，Kk表示k时刻的状态更新权值矩阵，由式(19)计算得到；
[0113] dk＝(1‑b)/(1‑bk+1) (14)
[0114] 其中，b∈[0,1]为遗忘因子。
[0115] 系统状态变化越剧烈，b应取值越大，以更加有效地捕捉系统变化动态。
[0116] II.3：状态与协方差更新。
[0117] 类似于式(6)，基于步骤II.2得到的 和Pk+1|k。进一步构造Sigma点集{ξi,k}，有：
[0118]
[0119]
[0120]
[0121]
[0122] Kk+1＝Ck+1/Sk+1 (19)
[0123]
[0124]
[0125] 其中， 和Pk+1分别为k+1时刻的状态向量及其对应的协方差矩阵。
[0126] 表示k时刻对k+1时刻系统第i个节点输出预测的估计值。
[0127] h(ξi,k)表示k时刻系统第i个节点输出量测值构造Sigma点集的输出函数。
[0128] rk表示为n×1维系统输出量测扰动向量。
[0129] 表示k时刻对k+1时刻系统第i个节点输出预测的Sigma点集的加权平均估计值。
[0130] Sk+1表示k+1时刻系统输出预测的Sigma点集的加权协方差矩阵。
[0131] Ck+1表示k+1时刻系统输出预测的Sigma点集与系统状态预测的Sigma点集的加权协方差矩阵。Kk+1表示k+1时刻的状态更新权值矩阵。
[0132] 由此重复进行步骤II.2和步骤II.3即可实现对动态电力系统状态向量的递推估计。
[0133] 系统节点电压原始残差矩阵为系统节点电压幅值与系统节点电压幅值递推估计值的差值矩阵 表示为：
[0134]
[0135] 其中，
[0136] xi,k为系统节点电压幅值矩阵的第i行、第k列元素。
[0137] 为系统节点电压幅值递推估计矩阵的第i行、第k列元素。
[0138] 为系统节点电压原始残差矩阵的第i行、第k列元素。
[0139] 差值矩阵 的列数为w，表示系统节点电压原始残差矩阵的采样窗口大小。
[0140] 随机扰动注入的系统节点电压残差随机矩阵Xk为：
[0141]
[0142] 其中，e为满足高斯分布的0均值随机注入扰动，其方差为
[0143] 图2所示为系统节点电压原始数据与残差随机数据的正态性检验对比图。
[0144] 其中：图2(a)表示在测试周期内系统节点电压状态真值的统计分布，图2(b)表示在测试周期内系统节点电压状态估计残差的统计分布。
[0145] 在图2(a)中，横坐标表示系统节点电压状态真值，纵坐标表示对应于系统节点电压状态真值区间的统计数量。在图2(b)中，横坐标表示系统节点电压状态的估计残差，纵坐
标表示对应于系统节点电压状态估计残差区间的统计数量。
[0146] 由图2中看出，基于非线性递推模型(3)到(5)以及优选的自适应无迹卡尔曼滤波状态估计算法，能够保证动态电力系统运行状态残差分布的随机特性，从而适用于基于残
差随机特性统计分析的系统态势感知性能指标及其异常检测控制限计算。
[0147] III.构造基于特征谱均值‑方差统计量的系统残差动态性能指标SAk，表达式为：
[0148]
[0149] 其中，λj,k为k时刻随机扰动注入的系统节点电压残差随机矩阵的协方差矩阵特征值。
[0150] 和 分别为系统正常状态下节点j对应特征值序列λj的均值和标准差，该系统正常状态节点特征值序列由满足高斯分布的0均值随机注入扰动e的均值和标准差近似替
代。
[0151] 指标SAt用于表征动态电力系统异常对系统节点电压残差矩阵特征值分布的影响。
[0152] 根据构造的自适应统计判定阈值，实现对动态电力系统的有效状态评估和异常检测。
[0153] 其中，自适应统计判定阈值δ的表达式为：
[0154]
[0155] 其中，
[0156] 其中，f为施加于单一节点电压的异常分量值。
[0157] λn为实非相关Wishart矩阵w‑1EET的最小特征值。
[0158] E为(2n‑1)×w维矩阵，其中每一列分别对应随机注入扰动向量e，其方差为
[0159] 表示异常分量值f的平方的数学期望。
[0160] 表示实非相关Wishart矩阵w‑1EET的最小特征值的数学期望。
[0161] 为实非相关Wishart矩阵w‑1EET的第i特征值的方差。
[0162] 自适应统计判定阈值δ为单一异常幅值下异常分量值f检测的充分条件，当系统中多个节点叠加异常分量大于异常分量值f时，自适应统计判定阈值δ均能够可靠检测。
[0163] 本实施例基于IEEE‑30节点测试系统以及节点典型日负荷曲线，对系统残差动态性能指标和自适应统计判定阈值的有效性进行验证，其结果如图3所示。在系统配置中，
1500s时刻对系统节点10施加初期异常信号f，其标准差为随机注入扰动标准差的2倍。由图
3中看出，系统残差动态性能指标和自适应统计判定阈值能够对初期异常信号f进行有效检
测。
[0164] 以上实验验证了本发明方法进行动态电力系统长时状态评估和异常检测的有效性。
[0165] 当然，以上说明仅仅为本发明的较佳实施例，本发明并不限于列举上述实施例，应当说明的是，任何熟悉本领域的技术人员在本说明书的教导下，所做出的所有等同替代、明
显变形形式，均落在本说明书的实质范围之内，理应受到本发明的保护。
