應用殘差網絡的微地震事件五分類檢測方法

關鍵詞： 微地震監測，事件檢測，數據增廣，殘差網絡，深度學習

0 引言

近年來，水力壓裂技術被廣泛用于頁巖氣、煤層氣等非常規油氣資源開發[1?5]。在水力壓裂過程中，通過檢測巖石破裂產生的微地震信號，可以反演震源位置、震源機制解、震級等參數，從而監測裂縫的發育過程并評估壓裂效果[6?11]。地面監測和井中監測是兩種常用的監測方式。與地面監測相比，井中監測的檢波器記錄具有較高的信噪比，可通過拾取P 波和S 波到時進行微地震定位[12]。野外記錄通常為連續的波形，對微地震信號的識別直接決定了微地震定位的準確度和可靠性。

目前，微地震事件自動檢測主要應用于天然地震領域，方法包括基于振幅的長短時窗比法（Short TermAverage/Long Term Average ， STA/LTA）[13?16]、基于Akaike 信息準則（Akaike Information Criteria，AIC）的方法[17?18]以及基于波形相似性的波形互相關法[19?22]等。目前，應用較為廣泛的為STA/LTA 算法，其基本原理是通過計算短、長時間窗內波形的平均能量值比值來評估信號的能量變化。該方法原理簡單、計算效率高，地震學家對該方法做了一系列的改進。例如，Withers 等[13]實現了一種自適應窗口的STA/LTA 算法，結合地震圖頻譜特性滿足事件檢測與定位的要求。宋維琪等[14]在STA/LTA 算法的基礎上，討論了基于特征值分析且抗噪能力較強的識別方法。邱磊等[15]通過建立閾值和背景噪聲之間的關系，提出了兩種基于參考閾值拾取地震波初至的方法，拾取結果的正確率比傳統STA/LTA 方法有所提高。劉晗等[16]對STA/LTA 算法和改進的加權系數法、多窗口算法和修正能量比法進行了對比，了解了時窗長度、觸發閾值和特征函數對算法性能的影響特征及其選取規律。STA/LTA 法具有較好的通用性且高效，但其拾取結果受閾值影響較大，依賴人工經驗，同時該方法對中強度背景噪聲十分敏感，容易出現漏檢事件。

近年來，隨著計算機領域人工智能的發展，深度學習相關技術在計算機視覺、自然語言處理、語音識別、醫學圖像處理等多個領域取得了巨大的成功[23?24]。機器學習算法也被廣泛地應用在地震勘探領域，尤其是地震信號的檢測和拾取，發展了一系列的智能算法[25?35]。例如，Mousavi 等[36]采用卷積神經網絡以及循環神經網絡搭建模型，形成了CRED 模型，用于噪聲與地震事件分類，該模型對于背景噪聲、非地震脈沖有著極高魯棒性，能夠高效地處理海量地震數據。趙明等[37]使用汶川地震余震事件構建訓練數據集，利用卷積神經網絡方法進行訓練，并用未參與模型訓練的汶川余震事件作為檢測數據集，其訓練和檢測準確率均達到95% 以上。Zhao 等[38]采用卷積神經網絡對連續數據進行時窗分割并分類，其輸入為一維的地震數據，驗證集精度達到93%。Zhang 等[39]構建了基于卷積神經網絡的微地震檢測網絡，通過比較分析多種微地震檢測分類方法，論證了深度學習算法在微地震事件檢測中具有較高的準確度，展示了在微地震實時監測領域中良好的應用前景。王維波等[40]利用卷積神經網絡算法構建微地震檢測二分類模型，通過水力壓裂實際監測的微地震信號構建訓練樣本集，訓練得到的網絡模型具有較好的抗噪性。隨后，王維波等[41]又提出利用殘差網絡對微地震監測數據的地震發射層析成像數據進行處理， 實現微地震事件自動識別。

深度學習算法可以直接從訓練數據中挖掘地震事件特征。預測新的數據時，只需要將微震事件分割成相互獨立的窗口，然后輸入到訓練好的模型中即可迅速進行事件檢測，節省大量的計算時間。現有的基于深度學習的微地震事件檢測方法通常把事件檢測看作圖像二分類問題，需要對連續監測的波形數據按照一定時窗長度進行分割。通常，相鄰時窗會有部分重疊以盡量減少漏檢。這種處理方法不僅增加了計算時間，降低了運算效率，同時如何選取重疊時窗大小也決定了分類結果的準確度。因此，本文提出一種基于深度學習的微震事件五分類檢測方法。首先，把野外采集的連續記錄波形等分，人工標定為噪聲、完整的P 波和S 波、只含有P 波、只含有S 波和多個微震事件五類； 其次，通過研究數據增廣方法，利用較少的原始標定樣本獲得多樣性的訓練集數據進行模型訓練； 最后，將訓練得到的五分類模型應用于實際微震數據檢測，通過對檢測得到的微震事件進行到時拾取和震源定位評估檢測精度，并與二分類算法結果進行對比驗證。

1 數據處理

本文使用了東方地球物理公司（BGP）提供的一套單井監測的水力壓裂數據，共有12 個三分量檢波器分布在深度為2800～3020 m 的監測井中，間隔為20 m，檢波器位于壓裂井段的上方。圖1 為本文使用的P 波和S 波速度模型，參數經測井數據提取并利用射孔事件校正。數據的時間采樣率為0. 5 ms。

1. 1 樣本標定

本文把微地震事件檢測視為五分類問題：圖2a為噪聲，被歸為第一類，簡稱C1；圖2b 為一個完整的微震事件，包含完整的P 波和S 波信息，被歸為第二類，簡稱C2；圖2c 只含有完整的P 波信息，被歸為第三類，簡稱C3；圖2d 只含有完整的S 波信息，被歸為第四類，簡稱C4；第五類較為復雜，定義為在一個樣本窗口中包含多個微震事件（圖2e）。此外，樣本中同時包含一個微震事件和另一個事件的P 波或者S 波也被歸為第五類，簡稱C5。

1. 2 數據增廣

數據增廣可以提高模型的泛化能力和魯棒性，有助于減少過擬合問題。在訓練神經網絡模型時，通常需要大量的訓練數據，但是大多數地球物理的數據集比較有限。通過數據增廣，能夠擴展數據集，增加訓練樣本的多樣性。常見的數據增廣方法包括旋轉、翻轉、裁剪、平移、縮放、添加噪聲等[23， 37]。針對井下微震事件檢測，本文參考Zhang 等[42]提出的方法進行數據增廣。其中，針對C2 樣本的增廣方法共有5 種。圖3 為原始樣本的三分量數據波形記錄，圖4 為對圖3 數據增廣后的結果。

圖4a 為顛倒檢波器順序：微震事件可能發生在檢波器陣列的上方或下方，顛倒檢波器順序可以模擬發生在不同位置的微震事件。

圖4b 為旋轉分量：微震事件可能發生在不同的方位角和不同的震源機制作用下。不同的方位角導致波形數據中的粒子運動變化。旋轉分量可以模擬具有不同偏振角產生的復雜波形的微地震事件[42]。

圖4c 為移動P 波或S 波相位：P 波與S 波的到時差與檢波器的幾何分布、地下速度模型以及事件發生的位置等因素有關。通過隨機移動P 波和S 波波形來模擬由各種物理因素引起的到時差。本文從? 15 ms 到15 ms 之間隨機選擇一個值進行移動P波到時或者S 波到時。此外，為了防止移動后S 波的到時早于P 波到時，選擇最小P、S 波到時差為25 ms。

圖4d 為調整P 波振幅：調整P 波振幅有助于提高模型的泛化能力。我們選擇對S 波之前的P 波波形隨機乘以值為0. 2～1. 2 的因子。當因子較小時，樣本主要顯示S 波波形，可模擬單震相事件。

圖4e 為添加噪聲：信噪比是影響微地震事件檢測精度的一個重要因素。噪聲類型的多樣性有利于提高網絡模型的檢測準確度。本文選擇了4 種不同類型的噪聲改變事件的信噪比。第一種為實際數據噪聲。首先從原始的監測數據中截取了50 個噪聲窗口，隨機選擇一個噪聲窗口疊加到微震事件波形中。

通過隨機設置信噪比模擬被不同程度噪聲污染的事件。第二種為高斯白噪聲。本文選擇的高斯白噪聲約為信號最大值的5×10?3～5×10?2 倍。第三種添加脈沖噪聲。脈沖噪聲絕對值為1，脈沖出現的概率范圍為10?2～10?3。第四種添加正弦波噪聲，正弦波頻率為3～6 Hz，振幅為0. 05～0. 30 個單位。

C3 與C4 樣本的數據增廣方法通過截取C2 樣本實現。以C2 樣本中所有檢波器中最后的P 波到時與最初的S 波到時的平均值為基點，隨機選擇? 25 ms～25 ms 的擾動，拆分窗口。分別得到一個C3 樣本以及一個C4 樣本。然而，拆分的樣本不滿足時窗長度要求，需要對P、S 波前后的噪聲進行插值或截斷，得到滿足時窗長度要求數據集。

針對C5 的數據增廣方法一共兩種：

（1） 兩個C2 樣本疊加。從C2 樣本中隨機選取兩個事件進行疊加。如果兩個選取的C2 樣本距離較近，則疊加后的樣本難以區分是C2 或是C5，因此本文規定兩個疊加事件的P 波到時差需要超過100 ms。

此外，通過乘以值為0. 5～1. 0 的一個因子改變振幅，更加貼合實際情況。

（2） 一個C2 樣本與一個C3 或者C4 樣本疊加。具體做法是將一個位于時窗前部S 波與一個任意位置的C2 樣本疊加，或者是一個位于時窗后部的P 波與一個任意位置的C2 樣本疊加。

本文從野外原始記錄的波形數據中分別截取5種數據樣本，并通過以上數據增廣方法進行數據增廣，數據集增廣前、后個數如表1 所示。數據集個數從91 增加到22431，5 種數據樣本的比例分別為19. 83%，20. 28%，19. 97%，20. 01% 以及19. 98%。最后將數據集按8：2 的比例劃分為訓練集以及驗證集，訓練集與驗證集個數分別為17944、4487。

2 方法與訓練結果

2. 1 殘差學習模塊

傳統卷積神經網絡層數增加時，易出現梯度消失的問題。He 等[24]在傳統卷積神經網絡的基礎上，增加跨層數據通路（圖5 黑色折線），加強了層與層之間的溝通，該神經網絡被稱為殘差網絡（Residual Net?work，簡稱ResNet）。該網絡主要由殘差塊堆疊而成。圖5 是本文所使用的殘差塊結構，輸入x 依次進入卷積層、批量規范化層（Batch NormalizationLayer， 簡稱BN 層）、ReLU 激活層、卷積層、BN 層以及ReLU 層，同時在最后一個ReLU 層之前，加入了跨層連接，得到映射函數F（x）+x。需要注意的是，殘差塊的輸入x 與最終輸出的數據形狀是相同的。

2. 2 殘差網絡結構

圖6 為本文殘差模型的架構。模型的輸入大小為 n×1001×36（ n 為訓練樣本個數，1001 為時間采樣點，36 為12 個三分量檢波器的道數）的波形數據。其中，第二個通道的排列順序為：第1 個檢波器的N、E、Z 分量；第2 個檢波器的N、E、Z 分量；以此類推。首先，模型第一層是通道數為64、步長為2 的7 × 7 的卷積層，隨后依次進入BN 層、RuLU 層、以及步長為2 的3 × 3 的最大池化層。此時的輸出將進入堆疊的殘差塊。按順序包括3 個通道數為64、3 個通道數為128、3 個通道數為256 和3 個通道數為512 的殘差塊。由于每個殘差塊包含2 個卷積層，再加上第一層7 × 7 的卷積層，該網絡一共包含了25 個卷積層。其中，第1、第8、第14 和第20 層的卷積層步長為2，其余步長均為1。隨后，輸出與全局池化層相連。最后，使用包含5 個神經元的全連接層并使用Softmax 分類器輸出5 個值，分別對應輸入樣本為C1、C2、C3、C4、C5 的概率，其中最大概率值對應的類別即為該網絡的預測結果。

深度學習的損失函數與地球物理反演中的目標函數十分相似，用于估計模型的預測值與真實值的差別程度。在分類問題中，多使用交叉熵損失函數衡量輸入的真實概率分布p 和預測概率分布q 之間的差異H （ p，q），即

式中： xi 表示第i 個數據; N 代表數據總數。損失函數可使用Adam 優化算法迭代求解，每次迭代后，卷積核和權重將會被更新。

2. 3 訓練過程和模型評估

模型在RTX3060 圖像處理器上進行訓練，訓練時長約為145 min。經過100 次迭代后，訓練集和驗證集的準確度損失曲線均得到收斂（圖7）。其中，訓練集準確度達到99. 96%，訓練集損失值約為0. 0017；驗證集準確度達到99. 82%，驗證集損失值約為0. 0082。

本文使用歸一化的混淆矩陣（Confusion Matrix）進一步分析了模型性能。混淆矩陣展示了一個分類器對樣本進行分類時產生的錯誤，即真實類別被分為其他類別的錯誤。圖8 為五分類模型驗證集的歸一化混淆矩陣，矩陣主對角線為預測標簽與真標簽吻合的概率值。在驗證集中C1 與C3 樣本全部預測正確，個數分別為890 與896，所以對應的主對角線的值為1；C2 個數為個910 個，其中被正確預測的個數為908，占比為99. 8%; 898 個C4 樣本中有897 個預測正確，占比99. 9%；893 個C5 樣本中有887 個預測正確，占比99. 4%，其中有6 個C5 樣本被錯誤預測為C2 樣本。

2. 4 數據降維可視化

降維可視化是將高維數據轉換為二維或三維的可視化展示，有助于更好地理解和分析數據，從而找到數據中的模式和規律。

神經網絡可以把高維度復雜的數據輸出為線性可分的語義特征，最后只需要一層線性分類器即可把數據進行分類，這說明倒數第二層的語義特征在高維度中線性可分。本文使用 T?SNE（ t?Distributed Sto?chastic Neighbor Embedding）降維方法。T ? SNE 是一種用于高維數據降維可視化的算法，其主要思想是將高維數據映射到低維空間中，使數據點之間的相對距離在低維空間中保持不變或最大化。

本文將驗證集在倒數第二層（全局池化層）輸出的語義特征進行降維（圖9）。雖然個別事件被錯誤分類，但五分類數據集確實線性可分。

3 五分類算法檢測流程與應用

3. 1 五分類事件檢測流程

本文提出的五分類算法將輸入樣本預測為噪聲、完整的微震事件、只含P 波、只含S 波、多個微震事件五類。但在實際應用中，最終需要獲得同時包含P 波和S 波的微地震事件，而五分類算法中的C3 與C4樣本并不是完整的微震事件。因此，需要通過移動時窗再次處理C3 與C4 樣本，從而得到完整的微震事件用于震源定位。

首先，把野外采集的連續波形記錄（圖10a）按照長度為1001 的時窗（ 約為500 ms）進行切分（ 圖10b），不同于二分類模型，五分類樣本的時窗沒有重疊。得到規定長度時窗后，進行歸一化處理，輸入到訓練好的模型中進行分類；然后利用訓練好的ResNet 模型把時窗預測為C1、C2、C3、C4 和C5 五類（圖10c）。C2 類型可直接用于微地震定位。通過分析可知，若相鄰兩個樣本預測結果為C3 和C4，則大概率是一個微地震事件被截斷成兩部分，將C3 樣本向后時移或者C4 樣本向前時移即可。本文選取200 ms 進行時窗移動能夠獲得完整的微震事件。如果C3 或C4 類型單獨出現，同樣向后或者向前移動200 ms，然后用于微地震定位。此外，C5 樣本情況較為復雜，可以對樣本進行人工處理后再用于定位。

3. 2 五分類模型應用

隨機選取未用于模型訓練的原始波形記錄作為測試數據進行模型測試， 數據長度為17517 ms（ 圖11）。該段數據的水力壓裂事件較為密集，信噪比范圍較廣，可充分測試五分類模型性能。首先，把該測試數據按1001 數據點的時窗大小截取，可以得到35 個樣本。所有樣本經過歸一化處理后，輸入到ResNet 模型中進行預測。通過人工檢查，測試數據中一共包括23 個微地震事件。

把預測結果按照順序并把不同類別標簽的時窗用不同顏色標注出來（圖12a）。據統計，五分類模型預測結果為：9 個C2、5 個C3、5 個C4 以及3 個C5。其中，相鄰兩個樣本預測結果為C3 和C4 的情況有3次，為了防止反復拾取，此時只移動C4，故從C3 和C4 中一共得到7 個微地震事件。經過人工檢查，3 個C5 事件一共包含7 個微地震事件，故成功拾取的微地震事件個數為23。

圖12b 與圖12c 是圖12a 中的兩個紅色虛線框位置的放大圖，分別是預測結果為C2 與C5 的低信噪比事件。該結果表明，ResNet 模型有很強的抗噪能力。

然后，將訓練好的五分類模型應用到該套水力壓裂數據的第14 個壓裂段的微地震事件檢測中。第14個壓裂段的壓裂時長約為135 min，把該數據按1001數據點（約為500 ms）的時窗大小截取，得到16204 個樣本。經過統計，該壓裂段樣本的預測結果為：10889個C1、2830 個C2、955 個C3、824 個C4 以及706 個C5。預測時間為82 s。本文還利用不重疊時窗的二分類模型對該16204 個樣本進行預測，預測結果為：12440 個C1、3764 個C2。

為了評估微地震檢測的準確度，利用Zhang 等[42]提出的微地震到時拾取泛化模型對檢測出的事件進行P、S 波到時拾取。該泛化模型的輸入為12 道三分量波形數據，輸出為12 道的P 波、S 波和噪聲的概率分布。對于輸入時窗的任意一道，可以將兩個最大概率值所在橫坐標作為拾取的P 或S 波到時。Zhang等[42]規定這兩個拾取點之間的距離需要大于高斯分布半徑（20 ms），且拾取點的概率值大于0. 7。由于很多檢測的事件信噪比較低，拾取模型[42]的泛化能力不足夠完全正確拾取P 波和S 波到時。為了保證定位的準確性，需要剔除拾取結果較差的事件。剔除標準為：①拾取網絡只拾取到單震相；②拾取的P 波或S波走時個數少于3 個；③時窗內包含多個事件。

經過剔除拾取結果較差的事件后，最終定位得到的五分類和二分類微震事件數分別為2844 個和1789 個。文中利用網格搜索的方法進行微地震定位。震源搜索深度為2900～3500 m ，水平方向為0～500 m，搜索間隔為2 m。定位結果如圖13 所示。定位結果分布在壓裂井的附近，較為合理，表明兩種模型的檢測結果均具有較高的準確度。但五分類檢測方法檢測出了更多數量的微震事件，進而可更好地監測壓裂過程和評估壓裂效果。

4 結論與認識

本文提出了一種基于殘差網絡的微地震事件五分類檢測方法。通過有針對性的數據增廣方法，用較少的標注數據獲得了較高精度的訓練模型，模型訓練集精度達到99. 96%，驗證集精度達到99. 82%。分別采用五分類與二分類算法對野外監測的連續數據進行檢測，并對檢測結果進行到時拾取，用于微地震定位。測試結果表明，五分類檢測方法檢測出了更多數量的微震事件。

值得注意的是，針對信噪比較差的實際數據，如果人工難以給出標簽，可以先通過繞射疊加算法（不需要拾取P 波和S 波到時）進行部分震源定位[43]，然后將能定位出事件的時窗作為訓練樣本，并通過數據增廣方法擴充數據集，從而用于模型的訓練和后續壓裂數據的處理。此外，本文提出的五分類檢測算法的泛化性有待進一步測試，未來可針對不同壓裂監測數據特征進行重新訓練或者使用遷移學習算法以提高網絡模型的泛化能力。處理新的工區數據時，只需獲得較少的微地震事件，通過文中采用的數據增廣方法，便可進行模型訓練。