基于快速動態(tài)時間彎曲和最小覆蓋球的多日負荷曲線聚類方法

劉曉峰，康 進，馬 翔，沃建棟，呂磊炎，吳 浩

（1. 浙江大學電氣工程學院，浙江杭州 310027；2. 國網(wǎng)浙江省電力有限公司金華供電公司，浙江金華 321016）

0 引言

“十四五”以來，以新能源為主體的新型電力系統(tǒng)得到大力發(fā)展，電力負荷作為終端，是新型電力系統(tǒng)建設(shè)的重點，也是源網(wǎng)荷儲的重要環(huán)節(jié)之一。為構(gòu)建清潔、低碳、安全、高效的能源體系，新能源分布式發(fā)電、柔性負荷、電動汽車等的比重大幅上升［1］。該類負荷受天氣、調(diào)控、用戶意愿的影響較大，使得總負荷在多日間的波動較明顯，給負荷曲線聚類研究帶來了難度。此外，隨著用電信息采集裝置的部署與數(shù)據(jù)采集技術(shù)的發(fā)展，電力系統(tǒng)負荷數(shù)據(jù)資源得到了快速積累［2］，給多日間負荷波動的研究帶來了可能性。因此，有必要且有條件進行多日負荷曲線聚類研究，探索負荷在多日間的波動特性，進而提升負荷調(diào)控能力以及深化電力體制改革［3-4］。

負荷曲線聚類是電力系統(tǒng)負荷特性分析的重要手段之一［5］。負荷曲線聚類研究通常包括負荷特征提取及降維、負荷間的相似度定義、聚類算法這3 個方面。近年來，眾多研究者針對上述3 個方面，從算法效率、聚類準確率等多個角度對傳統(tǒng)負荷曲線聚類研究進行了改進及深化。

在負荷特征提取及降維方面，現(xiàn)有研究主要從2 個角度開展。一部分研究針對特定的聚類方法，選取合適的負荷特征向量：文獻［6］通過離散小波變換提取負荷曲線的頻域特征，并與模糊K-modes 聚類方法結(jié)合得到聚類結(jié)果；文獻［7］提出一種譜聚類特征向量的優(yōu)化選取方法，提高了譜聚類算法的效率。另一部分研究采用傳統(tǒng)降維方法對負荷曲線進行降維：文獻［8］利用奇異值分解對日負荷曲線進行降維，從而提取特征；文獻［9］比較自組織映射、主成分分析等降維技術(shù)在電力系統(tǒng)負荷曲線聚類中的應(yīng)用效果。

在負荷間的相似度定義方面，現(xiàn)有研究主要從3 個角度開展，分別為采用歐氏距離等點對點的距離作為相似度定義［10］，采用動態(tài)時間彎曲等考慮負荷曲線時間滯后特性的相似度定義［11］，以及利用前2種方式綜合衡量負荷曲線間的相似度。

在聚類算法方面，常見算法有層次聚類算法、劃分聚類算法、密度聚類算法、模型聚類算法、譜聚類算法等［5］。由于譜聚類算法能夠劃分任意形狀的數(shù)據(jù)集，且對稀疏矩陣的處理效果優(yōu)于其他算法，因此本文采用該算法。

從聚類對象上而言，上述研究大多集中于單日負荷曲線。實際上多日間的負荷曲線可能差異較大，導致算法在不同日間應(yīng)用的穩(wěn)定性較差，因此有必要開展多日負荷曲線的降維與聚類研究。目前多日負荷曲線聚類方法主要是統(tǒng)計多日同一時刻的負荷分布情況以表征負荷的多日波動情況［12-13］。由于負荷曲線屬于時間序列，不同日的相同用電行為未必發(fā)生在同一時刻，因此應(yīng)考慮日負荷曲線的時間滯后特性。

針對現(xiàn)有研究的不足，本文提出一種基于快速動態(tài)時間彎曲FDTW（Fast Dynamic Time Warping）與最小覆蓋球的多日負荷曲線聚類方法。首先，基于FDTW 與多維尺度縮放MDS（Multi-Dimensional Scaling）提取負荷曲線的時間滯后特性，并將負荷曲線降至3 維；然后，為每個負荷迭代尋找最小覆蓋球，并基于最小覆蓋球定義負荷間的相似度計算公式；最后，將相似度矩陣應(yīng)用于基于k-最近鄰法的譜聚類算法。算例結(jié)果表明，本文所提方法在準確度和魯棒性上較傳統(tǒng)方法有一定優(yōu)勢。

1 基于FDTW和最小覆蓋球的多日負荷曲線聚類

1.1 多日負荷曲線的特點

為了說明多日負荷曲線的特點，圖1 給出了某地區(qū)一臺10 kV 專變2 周（2019 年8 月1 日至2019 年8月14日，其中8月1日為周四）內(nèi)的日負荷曲線圖。由圖1（a）所示的歸一化功率曲線可見：普通日負荷曲線、曲線3 和曲線4 均為雙峰曲線，其形狀大致相似，且這些曲線數(shù)量占85%以上，因此可以確定該專變下主要為雙峰型負荷；曲線1和曲線2分別具有晚高峰和錯峰曲線的特點，與雙峰曲線相似度極低，經(jīng)調(diào)查，曲線1 和曲線2 分別出現(xiàn)在8 月9 日（周五）和8月10日（周六），因此，異常曲線不完全是由節(jié)假日、非工作日產(chǎn)生的。由圖1（b）所示的實際功率曲線可見，曲線1 和曲線2 的實際功率較高，可排除數(shù)據(jù)采集、日內(nèi)負荷波動（此處負荷波動表示日內(nèi)負荷曲線的起伏，以日內(nèi)負荷波動命名，以與多日間的負荷波動相區(qū)別）等異常的影響，因此，推測此類曲線的出現(xiàn)可能是該專變下用戶加班、調(diào)休等原因所導致的。顯然，如果采用單日負荷曲線聚類的方法提取到曲線1或曲線2所對應(yīng)的負荷，則該負荷的聚類結(jié)果就會錯誤，因此，有必要考慮多日負荷曲線間的波動特性，研究多日負荷曲線的聚類方法。

圖1 某10 kV專變2周內(nèi)的日負荷曲線Fig.1 Daily load curves of a 10 kV specialized transformer within two weeks

曲線3 和曲線4 分別出現(xiàn)在8 月3 日（周六）和8月6 日（周二），二者均為雙峰曲線，但是曲線4 相對于曲線3 有約2 h 的時間滯后，因此，在多日負荷曲線的聚類過程中需考慮負荷曲線的時間滯后特性。

目前，從電網(wǎng)獲取有功功率數(shù)據(jù)的采樣間隔為15 min，相應(yīng)地，日負荷曲線的數(shù)據(jù)共有96 點，即數(shù)據(jù)維度為96 維?？紤]到多日負荷曲線的時間通常為2 周及以上，直接對負荷曲線進行聚類時數(shù)據(jù)維度超過1000維，可能導致聚類的準確度和魯棒性較低，因此，有必要在多日負荷曲線聚類的過程中對其進行降維。

1.2 算法思路

針對上述多日負荷曲線的特點，本文提出一種基于FDTW 和最小覆蓋球的多日負荷曲線聚類方法。該方法主要包含負荷曲線降維、基于最小覆蓋球的相似度定義、譜聚類這3個步驟。

1）負荷曲線降維。將多日負荷曲線按日分割，考慮負荷曲線的時間滯后特性，采用FDTW 形成日負荷曲線距離矩陣，并采用MDS 算法將該距離矩陣降至3維。

2）基于最小覆蓋球的相似度定義。在3 維空間中為每個負荷迭代尋找最小覆蓋球，并以3 種方式度量球間的相似度，以反映多日間的負荷波動大小，并減少負荷異常日的影響。

3）譜聚類。采用k-最近鄰法根據(jù)球間相似度矩陣構(gòu)造相似度圖，并將其應(yīng)用于譜聚類算法，得到多日負荷曲線的聚類結(jié)果。

所提方法的流程圖如附錄A圖A1所示。

2 基于FDTW和MDS的負荷曲線降維

2.1 降維思路

設(shè)多日負荷曲線的樣本數(shù)為n（即為10 kV 專變和公變數(shù)），持續(xù)時間為m日，每日采樣點為96 點，將負荷數(shù)據(jù)構(gòu)成的矩陣記為X∈Rn×96m。首先，對數(shù)據(jù)矩陣X進行曲線平滑及max-min 歸一化處理；其次，考慮負荷曲線的時間滯后特性，將n個負荷按日分割，得到重構(gòu)的負荷數(shù)據(jù)矩陣Y∈Rnm×96，并將矩陣Y中有數(shù)據(jù)采集異常的日負荷曲線剔除，得到Y(jié)1∈Rs×96(s≤nm)；然后，利用FDTW 計算日負荷曲線之間的距離，得到日負荷曲線距離矩陣D∈Rs×s；最后，利用MDS對矩陣D進行降維，得到3維負荷數(shù)據(jù)矩陣Z∈Rs×3。

2.2 FDTW

動態(tài)時間彎曲DTW（Dynamic Time Warping）在計算2 個時間序列的距離時自動扭曲時間序列（即點與點間的匹配在時間軸上進行局部的提前或滯后，但不打亂時間序列內(nèi)部各點的先后次序），使得2個序列的形態(tài)盡可能一致，得到最短距離［11］。

給定長度分別為x、y的時間序列，首先，計算2個序列各點之間的距離，形成距離矩陣d∈Rx×y，然后，設(shè)定起始條件L（1，1）=d（1，1），遞推規(guī)則如式（1）所示，構(gòu)造路徑矩陣L∈Rx×y，最終得到的L（x，y）即為2個序列的DTW距離。

在利用DTW 計算負荷曲線間的距離時算法復雜度較高，且可能會造成時間尺度過度扭曲，因此，本文采用FDTW［14］限制時間扭曲尺度，即算法搜索范圍，使時間復雜度由O(t2)降為O(τt)，其中t=max(x，y)，τ為限制算法搜索的范圍。

由于圖1（a）中曲線3 和曲線4 間的時間滯后約為2 h，因此，本文取τ=2 h，即8 個點。FDTW 提取負荷曲線時間滯后特性后形成距離矩陣D∈Rs×s，其第i行第j列的元素D( )i，j為第i條和第j條日負荷曲線間的FDTW距離。

2.3 MDS

MDS是一種能在降維過程中保持數(shù)據(jù)樣本間距離或者相似性的方法［15］。本文將MDS應(yīng)用于FDTW所得的距離矩陣D∈Rs×s，將其降至3維矩陣Z∈Rs×3，使得任意2 條日負荷曲線在3 維空間中的歐氏距離盡可能等于矩陣D所反映的距離，即：

式中：zi和zj分別為第i條和第j條日負荷曲線在3 維空間中的表示。在s個已知條件下求取3 維未知坐標，由于s?3，因此式（2）有解，本文采用梯度下降法對該式進行求解。

3 基于最小覆蓋球的相似度定義

3.1 最小覆蓋球的生成

本文的最小覆蓋球是指包含3 維空間中給定點的最小球［16］，其生成的具體流程如下。

1）以離均值點最近的3 個點構(gòu)造初始最小覆蓋球。

2）計算球外各點與球心的歐氏距離，對球外負荷點按照相應(yīng)歐氏距離從小到大的順序重新排序。

3）計算將球外距球心歐氏距離最小的3 個點包圍進來后的最小覆蓋球。

4）重復步驟2）、3），直至60%的點被包含進球內(nèi)。

5）重復步驟2），并計算將球外距球心歐氏距離最小的1個點包圍進來后的最小覆蓋球。

6）計算球內(nèi)負荷點數(shù)量和成本指標Icost，Icost為球體積增大的比例與球內(nèi)點數(shù)量增加的比例的比值，反映了在最小覆蓋球中添加某點所需要付出的成本，其表達式為：

式中：r-、n-和r+、n+分別為添加某點前、后的最小覆蓋球半徑和球內(nèi)負荷點數(shù)。

7）當Icost＞2 時，認定該點為異常負荷點，輸出包圍該點前的最小覆蓋球，迭代結(jié)束，否則，重復步驟5）、6）。

8）得到最小覆蓋球矩陣B∈Rn×4，其中前3 列為球心坐標，第4列為球半徑。

最小覆蓋球?qū)嵸|(zhì)上體現(xiàn)了剔除負荷異常日后各負荷點的最優(yōu)包絡(luò)，其球心及半徑反映了該負荷的日負荷曲線特征及波動特性。

3.2 相似度定義

2 個最小覆蓋球之間存在相離、相交、相含3 種情形，因此，須對3 種情形分別定義相似度計算方式。為此，先明確3種情形下的相似度定義要求。

1）兩球相離應(yīng)表示2個負荷的相似度為0。

2）兩球相交時的相似度應(yīng)與兩球相交的方向無關(guān)，而僅與相交部分的體積占比有關(guān)。

3）兩球相含時的相似度不僅與被含的體積相關(guān)，還與兩球半徑之差和球心間距離相關(guān)。具體地：當兩球半徑之差一定時，球心間距離越小則相似度越高；當球心間距離一定時，兩球半徑之差越小則相似度越高。

根據(jù)上述3 個要求，對給定半徑分別為r1和r2、體積分別為V1和V2的2個球，可定義相似度如下：

式中：Vol為兩球相交部分的體積，其計算公式如式（5）所示；rmin=min(r1，r2)；d12為球心間距離。

式中：h1、h2的含義如圖2所示。

圖2 兩球相交示意圖Fig.2 Schematic diagram for intersection of two balls

顯然，該相似度定義滿足上述3 個要求。當兩球相離時相似度最小，此時εsi=0；當兩球重合時，r1=r2且d12=0，此時εsi達到最大值1；由于式（6）、（7）恒成立，因此，εsi處于［0，1］之間，其值越大，則2 個負荷越相似。

式中：rmax=max(r1，r2)。式（6）反映了大小球相交的體積小于等于小球的體積，小于等于大小球體積的幾何平均值；式（7）反映了大小球半徑之和小于等于2 倍的大球半徑，小球半徑小于等于大小球半徑的幾何平均值。

下面論證該相似度定義滿足兩球相離、相交、相含時的邊界條件。兩球外切為兩球相離、相交的邊界條件，此時Vol=0，式（4）中相離和相交下的計算結(jié)果相同，因此邊界條件成立；兩球內(nèi)切為兩球相交、相含的邊界條件，此時Vol=4πr3min/3，r1+r2-d12=2rmin，式（4）中相交和相含下的計算結(jié)果相同，因此邊界條件成立；rmin=0 且d12=0 為兩球相離、相含的邊界條件，此時式（4）中相交和相含下的計算結(jié)果均為0，邊界條件成立。綜上，該相似度的定義具有一致性。考慮到本文的譜聚類算法將距離矩陣作為輸入，距離越小，則2個負荷越相似，因此定義ds=1/εsi，將相似度轉(zhuǎn)化為距離，以此得到對稱距離矩陣Ds∈Rn×n。

4 譜聚類及聚類評價指標

4.1 基于k-最近鄰法的譜聚類

譜聚類算法是一種基于譜圖理論的聚類算法，它將數(shù)據(jù)點映射為無向圖，將點與點之間的相似性表示為帶權(quán)重的邊，從而將聚類問題轉(zhuǎn)化為圖劃分問題［10］。譜聚類的關(guān)鍵是相似度圖的構(gòu)造，本文采用k-最近鄰法來構(gòu)造相似度圖。k-最近鄰法構(gòu)造相似度圖的思路為：對于包含n個點的數(shù)據(jù)集，令其中每個點xi對應(yīng)相似度圖上的一個節(jié)點vi，若vj屬于vi的前k個最近鄰節(jié)點集，則稱vj是vi的k近鄰，并連接vi和vj，其權(quán)重為距離矩陣Ds的第i行第j列元素。由于k近鄰關(guān)系不對稱（vj是vi的k近鄰，而vi未必是vj的k近鄰），因此，通過上述方法構(gòu)造的圖是有向的。通常有2 種方法將k-最近鄰法構(gòu)造的相似度有向圖轉(zhuǎn)化為無向圖：只有2 點互為對方的k近鄰，才將2點連接；只要2 點間有k近鄰關(guān)系，就將2 點連接。本文采用第2種方法構(gòu)造無向圖。

基于上述k-最近鄰法構(gòu)造的相似度無向圖，譜聚類的流程［7］如下。

1）由相似度圖構(gòu)建鄰接矩陣W，其第i行第j列元素wij為節(jié)點vi、vj間的權(quán)重，若兩節(jié)點間不連接，則權(quán)重為0。

2）由鄰接矩陣W構(gòu)建度矩陣Dd，其第i行第i列元素為和節(jié)點vi相連的所有邊的權(quán)重之和，即：

3）計算拉普拉斯矩陣Ld=Dd-W，進而構(gòu)建標準化后的拉普拉斯矩陣D-1/2dLdD-1/2d。

4）將D-1/2dLdD-1/2d最小的k個特征值對應(yīng)的特征向量f標準化，形成n×k維的特征矩陣F。

5）對特征矩陣F進行k-means聚類。

4.2 聚類評價指標

聚類評價指標是用來評價聚類質(zhì)量的指標，常用的聚類評價指標包括輪廓系數(shù)、CH 指標CHI（CalinskiHarabasz Index）和DB 指標DBI（Davies-Bouldin Index）［17］。

1）輪廓系數(shù)。

輪廓系數(shù)結(jié)合了聚類的凝聚度和分離度，用于評估聚類的效果，其值處于［-1，1］之間，其值越大表示聚類效果越好。樣本i的輪廓值s(i)定義為：

2）CHI。

CHI 定義為分離度與緊密度的比值，其值越大表示類自身越緊密，類與類之間越分散，即更優(yōu)的聚類結(jié)果。CHI的計算公式為：

式中：ICHI為CHI；SB為類中各點與類中心的距離平方和；SW為各類中心點與數(shù)據(jù)集中心點的距離平方和；kc為指定的聚類數(shù)。

3）DBI。

DBI 綜合考慮了類內(nèi)緊密性和類間分散性，其值越小表示聚類效果越好。DBI的計算公式為：

式中：IDBI為DBI；Si為類i中樣本點到其所屬類中心的距離平方和；Mij為類i與類j的類中心間距離。

5 算例分析

5.1 算例構(gòu)造

算例選用某地區(qū)2019 年8 月1 日至2019 年8 月14 日共2 周的10 kV 專變和公變共8 230 條負荷曲線。為驗證所提方法的合理性與優(yōu)越性，本文基于輪廓系數(shù)、CHI、DBI 這3 個指標，將本文方法與5 種對比方法的聚類結(jié)果進行比較。5 種對比方法分別為：基于歐氏距離的單日法（簡稱單日法），取一天的負荷數(shù)據(jù)，以歐氏距離衡量相似度并進行聚類［13］；基于歐氏距離的連續(xù)日法（簡稱連續(xù)日法），取多日負荷以歐氏距離衡量相似度并進行聚類；基于FDTW的連續(xù)日法（簡稱FDTW 法），取多日負荷以FDTW距離衡量相似度并進行聚類；基于歐氏距離的均值方差法（簡稱均值方差法），提取負荷曲線均值和方差以形成新的2 維數(shù)組，并采用歐氏距離衡量相似度以及進行聚類；基于歐氏距離的最小覆蓋球法（簡稱覆蓋球法），在第2 節(jié)降維步驟中以歐氏距離代替FDTW 計算日負荷曲線的距離矩陣，并將其代入最小覆蓋球算法。

5.2 成本指標閾值的確定

圖1 中曲線1 和曲線2 與其他曲線間的差異極大，因此，本文所提方法利用式（3）所示的Icost指標，在最小覆蓋球的生成過程中剔除此類異常日負荷曲線。若Icost閾值過小，則可能導致部分正常日負荷曲線被剔除；而若Icost閾值過大，則可能導致部分異常日負荷曲線無法被剔除。因此，本文需要確定Icost閾值。

為確定最優(yōu)的Icost閾值，本文取采樣比例為nch（采樣比例為所選數(shù)據(jù)集中數(shù)據(jù)總數(shù)占總數(shù)據(jù)集中數(shù)據(jù)總數(shù)的比例，本文采用隨機采樣，此處取nch=20%，40%，…，100%）的算例數(shù)據(jù)，在不同Icost閾值下，分別計算所提方法聚類結(jié)果的輪廓系數(shù)指標以評價聚類效果。輪廓系數(shù)越大，聚類效果越好，對應(yīng)的Icost閾值也越合理。不同情況下Icost的取值如下。

1）理想情況下，假設(shè)多日負荷曲線足夠多，且這些曲線均可在某一條日負荷曲線上疊加隨機高斯擾動得到，此時不存在異常日負荷曲線，則3 維空間中的各點將呈球狀均勻分布，球體積增大的比例與球內(nèi)點數(shù)量增加的比例一致，因此，Icost恒等于1。

2）實際情況下，多日負荷曲線中存在異常日負荷曲線，球內(nèi)異常點變得稀疏，此時球體積增大的比例超過球內(nèi)點數(shù)量增加的比例，Icost大于1。

綜上，Icost閾值大于等于1才有剔除異常日負荷曲線的意義，因此，本節(jié)將Icost閾值分別取為1.0、1.5、…、3.0，并設(shè)置不剔除異常日負荷曲線的結(jié)果進行對比。表1 為不同Icost和nch取值下重復采用本文方法聚類10 次得到的輪廓系數(shù)平均值。可以看出，不剔除異常日負荷曲線時輪廓系數(shù)指標最差，當Icost閾值為2.0時輪廓系數(shù)指標較優(yōu)，因此，設(shè)置Icost閾值可有效提高聚類質(zhì)量，本文取Icost閾值為2.0。

表1 不同Icost和nch取值下的輪廓系數(shù)平均值Table 1 Average value of silhouette coefficient under different values of Icost and nch

圖3 為在Icost閾值為2.0 下圖1 所示負荷曲線對應(yīng)的最小覆蓋球。由圖可見，除曲線1 和曲線2 外，其余曲線均在最小覆蓋球內(nèi)或球面上。該球的球心坐標為（-1.775 9，-0.062 5，0.127 9），半徑為0.982 6。曲線1、曲線2對應(yīng)的點坐標分別為（1.4086，1.4757，1.523 4）、（1.376 5，-3.748 4，-2.896 0），其與球心的歐氏距離分別為3.801 9 和5.715 5，相離較遠，因此2條曲線被剔除?？梢?，Icost閾值為2.0 時可有效剔除異常日負荷曲線。

圖3 圖1所示負荷曲線對應(yīng)的最小覆蓋球Fig.3 Minimum covering ball corresponding to load curves shown in Fig.1

5.3 負荷曲線-最小覆蓋球-相似度間的一致性

為了說明2 周負荷曲線、最小覆蓋球、相似度間的一致性關(guān)系，選取附錄A 圖A2 所示的2 周負荷曲線中的4 條曲線進行分析。由圖可看出：曲線1—3為典型的雙峰型負荷曲線，它們之間的相似度較高，且曲線1 和曲線3 的相似度最高；曲線4 為晚高峰型負荷曲線，與曲線1—3 的相似度均較低。因此，這4 條負荷曲線涵蓋了極不相似、較相似、十分相似這3種情況。

附錄A 圖A3 為附錄A 圖A2 中4 條負荷曲線對應(yīng)的最小覆蓋球。由圖可以看出：曲線1—3對應(yīng)的最小覆蓋球之間均相交，且曲線1和曲線3對應(yīng)的最小覆蓋球相交程度最高；曲線4 對應(yīng)的最小覆蓋球與曲線1—3對應(yīng)的最小覆蓋球之間均相離。

附錄A 表A1 為附錄A 圖A2 中4 條負荷曲線間的相似度εsi，其中曲線自身之間的相似度顯然為1。由表可以看出：曲線1—3間的εsi均大于0，且曲線1和曲線3 間的εsi最高，為0.842；曲線4 與曲線1—3間的εsi均為0。

綜上可知：2 周負荷曲線間相似度高的，其對應(yīng)的最小覆蓋球間的相交程度較高，對應(yīng)的相似度εsi也較高；2 周負荷曲線間相似度低的，其對應(yīng)的最小覆蓋球間的相交程度較低（甚至可能相離），對應(yīng)的相似度εsi也較低。因此，本文方法中2周負荷曲線、最小覆蓋球、相似度具有一致性，驗證了本文方法的合理性。

5.4 算例結(jié)果的比較分析

為確定最佳聚類數(shù)，本文采用輪廓系數(shù)、DBI、CHI作為評價指標，將3個指標最優(yōu)時對應(yīng)的聚類數(shù)作為各方法的最佳聚類數(shù)。經(jīng)過比較，各方法的最佳聚類數(shù)均為4。表2為不同方法的聚類性能比較。由表可見：均值方差法、覆蓋球法及本文方法在3 個指標上均有明顯優(yōu)勢，這說明無論是單日負荷曲線還是多日負荷曲線聚類，直接聚類的效果均較差，而均值方差法、覆蓋球法及本文方法這類通過提取負荷特征進行聚類的方法效果較好；本文方法比其他方法在輪廓系數(shù)和DBI 上均有所改善，在CHI 上僅小于覆蓋球法，但是差距較小。這說明本文方法的聚類效果有較明顯的優(yōu)勢。此外，由表中各方法的計算時間可見，本文方法耗時較長，具有一定劣勢，本文方法的時間損耗主要由FDTW 和最小覆蓋球算法產(chǎn)生，由于本文研究是基于MATLAB 平臺展開的，而MATLAB 對于FDTW 和最小覆蓋球算法的優(yōu)化十分有限，且無法使用多核運算，因而大幅影響了效率。雖然本文方法的計算時間約為1 h，但是聚類結(jié)果的應(yīng)用可以大幅縮短后續(xù)負荷預測、負荷建模的耗時，因此仍具有一定的實用性。

表2 6種方法的聚類性能比較Table 2 Comparison of clustering performance among six methods

圖4 為均值方差法、覆蓋球法及本文方法的負荷曲線聚類中心。

圖4（a）為均值方差法的負荷曲線聚類中心。其中，類1 曲線為典型居民負荷曲線，類2 曲線為典型商業(yè)負荷曲線，類3 曲線為典型工業(yè)負荷曲線，類4曲線為非典型負荷曲線。由圖可見，類4 曲線近似為類1和類3曲線的平均曲線，曲線特點不明確。

圖4（b）為覆蓋球法的負荷曲線聚類中心，4 類曲線含義與圖4（a）中相同。由于各負荷曲線在聚類前均經(jīng)過max-min 歸一化處理，因此，曲線數(shù)值分布在［0，1］之間，而類4 曲線的數(shù)值在各時刻均保持在0.5附近，這說明該類負荷曲線在各時刻的數(shù)值波動均較大，負荷曲線較雜亂，聚類結(jié)果不合理。

圖4 3種方法的負荷曲線聚類中心Fig.4 Clustering centers of load curves for three methods

表3 為均值方差法、覆蓋球法及本文方法下各類負荷曲線數(shù)量占比。由表可見，均值方差法和覆蓋球法的類4曲線數(shù)量占比分別為19.5%和26.3%，數(shù)量較多，可排除采樣異常的影響。因此，可判斷這2 類曲線的出現(xiàn)均是聚類算法的缺陷所導致的，其實際意義不明確。

表3 3種方法下各類負荷曲線數(shù)量占比Table 3 Percentage of each type of load curve under three methods

圖4（c）為本文方法的負荷曲線聚類中心。其中，類1 曲線為典型的居民負荷曲線，類2 曲線為典型的商業(yè)負荷曲線，類3 和類4 曲線為數(shù)值不同、數(shù)量大致相同的典型工業(yè)負荷曲線。由圖可見，本文方法聚類結(jié)果中各類曲線的定位和特點明確，具有更好的實用性。

5.5 算例魯棒性分析

為檢驗本文方法的魯棒性，在上述算例的基礎(chǔ)上疊加比例為±r（r=5%，10%，…，25%）的均勻分布擾動，以模擬在實際用戶負荷曲線采樣過程中因隨機因素造成的日內(nèi)負荷波動。

為說明不同程度的擾動對聚類結(jié)果的影響，本文設(shè)置Im指標描述不同干擾下各方法的聚類結(jié)果與無干擾下聚類結(jié)果的一致程度，其值越高則表示該方法的魯棒性越高。Im的計算公式為：

式中：nm為某方法在不同干擾下與在無干擾下聚類結(jié)果一致的樣本數(shù)量；nall為該方法的樣本總數(shù)。

圖5 為不同程度擾動下6 種方法的Im指標。由圖可知：隨著負荷曲線擾動程度的增加，各方法的Im指標均呈下降趨勢；隨著擾動程度的增加，單日法、連續(xù)日法、FDTW 法的Im指標下降較明顯，這說明這3 種方法的魯棒性明顯不足；均值方差法、覆蓋球法及本文方法由于考慮了多日負荷曲線并對其進行降維，因此，在各干擾下均有較高的Im指標，且本文方法的Im指標最高，這說明本文方法具有較高的魯棒性。

圖5 不同干擾下的Im指標Fig.5 Index of Im under different interferences

6 結(jié)論

本文針對目前負荷曲線聚類準確度和魯棒性的不足，提出一種基于FDTW 和最小覆蓋球的多日負荷曲線聚類方法。利用FDTW 和MDS 將負荷曲線降至3維，并在3維空間中對每個負荷尋找最小覆蓋球。在此基礎(chǔ)上，定義最小覆蓋球間的相似度，進而得到相似度矩陣。最后利用譜聚類算法得到聚類結(jié)果。所提方法具有以下特點：能同時考慮日負荷曲線內(nèi)的時間滯后特性和多日負荷曲線間的波動特性；在最小覆蓋球的迭代過程中自動剔除異常負荷。算例結(jié)果表明，該方法可有效提高聚類準確度及魯棒性，具有較好的實用價值。

附錄見本刊網(wǎng)絡(luò)版（http：//www.epae.cn）。