基于多方法耦合的礦井突水水源識別模型研究

宋立兵，董東林，王振榮，李 果，楊茂林

(1.神東煤炭集團(tuán)公司，陜西 神木 719315；2.中國礦業(yè)大學(xué)(北京)地球科學(xué)與測繪工程學(xué)院，北京 100083)

隨著煤礦的開挖深度的不斷增加，礦區(qū)內(nèi)水文地質(zhì)條件愈發(fā)復(fù)雜，礦井突水事故已成為威脅煤礦安全生產(chǎn)的5大災(zāi)害之一[1-4]，對礦井突水水源進(jìn)行快速、客觀、準(zhǔn)確的判別也成為了煤礦開采至關(guān)重要的環(huán)節(jié)[5]。保德礦地處黃河?xùn)|岸的山西省保德縣境內(nèi)，井田內(nèi)地形切割嚴(yán)重，總體地勢呈“V”字，中部低、南北高，礦區(qū)與佘家梁、張家溝、炭嶺溝、河石盤等小煤礦相臨。保德礦主采8#煤層采掘活動在奧灰水位以下，屬于帶壓開采，潛在威脅著8#煤層的安全開采，同時(shí)隨著礦井采空區(qū)面積增加，礦井老空水對保德煤礦未來深部煤層的采掘也具有一定威脅。因此，本文選用保德礦為研究示范區(qū)，以期在保德礦及其相似礦區(qū)礦井開采中發(fā)生突水時(shí)提供有效可靠的突水水源識別模型，保證未來開采工作安全、穩(wěn)定的進(jìn)行。

近年來，眾多學(xué)者基于數(shù)學(xué)函數(shù)理論，采用統(tǒng)計(jì)學(xué)、機(jī)器學(xué)習(xí)等相關(guān)模型對突水水源識別做出了大量研究[6-13]。其中，支持向量機(jī)(SVM)[8]、BP神經(jīng)網(wǎng)絡(luò)10]、極限學(xué)習(xí)機(jī)(ELM)[11]、Fisher判別分析法[12]、Bayes判別分析法[13]等方法為比較常用的方法，但單一的判別模型在應(yīng)用中均存在一定的局限性，如SVM核函數(shù)的不同會使數(shù)據(jù)呈現(xiàn)不同的分布結(jié)構(gòu)，進(jìn)而影響SVM模型的判別結(jié)果；BP神經(jīng)網(wǎng)絡(luò)收斂速度慢，且初始參數(shù)設(shè)置容易造成模型陷入局部最優(yōu)解；ELM雖通過模型迭代可獲得唯一的輸出權(quán)值，克服了傳統(tǒng)的BP神經(jīng)網(wǎng)絡(luò)的不足，但其隨機(jī)選取的初始連接權(quán)值和閾值可能造成輸出結(jié)果誤差、模型網(wǎng)絡(luò)不穩(wěn)定等問題；Fisher判別分析法對樣本結(jié)構(gòu)和樣本間關(guān)聯(lián)性要求較高；Bayes判別分析法在樣本各離子間相關(guān)性較小時(shí)才具有較好的判別效果。目前，已有學(xué)者對現(xiàn)有突水水源識別模型做出有效改進(jìn)，如陳志遠(yuǎn)等[14]建立了Fisher-SVM的線性降維與非線性組合模型，可有效地提取原始樣本數(shù)據(jù)特征，提升水源識別的準(zhǔn)確率；毛志勇等[15]在利用核主成分分析(KPCA)對原始水樣數(shù)據(jù)進(jìn)行降維后將粒子群算法(MPSO)與LM結(jié)合進(jìn)行突水水源判別，有效地提高的識別模型的整體性能；黃平華等[16]建立了Piper-PCA-Fisher識別模型并用于焦作煤礦，結(jié)果表明該方法在水源識別中具有較高的準(zhǔn)確性；鄧清海[17]將PCA與Bayes結(jié)合，對水化學(xué)參數(shù)進(jìn)行處理，確定了鶴壁礦區(qū)的水源樣本數(shù)據(jù)主成分，為該突水水源防治工作提供了依據(jù)。上述方法雖在實(shí)踐工作中的到了有效的實(shí)證運(yùn)用，但在綜合分析礦區(qū)水文地質(zhì)概況、各含水層水化學(xué)特征以及各離子間內(nèi)在聯(lián)系基礎(chǔ)上建立突水水源判別模型的研究卻相對較少。本文建立了耦合Piper-PCA-OT-Regression-Bayes的突水水源判別方法，不僅利用Piper三線圖和PCA對研究區(qū)水源樣本進(jìn)行了水化學(xué)特征分析和水樣主成分確定，并且引入OT-Regression對待測樣本的離群值進(jìn)行校正，以增加樣本間的區(qū)分度，最后采用Bayes判別分析法有效地對保德礦待測水樣進(jìn)行了分類識別。

1 數(shù)據(jù)獲取

表1 保德礦訓(xùn)練樣本水樣數(shù)據(jù)

表2 保德礦待測樣本水樣數(shù)據(jù)

2 研究方法

2.1 主成分分析法

主成分分析法即利用數(shù)學(xué)降維的思想，將多個(gè)指標(biāo)轉(zhuǎn)化為少數(shù)幾個(gè)可反映原始變量信息的綜合指標(biāo)的方法，降維后所得的綜合指標(biāo)即為主成分[18]。其相關(guān)理論及計(jì)算步驟如下[19]：

將原始數(shù)據(jù)矩陣X的p個(gè)向量X1、X2、…、Xp的線性組合為Y=AX，即：

式中，ai1+ai2+ai3+，…，+aip=1；Yi與Yj(i≠j；i、j=1，2，…，p)之間不相關(guān)。

1)將原始變量數(shù)據(jù)標(biāo)準(zhǔn)化處理，然后計(jì)算各變量之間的協(xié)方差矩陣Σ=Zij。

2)將計(jì)算好的協(xié)方差矩陣特征值按從大到小排列，即λ1≥λ2≥，…，≥λp，相應(yīng)的單位特征向量為T1、T2、…、Tp，則第i個(gè)主成分Yi的方差就等于Σ的第i大特征根λi。

3)第k個(gè)主成分的Yk的方差貢獻(xiàn)率為：

(2)

4)第m(m

5)若第m個(gè)主成分的累積貢獻(xiàn)率達(dá)到85%以上[20]，則可視為主成分的個(gè)數(shù)為m。

2.2 貝葉斯判別法

貝葉斯判別法多用于處理多元分布的數(shù)據(jù)，該方法通過計(jì)算假設(shè)概率的方法來達(dá)到判別的目的。首先根據(jù)貝葉斯公式，綜合未知參數(shù)的先驗(yàn)信息與已知的樣本信息，進(jìn)而得出后驗(yàn)信息，最后用于推斷出未知參數(shù)[17，20]。計(jì)算步驟如下[21]：

從G個(gè)樣本集中抽取了n個(gè)樣本，這n個(gè)樣本有p個(gè)變量，即可以假設(shè)存在一個(gè)p維空間，這n個(gè)樣本為空間中存在的n個(gè)離散的點(diǎn)。在對樣本進(jìn)行分類時(shí)，如果出現(xiàn)歸類錯(cuò)誤便產(chǎn)生損失量，如果當(dāng)某一未知樣本X劃分為任意類別時(shí)都會產(chǎn)生損失，但是未知樣本在類別為Ag時(shí)，產(chǎn)生的損失量最小，那么即可認(rèn)為這個(gè)未知樣本為Ag類。

如果已知g個(gè)母體的概率分布為fg(x)，設(shè)將一本應(yīng)該屬于母體Ag的未知樣本錯(cuò)誤的劃歸為母體Ah的概率是：P{h/g}，則下式成立：

此時(shí)產(chǎn)生的平均損失為：

設(shè)每個(gè)母體的先驗(yàn)概率ph為已知，則G個(gè)母體錯(cuò)誤歸類的平均損失為：

同理，設(shè)將一本應(yīng)該屬于母體Ah的未知樣本錯(cuò)誤的劃歸為母體Ag時(shí)產(chǎn)生的損失記為：L(g/h)，則：

當(dāng)Bayes推導(dǎo)出劃分空間{R}滿足的一下條件時(shí)，可將樣本錯(cuò)誤劃歸時(shí)產(chǎn)生的損失量降到最低：

因此，當(dāng)將未知樣本劃歸為Ag的后驗(yàn)概率比該樣本被劃歸為其他類別的后驗(yàn)概率大時(shí)，即可認(rèn)為該樣本屬于Ag。qgfg(x)最大即表明后驗(yàn)概率最大。故可以得出，判別任意一個(gè)樣本屬于某一類別時(shí)的判別函數(shù)表達(dá)式如下：

(13)

exp[-1/2(x-ag)′∑-1(x-ag)]

(14)

式中，x=(x1，x2，…，xp)′；ag為Ag的均值；∑為協(xié)方差矩陣；g=1，2，…，G。

進(jìn)一步計(jì)算可得正態(tài)母體多類線性判別函數(shù)：

Z(x)=b0g+b1gx1+…+bpgxp

(15)

式中，xp為樣本集中的第p個(gè)指標(biāo)的數(shù)值；bpg為判別系數(shù)。

最后比較Z1(x)至Zg(x)的數(shù)值大小，最大值所對應(yīng)的類別即為待測樣本所屬類別。

2.3 離群值檢驗(yàn)

離群值檢驗(yàn)是利用樣本數(shù)據(jù)的最小值、下四分位數(shù)、中位數(shù)、上四分位數(shù)、最大值五個(gè)統(tǒng)計(jì)量來繪制樣本的箱線圖，如圖1所示，進(jìn)而比較樣本的對稱性、分散程度等信息，達(dá)到描述樣本數(shù)據(jù)的目的。計(jì)算步驟如下：

圖1 箱線圖

1)將水樣中各離子的觀測值從小到大排列。

2)依據(jù)觀測數(shù)據(jù)，分別計(jì)算各離子的下四分位數(shù)(Q1)、中位數(shù)(Q2)、上四分位數(shù)(Q3)。

3)根據(jù)上四分位數(shù)與下四分位數(shù)，計(jì)算盒子長度(IQR)：

IQR=Q3-Q1

(16)

4)計(jì)算最大值max和最小值min：

max=Q3+1.5×IQR

(17)

min=Q1-1.5×IQR

(18)

若樣本中存在觀測值超過上四分位數(shù)加1.5倍四分位差，或者小于下四分位數(shù)減1.5倍四分位差，則視為離群值，在箱線圖中被單獨(dú)以點(diǎn)匯出。

2.4 回歸填補(bǔ)法

回歸填補(bǔ)法可通過水樣中各離子的相關(guān)系數(shù)矩陣，建立線性回歸模型，校正水樣中異常值，達(dá)到增加樣本間區(qū)分度的目的。具體步驟如下：

1)計(jì)算水樣的相關(guān)系數(shù)矩陣。

2)根據(jù)箱線圖和相關(guān)系數(shù)矩陣的確定自變量。

3)利用線性模型擬合函數(shù)，建立自變量與其余各指標(biāo)間線性回歸模型。

4)利用線性回歸模型求出各指標(biāo)的回歸值并填入原始數(shù)據(jù)，完成異常值的校正。

2.5 判別方法步驟

綜合判別模型步驟如下：

1)選取保德礦水樣實(shí)測數(shù)據(jù)，繪制Piper三線圖確定各含水層水樣水化學(xué)特征，剔除偏離水樣中心的異常水樣，剩余水樣即為各含水層代表性水樣。

2)利用PCA將樣本進(jìn)行降維處理，選取代表水樣的主成分?jǐn)?shù)，并進(jìn)行歸一化處理，將其作為綜合判別模型的訓(xùn)練集。

3)采用箱線圖離群值檢驗(yàn)法對待測水樣進(jìn)行異常值檢測，找出待測水樣中的異常值并利用回歸填補(bǔ)法建立測試水樣各主成分間線性回歸函數(shù)，進(jìn)而對測試水樣進(jìn)行異常值校正。

4)建立Bayes判別函數(shù)，分別對異常值處理前后的待測數(shù)據(jù)進(jìn)行歸類判別，得出模型判別結(jié)果并比較其優(yōu)越性。判別模型流程如圖2所示。

圖2 判別模型流程

3 結(jié)果與討論

3.1 水化學(xué)分析

圖3 各含水層水樣Piper三線圖

由圖3(a)可知，采空區(qū)水源樣本中第2、5、8三組水樣數(shù)據(jù)偏離含水層水樣中心，與其他水樣水質(zhì)類型差別較大，故視其為異常水樣；由圖3(b)可知，二疊系砂巖含水層水源樣本中第18組水樣偏離地層中心，與其他水樣水質(zhì)類型差別較大，故視其為異常水樣；由圖3(c)可知，石炭系砂巖含水層水源樣本中分布較為集中，可全部作為石炭系砂巖含水層代表水樣；由圖3(d)可知，奧灰含水層水源樣本中第37、39組水樣偏離地層中心，與其他水樣水質(zhì)類型差別較大，故視其為異常水樣。上述6組異常水樣排除后，剩余水樣即為個(gè)含水層代表水樣，可作為水源判別模型的訓(xùn)練樣本集進(jìn)行水源判別。

3.2 主成分分析

為減小測試水樣中的大量冗雜信息對模型判別準(zhǔn)度的影響，采用主成分分析法對水化學(xué)分析后所得各含水層代表水樣進(jìn)行主成分分析。測試樣本六大水化學(xué)指標(biāo)的協(xié)方差矩陣和解釋方差率分別見表3、表4。由表3可知，各水化學(xué)指標(biāo)間存在明顯的相似性，即各水化學(xué)指標(biāo)間存在多余重疊信息，對水源樣本進(jìn)行主成分分析，不僅可減少計(jì)算量，更能排除冗余信息對判別模型的影響。由表4可知，前四個(gè)水化學(xué)指標(biāo)的方差貢獻(xiàn)率較大，其累積方差貢獻(xiàn)率可達(dá)94.817%，即前四種水化學(xué)指標(biāo)對整體水樣的影響較大，可作為綜合指標(biāo)對整體樣本進(jìn)行有效概括。因此，在后續(xù)的模型中，將前四種水化學(xué)指標(biāo)作為水源樣本的主成分進(jìn)行水源判別。

表3 各水化學(xué)指標(biāo)協(xié)方差矩陣

表4 各指標(biāo)解釋方差率

3.3 離散值檢驗(yàn)與回歸填補(bǔ)

圖4 待測樣本各離子箱線圖

表5 待測水樣相關(guān)系數(shù)矩陣

表6 異常值校正后待測水樣

Ca2+=-0.1819(Na++K+)+79.6367

Mg2+=0.5427Ca2+-0.4767

3.4 Bayes判別結(jié)果與分析

為確保水源模型準(zhǔn)確性，首先將水化學(xué)特征分析以及主成分分析后確定的代表水樣歸一化處理后進(jìn)行模型回判訓(xùn)練，訓(xùn)練結(jié)果表明，所有訓(xùn)練樣本均能準(zhǔn)確識別。在此基礎(chǔ)上，對模型開展實(shí)證運(yùn)用，分別將保德礦20組校正前后的待測水樣分別進(jìn)行歸一化處理，然后代入Bayes判別模型中進(jìn)行水源判別，判別結(jié)果如圖5所示。由圖5可知，待測樣本未進(jìn)行離群值檢驗(yàn)和回歸值填充時(shí)，共3組待測樣本判別錯(cuò)誤，模型判別準(zhǔn)確率為85%，而校正后的測試結(jié)果表明，僅第4組采空區(qū)水樣判別錯(cuò)誤，其它水樣判別結(jié)果均與實(shí)際水樣類別相符合，判別準(zhǔn)確率高達(dá)95%。由此可見，Piper-PCA-OT-Regression-Bayes相結(jié)合的突水水源判別模型不僅可對大量的水化學(xué)數(shù)據(jù)進(jìn)行分析處理，而且在判別準(zhǔn)確率上也具有明顯的優(yōu)勢，可應(yīng)用于保德礦進(jìn)行突水水源識別。

圖5 模型判別結(jié)果

由上述兩個(gè)模型判別結(jié)果可知，A4組樣本在兩次判別中均出現(xiàn)誤判情況。由圖3(a)和圖3(c)可知，采空區(qū)訓(xùn)練樣本存在部分水樣與石炭系砂巖含水層水樣水化學(xué)離子濃度相近，在Piper三線圖菱形中分布區(qū)域相同，這可能是造成誤判的原因之一。此外，由表1、表4可知，校正前后該組待測水樣數(shù)據(jù)各離子濃度與采空區(qū)其余訓(xùn)練樣本各離子濃度相差較大，而與石炭系砂巖水待測水樣數(shù)據(jù)各離子濃度比較貼近，原因可能是該組水樣離子濃度測試結(jié)果有誤，或者取樣地點(diǎn)的含水層間產(chǎn)生了水力聯(lián)系，所取水樣為兩個(gè)含水層的混合水樣，故造成模型判別錯(cuò)誤。綜上，在今后的研究中，應(yīng)增加更多的水樣數(shù)據(jù)，確保數(shù)據(jù)的可靠性，同時(shí)可在模型中增加其他相關(guān)水化學(xué)離子信息以及pH值、總堿度、礦化度等相關(guān)指標(biāo)，進(jìn)而增加水樣的辨識度，保證模型判別結(jié)果的準(zhǔn)確率。

4 結(jié) 論

1)分析了保德礦各含水層水樣的水化學(xué)特征，確定采空區(qū)水質(zhì)類型為HCO3-Ca·Na·Mg型；二疊系砂巖含水層水質(zhì)類型為HCO3-Na型；石炭系砂巖含水層水質(zhì)類型同樣為HCO3-Na型；奧灰含水層水質(zhì)類型為HCO3-Ca型。同時(shí)，根據(jù)各含水層水質(zhì)類型排除了訓(xùn)練樣本中異常水樣，初步保證了訓(xùn)練樣本數(shù)據(jù)的可靠性。

4)利用Bayes判別模型對線性回歸模型校正前后的水樣數(shù)據(jù)分別進(jìn)行水源判別，判別結(jié)果準(zhǔn)確率分別為85%和95%，即本文建立的突水水源綜合判別模型判別結(jié)果準(zhǔn)確率高，具有一定應(yīng)用價(jià)值。