黃 悅，幸福堂，2，石癸鑫，吳孟龍

(1.武漢科技大學 資源與環(huán)境工程學院，湖北 武漢 430081；2.武漢科技大學 冶金礦產(chǎn)資源高效利用與造塊重點實驗室，湖北 武漢 430081)

近年來，隨著采礦業(yè)的發(fā)展，礦業(yè)安全生產(chǎn)事故數(shù)量和傷亡人數(shù)居高不下，給社會發(fā)展帶來了負面影響。據(jù)近年來我國礦山事故統(tǒng)計顯示，2005年至2019年，全國共發(fā)生礦業(yè)安全生產(chǎn)事故3 826起[1]，死亡人數(shù)達到14 251人[1]，事故起數(shù)和死亡人數(shù)超過發(fā)達國家，安全問題嚴重制約了我國礦山企業(yè)的發(fā)展。因此，準確預(yù)測礦業(yè)安全生產(chǎn)事故對制訂合理的礦業(yè)安全決策，以及提高礦業(yè)安全管理水平具有重要的意義。

礦業(yè)安全系統(tǒng)是一個復(fù)雜的循環(huán)系統(tǒng)[2]，礦業(yè)安全生產(chǎn)受自然條件、作業(yè)環(huán)境等諸多因素影響[3]。目前，國內(nèi)外眾多學者對礦業(yè)安全生產(chǎn)事故預(yù)測的方法展開了研究[4]。劉素兵等[5]將GM(1,1)與支持向量機融合，建立了灰色支持向量機組合預(yù)測模型；李懷俊等[6]通過參數(shù)尋優(yōu)改進了ARMA模型，提高了ARMA模型的預(yù)測精度；Barman等[7]通過螢火蟲算法(Firefly Algorithm，F(xiàn)A)優(yōu)化SVM的參數(shù)，提高了SVM模型的預(yù)測精度；劉俊娥等[8]提出了一種模糊信息?；cSVM相結(jié)合的方法，對瓦斯涌出量進行了預(yù)測；JANG Y J等[9]使用長短期記憶(LSTM)遞歸神經(jīng)網(wǎng)絡(luò)(RNN)的深度學習模型來預(yù)測業(yè)務(wù)失??；WU Menglong等[10]對經(jīng)驗?zāi)B(tài)分解后的序列單獨分析,建立預(yù)測模型來預(yù)測礦山安全生產(chǎn)形勢。然而上述方法預(yù)測效果的魯棒性較差，并且只能得到礦業(yè)安全生產(chǎn)事故的點預(yù)測結(jié)果，目前關(guān)于礦業(yè)安全區(qū)間預(yù)測的研究還較少。區(qū)間預(yù)測是一類既含有確定性又含有不確定性的預(yù)測[11]。為此，考慮將模糊信息?；氲降V業(yè)安全生產(chǎn)事故預(yù)測中，得到礦業(yè)安全生產(chǎn)事故時間序列的模糊區(qū)間預(yù)測結(jié)果。區(qū)間預(yù)測能夠提供對應(yīng)某一期望概率的區(qū)間預(yù)測結(jié)果，獲取比確定性預(yù)測結(jié)果更多的信息，有利于決策者評估礦業(yè)安全生產(chǎn)事故的不確定性和風險因素，做出更為合理的決策[12]。因此，準確預(yù)測礦業(yè)安全生產(chǎn)態(tài)勢的變化趨勢和變化空間，在實際應(yīng)用中具有重要的參考價值和指導意義。

筆者先將礦業(yè)安全生產(chǎn)態(tài)勢時間序列進行小波變換，然后采用樣本熵將子序列重組，將重組的分量按照一定窗口進行模糊信息?；玫侥軌蚍从吃紨?shù)據(jù)特點的信息粒；最后基于FCM對窗口化的 3個參數(shù)序列進行預(yù)測，以實現(xiàn)礦業(yè)安全生產(chǎn)事故時間序列的模糊區(qū)間預(yù)測。

1 模型的建立

通過小波變換和樣本熵(SE)算法將礦業(yè)安全生產(chǎn)事故的時間序列分解并重組，然后采用模糊信息?；姆椒?，將礦業(yè)安全生產(chǎn)事故時間序列模糊?；?，得到Low、R、Up共3個?；瘏?shù)時間序列；最后依據(jù)模糊C均值聚類算法對窗口化的時間序列進行預(yù)測，建立礦業(yè)安全生產(chǎn)事故時間序列預(yù)測模型。

1.1 時間序列的分解與重組

原始時間序列x(t)的小波變換可以表示為[13]：

(1)

式中：s為伸縮因子；t為時間；τ為平移因子；φ*為小波基函數(shù)φs,τ的復(fù)共軛。

小波基基于以下2個原則確定：

1)能—熵比準則

將能量與熵結(jié)合起來形成一個新準則，即能—熵比準則[14]，見下式：

R(j)=Eenergy(j)/Eentropy(j)

(2)

式中：Eenergy(j)為原始時間序列在第j分解水平時小波系數(shù)序列包含的能量；Eentropy(j)為原始時間序列在第j分解水平時小波系數(shù)序列的Shannon熵。

由式(2)可以看出，能量越大，熵越小，則能—熵比越大，小波基[15]越優(yōu)。

2)相關(guān)性準則

在處理和分析時間序列時，相關(guān)性分析[16]是一種在時域中對序列特性進行描述的重要方法。

礦業(yè)安全生產(chǎn)事故時間序列的細節(jié)和近似部分可通過低通和高通濾波器獲得，通過對第(j-1)分解水平的近似系數(shù)與低通濾波器系數(shù)進行卷積，可以獲得第j分解水平的近似系數(shù)：

(3)

式中：aj,k為近似系數(shù)；dj,k為細節(jié)系數(shù)；h、g分別為所選取的小波函數(shù)對應(yīng)的低通和高通濾波器的抽頭系數(shù)序列；m為維度；k=1,2,…,m。

經(jīng)過小波變換后的時間序列具有一定重復(fù)性和復(fù)雜性，故采用樣本熵算法對各分量時間序列進行重組。

樣本熵的計算過程如下[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ù)標準差的0.10～0.25；n為數(shù)據(jù)長度。

1.2 重組序列的模糊信息粒化

根據(jù)一定的窗口w將原始時間序列分割成一些子序列。原始時間序列窗口以1個季度為1個單元，然后進行模糊信息?；玫酱翱诨蟮?個?；瘏?shù)時間序列[18]。模糊信息?；竽：Ｗ拥碾`屬函數(shù)為[19]：

(7)

式中：x為論域中的變量；a、m、b為參數(shù)，對應(yīng)窗口模糊化后得到的3個變量Low、R、Up；Low為窗口變化的最小值，R為窗口變化的平均值，Up為窗口變化的最大值。

1.3 基于FCM的窗口化序列模糊預(yù)測

基于FCM對窗口化的時間序列進行預(yù)測，F(xiàn)CM的建模過程如下[20]：

1)確定聚類中心數(shù)。設(shè)原始礦業(yè)安全生產(chǎn)事故時間序列為Xt，t=1，2，…，N；待預(yù)測值為XN+1。首先對原始序列進行聚類，通過下式確定原始序列的聚類中心數(shù)：

(8)

式中：c為聚類中心數(shù)；Xmax、Xmin為原始序列中的最大值和最小值；[]為取整運算。

2)初始化聚類中心。根據(jù)式(8)計算得到聚類中心數(shù)，隨機選定初始聚類中心Ck。Ck是一個長度為c的向量，k=1，2…，c。

3)計算目標化函數(shù)。FCM的目標化函數(shù)為：

(9)

式中：Ck為模糊組i的聚類中心；uij為屬于第i個聚類的隸屬度，介于0～1；dij為第i個聚類中心與第j個數(shù)據(jù)點間的歐幾里得距離。

2 模型應(yīng)用

以2005—2019年礦業(yè)單位從業(yè)人員事故的萬人死亡率為研究對象,統(tǒng)計得到全國礦業(yè)事故萬人死亡率時間序列,如圖1所示。

圖1 2005—2019年全國礦業(yè)事故萬人死亡率時間序列

2.1 小波變換

采用能—熵比準則及相關(guān)性準則選取最優(yōu)小波基，并依據(jù)相關(guān)性準則確定臨界分解水平。

筆者采用較為常用的dbN(1～10)、coifN(1～5)、symN(2～10)3種小波系選取最優(yōu)小波基。

2.1.1 能—熵比準則的計算結(jié)果

基于能—熵比準則的計算結(jié)果如表1所示。

表1 基于能—熵比準則R(j)的計算結(jié)果

將表1中數(shù)據(jù)依據(jù)時間序列比值高低排序，結(jié)果如表2所示。

表2 基于能—熵比準則R(j)對小波基的排序

2.1.2 相關(guān)性準則的計算結(jié)果

基于相關(guān)性準則的計算結(jié)果如表3所示。

表3 基于相關(guān)性準則的計算結(jié)果

將表3中數(shù)據(jù)依據(jù)相關(guān)性系數(shù)大小進行排序，排序結(jié)果如表4所示。

表4 基于相關(guān)性準則對小波基的排序結(jié)果

2.1.3 小波基優(yōu)選結(jié)果

計算表2、表4中小波基的得分，根據(jù)其對應(yīng)的分值進行排序，排序結(jié)果如表5所示。

表5 小波基排序結(jié)果

由表5可知，小波基coif5得分最高，coif5為最優(yōu)小波基。

2.1.4 臨界分解水平的確定

基于臨界分解水平判別準則，對曲線的臨界分解水平進行分析判斷,得到各分解水平趨勢序列與原始礦業(yè)安全生產(chǎn)事故時間序列的相關(guān)性系數(shù),結(jié)果如圖2所示。

圖2 各分解水平的相關(guān)性系數(shù)

由圖2可知，礦業(yè)安全生產(chǎn)事故時間序列曲線分解水平大于5時，相關(guān)性系數(shù)發(fā)生顯著變化，故確定臨界分解水平為5。

2.1.5 時間序列的分解

在確定最優(yōu)小波基和臨界分解水平后，利用MATLAB小波工具箱對其進行小波變換。礦業(yè)安全生產(chǎn)事故時間序列5分解水平小波變換后的結(jié)果如圖3所示(圖中S代表原始序列，ai代表近似部分，di代表細節(jié)部分)。

(a)細節(jié)部分 (b)近似部分

由圖3可知，隨著分解水平的增加，細節(jié)部分逐漸被提取出來。經(jīng)過5分解水平小波變換后，可以滿足礦業(yè)安全生產(chǎn)事故時間序列分解的要求。

2.2 分量重組

基于礦業(yè)安全生產(chǎn)事故時間序列進行小波變換后的子序列，計算給定序列與子序列的樣本熵值，各分量的樣本熵如圖4所示。

圖4 各分量的樣本熵

根據(jù)子序列樣本熵的分布情況，取λ=0.2。d1為隨機分量；d2～d4重組為細節(jié)分量s；d5、a5重組為趨勢分量r。重組后各分量結(jié)果如圖5所示。

圖5 重組后各分量

2.3 重組分量的信息?；?/h3>

由圖5可知，樣本熵重組后的序列具有隨機性、波動性等特點。為此，先將原始數(shù)據(jù)按照一定窗口進行模糊信息?；瑢⒌V業(yè)安全生產(chǎn)態(tài)勢時間序列的萬人死亡率模糊信息?；癁長ow、R、Up共3個參數(shù)，3個子序列的模糊信息?；Y(jié)果如圖6所示。該例中Low具體為窗口化子序列每月變化的最小值，R具體為礦業(yè)安全生產(chǎn)態(tài)勢萬人死亡率每月變化的平均值，Up具體為礦業(yè)安全生產(chǎn)態(tài)勢萬人死亡率每月變化的最大值。

(a)d2序列?；?/p>

(b)d3序列?；?/p>

(c)d4序列粒化

2.4 基于FCM的窗口化子序列模糊預(yù)測

對于窗口化的時間序列，按照窗口模糊信息?；Y(jié)果將其分成6個數(shù)據(jù)集，根據(jù)每年數(shù)據(jù)情況確定討論區(qū)間的上下界計算模型的論域。基于樣本數(shù)據(jù)求得的隸屬度值見表6，樣本模糊化后的結(jié)果如表7 所示。

表6 樣本隸屬度

表7 樣本模糊關(guān)系

由表6和表7并結(jié)合式(9)計算，求得測試集預(yù)測結(jié)果，如表8所示。

表8 測試集預(yù)測結(jié)果

合并窗口化后，礦業(yè)安全生產(chǎn)事故時間序列的模糊區(qū)間預(yù)測結(jié)果如圖7所示。

圖7 窗口化區(qū)間預(yù)測結(jié)果

3 模型的驗證

利用訓練集150組樣本建立的模型，對2018年7月至2019年10月測試集的21組數(shù)據(jù)進行預(yù)測，并計算區(qū)間上下邊界值與中值的預(yù)測值，結(jié)果如表9 所示。

表9 測試集預(yù)測誤差對比

由表9可知：R值的平均相對誤差為17.737 1%，Up值的平均相對誤差為8.771 2%。Up值的預(yù)測精度明顯高于R值的預(yù)測精度，表明模糊信息粒化后的上下邊界值能夠較為準確地描述礦業(yè)安全生產(chǎn)事故的變化范圍。

將礦業(yè)安全生產(chǎn)事故時間序列的模糊區(qū)間預(yù)測模型與小波-SVM(支持向量機)、小波-ARMA(自回歸滑動平均模型)和小波神經(jīng)網(wǎng)絡(luò)的預(yù)測效果進行對比，各模型的預(yù)測誤差分析如表10所示。其中，MAE為預(yù)測值和實際值之間的平均絕對誤差；MRE為絕對誤差與真實值的比率，即平均相對誤差；RMSE為均方根誤差，其對異常值更為敏感。

表10 各模型預(yù)測誤差對比分析結(jié)果

由表10可知，本文預(yù)測模型的MAE、MRE、RMSE分別為0.445 7%、13.254 1%、0.372 3%，各項指標均低于對比模型，表明本文模型預(yù)測效果較為穩(wěn)定，精度較高。

4 結(jié)論

1)提出一種礦業(yè)安全生產(chǎn)事故時間序列的模糊區(qū)間預(yù)測模型，能夠?qū)⒌V業(yè)系統(tǒng)中的隨機因素等考慮進去，相較于定量值的點預(yù)測，區(qū)間預(yù)測的結(jié)果更為可靠，更具有實用價值。

2)R、Up值時間序列預(yù)測的平均相對誤差分別為17.737 1%、8.771 2%。邊界值的預(yù)測精度較高，表明模型的區(qū)間范圍較為可靠，也進一步地說明了模糊信息?；蟮倪吔缰到Y(jié)果的合理性。

3)相對于小波-支持向量機、小波-自回歸滑動平均模型和小波神經(jīng)網(wǎng)絡(luò)的預(yù)測結(jié)果，本文的預(yù)測模型平均絕對誤差(MAE)、平均相對誤差(MRE)、均方根誤差(RMSE)3個指標均為最小，表明訓練得到的模型具有較高的可靠性，對礦業(yè)安全決策具有一定的指導意義。