CN111273084B - 一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法 - Google Patents

Abstract

本发明提出了一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，通过分析残差序列在状态变量保持不变和突变两种状态下的特征构建判定变量，利用判定变量判定状态变量是否发生突变，进而自适应改变弱化因子的值，从而在保证屏蔽残差序列协方差的估计误差的基础上尽可能地减小弱化因子的值，从而加快强跟踪无迹卡尔曼滤波算法跟踪状态变量突变的速度，提高相量测量的准确性。

Description

一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法

技术领域

本发明属于广域同步相量测量领域，特别涉及一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法。

背景技术

由于实时性和测量精度上的优势，近年来卡尔曼滤波方法逐渐被国内外学者关注，使用该方法的同步相量测量方法，首先要建立电力信号的状态空间模型，然后利用依次得到的采样点对电力信号状态变量及量测量进行状态预测与滤波更新，从而达到对电力信号状态变量依次估计进而逐渐接近真实值的目的。由于电力信号与其频率的函数一定存在着非线性，致使只能解决线性问题的经典卡尔曼滤波方法无法解决电力信号参数测量问题，必须对电力信号状态空间模型进行线性化。由于无迹卡尔曼滤波方法能以二阶及以上的精度逼近非线性函数且能快速跟踪状态变量突变，近年来其在同步相量测量方法中应用越来越广。

无迹卡尔曼滤波方法采用前一时刻的状态变量误差协方差矩阵和UT变换预测新时刻的状态变量误差协方差，因此当状态变量发生突变时无迹卡尔曼方法无法检测出突变，致使预测的状态变量误差协方差仍然较小，进而导致方法对突变跟踪缓慢。强跟踪无迹卡尔曼方法根据正交性原理引入渐消因子来保证残差序列实时弱相关，增大突变时刻的预测状态变量误差协方差，进而达到跟踪状态变量突变的目的。计算渐消因子需要已知残差序列的协方差信息，但这种信息实际是未知的，目前存在的强跟踪无迹卡尔曼滤波方法采用对残差序列协方差进行估计的方法来解决这一问题，并引入弱化因子来屏蔽由残差序列协方差估计不准带来的方法发散影响。高噪声环境下，残差序列协方差估计精度急剧降低，这就需要设置较大的弱化因子来对更大的估计误差进行屏蔽，但过大的弱化因子会造成方法跟踪突变缓慢或无法跟踪突变，过小的弱化因子会由于无法屏蔽残差序列的估计误差造成方法发散。

发明内容

针对上述问题，本发明提出一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，包括如下步骤：

步骤1：按下式计算判定变量，判断判定变量是否在下分位点和上分为点之间，

式中：F_k为判定变量，ξ_i为i次迭代时的残差，N为时间窗长度；

步骤2：若所述判定变量在下分位点和上分为点之间，则按下式计算弱化因子的值，按无迹卡尔曼滤波方法的步骤更新状态变量误差协方差矩阵，而后进行量测更新和滤波更新，

式中：β_k为第k次迭代时弱化因子的取值，k₀为最近一次发生突变的迭代次数，k为目前迭代次数，V_i为第i次迭代时的残差序列协方差矩阵估计值，

为第i次迭代时量测量的后验协方差矩阵，R_i为第i次迭代时的噪声协方差矩阵；

步骤3：若所述判定变量不在下分位点和上分为点之间，则沿用上一时刻的弱化因子的值带入计算渐消因子，按强跟踪无迹卡尔曼方法的步骤更新状态变量误差协方差矩阵，而后进行量测更新和滤波更新。

进一步，若

成立，判定为稳定状态；反之，则判定为突变状态。

进一步，稳定状态时，系统与实际相匹配，渐消因子λ₀的值为1；当实际输出残差序列的协方差矩阵V_k预测不准时，稳定状态会出现渐消因子大于1的情况，此时应调节弱化因子的值将渐消因子强制为1，并在突变状态时延用此弱化因子。

进一步，上、下分位点依时间窗长度选择。

附图说明

图1是自适应强跟踪无迹卡尔曼滤波方法流程图；

图2是本发明的自适应强跟踪无迹卡尔曼滤波方法和强跟踪无迹卡尔曼滤波方法对每次采样的采样点进行滤波，得到的每次检测到系统突变的时刻的对比图；

图3是在突变状态时，强跟踪无迹卡尔曼滤波方法和本发明所提自适应强跟踪无迹卡尔曼滤波方法的幅值的RMSE差值图；

图4是在突变状态时，强跟踪无迹卡尔曼滤波方法和本发明所提自适应强跟踪无迹卡尔曼滤波方法的频率的RMSE差值图；

图5是在突变状态时，强跟踪无迹卡尔曼滤波方法和本发明所提自适应强跟踪无迹卡尔曼滤波方法的相角的RMSE差值图；

图6是利用强跟踪无迹卡尔曼滤波方法和本发明的自适应强跟踪无迹卡尔曼滤波方法对测试信号进行一次基波相量参数估计的幅值估计误差的绝对值对比图；

图7是利用强跟踪无迹卡尔曼滤波方法和本发明的自适应强跟踪无迹卡尔曼滤波方法对测试信号进行一次基波相量参数估计的频率估计误差的绝对值对比图；

图8是利用强跟踪无迹卡尔曼滤波方法和本发明的自适应强跟踪无迹卡尔曼滤波方法对测试信号进行一次基波相量参数估计的相角估计误差的绝对值对比图。

具体实施方式

下面结合附图，对本发明进行详细说明。

本发明提出了一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，通过分析残差序列在状态变量保持不变和突变两种状态下的特征构建判定变量，利用判定变量判定状态变量是否发生突变，进而自适应改变弱化因子的值，从而在保证屏蔽残差序列协方差的估计误差的基础上尽可能地减小弱化因子的值，从而加快强跟踪无迹卡尔曼滤波方法跟踪状态变量突变的速度。

1.电力系统基波相量状态空间模型

电力系统中电力信号在k时刻的采样值可表示为：

式中：y_k为信号的瞬时采样值；A₀、f₀、

分别为基波信号的幅值、频率和初相角；A_i、f_i、

(i≠0)分别为第(i+1)次谐波的幅值、频率和初相角；k为采样次数；T_s为采样间隔；v_k为与信号相互独立的高斯白噪声，其均值为0、方差为σ²。

且考虑到方法的实时性要求，系统的状态变量不应过多，PMU提供的同步相量信息在一般情况下更关注基波信息，于是本方法将谐波纳入噪声中滤除，定义电力信号的状态变量为：

易知，系统的状态转移方程为：

系统的观测方程为：

y_k＝x_k(1)sin[x_k(2)]+v_k (4)

由式(2)可知，迭代k次时的基波参量可按下式计算：

幅值：

A_0,k＝x_k(1) (5)

频率：

f_0,k＝x_k(3)/(2π) (6)

初相角：

2.强跟踪无迹卡尔曼滤波方法

以对称采样为例，介绍强跟踪无迹卡尔曼滤波方法如下：

考虑如下非线性高斯系统：

x_k＝f_k-1(x_k-1)+u_k-1

y_k＝h_k(x_k)+v_k (8)

式中：x_k∈R_n和y_k∈R_m分别为k时刻下n维状态变量与m维量测量；u_k-1～N(0,Q_k)和v_k～N(0,R_k)分别为系统过程噪声和量测噪声；f(·)为非线性状态转移函数；h(·)为非线性量测函数。

假设k-1时刻状态变量服从分布：

按照对称采样策略选取2n+1个sigma点：

对应于ξ_i(i＝0,1,...,2n)的权值为：

式中：c＝n+λ；λ为比例系数，取0；chol(·)为矩阵的Choleshy分解；[·]_i代表括号内矩阵的第i列。

(1)计算未引入渐消因子时的状态预测协方差阵、自协方差阵和互协方差阵：

式中：

为由

和

所构造的Sigma采样点。

(2)根据得到的

与

计算渐消因子：

V_k为实际输出残差序列的协方差阵，在实际中V_k是未知的，可由下式估算：

式中，β≥1为弱化因子，ξ_k为k时刻的残差序列，

为V_k的计算值，0＜ρ≤1为遗忘因子，取ρ＝0.95。

(3)在状态预测协方差阵

中引入渐消因子，重新计算

(4)量测更新

按

和

重新采样后，按照对称采样策略计算sigma点ξ_k|k-1,i(i＝0,1,...,2n)，通过非线性量测函数传播后，计算量测量均值、自协方差阵和互协方差阵：

获得新的量测值y_k后，进行滤波更新：

3.自适应强跟踪无迹卡尔曼滤波方法

在方法收敛后M_k的迹很小，造成

很小的估计误差传导至渐消因子中时被放大多倍，最终使得方法发散或无法正常跟踪突变。为解决V_k预测不准的问题，强跟踪无迹卡尔曼滤波在方法中引入了弱化因子β减小矩阵N_k的迹，屏蔽V_k预测不准造成的影响。

β的取值一般按经验选取，然而弱化因子过大会掩盖突变讯息导致方法对突变跟踪缓慢，而过小的弱化因子则无法屏蔽

估计误差使得方法发散。随着噪声含量增大，为防止方法发散，β的取值将会逐渐增大，方法对突变跟踪逐渐减缓。为此，本发明的方法中，通过自适应判定系统状态，尽可能减小弱化因子的取值。为便于说明，将状态变量未发生变化时的状态称为稳定状态，将状态变量突变时的状态称为突变状态。

(1)状态量突变检测

当方法收敛后，量测量真实值与估计值之间差值很小，噪声含量较大时其差值可近似忽略，收敛后在稳定状态时：

ξ_kξ_k ^T≈v_kv_k ^T (20)

为判定系统状态构建判定变量：

式中：N为时间窗长度。

由式(20)可知，收敛后稳定状态时：

此时，若噪声含量不变，随机变量F_k应近似服从自由度为N的卡方分布，

当α很小时，可认为

和

是小概率事件，几乎不可能发生。R_k为量测噪声协方差矩阵。

突变状态时，量测量真实值与估计值之间差值较大不可忽略，此时：

与原卡方分布差值明显，其取值大概率不在卡方分布下分位点和上分位点的区间内，根据此差异即可区分系统状态。

(2)弱化因子自适应计算

稳定状态时，系统与实际相匹配，渐消因子λ₀的值应当为1。当V_k预测不准时，稳定状态会出现渐消因子大于1的情况。此时应调节弱化因子的值将渐消因子强制为1，并在突变状态时延用此弱化因子。易知，稳定状态时弱化因子的值为：

为保证数学模型与实际匹配时任意时刻渐消因子均能为1，弱化因子的计算方法应变为：

(3)自适应强跟踪无迹卡尔曼滤波方法步骤，如图1所示。

基于以上理论，构建自适应强跟踪无迹卡尔曼滤波方法如下：

step1：选取适当的卡方分布上分位点和下分位点，设定状态量初始误差协方差、弱化因子、系统噪声协方差及状态变量初始值；

step2：按式(9)选取sigma点后，按式(11)、(12)计算未引入渐消因子时的状态预测协方差阵、自协方差阵和互协方差阵；

step3：获得新的采样点后，按式(14)、(15)、(16)计算V_k、N_k、M_k；

step4：构建时间窗，根据式(21)计算判定变量F_k；

step5：判断

是否成立，若成立，判定为稳定状态，反之则判定为突变状态；

step6：稳定状态时，按式(25)计算弱化因子的值β_k；突变状态时，保持前一时刻弱化因子的值不变，按式(13)计算渐消因子，并按式(17)更新状态预测协方差阵，重新采样后按式(18)进行量测更新；

step7：按式(19)进行滤波更新，根据更新后的状态变量值以及式(5)至式(7)输出此时的同步相量参数测量值；

step8：重复步骤2—7，直至采样结束或满足方法停止条件。

下面以具体的实施例来验证本发明方法的优越性。

本发明采用戴尔笔记本电脑(windows10系统、CPU型号为i5-8550U)作为仿真实验平台。通过将本发明所提方法与强跟踪无迹卡尔曼滤波方法对比验证所提方法有效性，测试信号如下：

其中，ω₁、ω₂分别为100.8π和101π，采样频率为9600Hz，即采样时间间隔T_s为1/9600秒，总采样时间为2秒，v_k为高斯白噪声，v_k～N(0,0.185²)，即信噪比为15dB。

两种方法参数设置如下：

强跟踪无迹卡尔曼滤波方法弱化因子β设置为16。自适应强跟踪无迹卡尔曼滤波方法选取时间窗长度为960个采样点，卡方分布上分位点选为1200，下分位点选为750。初始状态变量

系统量测噪声协方差阵R＝0.185²，初始误差协方差矩阵P₀＝I_3×3。

系统状态判定方法性能验证

为检验本发明所提系统状态判定方法有效性，本发明对式(26)的测试信号采样40次，分别使用自适应强跟踪无迹卡尔曼滤波方法和强跟踪无迹卡尔曼滤波方法对每次采样的采样点进行滤波，得到两种方法每次检测到系统突变的时刻。其中，强跟踪无迹卡尔曼滤波方法以得到较大渐消因子的时刻作为检测到系统突变的时刻。结果对比如图2所示。

由图2可知，强跟踪无迹卡尔曼滤波方法一般突变后1e-4s即可检测到突变，但最差需要5e-4s，且有6次未检测出系统发生了突变；本发明所提自适应强跟踪无迹卡尔曼滤波方法大概率在1e-4s及以下检测到突变，多数情况下比强跟踪无迹卡尔曼滤波方法快，且每次都可靠的检测出系统的突变。也即，本发明所提自适应强跟踪无迹卡尔曼滤波方法的系统状态判定方法检测到突变的时刻大概率先于强跟踪无迹卡尔曼滤波方法，且较强跟踪无迹卡尔曼滤波方法可靠性高。

突变状态自适应强跟踪无迹卡尔曼滤波方法性能分析

为评价方法性能，利用Matlab仿真软件工具，在突变状态时对强跟踪无迹卡尔曼滤波方法和本发明所提自适应强跟踪无迹卡尔曼滤波方法进行方法性能比较。为了评价方法的性能，定义幅值、频率、相角的均方根误差(root mean squared error，RMSE)为：

其中m为需要计算的采样点数，

表示状态变量估计值，x_k表示状态变量的实际值，其中电力信号的状态空间模型以及状态变量取值见6.1。

对式(26)的测试信号采样40次，截取了两种方法1.9s至2.0s的数据计算均方根误差。为便于比较，作图时已将强跟踪无迹卡尔曼滤波方法和本发明所提自适应强跟踪无迹卡尔曼滤波方法的RMSE作差，差值为正时强跟踪无迹卡尔曼滤波方法的RMSE比自适应强跟踪无迹卡尔曼滤波方法大。对比结果由图3-图5示出。由图3-图5可知，强跟踪无迹卡尔曼滤波方法的幅值RMSE大概率比自适应强跟踪无迹卡尔曼滤波方法的大，自适应强跟踪无迹卡尔曼滤波方法的频率和相角RMSE均基本比强跟踪无迹卡尔曼滤波方法小。

表2突变状态估计精度比较(RMSE)

由表2可知，本发明所提自适应强跟踪无迹卡尔曼滤波方法的幅值、频率、相角的RMSE平均值均比强跟踪无迹卡尔曼滤波方法的低，其性能优于强跟踪无迹卡尔曼滤波方法。

为较为直观地比较两种方法跟踪突变的速度和方法收敛后的精度，利用强跟踪无迹卡尔曼滤波方法和自适应强跟踪无迹卡尔曼滤波方法对式(26)的测试信号进行一次基波相量参数估计，为方便对比，输出结果被设置为估计误差的绝对值，对比结果由图6-图8示出。

由图6-图8可知，由于自适应强跟踪无迹卡尔曼滤波方法幅值、频率及相角均能较强跟踪无迹卡尔曼滤波方法更快的接近真实值。其中自适应强跟踪无迹卡尔曼滤波方法的幅值收敛过程中存在精度比强跟踪无迹卡尔曼滤波方法差的情况，但收敛后，自适应强跟踪无迹卡尔曼滤波方法的幅值估计精度明显优于强跟踪无迹卡尔曼滤波方法；自适应强跟踪无迹卡尔曼滤波方法频率收敛速度较强跟踪无迹卡尔曼滤波方法有明显的提升，且收敛后频率估计精度明显比强跟踪无迹卡尔曼滤波方法高；自适应强跟踪无迹卡尔曼滤波方法相角的估计误差也能很快的达到很小的值，收敛后，自适应强跟踪无迹卡尔曼滤波方法相角的估计误差比强跟踪无迹卡尔曼滤波方法小。

此实施例仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

Claims (4)

1.一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，包括如下步骤：

步骤1：按下式计算判定变量，判断判定变量是否在下分位点和上分位点之间，

式中：F_k为判定变量，ξ_i为i次迭代时的残差，N为时间窗长度，R_k为量测噪声协方差矩阵；

步骤2：若所述判定变量在下分位点和上分位点之间，则按下式计算弱化因子的值，按无迹卡尔曼滤波方法的步骤更新状态变量误差协方差矩阵，而后进行量测更新和滤波更新，

式中：β_k为第k次迭代时弱化因子的取值，k₀为最近一次发生突变的迭代次数，k为目前迭代次数，V_i为第i次迭代时的残差序列协方差矩阵估计值，

为第i次迭代时量测量的后验协方差矩阵，R_i为第i次迭代时的噪声协方差矩阵；

步骤3：若所述判定变量不在下分位点和上分位点之间，则沿用上一时刻的弱化因子的值带入计算渐消因子，按强跟踪无迹卡尔曼方法的步骤更新状态变量误差协方差矩阵，而后进行量测更新和滤波更新。

2.根据权利要求1所述的一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，其特征在于：若

成立，判定为稳定状态；反之，则判定为突变状态。

3.根据权利要求1所述的一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，其特征在于：稳定状态时，系统与实际相匹配，渐消因子λ0的值为1；当实际输出残差序列的协方差矩阵V_k预测不准时，稳定状态会出现渐消因子大于1的情况，此时应调节弱化因子的值将渐消因子强制为1，并在突变状态时延用此弱化因子。

4.根据权利要求1所述的一种自适应强跟踪无迹卡尔曼滤波同步相量参数估计方法，其特征在于：上、下分位点依时间窗长度选择。

Images