CN105322533A - 基于高斯-马尔科夫模型的自适应t型抗差状态估计方法 - Google Patents

Abstract

本发明涉及一种基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，属于电力系统调度自动化领域。该方法以t型分布对电力系统SCADA量测误差进行建模，利用自由度自适应调整t型估计的计算效率与抗差性。对模型进行简化处理，使得目标函数连续可微，利用与加权最小二乘法类似的牛顿法进行求解，程序兼容性好。该方法可克服部分杠杆量测，仅一次状态估计计算，即可完成不良数据的辨识，抗差性好，实现方便，可大幅提高电网抗差状态估计的运行效率，满足工程应用。

Description

基于高斯-马尔科夫模型的自适应t型抗差状态估计方法

技术领域

本发明涉及一种基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，属于电力系统调度自动化领域。

背景技术

自20世纪70年代状态估计引入调度中心后，已成为能量管理系统的核心与基石。状态估计能够对电力系统运行状态实时监测，从而为调度运行人员以及高级应用软件提供可靠的数据支撑。目前，实际电网采用的加权最小二乘法(weightedleastsquare，WLS)估计模型简单，计算速度快，但难于处理不良数据，从而对状态估计结果产生严重影响。

针对这一问题，研究较多的是抗差状态估计器，其中以M估计研究最多。M估计的抗差性和效率取决于状态量初值的可靠性、等价权函数及其临界值的合理性。电力系统状态估计中使用比较广泛的权函数有Huber权函数、Hampel权函数、Turkey权函数和含有相关量测的IGGIII权函数等。不同等价权函数的选取相当于建立了不同的误差分布模型。统计界近年来提出并有较深入理论研究支撑的t型估计，以t分布对量测误差建模，属于带有刻度参数的M估计。当t分布取较大的自由度时，趋于高斯分布，对应着最小二乘估计，在量测误差为高斯分布时具有较高的效率；当取较小的自由度时，为柯西分布，相应的估计具有很好的抗差性。恰当地选取自由度，可以使t型估计兼顾状态估计的抗差性和效率。目前，尚未见电力系统t型状态估计的报道。

发明内容

发明目的：本发明的目的是提出一种针对电力系统收集的量测数据进行的基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，用以克服不良数据对状态估计结果的影响。

本发明为基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，依次包括以下步骤：

(一)获得电力系统的网络参数和量测数据z；

(二)以量测数据z建立t型抗差状态估计的模型：

$(\hat{x},\widehat{\mathrm{\σ}})=\arg \underset{x,\mathrm{\σ}>0}{\min }\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\left\{f\right[\frac{z_i-h_i\left(x\right)}{\mathrm{\σ}}]+\log \left(\mathrm{\σ}\right)\}---\left(1\right)$

其中，z_i、h_i分别为m×1维量测向量、量测函数向量的第i个分量；w_i为第i个量测对应的权重；σ为一个未知的尺度参数，将其设置为常数1；f(u)＝log(1+u²/V)，V是自由度；x为n×1维的状态向量，包括节点电压幅值与相角；m，n分别为量测量及状态量的个数；

对式(1)以节点注入功率g(x)作为等式约束得到：

$\left\{\begin{array}{ccc}\underset{x}{\min }& J\left(x\right)=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)& \\ s.t.& g\left(x\right)=0& \\ & r_i=z_i-h_i\left(x\right)& i=\mathrm{1,2},...,m\end{array}\right.---\left(2\right)$

其中，r_i为m×1维残差向量的第i个分量；

(三)令l＝0，并设置单次迭代计数器k＝-0；

(四)循环迭代下，利用内点罚函数法对式(2)进行处理：

$\underset{x,\mathrm{\λ}}{\min }L(x,\mathrm{\λ})=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j{g}_j\left(x\right)---\left(3\right)$

其中，L(x，λ)为拉格朗日函数；λ_j为p×1维拉格朗日乘子向量的第j个分量；

上式的KKT条件为：

$\left\{\begin{array}{c}\frac{\∂L(x,\mathrm{\λ})}{\∂x}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\frac{{\∂f}_i}{\∂x}+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j\frac{{\∂g}_i}{\∂x}=H^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}=0\\ \∂L(x,\mathrm{\λ})/\∂L=g\left(x\right)\\ F^{\′}\left(x\right)=\mathrm{diag}\left\{F_{\mathrm{ii}}\left(x\right)\right\}=\mathrm{diag}\left\{\frac{{2w}_i}{V+r_i^2}\right\}\end{array}\right.---\left(4\right)$

其中，H为h(x)对x的雅可比矩阵，G为零注入节点等式约束g(x)对x的雅可比矩阵；

利用牛顿法对上式的非线性方程进行求解，求偏导得到增广拉格朗日函数的海森矩阵：

$\left\{\begin{array}{c}{\∂}^{2}L/{\∂x}^2=H^T{F}^{\′}[I-\mathrm{diag}\{{2r}_i^2/(V+r_i^2)\}H=Q\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂x\∂\mathrm{\λ})=G^T\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂\mathrm{\λ}\∂x)=G,{\∂}^{2}L(x,\mathrm{\λ})/{\∂\mathrm{\λ}}^2=0\end{array}\right.---\left(5\right)$

其中，I为单位对角阵；

修正公式：

$\left[\begin{array}{c}{\mathrm{\Δx}}^{\left(k\right)}\\ {\mathrm{\Δ\λ}}^{\left(k\right)}\end{array}\right]=\left[\begin{array}{cc}Q& G^T\\ G& 0\end{array}\right]=\left[\begin{array}{c}{H}^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}\\ g\left(x\right)\end{array}\right]---\left(6\right)$

(五)修正变量；

x^(k+1)＝x^(k)+Δx^(k+1)，λ^(k+1)＝λ^(k)+Δλ^(k+1)(7)

(六)设定单次循环迭代t型抗差状态估计收敛精度ε₁，当满足max|Δx^(k)，Δλ^(k)|＜ε₁或k＞10时，单次循环迭代收敛，转步骤(七)；否则k＝k+1，转步骤(四)；

(七)设置l＝l+1，并重新计算自由度：

$V=\frac{\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i|z_i-h_i\left(x\right)|}{m-n}---\left(8\right)$

当满足V＜0.01或l＞2时，程序收敛，输出状态估计结果。收敛精度ε₁为人工设定，通常取10^-4，可满足工程要求。

所述量测数据z包括：节点电压幅值、支路首端有功功率和无功功率、支路末端有功功率和无功功率；

所述网络参数包括：母线编号、名称、补偿电容，输电线路的支路号、首端节点和末端节点编号、串联电阻、串联电抗、并联电导、并联电纳、变压器变比和阻抗。上述网络参数用于形成电力系统的导纳矩阵，从而便于计算各量测数据的估计值。

本发明提出的基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，与以往方法相比，具有如下优点：

(1)传统WLS状态估计方法结果易受不良数据影响，本发明的方法以量测残差对数函数值作为目标函数，能够有效减小不良数据对状态估计结果的影响，抗差性强，估计结果合理性好。

(2)t分布当V→∞时，趋于正态分布，对应着估计方法为WLS估计，当量测系统中不含不良数据时有较高的计算效率；当V取较小值时，可取得较好的抗差性。本方法能够动态调整自由度，从而保证方法兼顾估计的抗差性和效率。

(3)该方法在循环迭代求解过程中通过对自由度动态调整，可有效克服局部最优解问题。

本发明拟将t型估计引入电力系统状态估计，针对t型估计求解较繁琐的特点，提出近似模型，并以真实潮流作为等式约束。随后，利用牛顿法对该模型进行求解，结合一个简单的电力系统算例分析该方法在实际应用中可能存在的局部最优点问题，通过分布调整自由度保证取得全局最优解。实际上，目前在实际现场的应用中，抗差估计的计算量往往较大，本发明提出的自适应自由度调整策略，可大幅减小程序迭代次数，具有优越的工程应用价值。

附图说明

图1：本发明方法流程图。

图2：本发明提出的基于高斯-马尔科夫模型的自适应t型抗差状态估计方法所应用的4节点系统参数与量测配置图。

图3：本发明采用的IEEE57节点标准系统结构图。

具体实施方式

实施例一

如图1流程图所示，基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，依次包括以下步骤：

(一)获得电力系统的网络参数和量测数据z；

(二)以量测数据z建立t型抗差状态估计的模型：

$(\hat{x},\widehat{\mathrm{\σ}})=\arg \underset{x,\mathrm{\σ}>0}{\min }\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\left\{f\right[\frac{z_i-h_i\left(x\right)}{\mathrm{\σ}}]+\log \left(\mathrm{\σ}\right)\}---\left(1\right)$

其中，z_i、h_i分别为m×1维量测向量、量测函数向量的第i个分量；w_i为第i个量测对应的权重；σ为一个未知的尺度参数，将其设置为常数1；f(u)＝l0g(1+u²/V)，V是自由度；x为n×1维的状态向量，包括节点电压幅值与相角；m，n分别为量测量及状态量的个数；

对式(1)以节点注入功率g(x)作为等式约束得到：

$\left\{\begin{array}{ccc}\underset{x}{\min }& J\left(x\right)=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)& \\ s.t.& g\left(x\right)=0& \\ & r_i=z_i-h_i\left(x\right)& i=\mathrm{1,2},...,m\end{array}\right.---\left(2\right)$

其中，r_i为m×1维残差向量的第i个分量；

(三)令l＝0，并设置单次迭代计数器k=0。

(四)循环迭代下，利用内点罚函数法对式(2)进行处理：

$\underset{x,\mathrm{\λ}}{\min }L(x,\mathrm{\λ})=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j{g}_j\left(x\right)---\left(3\right)$

其中，L(x，λ)为拉格朗日函数；λ_j为p×1维拉格朗日乘子向量的第j个分量；

上式的KKT条件为：

$\left\{\begin{array}{c}\frac{\∂L(x,\mathrm{\λ})}{\∂x}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\frac{{\∂f}_i}{\∂x}+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j\frac{{\∂g}_i}{\∂x}=H^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}=0\\ \∂L(x,\mathrm{\λ})/\∂L=g\left(x\right)\\ F^{\′}\left(x\right)=\mathrm{diag}\left\{F_{\mathrm{ii}}\left(x\right)\right\}=\mathrm{diag}\left\{\frac{{2w}_i}{V+r_i^2}\right\}\end{array}\right.---\left(4\right)$

其中，H为h(x)对x的雅可比矩阵，G为零注入节点等式约束g(x)对x的雅可比矩阵；

利用牛顿法对上式的非线性方程进行求解，求偏导得到增广拉格朗日函数的海森矩阵：

$\left\{\begin{array}{c}{\∂}^{2}L/{\∂x}^2=H^T{F}^{\′}[I-\mathrm{diag}\{{2r}_i^2/(V+r_i^2)\}H=Q\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂x\∂\mathrm{\λ})=G^T\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂\mathrm{\λ}\∂x)=G,{\∂}^{2}L(x,\mathrm{\λ})/{\∂\mathrm{\λ}}^2=0\end{array}\right.---\left(5\right)$

其中，I为单位对角阵；

修正公式：

$\left[\begin{array}{c}{\mathrm{\Δx}}^{\left(k\right)}\\ {\mathrm{\Δ\λ}}^{\left(k\right)}\end{array}\right]=\left[\begin{array}{cc}Q& G^T\\ G& 0\end{array}\right]=\left[\begin{array}{c}{H}^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}\\ g\left(x\right)\end{array}\right]---\left(6\right)$

(5)修正变量；

x^(k+1)＝x^(k)+Δx^(k+1)，λ^(k+1)＝λ^(k)+Δλ^(k+1)(7)

(6)设定单次循环迭代t型抗差状态估计收敛精度ε₁，当满足max|Δx^(k)，Δλ^(k)|＜ε₁或k＞10时，单次循环迭代收敛，转步骤(7)；否则k＝k+1，转步骤(4)；

(7)设置l＝l+1，并重新计算自由度：

$V=\frac{\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i|z_i-h_i\left(x\right)|}{m-n}---\left(8\right)$

当满足V＜0.01或l＞2时，程序收敛，输出状态估计结果。收敛精度ε₁为人工设定，通常取10^-4，可满足工程要求。

为了进一步验证本发明所述方法的有效性，与含不良数据辨识的WLS(WLS+BD)状态估计、二次-常数(quadratic-constant，QC)状态估计方法进行比较，其中，不良数据辨识方法采用零残差搜索法(参见E.Handschin，F.C.Schweppe，J.Kohlas，etal.Baddataanalysisforpowersystemstateestimation[J].IEEETransactionsonPAS，1975，94(1)：329-336.)，QC估计参见方法：郭伟，单渊达.M估计方法及其在电力系统状态估计中的应用[J].中国电机工程学报，2000，20(9)：26-31。

为了验证本发明所提方法的有效性，在VisualC++平台上编制了自适应t型抗差估计。对IEEE标准系统的数据进行试验，通过将潮流结果加2％高斯随机误差获得量测生数据，通过对生数据置0的方式得到试验用的不良数据，所有量测数据均为标幺值(pu)。下面进一步介绍本发明的二个实施例：

实施例二

电力系统中评价2个量测数据间的强相关性，由于残差r与量测误差v存在以下关系：

r＝Kv

其中，K为残差灵敏度矩阵，一般取：

K＝I-H(H^TWH)^-1H^TW

对于t型估计，若量测i的残差较大，则对角元素即对该量测降权处理，由此可以得出，t型估计模型中抑制不良数据的影响是可能的。

利用实施例一的方法对4节点系统进行分析，配置及系统参数如图2和表1所示。由图2可知，该系统共含4个母线、4条交流输电线路，并以1号母线作为平衡节点，共配置7对支路有功/无功功率量测数据、2号母线的节点注入有功/无功功率量测数据，以及2号母线电压幅值量测数据。在给定的量测系统和网络结构下，找出具有强相关性的量测对(CM)，并分析不良数据对该量测估计值的影响。相关性分析可通过统计K中的相关系数，量测i与量测j之间的相关度指标：

$P_{\mathrm{ij}}=\frac{|K_{\mathrm{ii}}{K}_{\mathrm{ij}}+K_{\mathrm{ij}}{K}_{\mathrm{jj}}|}{\sqrt{K_{\mathrm{ii}}^2+K_{\mathrm{ij}}^2}\·\sqrt{K_{\mathrm{ji}}^2+K_{\mathrm{jj}}^2}}$

当相关度指标满足0.7521＜P_ij＜0.9498时，为强相关量测对。该系统含4组强相关量测，各估计方法的强相关量测指标如表1所示。该4组量测对于WLS估计，均为强相关量测，而采用t型估计后，相关度指标显著减小。

表14节点系统强相关量测指标

对上述4组强相关量测添加不良数据，不同方法的估计结果如表2所示。以支路1-2有功功率(P_1-2)为例说明不同方法的抗差性能，由表2可知，该量测的真值为0.3883，WLS估计结果为0.2465，残差为0.1418，而t型估计值为0.3876，残差仅为0.0007，很好地辨识出不良数据。此外，对于该支路P_2-1功率量测，t型估计的误差也仅为-0.0006，表明t型估计成功辨识出了该组强相关不良数据。

表24节点系统强相关量测不同方法的估计结果比较

实施例三

进一步，对IEEE57节点系统进行分析，验证自适应调整自由度策略的有效性，系统结构图如图3所示。利用本发明的方法对IEEE-57节点系统进一步分析，采用定V的方式(以WLS状态估计结果作为系统状态量初值)，进行t型状态估计计算。同时，为了验证本发明动态调整自由度V策略的有效性，表3给出了上述不同自由度下状态估计的结果比较。其中，S₁及S₂指标如下：

$S_1=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}|x_i^{*}-x_i|$

$S_2=\max |x_i^{*}-x_i|$

其中，为第i个状态变量的估计值，x_i为状态变量的真实值。S₁代表状态估计的整体误差，S₂代表状态估计局部误差。

表3不同V值的t型估计结果

由表3可知，随着自由度V的减小，S₁及S₂指标逐渐减小，即较小的自由度具有更强的抗差能力。当自由度V＝1时，S₁指标为0.27415，平均误差为0.00243，最大误差为0.00481。实施例1的方法在第1次迭代后利用式(8)计算得到V＝0.0129，继续进行迭代计算。经过5次迭代后收敛后V＝0.00268，满足收敛条件。最终t型估计结果S1指标进一步减小至0.01638，仅为第1次迭代时的5.975％。同时，与定V方式比较可知，该方法并未陷入局部最优解。在计算效率上，定V策略经4轮迭代总时间为26.014ms，而本发明的方法仅为11.892ms。实验结果表明，本发明的混合策略在保证求解精度的前提下，可显著提高收敛速率。

Claims (2)

1.基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，其特征在于该方法依次包括以下步骤：

(一)获得电力系统的网络参数和量测数据z；

(二)以量测数据z建立t型抗差状态估计的模型：

$(\hat{x},\widehat{\mathrm{\σ}})=\arg \underset{x,\mathrm{\σ}>0}{\min }\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\left\{f\right[\frac{z_i-h_i\left(x\right)}{\mathrm{\σ}}]+\log \left(\mathrm{\σ}\right)\}---\left(1\right)$

其中，z_i、h_i分别为m×1维量测向量、量测函数向量的第i个分量；w_i为第i个量测对应的权重；σ为一个未知的尺度参数，将其设置为常数1；f(u)＝log(1+u²/V)，V是自由度；x为n×1维的状态向量，包括节点电压幅值与相角；m，n分别为量测量及状态量的个数；

对式(1)以节点注入功率g(x)作为等式约束得到：

$\left\{\begin{array}{ccc}\underset{x}{\min }& J\left(x\right)=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)& \\ s.t.& g\left(x\right)=0& \\ & r_i=z_i-h_i\left(x\right)& i=\mathrm{1,2},...,m\end{array}\right.---\left(2\right)$

其中，r_i为m×1维残差向量的第i个分量；

(三)令l＝0，并设置单次迭代计数器k＝0；

(四)循环迭代下，利用内点罚函数法对式(2)进行处理：

$\underset{x,\mathrm{\λ}}{\min }L(x,\mathrm{\λ})=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_{i}f\left(r_i\right)+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j{g}_j\left(x\right)---\left(3\right)$

其中，L(x，λ)为拉格朗日函数；λ_j为p×1维拉格朗日乘子向量的第j个分量；

上式的KKT条件为：

$\left\{\begin{array}{c}\frac{\∂L(x,\mathrm{\λ})}{\∂x}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i\frac{{\∂f}_i}{\∂x}+\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{\mathrm{\λ}}_j\frac{{\∂g}_i}{\∂x}=H^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}=0\\ \∂L(x,\mathrm{\λ})/\∂L=g\left(x\right)\\ F^{\′}\left(x\right)=\mathrm{diag}\left\{F_{\mathrm{ii}}\left(x\right)\right\}=\mathrm{diag}\left\{\frac{{2w}_i}{V+r_i^2}\right\}\end{array}\right.---\left(4\right)$

其中，H为h(x)对x的雅可比矩阵，G为零注入节点等式约束g(x)对x的雅可比矩阵；

利用牛顿法对上式的非线性方程进行求解，求偏导得到增广拉格朗日函数的海森矩阵：

$\left\{\begin{array}{c}{\∂}^{2}L/{\∂x}^2=H^T{F}^{\′}[I-\mathrm{diag}\{{2r}_i^2/(V+r_i^2)\}H=Q\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂x\∂\mathrm{\λ})=G^T\\ {\∂}^{2}L(x,\mathrm{\λ})/(\∂\mathrm{\λ}\∂x)=G,{\∂}^{2}L(x,\mathrm{\λ})/{\∂\mathrm{\λ}}^2=0\end{array}\right.---\left(5\right)$

其中，I为单位对角阵；

修正公式：

$\left[\begin{array}{c}{\mathrm{\Δx}}^{\left(k\right)}\\ {\mathrm{\Δ\λ}}^{\left(k\right)}\end{array}\right]=\left[\begin{array}{cc}Q& G^T\\ G& 0\end{array}\right]=\left[\begin{array}{c}{H}^T{F}^{\′}\left(x\right)r+G^T\mathrm{\λ}\\ g\left(x\right)\end{array}\right]---\left(6\right)$

(五)修正变量；

x^(k+1)＝x^(k)+Δx^(k+1)，λ^(k+1)＝λ^(k)+Δλ^(k+1)(7)

(六)设定单次循环迭代t型抗差状态估计收敛精度ε₁，当满足max|Δx^(k)，Δλ^(k)|＜ε₁或k＞10时，单次循环迭代收敛，转步骤(七)；否则k＝k+1，转步骤(四)；

(七)设置l＝l+1，并重新计算自由度：

$V=\frac{\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{w}_i|z_i-h_i\left(x\right)|}{m-n}---\left(8\right)$

当满足V＜0.01或l＞2时，程序收敛，输出状态估计结果。

2.基于高斯-马尔科夫模型的自适应t型抗差状态估计方法，其特征在于，所述量测数据z包括：节点电压幅值、支路首端有功功率和无功功率、支路末端有功功率和无功功率。