劉　棟，李夕兵，劉志祥，董隴軍，周勇勇，陳光輝

(中南大學(xué) 資源與安全工程學(xué)院，湖南 長沙 410083)

0　引言

采礦活動(dòng)向深部發(fā)展，破壞了深部的原巖應(yīng)力，由此引起的地質(zhì)災(zāi)害狀況頻發(fā)，微震監(jiān)測技術(shù)是目前監(jiān)測災(zāi)害的有效手段[1-2]。近幾年，微震監(jiān)測技術(shù)已經(jīng)在臺(tái)站布網(wǎng)優(yōu)化[3-4]、P波[5-6]、S波拾取[7]、微震定位[8-9]及微震爆破事件識(shí)別[10]等方面有了較大的進(jìn)展，并取得了一系列重要的研究成果。為了掌握巖體破壞的規(guī)律，對(duì)微震事件的聚集特征及活動(dòng)規(guī)律研究顯得尤為重要。在此方面，國內(nèi)外學(xué)者做了許多創(chuàng)新性的研究，例如吳愛祥[11]等利用最短距離聚類法對(duì)某礦山微地震活動(dòng)的時(shí)空分布進(jìn)行了分析，識(shí)別了井下微震活動(dòng)聚集區(qū)域；謝和平[12]等采用分形理論研究微震活動(dòng)的時(shí)空分布，發(fā)現(xiàn)分形維數(shù)在一定程度上反映了微震活動(dòng)的變化狀態(tài)；Martin[13]利用空間聚類方法對(duì)微震事件進(jìn)行聚類分析并解釋了微震活動(dòng)規(guī)律?？v觀這些研究仍存在一定的不足：在劃定微震活動(dòng)的聚集區(qū)域時(shí)只考慮微震事件的空間分布特征，沒有將時(shí)間特征考慮在內(nèi)。對(duì)此，本文采用時(shí)空共享近鄰聚類(STSNN)算法，將微震事件的時(shí)間空間特征相結(jié)合并展開聚類研究。2012年Liu等[14]將共享近鄰概念引入時(shí)空聚類，提出STSNN聚類算法，該算法采用實(shí)體的有效近鄰數(shù)目來定義其時(shí)空密度，更適于發(fā)現(xiàn)不同密度的時(shí)空簇，是時(shí)空聚類算法中較先進(jìn)的聚類算法。

開磷用沙壩礦存在嚴(yán)重的潛在地壓問題，于2013年安裝使用IMS微震監(jiān)測系統(tǒng)來監(jiān)測井下災(zāi)害。本文以用沙壩礦為工程背景，以用沙壩礦IMS微震監(jiān)測系統(tǒng)2014年收集到微震數(shù)據(jù)作為研究對(duì)象，通過P波、S波拾取后得到5 425個(gè)有效微震作為本文的數(shù)據(jù)集。采用STSNN算法并選擇合適的參數(shù)對(duì)其進(jìn)行聚類分析，計(jì)算其聚類有效性，并對(duì)最大類簇的微震事件活動(dòng)規(guī)律進(jìn)行探究，根據(jù)其變化規(guī)律發(fā)現(xiàn)產(chǎn)生較大事件的預(yù)警點(diǎn)。通過對(duì)微震事件活動(dòng)規(guī)律的研究可為大事件的產(chǎn)生提供有效的預(yù)判信息，為保證礦山安全生產(chǎn)發(fā)揮重要的作用。

1　時(shí)空共享近鄰聚類(STSNN算法)

1.1　基本定義

為了闡述時(shí)空近鄰共享聚類算法，首先介紹幾個(gè)重要的概念：

定義1時(shí)空事件：一個(gè)時(shí)空事件STi是具有空間坐標(biāo)和確定的時(shí)間標(biāo)志，表示為STi(Li,ti)。

定義2時(shí)空近鄰：時(shí)空事件STi的時(shí)空近鄰定義為以ΔD作為底面半徑、2ΔT作為高的圓柱體，ΔT為時(shí)間窗口。

定義3時(shí)空密度：時(shí)空事件STi的時(shí)空密度定義為STi時(shí)空直接可達(dá)的事件數(shù)量。

定義4時(shí)空可達(dá)：假如時(shí)空事件數(shù)據(jù)集ST1,ST2,…STi,…STn，STn是一個(gè)時(shí)空核心，如果STi-1是從STi時(shí)空直接可達(dá)，那么ST1也可視為從STn時(shí)空可達(dá)，表示為STn→ST1。

定義5噪聲：假設(shè)時(shí)空事件STi不是一個(gè)時(shí)空核心，在精確時(shí)空近鄰內(nèi)不存在時(shí)空核心，那么可將STi視為噪聲。

定義6時(shí)空類簇：對(duì)于一個(gè)時(shí)空數(shù)據(jù)集STD，時(shí)空類簇STC是滿足下列2個(gè)條件的非空集。

條件1：?STi∈STD,STj∈STD:STi→STj

條件2：?STi∈STD,STj∈STD:STi～STj

1.2　STSNN算法及參數(shù)選擇

STSNN算法需要用戶設(shè)置4個(gè)參數(shù)：k,kT,MinPts和ΔT，k是時(shí)空近鄰的事件數(shù)量，kT是判別兩個(gè)事件是否時(shí)空直接可達(dá)的閾值，MinPts用來識(shí)別時(shí)空核心，ΔT是定義時(shí)空近鄰的時(shí)間窗口。基于上述定義，STSNN算法的執(zhí)行可分為以下5個(gè)步驟：

Step1：對(duì)每個(gè)時(shí)空事件確定時(shí)空近鄰k。

Step2：對(duì)每個(gè)時(shí)空事件搜索其時(shí)空共享近鄰。

Step3：計(jì)算每個(gè)時(shí)空事件的時(shí)空密度以及進(jìn)一步探測時(shí)空核心。

Step4：擴(kuò)大時(shí)空類簇。選擇一個(gè)時(shí)空核心并使用遞歸算法將所有時(shí)空可達(dá)和時(shí)空相連的事件增加到時(shí)空類簇。

Step5：識(shí)別噪聲事件。

STSNN參數(shù)kT和MinPts通常設(shè)置為0.5k。另外2個(gè)參數(shù)，k和ΔT在本次研究中，ΔT初選值設(shè)定為[5,6]，k的初選值[5,6,7…17,18,19]，共有30種組合，為了方便后面表達(dá)，每一種組合的類別符合的如表1所示。

表1　30種組合的類別符號(hào)

2　聚類結(jié)果及有效性評(píng)價(jià)

分別利用STSNN算法對(duì)上述k和ΔT的30種組合進(jìn)行聚類運(yùn)算。當(dāng)k和ΔT分別選擇不同的值時(shí)，可以得到不同的聚類類簇，如圖1所示，當(dāng)ΔT一定k增大時(shí)，得到的類簇?cái)?shù)量越來越少，而最大類簇事件數(shù)沒有呈現(xiàn)規(guī)律性變化。圖2表示不同類別的噪聲率，除ID1(k=5，ΔT=5)和ID16(k=5，ΔT=6)外，其他類別的噪聲率相對(duì)較低。

圖1　不同類別的簇?cái)?shù)及最大類簇事件數(shù)Fig.1　The number of clusters of different classes and the maximum number of cluster events

圖2　不同類別的噪聲率Fig.2　Noise rates of different classes

聚類有效性評(píng)價(jià)是為了確定對(duì)數(shù)據(jù)劃分多少個(gè)類簇是最合適的。一個(gè)理想的聚類結(jié)果應(yīng)該滿足2方面要求[15]:(1)緊密性，即簇內(nèi)實(shí)體應(yīng)盡可能相似；(2)分離性，即不同簇間實(shí)體差異盡可能大。一個(gè)好的聚類結(jié)果應(yīng)該同時(shí)滿足緊密性與分離性。本文選擇RMSSDT指數(shù)[16]和PBM指數(shù)[17]作為有效性評(píng)價(jià)指標(biāo)。RMSSDT指數(shù)旨在度量實(shí)體間的離散程度，其值越小，則表示簇的均質(zhì)性越好；PBM指數(shù)采用融合多個(gè)度量指標(biāo)的策略來進(jìn)行聚類結(jié)果的有效性度量，其值越大，表示聚類結(jié)果越可靠，如圖3所示。

注：折線代表PBM；柱狀圖代表RMSSDT圖3　聚類有效性評(píng)價(jià)結(jié)果Fig.3　Clustering effectiveness evaluation results

利用RMSSDT指數(shù)和PBM指數(shù)對(duì)每一類別的有效性進(jìn)行評(píng)價(jià)，如圖3所示單純從PBM評(píng)價(jià)指標(biāo)(折線)來看，ID17(k=6，ΔT=6)、ID19(k=8，ΔT=6)和ID21(k=10，ΔT=6)聚類效果較好，優(yōu)勢較明顯，進(jìn)而從RMSSDT指數(shù)(柱狀圖)來比較，ID17(k=6，ΔT=6)聚類結(jié)果最好，并且其噪聲率不高，因此本文后續(xù)研究將利用ID17(k=6，ΔT=6)的聚類結(jié)果。該聚類共得到98個(gè)類簇，最大類簇的事件數(shù)為544，噪聲率為24.8%。最大的5個(gè)類簇事件數(shù)分別為544，246，234，161和100。圖4表示最大類簇在礦山中的主要聚集區(qū)域，從圖中可以發(fā)現(xiàn)，該類簇主要集中在用沙壩礦的斷層附近。

圖4　最大類簇微震事件的主要聚集區(qū)域Fig.4　The events aggregation area of the biggest cluster group

3　微震事件活動(dòng)特征

3.1　微震事件24 h分布

圖5　微震事件24 h分布(m表示震級(jí))Fig.5　The distribution of micro-seismic in 24 hours

圖5表示該區(qū)域聚類事件的24 h分布圖，微震事件主要發(fā)生在每天的11時(shí)到15時(shí)，大于0級(jí)的事件也主要發(fā)生在此時(shí)間段。用沙壩礦爆破時(shí)間主要集中在每天11時(shí)到14時(shí)。因此該區(qū)域的微震事件大部分分布在斷層區(qū)域，受采礦活動(dòng)誘發(fā)而產(chǎn)生。在非爆破時(shí)間，該區(qū)域的微震事件在時(shí)間分布上較為均勻，但偶爾也有震級(jí)大于-0.5級(jí)事件發(fā)生。

3.2　微震事件活動(dòng)率與視體積

微震事件活動(dòng)率反映了巖體內(nèi)部微裂紋的擴(kuò)展變化趨勢，表現(xiàn)了微裂紋的產(chǎn)生和發(fā)展速度，微震事件的突變特征體現(xiàn)了巖體開挖導(dǎo)致圍巖破壞的前兆信息。微震事件的視體積表征地震非彈性變形的體積與它釋放的能量和體變勢相關(guān)，累計(jì)視體積隨時(shí)間變化曲線的斜率通常被認(rèn)為是表示巖體應(yīng)變速率的重要指標(biāo)，曲線的斜率變化反映應(yīng)變率的微小變化，巖體失穩(wěn)變形的前兆是經(jīng)歷一段時(shí)間的加速變形過程。通過觀察圖6，可以按微震活動(dòng)率分為3個(gè)階段，AB階段微震事件活動(dòng)率降低，說明此階段可能有大事件發(fā)生，能量得到釋放；BC段微震活動(dòng)率變化緩慢，說明此階段巖體處于穩(wěn)定狀態(tài)，開始生成新的裂隙并呈穩(wěn)定發(fā)展階段；CD階段為巖體內(nèi)裂隙發(fā)展迅速，微震活動(dòng)率增加，因此可能發(fā)生較大事件。累積視體積曲線在B點(diǎn)(5月3日附近)和C(5月10日附近)點(diǎn)明顯上升，表明B，C點(diǎn)附近可能有破壞發(fā)生，導(dǎo)致巖體失穩(wěn)。在CD段開始急劇上升，表明此階段大量的微裂隙張開占據(jù)巖石空間導(dǎo)致巖石體積膨脹，視體積曲線突然增大的點(diǎn)可以作為發(fā)生大事件的前兆信息，因此在D(5月16日)點(diǎn)之后，很可能產(chǎn)生較大的破壞性事件。實(shí)際上在5月19日現(xiàn)場有一次震級(jí)為0.7級(jí)的微震事件發(fā)生。因此可以認(rèn)為微震事件率急劇下降并且累積視體積曲線忽然上升是巖石發(fā)生破壞的預(yù)兆現(xiàn)象。

圖6　微震活動(dòng)率與累積視體積事件序列曲線Fig.6　Time series curve of cumulative apparent volume and event activity

3.3　施密特?cái)?shù)與勁度系數(shù)

施密特?cái)?shù)(Schmidt Number)是運(yùn)動(dòng)學(xué)黏度與擴(kuò)散率之比，能測量震動(dòng)巖石流變的時(shí)空復(fù)雜程度，也是反映巖石不穩(wěn)定的參數(shù)，施密特?cái)?shù)越低，穩(wěn)定性越差。勁度系數(shù)(Stiffness Modulus)表征系統(tǒng)阻止地震變形和應(yīng)力增加的能力。如圖7所示施密特?cái)?shù)曲線與勁度系數(shù)曲線保持相同的變化趨勢，分別在5月2日、5月9和5月16日附近急劇上升然后迅速下降，降低程度幾乎相同，施密特?cái)?shù)處于較低水平時(shí)該區(qū)域的巖體穩(wěn)定性較差，產(chǎn)生破壞的可能性比較大。在這3個(gè)時(shí)間，從勁度系數(shù)曲線可以看出，巖體阻止系統(tǒng)變形的能力也瞬間下降，進(jìn)而也表明該區(qū)域在這3個(gè)時(shí)間較容易發(fā)生危險(xiǎn)事故。而實(shí)際現(xiàn)場監(jiān)測中，5月10日期間，發(fā)生了2次震級(jí)超過0.5級(jí)的微震事件，與5月9日附近施密特?cái)?shù)勁度系數(shù)先升后降的預(yù)測基本對(duì)應(yīng)。因此可以把施密特?cái)?shù)和勁度系數(shù)先升后降的點(diǎn)作為大事件將要發(fā)生的預(yù)警信息。

圖7　施密特?cái)?shù)與勁度系數(shù)事件序列曲線Fig.7　Time series curve of Schmidt Number and Stiffness Modulus

4　結(jié)論

1)利用時(shí)空共享近鄰(STSNN)聚類算法對(duì)用沙壩礦微震事件進(jìn)行聚類分析，可以發(fā)現(xiàn)微震事件的時(shí)空聚集類簇，并對(duì)其聚類有效性進(jìn)行分析，確定聚類參數(shù)選擇的合理性。對(duì)微震事件聚類分析來確定所研究的區(qū)域，比人為劃定的事件聚集區(qū)域更客觀有效。

2)對(duì)聚類最大類簇的微震事件活動(dòng)特征進(jìn)行分析，通過24 h分布圖可以發(fā)現(xiàn)用沙壩礦微震事件主要發(fā)生在每天的11時(shí)到15時(shí)，大于0級(jí)的時(shí)間也主要發(fā)生在此時(shí)間段，微震活動(dòng)率急劇下降并且累積視體積曲線忽然上升、施密特?cái)?shù)和勁度系數(shù)先升后降的點(diǎn)作為巖體失穩(wěn)發(fā)生破壞性事件的預(yù)警點(diǎn)。

[1]李夕兵. 巖石動(dòng)力學(xué)基礎(chǔ)及應(yīng)用[M]. 北京，科學(xué)出版社，2014.

[2]王春來,吳愛祥,劉曉暉.深井開采巖爆災(zāi)害微震監(jiān)測預(yù)警及控制技術(shù)[M].北京,冶金工業(yè)出版社,2013.

[3]周勇勇,李夕兵,劉志祥,等. 高精度三維地下微震臺(tái)網(wǎng)模糊優(yōu)選方法[J]. 中國安全生產(chǎn)科學(xué)技術(shù),2016, 12(7):82-86.

ZHOU Yongyong, LI Xibing, LIU Zhixiang, et al.A fuzzy optimum approach for high precision three-dimensional underground micro-seismic Network [J]. Journal of Safety Science and Technology, 2016, 12(7):82-86.

[4]張楚旋,李夕兵,董隴軍,等. 微震監(jiān)測傳感器布設(shè)方案評(píng)價(jià)模型及應(yīng)用[J]. 東北大學(xué)學(xué)報(bào)(自然科學(xué)版),2016(4):594-598,608.

ZHANG Chuxuan, LI Xibing, DONG Longjun, et al. Evaluation model of micro-seismic monitoring sensor layout scheme and its application[J].Journal of Northeastern University (Natural Science), 2016(4):594-598,608.

[5]LI X, SHANG X, WANG Z, et al. Identifying P-phase arrivals with noise: An improved Kurtosis method based on DWT and STA/LTA[J]. Journal of Applied Geophysics, 2016, 133: 50-61.

[6]LI X, SHANG X, Morales-Esteban A, et al. Identifying P phase arrival of weak events: the Akaike Information Criterion picking application based on the Empirical Mode Decomposition [J]. Computers & Geosciences, 2016.

[7]張楚旋,李夕兵,董隴軍,等. 三函數(shù)四指標(biāo)礦震信號(hào)S波到時(shí)拾取方法及應(yīng)用[J]. 巖石力學(xué)與工程學(xué)報(bào),2015(8):1650-1659.

ZHANG Chuxuan, LI Xibing, DONG Longjun, et al. A S-wave phase picking method with four indicators of three functions for micro-seismic signal in mines [J]. Chinese Journal of Rock Mechanics and Engineering, 2015(8):1650-1659.

[8]Linqi HUANG, Xibing LI, Longjun DONG, et al. Relocation method of micro-seismic source in deep mines [J]. Transactions of Nonferrous Metals Society of China, 2016(11):2988-2996.

[9]張飛, 劉德峰, 張衡,等. 基于IMS微震監(jiān)測系統(tǒng)的微震事件定位精度分析[J]. 中國安全生產(chǎn)科學(xué)技術(shù), 2013, 9(6):21-26.

ZHANG Fei, LIU Defeng, ZHANG Heng, et al. Location accuracy analysis of the seismic events based on the IMS micro-seismic monitoring system [J]. Journal of Safety Science and Technology, 2013, 9(6):21-26.

[10]GUOyan ZHAO, Ju MA, Long-jun DONG, et al. Classification of mine blasts and micro-seismic events using starting-up features in seismograms[J]. Transactions of Nonferrous Metals Society of China, 2015(10): 3410-3420.

[11]吳愛祥,武力聰,劉曉輝,等. 礦山微地震活動(dòng)時(shí)空分布[J]. 北京科技大學(xué)學(xué)報(bào),2012,34(6):609-613.

WU Aixiang，WU Licong，LIU Xiaohui，et al. Space-time distribution of micro-seismic activities in mines[J].Journal of University of Science and Technology Beijing，2012,34(6):609-613.

[12]謝和平,W. G. Pariseau. 巖爆的分形特征和機(jī)理 [J]. 巖石力學(xué)與工程學(xué)報(bào),1993,12(1):28-37.

XIE Heping，Pariseau W G. Fractal character and mechanism of rockbursts [J]. Chinese Journal of Rock Mechanics and Engineering,1993,12(1): 28-37.

[13]Hudyma MR. Analysis and interpretation of clusters of seismic events [D]. University of Western Australia, Perth, Australia, 2008.

[14]Liu Q., Deng M., Bi J., Yang, W. A Novel Method for Discovering Spatio-temporal Clusters of Different Sizes, Shapes and Densities in the Presence of Noise [J].International Journal of Digital Earth，2012, 7(2): 138-257.

[15]Berry M, Linoff G. Data Mining Techniques: For Marking, Sales and Customer Support[M]. New York: Wiley, 1997.

[16]Halkidi M, Batistakis Y, Vazirgiannis M. Clustering algorithm and validity measures[A]. Proceedings of 13th International Conference on Scientific and Statistics Database Management[C]. 2001, 3-22.

[17]Pakhira M K, Bandyopadbyay A, Maulik U. Validity index for crisp and fuzzy clustering[J]. Pattern Recognition, 2004, 37(3):487-501.