CN111693773A - 一种基于互近似熵数据筛选的谐波源责任划分方法 - Google Patents

Abstract

本发明针对如何解决背景谐波电压波动的各种仿真工况，分析比较所提方法与传统线性回归法的优劣的技术需求，本发明提供一种基于CAE数据筛选的谐波源责任划分方法，首先将采集的实测数据划分为若干区段，将各区段实测谐波电压、电流数据进行互近似熵值计算，保留满足CAE阈值要求的区段，以达到排除背景谐波干扰的目的；然后利用M估计稳健回归法对保留数据进行回归计算，最大限度避免异常值对回归计算的影响，得到系统谐波阻抗，从而实现准确的谐波责任划分。

Description

一种基于互近似熵数据筛选的谐波源责任划分方法

技术领域

本发明涉及一种基于互近似熵数据筛选的谐波源责任划分方法，属于电能质量控制技术领域。

背景技术

目前，随着大量的电力电子设备接入电网，其造成的谐波污染已成为电网电能质量的突出问题之一。为确保电网的供电可靠性，保证供电质量，对导致电网谐波超标的用电用户进行相应的惩罚势在必行。准确衡量用户谐波发射水平、明确划分谐波责任是对电力用户进行合理奖惩的前提。用户谐波发射水平的估计结果是由系统等值谐波阻抗的计算得到，因此评价用户谐波发射水平问题的关键在于能否对系统谐波阻抗进行准确的估计。

现有的谐波阻抗估计方法主要分为非干预式和干预式两类。干预式通过在系统中制造短时扰动，利用公共连接点(point of common coupling，PCC)处的暂态谐波电压和电流增量的计算来实现谐波阻抗的估计。主要包括谐波电流注入法、投切电容器法等。但干预式方法会对电网运行造成干扰，因此目前干预式法在谐波阻抗测量领域并不常用。

非干预式的方法主要包括波动量法、回归法等。采用波动量法是若不满足系统侧谐波电流源波动程度远远小于负荷侧波动程度时的前提条件，估计误差较大。对于线性回归法来说，其原理是根据谐波等效电路中PCC点处的电压、电流测量值，求解其对应方程回归系数，从而得到谐波阻抗。然而背景谐波的干扰是影响非干预式谐波估计准确性的关键因素，但传统线性回归法在背景谐波电压波动时拟合度较差，而实际中背景谐波电压往往是变化的，且具有随机性，因此传统线性回归法估计误差较大。

发明内容

本发明的目的在于，克服现有技术存在的技术缺陷，提出一种基于互近似熵数据筛选的谐波源责任划分方法，针对如何解决背景谐波电压波动的各种仿真工况，分析比较所提方法与传统线性回归法的优劣的技术需求。

本发明具体采用如下技术方案：一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，包括如下步骤：

步骤SS1：采集谐波数据电压U_pcc和谐波数据电流I_pcc，将所述谐波数据电压U_pcc和谐波数据电流I_pcc进行归一化处理；

步骤SS2：将归一化后的谐波数据划分为n个数据区段，每个数据区段包含L个数据点；

步骤SS3：计算每个数据区段的互近似熵值CAE，并将所述互近似熵值CAE与设定阈值比较，对满足设定阈值条件的数据进行保留，反之则进行舍弃；

步骤SS4：将所有数据区段的数据筛选完成后，利用M估计稳健回归法计算系统谐波阻抗，得到谐波源责任划分结果。

作为一种较佳的实施例，所述步骤SS1具体包括：

设I_ch为用户侧等效谐波源，I_sh为系统侧等效谐波源，Z_ch为用户侧等效谐波阻抗，Z_sh为系统侧等效谐波阻抗，h表示谐波次数；I_pcc-s为系统侧等效谐波源I_s贡献谐波电流，V_s为系统侧贡献谐波电压；I_pcc-c为用户侧等效谐波源I_c贡献谐波电流，V_c为用户侧贡献谐波电压；各参数表达如公式(1)～(4)所示：

I_pcc＝I_pcc-s+I_pcc-c (3)

V_pcc＝V_s+V_c (4)

由于电力系统中Z_c远大于Z_s，则上述关系式等效改写为：

联立式(5)和(6)得到:

V_pcc＝Z_s.I_pcc+V_s (7)

上述公式(1)～(7)中的元素皆为方向矢量；

当背景谐波波动很小时，Vs近似看作常数，因此由式(7)知，当系统谐波阻抗不变的情况下，V_pcc与I_pcc呈线性关系；V_pcc的波形趋势应与I_pcc的波形趋势一致；

式(7)模型的表达式如下：

V_pcc＝|Z_s|cosθ.I_pcc+|V_s|cosθ₁ (8)

其中θ为Z_s·I_pcc与V_pcc的夹角，θ1为V_s与V_pcc的夹角；则对于用户侧谐波责任H_c的计算由公式(9)表示：

当背景谐波保持平稳时，将式(8)中的|Z_s|cosθ和|V_s|cosθ₁近似看作常系数。

作为一种较佳的实施例，所述步骤SS1中的归一化处理的计算公式如下：

其中，x_N为归一化结果，x_min为全体数据中的最小值，x_max为全体数据中的最大值。

作为一种较佳的实施例，所述步骤SS2中的L取10。

作为一种较佳的实施例，所述步骤SS3中的所述互近似熵值CAE的计算方法如下：

设两个不同的时间序列分别为i(t)和j(t)，规定一个长度为m的窗口，如公式(10)和(11)所示，分别对i(t)和j(t)构造N-m+1个m维矢量X_p,和X_q，其中：

X_p＝[i(p),…,i(p+m-1)]

p＝1,2,…,N-m+1 (10)

X_q＝[j(q),…,j(q+m-1)]

q＝1,2,…,N-m+1 (11)

使用矢量的∞—范数描述X_p与X_q之间的距离d(X_p,X_q)＝||X_p-X_q||_∞；

给定相似容限r，r的表达式为：

r＝0.2×COV(X_P,X_q) (12)

对每一个p值统计X_p和所有X_q(q＝1,2,…,N-m+1)的矢量距离小于r的个数N_p,m,r，并计算N_p,m,r与总的矢量个数(N-m+1)的比值C_p,m,r：

由公式(13)可知，小于r的矢量个数越多，就意味着两数据段越相近；则两数据区段的互相关程度为：

增加窗口长度至m+1，重复上述公式(10)～公式(14)的运算过程，得到Φ_m+1,r；

最后，计算得到与m，r相关的互近似熵值：

CAE(m,r)＝Φ_m,r-Φ_m+1,r (15)

作为一种较佳的实施例，所述步骤SS4具体包括：

步骤SS41：利用最小二乘法计算出相应的回归系数a_i(a₁…a_n)，b_i(b₁…b_n)，并计算相应的误差量e_i；

步骤SS42：根据计算残差尺度s，计算相应的权重值ω；

步骤SS43：将权重计算结果ω_i(ω₁…ω_n)代入M估计稳健回归方程，求得新的回归计算结果；

比较前后两次计算的回归参数，若a_i,b_i全部满足阈值约束则保留，约束表达式为：

一般δ₀取值为10^-5，否则继续迭代计算，直到满足为止。

作为一种较佳的实施例，所述M估计稳健回归方程为：

其中，a和b均为回归系数，ω_i为权重计算结果。

作为一种较佳的实施例，所述残差尺度s的技术公式如下：

s＝Med|e_i-Med(e₁,...,e_n)|/0.6745 (17)

其中，Med表示一组数据的中位数。

作为一种较佳的实施例，所述权重值ω采用Huber法来计算：

式中：c为常数，其值一般取1.345，u_i为标准化的残差指标，u_i＝e_i/s。

作为一种较佳的实施例，所述步骤SS41具体包括：已知实测数据(x_i,y_i)(i＝1,2,…,n)，设线性关系式为y＝ax+b，则利用最小二乘法求取方程系数a和b；

残差e_i的目标函数表示为：

其中，a和b均为回归系数。

本发明所达到的有益效果：第一，本发明针对如何解决背景谐波电压波动的各种仿真工况，分析比较所提方法与传统线性回归法的优劣的技术需求，本发明提供一种基于CAE数据筛选的谐波源责任划分方法，首先将采集的实测数据划分为若干区段，将各区段实测谐波电压、电流数据进行互近似熵值计算，保留满足CAE阈值要求的区段，以达到排除背景谐波干扰的目的。然后利用M估计稳健回归法对保留数据进行回归计算，最大限度避免异常值对回归计算的影响，得到系统谐波阻抗，从而实现准确的谐波责任划分。第二，本发明与传统线性回归法相比，该方法抵抗背景谐波波动的能力更强，在各种背景谐波电压波动的情况下均能有较高的估计精度。第三，本发明利用M估计稳健回归法对保留的数据进行回归计算，最大限度避免异常值对回归计算的影响，进一步提高谐波阻抗估计精度。

附图说明

图1为诺顿等效电路图；

图2为责任划分等效原理图；

图3为本发明优选实施例的算法执行流程图；

图4为仿真电路图；

图5为PCC点5次谐波电压图；

图6为PCC点5次谐波电流图；

图7为回归误差随CAE阈值取值变化曲线图；

图8为筛选前谐波电压、电流数据点图；

图9为筛选后谐波电压、电流数据点图。

具体实施方式

下面结合附图对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案，而不能以此来限制本发明的保护范围。

实施例1：谐波责任划分的基本阐述；电力系统中谐波责任划分可用诺顿等值电路表示，如下图1所示，公共连接点的谐波电压、电流值由系统侧和用户侧的谐波源同时贡献。其中I_ch为用户侧等效谐波源，I_sh为系统侧等效谐波源，Z_ch为用户侧等效谐波阻抗，Z_sh为系统侧等效谐波阻抗，h表示谐波次数。

根据叠加定理，可将图1转换为系统侧和用户侧两部分叠加而成，如图2所示：

图中I_pcc-s为系统侧等效谐波源I_s贡献谐波电流，V_s为系统侧贡献谐波电压。同理，I_pcc-c为用户侧等效谐波源I_c贡献谐波电流，V_c为用户侧贡献谐波电压。各参数表达如公式(1)～(4)所示：

I_pcc＝I_pcc-s+I_pcc-c (3)

V_pcc＝V_s+V_c (4)

由于电力系统中Z_c远大于Z_s，，则上述关系式可等效改写为：

联立式(5),(6)可得，

V_pcc＝Z_s.I_pcc+V_s (7)

上述公式(1)～(7)中元素皆为方向矢量。

当背景谐波波动很小时，Vs可以近似看作常数，因此由式(7)可知，当系统谐波阻抗不变的情况下，V_pcc与I_pcc呈线性关系。V_pcc的波形趋势应与I_pcc的波形趋势一致。

式(7)模的表达式如下：

V_pcc＝|Z_s|cosθ.I_pcc+|V_s|cosθ₁ (8)

其中θ为Z_s·I_pcc与V_pcc的夹角，θ1为V_s与V_pcc的夹角。则对于用户侧谐波责任H_c的计算可由公式(9)表示：

当背景谐波保持平稳时，可将式(8)中的|Z_s|cosθ和|V_s|cosθ₁近似看作常系数，利用线性回归算法进行求取。

实施例2：互近似熵原理阐述；由上述原理可知，若排除背景谐波的影响，V_pcc与I_pcc的波形模式一致。由此计算而得的系统谐波阻抗的准确度会更高。因此，为了筛选出V_pcc与I_pcc波形趋势相似的数据区段，本文引入互近似熵的概念，利用波形段匹配的方法将V_pcc与I_pcc相似的波形数据段保留，不相似的数据段剔除，以实现排除背景谐波干扰的目的。

针对传统的熵存在需要大量采样数据、对噪声敏感和不易收敛等问题，StevenM.Pincus从衡量时间序列复杂性的角度，于20世纪90年代提出了近似熵(ApEn)。互近似熵是近似熵概念的拓展。令两个不同的时间序列分别为i(t)和j(t)，规定一个长度为m的窗口，如公式(10)和(11)所示。分别对i(t)和j(t)构造N-m+1个m维矢量X_p,X_q，其中

X_p＝[i(p),…,i(p+m-1)]

p＝1,2,…,N-m+1 (10)

X_q＝[j(q),…,j(q+m-1)]

q＝1,2,…,N-m+1 (11)

使用矢量的∞—范数描述Xp与X_q之间的距离d(X_p,X_q)＝||X_p-X_q||_∞

给定相似容限r，其表达式为：

r＝0.2×COV(X_P,X_q) (12)

对每一个p值统计X_p和所有X_q(q＝1,2,…,N-m+1)的矢量距离小于r的个数N_p,m,r，并计算N_p,m,r与总的矢量个数(N-m+1)的比值C_p,m,r：

由公式(13)可知，小于r的矢量个数越多，就意味着两数据段越相近。则两数据段的互相关程度为：

增加窗口长度至m+1，重复上述步骤(10)～(14)的运算过程，得到Φ_m+1,r。

计算得到与m，r相关的互近似熵值：

CAE(m,r)＝Φ_m,r-Φ_m+1,r (15)

在本发明中利用互近似熵值来衡量V_pcc与I_pcc波形的相似度，互近似熵值越小，说明比较数据区段内的波形越相似，也就意味着背景谐波的波动越小。以此阈值标准来对谐波数据进行筛选。

实施例3：M估计稳健回归法原理；传统线性回归法是用最小二乘法进行计算，其原理是使得残差的平方和最小。原始数据中存在异常值时，最小二乘法计算便会迁就远端数据，使得计算误差增大，因此传统线性回归法缺乏稳健性。本文利用M估计稳健回归来消除异常值对方法的影响[18]。以下以一元线性回归计算为例。已知实测数据(x_i,y_i)(i＝1,2,…,n)，设线性关系式为y＝ax+b，则可利用最小二乘法求取方程系数a和b。

由上述描述可知，残差e_i的目标函数可表示为：

传统的最小二乘法是估计使得随机误差e_i的值最小。为避免异常值对回归计算的影响，在M估计稳健回归中，对各点施加不同的权重：对残差小的点施加较大的权重，而对残差大的点加上较小的权重。则M估计稳健回归的形式为：

权重采用Huber法来计算：

式中：c为常数，其值一般取1.345，u_i为标准化的残差指标，u_i＝e_i/s，而s是残差尺度：

s＝Med|e_i-Med(e₁,...,e_n)|/0.6745 (17)

式中，Med表示一组数据的中位数。方法的具体计算过程如下：

(1)利用最小二乘法计算出相应的回归系数a_i(a₁…a_n)，b_i(b₁…b_n)，并计算相应的误差量e_i。

(2)根据计算残差尺度s，计算相应的权重值ω

(3)将权重计算结果ω_i(ω₁…ω_n)代入，求得新的回归计算结果。

比较前后两次计算的回归参数，若a_i,b_i全部满足阈值约束则保留，约束表达式为：

(一般δ₀取值为10-⁵)，否则继续迭代计算，直到满足为止。该迭代步骤保证了系统谐波阻抗回归计算的稳定性，避免了异常值的干扰，确保了计算结果的准确性。

实施例4：如图3所示，本发明提出一种基于互近似熵数据筛选的谐波源责任划分方法，包括以下四个步骤：

步骤一：采集V_pcc，I_pcc谐波数据；为避免数据幅值、量纲对互近似熵计算结果的影响，将谐波数据进行归一化，归一化的计算公式如下所示：

x_N为归一化结果，x_min为全体数据中的最小值，x_max为全体数据中的最大值。

步骤二：将归一化后的谐波数据划分为n个数据区段，每个数据区段包含L个数据点，为避免截断有效数据，L的选择不宜过大，文中L＝10。

步骤三：计算每个区段的CAE值，并将其与设定阈值比较。满足条件的数据保留，反之则舍弃。

步骤四：将所有区段数据筛选完成后即可利用M估计稳健回归法计算系统谐波阻抗，从而得到谐波源责任划分结果。

算例效果评估:以图4中仿真电路对算法的有效性和优越性进行验证。

电路参数如下表所示，图4中系统电源为10kV，其中背景谐波为谐波源I0，系统谐波阻抗为Zs，负荷为Zc，用户侧谐波源为Ic。

具体参数如下表所示：

表1仿真电路参数

系统仿真时间设置为40s，每0.02s采集一组5次谐波电压与电流数据作为样本。PCC点5次谐波电压波形及5次谐波电流波形如下图5所示：

为不遗漏相似波形段的同时，保证有效数据不被截断，本实施例4的实验中互近似熵计算的参数设置为L＝10，m＝2。每10个数据点为一个数据区段，将仿真数据分为200组数据区段，并将其中受背景谐波影响较小的区段筛选出来。小于CAE阈值的数据段保留，反之则舍弃。本实施例4采用枚举法对CAE阈值参数的选择进行了分析，互近似熵的阈值计算选取范围设置为0.03～0.1。图6为CAE阈值变化对应的系统谐波阻抗回归误差。其中误差的计算方法为：

其中E表示估计值，T表示理论值。

由图7中所示，当CAE阈值大于0.06时，误差结果处于0.15左右，当CAE阈值小于0.06时，误差结果急剧下降。这意味着所提算法筛选掉了一些对结果有干扰的数据段，当阈值取为0.045时，谐波阻抗计算的误差可以稳定在0.026左右。因此，CAE阈值选为0.045。

以经过筛选后的数据为基准进行M估计稳健回归计算。数据筛选结果如图8所示，由于背景谐波波动较为严重，筛选前的谐波数据点杂乱无序。经过筛选后，原有的2000个谐波数据点仅保留140个有效数据点(由于随着采集时间的延长，采集点数的增加，有效数据点也会随之增加，因此实际情况中不存在有效数据点不足的问题)。由图9可知，有效数据点呈现出了清晰的线性关系，证明了本发明所提方法的有效性。

本发明分别利用随机矢量法(方法一)，主导波动量法(方法二)，以及本文所提互近似熵与最小二乘线性回归估计组合方法(方法三)，与本文计算结果作对比。结果如表2所示：

表2谐波阻抗计算结果

由表中结果可知，方法一和方法二本身虽然考虑到了背景谐波的影响，分别采用随机矢量的独立性和主导波动量的方法规避背景谐波的影响，但最后谐波阻抗计算的结果依然与理论值偏差较大。方法三与本文所提方法均采用CAE算法将受背景谐波干扰严重的点进行了舍弃，提高了谐波阻抗计算的准确率。方法三利用传统的最小二乘回归法计算了谐波阻抗，其计算准确性要略逊于本文所提方法。但最终误差都控制在了5％以内。

系统谐波阻抗的估算结果直接决定了谐波责任划分的准确性，本发明对每个数据点责任划分结果取平均值，则上述四种方法谐波的谐波责任划分结果如下表3所示，

表3谐波责任划分结果

本发明的创新点在于，提出了一种基于CAE数据筛选的谐波源责任划分方法，首先利用CAE对谐波电压、电流数据进行了预处理；CAE数据筛选方法有效规避了背景谐波波动对谐波阻抗计算的影响，提高了谐波责任划分的准确性。然后通过M估计稳健回归方法计算得到了谐波阻抗值，M估计稳健回归方法有效提高了回归计算的稳定性。通过该方法与CAE数据筛选法的有效结合，谐波阻抗计算准确率可达到97％以上。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。

Claims (10)

1.一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，包括如下步骤：

步骤SS1：采集谐波数据电压U_pcc和谐波数据电流I_pcc，将所述谐波数据电压U_pcc和谐波数据电流I_pcc进行归一化处理；

步骤SS2：将归一化后的谐波数据划分为n个数据区段，每个数据区段包含L个数据点；

步骤SS3：计算每个数据区段的互近似熵值CAE，并将所述互近似熵值CAE与设定阈值比较，对满足设定阈值条件的数据进行保留，反之则进行舍弃；

步骤SS4：将所有数据区段的数据筛选完成后，利用M估计稳健回归法计算系统谐波阻抗，得到谐波源责任划分结果。

2.根据权利要求1所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS1具体包括：

设I_ch为用户侧等效谐波源，I_sh为系统侧等效谐波源，Z_ch为用户侧等效谐波阻抗，Z_sh为系统侧等效谐波阻抗，h表示谐波次数；I_pcc-s为系统侧等效谐波源I_s贡献谐波电流，V_s为系统侧贡献谐波电压；I_pcc-c为用户侧等效谐波源I_c贡献谐波电流，V_c为用户侧贡献谐波电压；各参数表达如公式(1)～(4)所示：

I_pcc＝I_pcc-s+I_pcc-c (3)

V_pcc＝V_s+V_c (4)

由于电力系统中Z_c远大于Z_s，则上述关系式等效改写为：

联立式(5)和(6)得到:

V_pcc＝Z_s.I_pcc+V_s (7)

上述公式(1)～(7)中的元素皆为方向矢量，；

当背景谐波波动很小时，Vs近似看作常数，因此由式(7)知，当系统谐波阻抗不变的情况下，V_pcc与I_pcc呈线性关系；V_pcc的波形趋势应与I_pcc的波形趋势一致；

式(7)模型的表达式如下：

V_pcc＝|Z_s|cosθ.I_pcc+|V_s|cosθ₁ (8)

其中θ为Z_s·I_pcc与V_pcc的夹角，θ₁为V_s与V_pcc的夹角；则对于用户侧谐波责任H_c的计算由公式(9)表示：

当背景谐波保持平稳时，将式(8)中的|Z_s|cosθ和|V_s|cosθ₁近似看作常系数。

3.根据权利要求1所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS1中的归一化处理的计算公式如下：

其中，x_N为归一化结果，x_min为全体数据中的最小值，x_max为全体数据中的最大值。

4.根据权利要求1所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS2中的L取10。

5.根据权利要求1所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS3中的所述互近似熵值CAE的计算方法如下：

设两个不同的时间序列分别为i(t)和j(t)，规定一个长度为m的窗口，如公式(10)和(11)所示，分别对i(t)和j(t)构造N-m+1个m维矢量X_p,和X_q，其中：

X_p＝[i(p),…,i(p+m-1)]

p＝1,2,…,N-m+1 (10)

X_q＝[j(q),…,j(q+m-1)]

q＝1,2,…,N-m+1 (11)

使用矢量的∞—范数描述X_p与X_q之间的距离d(X_p,X_q)＝||X_p-X_q||_∞；

给定相似容限r，r的表达式为：

r＝0.2×COV(X_P,X_q) (12)

对每一个p值统计X_p和所有X_q(q＝1,2,…,N-m+1)的矢量距离小于r的个数N_p,m,r，并计算N_p,m,r与总的矢量个数(N-m+1)的比值C_p,m,r：

由公式(13)可知，小于r的矢量个数越多，就意味着两数据段越相近；则两数据区段的互相关程度为：

增加窗口长度至m+1，重复上述公式(10)～公式(14)的运算过程，得到Φ_m+1,r；

最后，计算得到与m，r相关的互近似熵值：

CAE(m,r)＝Φ_m,r-Φ_m+1,r (15)

6.根据权利要求1所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS4具体包括：

步骤SS41：利用最小二乘法计算出相应的回归系数a_i(a₁…a_n)，b_i(b₁…b_n)，并计算相应的误差量e_i；

步骤SS42：根据计算残差尺度s，计算相应的权重值ω；

步骤SS43：将权重计算结果ω_i(ω₁…ω_n)代入M估计稳健回归方程，求得新的回归计算结果；

比较前后两次计算的回归参数，若a_i,b_i全部满足阈值约束则保留，约束表达式为：

否则继续迭代计算，直到满足为止。

7.根据权利要求6所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述M估计稳健回归方程为：

其中，a和b均为回归系数，ω_i为权重计算结果。

8.根据权利要求6所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述残差尺度s的技术公式如下：

s＝Med|e_i-Med(e₁,...,e_n)|/0.6745 (17)

其中，Med表示一组数据的中位数。

9.根据权利要求6所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述权重值ω采用Huber法来计算：

式中：c为常数，其值一般取1.345，u_i为标准化的残差指标，u_i＝e_i/s。

10.根据权利要求6所述的一种基于互近似熵数据筛选的谐波源责任划分方法，其特征在于，所述步骤SS41具体包括：已知实测数据(x_i,y_i)(i＝1,2,…,n)，设线性关系式为y＝ax+b，则利用最小二乘法求取方程系数a和b；

残差e_i的目标函数表示为：

其中，a和b均为回归系数。

Images