CN106845010B - 基于改进SVD降噪和Prony的低频振荡主导模式辨识方法 - Google Patents

Abstract

本发明公开了一种基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其包括根据输入信号和基本不等式原理，构造出SDV算法中矩阵行数和矩阵列数具有最大乘积的Hankel矩阵；根据输入信号，绘制其信噪比曲线，并对信噪比曲线进行分析，确定最佳有效奇异值阶次；根据最佳有效奇异值阶次对Prony算法中的辨识阶次进行选择，确立最佳辨识阶次；利用具有Hankel矩阵和最佳有效奇异值阶次的SVD算法对输入信号进行处理，得到降噪信号；通过具有最佳辨识阶次的Prony算法对降噪信号进行分析，辨识低频振荡主导模式；该基于改进SVD降噪和Prony的低频振荡主导模式辨识方法具有噪声抑制能力强，辨识精度和准确度高等优点。

Description

基于改进SVD降噪和Prony的低频振荡主导模式辨识方法

技术领域

本发明涉及电力***领域，具体涉及一种基于改进SVD降噪和Prony的低频振荡主导模式辨识方法。

背景技术

随着大电网互联的推进和电力***规模不断扩大，在提高电网运行的可靠性和经济性的同时，也带来了新的安全隐患；近年来多次发生低频振荡严重的危及了电网的安全稳定运行，引起了工业界和学术界的广泛关注；因此，正确分析低频振荡特征参数是有效抑制电力***低频振荡现象的重要基础。

基于受扰轨迹的低频振荡分析可直接对***输出响应进行分析，无需详细的***模型和大规模特征值计算，能够适应***运行方式和参数的变化，并能反映***在扰动之后的动态过程并且计及各种非线性因素的影响；在各种信号分析的方法中，傅里叶变换、小波变换等应用非常广泛，但很难提取出信号的衰减系数，即难以求出阻尼比这一重要特征。

近年来，Prony算法在电力***分析与控制领域得到了广泛的应用，在电力***的研究中，Prony分析方法有广泛的适用性，特别是在小信号稳定控制领域的***辨识中的优势十分明显；利用Prony分析实测数据得到的信息比小扰动分析得到的信息更准确，通过对实际***的Prony分析可直接得到***中主导特征根及其传递函数留数等信息，用一个最优的***降阶模型来逼近原高阶模型。

然而，Prony算法对输入信号的要求较高，噪声干扰会严重影响Prony极点公式的估计的精度，从而使计算的结果出现较大误差；在理想情况下，Prony算法的求解并不复杂，但在白噪声背景下，该复指数模型的最优求解是一个困难的非线性最小二乘问题，而SVD降噪中又存在Hankel矩阵阶数和有效奇异值阶次的难以确定等问题。

发明内容

针对现有技术中的上述不足，本发明提供的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法具有噪声抑制能力强，辨识精度和准确度高等优点。

为了达到上述发明目的，本发明采用的技术方案为：提供一种基于改进SVD 降噪和Prony的低频振荡主导模式辨识方法，其包括根据输入信号和基本不等式原理，构造出SDV算法中矩阵行数和矩阵列数具有最大乘积的Hankel矩阵；根据输入信号，绘制其信噪比曲线，并对信噪比曲线进行分析，确定最佳有效奇异值阶次；根据最佳有效奇异值阶次对Prony算法中的辨识阶次进行选择，确立最佳辨识阶次；利用具有Hankel矩阵和最佳有效奇异值阶次的SVD算法对输入信号进行处理，得到降噪信号；通过具有最佳辨识阶次的Prony算法对降噪信号进行分析，辨识低频振荡主导模式。

进一步地，S1的具体步骤为：设Hankel矩阵行数为m，矩阵列数为n，输入信号为X(N)＝{x₁，x₂，....，x_N}，输入信号中的信号点数为N；根据不等式原理中当m和n相等或最接近时两者的乘积最大，确立矩阵行数m和矩阵列数n 的值，并通过相空间重构构造m×n阶的Hankel矩阵H。

进一步地，Hankel矩阵H的矩阵行数m为：

Hankel矩阵H的矩阵列数n为：

n＝N+1-m：

其中，N为输入信号中的信号点数；

Hankel矩阵H为：

其中，N＝m+n-1；D_m×n为无噪干扰的信号子空间；W_m×n为噪声信号子空间， {x₁，x₂，....，x_N}为输入信号。

进一步地，S2的具体步骤为：根据SDV算法中有效奇异值阶次不同和得到输入信号降噪后的信噪比不同，绘制出信噪比曲线；并选择输入信号中信噪比最大时对应的奇异值阶次，即为最佳有效奇异值阶次。

进一步地，最佳有效奇异值阶次即为最佳辨识阶次的中拟合的最优子集个数。

进一步地，S4的具体步骤为：对Hankel矩阵进行奇异值分解，得到分解后的Hankel矩阵和其矩阵的秩和奇异值；对前K个奇异值进行保存，将剩余的奇异值置零，再利用奇异值分解的逆过程得到重构矩阵，将重构矩阵依据相空间重构的方法进行逆变换，得到降噪信号；其中，K的数值等于最佳有效奇异值阶次的数值。

进一步地，分解后的Hankel矩阵为：

其中，m为矩阵行数，n为矩阵列数，U、V均为正交矩阵，∑为非负对角阵，即：

其中，r为Hankel矩阵H的秩，σ_i为Hankel矩阵H的奇异值。

进一步地，S5的具体步骤为：设低频振荡模式为具有任意振幅、相位、频率和衰减因子的P个指数函数的线性组合，其离散时间的函数形式为：

其中，A_i为幅值，θ_i为初相，f_i为频率，σ_i为衰减因子，p_i为拟合的指数函数的个数，N是采样个数，Δt是采样时间间隔；将

作为实际采样点y(n)的近似，构建代价函数，并令代价函数的值最小，获得离散时间函数；根据离散时间函数和Prony算法的法方程，得到主导模式的振幅、相位、频率和衰减因子。

进一步地，Prony算法的法方程为：

其中，

i，j＝0，1，...，p，x^*(n-i)是x(n-i)的共轭，

p为指数函数的个数，a₁，a₂，...，a_p为待求解系数。

离散时间函数为：

其中，

(n＝0，1，...，N-1)，

e(n)为定义实际测量值y(n)和估计值

的误差，

b₁，b₂，...，b_p为待求解系数。

进一步地，低频振荡主导模式的振幅、相位、频率和衰减因子为：

其中，Re表示取实部，Im表示取虚部，A_i为振幅，θ_i为相位，f_i为频率，σ_i为衰减因子。

本发明的有益效果为：该基于改进SVD降噪和Prony的低频振荡主导模式辨识方法采用基本不等式确定最佳Hankel矩阵阶次并提出利用信噪比来解决奇异值阶次选择问题，利用改进SVD去噪技术对数据进行预处理，提高了信号的信噪比，减小了噪声对Prony分析结果的影响；且进一步通过算例仿真验证了该算法具有噪声抑制能力强、辨识出的主导振荡模式精度高等优点，能够较为准确的辨识电力***低频振荡主导模式。

附图说明

图1给出了在MATLAB环境中振荡信号加噪后的波形图。

图2示意性的给出了对Hankel矩阵进行奇异值分解后，保留奇异值阶次分别为1，2...10的信噪比曲线图。

图3示意性的给出了在MATLAB环境中振荡信号降噪后曲线与原信号曲线图。

图4示意性的给出了Prony38阶的拟合曲线图。

图5示意性的给出了Prony的6阶最优子集拟合曲线图。

图6示意性的给出了6阶拟合平方误差曲线图。

图7示意性的给出了功角曲线的示意图。

图8示意性的给出了对功角曲线进行奇异值分解后，保留奇异值阶次分别为1，2，...，12的信噪比曲线图。

图9示意性的给出了WSCC3机9节点***中的降噪曲线与原信号曲线的对比图。

图10示意性的给出了Prony50阶的拟合曲线图。

图11示意性的给出了Prony的7阶最优子集拟合曲线图。

图12示意性的给出了7阶拟合平方误差曲线图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一种实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明的保护范围。

为简单起见，以下内容省略了该技术领域技术人员所公知的技术常识。

该基于改进SVD降噪和Prony的低频振荡主导模式辨识方法包括：

S1、根据输入信号和基本不等式原理，构造出SDV算法中矩阵行数和矩阵列数具有最大乘积的Hankel矩阵；在具体实施中，对于同一个信号而言，可以重构出不同结构的Hankel矩阵，而不同结构的Hankel矩阵之间会使信号的SVD 分离结果产生很大差别，直接影响信号的降噪效果；为了实现信号和噪声的充分分离，因构造Hankel矩阵的行数和列数的乘积尽可能的最大。

在实际操作中，设Hankel矩阵行数为m，矩阵列数为n，输入信号为 X(N)＝{x₁，x₂，....，x_N}，输入信号中的信号点数为N；根据不等式原理中当m 和n相等或最接近时两者的乘积最大，确立矩阵行数m和矩阵列数n的值，进而保证信号和噪声可得到充分的分离，并通过相空间重构构造m×n阶的Hankel 矩阵H。

满足Hankel矩阵行数m和Hankel矩阵列数n的乘积最大主要取决于信号点数N的奇偶性；在具体实施中，结合信号点数N的奇偶性确定Hankel矩阵行数m和列数n，即Hankel矩阵H的矩阵行数m为：

Hankel矩阵H的矩阵列数n为：

n＝N+1-m；

其中，N为输入信号中的信号点数；在确定出Hankel矩阵行数m和列数n 之后，通过相空间重构构造m×n阶的Hankel矩阵H，构造的Hankel矩阵H为：

其中，N＝m+n-1；D_m×n为无噪干扰的信号子空间；W_m×n为噪声信号子空间， {x₁，x₂，....，x_N}为输入信号。

S2、根据输入信号，绘制其信噪比曲线，并对信噪比曲线进行分析，确定最佳有效奇异值阶次；在具体实施中，根据SDV算法中有效奇异值阶次不同和得到输入信号降噪后的信噪比不同，绘制出信噪比曲线；并选择输入信号中信噪比最大时对应的奇异值阶次，即为最佳有效奇异值阶次；其有效地解决了现有技术中奇异值阈值难以确定的问题；其选出信噪比最高时的信号重构，可达到信噪比最高和降噪效果最明显，进而减小了噪声对Prony分析结果的影响。

S3、根据最佳有效奇异值阶次对Prony算法中的辨识阶次进行选择，确立最佳辨识阶次；其根据最佳有效奇异值阶次进行Prony辨识的阶次确定，有效地解决了Prony辨识阶次选择的困难；在具体实施中，信噪比的定义为：

SNR＝10log₁₀(P_s/P_n)；

其中，P_s为原信号能量，P_n噪声能量；且最佳有效奇异值阶次即为最佳辨识阶次的中拟合的最优子集个数。

S4、利用具有Hankel矩阵和最佳有效奇异值阶次的SVD算法对输入信号进行处理，得到降噪信号；在具体实施中，对Hankel矩阵进行奇异值分解，得到分解后的Hankel矩阵和其矩阵的秩和奇异值；对前K个奇异值进行保存，将剩余的奇异值置零，再利用奇异值分解的逆过程得到重构矩阵，将重构矩阵依据相空间重构的方法进行逆变换，得到降噪信号；其中，K的数值等于最佳有效奇异值阶次的数值，分解后的Hankel矩阵为：

其中，m为矩阵行数，n为矩阵列数，U、V均为正交矩阵，∑为非负对角阵，即：

其中，r为Hankel矩阵H的秩，σ_i为Hankel矩阵H的奇异值。

S5、通过具有最佳辨识阶次的Prony算法对降噪信号进行分析，辨识低频振荡主导模式；在具体实施中，Prony算法是提取平稳振荡模式的常用算法，它针对等间距采样点；设低频振荡模式为具有任意振幅、相位、频率和衰减因子的P个指数函数的线性组合，其离散时间的函数形式为：

其中，A_i为幅值，θ_i为初相，f_i为频率，σ_i为衰减因子，p_i为拟合的指数函数的个数，N是采样个数，Δt是采样时间间隔；将

作为实际采样点y(n)的近似，构建代价函数，并令代价函数的值最小，获得离散时间函数；根据离散时间函数和Prony算法的法方程，得到主导模式的振幅、相位、频率和衰减因子。

在实际操作中，将

作为实际采样点y(n)的近似，其参数辨识的方法是构造代价函数ε，令

为使ε达到最小，从而获得

中的各参数，这需求解非线性方程组，通过现有技术中的一系列的数学变化，可推出差分方程式如下：

为了建立Prony算法，定义实际测量值y(n)和估计值

的误差为e(n)，即

将

带入至差分方程式中，进而得到方程

(n＝o，1，...，N-1)；其中，

因此如果把目标函数修订为使得

最小，则可以找到一组线性方程：

为使目标函数

为最小值，令

则有

其中，x^*(n-i)是x(n-i)的共轭；此时，定义

i，j＝0，1，...，p，即可得到Prony算法的法方程为：

其中，

i，j＝0，1，...，p，x^*(n-i)是x(n-i)的共轭，

p为指数函数的个数，a₁，a₂，...，a_p为待求解系数；可得到系数a₁，a₂，...，a_p，进一步求解特征多项式1+a₁z^-1+...+a_pz^-p＝0得到特征根Zi， i＝1，2，...，p，并启用其简化

(n＝1，2，...，N-1)，可得到离散时间函数，其离散时间函数为：

其中，

(n＝0，1，...，N-1)，

为定义实际测量值y(n)和估计值

的误差，

b1，b2，...，b_p为待求解系数；最后对离散时间函数和Prony算法的法方程进行求解，可算出低频振荡主导模式的振幅、相位、频率和衰减因子，进而确立低频振荡的主导模式；在具体实施中，其低频振荡主导模式的振幅、相位、频率和衰减因子为：

其中，Re表示取实部，Im表示取虚部，A_i为振幅，θ_i为相位，f_i为频率，σ_i为衰减因子。

该基于改进SVD降噪和Prony的低频振荡主导模式辨识方法根据信噪比进行奇异值有效阶次的选择，很好的解决了SVD在重构降噪时关于阈值选取不恰导致降噪效果不明显的问题；根据SVD最佳有效奇异值阶次确定Prony阶次的选择，避免了Prony阶次选择不当，导致辨识结果误差大的问题；且改进的SVD 降噪效果明显，提高了信号的信噪比，减小了噪声对Prony分析结果的影响。

在实际操作中，为了验证本发明所提出的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法的有效性，借助Matlab/Simulink仿真平台对其进行仿真验证。

在具体实施中，在MATLAB环境中产生一组振荡信号的波形，设该组振荡信号的波形为：

如图1所示，图1给出了加噪后的波形图；其中，取采样频率为1kHz，观测窗口长度为5s，向X信号中加入均值为0，方差为1的高斯白噪声，得到加噪后的波形。

如图2所示，图2示意性的给出了对Hankel矩阵进行奇异值分解后，保留奇异值阶次分别为1，2...10的信噪比曲线图；由该信噪比曲线图可知当保留的奇异值个数为6的时候，降噪后的信号信噪比最高，最适合用来进行Prony分析；此时，选取保留最佳有效奇异值阶次为6，并进行信号重构降噪。

如图3所示，图3示意性的给出了降噪后曲线与原信号曲线图，由降噪后曲线与原信号曲线可看出拟合出的曲线具有良好的精度，降噪效果明显，适合用于Prony分析。

如图4所示，图4示意性的给出了Prony38阶的拟合曲线图，其中，拟合曲线与降噪后曲线基本吻合，辨识精度高；当采样数据过密将会影响Prony辨识的精度，因此将降噪后的信号进行间隔30个采样，新序列为0.03s的步长，总时长为5s，将新序列进行Prony 38阶拟合。

如图5和图6所示，图5示意性的给出了Prony的6阶最优子集拟合曲线图，图6示意性的给出了6阶拟合平方误差曲线图；根据SDV算法中保留的最佳有效奇异值阶次为6，在Prony 38阶拟合中选取6阶最优子集进行拟合，其拟合精度高；所以利用信噪比确定奇异值阶次进行重构的同时，也解决了Prony定阶的困难问题，其拟合平方差小于0.01，符合Prony拟合精度要求。

理想值、SVD+Prony和传统Prony在频率和阻尼两方面的对比如下表所示：

由上表可以看出，本发明提出的改进SVD降噪和Prony与传统的Prony在频率和阻尼两方面的对比，在同样为38阶拟合，选取6阶最优子集进行拟合的情况下，传统Prony只能辨识出模式2的频率和阻尼，并不能辨识出模式3的阻尼和频率，对于模式1的频率和阻尼辨识误差达到60％以上，而利用改进的SVD 降噪和Prony算法进行辨识时3个模式都能很好的进行辨识，特别是在频率辨识误差低于2.5％，阻尼辨识误差低于4％，可以看出改进的SVD降噪和Pronyy 算法可以很好的对噪声信号进行辨识。

在实际操作中，选用WSCC3机9节点***作算例在电力***分析综合程序 (PSASP)中进行仿真分析，***发电机总容量为567.5MW，有功负荷为315MW，发电机采用3阶E′q变化模型，励磁***选取PSASP中的1型励磁***，负荷采用恒阻抗负荷。

对功角曲线进行改进SVD和Prony算法分析，将得到的辨识结果与PSASP 中的小干扰稳定分析得到的***主导模式进行对比，验证SVD和Prony算法在处理噪声信号方面的有良好的精度。

考虑到在实际测量信号中含有高频率的噪声，为了更真实的再现实际信号，验证本发明提出的方法对现场实际信号辨识的有效性，对算例中发生突然三相短路后的功角波形加入高次谐波和白噪声，使其得到的曲线更接近于现场实际测量信号。

如图7和图8所示，图7示意性的给出了功角曲线的示意图，图8示意性的给出了对功角曲线进行奇异值分解后，保留奇异值阶次分别为1，2，...，12 的信噪比曲线图；由信噪比曲线可看出当有效奇异值阶次保留7个时，信噪比最高，选取保留最佳有效奇异值阶次为7。

如图9所示，图9示意性的给出了降噪曲线与原信号曲线的对比图，可看出拟合出的曲线具有良好的精度，降噪效果明显，适合用于Prony分析。

如图10和图11所示，图10示意性的给出了Prony50阶的拟合曲线图，图 11示意性的给出了Prony的7阶最优子集拟合曲线图；将降噪后的曲线进行 Prony50阶拟合，拟合曲线与降噪后曲线基本吻合，辨识精度高；根据SVD保留的奇异值个数为7个，则在Prony 50阶拟合中选取7阶最优子集进行拟合，拟合精度高，能够很好的辨识出***中的主导振荡模式。

如图12所示，图12示意性的给出了7阶拟合平方误差曲线图，从中可知拟合平方差小于0.01，符合Prony拟合精度要求。

理想值、SVD+Prony和传统Prony在频率和阻尼两方面的对比如下表所示：

在同样为50阶拟合，选取7阶最优子集进行拟合的情况下，传统Prony在辨识阻尼方面误差达到38％以上，而利用改进的SVD降噪和Prony算法进行辨识时，对***中的2个主导振荡模式都能很好的进行辨识，特别是在频率辨识误差都低于1.7％，阻尼辨识误差低于9％，较传统Prony的辨识精度有较大的提高，所以本发明采用的改进的SVD降噪和Prony算法可以很好的对噪声信号进行辨识。

在具体实施中，上述3机9节点算例仿真验证了该算法具有噪声抑制能力强、辨识出的主导振荡模式精度高等优点，能够较为准确的辨识电力***低频振荡主导模式。

该基于改进SVD降噪和Prony的低频振荡主导模式辨识方法采用基本不等式确定最佳Hankel矩阵阶次并提出利用信噪比来解决奇异值阶次选择问题，利用改进SVD去噪技术对数据进行预处理，提高了信号的信噪比，减小了噪声对 Prony分析结果的影响；且进一步通过算例仿真验证了该算法具有噪声抑制能力强、辨识出的主导振荡模式精度高等优点，能够较为准确的辨识电力***低频振荡主导模式。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将使显而易见的，本文所定义的一般原理可以在不脱离发明的精神或范围的情况下，在其他实施例中实现。因此，本发明将不会被限制与本文所示的这些实施例，而是要符合与本文所公开的原理和新颖性特点相一致的最宽的范围。

Claims (7)

1.一种基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于，包括：

S1、根据输入信号和基本不等式原理，构造出SVD算法中矩阵行数和矩阵列数具有最大乘积的Hankel矩阵；

S2、根据有效奇异值阶次不同，得到输入信号降噪后的信噪比不同，从而绘制出信噪比曲线，进而根据信噪比曲线确定最佳有效奇异值阶次；具体步骤为：选择输入信号中信噪比最大时对应的奇异值阶次，即为最佳有效奇异值阶次；

S3、根据所述最佳有效奇异值阶次对Prony算法中的辨识阶次进行选择，确立最佳辨识阶次，最佳有效奇异值阶次即为最佳辨识阶次中拟合的最优子集个数；

S4、利用具有所述Hankel矩阵和最佳有效奇异值阶次的SVD算法对输入信号进行处理，得到降噪信号；

S5、通过具有最佳辨识阶次的Prony算法对降噪信号进行分析，辨识低频振荡主导模式。

2.根据权利要求1所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于，所述S1的具体步骤为：设Hankel矩阵行数为m，矩阵列数为n，输入信号为X(N)＝{x₁，x₂，….，x_N}，输入信号中的信号点数为N；根据不等式原理中当m和n相等或最接近时两者的乘积最大，确立矩阵行数m和矩阵列数n的值，并通过相空间重构构造m×n阶的Hankel矩阵H。

3.根据权利要求1所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于：所述Hankel矩阵H的矩阵行数m为：

所述Hankel矩阵H的矩阵列数n为：n＝N+1-m；其中，N为输入信号中的信号点数；

所述Hankel矩阵H为：

其中，N＝m+n-1；D_m×n为无噪干扰的信号子空间；W_m×n为噪声信号子空间，{x₁，x₂，….，x_N}为输入信号。

4.根据权利要求1所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于，所述S4的具体步骤为：对所述Hankel矩阵进行奇异值分解，得到分解后的Hankel矩阵和其矩阵的秩和奇异值；对前K个奇异值进行保存，将剩余的奇异值置零，再利用奇异值分解的逆过程得到重构矩阵，将重构矩阵依据相空间重构的方法进行逆变换，得到降噪信号；其中，K为有效奇异值阶次。

5.根据权利要求4所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于：所述分解后的Hankel矩阵为：

其中，m为矩阵行数，n为矩阵列数，U、V均为正交矩阵，Σ为非负对角阵，即：

其中，r为Hankel矩阵H的秩，σ_i为Hankel矩阵H的奇异值。

6.根据权利要求1所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于，所述S5的具体步骤为：设低频振荡模式为具有任意振幅、相位、频率和衰减因子的p个指数函数的线性组合，其离散时间的函数形式为：

其中，A_i为幅值，θ_i为相位，f_i为频率，σ_i为衰减因子，p为拟合的指数函数的个数，N是采样个数，Δt是采样时间间隔；将

作为实际采样点y(n)的近似，构建代价函数，并令代价函数的值最小，获得离散时间函数；根据所述离散时间函数和Prony算法的法方程，得到主导模式的振幅、相位、频率和衰减因子。

7.根据权利要求6所述的基于改进SVD降噪和Prony的低频振荡主导模式辨识方法，其特征在于：所述低频振荡主导模式的振幅、相位、频率和衰减因子为：

其中，Re表示取实部，Im表示取虚部，A_i为振幅，θ_i为相位，f_i为频率，σ_i为衰减因子。

Images