基于主成分回歸的自然發(fā)火程度判識模型

杜 斌

（1.煤炭科學(xué)研究總院，北京100013；2.中煤科工集團(tuán)沈陽研究院有限公司，遼寧 撫順113122；3.煤礦安全技術(shù)國家重點實驗室，遼寧 撫順113122）

采空區(qū)遺煤自然發(fā)火是影響我國礦井安全生產(chǎn)的重大危險源之一。當(dāng)工作面上隅角及回風(fēng)流中的CO 體積分?jǐn)?shù)出現(xiàn)上升趨勢時，需要對自然發(fā)火程度進(jìn)行辨識及預(yù)測，以便及時采取有效應(yīng)對措施；同樣，封閉火區(qū)啟封前也需要對火區(qū)內(nèi)遺煤的自然發(fā)火程度進(jìn)行判識。目前主要是建立自然發(fā)火標(biāo)志性氣體預(yù)測預(yù)報指標(biāo)體系[1]，以CO 及其派生指標(biāo)、C2H4、C2H2為主指標(biāo)，以鏈烷比和烯烷比以及溫度為輔助指標(biāo)，通過氣體分析法對煤自然發(fā)火程度進(jìn)行判識，賈寶山[2]和董紹樸[3]等通過程序升溫實驗和主成分分析的方法對自然發(fā)火標(biāo)志性氣體進(jìn)行了優(yōu)選。這些研究成果都表明了自然發(fā)火標(biāo)志性氣體可以用于判識煤的氧化程度。林棉金[4]研究了自然發(fā)火過程中生成氣體的分階段宏觀特性，楊波等[5-6]研究了不同氧化程度煤在氧化蓄熱升溫過程中的細(xì)觀孔隙-裂隙結(jié)構(gòu)演變規(guī)律，王彩萍等[7-8]研究了煤氧化過程中微觀分子結(jié)構(gòu)和官能團(tuán)的變化規(guī)律。這些研究都指出了煤氧化過程中的分階段特性，但是對劃分各階段的臨界點及劃分方法尚未提出。尤其是目前尚未提出煤氧化程度的定量化評價方法。為此通過對煤樣進(jìn)行程序升溫實驗，借助于多元統(tǒng)計方法和SPSS 統(tǒng)計分析軟件，對實驗數(shù)據(jù)進(jìn)行相關(guān)性分析、主成分分析及回歸，構(gòu)建自然發(fā)火程度預(yù)測模型。

1 煤樣程序升溫實驗

1.1 煤樣的制備與實驗設(shè)計

煤樣采用變質(zhì)程度較低的褐煤，采集自麥垛山煤礦6#煤層工作面。按照國標(biāo)[9-10]的要求，對煤樣進(jìn)行粉碎、研磨、篩分，篩選出粒徑小于100 目（0.15 mm）的煤樣封存于密封袋中。對煤樣進(jìn)行工業(yè)分析，煤樣工業(yè)分析結(jié)果見表1。

表1 煤樣工業(yè)分析結(jié)果Table 1 Industrial analysis results of the coal sample%

實驗過程為：①取1 g 制備好的煤樣，放入程序升溫實驗裝置的煤樣罐中，連接好進(jìn)氣和出氣管路；②升溫速率設(shè)置：25～80 ℃為0.5 ℃/min，80～200 ℃為1.0 ℃/min，200～300 ℃為2.0 ℃/min；③供氣流量100 mL/min，用氣相色譜儀分析實驗生成的氣體組分及體積分?jǐn)?shù)，每隔20 min 記錄1 次實驗結(jié)果。

1.2 實驗結(jié)果

程序升溫實驗氣體產(chǎn)物分析結(jié)果見表2。

表2 程序升溫實驗氣體產(chǎn)物分析結(jié)果Table 2 Analysis results of gases generated in the programmed heating experiment

從實驗結(jié)果可以發(fā)現(xiàn)，從448 ℃開始，O2體積分?jǐn)?shù)出現(xiàn)斷崖式下降，說明此時已經(jīng)進(jìn)入燃燒期，煤樣明火燃燒時消耗了大量O2；同時煤體溫度先上升至531 ℃后逐漸下降，說明煤樣明火燃燒后進(jìn)入陰燃階段，此時耗氧量逐漸減小，氧體積分?jǐn)?shù)逐漸升高。待O2體積分?jǐn)?shù)上升回復(fù)至與空氣中O2體積分?jǐn)?shù)相同時，說明煤樣完全熄滅。

N2體積分?jǐn)?shù)在實驗過程中基本無變化，這主要是因為N2主要來自于外部環(huán)境，其在實驗過程中既不是反應(yīng)物也不是生成物，所以N2體積分?jǐn)?shù)主要是受氣體體積相對比的影響而上下浮動。

除了O2、N2之外，CxHy和CxOy的變化趨勢大致相同，具有較強的相關(guān)性。

2 相關(guān)性分析

目前對自然發(fā)火程度進(jìn)行預(yù)測主要是通過氣體分析法和測溫法，由于測溫法主要是借助溫度傳感器等，無法進(jìn)行非接觸式測量，因此氣體分析法使用較為廣泛。

但是在程序升溫過程中產(chǎn)生的其他CxHy和CxOy氣體與溫度的相關(guān)程度大小尚未被定量描述，因此可以將各種氣體看作隨機(jī)向量，通過相關(guān)性分析，篩選出與溫度相關(guān)系數(shù)大的氣體作為描述自然發(fā)火進(jìn)程的標(biāo)志性氣體。

標(biāo)志性氣體主要考慮從與溫度具有較強相關(guān)性的CxHy和CxOy中選取，但是CxHy和CxOy的體積分?jǐn)?shù)量級變化范圍較大，涵蓋了10-1～10-6，需要將實驗結(jié)果進(jìn)行標(biāo)準(zhǔn)化變換。標(biāo)準(zhǔn)化變換如下[11]：

計算相關(guān)系數(shù)矩陣R，由此繪制的相關(guān)系數(shù)矩陣圖如圖1。

圖1 相關(guān)系數(shù)矩陣圖Fig.1 Correlation coefficient matrix

圖1 中藍(lán)色代表負(fù)相關(guān)，紅色代表正相關(guān)。圓圈面積的大小表示了相關(guān)性的強弱。與溫度相關(guān)系數(shù)最大的是CO，溫度與O2、N2呈負(fù)相關(guān)，與其他氣體具有較強的正相關(guān)性。

3 主成分分析

CO 在自然發(fā)火過程中出現(xiàn)時間最早且貫穿整個自然發(fā)火進(jìn)程，便于檢測，因此常常作為預(yù)測預(yù)報自然發(fā)火的標(biāo)志性氣體。但是由于自然發(fā)火過程中產(chǎn)生的各氣體與溫度T 都具有較強的相關(guān)性，使用單一氣體難以全面刻畫自然發(fā)火進(jìn)程，因此需要借助主成分分析，在盡量不丟失描述自然發(fā)火進(jìn)程信息的情況下，將具有相關(guān)性的各氣體轉(zhuǎn)化為線性無關(guān)的主成分。

3.1 主成分分析的基本思想

主成分分析就是設(shè)法將原來具有相關(guān)性的指標(biāo)重新組合成一組新的互相無關(guān)的幾個綜合指標(biāo)來代替原來指標(biāo)，同時根據(jù)實際需要，從中選取幾個較少的綜合指標(biāo)盡可能多地反映原來指標(biāo)的信息[12]。通常數(shù)學(xué)上的處理就是將原來p 個指標(biāo)做線性組合，作為新的綜合指標(biāo)。

式中：Xi（i=1，2，…，p）為隨機(jī)向量，分別代表程序升溫實驗中出現(xiàn)的氣體；Fi（i=1，2，…，p）為主成分；ai1、ai2、…、aip為系數(shù)。

求取系數(shù)向量（ai1，ai2，…，aip）主要是通過相關(guān)系數(shù)矩陣求特征值和特征向量，單位特征向量即為系數(shù)向量。按照特征值取值大小確定主成分的排名。如果前k 個主成分的貢獻(xiàn)率達(dá)到了85%，就認(rèn)為前k個主成分基本包含了全部測量指標(biāo)所具有的信息。

3.2 主成分分析的結(jié)果

借助于SPSS 軟件[13-14]的主成分分析功能，計算得到相關(guān)系數(shù)矩陣的特征值為λ1=6.373，λ2=0.562，λ3=0.049，λ4=0.012，λ5=0.004，λ6=7.837 ×10-5，λ7=3.030×10-5。其中λ1/Σλi=91.038%，因此可以認(rèn)為第1主成分基本包含了描述自然發(fā)火進(jìn)程的全部信息。第1 主成分的系數(shù)向量為特征值λ1對應(yīng)的單位特征向量，記為a1＝［0.145，0.141，0.154，0.156，0.150，0.151，0.151］。

所以，第1 主成分F1為：

比較系數(shù)向量的p 個分量，可以發(fā)現(xiàn)a1p都在0.150 附近浮動，說明各氣體Xi對第1 主成分的貢獻(xiàn)相近。

4 主成分回歸

4.1 回歸模型的建立

為了對自然發(fā)火程度進(jìn)行預(yù)測，還需要將溫度T 作為因變量，各氣體Xi作為自變量，通過主成分回歸分析，構(gòu)建溫度T 和各氣體Xi之間的回歸預(yù)測模型。之前通過主成分分析，已經(jīng)在盡量保留原來指標(biāo)信息的情況下將具有相關(guān)性的p（p=7）維氣體Xi降維，得到了m（m=1）維線性無關(guān)的主成分F。這樣就簡化了回歸方程的結(jié)構(gòu)，只需要構(gòu)建溫度T 和主成分F 之間的回歸預(yù)測模型，然后將主成分F 與各氣體Xi的表達(dá)式進(jìn)行回代，就可以得到煤的自然發(fā)火程度預(yù)測模型。

首先，在XOY 坐標(biāo)系中繪制散點圖，觀察溫度T 與主成分F 之間的關(guān)系，對T 與F 可能符合的回歸模型進(jìn)行初步辨識。主成分F1與因變量T 的散點圖如圖2。

圖2 主成分F1 與因變量T 的散點圖Fig.2 Scatter diagram of principal component F1 and dependent variable T

擬合優(yōu)度R2=0.886，說明擬合效果良好。

然后將式（7）代入上述對數(shù)型回歸預(yù)測模型，可以得到：

式中：X1為CO；X2為CO2；X3為CH4；X4為C2H6；X5為C2H4；X6為C3H8；X7為C3H6。

這樣就得到了煤體溫度T 和各氣體Xi之間的自然發(fā)火預(yù)測模型。

4.2 回歸模型的驗證

為了驗證對數(shù)型模型的有效性，借助SPSS 統(tǒng)計分析軟件的曲線估計功能，分別選取線性、二次、三次、S 型、指數(shù)型等模型，計算模型的擬合優(yōu)度R2并和對數(shù)型模型的擬合優(yōu)度R2見表3。

表3 預(yù)測模型擬合優(yōu)度R2Table 3 Goodness of fit（R2）in prediction models

通過對比可以發(fā)現(xiàn)，曲線估計功能擬合效果最好的為三次曲線，但其擬合優(yōu)度依然小于對數(shù)型曲線。所以可以將對數(shù)型曲線作為自然發(fā)火程度的預(yù)測模型。

5 結(jié) 語

通過開展程序升溫實驗，分析了自然發(fā)火生成的各氣體與溫度之間的相關(guān)性，發(fā)現(xiàn)單一氣體難以全面刻畫自然發(fā)火的程度。借助SPSS 統(tǒng)計分析軟件，分析了各氣體之間的相關(guān)性大小，通過主成分分析、主成分回歸等方法，發(fā)現(xiàn)溫度T 和主成分F 之間大致呈對數(shù)型函數(shù)關(guān)系，構(gòu)建了以溫度T 為因變量、氣體Xi為自變量的自然發(fā)火程度預(yù)測數(shù)學(xué)模型，可以為今后自然發(fā)火程度的預(yù)測起到一定的借鑒意義。