梁偉平, 牛博通
(華北電力大學(xué) 控制與計(jì)算機(jī)工程學(xué)院,河北 保定 071003)
煤的發(fā)熱量是電站系統(tǒng)的重要參數(shù),尤其在電站鍋爐熱平衡、熱效率計(jì)算、確定最佳的風(fēng)煤比、以及估計(jì)燃燒是否達(dá)到理論溫度等過(guò)程中,煤質(zhì)發(fā)熱量更是不可缺少[1-2]。隨著智能算法的發(fā)展,神經(jīng)網(wǎng)絡(luò)開(kāi)始用于建立煤質(zhì)發(fā)熱量的軟測(cè)量模型[3~5],但是神經(jīng)網(wǎng)絡(luò)算法訓(xùn)練時(shí)間長(zhǎng)、易陷入局部最優(yōu)點(diǎn)等缺點(diǎn)限制了其效果。隨后基于結(jié)構(gòu)風(fēng)險(xiǎn)最小化準(zhǔn)則的支持向量機(jī)算法被引入[6~9],解決了神經(jīng)網(wǎng)絡(luò)算法的問(wèn)題,但隨著DCS和存儲(chǔ)技術(shù)的發(fā)展,數(shù)據(jù)規(guī)模不斷增大,一般算法直接處理大量數(shù)據(jù)會(huì)導(dǎo)致建模速度慢,若從大規(guī)模數(shù)據(jù)中選取小樣本,會(huì)因憑借人為因素和經(jīng)驗(yàn)主義影響模型的精度和泛化能力。另外模型的特征變量是決定一個(gè)模型是否準(zhǔn)確和簡(jiǎn)單的重要因素,對(duì)特征變量進(jìn)行分析,選出最優(yōu)變量集會(huì)使模型更加準(zhǔn)確和簡(jiǎn)單,有利于提高泛化能力[10]。
本文采用Tsang等提出的CVR算法[11],該方法利用一種計(jì)算幾何中的最小閉包球理論(Minimum Enclosing Ball,MEB)與支持向量機(jī)(Support Vector Machine,SVM)相結(jié)合,使算法的計(jì)算復(fù)雜度和樣本容量成正比,空間復(fù)雜度與樣本容量無(wú)關(guān)[12],相比較其它SVM方法,當(dāng)建模數(shù)據(jù)增加時(shí),因受二次規(guī)劃(Quadratic Programming,QP)問(wèn)題的影響,模型的計(jì)算時(shí)間會(huì)隨之增加,而CVR隨著建模數(shù)據(jù)的增加,計(jì)算精度會(huì)有所改善(與大樣本中是否包含更多的信息有關(guān)),但計(jì)算時(shí)間的增長(zhǎng)遠(yuǎn)小于其它算法,解決了大規(guī)模數(shù)據(jù)下建模耗時(shí)的問(wèn)題,彌補(bǔ)了建模需要根據(jù)人為經(jīng)驗(yàn)選取小樣本導(dǎo)致模型不準(zhǔn)確的缺點(diǎn);在特征變量選擇上,利用R. J. May等提出的PMI方法分析特征變量[13],該方法基于信息熵理論中的互信息(Mutual Information,MI),通過(guò)條件期望消除了變量之間的聯(lián)系,消除了輸入變量的耦合關(guān)系對(duì)MI算法的影響[14],可以建立最經(jīng)濟(jì)簡(jiǎn)單、最大限度快速的軟測(cè)量模型。
信息熵的概念由Shannon在1948年提出,是用來(lái)量化描述變量所帶信息量的方法,MI又在此基礎(chǔ)上進(jìn)一步的度量了一個(gè)變量中含有的關(guān)于另一個(gè)變量的信息[15],R.Battiti等就基于互信息對(duì)神經(jīng)網(wǎng)絡(luò)的特征變量進(jìn)行了篩選[16],變量之間可能的耦合關(guān)系使MI計(jì)算受到了影響。PMI方法引入條件期望,將處理后的變量再得到互信息值,解決了這一問(wèn)題。
信息熵的具體計(jì)算公式為:
(1)
式中:pi為在各個(gè)取值下概率分布,n為數(shù)據(jù)個(gè)數(shù)。如果存在有向量X與Y相關(guān),pij為X和Y的聯(lián)合分布概率,則二維聯(lián)合熵定義為:
(2)
那么互信息的公式為:
(3)
由于一般樣本數(shù)據(jù)的概率分布未知,所以需要采用概率密度估計(jì)的方法來(lái)解決,因此實(shí)際情況下的公式變?yōu)槿缦拢?/p>
(4)
式中:xi,yi分別為X,Y的第i個(gè)取值,f為基于n個(gè)樣本數(shù)據(jù)的概率估計(jì)密度函數(shù)。
在MI的計(jì)算中,若輸入X,Z之間有耦合,那么I(Y,Z) 的值會(huì)大于實(shí)際,所以用條件期望mX(Z) 和mY(Z)消除Z的影響,如下:
(5)
U=X-mX(Z)
(6)
V=Y-mY(Z)
(7)
式中:zi為Z的第i個(gè)取值。
那么X,Y的偏互信息可記為:
IPMI(X,Y)=IPMI(U,V)
(8)
概率密度估計(jì)方法采用核密度估計(jì),這里核函數(shù)采取高斯核函數(shù),那么概率密度函數(shù)的估計(jì)公式為[17]:
(9)
式中:d為X的維數(shù);∑為X的協(xié)方差;h為帶寬;‖x-xi‖為馬氏距離。
馬氏距離公式為:
‖x-xi‖=(x-xi)TΣ-1(x-xi)
(10)
帶寬的選擇利用高斯帶寬,公式為:
(11)
式中:σ為樣本標(biāo)準(zhǔn)差。
結(jié)束條件采用赤池信息量準(zhǔn)則(Akaike Information Criterion, AIC),其公式為[18]:
(12)
式中:ri為根據(jù)已選變量計(jì)算的Y回歸殘差;p為已選變量個(gè)數(shù)。隨選出變量的變化,TAIC逐漸減小,當(dāng)不再減少時(shí)終止篩選。
最終的PMI變量選擇的算法流程如下:
1)初始化S為空集,C不為空集;
2)計(jì)算C中每個(gè)變量與Y的互信息,將互信息最大的Cs移入S,并更新C;
3)若C不為空,對(duì)每一個(gè)Cj∈C,計(jì)算v=Cj-mCj(S)和計(jì)算u=Cj-mCj(S);
4)計(jì)算I(u,v),選取使I(u,v)最大的Cs,用Cs計(jì)算AIC;
5)若AIC減小,則將Cs移入S,返回步驟3,否則終止。
(1)近似MEB與CVR
給定一點(diǎn)集S={x1,…,xm},xi∈Rd,那么集合S的最小閉包球是指包含集合S中所有數(shù)據(jù)點(diǎn)的最小球,表示為MEB(S)。傳統(tǒng)最小閉包球算法不能有效的解決d大于30的問(wèn)題,因此提出了一種快速近似算法,任意給定ε>0,如果存在R≤rMEB(S)并且S∈B(c,(1+ε)),則稱B(c,(1+ε))為B(c,R)的(1+ε)-近似,如下圖1所示。
圖1 近似閉包球
解決這個(gè)子集的問(wèn)題也能得到近似的估計(jì)和正確的結(jié)果[19],并且最終的核心集中的數(shù)量與迭代的次數(shù)無(wú)關(guān),主要與ε有關(guān)。
對(duì)于最小閉包球我們可以得到式(13):
minR2:‖φ(xi)-c‖2≤R2,i=1,…,m
(13)
φ(x)是由核函數(shù)誘導(dǎo)的特征映射函數(shù),得到上式的對(duì)偶形式如式(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)]是核矩陣。
在支持向量回歸中,有訓(xùn)練集合{zi=(xi,yi)},其中xi為輸入變量,yi為輸出變量。經(jīng)過(guò)核函數(shù)處理后的線性擬合方為:f(x)=ωTφ(xi)+b,采用ε-敏感損失函數(shù),那么最優(yōu)化的公式可以表示為式(15)[20]:
(15)
這里μ是控制損失函數(shù)參數(shù)的,得到對(duì)應(yīng)的對(duì)偶的矩陣形式如式(16):
(16)
(17)
但是得到的如式(16)的二次規(guī)劃問(wèn)題與最小閉包球問(wèn)題的式(14)的形式不同。
(2)中心約束MEB問(wèn)題
現(xiàn)在我們對(duì)每一個(gè)φ(xi)增加一個(gè)額外的項(xiàng)δi∈R,如式:φ(xi)~[φ(xi)T,δi]同時(shí)約束球心的最后一維變?yōu)椋篶~[cT,0]T那么問(wèn)題原始式變?yōu)槿缦?
minR2:‖φ(xi)-c‖2+δ2≤R2
(18)
maxαT(diag(K)+δ)-αTKα
s.t.α≥0,αT1=1
(19)
那么只要有最優(yōu)的,就可以計(jì)算出
(20)
(21)
因?yàn)棣罷1=1,那么對(duì)于任意η∈R加入不會(huì)影響α的結(jié)果,所以有[21]:
當(dāng)定義:
(22)
(23)
(24)
選擇適當(dāng)?shù)摩潜WC式(11)中的δ≥0,那么式(24)變?yōu)椋?/p>
(25)
式(25)與式(14)有著相同的形式,因此利用中心約束最小閉包球問(wèn)題,就可以及將ε-不敏感損失函數(shù)的支持向量機(jī)回歸問(wèn)題轉(zhuǎn)化為最小閉包球問(wèn)題。
(3) CVR算法流程
CVR算法步驟如下表所示:
①初始化選擇S0,C0,R0;
②在第t次迭代中,如果沒(méi)有φ(zi)在最小閉包球B(Ct,(1+ε)Rt)外,就終止訓(xùn)練;
③找到距離Ct最遠(yuǎn)的φ(zi),使St+1=St∪φ(zi);
④得到最新的閉包球集,計(jì)算新的Ct,Rt;
⑤迭代次數(shù)t=t+1,并返回第二步。
核心集中的向量為核心向量,通過(guò)這些核心向量可以得到原SVM問(wèn)題的解。
模型所用數(shù)據(jù)來(lái)自某摻煤燃燒電廠的煤質(zhì)分析報(bào)告,數(shù)據(jù)包括:全水分Mt(%)、空氣干燥基水分Mad(%),空氣干燥基灰分Aad(%)、空氣干燥基揮發(fā)分Vad(%)、收到基低位發(fā)熱量Qnet, ar(MJ/kg)、收到基全硫St,ar(%)、空干基全硫St, ad(%)等7維特征變量,為了去除數(shù)據(jù)噪聲,利用五點(diǎn)三次滑動(dòng)平均法對(duì)數(shù)據(jù)進(jìn)行降噪處理,部分低位發(fā)熱量去噪前后的對(duì)比分布如圖2所示:
圖2 低位發(fā)熱量去噪前后對(duì)比分布
圖2表明五點(diǎn)三次滑動(dòng)平均法有效的濾過(guò)了噪聲的干擾,另外不難發(fā)現(xiàn)摻配煤下的煤質(zhì)波動(dòng)非常頻繁,且波動(dòng)較大。然后再利用MATLAB中mapstd工具箱對(duì)數(shù)據(jù)進(jìn)行歸一化處理,消除數(shù)據(jù)單位的影響。
建立 CVR軟測(cè)量模型,特征變量的選取對(duì)于模型的計(jì)算復(fù)雜度、精度以及泛化能力都影響重大,特征變量之間有相關(guān)性和耦合性關(guān)系,可以利用PMI的方法分析預(yù)測(cè)模型的最優(yōu)變量子集。根據(jù)PMI變量篩選的步驟,得到TAIC的變化曲線如圖3所示。
特征變量初始特征集是根據(jù)最大互信息得出,各個(gè)特征變量與低位發(fā)熱量的互信值如表1所示。
圖3 PMI篩選TAIC曲線
%
全水分信息與低位發(fā)熱量互信值為0.334,是與低位發(fā)熱量互信值最大的變量,所以初始特征集中的第1個(gè)變量是全水分。在初始特征集的基礎(chǔ)上PMI計(jì)算得到的第2個(gè)變量是空干基全硫,第3個(gè)變量是收到基全硫,第4個(gè)變量是空干基揮發(fā)分,當(dāng)特征集加入第5個(gè)變量也就是空干基灰分時(shí)TAIC曲線由下降趨勢(shì)變?yōu)樯仙?,說(shuō)明雖然空干基灰分和空干基水分雖然與低位發(fā)熱量有互信息關(guān)系,但是在特征變量已經(jīng)包含前4個(gè)特征的情況下,空干基灰分和空干基水分與已選特征集之間的耦合性和相關(guān)性會(huì)導(dǎo)致將其加入特征集之后使預(yù)測(cè)效果變差。所以選取使特征變量與預(yù)測(cè)變量互信息最強(qiáng)并且特征變量間耦合性最小的全水分、空干基全硫、收到基全硫、空干基灰分4個(gè)變量作為模型的輸入特征變量。
因?yàn)楹诵闹С窒蛄繖C(jī)在處理大數(shù)據(jù)建模時(shí)對(duì)樣本沒(méi)有容量限制,可以利用最小閉包球的特點(diǎn)快速篩選支持向量,所以無(wú)需選取最優(yōu)樣本,因此選用去噪后的全部6 180 組煤質(zhì)數(shù)據(jù),以PMI選擇后的變量為特征變量,進(jìn)行仿真。
根據(jù)文獻(xiàn)[12],模型初始參數(shù)C對(duì)預(yù)測(cè)精度與預(yù)測(cè)時(shí)間影響較小,一般可取105,膨脹系數(shù)mu則對(duì)模型影響較大,對(duì)mu利用量子遺傳算法進(jìn)行尋優(yōu),選擇初始種群為30,轉(zhuǎn)角步長(zhǎng)為0.02,變異概率為0.1,迭代次數(shù)為30,得到當(dāng)為全部輸入變量時(shí)最優(yōu)參數(shù)為mu=2×10-6和當(dāng)變量為選擇之后時(shí)的最優(yōu)參數(shù)為mu=5×10-7。
初始拉格朗日參數(shù)選為1/2,剩下為0。初始核心集是先找到任意一點(diǎn)z0,再找到距離z0最遠(yuǎn)的點(diǎn)za, 距離za最遠(yuǎn)的點(diǎn)zb,za,zb為核心集,在第二步計(jì)算點(diǎn)到球心的距離時(shí)利用核函數(shù)避免了進(jìn)行φ(zl)的顯示計(jì)算如下式:
(26)
計(jì)算樣本與球心距離時(shí)采用Smola與Scholkopf提出的一種加速方法,文中提到當(dāng)樣本中隨機(jī)選取59個(gè)構(gòu)成子集時(shí),最遠(yuǎn)點(diǎn)的在其中的可能性為95%,可以大大降低時(shí)間復(fù)雜度[22],利用此方法使得隨機(jī)選中點(diǎn)不同,所以可以利用多次計(jì)算取最優(yōu)結(jié)果的方法。
為了評(píng)價(jià)模型的預(yù)測(cè)能力,選取均方差(MSE)、平均相對(duì)誤差(MRE)兩個(gè)指標(biāo)反應(yīng)模型的預(yù)測(cè)能力。
PMI-CVR算法的流程如圖4所示:
圖4 PMI-CVR算法流程圖
算法中TAIC是一個(gè)數(shù),將其初值設(shè)置為一個(gè)較大的數(shù)值,計(jì)算邊緣概率密度時(shí)的高斯帶寬取1,計(jì)算聯(lián)合概率密度時(shí)d取2,判斷是否為最優(yōu)結(jié)果時(shí),是利用遺傳算法滿足一定精度要求的條件,或著是當(dāng)遺傳算法的迭代次數(shù)達(dá)到設(shè)定的最大迭代次數(shù)時(shí)的最優(yōu)值。
將處理后的數(shù)據(jù)隨機(jī)選取6 130組作為訓(xùn)練樣本,剩下50組為測(cè)試樣本,按照上述方法建立5組隨機(jī)訓(xùn)練樣本和對(duì)應(yīng)的測(cè)試樣本,分別利用5組數(shù)據(jù)對(duì)CVR模型和PMI-CVR模型進(jìn)行預(yù)測(cè),模型的均方差、平均相對(duì)誤差以及計(jì)算時(shí)間取 5次結(jié)果的平均值,如表2所示:
表2 PMI-CVR與CVR對(duì)比
根據(jù)表2數(shù)據(jù)發(fā)現(xiàn),采用PMI的方法進(jìn)行變量篩選后,預(yù)測(cè)更加平穩(wěn), PMI-CVR模型預(yù)測(cè)的平均相對(duì)誤差有所降低,說(shuō)明經(jīng)過(guò)PMI對(duì)有耦合性的特征變量進(jìn)行剔除后,模型更加準(zhǔn)確,泛化能力強(qiáng),同時(shí)輸入變量的減少使模型更加簡(jiǎn)單,計(jì)算時(shí)間變短,計(jì)算更快。
兩種方法中預(yù)測(cè)最好的一組預(yù)測(cè)結(jié)果分別如圖5、圖6所示。
圖5 CVR低位發(fā)熱量預(yù)測(cè)
圖6 PMI-CVR低位發(fā)熱量預(yù)測(cè)
利用經(jīng)過(guò)PMI選擇后的特征變量為輸入特征變量,分別采用CVM和LSSVM兩種算法進(jìn)行建模。建模的樣本數(shù)量以1 000為步距,分別建立 1 000~6 000 6組數(shù)據(jù)規(guī)模下的預(yù)測(cè)模型,每種樣本規(guī)模下分別進(jìn)行5次建模,得到的5組數(shù)據(jù)的平均值作為該樣本規(guī)模下的結(jié)果, 兩種方法的誤差趨勢(shì)和建模時(shí)間上的趨勢(shì),隨著樣本的增加兩種方法建模誤差對(duì)比圖如圖7所示。
圖7 CVR與LSSVM建模誤差對(duì)比
隨著樣本的增加兩種方法建模計(jì)算時(shí)間如圖9所示:
圖8 CVR與LSSVM建模時(shí)間對(duì)比
對(duì)比圖7與圖8,當(dāng)建模樣本規(guī)模增加時(shí),PMI-CVR與LSSVM的建模相對(duì)誤差相差不到0.005,但是建模時(shí)間上LSSVM耗時(shí)增加明顯,而PMI-CVR建模的時(shí)間變化很小。
本文利用偏互信息與核心支持向量機(jī)相結(jié)合,在大規(guī)模煤質(zhì)數(shù)據(jù)下建立了PMI-CVR模型,相比未經(jīng)過(guò)PMI處理的CVR模型,PMI-CVR更加簡(jiǎn)單,計(jì)算結(jié)果泛化能力更強(qiáng);在相同的輸入特征變量下,隨著樣本數(shù)據(jù)量增加的PMI-CVR模型比LSSVM模型計(jì)算更快,同時(shí)在準(zhǔn)確度上二者相差不大,所以PMI-CVR對(duì)大數(shù)據(jù)的處理能力更強(qiáng)。最終PMI-CVR的低位發(fā)熱量預(yù)測(cè)模型相對(duì)誤差為0.025,計(jì)算時(shí)間為0.272 s,證明該模型在大規(guī)模數(shù)據(jù)下更加有優(yōu)勢(shì)。