基于K均值聚類及高斯過程回歸集成的鉛酸電池荷電狀態(tài)預(yù)測

田明光，翟 旭，江穎潔，田保鵬，劉曉晨

（國網(wǎng)山東省電力公司信息通信公司，山東 濟(jì)南 250001）

0 引言

鉛酸電池以其電壓穩(wěn)定、貯存性能好、造價(jià)較低等優(yōu)勢被廣泛地應(yīng)用于各類供電系統(tǒng)中[1]。在鉛酸電池的多種內(nèi)部狀態(tài)參數(shù)中，電池荷電狀態(tài)（State of Charge，SOC）是最重要的，它反映著電池當(dāng)前的剩余電量。準(zhǔn)確地估計(jì)鉛酸電池的SOC有助于人們及時(shí)對電池進(jìn)行充電，從而確保供電系統(tǒng)順利運(yùn)行。因此，鉛酸電池SOC的預(yù)測成為當(dāng)前電池研究領(lǐng)域的一大熱點(diǎn)問題[2]。

現(xiàn)階段，鉛酸電池SOC的預(yù)測方法有很多種，其中最常見的兩種方法是基于濾波的方法[3]和基于機(jī)器學(xué)習(xí)的方法[4]。基于濾波的方法將經(jīng)典的濾波理論（如粒子濾波、卡爾曼濾波等）與常用的數(shù)學(xué)模型（如多項(xiàng)式模型、指數(shù)衰減模型等）相結(jié)合來擬合鉛酸電池SOC的變化曲線。然而，基于濾波的方法只能在數(shù)據(jù)分布已知的條件下取得令人滿意的結(jié)果，具有嚴(yán)重的局限性。隨著機(jī)器學(xué)習(xí)理論的日漸成熟，第二種方法越來越受到重視。基于機(jī)器學(xué)習(xí)的方法通過各種機(jī)器學(xué)習(xí)算法（如高斯過程回歸算法（Gaussian Process Regression，GPR）[5]、支持向量機(jī)算法（Support Vector Machine，SVM）[6]、人工神經(jīng)網(wǎng)絡(luò)算法（Artificial Neural Network，ANN）[7]等）對采集的歷史數(shù)據(jù)（這些數(shù)據(jù)通常被稱為原始訓(xùn)練集）進(jìn)行建模，從而有效地預(yù)測鉛酸電池的SOC。基于機(jī)器學(xué)習(xí)的方法可以進(jìn)一步分為基于全局信息的方法和基于局部信息的方法。前者是直接在原始訓(xùn)練集上訓(xùn)練模型作為最終的預(yù)測模型，這種方法看似簡單合理，卻容易造成信息的冗余。后者是通過某種策略（例如：聚類、信息匹配）提取原始訓(xùn)練集中部分有效數(shù)據(jù)作為訓(xùn)練集，并在新的訓(xùn)練集上訓(xùn)練預(yù)測模型。此方法雖然能夠去除冗余的信息，卻有可能造成信息的丟失。

為了解決上述信息的丟失或冗余問題，本文以高斯過程回歸為基礎(chǔ)模型，并結(jié)合聚類算法和集成學(xué)習(xí)理論，提出一種基于K均值聚類的高斯過程回歸集成算法（K-means Cluster with Ensemble Gaussian Process Regression，KC-EGPR）來預(yù)測鉛酸電池的SOC。該算法的核心思想是采用聚類算法生成原始訓(xùn)練集的多個(gè)包含局部數(shù)據(jù)信息的子集，再利用集成學(xué)習(xí)理論將數(shù)據(jù)層融合問題轉(zhuǎn)化為決策層融合問題。KC-EGPR的步驟如下：首先通過K-means算法對原始訓(xùn)練集進(jìn)行聚類生成若干子訓(xùn)練集，每個(gè)子集包含了原始訓(xùn)練集的某種局部信息；再在每個(gè)子集上訓(xùn)練高斯過程回歸模型（GPR），多個(gè)GPR構(gòu)成一個(gè) GPR的集合；最后利用自適應(yīng)提升算法（AdaBoost）對多個(gè)訓(xùn)練的GPR進(jìn)行集成，得到最終的預(yù)測模型。在三組鉛酸電池?cái)?shù)據(jù)集上的實(shí)驗(yàn)結(jié)果表明，所提方法的性能優(yōu)于對比模型，且具有廣闊的應(yīng)用前景。

1 基礎(chǔ)算法分析

1.1 K均值聚類算法

K均值（K-means）聚類算法[8]是聚類分析中最常用的算法之一，它最大的特點(diǎn)是簡單、易理解、運(yùn)算速度快。K-means聚類算法需要預(yù)先設(shè)定簇（類別）的數(shù)目K值及每個(gè)簇的初始中心點(diǎn)，然后通過迭代的方式不斷地更新中心點(diǎn)和優(yōu)化聚類結(jié)果。K-means聚類算法的步驟描述如下：

輸入：訓(xùn)練樣本集 S = { x1, x2,… , xp, … , xm}，其中x∈Rn，m表示樣本的數(shù)目；簇的數(shù)目k；迭代p終止的閾值ε；

輸出：k個(gè)簇（ S D1， SD2,… , S Dk）；

步驟 1 選取k個(gè)聚類中心點(diǎn)，依次表示為u1, u2,...,uk。

步驟2 將k個(gè)簇依次置空：

對于訓(xùn)練集S中的每個(gè)樣本px，計(jì)算它到k個(gè)中心點(diǎn)的距離并將其分到距離最小的中心點(diǎn)所對應(yīng)的簇中，用公式表示為：

其中pC表示樣本px所屬的簇。

步驟 3 對于每個(gè)簇i，根據(jù)下面公式重新計(jì)算它的中心點(diǎn)iu：

其中， su m _ feature(i)表示簇i中所有樣本的特征和， s um _ n umber(i)表示簇i中所有樣本的數(shù)目。

步驟 4 不斷重復(fù)步驟 2和 3，直到代價(jià)函數(shù)J(u,c)小于給定的閾值ε。J(u,c)表示樣本集S中所有樣本到其對應(yīng)簇的中心點(diǎn)的歐式距離之和，它的計(jì)算公式如下：

1.2 高斯過程及其回歸算法

高斯過程（Gaussian Process，GP）[9]也叫正態(tài)隨機(jī)過程，它對處理小樣本、非線性、高維數(shù)等復(fù)雜問題具有很好的適應(yīng)性。給定一個(gè)有限的樣本集S = { xi, yi和隱函數(shù) f :Rn→R，其中xi為第i個(gè)樣本的n維特征向量（ x ∈Rn），y為第i個(gè)樣本的i觀察目標(biāo)值(yi∈R)，M表示樣本的數(shù)目。若f(x1),f(x2),… ,f(xM)可構(gòu)成隨機(jī)變量的一個(gè)集合，其性質(zhì)完全由均值函數(shù) m (x)和核函數(shù) k (xi, xj)確定，且具有聯(lián)合高斯分布，則它們形成的隨機(jī)過程就稱之為高斯過程，表示為：

其中 []EX表示期望。

將噪聲ε考慮到觀測目標(biāo)值y時(shí)，可建立高斯過程回歸模型，定義如下：

其中ε為與 f (x)不相關(guān)的獨(dú)立高斯白噪聲，服從均值為0，方差為 σM2的正態(tài)分布。因?yàn)?f (x)服從高斯分布，所以y同樣服從高斯分布。因此樣本集S中所有的觀察目標(biāo)值y聯(lián)合分布的集合可以形成一個(gè)高斯過程，表示為：

給定測試樣本 x*，其預(yù)測值為 f (x*)。根據(jù)貝葉斯原理，樣本集S中的觀察目標(biāo)值y與 f (x*)的聯(lián)合分布表示為：

其中， K (x*, X ) = K(X, x*)為測試樣本 x*與訓(xùn)練集X之間的 M ×1階協(xié)方差矩陣。

因此，可以計(jì)算出測試樣本 x*的預(yù)測值 f (x*)所滿足的后驗(yàn)條件分布，即GP模型回歸方程為：

其中，則 m (f(x*))是測試點(diǎn) x*對應(yīng)的預(yù)測值f(x*)的均值，c o v(f(x*))是測試點(diǎn) x*對應(yīng)的預(yù)測值f(x*)的方差。

在高斯過程回歸中，核函數(shù) k (xi, xj)起到了重要的作用，它定義如下：如果原始特征向量的內(nèi)積是＜ xi, xj＞ ，映射后的內(nèi)積為 ＜ φ( xi) ,φ(xj)＞，那么定義核函數(shù)為 k (x,x ) = φ( x )Tφ(x )。常用的核函數(shù)有線ijij性核函數(shù)（Linear Kernel）、多項(xiàng)式核函數(shù)（Polynomial Kernel）和徑向基函數(shù)(Radial Basis Function, RBF)。本文選用RBF核函數(shù)，它的定義如下：

其中，σ表示RBF核的半徑。

1.3 自適應(yīng)提升算法

自適應(yīng)提升（Adaptive Boosting，AdaBoost）算法[10]因其泛化能力強(qiáng)、無需參數(shù)調(diào)整等優(yōu)點(diǎn)成為集成學(xué)習(xí)領(lǐng)域最常用的算法之一。AdaBoost算法的核心思想如下：首先在同一訓(xùn)練集上訓(xùn)練若干不同的預(yù)測模型，然后通過迭代的方式計(jì)算每個(gè)模型的錯(cuò)誤率，再根據(jù)錯(cuò)誤率計(jì)算出每個(gè)模型對應(yīng)的權(quán)重，最后根據(jù)權(quán)重將所有的模型線性組合得到最終的預(yù)測模型。

AdaBoost 算法既可用于分類模型的集成，也可以用于回歸模型的集成。在回歸問題中，AdaBoost算法的基本步驟描述如下：

輸入：訓(xùn)練樣本集 S = { xi, yi}，其中 xi為第i個(gè)樣本的特征向量，yi為第i個(gè)樣本的回歸值，N表示樣本的數(shù)目；迭代次數(shù)T；判斷樣本是否被正確預(yù)測的閾值ε；

輸出：最終的預(yù)測模型 Mf(x)；

步驟 1 初始化訓(xùn)練樣本集S的權(quán)值分布 D1，使得每個(gè)樣本的權(quán)重相同，具體方式如下：

步驟 2 執(zhí)行T輪學(xué)習(xí)，第t（ t = 1,2,… ,T）輪學(xué)習(xí)過程如下：

（1）在權(quán)值分布 Dt的訓(xùn)練集S上訓(xùn)練模型 Mt(x)；

（2）評(píng)估 Mt(x)在訓(xùn)練集S上的錯(cuò)誤率 et：

其中ε為預(yù)先設(shè)置的閾值，如果樣本 xi的預(yù)測值與真實(shí)值之間的誤差小于ε，則認(rèn)為樣本 xi被正確預(yù)測；否則，樣本 xi被錯(cuò)誤預(yù)測；

（3）根據(jù)錯(cuò)誤率 et計(jì)算 Mt(x)的權(quán)重βt：

如果模型 Mt(x)的錯(cuò)誤率 et大于或等于0.5，則舍棄該模型；

（4）更新訓(xùn)練集S上的權(quán)值分布 Dt+1：

其中，Zt為標(biāo)準(zhǔn)化因子。若樣本 xi被正確預(yù)測，則減小 xi對應(yīng)的樣本權(quán)重；否則，增大 xi對應(yīng)的樣本權(quán)重；

盡管經(jīng)典的 AdaBoost算法在集成學(xué)習(xí)領(lǐng)域取得了巨大的成功，但它仍然存在一個(gè)缺點(diǎn)——容易導(dǎo)致過擬合。在上述過程中，訓(xùn)練集既被用來訓(xùn)練模型，又被當(dāng)作評(píng)估數(shù)據(jù)集來計(jì)算模型的錯(cuò)誤率。換言之，評(píng)估數(shù)據(jù)集和訓(xùn)練集實(shí)際上來源于同一數(shù)據(jù)集，從而導(dǎo)致最終的模型在訓(xùn)練集上能取得較好的預(yù)測性能，但是在測試集上結(jié)果較差，具有較差的泛化能力。為了避免過擬合現(xiàn)象，本文在AdaBoost算法中引入一個(gè)與訓(xùn)練集沒有相同樣本的獨(dú)立評(píng)估數(shù)據(jù)集（詳見算法1），利用該數(shù)據(jù)集評(píng)估每個(gè)模型的錯(cuò)誤率。

2 基于K均值聚類的高斯過程回歸集成算法

2.1 算法原理和基本步驟

基于K均值聚類的高斯過程回歸集成算法（K-means Cluster with Ensemble Gaussian Process Regression，KC-EGPR）的核心思想是采用聚類算法生成原始訓(xùn)練集的多個(gè)包含局部數(shù)據(jù)信息的子集，再利用集成學(xué)習(xí)理論將數(shù)據(jù)層融合問題轉(zhuǎn)化為決策層融合問題，該算法的流程圖如圖1所示。

根據(jù)圖1，KC-EGPR算法可以分為訓(xùn)練階段和測試階段。

訓(xùn)練階段的步驟描述如下：

輸出：集成模型 A da_GPR；

圖1 基于K均值聚類的高斯過程回歸集成算法的流程圖Fig.1 Flowchart of the ensemble algorithm based on K-means clustering and Gaussian process regression

步驟1 對訓(xùn)練數(shù)據(jù)集TRD進(jìn)行K-means聚類，生成K個(gè)簇，每個(gè)簇作為TRD的子訓(xùn)練集，依次表示為 S D1, S D2,… , S DK；

步驟2 在每個(gè)子訓(xùn)練集 S Dk上訓(xùn)練一個(gè)高斯過程回歸模型，表示為 G PRk；

步驟3 結(jié)合獨(dú)立評(píng)估數(shù)據(jù)集，利用AdaBoost算法對K個(gè)GPR模型進(jìn)行集成，得到集成的GPR模型，表示為：

其中，kβ表示第k個(gè)GPR模型的權(quán)重，且

測試階段的步驟描述如下：

輸入：測試樣本TS；集成的GPR模型 A da_GPR；

輸出：預(yù)測結(jié)果PR；

步驟1分別用 A da_GPR中的每個(gè)GPR模型對TS進(jìn)行預(yù)測，預(yù)測結(jié)果依次表示為 Y1, Y2,… ,YK；

步驟2根據(jù) A da_GPR中每個(gè)GPR模型的權(quán)重，最終的預(yù)測結(jié)果PR可以按照如下公式計(jì)算：

2.2 算法評(píng)價(jià)指標(biāo)

為了評(píng)估 KC-EGPR算法的性能，本文采用絕對百分比均值誤差（Mean Absolute Percentage Error，MAPE）和均方根誤差（Root Mean Square Error，RMSE）兩個(gè)評(píng)價(jià)指標(biāo)，它們的計(jì)算公式如下：

其中，iY表示第i個(gè)樣本的真實(shí)值，iY′表示第i個(gè)樣本的預(yù)測值，N表示樣本的數(shù)目。從上述公式中可以看出，MAPE表示模型相對誤差的均值，RMSE表示模型誤差的標(biāo)準(zhǔn)差，它們的值越小，模型的預(yù)測性能越好。

3 實(shí)驗(yàn)結(jié)果分析

3.1 實(shí)驗(yàn)數(shù)據(jù)集

本文采用恒流放電試驗(yàn)法來獲取鉛酸電池的實(shí)驗(yàn)數(shù)據(jù)，實(shí)驗(yàn)對象為HAB類型的鉛酸蓄電池，額定容量為100 Ah，額定電壓為2 V。電池的放電實(shí)驗(yàn)描述如下：開始電壓為2.2 V，截止電壓為1.9 V，放電電流32 A。在上述過程中，記錄電池的電壓和溫度作為樣本的特征，同時(shí)估算出電池的SOC作為樣本的回歸量，共200個(gè)樣本數(shù)據(jù)。為了驗(yàn)證KCEGPR算法的有效性，本文分別將前140、150、160個(gè)數(shù)據(jù)作為訓(xùn)練數(shù)據(jù)集（TRD），剩余的60、50、40個(gè)數(shù)據(jù)作為測試數(shù)據(jù)集（TED），上述三組數(shù)據(jù)集依次表示為DS1、DS2、DS3。考慮到KC-EGPR算法需要獨(dú)立評(píng)估數(shù)據(jù)集，故從訓(xùn)練數(shù)據(jù)集中分出后10%的數(shù)據(jù)作為獨(dú)立評(píng)估數(shù)據(jù)集（IED）。因此，在DS1、DS2和DS3下，TRD、IED和TED的樣本數(shù)目如表1所示。

表1 TRD、IED和TED在三組數(shù)據(jù)集下的樣本數(shù)目Tab.1 The number of samples of TRD, IED, and TED

3.2 KC-EGPR算法的參數(shù)優(yōu)化

KC-EGPR預(yù)測性能受眾多參數(shù)的影響，其中最重要的參數(shù)是聚類分析中的簇的數(shù)目K。在本文中，我們以1為步長，在區(qū)間[2,8]內(nèi)逐次調(diào)整K值，選擇使 KC-EGPR性能達(dá)到最佳的K為最終的參數(shù)。圖2 描述了KC-EGPR算法在三組數(shù)據(jù)集下的預(yù)測性能（RMSE值）隨簇的數(shù)目K的變化曲線。

從圖2中可以觀察到，隨著簇的數(shù)目K的增加，KC-EGPR的 RMSE值呈先減小后增加的趨勢（該趨勢在 DS3數(shù)據(jù)集上最為明顯），這意味著 KCEGPR的預(yù)測性能先提升后降低。上述現(xiàn)象可以按照如下方式解釋：初始時(shí)（ =2K ），K的值較低，樣本集經(jīng)過聚類后形成的簇的數(shù)目較少，每個(gè)簇包含大量的樣本并且仍然存在冗余信息，因此 KCEGPR的預(yù)測性能偏低；隨著K值的增加（ =3,4K ），樣本集經(jīng)過聚類后形成的簇的數(shù)目增加，這意味著樣本集被劃分的更細(xì)，每個(gè)簇中的冗余信息明顯減少，KC-EGPR的預(yù)測性能明顯提高。隨著K值的持續(xù)增加（( 5)K≥ ），樣本集經(jīng)過聚類后形成過多的簇，每個(gè)簇中包含的信息量嚴(yán)重不足，反而導(dǎo)致KC-EGPR的預(yù)測性能降低。

圖2 KC-EGPR算法在數(shù)據(jù)集DS1、DS2和DS3下的RMSE值隨簇的數(shù)目K的變化曲線Fig.2 The variations of RMSE with different values of K for KC-EGPR algorithm on DS1, DS2, and DS3

3.3 KC-EGPR與基于全局信息的 GPR模型性能比較

前文提到，原始訓(xùn)練集中可能存在冗余信息，直接在原始訓(xùn)練集上訓(xùn)練模型（即基于全局信息的方法）往往不能取得令人滿意的結(jié)果。為了驗(yàn)證KC-EGPR能夠減輕原始訓(xùn)練集中冗余信息帶來的負(fù)面影響，在本節(jié)中，我們將 KC-EGPR的預(yù)測性能與基于全局信息的 GPR模型的預(yù)測性能進(jìn)行比較。為了便于描述，基于全局信息的GPR被命名為GL-GPR。此外，需要注意的是KC-EGPR預(yù)測性能與簇的數(shù)目K有著直接的關(guān)聯(lián)。在本文中，我們以1為步長從2到8調(diào)整K值，選擇使KC-EGPR性能達(dá)到最佳的K為最終的參數(shù)。表 2描述了 KCEGPR和GL-GPR在三組數(shù)據(jù)集下的預(yù)測性能。

表2 KC-EGPR和GL-GPR在DS1、DS2和DS3下的預(yù)測性能Tab. 2 Performance comparisons between KC-EGPR and GL-GPR on DS1, DS2, and DS3

從表2中可以觀察到，對于每組鉛酸電池?cái)?shù)據(jù)集，KC-EGPR的 MAPE和 RMSE的值均低于GL-GPR 對應(yīng)的值，這意味著 KC-EGPR的預(yù)測性能均優(yōu)于GL-GPR的預(yù)測性能。具體而言，從MAPE的角度分析，KC-EGPR的MAPE值相比于GL-GPR在 DS1、DS2和 DS3上分別降低了 2.95%、2.29%和1.98%。從RMSE的角度分析，KC-EGPR的RMSE值相比于GL-GPR在三組鉛酸電池?cái)?shù)據(jù)集上平均降低了0.7900，其中最大的降幅發(fā)生在DS3數(shù)據(jù)集上，達(dá)到了0.8712。

4 結(jié)束語

準(zhǔn)確地預(yù)測鉛酸電池的荷電狀態(tài)（SOC）能夠確保電池正常供電，對其所在系統(tǒng)的安全性和可靠性有著重要的意義。為了提高預(yù)測準(zhǔn)確率，本文提出了一種基于 K均值聚類的高斯過程回歸集成（KC-EGPR）算法。該算法的核心思想是采用聚類算法生成原始訓(xùn)練集的多個(gè)包含局部數(shù)據(jù)信息的子集，再利用集成學(xué)習(xí)理論將數(shù)據(jù)層融合問題轉(zhuǎn)化為決策層融合問題。相比于現(xiàn)有的模型，本文提出的算法能夠更有效地預(yù)測鉛酸電池的SOC，具有廣闊的應(yīng)用前景。

[1] Hua Chih Chiang, Lin Meng Yu. A study of charging control of lead-acid battery for electric vehicles[C]. Proceedings of the IEEE International Symposium on Industrial Electronics,2000, 1(1): 135-140.

[2] Vasebi Amir, Partovibakhsh Maral, Bathaee S. Mohammad Taghi. A novel combined battery model for state-of-charge estimation in lead-acid batteries based on extended Kalman filter for hybrid electric vehicle applications[J]. Journal of Power Sources, 2007, 174(1): 30-40.

[3] Chen Su, Kang Chunjian, Zhang Zhang, et al. A method for SOC estimation for lead-acid battery based on multi-model adaptive extended Kalman filtering estimation[C]. Proceedings of the Industrial Electronics Society, 2016, 1(1): 18-24.

[4] Gong Qing Zhu, Zhang Yu Ping. Application of LS-SVM in battery's SOC estimation in telecommunication power systems[J].Electrical Measurement & Instrumentation, 2010, 23(2): 145-152.

[5] Rasmussen Carl Edward. A unifying view of sparse approximate Gaussian process regression[M]. JMLR.org, 2005.

[6] 克里斯蒂亞尼尼. 支持向量機(jī)導(dǎo)論[M]. 電子工業(yè)出版社,2004.

[7] 毛健, 趙紅東, 姚婧婧. 人工神經(jīng)網(wǎng)絡(luò)的發(fā)展及應(yīng)用[J].電子設(shè)計(jì)工程, 2011, 19(24): 62-65.

[8] Hartigan J. A, Wong M. A. Algorithm as 136: a k-means clustering algorithm[J]. Journal of the Royal Statistical Society, 1979, 28(1): 100-108.

[9] 孫斌, 姚海濤, 劉婷. 基于高斯過程回歸的短期風(fēng)速預(yù)測[J]. 中國電機(jī)工程學(xué)報(bào), 2012, 32(29): 104-109.

[10] Freund Yoav, Schapire Robert E. Experiments with a new boosting algorithm[J]. 1996, 13(1): 148-156.