改進多目標螢火蟲優化的軟子空間聚類算法及在短期負荷預測中的應用

張 曦 康 平 付雪峰,2 葉 軍,2 趙 嘉,2,3*

1(南昌工程學院信息工程學院 江西 南昌 330000) 2(南昌工程學院江西省水信息協同感知與智能處理重點實驗室 江西 南昌 330000) 3(南昌工程學院鄱陽湖流域水工程安全與資源高效利用國家地方聯合工程實驗室 江西 南昌 330000)

0 引 言

在大數據時代，高維數據聚類是一項重要且具挑戰性的工作。高維數據中大量的無關屬性使目標簇所在的子空間只與某些維度有關；稀疏的數據分布使數據之間的距離幾乎相等，導致傳統聚類方法中的全維空間距離度量失效。為解決以上問題，軟子空間聚類(Soft Subspace Clustering,SSC)為每個維度分配一個0到1之間的權值，以權重來表示維度對目標簇的貢獻程度[1]。

許多傳統的軟子空間聚類將若干標準合并成一個目標函數進行優化，由于這些標準之間可能會相互沖突，其在目標函數中引入加權參數來平衡這些標準以提高性能，但是如何設置適當的加權參數仍然是一個未解決的難題，而且聚類效果對加權參數的設置極為依賴，同一參數無法適用于不同的數據集。因此，相關學者將軟子空間聚類問題轉化為多目標優化問題(Multi-objective Optimization Problem,MOP)以減少這種加權參數，并運用多目標進化算法對多個目標函數同時進行優化，以提高軟子空間聚類的性能[2]。Zhu等[3]提出了一種基于多目標進化算法的軟子空間聚類(Soft subspace clustering with a multi-objective evolutionary approach,MOSSC)。Xia等[4]考慮到權重的分布，將負權值熵和簇間分離度整合成一項，并將其與加權簇內緊湊度視為兩個目標同時進行優化，提出了一種基于多目標進化方法的軟子空間聚類(A novel soft subspace clustering with a multi-objective evolutionary approach,MOEASSC)。在MOEASSC的基礎上，Liu等[5]采用基于分解的多目標進化算法(A multi-objective evolutionary algorithm based on decomposition,MOEA/D)[6]對簇內緊湊度、簇間分離度、負權值熵三個目標同時進行優化，提出了一種對高維數據聚類的改進多目標進化方法(An improved multi-objective evolutionary approach for clustering high dimensional data,MOEC)。雖然MOEC比一些軟子空間聚類算法的性能要好，但是其采用的MOEA/D算法尋優精度不高，獲得的最優解集可能分布不均勻[7]，影響了聚類效果。

2013年劍橋學者Yang[8]提出了一種多目標螢火蟲算法(Multi-Objective Firefly Algorithm,MOFA)，相對于幾種常用的多目標進化算法而言能得到較均勻的最優解集。由于MOFA易提前收斂、尋優精度較低，許多MOFA的改進算法相繼被提出。謝承旺等[7]采用正交實驗方法產生初始群體，利用外部檔案中的精英解指導螢火蟲移動，并運用3點最短路徑方法保持檔案群體多樣性，提出了一種混合型多目標螢火蟲算法。Lv等[9]在螢火蟲學習公式中加入補償因子進行迭代，并隨機選取一個外部檔案粒子作為精英粒子參與種群的進化，提出了一種基于補償因子與精英學習的多目標螢火蟲算法(Multi-objective firefly algorithm based on compensation factor and elite learning,CFMOFA)。現有的MOFA能較好地解決一些MOP問題，但其性能仍有很大的提升空間[9]。

基于上述思想，本文采用多目標螢火蟲算法優化簇內緊湊度、簇間分離度及負權值熵三個目標函數，同時，為彌補MOFA易提前收斂、尋優精度較低的缺陷，將螢火蟲位置更新公式中的步長α、吸引力β0進行動態定義，并設計一種螢火蟲單行隨機學習的機制來提高最優解集分布的均勻性，提出一種改進多目標螢火蟲優化的軟子空間聚類算法(IMOFASSC)。IMOFASSC將簇心和權重擬化成螢火蟲種群，通過螢火蟲位置的更新來搜索最優簇心和權重，并發掘子空間中隱藏的簇。在UCI標準數據集和癌癥基因表達數據集上的實驗結果表明，IMOFASSC算法在處理低維和高維數據時較其他算法具有更好的聚類效果。

短期負荷預測對于電力系統的可靠性至關重要，因為負荷預測結果是為確定系統是否易受攻擊而進行離線網絡分析的基礎。如果預測不充分，則會導致分配不足的備用容量和使用昂貴的調峰裝置，或者導致不必要的大備用容量，這兩者都與增加系統運行成本有關。此外，在公開市場環境中，短期負荷的準確預測是電能交易和現貨價格確定的基礎，旨在獲得最低的購電成本。因此，開展短期電力負荷預測研究工作，提高其預測精度，具有重大的現實意義[10]。基于該背景，本文將IMOFASSC算法應用到短期電力負荷預測中，仿真實驗結果表明IMOFASSC算法在短期電力負荷預測中具有一定的推廣應用價值。

1 相關工作

1.1 軟子空間聚類

軟子空間聚類算法的最小化目標函數可以概括為：

J(V,U,W)=Jwc(V,U,W)+γ·Jadd(V,U,W)

(1)

Jwc通常表示為：

(2)

式中：Jwc表示加權簇內緊湊度，其值越小簇內樣本點越緊湊；Jadd是一個從其他角度評估聚類結果的附加項；γ是用戶設定的參數，用于平衡Jwc和Jadd；V=[vik]C×D為簇心矩陣；U=[uij]C×N是劃分矩陣，表示N個樣本分別屬于C個簇的程度；W=[wik]C×D為權重矩陣，表示C個簇分別在D個維度上的權重；X=[xjk]N×D是D維空間中N個樣本的數據矩陣；m(m≥1)和β(β≥1)是兩個參數。

最早提出的軟子空間聚類算法不包含Jadd附加項[11-12]，為進一步提高算法性能，許多學者在目標函數中添加一些附加項。為避免權重無限大，Jing等[13]提出了一種模糊加權K均值(Fuzzy weighting k-Means,FWKM)算法，其中：Jadd是全部樣本所在簇的模糊權重的總和；γ是整個樣本的平均離散度。Gan等[14]提出了模糊子空間聚類(Fuzzy Subspace Clustering,FSC)算法，其中：Jadd是模糊權重總和；γ是一個恒定參數。為使算法選擇更多的維度來搜索子空間，熵加權K均值(Entropy weighting k-Means,EWKM)算法[15]和局部自適應聚類(Local Adaptive Clustering,LAC)算法[16]都采用負權值熵作為附加項。Deng等[17]在最小化EWKM目標函數的同時最大化加權簇間分離度，提出了增強的軟子空間聚類(Enhanced Soft Subspace Clustering,ESSC)算法，其目標函數如下：

J(V,U,W)=Jwc(V,U,W)+γ·Jadd(V,U,W)-

η·Jbs(V,U,W)

(3)

式中：Jbs表示加權簇間分離度；η是用戶預設的控制Jbs影響的參數。

在式(1)和式(3)對應的軟子空間聚類算法中，γ和η值難以設置，因為式中項與項之間可能會沖突，沒有確切的方法定義適當的γ和η值來平衡這些項[4,17]。例如，對于Jadd項，EWKM、LAC和ESSC中的負權值熵應被最小化以在更多維度上搜索子空間，可能導致某些簇的子空間被一些不太相關的特征錯誤地構造，雖然最后也能得到較小的Jwc，但是結果不準確；對于Jbs項，在ESSC中它也可能與Jwc沖突，例如，具有較小Jwc的某些簇可能在樣本點密集的區域中擠在一起，使它們之間的分離度較小。

1.2 多目標螢火蟲算法

自然界中，螢火蟲利用自身的發光特性來相互吸引、傳遞信息，實現求偶交配和繁殖的目的。MOFA包含亮度和吸引力兩個尋優要素，亮度決定螢火蟲的移動方向，吸引力體現移動距離，迭代中螢火蟲不斷被較亮螢火蟲吸引并向其移動。在M維目標優化問題中，螢火蟲i適應值函數表示為：

f(xi)=(f1(xi),f2(xi),…,fM(xi))

(4)

其亮度由f(xi)決定，xi表示螢火蟲i的位置。螢火蟲的吸引力是相對的，它隨螢火蟲i與螢火蟲j的距離的變化而變化，同時吸引力也與吸收因子有關，由于傳播媒介對光的吸收，光的亮度隨著與光源的距離增大而變小。吸引力β可定義為：

(5)

(6)

式中：β0為初始吸引力，即在rij=0處的吸引力，一般取為1；γ為光吸收系數，標志吸引力的變化，一般取γ=[0.01,100]；rij為螢火蟲i到螢火蟲j的歐氏距離；d是數據維度。

MOFA以Pareto支配[18]的概念確定螢火蟲之間的吸引關系，即：如果螢火蟲j支配螢火蟲i(記為j?i)，螢火蟲i就被更亮的螢火蟲j吸引，位置更新公式如下：

xi(t+1)=xi(t)+βji(xj(t)-xi(t))+α·εi

(7)

式中：t為當前迭代次數；xi(t)、xj(t)為螢火蟲i和j的第t代空間位置；βji表示螢火蟲j對螢火蟲i的吸引力；α為步長因子，一般取[0,1]上的常數；εi是服從均勻分布或高斯分布或者其他分布的隨機數向量，一般用rand-0.5來表示，rand為[0,1]上服從均勻分布的隨機數。

若螢火蟲不被其他任何螢火蟲支配則采用式(8)更新它的位置。

xi(t+1)=g*+α·εi

(8)

標準MOFA通常通過式(9)將多個目標函數以隨機加權和的方式轉換為單目標函數，g*是ψ(x)值最優時螢火蟲的位置。

(9)

式中：M為目標函數個數；ωm∈[0,1]，并且每一代ωm都不相同。

每一代螢火蟲完成進化后，將所有最優解保存在外部檔案(External Archive,EA)中，并對外部檔案進行維護，以確保外部檔案中所有的解互不支配。為節省計算資源，一般會設定外部檔案的最大容量。當最優解的個數超過外部檔案數時采用自適應網格刪除法[19]確保算法所得的Pareto最優解分布更均勻。

2 改進多目標螢火蟲算法優化的軟子空間聚類(IMOFASSC)

2.1 目標函數

本文采用文獻[5]提出的一種對高維數據聚類的改進多目標進化方法MOEC的目標函數作為IMOFASSSC的目標函數，見式(10)。

(10)

(11)

式中：xjk為第j個數據點第k維的值；Awi表示第i個簇重要權重的平均值，即該簇中值大于1/D的權重的平均值，Awi的最小化可以抑制一些非常大的權重；Sepi表示第i個簇心與其他簇心之間的歐幾里得距離的總和；σ一般設為0.000 1，可避免分母為零。

2.2 多目標螢火蟲算法改進

根據我們之前的研究[20]，螢火蟲算法的性能很大程度上取決于其參數，如式(7)中固定步長α和初始吸引力β0。當螢火蟲算法收斂時，α應滿足當t→∞時limα(t)=0。根據式(5)，β由兩只螢火蟲之間的距離決定，隨著迭代的進行，算法不斷收斂，螢火蟲之間的距離逐漸減小到零，在后期的迭代中吸引力β恒等于β0，對局部搜索沒有幫助。基于以上思路，本文對MOFA中的α和β0做出同樣的修改。

α(t)=(1-t/Maxiter)α0

(12)

(13)

式中：rand1為第rand1代螢火蟲的步長；rand1為初始步長；rand1為最大迭代次數；rand1是第rand1代螢火蟲的初始吸引力；rand1是第rand1代螢火蟲rand1和rand1之間的吸引力；rand1和rand2表示范圍在j服從均勻分布的兩個隨機因子。

MOFA外部檔案中的最優解又稱精英解，精英解包含種群的優質信息，而在種群進化時精英解沒有得到利用，浪費了計算代價。此外，標準MOFA沒有考慮到螢火蟲彼此不存在支配關系時該如何移動的問題，以上原因導致了螢火蟲種群進化緩慢，多樣性差。針對這些問題，本文設計一種螢火蟲單行隨機學習的機制。

(1) 當螢火蟲j支配螢火蟲i時，即j?i時，則i利用式(14)向j移動。

xi(t+1)=xi(t)+β(t)·(xj(t)-xi(t))+α·εi

(14)

(2) 當螢火蟲i不被螢火蟲j支配時，則以螢火蟲i的每一行為研究對象，在外部檔案EA中隨機選取一個精英解g*，利用式(15)對i的當前行進行位置更新并用一個新粒子i′來記錄更新結果。所有行全部更新之后得到xi′，判斷此時的i′和i之間的支配關系，若i′?i，則xi(t+1)=xi′，否則不作處理。

(15)

2.3 編碼方式及初始化

螢火蟲個體采用混合實數編碼，如圖1所示，其中i=(1,2,…,C)，C和D的定義同式(2)。例如：在螢火蟲第i行的位置中，前半部分代表第i個簇的簇心，后半部分代表第i個簇的權重。在計算目標函數之前，應將各個權重除以它們所在簇的權重總和，以確保每個簇的權重和等于1。

圖1 螢火蟲個體編碼

C

/D

V

W

U

(16)

式中：uij表示第j個樣本屬于第i個簇的程度，其確定方法是計算第j個樣本到各個簇心的加權距離，若j到某個簇心的距離最小，則u=1，否則u=0；1≤p≤C。

2.4 IMOFASSC算法步驟

IMOFASSC的流程如算法1所示。

算法1IMOFASSC算法

1.初始化螢火蟲種群并設置參數，如：螢火蟲數量Npop，外部檔案最大容量Epop，最大迭代次數Maxiter；

2.利用式(10)計算每個螢火蟲的三個目標函數值；

3.while(t

4.fori=1:Npop

5.forj=1:Npopd

6.ifxj?xi

7.螢火蟲j利用式(14)向j移動；

8.else

9.forh=1:C

11.end for

12.ifxi′?xi

13.xi=xi′；

14.end if

15.end if

16.利用式(10)更新螢火蟲i的三個目標函數值；

17.end for

18.end for

19.更新并保存外部檔案中的Pareto最優解；

20.t=t+1；

21.end while

22.輸出外部檔案中最優的V和W

3 實驗與結果分析

為檢驗IMOFASSC算法的優良性能，本文使用6個UCI標準數據集和4個癌癥基因表達數據集[21-22]進行測試，并將IMOFASSC與FSC[14]、EWKM[15]、ESSC[16]、MOEASSC[4]和MOEC[5]等5種先進的軟子空間聚類算法進行比較。數據集基本信息在表1中給出。

表1 數據集信息

本文采用Rand Index(RI)[23]和Normalized Mutual Information(NMI)[24]兩個指標評價聚類效果。

(17)

(18)

式中：K為標簽數量；C為聚類簇數；N為樣本總數；f00表示具有不同標簽的點被分配到不同簇的樣本對數；f11表示具有相同標簽的點被分配到同簇的樣本對數；ni表示具有第i個標簽的樣本數；nj表示簇j中的樣本數；nij表示具有第i個標簽的樣本點被分配到簇j的樣本對數。RI和NMI的值在[0,1]區間內，兩指標越接近1，聚類效果越好。

實驗所有數據集歸一化處理，參數設置見表2。為公平比較，設置所有算法的最大迭代次數為5 000，每種算法運行30次，記錄最佳結果見表3、表4。由表3、表4中的RI、NMI評價指標可知，IMOFASSC算法在6個UCI數據集上的RI、NMI最好值均明顯優于其他5種算法，表明IMOFASSC在處理低維數據集時具有較高的精度。其次，在CNS_tumors、DLBCL_B、StJude_Leukemia三個癌癥基因表達數據集上，IMOFASSC取得的RI、NMI最好值均最優，這說明在高維數據維度災難的影響下，IMOFASSC不僅性能更好而且具有一定的穩定性。最后，從所有數據集上來看，IMOFASSC一共在9個數據集上取得了最優的RI、NMI值。此外，還可發現MOEASSC和MOEC的整體性能也比單目標的FSC、EWKM和ESSC要好，這印證了第1節中的觀點，即式(1)和式(3)中的項彼此沖突，γ和η參數的設置對聚類效果影響極大，多目標進化算法實現的軟子空間聚類可以同時優化這些項有效地解決這個問題。

表2 算法參數設置

表3 每種算法運行30次的最佳結果(RI)

表4 每種算法運行30次的最佳結果(NMI)

綜上所述，對于低維UCI數據集，IMOFASSC算法能夠有效地提高聚類精度，同時在高維癌癥基因表達數據集上也取得了較好的聚類效果。

4 短期負荷預測模擬實驗

4.1 IMOFASSC-SVR模型

短期負荷預測是一個受時間、經濟、季節、天氣和隨機效應多因素影響的復雜的非線性問題。在眾多預測方法中，支持向量回歸(Support vector machine for regression,SVR)在解決有限樣本、凸二次規劃及非線性高維的問題時具有優異的性能，已成為短期負荷預測領域的研究熱點問題之一[25]。但是，SVR預測性能對訓練樣本數據的相關性極為敏感。傳統SVR預測方法未對訓練樣本進行篩選，訓練樣本數據量大，相關性較弱，導致預測性能不佳。針對上述問題，本文采用IMOFASSC實現SVR訓練樣本的自動選取，建立IMOFASSC-SVR模型。

本模型利用IMOFASSC對歷史日和預測日的影響因子數據進行聚類，并找到預測日所在的類，則該類中的歷史日即為預測日的相似日，這些歷史日的負荷數據與相應的影響因子數據共同構成SVR預測模型的訓練樣本。為檢驗模型預測性能，采用平均絕對百分誤差(Mean Absolute Percentage Error,MAPE)作為評價指標，計算公式如下：

(19)

具體預測步驟如下：

(1) 尋找歷史負荷數據，確定影響因子，并進行歸一化處理。

(2) 利用改進的軟子空間聚類選取與預測日相似度較高的歷史負荷數據，與影響因子一起作為訓練樣本。

(3) 將訓練樣本的影響因子作為SVR的輸入向量，對應的負荷數據作為輸出向量，訓練SVR預測模型。

(4) 在訓練好的SVR預測模型中輸入預測日的影響因子，輸出值反歸一化后即為相應的預測結果。

(5) 計算預測誤差MAPE。

4.2 仿真實驗

本節實驗數據引自第九屆電工數學建模競賽試題A，即以某地區2014年10月25日至12月24日每天96個時刻點的負荷數據和氣象數據作為的訓練樣本，2014年12月25日至12月31日每天96個時刻點的實際負荷值作為測試樣本進行實驗，并將IMOFASSC-SVR模型與SVR模型進行比較。

為公平起見，IMOFASSC算法中的聚類簇數為3，種群規模為20，最大迭代次數為100，其他參數與第3節相同。IMOFASSC-SVR與SVR模型的懲罰系數C、不敏感損失系數ε和核參數σ的取值分別設為默認值，即[0.01,1 000]、[0,1]和[0.1,100]，SVR核函數設為徑向基核函數。實驗中IMOFASSC-SVR和SVR模型獨立運行20次，記錄兩種模型得到的最佳預測結果，繪制出兩種模型的預測結果，見圖2。

圖2 兩種模型的預測結果

圖2表示這兩種模型7天672個采樣點的負荷預測值與實際負荷值的比較。可以看出，兩種模型的預測負荷曲線與實際負荷曲線的走勢基本相同，但在前3天和第4至第5天期間，兩種模型的預測負荷曲線與實際負荷曲線基本重合，說明兩種模型的預測效果相當；在第4天預測階段，兩種模型的預測負荷曲線均與實際負荷曲線存在較大的間距，表明兩種模型的預測效果不佳；在第7天的預測中，SVR模型的預測曲線與實際負荷曲線存在明顯的間距，IMOFASSC-SVR模型的預測曲線與實際負荷曲線間距較小，說明IMOFASSC-SVR模型的預測效果較好。

為進一步比較模型性能，將IMOFASSC-SVR模型的預測誤差MAPE與標準SVR預測模型進行對比，結果見表5。

表5 四種模型的預測誤差(MAPE)(%)

續表5

表5顯示了兩種模型分別對12月25日至12月31日的電力負荷進行預測時產生的誤差MAPE值，以及兩種模型在7天預測中的MAPE誤差的平均值。可以看出，在12月25日、26日、29日、31日的預測中，由于IMOFASSC-SVR采用IMOFASSC自動選取預測日的相似日訓練樣本使訓練出的SVR模型更準確，有利于后續預測，因此IMOFASSC-SVR模型的MAPE值均小于SVR模型，尤其31日的誤差較SVR模型明顯降低了4.37百分點，表明IMOFASSC-SVR模型在這幾天的預測效果優于SVR模型；對于27日、28日、30日，SVR模型取得了最優的MAPE值，這是因為當預測日的氣象情況差別很大時，在氣象數據的聚類過程中僅設置一個固定的聚類簇數變得不準確；從兩種模型7天的平均MAPE值來看，IMOFASSC-SVR模型7天的平均MAPE值為2.49%，SVR模型為3.34%，表明IMOFASSC-SVR模型在這7天的整體預測效果最優；從全局來看，IMOFASSC-SVR模型在7天的預測中有4天取得了最優的預測效果，并且7天的整體預測效果也是最優的，而SVR模型只在3天中的預測效果較好，以上分析表明IMOFASSC-SVR模型的預測性能優于SVR模型。

5 結 語

本文提出一種改進多目標螢火蟲優化的軟子空間聚類算法。IMOFASSC首先對MOFA的步長α、吸引力β0進行了動態定義；設計一種螢火蟲單行隨機學習機制來提高最優解集分布的均勻性；將改進的MOFA運用到軟子空間聚類問題中，同時優化簇內緊湊度、簇間分離度及負權值熵三個目標函數。IMOFASSC將簇心矩陣和權重矩陣擬化成螢火蟲種群，通過螢火蟲位置的更新來搜索最優簇心和權重，并發掘子空間中隱藏的簇。在UCI標準數據集、癌癥基因表達數據集和短期電力負荷預測的實驗結果表明，IMOFASSC不僅在低維數據聚類中有較高的精度，而且在處理高維數據時也表現出良好的性能，可以有效地解決軟子空間聚類中的MOP問題。在短期電力負荷預測模擬實驗中，IMOFASSC有效地去掉訓練樣本中一些無關數據，加強數據間的相關性，使IMOFASSC-SVR模型預測精度較SVR更高，表明了IMOFASSC在短期負荷預測中具有一定推廣應用價值。下一步的研究方向是繼續研究改進IMOFASSC，使其能同時識別子空間和簇的數量。