基于改進K-means聚類k值選擇算法的配網電壓數據異常檢測

劉明群，何 鑫，覃日升，姜 訸，孟 賢

(云南電網有限責任公司電力科學研究院，云南 昆明 650217)

近年來，配電網電壓監測系統建設加強了對配電網電壓、電流等數據的管理與計算分析，實現了電網故障預警、電壓準實時在線監測等功能[1]。然而，由于配電網規模和配電網監測數據日益增大，傳統在線監測方法已無法滿足監測系統對數據挖掘的快速性和準確性要求。因此，為保障配電網安全可靠運行，配電網電壓監測系統數據異常檢測方法研究具有重要意義。

聚類算法通過對配電網電壓數據聚類，挖掘數據特征、區分數據類型、實現實時故障預警。聚類算法可分為劃分聚類法、層次聚類法和密度聚類法等，其中劃分聚類法計算效率較高，但需事先假定聚類數值[2]，包括K-means、K-medoids和CLARA算法。特別地，K-means算法原理簡單且效率高，非常適合配電網電壓數據在線監測；K-means算法僅在多項式時間內收斂到局部最優解[3]，但合理地選擇初值將有利于收斂到全局最優解。文獻[4]采用K-means++算法選取K-means算法初值，但運用K-means算法要事先假定聚類數。聚類數可根據樣本情況判斷[5]，然而，在相關經驗不足或數據量過大等情況下無法給出最佳聚類數。對于靜態數據庫而言，聚類數不會改變，但當在線監測系統數據庫是動態時，聚類數會隨著配電網故障等問題的產生而動態變化，導致K-means算法的聚類效果變差[6]，因而需要其他算法輔助選擇聚類數。

為解決聚類數選擇問題，目前已有學者提出多種聚類數選擇算法。文獻[7]提出輪廓系數法，算法原理簡單且只需給定聚類數上限，對最佳聚類數估計效果很好，但計算速度較慢，不適合在線監測配電網電壓；文獻[8]運用DBI算法為K-means算法選取聚類數，DBI算法對最佳聚類數估計效果較好，但計算速度稍慢；文獻[9]運用Canopy算法快速自適應選取聚類數，但該算法需根據交叉驗證法或先驗知識設定松閾值和緊閾值，且緊閾值嚴重影響聚類數的選取，因此算法自適應能力不夠強；文獻[10-11]運用elbow method選取聚類數，該方法因簡單直觀且計算速度快而被廣泛應用，但常常因“肘部點”不明顯而無法估計最佳聚類數。上述常用算法存在自適應能力不強、計算速度慢和準確率不夠高等問題，不適合異常數據在線檢測。

針對上述常用算法問題，本文提出一種快速選取聚類數的自適應算法。所提算法基于elbow method和輪廓系數法，首先利用自適應變化閾值求解聚類數下限，接著在聚類數上、下限內計算輪廓系數。為提高算法速度，提出“一個極大值”規則，避免計算所有輪廓系數。該算法充分考慮elbow method的快速性和輪廓系數法的高準確率特點，為K-means算法自動選取聚類數，使K-means算法在線監測成為可能。最后，為評價所提算法，以2個實際配網電壓數據為例，通過仿真對比其他聚類數選擇算法。結果表明，相比于所對比算法，所提算法能以最高準確率和最快計算速度自適應選取最佳聚類數。

1 K-means聚類算法

K-means是一種基于劃分的無監督聚類算法，能將數據集分成k類，其中k是事先假定的。K-means算法隨機產生k個聚類中心，根據最近鄰原則將數據點歸類離其最近的聚類中心，形成k個類，并重新計算各類的聚類中心，重復上述步驟直到聚類中心不再改變位置或達到規定的迭代次數。

K-means算法聚類的目標是使各類誤差平方和(sum of the squared errors, SSE)最小，即

(1)

式中Ci為第i個聚類；p為Ci中的樣本點；mi為Ci的聚類中心，即Ci中所有樣本的均值；eSSE為所有樣本的聚類誤差，代表聚類效果的好壞。

根據拉格朗日定理和最小二乘法原理，確保SSE最小的聚類中心[3]應滿足：

(2)

由式(2)可得：

(3)

其中，ni為第i聚類的樣本總量。因此，聚類中心是該聚類數據的平均值。K-means算法每一次迭代將聚類中心取為該聚類樣本的平均值，確保SSE在本次迭代內達到最小，交替采用最近鄰原則和以均值計算聚類中心，使SSE不斷下降，直到平衡收斂。

2 基于elbow method和輪廓系數的聚類數自適應確定

結合elbow method和輪廓系數，提出改進的elbow method和輪廓系數算法(improved elbow method and silhouette coefficient，IES)，用于自適應確定聚類數，從而與K-means算法結合為自適應K-means算法。對于假定的聚類數上限kmax，首先，IES算法基于elbow method的聚類評價指標SSE值確定聚類數下限kmin；然后，在聚類數搜索范圍[kmin,kmax]內基于輪廓系數搜尋最佳聚類數k*，并利用提出的“一個極大值”規則避免計算[kmin,kmax]內所有聚類數對應的輪廓系數。當最佳聚類數確定后，可用K-means算法進行聚類。

2.1 改進elbow method

IES算法的核心是在一定的定義域內利用輪廓系數尋找最佳聚類數。使輪廓系數最大的聚類數為最佳聚類數k*。輪廓系數法一般在定義域[2,kmax]內計算每一個聚類數k對應的輪廓系數，但輪廓系數計算速度慢，不滿足配電網在線監測的快速性要求。因此，IES算法利用計算速度更快的SSE確定聚類數下限kmin，將定義域范圍縮小為[kmin,kmax]。

聚類數為k時elbow method的誤差平方和記為eSSE(k)。k=i時的相對SSE定義為

(4)

式中eSSE(1)為k=1時eSSE(k)的值，也是其最大值。

由于eSSE(k)、相對SSE單調遞減且離散不可導，無法用極值點和拐點估計最佳聚類數，因此考慮設置閾值來估計聚類數下限。同時，由于不同數據集的手肘曲線不同，因而所設置閾值應隨之自適應變化。

定義取值范圍為(0,1)的常數Kelbow，稱為手肘系數。隨著k的增大，相對SSE從eSSE(1)%開始下降，當k=kmax時，相對SSE最大下降量為eSSE(1)%-eSSE(kmax)%，而當相對SSE下降量達到最大下降量的Kelbow倍時，對應的相對SSE定義為手肘閾值eSSE，elbow%，可計算如下：

eSSE,elbow%=eSSE(1)%-Kelbow(eSSE(1)%-

eSSE(kmax)%)=100%-Kelbow·

(100%-eSSE(kmax)%)

(5)

根據大量仿真經驗，取Kelbow=0.5，則

(6)

利用式(6)，將聚類數下限kmin取為使eSSE(k)%≤eSSE，elbow%成立的最小k值。

在聚類數上限為kmax以及Kelbow=0.5時，kmax對應人為限定的相對SSE最大下降量，所選的kmin使相對SSE下降量達到最大下降量的一半以上。應指出，IES算法解出的kmin≥2。

手肘閾值具有自適應性，因為當k>k*后eSSE(k)%變化緩慢，即使kmax較大，eSSE(kmax)%仍接近于eSSE(k*)%，而不同數據集的eSSE(k*)%不同，因此，由eSSE(kmax)%計算的手肘閾值也隨之自適應變化。對不同數據集，當kmax等于樣本總數時，均有eSSE(kmax)=eSSE(kmax)%=0[3]，由式(6)可知，手肘閾值達到其最小值eSSE，elbow%=50%，手肘閾值恒為常數且不隨數據集不同而變化，因此kmax可取值較大，但不應過大。

不同Kelbow取值影響IES算法計算時間和準確率(能否解出k*)。若降低Kelbow，手肘閾值將增大，使得解出的kmin變小，有利于IES算法解出k*(若kmin>k*，則IES算法無法解出k*)，但kmin變小又會使得輪廓系數法的定義域[kmin,kmax]范圍變大，不利于縮短輪廓系數法計算時間。同樣分析Kelbow增大的情況，可知Kelbow不應過大或過小。由此也能看出，手肘閾值隨數據集不同而自適應增大或減小，是在自適應兼顧IES算法對不同數據集的準確率和計算速度。實際應用中可根據配網電壓歷史數據進行測試，適當調整Kelbow取值。改進elbow method確定聚類數下限流程如圖1所示。

圖1 改進elbow method確定聚類數下限流程Figure 1 The flowchart of determinizing the lower limit of clustering number with the improved elbow method

具體步驟如下：

1)讀取人為設定的kmax和數據集；2)對k=1、kmax分別進行K-means聚類，對聚類結果應用式(1)計算SSE，即eSSE(1)、eSSE(kmax)；3)結合式(4)、(6)計算eSSE，elbow%，令k=2；4)對當前k值進行K-means聚類并計算eSSE(k)；5)采用式(4)計算eSSE(k)%；6)判斷eSSE(k)%≤eSSE，elbow%是否成立，不成立則令k自增1并返回步驟4)，成立則跳出循環進入步驟7)；7)記錄此時的k為kmin，該流程結束。

2.2 “一個極大值”規則

聚類數為k時相應的輪廓系數記為S(k)?；谳喞禂邓褜ぷ罴丫垲悢?，在[kmin,kmax]區間內利用輪廓系數算法搜尋最佳聚類數。然而，為確保大于最佳聚類數k*，聚類數上限kmax的設置可能會過大，而輪廓系數計算速度慢，對[kmin,kmax]區間內每一個聚類數k計算輪廓系數將消耗大量時間。

為提高算法速度，IES算法借鑒gap statistic算法中“一個標準錯誤”(1-standard-error)的規則[12](文獻[13]也在其他算法中使用該規則)，提出“一個極大值”規則，即令聚類數k在[kmin,kmax]區間內每次增加1，依次計算輪廓系數S(k)，當S(k)首次出現極大值時停止計算S(k)。使用該規則得到多個輪廓系數，選其中最大值對應的k為最佳聚類數k*。當S(k)在定義域[kmin,kmax]內不存在極大值時，“一個極大值”規則失效，需計算定義域內所有S(k)，選最大值對應的k為k*。本文中“一個極大值”規則在k=K生效是指：對于K>kmin，當k增大到K+1時，出現S(k)的極大值S(K)，IES算法停止計算S(k)。S(K)為極大值是指：S(K)>S(K-1)且S(K)>S(K+1)。

“一個極大值”規則避免計算所有輪廓系數，相當于降低了實際假定的kmax，從而提高算法速度。

2.3 基于輪廓系數自適應確定聚類數

聚類數下限確定后應用“一個極大值”規則，在定義域[kmin,kmax]內利用輪廓系數法求解最佳聚類數。在已計算不同聚類數對應的輪廓系數中自動尋找最大輪廓系數，所對應聚類數為最佳聚類數k*，從而實現自適應確定聚類數?；谳喞禂底赃m應確定聚類數流程如圖2所示。

圖2 基于輪廓系數自適應確定聚類數流程Figure 2 The flowchart of determinizing the clustering number with the improved elbow method based on the silhouette coefficient adaptive determination

具體步驟如下：

1)建立空數組{S}，令k=kmin；2)對當前k值進行K-means聚類，接著對聚類結果計算輪廓系數并放入數組{S}；3)若數組{S}中的元素已達3個及以上，則說明可以判斷是否出現極大值，進入步驟4)，否則令k自增1并回到步驟2)；4)判斷是否出現極大值，即S[-2]>S[-3]、S[-2]>S[-1]是否同時成立，其中S[-1]是數組{S}倒數第1個元素，即本次循環計算得到的輪廓系數；S[-2]、S[-3]分別是倒數第2、3個元素，若不出現極大值，令k自增1并回到步驟2)，若出現極大值則跳出循環進入步驟5)；5)在數組{S}中尋找最大輪廓系數，并記錄對應的聚類數為最佳聚類數k*，IES算法結束。

2.4 基于自適應K-means的實時異常檢測模型

IES算法從圖1流程開始至圖2流程結束。自適應確定聚類數的IES算法與K-means算法結合為自適應K-means算法。正常運行時配電網電壓數據波動范圍較穩定，因此，可利用K-means算法對正常運行數據聚類并得到聚類中心，通過判斷新輸入數據到聚類中心距離是否超過距離閾值H，從而判斷數據是否異常。

H=(h1,h2,…,hk)表示各聚類的閾值，其中k是聚類數，聚類中數據到聚類中心距離的最大值乘以常數D作為H，綜合考慮文獻[14]、[15]的實驗結果，取D=1.04。若某數據X到k個聚類中心Ci距離均超過相應閾值，則判定為異常數據，即異常數據滿足：

|X-Ci|>hi,i=1,2,…,k

(7)

IES算法能在異常檢測中更新正常數據最佳聚類數，并能在發生異常時幫助挖掘異常數據特征。隨著歷史正常運行數據的不斷增多，正常數據的最佳聚類數可能改變，因此，每隔一段時間需用IES算法自適應求解并更新最佳聚類數。當發生異常時，在分析數據異常模式之前，為充分利用當前所有異常數據，可通過IES算法對當前所有正常和異常數據的最佳聚類數自適應快速求解，然后利用K-means算法將異常與正常數據一起聚類，為挖掘異常數據特征和探測異常來源提供信息。除上述基于自適應K-means聚類的方法，文獻[14]還利用了其他方法分析數據異常模式。

基于自適應K-means的實時異常檢測總流程如圖3所示，具體步驟如下：

1)對配電網電壓歷史正常運行數據進行K-means聚類，并根據聚類得到的最優聚類中心和聚類結果更新距離閾值H；2)計算新輸入數據到各個聚類中心的距離并與距離閾值比較；3)若新輸入數據屬于異常數據則標記為異常，否則將其加入歷史正常運行數據；當歷史正常運行數據新增數量達到一定量時，利用IES算法求解并更新最佳聚類數；4)將異常數據與歷史正常運行數據共同作為新的數據集DS；5)IES算法根據事先假定聚類數上限計算數據集DS的最佳聚類數k*；6)用K-means將數據集DS分為k*個聚類；7)利用K-means聚類結果分析數據異常模式。

3 算例分析

3.1 數據集

以2個實際配電網電壓數據集為例(記為D1和D2)，與DBI算法和輪廓系數法進行仿真比較，驗證所提IES算法的有效性。

D1有1 000個樣本點，每個點對應一個時刻三相電壓有效值。為體現所提IES算法的普適性，對A、B、C三相電壓分別加入異常數據進行實驗。對于A相電壓，隨機抽取10%，即100個正常數據，并加上4%～10%噪聲，生成1個數據集，重復進行50次，生成50個數據集，記為A組數據(注意：每個數據集只含10%異常數據，其中50個正常數據加上正噪聲4%～10%，50個正常數據加上負噪聲-10%～-4%)。用同樣方法對B、C相電壓各生成50個數據集，分別記為B、C組數據。

D2有4 000個樣本點，每個點對應一個時刻三相電壓有效值。為說明選取不同聚類數對聚類效果的影響，對A、B、C三相電壓均加入異常數據進行實驗。對于A相，隨機抽取15%，即600個正常數據，在±5%處加高斯分布隨機函數G作為噪聲。設正常數據值為Ndata，則加噪聲后的值為Ndata·(1±0.05)+G。用同樣方法處理B、C相電壓，最終得到的數據集記為D2noise。

3.2 評價指標

各算法對最佳聚類數估計的準確率定義為

(8)

式中N為實驗次數，對A、B、C三相電壓的每一組數據實驗50次，故N=50；NT為對最佳聚類數估計正確的次數。

對于A、B、C三相電壓的3組數據，最佳聚類數為3，即分為1個正常數據聚類和2個異常數據聚類。用DBI算法、輪廓系數法和IES算法估計所有數據集的最佳聚類數k*，若k*=3則估計正確，若k*≠3則估計錯誤。

由于出現異常時需用自適應K-means算法對正常和異常數據聚類，聚類效果影響下一步分析數據異常模式，因此，需說明K-means對正常和異常數據聚類時不同最佳聚類數估計值對聚類效果的影響。為方便說明，簡化為分析二分類異常檢測效果，再給出異常檢測效果評價指標。對于二分類異常檢測，可根據真實和檢測情況將檢測結果分為4類，如表1所示，TN為實際正常且被檢測為正常的樣本，FP為實際正常但被檢測為異常的樣本，FN為實際異常但被檢測為正常的樣本，TP為實際異常且被檢測為異常的樣本。

表1 檢測結果分類Table 1 Classification of test results

對于異常檢測，相比于不將正常樣本判定為異常，更重要的是檢測到更多的異常點[16]。因此，采用召回率評估異常檢測效果，其值越高意味著檢測到越多真實異常點，其最大值為1。召回率：

(9)

3.3 實驗結果與分析

計算環境如下：計算機CPU為Core i7-10700，內存16 GB，主頻2.90 GHz，操作系統為Windows 10(64 bit)，數據分析工具為Python 3、Jupyter NoteBook。

對A、B、C三相電壓的3組數據分別用3種算法進行最佳聚類數估計，為滿足在線監測的自適應性，應取足夠大的聚類數上限kmax，以保證其大于最佳聚類數，因此取kmax=20，實驗結果如表2～4所示。

表2 對A組數據應用3種算法Table 2 Three algorithms are applied to group A data

表3 對B組數據應用3種算法Table 3 Three algorithms are applied to group B data

表4 對C組數據應用3種算法Table 4 Three algorithms are applied to group C data

由表2可以看出，輪廓系數法和IES算法對最佳聚類數的估計值穩定為3；DBI算法的估計值僅在2、3之間波動，其中，對于A組9個數據集，DBI算法解出k*=2，其余41個數據集解出k*=3。對表3、4不再贅述。根據表2～4，計算3種算法對A、B、C三相電壓的3組數據最佳聚類數估計的準確率ω，如表5所示。

表5 3種算法準確率Table 5 Accuracy of three algorithms accuracy %

由于3種算法對k*的估計值只有2和3，因此分別取k*=2、3，對A、B、C三相電壓的3組數據進行K-means聚類，其中，k*=3時將聚類后的2個異常數據聚類合并，從而得到正常和異常數據二分類(注意：實際中不將異常數據聚類合并，因為會損失異常數據特征信息，此處僅是為了計算召回率)。計算不同k*取值下各組數據召回率均值，如表6所示。

表6 不同k*取值下各組數據召回率均值Table 6 Mean recall rates of data in each group under different k* values

由表6可見，選取合適的聚類數能大幅提升異常檢測效果。對于A、B、C三相電壓的3組數據取k*=2顯然不合適，而在表2～4中，DBI算法對k*的估計值多次為2，結合表5可知輪廓系數法對最佳聚類數估計的準確率比DBI算法更高，而IES算法保持了輪廓系數法的高準確率。

為進一步說明輪廓系數法和IES算法在準確率方面比DBI算法更適合為K-means選擇最佳聚類數，對A組數據中某一數據集進行K-means聚類，取k*=3，聚類結果如圖4所示，可明顯區分3類數據，但對于該數據集，DBI算法解出k*=2，而輪廓系數法和IES算法解出k*=3。

圖4 k*=3時K-means算法聚類結果Figure 4 Clustering results of K-means algorithm when k* =3

圖4中聚類2、3為異常數據類，聚類1為正常數據類，此時召回率的值為1，表明K-means算法適合對該數據集聚類，而輪廓系數法和IES算法比DBI算法更適合為K-means選擇最佳聚類數。

為說明選取不同聚類數對聚類效果的影響，分別取聚類數為2～9，對數據集D2noise進行K-means聚類；為方便計算各聚類結果的召回率，人為將異常數據聚類合并，從而得到正常和異常數據二分類。數據集D2noise召回率隨聚類數變化曲線如圖5所示,可見召回率隨著聚類數增大而增大，說明聚類效果越來越好，當聚類數為7、8、9時達到最大值1。應指出，3種聚類數選擇算法對于D2noise的最佳聚類數估計值均為7。

圖5 數據集D2noise召回率隨聚類數變化曲線Figure 5 The curve of recall rate of D2noise changing with clustering number

進一步分析發現，聚類數大于7時會發生模型過擬合。聚類數分別為7、8時數據集D2noise的K-means聚類結果如圖6、7所示。對比圖6、7可知，圖6為最佳聚類，而圖7中將正常數據過擬合為2個聚類(聚類4、7)，不利于對數據進行分析。

圖6 聚類數為7時數據集D2noise的K-means聚類結果Figure 6 K-means clustering results of D2noise when the clustering number is 7

用3種算法對A組50個數據集估計最佳聚類數，記錄平均運行時間和最小、最大運行時間，如表7所示，可見輪廓系數法的最小運行時間大于DBI算法最大運行時間，從運行時間均值也能看出DBI算法運行速度更快。IES算法的運行時間均值小于其他2個算法，計算速度最快，最符合在線監測的快速性。IES算法運行時間波動范圍大于其他2個算法，是因為對于50個數據集，IES算法均解出kmin=3，但“一個極大值”規則在不同的聚類數k值(k>3)處生效，因此，不同數據集的計算量不同，導致運行時間波動。

表7 3種算法運行時間對比Table 7 Running time comparison of three algorithms

綜上所述，盡管DBI算法計算速度稍快于輪廓系數法，但DBI算法準確率是3種算法中最低的。與DBI算法相比，IES算法不僅計算速度更快，而且準確率更高；與輪廓系數法相比，IES算法不僅保持相同的準確率，而且計算速度更快。因此，IES算法兼顧準確率和計算速度，在保證高準確率的前提下縮短了計算時間，提高了K-means算法在線監測的準確率和高效性。

4 結語

K-means聚類算法計算速度快、準確率高，適合配電網在線監測，但當假定聚類數不合適時，可能導致聚類結果不理想。本文提出了一種快速選取聚類數的自適應IES算法，為K-means算法自動選取聚類數，使K-means算法在線監測配電網成為可能。以召回率評價二分類異常檢測效果為例，說明為K-means選取合適聚類數對異常檢測的重要性。IES算法首先利用自適應變化閾值求解聚類數下限，接著在聚類數上、下限內計算輪廓系數。為提高算法速度，提出“一個極大值”規則，避免計算所有輪廓系數。所提IES算法有如下優點。

1)自適應能力強。IES算法只需給定聚類數上限這一參數，且該上限允許較大，即使動態數據庫的最佳聚類數發生一定的改變，也能保證大于最佳聚類數。所提出用于確定聚類數下限的閾值可隨數據集不同而自適應變化，從而自適應兼顧IES算法準確率和計算速度。

2)計算速度快。IES算法利用計算迅速的SSE求解聚類數下限，縮小了最佳聚類數的搜尋范圍，又利用所提出的“一個極大值”規則減少計算量，提高了計算速度。

3)準確率高。IES算法充分利用了輪廓系數高準確率的特點。

算例表明，所提IES算法能自適應快速選取最佳聚類數，與輪廓系數法相比，IES算法準確率相同而計算速度更快，與DBI算法相比，IES算法不僅準確率更高，而且計算速度更快。因此，IES算法兼顧準確率和計算速度，更有利于應用于配電網在線監測。