黃 悅,幸福堂,2,石癸鑫,吳孟龍
(1.武漢科技大學(xué) 資源與環(huán)境工程學(xué)院,湖北 武漢 430081;2.武漢科技大學(xué) 冶金礦產(chǎn)資源高效利用與造塊重點(diǎn)實(shí)驗(yàn)室,湖北 武漢 430081)
近年來,隨著采礦業(yè)的發(fā)展,礦業(yè)安全生產(chǎn)事故數(shù)量和傷亡人數(shù)居高不下,給社會(huì)發(fā)展帶來了負(fù)面影響。據(jù)近年來我國礦山事故統(tǒng)計(jì)顯示,2005年至2019年,全國共發(fā)生礦業(yè)安全生產(chǎn)事故3 826起[1],死亡人數(shù)達(dá)到14 251人[1],事故起數(shù)和死亡人數(shù)超過發(fā)達(dá)國家,安全問題嚴(yán)重制約了我國礦山企業(yè)的發(fā)展。因此,準(zhǔn)確預(yù)測礦業(yè)安全生產(chǎn)事故對制訂合理的礦業(yè)安全決策,以及提高礦業(yè)安全管理水平具有重要的意義。
礦業(yè)安全系統(tǒng)是一個(gè)復(fù)雜的循環(huán)系統(tǒng)[2],礦業(yè)安全生產(chǎn)受自然條件、作業(yè)環(huán)境等諸多因素影響[3]。目前,國內(nèi)外眾多學(xué)者對礦業(yè)安全生產(chǎn)事故預(yù)測的方法展開了研究[4]。劉素兵等[5]將GM(1,1)與支持向量機(jī)融合,建立了灰色支持向量機(jī)組合預(yù)測模型;李懷俊等[6]通過參數(shù)尋優(yōu)改進(jìn)了ARMA模型,提高了ARMA模型的預(yù)測精度;Barman等[7]通過螢火蟲算法(Firefly Algorithm,F(xiàn)A)優(yōu)化SVM的參數(shù),提高了SVM模型的預(yù)測精度;劉俊娥等[8]提出了一種模糊信息?;cSVM相結(jié)合的方法,對瓦斯涌出量進(jìn)行了預(yù)測;JANG Y J等[9]使用長短期記憶(LSTM)遞歸神經(jīng)網(wǎng)絡(luò)(RNN)的深度學(xué)習(xí)模型來預(yù)測業(yè)務(wù)失??;WU Menglong等[10]對經(jīng)驗(yàn)?zāi)B(tài)分解后的序列單獨(dú)分析,建立預(yù)測模型來預(yù)測礦山安全生產(chǎn)形勢。然而上述方法預(yù)測效果的魯棒性較差,并且只能得到礦業(yè)安全生產(chǎn)事故的點(diǎn)預(yù)測結(jié)果,目前關(guān)于礦業(yè)安全區(qū)間預(yù)測的研究還較少。區(qū)間預(yù)測是一類既含有確定性又含有不確定性的預(yù)測[11]。為此,考慮將模糊信息?;氲降V業(yè)安全生產(chǎn)事故預(yù)測中,得到礦業(yè)安全生產(chǎn)事故時(shí)間序列的模糊區(qū)間預(yù)測結(jié)果。區(qū)間預(yù)測能夠提供對應(yīng)某一期望概率的區(qū)間預(yù)測結(jié)果,獲取比確定性預(yù)測結(jié)果更多的信息,有利于決策者評估礦業(yè)安全生產(chǎn)事故的不確定性和風(fēng)險(xiǎn)因素,做出更為合理的決策[12]。因此,準(zhǔn)確預(yù)測礦業(yè)安全生產(chǎn)態(tài)勢的變化趨勢和變化空間,在實(shí)際應(yīng)用中具有重要的參考價(jià)值和指導(dǎo)意義。
筆者先將礦業(yè)安全生產(chǎn)態(tài)勢時(shí)間序列進(jìn)行小波變換,然后采用樣本熵將子序列重組,將重組的分量按照一定窗口進(jìn)行模糊信息?;?,得到能夠反映原始數(shù)據(jù)特點(diǎn)的信息粒;最后基于FCM對窗口化的 3個(gè)參數(shù)序列進(jìn)行預(yù)測,以實(shí)現(xiàn)礦業(yè)安全生產(chǎn)事故時(shí)間序列的模糊區(qū)間預(yù)測。
通過小波變換和樣本熵(SE)算法將礦業(yè)安全生產(chǎn)事故的時(shí)間序列分解并重組,然后采用模糊信息?;姆椒ǎ瑢⒌V業(yè)安全生產(chǎn)事故時(shí)間序列模糊?;玫絃ow、R、Up共3個(gè)?;瘏?shù)時(shí)間序列;最后依據(jù)模糊C均值聚類算法對窗口化的時(shí)間序列進(jìn)行預(yù)測,建立礦業(yè)安全生產(chǎn)事故時(shí)間序列預(yù)測模型。
原始時(shí)間序列x(t)的小波變換可以表示為[13]:
(1)
式中:s為伸縮因子;t為時(shí)間;τ為平移因子;φ*為小波基函數(shù)φs,τ的復(fù)共軛。
小波基基于以下2個(gè)原則確定:
1)能—熵比準(zhǔn)則
將能量與熵結(jié)合起來形成一個(gè)新準(zhǔn)則,即能—熵比準(zhǔn)則[14],見下式:
R(j)=Eenergy(j)/Eentropy(j)
(2)
式中:Eenergy(j)為原始時(shí)間序列在第j分解水平時(shí)小波系數(shù)序列包含的能量;Eentropy(j)為原始時(shí)間序列在第j分解水平時(shí)小波系數(shù)序列的Shannon熵。
由式(2)可以看出,能量越大,熵越小,則能—熵比越大,小波基[15]越優(yōu)。
2)相關(guān)性準(zhǔn)則
在處理和分析時(shí)間序列時(shí),相關(guān)性分析[16]是一種在時(shí)域中對序列特性進(jìn)行描述的重要方法。
礦業(yè)安全生產(chǎn)事故時(shí)間序列的細(xì)節(jié)和近似部分可通過低通和高通濾波器獲得,通過對第(j-1)分解水平的近似系數(shù)與低通濾波器系數(shù)進(jìn)行卷積,可以獲得第j分解水平的近似系數(shù):
(3)
式中:aj,k為近似系數(shù);dj,k為細(xì)節(jié)系數(shù);h、g分別為所選取的小波函數(shù)對應(yīng)的低通和高通濾波器的抽頭系數(shù)序列;m為維度;k=1,2,…,m。
經(jīng)過小波變換后的時(shí)間序列具有一定重復(fù)性和復(fù)雜性,故采用樣本熵算法對各分量時(shí)間序列進(jìn)行重組。
樣本熵的計(jì)算過程如下[17]:
1)將小波變換后的序列按序號形成n維的向量組Xm(1),…,Xm(N-m+1),其中Xm(i)={x(i),x(i+1),…,x(i-m+1)},1≤i≤N-m+1。
(4)
3)定義B(m)(r)為:
(5)
4)將維數(shù)m增加至m+1維,對m+1維向量重復(fù)步驟1)~3)得到B(m+1)(r)。
5)該序列的樣本熵值為:
(6)
式中:m為重構(gòu)維數(shù),一般取m=2;r為相似容限,一般為原始數(shù)據(jù)標(biāo)準(zhǔn)差的0.10~0.25;n為數(shù)據(jù)長度。
根據(jù)一定的窗口w將原始時(shí)間序列分割成一些子序列。原始時(shí)間序列窗口以1個(gè)季度為1個(gè)單元,然后進(jìn)行模糊信息粒化,得到窗口化后的3個(gè)?;瘏?shù)時(shí)間序列[18]。模糊信息粒化后模糊粒子的隸屬函數(shù)為[19]:
(7)
式中:x為論域中的變量;a、m、b為參數(shù),對應(yīng)窗口模糊化后得到的3個(gè)變量Low、R、Up;Low為窗口變化的最小值,R為窗口變化的平均值,Up為窗口變化的最大值。
基于FCM對窗口化的時(shí)間序列進(jìn)行預(yù)測,F(xiàn)CM的建模過程如下[20]:
1)確定聚類中心數(shù)。設(shè)原始礦業(yè)安全生產(chǎn)事故時(shí)間序列為Xt,t=1,2,…,N;待預(yù)測值為XN+1。首先對原始序列進(jìn)行聚類,通過下式確定原始序列的聚類中心數(shù):
(8)
式中:c為聚類中心數(shù);Xmax、Xmin為原始序列中的最大值和最小值;[]為取整運(yùn)算。
2)初始化聚類中心。根據(jù)式(8)計(jì)算得到聚類中心數(shù),隨機(jī)選定初始聚類中心Ck。Ck是一個(gè)長度為c的向量,k=1,2…,c。
3)計(jì)算目標(biāo)化函數(shù)。FCM的目標(biāo)化函數(shù)為:
(9)
式中:Ck為模糊組i的聚類中心;uij為屬于第i個(gè)聚類的隸屬度,介于0~1;dij為第i個(gè)聚類中心與第j個(gè)數(shù)據(jù)點(diǎn)間的歐幾里得距離。
以2005—2019年礦業(yè)單位從業(yè)人員事故的萬人死亡率為研究對象,統(tǒng)計(jì)得到全國礦業(yè)事故萬人死亡率時(shí)間序列,如圖1所示。
圖1 2005—2019年全國礦業(yè)事故萬人死亡率時(shí)間序列
采用能—熵比準(zhǔn)則及相關(guān)性準(zhǔn)則選取最優(yōu)小波基,并依據(jù)相關(guān)性準(zhǔn)則確定臨界分解水平。
筆者采用較為常用的dbN(1~10)、coifN(1~5)、symN(2~10)3種小波系選取最優(yōu)小波基。
2.1.1 能—熵比準(zhǔn)則的計(jì)算結(jié)果
基于能—熵比準(zhǔn)則的計(jì)算結(jié)果如表1所示。
表1 基于能—熵比準(zhǔn)則R(j)的計(jì)算結(jié)果
將表1中數(shù)據(jù)依據(jù)時(shí)間序列比值高低排序,結(jié)果如表2所示。
表2 基于能—熵比準(zhǔn)則R(j)對小波基的排序
2.1.2 相關(guān)性準(zhǔn)則的計(jì)算結(jié)果
基于相關(guān)性準(zhǔn)則的計(jì)算結(jié)果如表3所示。
表3 基于相關(guān)性準(zhǔn)則的計(jì)算結(jié)果
將表3中數(shù)據(jù)依據(jù)相關(guān)性系數(shù)大小進(jìn)行排序,排序結(jié)果如表4所示。
表4 基于相關(guān)性準(zhǔn)則對小波基的排序結(jié)果
2.1.3 小波基優(yōu)選結(jié)果
計(jì)算表2、表4中小波基的得分,根據(jù)其對應(yīng)的分值進(jìn)行排序,排序結(jié)果如表5所示。
表5 小波基排序結(jié)果
由表5可知,小波基coif5得分最高,coif5為最優(yōu)小波基。
2.1.4 臨界分解水平的確定
基于臨界分解水平判別準(zhǔn)則,對曲線的臨界分解水平進(jìn)行分析判斷,得到各分解水平趨勢序列與原始礦業(yè)安全生產(chǎn)事故時(shí)間序列的相關(guān)性系數(shù),結(jié)果如圖2所示。
圖2 各分解水平的相關(guān)性系數(shù)
由圖2可知,礦業(yè)安全生產(chǎn)事故時(shí)間序列曲線分解水平大于5時(shí),相關(guān)性系數(shù)發(fā)生顯著變化,故確定臨界分解水平為5。
2.1.5 時(shí)間序列的分解
在確定最優(yōu)小波基和臨界分解水平后,利用MATLAB小波工具箱對其進(jìn)行小波變換。礦業(yè)安全生產(chǎn)事故時(shí)間序列5分解水平小波變換后的結(jié)果如圖3所示(圖中S代表原始序列,ai代表近似部分,di代表細(xì)節(jié)部分)。
(a)細(xì)節(jié)部分 (b)近似部分
由圖3可知,隨著分解水平的增加,細(xì)節(jié)部分逐漸被提取出來。經(jīng)過5分解水平小波變換后,可以滿足礦業(yè)安全生產(chǎn)事故時(shí)間序列分解的要求。
基于礦業(yè)安全生產(chǎn)事故時(shí)間序列進(jìn)行小波變換后的子序列,計(jì)算給定序列與子序列的樣本熵值,各分量的樣本熵如圖4所示。
圖4 各分量的樣本熵
根據(jù)子序列樣本熵的分布情況,取λ=0.2。d1為隨機(jī)分量;d2~d4重組為細(xì)節(jié)分量s;d5、a5重組為趨勢分量r。重組后各分量結(jié)果如圖5所示。
圖5 重組后各分量
由圖5可知,樣本熵重組后的序列具有隨機(jī)性、波動(dòng)性等特點(diǎn)。為此,先將原始數(shù)據(jù)按照一定窗口進(jìn)行模糊信息?;?,將礦業(yè)安全生產(chǎn)態(tài)勢時(shí)間序列的萬人死亡率模糊信息粒化為Low、R、Up共3個(gè)參數(shù),3個(gè)子序列的模糊信息?;Y(jié)果如圖6所示。該例中Low具體為窗口化子序列每月變化的最小值,R具體為礦業(yè)安全生產(chǎn)態(tài)勢萬人死亡率每月變化的平均值,Up具體為礦業(yè)安全生產(chǎn)態(tài)勢萬人死亡率每月變化的最大值。
(a)d2序列粒化
(b)d3序列?;?/p>
(c)d4序列?;?/p>
對于窗口化的時(shí)間序列,按照窗口模糊信息?;Y(jié)果將其分成6個(gè)數(shù)據(jù)集,根據(jù)每年數(shù)據(jù)情況確定討論區(qū)間的上下界計(jì)算模型的論域。基于樣本數(shù)據(jù)求得的隸屬度值見表6,樣本模糊化后的結(jié)果如表7 所示。
表6 樣本隸屬度
表7 樣本模糊關(guān)系
由表6和表7并結(jié)合式(9)計(jì)算,求得測試集預(yù)測結(jié)果,如表8所示。
表8 測試集預(yù)測結(jié)果
合并窗口化后,礦業(yè)安全生產(chǎn)事故時(shí)間序列的模糊區(qū)間預(yù)測結(jié)果如圖7所示。
圖7 窗口化區(qū)間預(yù)測結(jié)果
利用訓(xùn)練集150組樣本建立的模型,對2018年7月至2019年10月測試集的21組數(shù)據(jù)進(jìn)行預(yù)測,并計(jì)算區(qū)間上下邊界值與中值的預(yù)測值,結(jié)果如表9 所示。
表9 測試集預(yù)測誤差對比
由表9可知:R值的平均相對誤差為17.737 1%,Up值的平均相對誤差為8.771 2%。Up值的預(yù)測精度明顯高于R值的預(yù)測精度,表明模糊信息?;蟮纳舷逻吔缰的軌蜉^為準(zhǔn)確地描述礦業(yè)安全生產(chǎn)事故的變化范圍。
將礦業(yè)安全生產(chǎn)事故時(shí)間序列的模糊區(qū)間預(yù)測模型與小波-SVM(支持向量機(jī))、小波-ARMA(自回歸滑動(dòng)平均模型)和小波神經(jīng)網(wǎng)絡(luò)的預(yù)測效果進(jìn)行對比,各模型的預(yù)測誤差分析如表10所示。其中,MAE為預(yù)測值和實(shí)際值之間的平均絕對誤差;MRE為絕對誤差與真實(shí)值的比率,即平均相對誤差;RMSE為均方根誤差,其對異常值更為敏感。
表10 各模型預(yù)測誤差對比分析結(jié)果
由表10可知,本文預(yù)測模型的MAE、MRE、RMSE分別為0.445 7%、13.254 1%、0.372 3%,各項(xiàng)指標(biāo)均低于對比模型,表明本文模型預(yù)測效果較為穩(wěn)定,精度較高。
1)提出一種礦業(yè)安全生產(chǎn)事故時(shí)間序列的模糊區(qū)間預(yù)測模型,能夠?qū)⒌V業(yè)系統(tǒng)中的隨機(jī)因素等考慮進(jìn)去,相較于定量值的點(diǎn)預(yù)測,區(qū)間預(yù)測的結(jié)果更為可靠,更具有實(shí)用價(jià)值。
2)R、Up值時(shí)間序列預(yù)測的平均相對誤差分別為17.737 1%、8.771 2%。邊界值的預(yù)測精度較高,表明模型的區(qū)間范圍較為可靠,也進(jìn)一步地說明了模糊信息粒化后的邊界值結(jié)果的合理性。
3)相對于小波-支持向量機(jī)、小波-自回歸滑動(dòng)平均模型和小波神經(jīng)網(wǎng)絡(luò)的預(yù)測結(jié)果,本文的預(yù)測模型平均絕對誤差(MAE)、平均相對誤差(MRE)、均方根誤差(RMSE)3個(gè)指標(biāo)均為最小,表明訓(xùn)練得到的模型具有較高的可靠性,對礦業(yè)安全決策具有一定的指導(dǎo)意義。