本发明涉及电力系统分析技术领域，特别是涉及一种基于电力系统多维特征的功率振荡类型判别方法。
背景技术：
随着我国电力系统规模的不断扩大，低频振荡发生的风险日益增大，而且表现出许多新的特征。电力系统主要存在两种类型的低频振荡，一种是由于系统阻尼不足而引发的负阻尼振荡，另一种是由于周期性功率扰动引起的强迫功率振荡。在实际电力系统中，负阻尼振荡和强迫功率振荡由于产生机理不同需要采取不同的抑制措施。但是该两种振荡波形相似，有时难以区分其振荡类型，所以研究能够快速有效识别振荡类型的方法具有重要意义。
目前提出的方法包括基于波形的判别方法、基于能量的判别方法等，其主要判据来源均为基于数学模型进行理论推导，计算时域或频域下的某一种指标，认为该指标为区别不同机理的振荡的本质特征。但随着研究的深入，多个不同的本质特征被提出，如起振阶段的频率响应分量数量、包络线的变化趋势、端口能量变化等。在目前复杂大电网的情况下，这些提取的特征是否为不同振荡机理的本质特征、单一某种指标是否足以判别(即判据的充分性和必要性)，有待研究。实际上，已有文献对某些判据的充要性提出了质疑，并给出了反例，证明了某些判据为非充分必要条件。例如强迫振荡为拍频振荡时，其起振波形与负阻尼振荡的起振波形特征类似，此时基于起振段波形的判别方法很有可能误判。此外，振荡发生后系统中收集到大量的数据，仅就局部的某一特征给出判别结果，如波形包络线形状、频率响应分量数量等，大量的其他有效信息被忽略。而且，随着电网规模越来越大，电网的特性也越来越复杂，仅靠人工经验难以全面把握电网的安全特征及规律，容易造成信息遗漏，且难以发现电网中潜在的耦合关系，特征选择方法对特征间的协同效应考虑不足，可靠性差，准确率低。因此亟需研究一种新的低频振荡类型识别方法。
技术实现要素：
发明目的：本发明的目的是提供一种能够解决现有技术中存在的缺陷的基于电力系统多维特征的功率振荡类型判别方法。
技术方案：为达到此目的，本发明采用以下技术方案：
本发明所述的基于电力系统多维特征的功率振荡类型判别方法，包括以下步骤：
s1：建立电力系统仿真模型，通过调节发电机励磁、电力系统负荷或施加短路故障使得电力系统呈现弱阻尼特性来进行负阻尼振荡的批量仿真，通过在原动机转矩、励磁或者负荷上施加扰动源来进行强迫功率振荡的批量仿真，从而获取批量数据样本；所述扰动源为周期性正弦波或者方波；
s2：对振荡的数据样本计算低频振荡信号时域、频域、能量、相关性、复杂度和模态共六方面的特征指标集；所述复杂度为样本熵；
s3：使用互信息特征选择方法对数据样本的特征指标集进行特征选择，得到经过特征选择后的指标集；
s4：将特征选择后的指标集使用机器学习分类器进行监督学习，得到功率振荡事件类型的识别模型；
s5：使用批量数据样本对功率振荡事件类型的识别模型进行交叉验证；
s6：对电力系统采集得到的pmu信号计算特征指标集，将特征指标集输入到功率振荡事件类型的识别模型，从而判别实际系统发生的振荡类型。
进一步，所述步骤s1中的批量数据样本包括发电机输出有功功率信号、发电机输出无功功率信号、发电机转子角速度信号和发电机端电压信号。
进一步，所述步骤s2中的时域方面的特征指标集包括发电机有功功率信号的均值、样本标准差、方根幅值、均方根值、峰值、歪度指标、峭度指标、峰值系数、裕度指标、波形指标和脉冲指标；频域方面的特征指标集包括中心频率、方差、偏斜度、峰度、频度中心、频率标准差、均方根频率、波形稳定系数、变异系数、歪度、峭度和均方根比率。
进一步，所述步骤s2中的能量方面的特征指标集包括低频振荡能量函数的时域指标、频域指标和能量时空分布熵；其中，低频振荡能量函数通过式(1)得到，能量时空分布熵通过式(2)得到：
egi＝∫δpgi2πδfidt+∫δqgid(δlnui)(1)
式(1)中，egi为第i台发电机的低频振荡能量函数；δpgi为第i台发电机输出的有功功率相对稳态值的变化量；δfi为第i台发电机的频率偏移量；δqgi为第i台发电机输出的无功功率相对稳态值的变化量；δlnui为第i台发电机母线电压的自然对数值的相对稳态值的变化量；
式(2)中，eσ为系统振荡能量之和；n为发电机总数；soe为能量时空分布熵。
进一步，所述步骤s2中的相关性方面的特征指标集包括互相关函数和自相关函数；其中，互相关函数r12通过式(3)得到，自相关函数r(τ)通过式(4)得到：
式(3)中，f1(t)为一般节点有功功率信号关于时间t的函数，f2(t+τ)为参考节点有功功率信号关于时间t+τ的函数，参考节点是指低频振荡信号电压方差最大的节点，一般节点是指除参考节点以外的其他节点；
式(4)中，xt为有功功率信号关于时间t的函数，xt+τ为有功功率信号关于时间t+τ的函数，μ为有功功率信号的期望，σ为有功功率信号的标准差。
进一步，所述步骤s2中的复杂度方面的特征指标集包括样本熵，样本熵为等时间间隔采样的发电机有功功率的样本熵。
进一步，所述步骤s2中的模态方面的特征指标集包括频率和阻尼比，通过使用总体最小二乘-旋转不变算法对振荡信号进行模态分析得到。
进一步，所述步骤s3中的互信息特征选择方法是基于互信息对特征指标进行评价的特征选择方法；互信息i(x；y)通过式(5)得到，对互信息进行评价通过式(6)所示的互信息评价函数j实现；
式(5)中，p(x)为随机变量x的边际分布，p(y)为随机变量y的边际分布，p(x,y)为随机变量(x,y)的联合分布，x为特征指标变量，y为分类标签变量；
式(6)中，
进一步，所述步骤s4中的机器学习分类器为svm支持向量机、决策树、线性判别分析和最邻近分类中的任意一种。
进一步，所述步骤s5中的交叉验证为k折交叉验证。
有益效果：本发明公开了一种基于电力系统多维特征的功率振荡类型判别方法，与现有技术相比，具有以下有益效果：
1)本发明通过对低频振荡信号计算时域指标，频域指标、能量指标、互相关指标、自相关指标、样本熵指标和模态指标，建立了较完备的指标集，能够较完整地描述电力系统振荡的特征信息；
2)本发明使用了互信息特征选择方法，相比广泛使用的fisher判别法，互信息特征选择可以度量变量之间的非线性关系。使用互信息特征选择方法得到的特征进行模型训练，有助于提高训练模型的泛化能力和降低训练模型的复杂度,从而有效防止过拟合现象产生；
3)本发明使用了机器学习分类器对电力系统功率振荡事件类型进行识别，相比传统分类方法，可以有效提高分类的精度和训练模型的泛化能力。
附图说明
图1为本发明具体实施方式中方法的流程图；
图2为本发明具体实施方式中典型的负阻尼机理低频振荡有功功率波形图；
图3为本发明具体实施方式中典型的强迫振荡有功功率波形图；
图4为本发明具体实施方式中形成特征指标集的流程图；
图5为本发明具体实施方式中互信息特征选择方法流程图；
图6为本发明具体实施方式中识别振荡类型的整体实现流程图。
具体实施方式
下面结合具体实施方式和附图对本发明的技术方案作进一步的介绍。
本具体实施方式公开了一种基于电力系统多维特征的功率振荡类型判别方法，如图1和图6所示，包括以下步骤：
s1：建立电力系统仿真模型，通过调节发电机励磁、电力系统负荷或施加短路故障使得电力系统呈现弱阻尼特性来进行负阻尼振荡的批量仿真，通过在原动机转矩、励磁或者负荷上施加扰动源来进行强迫功率振荡的批量仿真，从而获取批量数据样本；所述扰动源为周期性正弦波或者方波。步骤s1中的批量数据样本包括发电机输出有功功率信号、发电机输出无功功率信号、发电机转子角速度信号和发电机端电压信号。步骤s1的具体步骤如下：
步骤1.1在matlab/simulink里搭建四机两区模型，设置总额定负荷为2734mw，区域振荡频率为0.64hz，进行仿真，将仿真运行至50s，使仿真达到稳定状态；
步骤1.2在仿真达到稳态的条件下，改变四机两区的负荷从90％-103％的额定负荷变化，通过调节发电机励磁、电力系统负荷或在两区之间的联络线上施加三相短路故障的方法，使得系统阻尼特性为负阻尼，每隔0.5％的负荷变化记录一组数据，仿真时间15s，为获取与强迫振荡波形相似的负阻尼波形，截取1.5s以后的数据段。为模拟电力系统中pmu的工作状态，数据采样频率为25hz，获得负阻尼振荡的数据样本。图2为典型的负阻尼机理低频振荡有功功率波形图；
步骤1.3在仿真达到稳态的条件下，改变四机两区的负荷从90％-110％的额定负荷变化，在原动机转矩、励磁或者负荷上施加周期性正弦波或方波等扰动源，设置仿真时间15s，每隔0.5％的负荷变化记录一组数据。为模拟电力系统中pmu的工作状态，数据采样频率为25hz，获得强迫振荡的数据样本。图3为典型的强迫振荡有功功率波形图。
s2：对振荡的数据样本计算低频振荡信号时域、频域、能量、相关性、复杂度和模态共六方面的特征指标集；时域方面的特征指标集包括发电机有功功率信号的均值、样本标准差、方根幅值、均方根值、峰值、歪度指标、峭度指标、峰值系数、裕度指标、波形指标和脉冲指标；频域方面的特征指标集包括中心频率、方差、偏斜度、峰度、频度中心、频率标准差、均方根频率、波形稳定系数、变异系数、歪度、峭度和均方根比率。复杂度方面的特征指标集包括样本熵，样本熵为等时间间隔采样的发电机有功功率的样本熵。模态方面的特征指标集包括频率和阻尼比，通过使用总体最小二乘-旋转不变算法对振荡信号进行模态分析得到。所述复杂度为样本熵。具体步骤如下：
步骤2.1计算时域的各个统计特征指标，各特征指标如下：
其中，x(n)是n＝1,2,...,n的信号序列，n是数据点的个数。依照上式，计算发电机有功功率的时域指标。
步骤2.2计算频域的各个统计特征指标，各特征指标如下：
其中，s(k)是k＝1,2,...,k时的频谱,k是谱线的数量，fk第k个谱线的频率。依照上式，计算发电机有功功率的频域指标。
步骤2.3计算系统发生振荡时的能量函数，能量函数的具体计算方法如下：
egi＝∫δpgi2πδfidt+∫δqgid(δlnui)(1)
式(1)中，egi为第i台发电机的低频振荡能量函数；δpgi为第i台发电机输出的有功功率相对稳态值的变化量；δfi为第i台发电机的频率偏移量；δqgi为第i台发电机输出的无功功率相对稳态值的变化量；δlnui为第i台发电机母线电压的自然对数值的相对稳态值的变化量；
通过计算得到的能量函数，计算其时域指标、频域指标和能量时空分布熵作为能量指标。其中，能量时空分布熵的计算方法如下：
式(2)中，eσ为系统振荡能量之和；n为发电机总数；soe为能量时空分布熵。
步骤2.4计算互相关指标,互相关函数计算方法如下：
式(3)中，f1(t)为一般节点有功功率信号关于时间t的函数，f2(t+τ)为参考节点有功功率信号关于时间t+τ的函数，参考节点是指低频振荡信号电压方差最大的节点，一般节点是指除参考节点以外的其他节点。依照上式，计算发电机有功功率时间序列的互相关函数，选取互相关函数取值最大时的延时作为互相关指标；
步骤2.5计算自相关指标,自相关函数计算方法如下：
式(4)中，xt为有功功率信号关于时间t的函数，xt+τ为有功功率信号关于时间t+τ的函数，μ为有功功率信号的期望，σ为有功功率信号的标准差。依照上式，计算发电机有功功率时间序列的互相关函数，选取当自相关函数延时不等于0，且函数取值最大时的延时作为自相关指标；
步骤2.6计算振荡信号的样本熵,样本熵计算方法如下：
(1)将等时间间隔采样的发电机有功功率作为一个待处理的时间序列u,定义算法相关参数m和r，重构m维向量xm(1),xm(2),...,xm(n-m+1),其中xm(i)＝[ui(1),ui(2),...,ui(n-m+1)]；
(2)对于1≤i≤n-m+1,统计满足以下条件的个数：bim(r)＝(满足max|ui(a)-uj(a)|≤r的xm(j)的数量)/(n-m)，i≠j),其中ui(a)为xm(i)的第i个元素，uj(a)为xm(j)的第j个元素，记bim(r)对所有i值的平均值为bm(r)；
(3)取k＝m+1,用相同方法计算bk(r),则样本熵为：-ln[bk(r)/bm(r)]。
计算发电机有功功率时间序列的标准差std，r选取为0.2*std,m为2。依照以上方法，计算样本熵指标。
步骤2.7使用总体最小二乘-旋转不变技术(tls-esprit)算法对振荡信号进行模态分析，将频率和阻尼比作为模态指标：tls-esprit基于子空间技术，把待估计信号分解成信号子空间和噪声子空间，通过信号空间估计出信号参数。选取阶数值为10。
步骤2.8对每个样本进行上述特征指标的采集，获得描述电力系统振荡的特征信息，特征指标集的形成流程图如图4所示。
s3：使用互信息特征选择方法对数据样本的特征指标集进行特征选择，得到经过特征选择后的指标集。互信息特征选择方法是基于互信息对特征指标进行评价的特征选择方法，如图5所示；互信息i(x；y)通过式(5)得到，对互信息进行评价通过式(6)所示的互信息评价函数j实现；
式(5)中，p(x)为随机变量x的边际分布，p(y)为随机变量y的边际分布，p(x,y)为随机变量(x,y)的联合分布，x为特征指标变量，y为分类标签变量；
式(6)中，
本具体实施方式中设置特征选择的个数为10，即可选出最能够表征功率振荡类型的10个特征指标。
s4：将特征选择后的指标集使用机器学习分类器进行监督学习，得到功率振荡事件类型的识别模型。机器学习分类器为svm支持向量机、决策树、线性判别分析和最邻近分类中的任意一种。
s5：使用批量数据样本对功率振荡事件类型的识别模型进行交叉验证。交叉验证为k折交叉验证。将初始采样风格为k个子样本，一个单独的子样本被保留作为验证模型的数据，其他k-1个样本用来训练。交叉验证重复k次，每个子样本验证一次。选择参数k的值为10，进行10折交叉验证，来验证训练模型的分类正确率。
s6：对电力系统采集得到的pmu信号计算特征指标集，将特征指标集输入到功率振荡事件类型的识别模型，从而判别实际系统发生的振荡类型。