基于偏互信息與核心向量機的煤質大數據預測

梁偉平， 牛博通

(華北電力大學 控制與計算機工程學院，河北 保定 071003)

0 引 言

煤的發熱量是電站系統的重要參數，尤其在電站鍋爐熱平衡、熱效率計算、確定最佳的風煤比、以及估計燃燒是否達到理論溫度等過程中，煤質發熱量更是不可缺少[1-2]。隨著智能算法的發展，神經網絡開始用于建立煤質發熱量的軟測量模型[3～5]，但是神經網絡算法訓練時間長、易陷入局部最優點等缺點限制了其效果。隨后基于結構風險最小化準則的支持向量機算法被引入[6～9]，解決了神經網絡算法的問題，但隨著DCS和存儲技術的發展，數據規模不斷增大，一般算法直接處理大量數據會導致建模速度慢，若從大規模數據中選取小樣本，會因憑借人為因素和經驗主義影響模型的精度和泛化能力。另外模型的特征變量是決定一個模型是否準確和簡單的重要因素，對特征變量進行分析，選出最優變量集會使模型更加準確和簡單，有利于提高泛化能力[10]。

本文采用Tsang等提出的CVR算法[11]，該方法利用一種計算幾何中的最小閉包球理論(Minimum Enclosing Ball，MEB)與支持向量機(Support Vector Machine，SVM)相結合，使算法的計算復雜度和樣本容量成正比，空間復雜度與樣本容量無關[12]，相比較其它SVM方法，當建模數據增加時，因受二次規劃(Quadratic Programming，QP)問題的影響，模型的計算時間會隨之增加，而CVR隨著建模數據的增加，計算精度會有所改善(與大樣本中是否包含更多的信息有關)，但計算時間的增長遠小于其它算法，解決了大規模數據下建模耗時的問題，彌補了建模需要根據人為經驗選取小樣本導致模型不準確的缺點；在特征變量選擇上，利用R. J. May等提出的PMI方法分析特征變量[13]，該方法基于信息熵理論中的互信息(Mutual Information，MI),通過條件期望消除了變量之間的聯系,消除了輸入變量的耦合關系對MI算法的影響[14]，可以建立最經濟簡單、最大限度快速的軟測量模型。

1 PMI與CVR原理

1.1 PMI

信息熵的概念由Shannon在1948年提出，是用來量化描述變量所帶信息量的方法,MI又在此基礎上進一步的度量了一個變量中含有的關于另一個變量的信息[15]，R.Battiti等就基于互信息對神經網絡的特征變量進行了篩選[16]，變量之間可能的耦合關系使MI計算受到了影響。PMI方法引入條件期望，將處理后的變量再得到互信息值，解決了這一問題。

信息熵的具體計算公式為：

(1)

式中：pi為在各個取值下概率分布，n為數據個數。如果存在有向量X與Y相關，pij為X和Y的聯合分布概率，則二維聯合熵定義為:

(2)

那么互信息的公式為：

(3)

由于一般樣本數據的概率分布未知，所以需要采用概率密度估計的方法來解決，因此實際情況下的公式變為如下：

(4)

式中：xi,yi分別為X,Y的第i個取值，f為基于n個樣本數據的概率估計密度函數。

在MI的計算中，若輸入X，Z之間有耦合，那么I(Y,Z) 的值會大于實際，所以用條件期望mX(Z) 和mY(Z)消除Z的影響，如下:

(5)

U=X-mX(Z)

(6)

V=Y-mY(Z)

(7)

式中：zi為Z的第i個取值。

那么X,Y的偏互信息可記為：

IPMI(X,Y)=IPMI(U,V)

(8)

概率密度估計方法采用核密度估計，這里核函數采取高斯核函數，那么概率密度函數的估計公式為[17]:

(9)

式中：d為X的維數；∑為X的協方差；h為帶寬；‖x-xi‖為馬氏距離。

馬氏距離公式為：

‖x-xi‖=(x-xi)TΣ-1(x-xi)

(10)

帶寬的選擇利用高斯帶寬，公式為：

(11)

式中：σ為樣本標準差。

結束條件采用赤池信息量準則(Akaike Information Criterion, AIC),其公式為[18]：

(12)

式中：ri為根據已選變量計算的Y回歸殘差；p為已選變量個數。隨選出變量的變化，TAIC逐漸減小，當不再減少時終止篩選。

最終的PMI變量選擇的算法流程如下：

1)初始化S為空集,C不為空集；

2)計算C中每個變量與Y的互信息，將互信息最大的Cs移入S，并更新C；

3)若C不為空，對每一個Cj∈C，計算v=Cj-mCj(S)和計算u=Cj-mCj(S)；

4)計算I(u,v)，選取使I(u,v)最大的Cs，用Cs計算AIC；

5)若AIC減小，則將Cs移入S,返回步驟3，否則終止。

1.2 CVR

(1)近似MEB與CVR

給定一點集S={x1,…,xm},xi∈Rd，那么集合S的最小閉包球是指包含集合S中所有數據點的最小球，表示為MEB(S)。傳統最小閉包球算法不能有效的解決d大于30的問題，因此提出了一種快速近似算法，任意給定ε>0，如果存在R≤rMEB(S)并且S∈B(c,(1+ε)),則稱B(c,(1+ε))為B(c,R)的(1+ε)-近似，如下圖1所示。

圖1 近似閉包球

解決這個子集的問題也能得到近似的估計和正確的結果[19]，并且最終的核心集中的數量與迭代的次數無關，主要與ε有關。

對于最小閉包球我們可以得到式(13)：

minR2:‖φ(xi)-c‖2≤R2，i=1,…,m

(13)

φ(x)是由核函數誘導的特征映射函數，得到上式的對偶形式如式(14)：

maxαTdiag(K)-αTKαα≥0,αT1=1

(14)

(14)式中：α=[αi,…，αm]T是拉格朗日乘子,0=[0,…,0]T,1=[1,…,1]T，Km×m=[K(xi,xi)]=[φ(xi)Tφ(xi)]是核矩陣。

在支持向量回歸中，有訓練集合{zi=(xi,yi)}，其中xi為輸入變量，yi為輸出變量。經過核函數處理后的線性擬合方為：f(x)=ωTφ(xi)+b，采用ε-敏感損失函數，那么最優化的公式可以表示為式(15)[20]：

(15)

這里μ是控制損失函數參數的，得到對應的對偶的矩陣形式如式(16)：

(16)

(17)

但是得到的如式(16)的二次規劃問題與最小閉包球問題的式(14)的形式不同。

(2)中心約束MEB問題

現在我們對每一個φ(xi)增加一個額外的項δi∈R，如式：φ(xi)～[φ(xi)T,δi]同時約束球心的最后一維變為：c～[cT,0]T那么問題原始式變為如下:

minR2:‖φ(xi)-c‖2+δ2≤R2

(18)

maxαT(diag(K)+δ)-αTKα

s.t.α≥0,αT1=1

(19)

那么只要有最優的，就可以計算出

(20)

(21)

因為αT1=1，那么對于任意η∈R加入不會影響α的結果，所以有[21]：

當定義：

(22)

(23)

(24)

選擇適當的η保證式(11)中的δ≥0，那么式(24)變為：

(25)

式(25)與式(14)有著相同的形式，因此利用中心約束最小閉包球問題，就可以及將ε-不敏感損失函數的支持向量機回歸問題轉化為最小閉包球問題。

(3) CVR算法流程

CVR算法步驟如下表所示：

①初始化選擇S0,C0,R0;

②在第t次迭代中,如果沒有φ(zi)在最小閉包球B(Ct,(1+ε)Rt)外，就終止訓練;

③找到距離Ct最遠的φ(zi)，使St+1=St∪φ(zi);

④得到最新的閉包球集，計算新的Ct,Rt;

⑤迭代次數t=t+1，并返回第二步。

核心集中的向量為核心向量，通過這些核心向量可以得到原SVM問題的解。

2 數據預處理及變量分析

2.1 數據預處理

模型所用數據來自某摻煤燃燒電廠的煤質分析報告，數據包括：全水分Mt(%)、空氣干燥基水分Mad(%)，空氣干燥基灰分Aad(%)、空氣干燥基揮發分Vad(%)、收到基低位發熱量Qnet, ar(MJ/kg)、收到基全硫St,ar(%)、空干基全硫St, ad(%)等7維特征變量，為了去除數據噪聲，利用五點三次滑動平均法對數據進行降噪處理，部分低位發熱量去噪前后的對比分布如圖2所示：

圖2 低位發熱量去噪前后對比分布

圖2表明五點三次滑動平均法有效的濾過了噪聲的干擾，另外不難發現摻配煤下的煤質波動非常頻繁，且波動較大。然后再利用MATLAB中mapstd工具箱對數據進行歸一化處理，消除數據單位的影響。

2.2 PMI變量分析

建立 CVR軟測量模型，特征變量的選取對于模型的計算復雜度、精度以及泛化能力都影響重大，特征變量之間有相關性和耦合性關系，可以利用PMI的方法分析預測模型的最優變量子集。根據PMI變量篩選的步驟，得到TAIC的變化曲線如圖3所示。

特征變量初始特征集是根據最大互信息得出，各個特征變量與低位發熱量的互信值如表1所示。

圖3 PMI篩選TAIC曲線

%

全水分信息與低位發熱量互信值為0.334，是與低位發熱量互信值最大的變量，所以初始特征集中的第1個變量是全水分。在初始特征集的基礎上PMI計算得到的第2個變量是空干基全硫，第3個變量是收到基全硫，第4個變量是空干基揮發分，當特征集加入第5個變量也就是空干基灰分時TAIC曲線由下降趨勢變為上升，說明雖然空干基灰分和空干基水分雖然與低位發熱量有互信息關系，但是在特征變量已經包含前4個特征的情況下，空干基灰分和空干基水分與已選特征集之間的耦合性和相關性會導致將其加入特征集之后使預測效果變差。所以選取使特征變量與預測變量互信息最強并且特征變量間耦合性最小的全水分、空干基全硫、收到基全硫、空干基灰分4個變量作為模型的輸入特征變量。

3 實驗仿真

3.1 模型初始參數及構建

因為核心支持向量機在處理大數據建模時對樣本沒有容量限制，可以利用最小閉包球的特點快速篩選支持向量，所以無需選取最優樣本，因此選用去噪后的全部6 180 組煤質數據，以PMI選擇后的變量為特征變量，進行仿真。

根據文獻[12],模型初始參數C對預測精度與預測時間影響較小，一般可取105，膨脹系數mu則對模型影響較大，對mu利用量子遺傳算法進行尋優，選擇初始種群為30，轉角步長為0.02,變異概率為0.1，迭代次數為30，得到當為全部輸入變量時最優參數為mu=2×10-6和當變量為選擇之后時的最優參數為mu=5×10-7。

初始拉格朗日參數選為1/2，剩下為0。初始核心集是先找到任意一點z0，再找到距離z0最遠的點za, 距離za最遠的點zb，za,zb為核心集，在第二步計算點到球心的距離時利用核函數避免了進行φ(zl)的顯示計算如下式：

(26)

計算樣本與球心距離時采用Smola與Scholkopf提出的一種加速方法，文中提到當樣本中隨機選取59個構成子集時，最遠點的在其中的可能性為95%，可以大大降低時間復雜度[22]，利用此方法使得隨機選中點不同，所以可以利用多次計算取最優結果的方法。

為了評價模型的預測能力，選取均方差(MSE)、平均相對誤差(MRE)兩個指標反應模型的預測能力。

PMI-CVR算法的流程如圖4所示：

圖4 PMI-CVR算法流程圖

算法中TAIC是一個數，將其初值設置為一個較大的數值，計算邊緣概率密度時的高斯帶寬取1，計算聯合概率密度時d取2，判斷是否為最優結果時，是利用遺傳算法滿足一定精度要求的條件，或著是當遺傳算法的迭代次數達到設定的最大迭代次數時的最優值。

3.2 實驗與結果分析

將處理后的數據隨機選取6 130組作為訓練樣本，剩下50組為測試樣本，按照上述方法建立5組隨機訓練樣本和對應的測試樣本，分別利用5組數據對CVR模型和PMI-CVR模型進行預測，模型的均方差、平均相對誤差以及計算時間取 5次結果的平均值，如表2所示：

表2 PMI-CVR與CVR對比

根據表2數據發現，采用PMI的方法進行變量篩選后，預測更加平穩， PMI-CVR模型預測的平均相對誤差有所降低，說明經過PMI對有耦合性的特征變量進行剔除后，模型更加準確，泛化能力強，同時輸入變量的減少使模型更加簡單，計算時間變短，計算更快。

兩種方法中預測最好的一組預測結果分別如圖5、圖6所示。

圖5 CVR低位發熱量預測

圖6 PMI-CVR低位發熱量預測

3.3 不同算法對比

利用經過PMI選擇后的特征變量為輸入特征變量，分別采用CVM和LSSVM兩種算法進行建模。建模的樣本數量以1 000為步距，分別建立 1 000～6 000 6組數據規模下的預測模型，每種樣本規模下分別進行5次建模，得到的5組數據的平均值作為該樣本規模下的結果， 兩種方法的誤差趨勢和建模時間上的趨勢，隨著樣本的增加兩種方法建模誤差對比圖如圖7所示。

圖7 CVR與LSSVM建模誤差對比

隨著樣本的增加兩種方法建模計算時間如圖9所示：

圖8 CVR與LSSVM建模時間對比

對比圖7與圖8，當建模樣本規模增加時，PMI-CVR與LSSVM的建模相對誤差相差不到0.005，但是建模時間上LSSVM耗時增加明顯，而PMI-CVR建模的時間變化很小。

4 結論

本文利用偏互信息與核心支持向量機相結合，在大規模煤質數據下建立了PMI-CVR模型，相比未經過PMI處理的CVR模型，PMI-CVR更加簡單，計算結果泛化能力更強；在相同的輸入特征變量下，隨著樣本數據量增加的PMI-CVR模型比LSSVM模型計算更快，同時在準確度上二者相差不大，所以PMI-CVR對大數據的處理能力更強。最終PMI-CVR的低位發熱量預測模型相對誤差為0.025，計算時間為0.272 s，證明該模型在大規模數據下更加有優勢。