基于Python的瓦斯?jié)舛華RIMA預(yù)測(cè)模型構(gòu)建及其應(yīng)用

張 震，朱權(quán)潔,2，李青松，張爾輝，劉鴻偉

(1. 華北科技學(xué)院 安全工程學(xué)院，北京 東燕郊 065201；2.華北科技學(xué)院 河北省礦井災(zāi)害防治重點(diǎn)實(shí)驗(yàn)室，北京 東燕郊 065201；3.貴州省礦山安全科學(xué)研究院，貴陽(yáng) 550025)

0 引言

瓦斯災(zāi)害是我國(guó)煤礦生產(chǎn)中所面臨的第一大災(zāi)害，其引起的傷亡人數(shù)遠(yuǎn)多于其它災(zāi)害[1,2]。提前獲知瓦斯?jié)舛茸兓厔?shì)，能夠及時(shí)、有效的防治災(zāi)害發(fā)生。因此，開(kāi)展礦井瓦斯?jié)舛茸兓?guī)律的研究勢(shì)在必行。

近年來(lái)，隨著瓦斯災(zāi)害防治領(lǐng)域研究的不斷深化，我國(guó)在瓦斯?jié)舛劝l(fā)展變化規(guī)律研究方面取得了長(zhǎng)足進(jìn)步[3，4]。其中，黃東[5]等基于灰色系統(tǒng)理論的基礎(chǔ)上有效結(jié)合了自記憶原理，構(gòu)建了新型斯?jié)舛阮A(yù)測(cè)模型，精準(zhǔn)的預(yù)測(cè)出了瓦斯?jié)舛任磥?lái)變化規(guī)律的總體趨勢(shì)。趙亮[6]基于混沌理論重構(gòu)了瓦斯?jié)舛葧r(shí)間序列的相空間，建立了適用工作面的瓦斯?jié)舛阮A(yù)測(cè)模型。通過(guò)實(shí)際應(yīng)用表明，實(shí)際監(jiān)測(cè)的瓦斯?jié)舛扰c該模型預(yù)測(cè)的瓦斯?jié)舛染哂休^高的吻合度。姜雷[7]開(kāi)展了各種瓦斯?jié)舛扔绊懸蛩胤蔷€性特征研究，并以BP神經(jīng)網(wǎng)絡(luò)算法為核心算法建立了瓦斯?jié)舛阮A(yù)測(cè)模型，取得了較好的預(yù)測(cè)效果。由于單一的BP神經(jīng)網(wǎng)絡(luò)模型存在收斂速度慢、部分節(jié)點(diǎn)確定缺乏科學(xué)性的缺陷。為了提高預(yù)測(cè)精度，韓婷婷[8]等提出了一種采用馬爾科夫模型修正神經(jīng)網(wǎng)絡(luò)模型瓦斯?jié)舛任磥?lái)趨勢(shì)預(yù)測(cè)結(jié)果的方法，將瓦斯?jié)舛阮A(yù)測(cè)值與實(shí)測(cè)值的最大誤差降低到6%以下。但馬爾科夫算法存在計(jì)算量大等不足，降低了預(yù)測(cè)效率。為此，姚青華[9]等通過(guò)將遺傳算法與BP神經(jīng)網(wǎng)絡(luò)算法深度結(jié)合的方式，改進(jìn)了傳統(tǒng)BP神經(jīng)網(wǎng)絡(luò)系統(tǒng)預(yù)測(cè)模型迭代次數(shù)多等缺陷，降低了瓦斯?jié)舛阮A(yù)測(cè)模型的絕對(duì)誤差。由于瓦斯災(zāi)害具有破壞性、突發(fā)性，及時(shí)預(yù)測(cè)瓦斯?jié)舛茸兓厔?shì)可為煤礦提前采取安全防范措施提供依據(jù)。因此，張昭昭等基于動(dòng)態(tài)神經(jīng)網(wǎng)絡(luò)算法建立了瓦斯?jié)舛茸兓?guī)律實(shí)時(shí)預(yù)測(cè)模型，該模型根據(jù)瓦斯?jié)舛鹊膭?dòng)力學(xué)特性及時(shí)更新[10]，提高了瓦斯?jié)舛茸兓?guī)律的預(yù)測(cè)精度。由于實(shí)時(shí)采集的瓦斯?jié)舛葦?shù)據(jù)海量，對(duì)算法處理數(shù)據(jù)的要求高。劉曉悅[11]等針對(duì)海量的瓦斯?jié)舛葦?shù)據(jù)處理要求高的問(wèn)題，基于云計(jì)算，利用Elman神經(jīng)網(wǎng)絡(luò)模型開(kāi)展了煤礦瓦斯?jié)舛阮A(yù)測(cè)的可行性研究，實(shí)現(xiàn)了瓦斯數(shù)據(jù)的高效處理，快速的預(yù)測(cè)了瓦斯?jié)舛茸兓厔?shì)。

綜上所述，針對(duì)現(xiàn)有隨機(jī)過(guò)程預(yù)測(cè)算法復(fù)雜、難以掌握等問(wèn)題，基于Python語(yǔ)言開(kāi)展了瓦斯?jié)舛葧r(shí)間序列預(yù)測(cè)研究。本研究在前人研究的基礎(chǔ)上，探討了時(shí)間序列分析基本原理，建立了基于差分自回歸移動(dòng)平均(ARIMA)模型的瓦斯?jié)舛阮A(yù)測(cè)模型。通過(guò)ARIMA瓦斯?jié)舛阮A(yù)測(cè)模型分析、捕獲瓦斯?jié)舛葧r(shí)間序列的變化規(guī)律，將該時(shí)間序列的變化規(guī)律向前延伸，從而達(dá)到精準(zhǔn)預(yù)測(cè)瓦斯?jié)舛鹊哪康摹?/p>

1 Python語(yǔ)言簡(jiǎn)介

Python語(yǔ)言作為一款功能完備、運(yùn)行穩(wěn)定的計(jì)算機(jī)程序設(shè)計(jì)語(yǔ)言，其被廣泛應(yīng)用于Web和Internet開(kāi)發(fā)、人工智能及網(wǎng)絡(luò)爬蟲(chóng)等眾多領(lǐng)域。在科學(xué)計(jì)算和統(tǒng)計(jì)學(xué)方面，Python設(shè)計(jì)開(kāi)發(fā)有大量公開(kāi)源代碼的分析計(jì)算模塊，提供了分析計(jì)算結(jié)果可視化展現(xiàn)的命令窗口。Python還被稱為“膠水語(yǔ)言”，其能輕松結(jié)合其他語(yǔ)言的功能模塊協(xié)同工作，可擴(kuò)展性強(qiáng)。通過(guò)Python計(jì)算機(jī)語(yǔ)言建立的ARIMA(p,d,q)模型簡(jiǎn)單易懂、使用靈活，因此，選用Python語(yǔ)言作為工具建立礦井瓦斯?jié)舛茸兓?guī)律預(yù)測(cè)模型。

如圖1為構(gòu)建預(yù)測(cè)模型所用到的Python工具包。其中：

圖1 Python工具包示意圖

Pandas含有豐富的數(shù)據(jù)模型，可用于調(diào)用采集的瓦斯?jié)舛缺O(jiān)測(cè)數(shù)據(jù)，包括瓦斯?jié)舛戎导皩?duì)應(yīng)的監(jiān)測(cè)時(shí)間。

Numpy庫(kù)納入了大量?jī)?yōu)化的內(nèi)置數(shù)學(xué)函數(shù)，為非穩(wěn)定時(shí)間序列的平滑處理及均方根誤差的計(jì)算提供了模塊函數(shù)。

Matplotlib主要用于繪制瓦斯?jié)舛葦?shù)據(jù)、自相關(guān)及偏自相關(guān)圖等。其繪圖腳本簡(jiǎn)單易懂，支持滿足論文要求的二維平面圖、三維立體圖的繪制。

Statsmodels為瓦斯?jié)舛阮A(yù)測(cè)模型建立提供了ARIMA、ARMA、plot_acf等功能強(qiáng)大的模塊函數(shù)。擁有時(shí)間序列分析、非參數(shù)模型和估計(jì)以及空間面板模型等各種常見(jiàn)的模型。

2 建模原理

時(shí)間序列預(yù)測(cè)分析主要是對(duì)具有時(shí)序特性的一組數(shù)據(jù)進(jìn)行定量分析，挖掘時(shí)間序列數(shù)據(jù)潛在的趨勢(shì)、周期性、循環(huán)變化規(guī)律，并利用這種規(guī)律推斷出時(shí)間序列發(fā)展變化趨勢(shì)。ARIMA作為時(shí)間序列預(yù)測(cè)分析方法之一，其對(duì)時(shí)間序列預(yù)測(cè)主要包括時(shí)間序列平穩(wěn)性檢驗(yàn)與處理、模型定階、模型檢驗(yàn)及模型預(yù)測(cè)四個(gè)步驟。具體如下：

(1) 時(shí)間序列平穩(wěn)性檢驗(yàn)與處理

由于ARIMA模型不能捕獲非穩(wěn)定時(shí)間序列的變化規(guī)律，因此，首先需對(duì)時(shí)間序列進(jìn)行平穩(wěn)性檢驗(yàn)。時(shí)間序列平穩(wěn)性檢驗(yàn)的途徑一般有兩種：看圖法，通過(guò)繪制時(shí)序數(shù)據(jù)的均值與方差曲線，觀察其變化趨勢(shì)評(píng)斷時(shí)間序列是否平穩(wěn)；單位根檢驗(yàn)法，利用解算出的ADF值判斷時(shí)序數(shù)據(jù)是否在給定的置信區(qū)間上顯著，若顯著則證明該時(shí)序數(shù)據(jù)穩(wěn)定。如式(1)為隨機(jī)序列均值的表達(dá)式：

μt=E(yt),t=0,±1,±2,…n

(1)

其中，E(yt)表示隨機(jī)序列yt在t時(shí)刻的期望值。

若判斷出時(shí)間序列為非平穩(wěn)序列，就須利用平滑、差分等方法對(duì)時(shí)間序列進(jìn)行平穩(wěn)化處理。

① 平滑處理

平滑處理能夠消除或降低時(shí)序數(shù)據(jù)中較為明顯的不規(guī)則變動(dòng)、隨機(jī)誤差等成分的影響。公式(2)為簡(jiǎn)單移動(dòng)平均方程。

(2)

式中，St為時(shí)間點(diǎn)t的平滑值；2q+1是每次用來(lái)平均觀測(cè)值的個(gè)數(shù)。

② 差分處理

差分處理是通過(guò)對(duì)周期間隔的時(shí)序數(shù)據(jù)進(jìn)行求減，從而達(dá)到消除周期性因素對(duì)時(shí)間序列平穩(wěn)性影響的目的。公式(3)為d次差分方程：

wt=(1-Β)dyt

(3)

式中，wt表示差分后的平穩(wěn)數(shù)列；B為延遲算子，表示時(shí)間序列的時(shí)間指標(biāo)。

(2) 模型定階

自相關(guān)檢驗(yàn)(ACF)與偏自相關(guān)檢驗(yàn)(PACF)是獲取ARIMA階層p與階數(shù)q的常用方法。通過(guò)觀察繪制的ACF與PACF圖像，將ACF與PACF曲線第一次穿過(guò)上置信區(qū)間的橫坐標(biāo)值確定為p、q的值。式(4)為自相關(guān)函數(shù)的數(shù)學(xué)表達(dá)式：

(4)

如下式為偏自相關(guān)函數(shù)的數(shù)學(xué)表達(dá)式：

(5)

式中，φkk=φk-1,φkkφk-j；ρk表示隨機(jī)序列的自相關(guān)函數(shù)。

(3) 模型檢驗(yàn)

在利用建立的ARIMA模型正式對(duì)瓦斯?jié)舛茸兓?guī)律預(yù)測(cè)前，需對(duì)模型開(kāi)展合理性評(píng)價(jià)。模型評(píng)價(jià)是為了避免利用觀看圖像確定的p、q值存在偏差，影響后續(xù)瓦斯?jié)舛阮A(yù)測(cè)的精準(zhǔn)度。一般的，主要利用均方根誤差(RMSE)與正態(tài)性檢驗(yàn)對(duì)模型擬合優(yōu)度進(jìn)行評(píng)估。均方根誤差公式為：

(6)

式中，N為觀測(cè)次數(shù)。

正態(tài)性檢驗(yàn)是通過(guò)獲取的監(jiān)測(cè)數(shù)據(jù)判斷總體是否服從正態(tài)分布的檢驗(yàn)。通過(guò)繪制ARIMA模型標(biāo)準(zhǔn)殘差Q-Q圖來(lái)評(píng)定建立的預(yù)測(cè)模型是否合理。

(4) 模型預(yù)測(cè)

時(shí)間序列分析的主要作用是利用數(shù)學(xué)模型描述一組時(shí)間序列的自相關(guān)性，進(jìn)而預(yù)測(cè)時(shí)間序列未來(lái)發(fā)展的變化規(guī)律。通過(guò)上述內(nèi)容確定的最優(yōu)預(yù)測(cè)模型對(duì)時(shí)間序列進(jìn)行預(yù)測(cè)，并觀察預(yù)測(cè)數(shù)據(jù)折線是否在設(shè)定的置信區(qū)間內(nèi)，來(lái)評(píng)估預(yù)測(cè)效果的優(yōu)劣。

3 實(shí)例分析

貴州某礦屬于高瓦斯礦井。本研究利用瓦斯傳感器，采集該礦井1302綜采工作面瓦斯?jié)舛?，并將采集的?shù)據(jù)存儲(chǔ)于地面主機(jī)硬盤(pán)。選用2018年3月5日 0點(diǎn)至2018年3月7日 3點(diǎn)時(shí)間段內(nèi)采集的瓦斯數(shù)據(jù)作為樣本數(shù)據(jù)，通過(guò)建立的瓦斯?jié)舛茸兓?guī)律分析模型，預(yù)測(cè)未來(lái)5小時(shí)內(nèi)的瓦斯?jié)舛戎怠?/p>

在實(shí)驗(yàn)中，通過(guò)采集的礦井瓦斯?jié)舛葦?shù)據(jù)構(gòu)成的隨機(jī)變量序列并非平穩(wěn)的。由于ARIMA模型無(wú)法捕捉不平穩(wěn)時(shí)序數(shù)據(jù)的變化規(guī)律，因此，在研究瓦斯?jié)舛入S機(jī)變量序列時(shí)，首先需對(duì)其進(jìn)行平穩(wěn)性檢驗(yàn)與處理。

通過(guò)調(diào)用的瓦斯?jié)舛葰v史監(jiān)測(cè)數(shù)據(jù)組成時(shí)間序列，并利用均值、方差判定時(shí)間序列的穩(wěn)定性。圖2為利用該礦2018年3月5日 0點(diǎn)至2018年3月7日 3點(diǎn)內(nèi)，監(jiān)測(cè)的瓦斯?jié)舛葦?shù)據(jù)繪制的曲線圖。如圖2所示，該時(shí)間段內(nèi)，瓦斯?jié)舛葦?shù)據(jù)均值整體趨勢(shì)緩慢下降，且不同時(shí)間點(diǎn)的方差值也不同，因此，該瓦斯?jié)舛葧r(shí)間序列為非穩(wěn)定序列。本文采用平滑、差分法剔除周期性、不規(guī)則變動(dòng)、隨機(jī)誤差等因素，使時(shí)間序列平穩(wěn)。

① 瓦斯?jié)舛葧r(shí)間序列平滑處理

如圖3(a)為對(duì)瓦斯?jié)舛葧r(shí)間序列進(jìn)行平滑處理后的效果圖。由圖3(a)可知，平滑處理后的時(shí)間序列變得較為平穩(wěn)，但仍屬于不平穩(wěn)序列。

② 瓦斯?jié)舛葧r(shí)間序列差分處理

由圖3(b)可知，差分處理后的時(shí)間序列均值與方差曲線波動(dòng)較小趨于常數(shù)，可直觀的判定其為穩(wěn)定序列。

③ 瓦斯?jié)舛葧r(shí)間序列平穩(wěn)性檢驗(yàn)

由表1得，1%的Critical value值遠(yuǎn)大于Test Statistic值，由此表明，在99%的置信度下處理后的序列穩(wěn)定，滿足ARIMA預(yù)測(cè)的基本要求。

圖2 瓦斯?jié)舛仍急O(jiān)測(cè)數(shù)據(jù)圖

圖3 瓦斯?jié)舛绕椒€(wěn)處理效果圖

表1 ADF test

Test Statisticp-value#Lags UsedNumber of Observations UsedCritical value (1%)Critical value (5%)Critical value (10%)-5.702883e+007.604219e-076.000000e+001.510000e+02-3.474416e+00-2.880878e+00-2.577081e+00

當(dāng)瓦斯?jié)舛葧r(shí)間序列變?yōu)槠椒€(wěn)序列后，就需利用自相關(guān)圖與偏自相關(guān)圖確定參數(shù)p、q。圖4、圖5分別為瓦斯?jié)舛葧r(shí)間序列的自相關(guān)圖與偏自相關(guān)圖，其中，置信區(qū)間為上下灰色虛線之間的區(qū)域。通過(guò)觀察ACF、PACF第一次穿過(guò)上置信區(qū)間的橫坐標(biāo)值，p取1，q取1。

圖4 一階差分瓦斯?jié)舛刃蛄蠥CF

圖5 一階差分瓦斯?jié)舛刃蛄蠵ACF

通過(guò)調(diào)用Python的ARIMA模塊函數(shù)，并結(jié)合獲取的p、d、q值，建立瓦斯?jié)舛阮A(yù)測(cè)模型ARIMA(1,1,1)。建模代碼如下：

from statsmodels.tsa.arima_model impo-rt ARIMA

model = ARIMA(ts_diff_1, order=(1,1,1))

diff_recover_1 = predict_ts.add(diff_shift_ts)

rol_recover = diff_recover*12 - rol_sum.shift(1)

log_recover = np.exp(rol_recover)

在利用ARIMA(1,1,1)模型對(duì)瓦斯?jié)舛阮A(yù)測(cè)前，需對(duì)其進(jìn)行評(píng)價(jià)。模型評(píng)價(jià)是對(duì)模型的擬合優(yōu)度進(jìn)行評(píng)價(jià)，假如擬合效果差，則必須重新對(duì)瓦斯?jié)舛葧r(shí)間序列進(jìn)行差分處理，并確定模型階數(shù)。本文通過(guò)對(duì)2018年3月5日 0點(diǎn)至2018年3月7日 3點(diǎn)時(shí)間段內(nèi)瓦斯?jié)舛葰v史監(jiān)測(cè)數(shù)據(jù)與預(yù)測(cè)數(shù)進(jìn)行擬合，擬合效果如圖6。

通過(guò)圖6(a)(紅、綠、藍(lán)色曲線分別為瓦斯?jié)舛葰v史數(shù)據(jù)、ARIMA模型預(yù)測(cè)數(shù)據(jù)及ARMA模型預(yù)測(cè)數(shù)據(jù))對(duì)瓦斯?jié)舛阮A(yù)測(cè)數(shù)據(jù)與監(jiān)測(cè)數(shù)據(jù)對(duì)比，可直觀的判斷出瓦斯?jié)舛阮A(yù)測(cè)結(jié)果理想，且ARIMA(1,1,1)模型預(yù)測(cè)效果好于ARMA(1,1)模型。由圖6(b)(紅、綠、黃、藍(lán)曲線分別代表ARIMA(1,1,1)、ARIMA(0,1,1)、ARIMA(1,1,0)瓦斯?jié)舛阮A(yù)測(cè)數(shù)據(jù)及歷史數(shù)據(jù))可知，利用ARIMA(1,1,1)、ARIMA(0,1,1)、ARIMA(1,1,0)模型預(yù)測(cè)出的瓦斯?jié)舛扰c原始瓦斯?jié)舛鹊臄M合優(yōu)度較好。其中，ARIMA(1,1,1)模型預(yù)測(cè)精度高于ARIMA(0,1,1)、ARIMA(1,1,0)模型。表2為瓦斯?jié)舛阮A(yù)測(cè)效果表。由表2得，ARMA(1,1)、ARIMA(1,1,1)、ARIMA(0,1,1)及ARIMA(1,1,0)模型預(yù)測(cè)的瓦斯?jié)舛阮A(yù)測(cè)誤差均在允許范圍內(nèi)，預(yù)測(cè)精度都較為理想。其中，ARIMA(1,1,1)的均方根誤差最趨于0，預(yù)測(cè)精度最高。由圖7可知，中間部分點(diǎn)幾乎與直線重合，由此證明瓦斯?jié)舛阮A(yù)測(cè)模型ARIMA(1,1,1)的殘差相服從正態(tài)分布，即所建立的模型是合理的。

表2 瓦斯?jié)舛阮A(yù)測(cè)效果

通過(guò)上述分析驗(yàn)證，最終確定選用ARIMA(1,1,1)模型對(duì)瓦斯?jié)舛缺O(jiān)測(cè)數(shù)據(jù)的發(fā)展變化規(guī)律進(jìn)行預(yù)測(cè)。通過(guò)觀察圖8可發(fā)現(xiàn)，未來(lái)5小時(shí)內(nèi)(2018年3月7日 3點(diǎn)至2018年3月7日 8點(diǎn))瓦斯?jié)舛阮A(yù)測(cè)曲線均在95%的置信區(qū)間里，因此，瓦斯?jié)舛阮A(yù)測(cè)精度較高。

圖6 瓦斯?jié)舛葦M合效果圖

圖7 ARIMA模型的標(biāo)準(zhǔn)殘差Q-Q圖

圖8 瓦斯?jié)舛任磥?lái)變化規(guī)律預(yù)測(cè)效果圖

4 結(jié)論

(1) 利用建立的預(yù)測(cè)模型對(duì)某礦2018年3月5日 0點(diǎn)至2018年3月7日 3點(diǎn)時(shí)間段內(nèi)的瓦斯?jié)舛冗M(jìn)行預(yù)測(cè)。通過(guò)對(duì)該段時(shí)間內(nèi)的瓦斯?jié)舛阮A(yù)測(cè)數(shù)據(jù)與歷史數(shù)據(jù)對(duì)比分析，發(fā)現(xiàn)二者之間的均方根誤差為2.34%，證明該預(yù)測(cè)模型預(yù)測(cè)精度較高。

(2) 針對(duì)預(yù)測(cè)模型算法復(fù)雜、應(yīng)用不靈活的問(wèn)題，提出基于Python的瓦斯?jié)舛葧r(shí)間序列預(yù)測(cè)方法。實(shí)例分析表明，該預(yù)測(cè)算法易于理解、開(kāi)發(fā)要求低，實(shí)現(xiàn)了瓦斯?jié)舛阮A(yù)測(cè)的可視化，預(yù)測(cè)效果較好，為降低瓦斯災(zāi)害提供了一定的技術(shù)支撐。

(3) 通過(guò)對(duì)比分析ARIMA與ARMA模型瓦斯?jié)舛阮A(yù)測(cè)效果可知，ARIMA模型提高了瓦斯?jié)舛阮A(yù)測(cè)精準(zhǔn)度，可為煤礦提前獲知瓦斯?jié)舛茸兓?guī)律提供科學(xué)依據(jù)。