基于DBSCAN聚类和区间回归的多谐波责任划分

01监测数据采集与谐波责任划分分析1.1

谐波监测数据采集与样本集构建电能质量监测装置监测点如图1所示，监测装置一般部署在10kV母线公共连接点上，数据采样间隔为3min，可获得母线谐波电压数据和多条馈线的谐波电流数据。图1

电能质量监测系统与数据采集示意

Fig.1

Schematicdiagramofpowerqualitymonitoringsystemanddataacquisition

目前电网公司所使用的谐波监测装置输出的测量数据一般为一段监测周期内的统计值，如最大值、最小值、平均值和95%概率大值。注入电力系统关注节点的谐波监测数据可用该监测周期内的最小值与最大值组成的区间来表征，取其构成的区间来表示谐波电流的范围，设存在区间谐波电流样本集X为式中：[xij]为j次谐波在第i个监测时段内监测数据所属的区间，[xij]=[xij(min),xij(max)]，xij(max)、xij(min)分别为监测数据的最大值、最小值。1.2

多谐波源谐波责任当电力系统中存在多个谐波源分散分布时，任一关注节点上的谐波电压畸变都是由所有谐波源注入谐波电流而引起的共同结果。对某一关注节点而言，为对谐波实现科学准确的管控，须定量分析电力系统中每个谐波源在该节点谐波电压畸变中应准确承担的谐波责任。假设该电力系统中含有n个谐波源，其系统谐波电压方程为式中：U为节点的谐波电压向量；I为谐波电流向量；Z为谐波转移阻抗矩阵。展开可表示为式中：第m行对应关注节点m；其谐波电压Um等于系统中n个谐波源共同在该节点贡献之和，即式中：Zms为关注节点m和谐波源s(s=1,2,···,n)之间的谐波阻抗（当s≠m时为互阻抗，当s=m时为自阻抗）；Ums为谐波源s（s=1,2,···,n）对关注节点m的谐波电压贡献值。设谐波源s在关注节点m单独作用时产生的h次谐波电压为Ush，除谐波源s以外的所有谐波源在关注节点m共同作用时产生的h次谐波电压为U0h，U0h也可理解为谐波源s在关注节点m作用时的h次背景谐波电压。根据叠加原理，则式（3）可表示为式中：Zsh为谐波源s与关注节点m之间的h次谐波转移阻抗；Ish为谐波源s注入关注节点m的h次谐波电流。Umh、Ush和U0h满足如图2所示的向量关系，图2中θsm为Ush和Umh的夹角，θ0m为U0h和Umh的夹角，谐波源s在关注节点m上产生的h次谐波电压责任可由Ush在Umh方向上的投影长度占Umh=|Umh|的百分比表示。图2

谐波电压向量关系

Fig.2

Harmonicvoltagevectorrelationship以图2的Umh为参考向量，根据式（4）其幅值Umh为式中：Ush=|Ush|；U0h=|U0h|；Zsh=|Zsh|；Ish=|Ish|。设谐波源s的谐波电流系数Zshcos⁡θsm为α1，背景谐波电压U0h在Umh方向上的投影U0hcos⁡θ0m为α0，即谐波源s对关注节点的h次谐波电压贡献可通过谐波电流幅值与其系数的乘积表示，则式（5）可表示为Umh和Ish由测量获取，由线性回归计算得到α1和α0的值。因此，依照谐波责任划分定义，可得谐波源s对关注节点m的h次谐波电压责任Fsh为1.3

区间谐波责任如果α1和α0的值在谐波责任划分时间段内保持稳定，则可由监测数据Umh和Ish回归计算得到α1和α0的值，进而由式（7）计算得到谐波源s应承担的谐波责任。但在实际电力系统中，常发生负荷投切、无功补偿装置切换等情况，导致系统谐波阻抗发生变化，且系统背景谐波电压往往存在剧烈波动，即α1和α0处于动态变化中。如图3所示，当谐波阻抗Zsh和背景谐波电压U0h发生变化时将对谐波责任划分造成较大影响。图3

谐波阻抗变化与谐波电压波动示意

Fig.3

Schematicdiagramofharmonicimpedancevariationandharmonicvoltagefluctuationα1和α0虽然处于动态变化中，但其变化是有界的，α1和α0在最小值与最大值区间内变化，区间变量[α1]和[α0]分别为在任意一个监测周期内，监测数据谐波电压Umh和谐波电流Ish介于最小值和最大值之间，同理可构造区间谐波电压数据[Umh]和区间谐波电流数据[Ish]分别为结合式（8）和式（9），可将式（6）转变为区间划分形式，即在划分谐波责任的关注时段内，谐波阻抗Zsh和背景谐波电压U0h可能发生变化，若将关注时段的数据集划分为k个时间段，且每个时间段内的谐波源运行保持相对平稳，则将每个时间段称为同场景时段。02背景谐波变化下的谐波责任划分当谐波监测装置数据采集周期较长时，负荷投切、无功补偿装置切换和谐波电压剧烈波动等情况将对如何准确划分谐波源s应承担的谐波责任产生较大影响，不同时段谐波源s的应承担谐波责任将差别较大。因此，本文提出以[α1]和[α0]为划分依据将谐波样本集进行运行场景划分，将同场景时段监测数据集划分为同一簇，寻找最优的样本划分。2.1

基于DBSCAN聚类的场景划分令样本集E表示为{[X],[Y]}，利用区间动态聚类算法将样本集E划分为k个同场景时段。利用监测到的谐波样本数据集E进行上述的区间线性回归估算，首先需要将该监测数据划分为不同的场景，进行模型参数估算。因为电力系统的波动性和复杂性且受负荷暂态过程的影响，监测数据并不是理想中整齐排列的不同区间，存在较多干扰数据。这就要求聚类算法具有较强的抗干扰能力和较好的聚类效果。相比于经典的k-means，DBSCAN具有更好的抗干扰能力。因此，本文利用DBSCAN进行动态聚类区间划分。DBSCAN在给定对象半径Eps和最小数目min

Pts的条件下，通过判断样本集E的数据点是否满足条件，将[X]和[Y]数据划分为不同的线性数据簇，有效地分离出不同谐波阻抗和背景谐波电压的数据簇。DBSCAN算法流程为：1）将样本集E的全部数据标记为待处理状态；2）以样本集E中的每个数据P所处位置为圆心，Eps为半径作圆，得到包含不同对象数的圆域；3）判断数据P的圆域所包含的对象数是否小于最小数目min

Pts，若小于则将数据P标记为边界点或噪声点，若大于等于则进入步骤4）；4）将数据P标记为核心数据点，并建立新簇Ci，将数据P圆域内未归类为任意一簇的数据点并入簇Ci；5）检查样本集E中的所有数据点是否已全部并入某簇或标记为边界点或噪声点，若否则回到步骤2），反复循环以将每个数据点都进行处理，若是则输出结果。误差平方和（SSE）是考虑欧几里得距离的聚类效果评价指标，即计算所有子类到对应类簇聚类中心的欧氏距离dSSE为式中：L为聚类分区数；ci为类簇x的聚类中心；d(ci,x)为向量间的欧式距离。聚类效果与L密切相关，随着L从1动态向真实同场景时段数k逐步增加，dSSE将随之减小；当L值越过k值再增大时，dSSE将趋于平缓并可能略有增大，将拐点处作为最佳聚类处。2.2

滑窗动态相关性分析聚类完成后，监测数据集被划分为不同的簇，但只完成了区间划分，并未剔除各簇内的干扰数据。监测数据的背景谐波电压波动较大时，簇内部分区间数据与其他样本数据分别对应的背景谐波电压可能相差较大，如不经处理直接进行区间线性回归估算将破坏Umh与Ish的线性关系。本文采用滑窗动态相关性分析方法表征母线谐波电压和馈线谐波电流的相关性，对于谐波样本集X和Y，其相关性可以用系数r(x,y)表示为式中：xi和yi分别为样本集X和Y中的第i个元素；分别为X和Y中全部元素的均值，r(x,y)∈[−1,1]。在长时段监测周期内，各馈线用户用电行为对关注节点m处的谐波电压影响并非为恒定线性关系，易受干扰因素破坏线性关系，因此需分析不同时刻数据间的相关程度，筛选相关性较高的数据，以便精确划分馈线谐波责任。本文通过计算给定窗宽下序列相关系数，滑动窗口获取全时段内动态相关性数列作为筛选依据。设给定窗宽为C，动态相关性数列为其数列所含元素计算方法为式中：1≤i≤N–C，N为样本集采样数；x(i,i+C)和y(i,i+C)分别为样本集X和Y中第i个元素到第i+C个元素。当的值趋向于1或者–1时，表明Umh与Ish数据段具有良好的线性关系。在实际运用中通过反复观察筛选后以Ish为横坐标以Umh为纵坐标的散点图是否具有相当好的线性分布来选取的阈值，经大量比对，本文设定窗宽C为5，的阈值为0.85。筛选出线性关系良好的区间数据[X]和[Y]，继而对组内样本数据进行区间线性回归估算与样本分派。2.3

区间回归分析区间回归分析通过构建区间自变量和区间因变量的函数模型进而借助区间自变量预测未知的区间因变量。传统的回归分析大都基于点数据开展，针对模糊区间数据本文采用能自动选择最优参考点的算法PM来进行区间线性回归以估算区间方程参数，并在实测数据上验证该方法所构建的区间回归模型的参数估算准确性。样本集X的任意区间[xij](i=1,2,···,n;j=1,2,···,p)上的任一点xij可表示如式（14）所示的形式，其中，0≤γj≤1，确定一个γj值就可以确定区间[xij]上一点xij的值，即对于任意同场景时段区间，构建该区间的上下限回归模型为式中：αj(max)、αj(min)分别为上下限回归系数；εi(max)、εi(min)分别为上下限误差。将式（14）代入式（15），可得如区间线性回归上下限模型为计算区间线性回归下限模型的回归系数，并定义系数δj(min)和ωj(min)，使其满足式（17）所示的函数关系。将式（17）代入式（16），区间线性回归下限模型可进一步表示为式（18）可表示为矩阵形式，即根据最小二乘法，可估算α(min)和γj(min)为同理，对式（16）所示的区间线性回归上限模型运用式（17）~（20）的分析方法，易得到如式（21）所示区间线性回归上限模型的α(max)和γj(max)为综上所述，样本yi的预测区间为[X(min)α(min),X(max)α(max)]。由以上分析可知，γj无须提前设定，在区间线性回归模型的参数估算过程中自动计算了γj值。传统的区间回归分析方法需要人为设置参数点，PM算法能够根据区间样本集的特征，在区间自变量和因变量上自动选择线性程度最高的参数点xij。2.4

样本归属根据式（20）（21）分别计算各簇区间线性回归下限模型的回归系数α(min)和区间线性回归上限模型的回归系数α(max)，继而对所有样本数据根据式（20）（21）计算第i个区间样本在第a(a=1,2,···,k)个同场景时段的预测区间为[Xa(min)α(min),Xa(max)α(max)]，i=1,2,···,n。在得到预测区间[Xa(min)α(min),Xa(max)α(max)]后，可计算第i个样本属于第a个同场景时段的区间估算的误差εia为以获得区间样本归属某簇误差最小为目标，设最佳样本区间划分方案Mk通过如式（23）所示的准则获得，直至第i个区间样本归属某一簇的区间估算误差比归属任一簇的区间估算误差都要小。式中：Oa表示第i个样本归属的第a个同场景时段簇。03谐波责任指标在同场景时段，谐波阻抗和背景谐波电压基本不变，可认为[α1]和[α0]不再为区间值，而为定值。在第a个同场景时段（a=1,2,···,

k），可将式（8）的区间值[α1]和[α0]替换为利用同场景时段的区间监测数据由式（10）回归计算得到区间方程系数α1(a)，进一步由式（7）得到每一个同场景时段的谐波责任为式中：[Fsh(a)]为第a个同场景时段下谐波源s的h次谐波电压责任区间。按上述区间谐波责任划分方法，对监测数据进行场景划分，将区间值[α1]和[α0]分别趋向某一稳定值的时段划分为同场景时段。估算得到各同场景时段的区间方程系数α1(a)(a=1,2,···,

k)，再由式（25）确定谐波责任区间。若关注时段内含有k个同场景时段，对于其中的第a个同场景时段(a=1,2,···,

k)，假设该同场景时段由N(a)个测量周期组成。则式（25）可进一步表示为式中：[Fsh]a为第a个同场景时段下谐波源s能直接用于谐波责任划分的h次谐波电压责任具体值；区间上划线表示取平均值。为了评估关注时段内谐波源s的应承担的谐波责任，累计算出第1至第a个同场景时段的区间谐波责任为本文方法流程如图4所示。图4

谐波责任划分方法流程

Fig.4

Harmonicresponsibilitydivisionmethodflowchart

04实例分析为验证本文方法在实际电力系统谐波责任划分中的有效性，本文以含背景谐波变化的某110kV实际变电站10kV母线及该母线所包含的3条馈线为例进行谐波责任划分。电能质量监测装置监测点如图1所示，监测周期为7天，数据采样间隔为3min，对10kV母线PCC处采集母线谐波电压数据，对3条馈线采集谐波电流数据。在监测周期内，得到母线各次谐波电压和3条馈线各次谐波电流的3360组连续采样数据。由于数据体量较大，为清晰展示并验证本文方法对多次谐波责任划分的有效性，本文列出9次谐波和7次谐波责任计算结果。4.1

9次谐波责任划分为获得最佳聚类效果，采用本文方法进行多次聚类，记录聚类分区数L和每次聚类后的dSSE，dSSE和L对应的关系如图5所示。从图5可看出，聚类分区数L=7是3条曲线共同明显的拐点，此时dSSE下降趋势变得缓慢。因此，将采样周期内采样得到数据集划分为7个同场景时段。图5

dSSE与L关系

Fig.5

RelationshipbetweendSSEandL对划分后的各同场景时段数据集，分别利用其3条馈线9次谐波电流样本和母线9次谐波电压样本计算如式（8）所示的各谐波源的谐波电流区间参数[α1]，结果如表1所示。表1

9次谐波各场景区间参数[α1]计算结果Table1

Calculationresultsoftheparameter[α1]fortheintervalofeachscenarioofthe9thharmonic此后，运用本文方法计算7个场景各336个区间样本的动态9次谐波责任指标。以所有谐波源在母线共同造成的结果为基准值，对3条馈线的动态9次谐波责任计算结果作箱线图，如图6所示。其中，负动态谐波责任值表示该馈线被迫吸收谐波，与母线共同承担了由谐波源造成的后果，是电能质量问题中的受害者。因馈线2并未产生谐波反而被迫吸收了部分谐波，在以母线谐波电压畸变为基准的前提下部分馈线的动态谐波责任指标会大于1。图6

动态9次谐波责任箱线图

Fig.6

Dynamic9thharmonicresponsibilityboxlinediagram

结合图6可知，3条馈线的谐波责任在一定范围内波动，这与3条馈线的运行场景和所含用户的用电行为关系紧密。馈线3的9次谐波责任波动范围较大且应承担主要责任，馈线19次谐波责任较为集中且应负次要责任，馈线2被迫吸收了其他两条馈线产生的谐波且不同时刻吸收程度存在较大差异，其规律反映进行场景划分的必要性。计算各同场景时段的9次谐波责任的标幺值，其结果如表2所示。以场景为观察轴，同一馈线不同场景下的9次谐波责任存在明显波动，不同馈线波动剧烈程度存在差异，但基本在某一值附近波动，3条馈线的9次谐波责任随场景改变而变化的趋势基本一致，以上变化规律与该变电站存在的运行场景变化规律基本一致。表2

各场景下9次谐波责任标幺值的划分Table2

Divisionofresponsibilityp.u.forthe9thharmonicineachscenario由式（27）计算监测周期内的全时段9次总谐波责任，并分析各馈线的9次谐波责任在不同场景下差异，如表3所示。结合表1和表3可知，3条馈线在不同场景下9次谐波责任计算系数α1具有一定的波动，且馈线2波动较大，如未进行场景划分按固定谐波责任计算系数直接进行监测周期内的谐波责任划分，其结果准确度和可信度将存疑。监测周期内3条馈线的最大责任值和最小责任值具有一定差值，且馈线2差值较大，与表1中馈线2的α1值波动最大存在对应关系，其差值说明进行场景划分的必要性，对于谐波责任划分须分情况分场景进行讨论，建立从局部评价到整体评价的谐波评价体系。表3

9次谐波责任标幺值的划分统计分析Table3

Statisticalanalysisofthedivisionofresponsibilityp.u.forthe9thharmonic对监测数据采用概率性谐波分析方法，做出3条馈线9次谐波电流和母线9次谐波电压的概率密度函数曲线，如图7所示。从图7中可看出馈线3的9次谐波电流概率密度曲线和母线的9次谐波电压概率密度曲线形状和走势高度类似，馈线1和馈线2的9次谐波电流概率密度曲线位置和形状高度类似。这与图6和表3可得出的结论“馈线3应承担主要9次谐波责任，馈线1和馈线2在9次绝对谐波责任数值上接近”具有一致性。图7

馈线9次谐波电流和母线9次谐波电压概率密度曲线

Fig.7

Feeder9thharmoniccurrentandbusbar9thharmonicvoltageprobabilitydensitycurve

4.2

7次谐波责任划分同理，按本文4.1节的谐波责任划分方法和流程步骤对7次谐波进行谐波责任划分。计算各“同场景时段”的7次谐波责任，其结果如表4所示。计算监测周期内的全时段7次总谐波责任，并分析各馈线的7次谐波责任在不同场景下差异，如表5所示。表4

各场景下的7次谐波责任标幺值的划分Table4

Divisionoftheresponsibilityp.u.forthe7thharmonicineachscenario表5

7次谐波责任标幺值的划分统计分析Table5

Statisticalanalysisofthedivisionofresponsibilityp.u.forthe7thharmonic从表4和表5可看出，3条馈线的7次谐波责任值和9次谐波存在类似