卷積神經網絡在遠-近地震震相拾取中的應用及模型解釋＊

申中寅吳慶舉

（中國北京 100081 中國地震局地球物理研究所）

引言

隨著地震觀測規模與密度的提高，大量涌現的地震數據對震相拾取算法的性能及可靠性提出了更高要求.基于地震信號對特定物理量的擾動，震相拾取包括震相識別與到時測量兩部分.這些物理量包括能量（如長短窗均值比）（Allen，1982）、峰度（Saragiotiset al，2004）、赤池信息準則（Sleeman，van Eck，1999）等指標.相較于傳統方法，卷積神經網絡（convolutional neural network，縮寫為CNN）能充分利用地震波形的全局特征，具有良好的魯棒性和泛化能力（ 趙明等，2019a，b；Zhuet al，2019）.根據 CNN 在走時測量的參與程度，相關算法依次為：① CNN識別震相，利用傳統方法測量到時（如CNN+dbshear，趙明等，2019a）；② CNN識別震相，利用得分函數測量到時（Zhuet al，2019）；③ 地震波形經CNN直接輸出震相類型及其到時，如 U-net （Ronnebergeret al，2015；趙明等，2019b）和多任務 CNN （李健等，2020）.從①到③，整個震相拾取流程趨近端對端模式（輸入原始波形，直接輸出震相到時），并伴隨模型和標簽復雜度的增加.與此同時，對CNN結構與訓練的探索（于子葉等，2020；周本偉等，2020）為方法優化及其內在原理闡釋開辟了道路.

另一方面，對CNN透明化的嘗試推進了相關解釋算法在查漏優化、可靠性評估及逆向學習等環節的應用（Selvarajuet al，2016）.在震相識別領域，包括CNN在內的機器學習方法，當前大致處于復現人工挑揀的水平，因此有必要深入了解模型的潛在缺陷及其震相判別機制.在眾多解釋算法中，平滑GradCAM++和類模型可視化（class model visualisation，縮寫為CMV）（Simonyanet al，2014）分別從具體個案和模型整體出發，定量評估了CNN的決策敏感區域以及特征復現能力.其中 GradCAM 系列（Selvarajuet al，2016；Chattopadhyayet al，2018；Omeizaet al，2019）整合了類激活映射（class activation mapping，縮寫為 CAM ）的類別敏感性（Zhouet al，2016）與導向背傳播的像素級分辨率（guided back-propagation，縮寫為 GBP）（Springenberget al，2014），而平滑 GradCAM++ （Omeizaet al，2019）又進一步優化了最末卷積層的梯度權重及其導向背傳播過程.CMV通過反向傳播損失函數合成各類別的高分輸入圖形，是評估模型特征提取能力的重要依據.在不觸及模型及數據結構的前提下，上述算法已成功用于CNN圖像處理模型的解釋.然而在震相識別領域，CNN解釋算法的應用與實踐尚處起步階段.為此，本研究采用CMV與平滑GradCAM++分析訓練得到的CNN模型，以考察其震相識別的可靠程度與內在機制.

鑒于此，本文擬采用CNN識別震相的得分函數測量到時來拾取北京國家地球觀象臺（以下簡稱北京臺）記錄的遠震（P，S）和近震（Pg，Sg）波形，以期考察輕量CNN（參數不超過1萬）在小樣本集（少于1萬）的表現.此外，模型解釋算法也被用于評估訓練所得CNN的性能及其決策機制.

通過有放回隨機抽樣構造訓練樣本：隨機抽取固定數目的樣本并確保各震相數量相同（每個震相取2 000個，樣本子集大小為5×2000＝10 000）；重復上述抽樣過程，生成子樣本（10次采樣生成10個子樣本）.子樣本集合作為整體構成一個完整的輪次參與訓練.統計結果表明，子樣本集合囊括了幾乎全部震相及大部分“噪聲”.較小規模（5×500＝2500）的測試樣本同樣由10個子集構成.上述的隨機采樣方法在充分利用現有數據的同時，還有助于訓練過程中局部最優點的遍歷.

2 模型搭建與訓練

2.1 CNN模型結構

鑒于不同震相持續時間的差異，為模型輸入長（40 s）?短（10 s）窗三分量記錄，以應對不同時間尺度的波形特征.模型的輸出向量依次由噪聲N和震相P，S，Pg和Sg的概率組成.隱藏層結構參考Zhu和Beroza （2019），隨卷積-池化依次倍縮（倍增）特征層的維度（數目）.由于激活函數選用ReLU，卷積層權重的隨機初始化采用he_normal分布，以避免梯度傳播的不穩定現象（Heet al，2015）.相關論述可參閱附件2.研究還進一步考察殘差塊添加對CNN性能的影響.殘差塊采用“瓶頸（bottle-neck）”結構以控制參數規模（Heet al，2016）.

CNN和瓶頸式殘差塊結構如圖2所示.圖2a展示了不同深度的卷積-池化序列.以4層CNN 為例，輸入張量（1000采樣點，6道）在雙邊鏡像延長至（1024，6）后，歷經 4次卷積-池化成為（8，64），最后展平并全連接至輸出層.其它深度的網絡結構依次類推.由于卷積層尺寸在深度為7時倍縮至1，模型深度搜索止于7層.圖2b展示了瓶頸式殘差塊的內部構造：數據流在依次卷積減-增特征層后恢復初始形狀，并與自身對位相加輸出結果.

圖2 卷積神經網絡（a）和瓶頸式殘差塊（b）結構示意圖紅色虛線對應于層深為4的CNN結構Fig.2 The structure of the convolutional neural network （a） and a bottle-neck residual block （b）The red dotted line corresponds to the CNN structure with a depth of 4 layers

2.2 模型訓練

模型訓練采用早停及L2正則化控制過擬合.其中早停用于動態判定訓練迭代的終止：以驗證集的準確率（損失函數）為準，當驗證集準確率（損失函數）連續3次低于（高于）最優值時，終止訓練并輸出最大準確率（最小損失）模型.研究主要考察最大準確率模型.L2正則化則懲罰模型復雜程度，最終損失函數為

式中：wm和M為模型參數及其總數；λ是L2因子的權重，可由數值實驗確定合理范圍.考慮到λ搜索和殘差塊添加的運算成本，模型最優結構的搜索將分步進行（結果見第3節）：

1） 從1至7依次增加卷積-池化層數目，考察網絡深度對震相識別能力的影響.考慮到模型初始化的隨機效應，訓練將重復10次進行.不考慮L2正則化（λ＝0）.

2） 在1）最優模型的基礎上，搜索合適的λ取值λ*，所得模型及λ*將用于后續測試.λ搜索由疏到密分兩步進行.

3） 在 2）最優模型的基礎上，分層依次添加殘差模塊（Heet al，2016）.當模型準確率連續三次小于同深度最優結果時，終止當前層位的殘差塊添加.

在步驟2）和3）的訓練中，在既有最優模型的基礎上進行了參數微調.這不僅節省訓練時間，還充分利用了步驟1）的隨機搜索結果，有助于模型參數的客觀比較.

2.3 性能指標

混淆矩陣（confusion matrix）記錄了 CNN 對測試集樣本不同震相的識別能力，是模型性能評估的重要手段.在圖3中，混淆矩陣的行和列序號分別對應樣本的真實標簽和預測結果：元素Cij即標簽為i的樣本被分類為j的數目.一般而言，模型性能越好，對角線元素越大.為全面客觀地描述模型表現，混淆矩陣派生出包括準確率、震相精度、召回率和F1得分等指標.

圖3 混淆矩陣Fig.3 Confusion matrix

準確率A反映了模型的總體震相識別水平，是混淆矩陣對角線元素之和占樣本總數的比例.精度P和召回率R描述了模型對特定震相的識別能力.其中震相精度反映了預測類別的正確比例.召回率則關注指定標簽的成功檢出率.模型準確率、震相i的精度、震相i的召回率表達式分別為：

F1得分被用于考察模型的整體性能，該準則傾向于為指標均衡的模型賦予高分.當精度與召回率同等權重時，震相i的F1得分定義為

而整體F1得分

則作為模型性能的最終衡量指標，衡量震相識別的總體表現.

3 訓練結果

卷積層深度.未清洗樣本訓練時不同模型深度N的準確率A、損失函數和F1t如圖4a，b和附件1的表1所示.受模型初始化等隨機因素影響，各輪次結果波動顯著.盡管如此，依舊可見以下規律：單卷積層顯著劣于其它結構；當卷積層深度大于等于2時，CNN準確率維持較高水平；中等深度（4，5層）模型表現較好.其中4層CNN準確率最高（A＝0.826，F1t＝82.4），而最低損失值則出現在5層網絡.單層網絡表現較差主要源于過簡結構的欠擬合，而深層（特別是7層）準確率低的結構成因有待論證.具體數值及統計結果可參見附件1的表1.

圖4 卷積層深度對模型性能的影響（a） 未清洗數據，第4輪訓練中不同卷積層深度（線上序號）模型的表現；（b） 未清洗數據，各輪次訓練中不同卷積層數模型的最高準確率和最低損失函數；（c） 已清洗數據，第10輪訓練中不同卷積層深度（線上序號）模型的表現；（d） 已清洗數據，各輪次訓練中不同卷積層數模型的最高準確率及最低損失值Fig.4 The influence of convolutional layer depth on model performance（a） Data unwashed，the model perfomance for different depths （marked by numbers） during the 4th training round；（b） Data unwashed，the maximum accuray and minimum loss with different depths in the total 10 training rounds；（c） Data washed，the model perfomance for different depths （marked by numbers） during the 10th training round；（d） Data washed，the maximum accuray and minimum loss with different depths in the total 10 training rounds

L2正則化.在最優深度的基礎上（N＝4，A＝0.826），測試了L2正則化系數（λ）的影響，采樣間隔按疏密依次進行（附件1中表2）.疏測試平行進行3次，比較λ的不同量級：0.001，0.1，1，10，100.當λ為10，100時，模型準確率分別為0.842和0.838，明顯高于未正則化的結果.為獲得λ的最優范圍，密測試從10到100進行（步長為10）.由于結果穩定性良好（最大偏差不超過0.005），密測試只進行1輪.結果表明，當λ＝30時模型準確率最高（0.860）.

殘差塊添加.在L2正則化最優模型（N＝4，λ＝30，A＝0.860）的基礎上逐層添加殘差塊.相比卷積層深度及L2正則化，殘差塊帶來的改善十分有限（模型準確率＋0.002至＋0.005不等）.因此，本文不再考察殘差塊對CNN的影響.

數據清洗及重新訓練.基于最優模型（N＝4，λ＝30，A＝0.860），逐一審查訓練集和測試集的錯判案例，并清洗實際波形與標簽對應不佳的樣本.被清洗的數據主要來自：噪聲疊加，信噪比過低，震相分析標注的偏差，噪聲自動截取時相鄰震相的混入等.由于剔除了不合適的震相和噪聲疊加樣本（特別是P與S），單個訓練集與驗證集分別收縮至5×800和5×175.

給定λ＝30，卷積層深度N取1，2，···，7訓練10次，結果見圖4c，d和附件1表1.得益于訓練樣本波形的改善以及驗證集質量的提高，數據清洗使模型準確率大幅提高.雖然最高準確率出現在5層CNN （A＝0.971，F1t＝97.1），模型性能隨網絡深度變化的趨勢依然不變：適中的深度（N＝4，5，6）普遍優于過淺、過深情形（N＝1，2，7）.與此同時，卷積層深度為5的CNN依然具有最低的損失值.

4 模型解釋

本節使用兩種模型解釋算法評估最優CNN模型（清洗樣本集，N＝5，λ＝30，A＝0.971）.其中類模型可視化CMV旨在反演各震相的最高得分波形，平滑GradCAM++ （Omeizaet al，2019）則用于勾勒輸入波形的決策敏感區域.這些梯度算法能在不觸及模型結構的前提下，定量評估目標震相的特征提取及其判別機理.

4.1 類模型可視化（CMV）

作為基于模型本身的解釋方法，CMV根據得分函數梯度修正空白輸入（I0＝def0），以獲取使震相c得分yc最大的特征波形.為減少類別混染，yc取“softmax”激活之前的數值.為控制幅值及波動水平，目標函數具有形式

圖5 類模型可視化算法的偽代碼（Nguass＝20）Fig.5 Pseudo code of Class Model Visualisation（Ngauss＝20）

完整的CMV反演結果可見附件3，圖7a選取了其中的長窗分量，可見如下特征：

1） 短窗振幅顯著小于長窗；

2） P和Pg垂直分量優勢明顯，Sg水平分量有所發育，S水平分量占優；

3） 噪聲的較大振幅出現在窗口前半部分；

4） 地震震相波形各異：P與Pg振幅突增明顯，存在明顯后續能量；S受前驅干擾影響顯著，后半窗波形模糊；Sg則以清晰的Pg前驅能量為特征.

其中2）——4）表明訓練所得CNN模型能捕捉到目標震相的基本特征，但1）暗示短窗數據未能有效參與預測.此外，震相S呈現出較低的波形質量，震相Sg波形明顯偏離窗口中心位置.這些都是CNN應用中必須注意的問題，也是未來需要改進的方向.

4.2 平滑GradCAM++

GradCAM++旨在評估波形不同部位對模型預測的影響，算法可分解為“Grad”和“CAM”兩部分.其中Grad通過誤差量的導向背傳播（guided back-propagation，縮寫為GBP）描繪圖像細節對決策函數的貢獻.而CAM則反映最末卷積層不同位置對分類結果的影響.CAM結果經上采樣與Grad對位相乘，得到輸入波形對模型決策的敏感度分布.相關公式可見Chattopadhyay等（2018），圖6為平滑 GradCAM++算法的偽代碼.Omeiza等（2019）的平滑GradCAM++藉由加噪平均進一步提高了模型解釋水平，也是本文采用的算法.

圖6 平滑 GradCAM++的偽代碼Fig.6 Pseudo code of smooth GradCAM++

為獲得震相敏感區域，平滑GradCAM++只選取高分案例（得分＞0.99）.相比二維圖像處理，震相識別更易出現接近1.0的得分.此時過小的損失函數會隨GBP消失為零，導致空白結果.噪聲疊加可有效防止過高得分引起的GBP消失，這也是選用平滑GradCAM++的重要理由.研究向輸入波形添加高斯噪聲（1/2最大振幅）并平均所得GBP （10次），并與CMV一同構成了CNN解釋的依據.

圖7b展示了震相波形（綠色）及其平滑GradCAM++結果（灰色），完整內容可附件3.與CMV類似，平滑GradCAM++的短窗振幅普遍偏低.長窗垂直分量在模型預測中起著決定性的作用.P與Pg對窗口中心兩側較為敏感：前側的平靜與后側的振動均對震相判別不可或缺.Sg也具有前側（Pg振動）和后側（Sg振動）兩個敏感區域，只是后者幅度相對較小.對于震相S，窗口后半段并未出現可見的敏感區域.結合CMV （4.1節）結果可以看出，模型在波形分量平衡及震相S識別等方面尚有較大改進空間.

圖7 模型解釋結果（a） 類模型可視化（CMV）；（b） 平滑 G radCAM++.綠色背景為原始波形Fig.7 Results of model interpretation（a） Class model visualization （CMV）； （b） Smooth GradCAM++.The green background stands for the original waveforms

5 連續波形下的模型部署

參照Zhu等 （2019）將2020年4月的連續波形用于震相拾取實測，以規避訓練集和測試集的影響.只有當震相（P，S，Pg和Sg）概率大于閾值時才認定為有效觸發，得分函數最大時刻即為震相測量到時.為平衡掃描效率和精度，研究采用變步長滑動窗口進行震相的搜索和測量，窗口長度以20 s 長窗為準.特征函數選用模型的輸出得分.具體流程如下：

1） 長步長（1 s）搜索震相，合并相鄰同震相窗口；

2） 短步長（0.1 s）計算目標震相的特征函數，最大值對應時刻即震相的測量到時.

在連續波形掃描之前（5.2節），有必要考察震相片段的局部掃描結果（5.1節），以便合理選定閾值.

5.1 震相片段的掃描

為確定不同震相閾值，首先對波形片段（基于測震目錄）進行短步長掃描，概率得分及走時偏差如圖8和附件1表3 （前4行）所示.不同震相表現為：震相P，具有大量較低得分，走時偏差集中在0 s附近；震相S，得分與走時偏差均顯著分散；震相Pg，得分較高（＞0.999），除過個別情況（走時偏差＞10 s，與Sg混淆），走時測量略有提前；震相Sg，得分較高（＞0.999），走時整體延遲.其中Pg-Sg走時偏差印證了CMV反演波形中震相初動相對窗口中心位置的偏離.為此，我們以0.999作為短步長測量的得分閾值，以平衡震相拾取的準確性與召回率.考慮到模型解釋中震相S的較差表現（圖7），震相S的得分閾值被調高至0.999 999以屏蔽干擾.為避免稀疏采樣引起的漏檢，長步長（1 s）搜索采用了更為寬松的閾值（P＞0.9，S＞0.99，Pg＞0.9，Sg＞0.9）.本研究設定的閾值（0.999）明顯高于 Zhu等 （2019）的0.5和Ross等（2016）的0.98.這可能來自樣本及訓練等環節，具體可見6.3節.

圖8 連續觀測震相片段的走時測量結果（a） 各震相走時偏差及模型得分的分布；（b） 各閾值下諸震相的識別數目（藍色方點）與走時偏差（紅色圓點）Fig.8 Traveltime mearurements for seismic phase clips cut from continuous observation（a） Travetime residuals and CNN-predicted phase scores；（b） Detection number （blue square） and traveltime residuals （red dot） of the phases

5.2 連續波形掃描

根據預設閾值掃描了2020年4月的連續波形.鑒于Pg-Sg總成對出現且到時間隔有限（小于 6 0 s），事先修剪了未成對（或間隔超過 6 0 s）的Pg，Sg拾取記錄.當震相拾取到時與目錄到時相差小于 5 s時才定為檢出（S 震相放寬至 1 0 s），結果見圖9和附件1表3 （最末行）.可以看出，震相間的干擾十分有限（僅有1例Sg誤識別為Pg）.識別偏差主要來自CNN的噪聲過敏感（噪聲被誤判為震相）與震相欠敏感（震相被誤判為噪聲）.正確識別結果的得分函數分布如圖10所示.為突出相對變化 趨勢，得分函數經過了單調變換x→ ? lg（1－x）、去均值和最大值歸一等處理.其中震相P僅表現為高頻成分的陡增，卻依然被準確拾取.而波形清晰的Pg，Sg卻呈現出明顯的走時偏差.對于錯誤結果，本節將從召回率和準確率兩方面分別討論.

圖9 連續波形掃描結果的混淆矩陣Fig.9 The confusion matrix for the continuous waveform scanning

圖10 正確的震相識別水平虛線由低到高依次對應長、短步長的得分閾值，原始波形經歸一化處理，下同Fig.10 Seismic phases detected correctlyHorizontal dotted lines correspond to thresholds for long（lower） and short （upper） scanning steps，the waveform is normalized，the same below

震相召回率反映了CNN的靈敏程度.震相未能召回的原因主要有：震相概率得分過低，走時測量偏差過大.不同震相的召回率分別為：P （66.7%），S （10%），Pg （40%），Sg （60%）.對于震相P，較低的信噪比壓低了概率得分，是漏檢的主要原因（大于20例，圖11b）.而波形清晰的漏檢案例常具有豐富的低頻成分，并與較高的得分閾值相關（圖11a）.除了較高的閾值，震相S的召回率也受波形質量與模型性能制約（圖11c?d）.噪聲擾動（導致過低得分）及走時拾取的系統偏差是Pg-Sg漏檢的重要原因（圖11e?f）.

震相精度反映了CNN的抗噪能力.不同震相精度分別為：P （29.3%），S （41.7%），Pg （6.8%），Sg （10.2%）.較低的震相精度表明大量噪聲波形被誤判為地震震相.對于震相P，除過個別波形失真、疑似地震（圖11g）、震相誤讀（圖11i）外，誤判案例主要來自周期為4——5 s的背景噪聲的擾動（總數不小于150，圖11j）.此外，得分函數振蕩造成的搜索窗口分裂也是導致誤判的原因之一（圖11h）.震相S較高的得分閾值在造成大量漏檢的同時，也有效降低了噪聲混入的概率.除過個別震相誤判（圖11k），震相S的假正例主要來自模型對背景噪聲的過敏感（圖11l）.后者也是Pg-Sg假正例的主要來源（圖11n）.與此同時，模型亦成功拾取了未錄入測震目錄的Pg-Sg事件（圖11m）.

圖11 錯誤震相的識別（e） Pg 走時測量偏差過大；（f） Pg-Sg 波形不清晰；（g） 疑似地震事件；（h） 搜索窗口分裂；（i） 后續震相（ScP）干擾；（j） 常見 P震相誤判波形；（k） 誤判為 S 的 S cS 震相；（l） S 震相的常見誤判波形；（m） 未錄入地震目錄的 P g-Sg；（n） Pg-Sg 常 見誤判波形Fig.11 Cases of the wrong detections（e） Too enormous Pg travetime residual；（f） Unclear Pg-Sg waveform；（g） Suspisious earthquake；（h） Splitting in searching window；（i） Later-coming phase （ScP）；（j） Common false P detection；（k） ScS identified as S；（l） Common false S detection； （m） Pg-Sg event not in catalogue；（n） Common false Pg-Sg detection

圖11 錯誤震相的識別（a） P 震 相波形清晰；（b） P 震 相波形不清晰；（c） S 震 相波形清晰；（d） S 震 相波形不清晰Fig.11 Cases of the wrong detections（a） Clear P waveform；（b） Unclear P waveform；（c） Clear S waveform；（d） Unclear S waveform

6 討論

上文介紹了樣本制備、訓練以及模型解釋、部署的整個流程，此節將重點說明其中涉及的一些問題.

6.1 訓練樣本規模

雖然模型的參數量被控制在較低的量級，其復雜度依然高于訓練樣本數目.為此我們對比了震相樣本規模從200到800時各深度模型的訓練結果，訓練重復5次.

圖12a展示了CNN準確率隨訓練樣本規模的變化趨勢.4層CNN的準確率普遍較高，5層次之.隨著樣本規模的增加，深層（6，7）模型的表現逐漸優于淺層模型（2，3）.圖12b對應樣本實際使用量與訓練樣本規模之間的關系.受樣本池規模所限，樣本使用量的增速隨樣本規模的增加持續減小.這雖然會影響模型性能的提高，但無法解釋圖12a呈現的復雜趨勢.

圖12 樣本規模對模型準確率及樣本利用情況的影響（a）模型準確率與訓練樣本規模的關系；（b）訓練樣本規模與樣本使用量的關系Fig.12 Relationship among train data size，model accuracy，and samples used（a） Relationship between train data size and model accuracy （b） Relationship between train data size and the number of samples used

圖13列舉了模型準確率隨深度的變化，是圖12a的詳細補充.當訓練樣本規模小于500時，模型準確率與樣本規模正向相關.而當樣本規模超過500時，模型準確率增速明顯放緩（3，4，6層）乃至下降（2，5，7層）.除淺層（2，3層）外，模型最高準確率均未對應最大樣本規模.此外，即使樣本規模僅為200，4層CNN的準確率依然超過了0.96，反映了輕量CNN對小規模訓練樣本的較強適應能力（Zhuet al，2019）.

圖13 各深度模型準確率與訓練樣本規模的關系Fig.13 Relationship between the raining data size and model accuracy

由此可見，擴充樣本有助于模型準確率的提高，當前樣本規模尚足以滿足模型訓練需要.同時，訓練樣本規模亦無法解釋深層（6，7層）CNN準確率的下降.

6.2 正則化系數

在CNN的深度探索階段，研究采用了固定的正則化系數λ.為考察λ取值對結果的影響，本節將系統考察最優λ隨CNN深度的分布趨勢.

圖14a展示了模型準確率隨λ的變化趨勢，訓練進行3次并選取最高準確率結果.隨深度增加，最優λ也相應增大.即便考慮高λ值對深層CNN準確率的提高，5層模型依舊具有最好的表現.圖14b展示了λ與（模型系數的平方均值）的關系.由于正則項約束的增強，模型均方隨λ增加持續下降.而對于給定λ，模型系數的平方均值則與層數反向相關，體現了模型總體激活程度的減弱.最小二乘擬合結果揭示了最優λ與模型系數均方滿足

表明不同深度的最優模型具有大致相同的正則項（圖14c），其內在機制有待進一步澄清.

最優λ隨模型深度的分布可解釋兩次深度搜索結果的不同（3.1節與3.4節）.根據圖14a，當λ＜1時4層CNN具有明顯優勢，對應了第一輪深度搜索的結果（無正則項，相當于λ＝0）.而當進行第二輪搜索時，λ＝30 （第一輪4層CNN的λ搜索結果）使5層CNN的性能大幅提升，并使后者成為最優搜索模型.據此可知，當λ＝80時（當前取值30），5層網絡的性能還有望進一步提升.

圖14 正則化對各深度模型的影響（a） 正則化系數λ、模型層數與模型準確率之間的關系，圓點對應最高模型準確率；（b） 正則化系數λλ與模型均方之間的關系，其中空心圖形對應最高準確率模型（注意2層數據點疊覆于3層之下）；（c） 模型損失函數的構成Fig.14 The effects of regularization on CNN with different depths（a） Relationship among regularization factor λ，CNN depth，and accuracy；（b） Relationship between regularization factor λ and the squaremean of model weights，with hollow patterns standing for CNN with highest accuracy （2-layer dot is overlapped by that of 3-layer）；（c） Variation of loss function with CNN depth

6.3 連續波形掃描的過高閾值

本文的震相得分閾值明顯高于同類研究的結果（Rosset al，2016；Zhuet al，2019），這雖然出于屏蔽噪聲干擾的需要，但也與震相本身的過高得分密不可分.對于后者，本節將從數據泄露和模型過擬合兩方面進行討論.

數據泄露會造成虛高的訓練精度，通常來自模型對個別特征過度依賴以及驗證集同訓練集的混染.前者常見于聚類分析的不當參數選取，表現為個別指標完全控制模型分類.作為基于波形的分類算法，CNN對震相特定結構的過度敏感會導致CMV和平滑GradCAM++結果的異常，這已為模型解釋（第4節）所排除.此外，訓練集與驗證集截取時段的嚴格分隔也確保了二者的有效分離.然而，調參與訓練階段對樣本的共享作為潛在的數據泄露機制，有待于樣本的擴充與細化.

過擬合也會導致模型虛高的準確率，具體表現為測試集的準確率相對訓練集顯著偏低.訓練樣本欠缺和模型過高的復雜程度是過擬合的重要誘因.研究采用L2正則化、最大池化，以及早停降低過擬合風險.本節將從模型選擇、訓練樣本數量、L2正則化等三方面進行討論.

受隨機性影響，模型選取會以多次訓練的最高準確率結果為準.由于CNN訓練往往收斂于局部最優點而全局最優點常為過擬合點，有必要考察最高準確率模型（更接近全局最優點）的過擬合情況.圖15a對比了最高準確率模型相對其它模型（5層第4輪訓練）的震相得分，可見前者的P，Sg得分反而偏低，暗示模型收斂點對過擬合的影響十分有限.

訓練樣本數量也是影響模型過擬合的重要因素.為方便對比，選取準確率相近的4層模型（A=0.962±0.03）.圖15b展示了不同訓練樣本規模下的震相得分，未見震相得分的系統變化.

L2正則化通過平滑CNN參數來控制模型的復雜程度，其最優系數隨模型深度增加（6.2節）.連續波形掃描所用5層模型在訓練過程中用到的正則系數（30）顯著小于該深度的最優取值（80）.圖15c對比了不同正則化系數訓練下5層模型的震相得分（模型準確率均為0.971），可見L2正則化對震相得分的抑制作用.

圖15 不同模型的震相得分對比（a） 模型選擇；（b） 樣本規模；（c） 正則化系數.黑色實線對應y ＝ x ，震相得分經單調變換 x → ? lg（1－x）以便于展示Fig.15 Phase scores for the selected models（a） Model selection；（b） Training sample size；（c） L2 regularization factor.The black solid line corresponds to y＝ x ，with phase scores transformed monotonously by x → ? lg（1－x） for better view

綜合上述，正則化系數的不當選取是導致震相得分虛高的重要原因，不排除數據泄露的可能影響.盡管研究涉及的樣本規模不足以造成過擬合程度的系統變化，高質量訓練樣本的積累依然是CNN訓練質量的前提和保證.

7 結論

本文探討了卷積層深度、正則化、殘差塊添加和數據清洗等因素對CNN震相（P，S，Pg，Sg）拾取能力的影響.結果表明，適中的卷積層深往往具有更優的性能.深度為4，5的CNN 網絡分別具有（8，64），（4，64）的最末卷積層（圖2a）.在提取足量特征（64）的同時，上述層位依舊具有一定的局部分辨能力，從而能有效地處理震相序列沿時間軸的分布.與此同時，合適的正則化系數及數據清洗也能有效提高震相識別的準確率.相比而言，殘差塊的有限作用可能與地震數據（相對二維圖像識別）較低的復雜度（于子葉等，2020）有關.

為了客觀評估模型的震相識別機制，采用CMV和平滑GradCAM++解釋所得CNN.其中CMV反演波形基本復現了各震相的主要特征，包括不同震相垂直與水平分量的相對大小.而平滑GradCAM++則揭示了前方平靜段（P與Pg）或前驅能量（Sg）對震相識別的重要作用.與此同時，模型解釋也揭示了模型及訓練的問題.對于CMV，S震相波形噪聲過大，Sg波形初動偏離過大；對于平滑GradCAM++，S震相的波形區段對得分幾乎不作任何貢獻.而短窗分量的過低參與度導致長窗垂直分量主導CNN預測，不利于多分量地震圖的綜合利用.

最后，將CNN模型用于連續波形的掃描.掃描按長、短步長依次識別、拾取震相.結果表明模型具有初步的遠-近震相拾取能力，但靈敏度和抗噪性亟需提升.在討論章節，我們初步分析了訓練樣本規模、正則化參數、數據泄露和過擬合等因素對所得結果的影響，為后續深入研究提供了一定參考.

切實提高訓練樣本的質量與數量，改善S震相訓練水平并平衡各分量的模型參與度，是未來深入研究的重要方向.為此，后續研究可從以下幾個方面著手：

1） 適當延長震相S的選取時段，或從理論地震圖生成訓練樣本，從源頭緩解震相的不足.

2） 在模型架構上，輸入的長短窗分量流入各自的卷積?池化序列，以壓制長窗垂直分量的優勢地位（附件1中圖1a）.

3） 采用軟標簽編碼訓練樣本，適當平移震相波形，在豐富樣本的同時優化得分函數對走時的敏感（附件1中圖1b）.考慮到短窗長度，震相得分的非零區間半徑取5秒.對于Pg與Sg的中間時段，該值參考Pg-Sg的到時差值（中點處Pg，Sg標簽均為0.5）.

4） 隨機抽樣生成的訓練樣本集合單獨參與CNN的訓練，得到各自對應的子模型，最終結果取決于每個子模型得分的“投票”結果.這種集成學習策略有望利用隨機性消除局部極值點的影響，壓制過擬合造成的“極端”預測結果.

本文模型訓練采用的操作系統為Fedora 30，中央處理器為Intel Core i7-6700CPU@3.40 GHz×8，硬盤為 ST1000LM044 HN-M101SAD （1 TB，7200 轉），1 個 epoch （以清洗后樣本集，深度為5的CNN為例）耗時約10分鐘.模型的搭建與訓練在python3.7下的tensorflow2.2環境完成，震相數據處理采用obspy包，圖件繪制采用matplotlib和cartopy包.中國地震局地球物理研究所李麗研究員為本文提出了寶貴建議.北京國家觀象臺朱戰斌高級工程師與周江林高級工程師為本文提供了數據支持，審稿專家為本文提出的寶貴修改意見，作者在此一并表感謝.