震源參數反演及精度評定的Bootstrap方法

王樂洋， 李志強

1 東華理工大學測繪工程學院， 南昌 330013 2 山東建筑大學測繪地理信息學院， 濟南 250101

0 引言

地震是地殼巖層受力突破彈性極限(抗剪切強度)后快速破裂錯動引起的地表振動或破壞，是對人類文明威脅最大的自然災害之一.震源參數反演是研究地震斷層錯動特征、震源機制、觸發關系等的重要方式，也是地震危險性評價及地震預警的基礎.地表形變是地震發生時最直接的表現，是研究地震斷層構造、破壞程度及地震機理等的重要資料(王樂洋等，2017).在獲取地表形變數據的過程中，由于各種噪聲的干擾導致形變數據的準確性和可靠性受到影響，進而導致震源參數的精度損失.因此，獲取震源參數的精度信息、誤差水平及可信程度顯得尤為重要，對震源參數進行精度評定研究成為了提升成果質量的重要內容，對準確、深刻的理解地震也具有一定的實際意義.

由于依據彈性半空間矩形位錯模型(Okada，1985)建立的關于地表形變與震源參數之間的關系為復雜多維的非線性函數，因此需要依靠搜索算法進行非線性反演.常用的搜索算法有模擬退火法(Kirkpatrick et al.，1983)、粒子群算法(Kennedy and Eberhart，1995；馮萬鵬和李振洪，2010)、神經網絡法(Barradas and Vieira，2000)、遺傳算法(Genetic Algorithm，GA)(Holland，1975)等，該類方法無需獲取參數的解析表達式，也無需求導計算，因此常用于非線性反演獲取參數最優值.遺傳算法不必明確地描述優化問題的全部特征，且不存在連續性、可導或單峰等特殊要求，因此具有通用性.此外，該搜索算法全局性性強，不易陷入局部最優，并且適合并行計算(王新洲，2002)，因此較適用于數據量大的大地測量震源參數非線性反演問題研究(Nolte and Frazer，1994；萬永革等，1997；王文萍和王慶良，1999).

將震源參數反演過程視為一個非線性函數，目前已有的非線性精度評定方法有近似函數法(徐培亮，1986；Xue et al.，2015；Wang and Zhao，2019)和近似非線性函數的概率密度分布方法(王樂洋和趙英文，2019；Wang and Yu，2019；Wang and Zou，2019).近似函數法需要利用泰勒級數公式對非線性模型展開，根據誤差傳播定律獲取參數估值的方差或協方差.但由于該方法無法避免求導計算，因此在震源參數精度評定研究中適用性較差.通過采樣來代替復雜求導計算的非線性精度評定方法包括蒙特卡羅(Monte Carlo, MC)方法、自適應蒙特卡羅(Adaptive Monte Carlo, AMC)方法、Jackknife方法、Sterling插值方法等.該類方法一般通過重復采樣獲取隨機樣本來近似非線性函數的概率密度分布，進而輸出待統計量的精度信息，這些采樣方法在一般的非線性模型精度評定問題中能夠取得較好的效果.在震源參數反演精度評定方面，Monte Carlo方法得到了較為廣泛的應用(Wright et al.，1999；Funning et al.，2005；Dawson and Tregoning，2007；Xu et al.，2017)，亦有學者將Jackknife方法、SUT(Scaled Unscented Transformation)方法、Sterling插值方法用于地震大地測量反演精度評估(Ichinose et al.，2003；Prieto et al.，2007；Wang and Zou，2021；Wang and Ding，2020).

Bootstrap方法(Johnson，2001；Chen，2018)也稱為自助法，是一種用于估計和修正統計估計值的偏差或方差信息的一種統計推斷方法.Efron(1979)首次提出了通過重復采樣得到的自助世界的經驗分布來逼近總體真實分布的思想，并給出了其基本采樣策略和相關理論證明.由于Bootstrap方法通過重采樣代替求導計算，只需選擇采樣次數并結合有放回抽樣策略獲取自助樣本，最后解算每個樣本便可求得未知統計量的均值和協方差陣，因此歷年來有較多學者將該方法用于偏差分析和方差計算(Sahinler and Topuz，2007；Kitagawa and Konishi，2010；O′Hagan et al.，2019).目前，關于Bootstrap方法的研究成果較為豐富，主要集中于該方法在統計學上的性質(Martínez-Camblor and Corral，2012；Pan and Politis，2016).目前，已有部分學者將Bootstrap方法用于地球科學反演精度評估問題，Wéber(2005)利用貝葉斯方法和Bootstrap重采樣對震源參數進行了誤差分析；Cesca等(2014)通過Bootstrap采樣來評估發震斷層位置的不確定性；Jia等(2017)采用Bootstrap方法來估計質心深度的不確定性.但Bootstrap采樣方法在大地測量震源參數反演精度評定方面的性能未得到較為全面的研究，尚未給出將該方法用于獲取震源參數精度信息的詳細采樣步驟及具體實現流程，并且在震源參數中誤差計算及區間估計方面的對比評價也未得到重視.

針對以上問題，本文從為了獲取參數估值更為合理的精度信息出發，將Bootstrap方法用于震源參數非線性反演的精度評定研究中，通過對地表形變數據實施自助重采樣并結合遺傳算法，進而設計得到震源參數反演精度評定的Bootstrap方法計算流程.最后通過6個模擬地震、Amatrice地震和Visso地震，對本文精度評定方法的有效性和可靠性進行驗證.

1 震源參數反演與遺傳算法

1.1 震源參數反演

均勻彈性半空間矩形位錯模型是反演震源參數的常用有效模型(Okada，1985)，本文采用Okada矩形位錯模型并結合GPS臺站地表形變觀測資料進行震源參數反演.地表形變量觀測值與震源參數之間為高度非線性關系，其函數模型可表示為

d=G(ξ)+e,

(1)

式中，d為GPS觀測到的地表形變量，它包含dx、dy及dz三個分量，分別表示東西、南北和垂直三個不同方向的地表形變.G(·)為利用Okada彈性位錯理論模型建立的關于GPS地表形變量與震源參數的復雜非線性關系.e為GPS地表形變量觀測誤差.ξ為震源參數，包括發震斷層位置(經度、緯度、深度)，斷層尺寸(長度、寬度)，斷層幾何(走向角、傾角、滑動角)以及滑動量；本文以頂深、底深、走向角、傾角、長度、滑動角、滑動量以及斷層中心位置為參數進行反演和精度評定，其中通過頂深和底深可計算獲得深度，長度可聯立頂深、底深及傾角獲取.

1.2 遺傳算法在震源參數反演中的應用

由于震源參數未知向量與地表形變觀測值向量之間是復雜的多維非線性關系，本文將遺傳算法運用于震源參數反演問題.遺傳算法最早由美國教授Holland(1975)提出，該算法是一種啟發式高級全局搜索算法，它所依據的是自然界生物優勝劣汰的生存法則，是基于達爾文的進化論和孟德爾的遺傳學說建立的.遺傳算法僅通過目標函數而不需要其他附加條件對種群(將待求參數稱之為種群，每個參數稱之為個體)進行全局搜索.對提高算法的搜索效率并保證反演得到的震源參數具有物理意義，反演前通常需要給定震源參數的搜索區間作為約束.本文在使用遺傳算法進行震源參數反演時，首先在全參數空間內進行全局搜索，當獲取初步的反演結果后再縮減搜索空間并通過反復迭代計算最后求得穩定的最優解.

將遺傳算法用于震源參數反演的大致步驟是首先使用二進制對區間內的數進行編碼；其次，隨機生成一組任意排列的字符串，即產生初始種群；然后根據提前設置的適應度函數判斷個體的優劣性，淘汰適應度低的個體，選擇復制適應度較高的個體作為父代并且進行交叉和變異操作進而生成下一代個體；最后通過不斷繁殖、選擇、交叉及變異，便可朝最優解的方向進化獲得最優的震源參數.遺傳算法實現非線性反演的具體過程總結為

1)選擇初始群體大小pop及最大迭代次數maxgen.

2)編碼：編碼是將震源參數解的搜索范圍轉換到遺傳算法所能認知的搜索空間的一種轉換方式，通過編碼將具體的非線性反演問題與自然界生物進化建立聯系.

3)生成初始種群：在編碼區間內隨機生成數量為pop的初始參數，作為遺傳算法求解優化問題時的初始解的集合.

4)個體適應度函數值的計算：適應度函數是用來評價個體優劣的指標.震源參數反演中常將殘差的均方根即式(2)作為目標函數，計算每個個體的適應度函數值，區分出每個個體的優劣性.

5)選擇復制：將步驟4)所計算出的個體適應度函數值作為評價標準，選擇出適應度較好的個體作為父代來繁衍后代種群.

6)交叉：交叉是個體兩兩配對的過程即基因重組過程，具體實現過程為從父代中選取部分基因，并以某個交叉概率交換部分染色體進而生成新的個體.該過程的目的是為了生成新一代個體以提升算法搜索能力.

7)變異：基因遺傳算法中的變異過程是指等位基因替換引起的基因結構的變化，是自然界生物進化的必然要求.變異的目的是保證種群多樣性，防止陷入局部最優.

8)更新：將步驟7)獲取到的子代作為更適應環境的父代種群，重復步驟4)～步驟7)，重新計算個體適應度函數值、選擇、交叉、變異，直到迭代次數大于設定的最大迭代次數即Itei>maxgen則跳出迭代.那么，最新一代中適應度函數值最小所對應的個體即為震源參數的最優參數估值.

在使用基因遺傳算法反演震源參數的過程中，本文采用地表形變觀測值與Okada模型預測值差值的均方根(Root Mean Square, RMS)誤差作為目標函數來評估反演結果，其表達式為

(2)

式中，n為地表GPS觀測點個數，Pi表示每個觀測值的權值，di為形變實際觀測值，ci為通過反演模型得到的地表形變預測值.

2 震源參數反演及精度評定的Bootstrap方法

Bootstrap方法將原始觀測值樣本作為總體，利用有放回隨機抽樣法從總體分布函數中得到一系列獨立樣本(稱為自助樣本)，并通過計算每個自助樣本來獲取未知統計量抽樣分布的經驗估計(Efron，1979；Sahinler and Topuz，2007).該方法的優勢性主要表現在它不需要對未知模型及分布做任何假設，也無需推導估計量的精確表達式，只需通過檢驗樣本內統計量的變化便能夠估計未知參數的均值和方差.Bootstrap方法通過重采樣獲取自助樣本的方式來代替復雜的求導運算，因此本文將Bootstrap方法引入到震源參數精度評定問題，試圖通過在避免求導計算的情形下依然能夠獲取合理的精度評定結果.

2.1 精度評定的Bootstrap方法

(3)

震源參數估值的方差-協方差陣計算公式為

(5)

=(1-α),

(6)

式中，pr表示概率，α為顯著性水平，Zα/2為對應的標準分數.

2.2 精度評定計算步驟

為充分利用總體包含在樣本中的先驗信息和數據性質，Bootstrap方法通過有放回隨機抽樣將GPS觀測點三個方向形變量數據采樣成多個自助樣本，原始樣本集中的每個樣本點被采樣到的概率相同，每個自助樣本的樣本容量與原始樣本相同.顧及GPS站點形變量含有誤差，故對三個方向形變量數據按照其對應的觀測精度進行加權處理，通常將三個方向形變量所對應的觀測中誤差作為定權的依據.在重采樣獲得自助樣本的過程中，需保持自助樣本中的形變觀測值與其權值Pi相對應，具體表現為重采樣三個方向形變量數據的中誤差.形變值與對應中誤差兩者的采樣方式一致.根據以上分析，總結得到震源參數精度評定的Bootstrap方法的完整計算步驟：

1)假設n個GPS觀測點的三個方向形變量為隨機項樣本N=(ux,uy,uz)，其對應的觀測中誤差為ε=(σx,σy,σz)；

2)將1/n的概率分別設置在N中的每個樣本值上；

3)產生服從均勻分布的n個隨機數(i1,i2,…in)，其中1≤ij≤n,j=1,2,3,…,n；

根據以上計算步驟，可得Bootstrap方法評定震源參數精度的計算流程圖，如圖1所示.

圖1 震源參數精度評定的Bootstrap方法迭代流程圖Fig.1 The iteration flow chart of the Bootstrap resampling method for precision estimation of fault source parameters

由以上的震源參數精度評定的Bootstrap方法計算流程可知，本文精度評定方法獲得的所有自助樣本中的數據均來源于原始固定的GPS觀測點三個方向形變量數據，并未根據更多的觀測形變信息進行計算.Bootstrap方法的采樣過程不會產生額外的模型誤差，也不會改變模型態性.盡管每個自助樣本的樣本容量與原始樣本相同，但不是改變GPS觀測點三個方向形變量數據的先后排列順序，有放回隨機抽樣過程使得每個自助樣本中可能存在重復的數據點，而另外一些樣本點則沒有出現.因此，每個自助樣本將隨機地異于原始樣本，導致每個自助樣本反演得到的震源參數估值存在細微差異.

3 實驗與分析

為驗證Bootstrap方法在震源參數非線性反演精度評定中的可靠性和有效性以及探討該方法在實際震例中的性能，本文共設計了6個模擬實驗，同時也將本文精度評定方法用于獲取Amatrice地震和Visso地震的震源參數精度信息，并與Jackknife方法和Monte Carlo方法進行對比分析.當使用到遺傳算法時，為同時保證較高的反演精度和反演效率，初始群體大小取值為450，最大迭代次數選擇200次，交叉概率為0.85，變異概率選擇0.1.Bootstrap采樣次數的選擇取決于統計量的具體研究問題并且依賴于需要做的檢驗(例如方差或置信區間)(Efron，1981)，諸多統計學家給出了抽樣次數的經驗取值，若要獲取未知參數較為合理的方差或中誤差，采樣次數應為50～200之間，當采樣次數超過200時，方差近似值的改進很小(Efron，1981；Efron and Tibshirani，1986；Léger et al.，1992).經實驗，通過增加自助采樣次數在一定程度上能夠提升精度評定結果質量，但隨著采樣次數的增加，中誤差結果趨于較為穩定的值，而其計算量將成倍遞增.本文依據已有研究的經驗取值并顧及復雜多維的地震非線性反演實際問題，在保證較高的反演精度條件下，同時考慮較高的反演效率，因此本文將Bootstrap方法的采樣次數設置為300.

圖2 仿真設置1的斷層平面位置及觀測點位分布Fig.2 The fault plane position and observation point positions of the first simulation setting

3.1 模擬地震

通過仿真獲取某一地震區域GPS觀測點三個方向的形變量資料(仿真設置1)：斷層面的幾何中心位于坐標系原點，37個GPS觀測點分布在斷層面幾何中心的40 km范圍內，模擬斷層平面位置及GPS觀測點分布如圖2所示(其中橫坐標E為東方向，縱坐標N為北方向，矩形框表示斷層投影平面)；斷層面幾何參數的仿真設置為：斷層頂深為2 km、底深為18 km、走向角為45°、傾角為60°、斷層長度為40 km，滑動角為70°、滑動量為1.5 m；將GPS觀測點坐標及參數模擬真值代入Okada模型正演獲取形變量真值；最后在形變值真值上施加服從均值為0、單位權中誤差為3 mm的隨機正態分布噪聲進而獲取模擬形變量.仿真設置1的GPS同震形變場在水平方向和垂直方向的模擬形變值如圖3所示.

在獲取GPS形變量數據后，分別采用Jackknife方法和Bootstrap方法對地表形變值進行采樣并獲取震源參數均值和中誤差；為推斷總體均值，利用Bootstrap方法計算震源參數的置信區間，并與Jackknife方法的結果進行對照.各方法獲取仿真設置1的震源參數結果、置信上下限及區間間隔列于表1.

由表1中的震源參數估值可以看出，使用遺傳算法搜索到的震源參數估值與參數模擬真值基本一致，表明本文采用遺傳算法進行非線性反演是正確且可靠的.Jackknife方法和Bootstrap方法得到的震源參數在符號和量級上與遺傳算法的結果相同，且Bootstrap方法的結果與真值的差值二范數小于遺傳搜索算法和Jackknife方法的結果；根據遺傳算法、Jackknife方法和Bootstrap方法獲取的震源參數得到的矩震級均為MW6.97，與模擬的矩震級MW6.98的差異僅在小數點后第2位.這表明基于重采樣的Bootstrap方法能夠適用于震源參數非線性反演，也說明該方法通過對形變量數據實施自助采樣在一定程度上能夠降低反演估值的偏差，因此其震源參數反演值比遺傳算法和Jackknife方法的結果更加靠近模擬真值.

圖3 仿真設置1的GPS同震形變場(a) 水平位移場； (b) 為垂直位移場.Fig.3 The coseismic displacement fields of the first simulation setting obtained by GPS(a) The horizontal displacement field； (b) The vertical displacement field.

表1 仿真設置1的震源參數搜索區間、參數估值及置信水平為95%情況下的置信區間計算結果Table 1 Search ranges, inversion values, and confidence intervals with 95% confidence level of source parameters in the first simulation setting

根據表1中各參數的置信上限、置信下限及區間間隔可以看出，在95%置信水平下，Jackknife方法獲取的震源參數置信區間均比Bootstrap方法獲取的置信區間更寬，說明Jackknife方法反演得到的震源參數的可靠性較弱，因此Bootstrap方法的反演精度較高.這表明在震源參數非線性反演問題中，Bootstrap方法獲取的參數估值比Jackknife方法的結果更具可信度.由于置信區間的計算公式取決于所用的統計量，Jackknife方法的樣本量為總觀測個數，而Bootstrap方法的樣本量可根據模型復雜度或樣本容量進行調整.因此，與Jackknife方法相比，Bootstrap方法在震源參數反演中更具可靠性和適用性.

圖4 仿真設置6的GPS同震形變場(a) 水平位移場； (b) 垂直位移場.Fig.4 The coseismic displacement fields of the sixth simulation setting obtained by GPS(a) The horizontal displacement field； (b) The vertical displacement field.

表2 六種實驗方案及其說明Table 2 The six experimental schemes and their explanations

在獲取以上6個仿真設置的形變量數據后，分別利用Bootstrap方法、Jackknife方法和Monte Carlo方法獲取震源參數的中誤差，3種精度評定方法的計算結果列于表3，并將6個不同仿真設置的中誤差可視化繪于圖5((a—f)所繪中誤差分別表示仿真設置1～6的結果).

通過表3中的參數中誤差結果可以看出，Bootstrap方法在6個不同的仿真設置下所獲得的震源參數中誤差與Monte Carlo方法的結果在符號和量級上相同、在數值上相近，表明基于重采樣的Bootstrap方法能夠對反演得到震源參數進行有效的精度評定.從計算結果上看，相比于Monte Carlo方法所獲得的震源參數中誤差，Jackknife方法得到的結果偏離較大，是由于Jackknife法易受原始樣本容量大小的影響且其重采樣過程獲取的刀切樣本較少所致.表明Jackknife方法的結果可能會低估震源參數反演估值的精度，因此無法反映參數估值的真實精度.Bootstrap方法獲取的震源參數中誤差小于Jackknife法的結果且與Monte Carlo方法的結果較為接近.從圖5也可直觀看出，在不同仿真設置情況下，Jackknife方法獲取的參數中誤差總是表現出較多深色，而Bootstrap方法的結果與Monte Carlo方法所得中誤差顏色較為相近且淺色居多，這表明將Bootstrap方法用于解決震源參數精度評定問題是可行的，也表明Bootstrap方法在獲取震源參數精度信息方面較Jackknife方法更加精確有效.

表3 六種不同仿真設置下震源參數的中誤差計算結果Table 3 The standard errors of source parameters under the six different simulation settings

圖5 不同仿真設置情況下的震源參數中誤差圖圖中，“TD”表示頂深(Top Depth)、“BD”表示底深(Bottom Depth)、“Strike”表示走向角(Strike)、“Dip”表示傾角(Dip)、“Length”表示長度(Length)、“Rake”表示滑動角(Rake)、“Slip”表示滑動量(Slip)Fig.5 The standard error of source parameters under the six different simulation settingsLet “TD” denote the Top Depth, and “BD” indicate the Bottom Depth. “Strike”, “Dip”, “Length”, “Rake” and “Slip” corresponding to the source parameters Strike, Dip, Length, Rake and Slip, respectively.

通過表3中仿真設置1和仿真設置2的結果對比、以及仿真設置3和仿真設置4的結果對比可以看出，在GPS觀測點個數相同的情況下，正態分布隨機噪聲越大，通過執行不同的采樣方法獲取的震源參數的中誤差也均相應增大.對比圖5a與5b、圖5c與5d也可以直觀看出，單位權中誤差為5 mm的隨機噪聲情況下得到的震源參數中誤差比3 mm情況下誤差顏色更深，表明不同誤差水平對震源參數精度的影響不同：隨機誤差越大，震源參數的精度越低.

對比表3中仿真設置1、3、5的結果可以看出，當GPS觀測點三個方向形變量數據的隨機噪聲為同分布情況下，將仿真設置1中的GPS觀測點個數增加一倍，并未對仿真設置5的參數中誤差結果造成較大的影響.并且隨著觀測點個數的增加，參數中誤差并未呈現出較明顯的遞減趨勢.通過對比仿真設置2、4的結果，圖5a、5c、5e，以及圖5b與5d的中誤差也可發現，當誤差水平相同條件下，觀測點的個數對震源參數的中誤差影響較小.因此，若試圖通過增加GPS觀測點數量來提升震源參數的精度可能并不十分顯著.

從表3中仿真設置6的結果可以看出，在GPS站點隨機分布情況下，Bootstrap方法依然能夠估出合理的震源參數中誤差，而Jackknife方法的結果與Monte Carlo方法結果偏離較大.這表明Bootstrap方法同樣適用于站點分布不均勻的普遍情況，并且能夠獲取比Jackknife方法更加合理的震源參數精度信息.對比仿真設置1與仿真設置6的結果可以看出，在相同的隨機誤差水平及GPS觀測點個數的條件下，通過執行不同采樣方法獲取的不均勻分布站點情況下的震源參數中誤差較均勻分布情況下的結果更大，這可能是因為不均勻分布站點對斷層面的約束力不及均勻分布情況，進而導致更大的震源參數中誤差.

3.2 Amatrice地震

2016年8月24日，意大利中部阿馬特里斯地區發生一次MW6.2地震，造成約300人傷亡且震中附近大量房屋被毀.美國地質調查局(United States Geological Survey, USGS)發布該地震的震中位置為(13.188°E，42.723°N)，震源深度為4.4 km，屬于淺層地震.震后，許多研究機構和學者快速對該地震進行了研究，分析了地震的震源機制、觸發關系以及庫侖應力(Wang et al., 2018)；亦有學者對該地震進行了震源參數反演和精度評估(Xu et al.，2017)，但針對不同的重采樣精度評定方法的性能評估對比、以及Jackknife方法和Bootstrap方法區間估計的對比分析未得到深入的研究.

本文獲取了Amatrice地震的106個GPS觀測臺站的形變量數據(https:∥ring.gm.ingv.it/)，將Amatrice地震GPS同震形變場在水平方向和垂直方向的實際形變值、利用Okada模型正演得到的預測值及其對應殘差分別繪于圖6(a、b、c).GPS形變整體擬合較好，水平方向的擬合殘差RMS為0.11 cm，垂直方向的擬合殘差RMS為0.36 cm，略小于Wang等(2018)得到的結果.水平方向的擬合精度整體上高于垂直方向，僅在主要形變區域的南側及西北側存在少量較大垂直殘差.

在獲取形變量數據后，經過坐標轉換將原始數據中的經緯度坐標轉換為以USGS發布的震中位置為原點的大地坐標；在均勻位錯模型和遺傳算法的基礎上，利用Bootstrap方法對該地震的形變量數據進行重采樣，通過遺傳算法搜索震源參數并獲取斷層參數的置信區間；再采用Bootstrap方法獲取參數的中誤差并與Jackknife方法和Monte Carlo方法的計算結果進行對比評估.參數的搜索區間、部分已有文獻研究結果、Monte Carlo方法獲取的中誤差，以及Jackknife方法和Bootstrap方法的反演估值、置信上下限、區間間隔和參數中誤差均列于表4.

表4 Amatrice地震震源參數搜索區間、反演估值、95%置信區間及參數中誤差計算結果Table 4 Search ranges, inversion values, and confidence intervals with 95% confidence level of source parameters for the Amatrice earthquake

圖6 Amatrice地震和Visso地震GPS同震形變場及殘差(a) Amatrice地震水平位移場； (b) Amatrice地震垂直位移場； (c) Amatrice地震的殘差圖； (d) Visso地震水平位移場； (e) Visso地震垂直位移場； (f) Visso地震的殘差圖；“Dis.”為實際形變值(Displacement)，“Pre.”為模型預測值(Prediction)，“Hor.”為水平方向殘差，“Ver.”為垂直方向殘差.Fig.6 The coseismic displacement fields and the corresponding residuals of Amatrice earthquake and Visso earthquake obtained by GPS(a) The horizontal displacement field of Amatrice earthquake； (b) The vertical displacement field of Amatrice earthquake； (c) The corresponding residuals of Amatrice earthquake; (d) The horizontal displacement field of Visso earthquake； (e) The vertical displacement field of Visso earthquake； (f) The corresponding residuals of Visso earthquake; “Dis.” is the actual displacement, “Pre.” is the model prediction; “Hor.” is the horizontal residual, and “Ver.” is the vertical residual.

根據表4中的震源參數估值可以看出，遺傳算法、Jackknife方法及Bootstrap方法獲取的Amatrice地震震源參數較為接近.由于各研究者采用的反演方法和數據的不同，這3種方法得到的震源參數與已有文獻的研究成果有所差異.通過對比各研究所得的矩震級可以發現，Bootstrap方法獲取的震源參數得到的矩震級與已有研究成果較為接近，之間的差異在僅在小數點后第2位.這證明了將Bootstrap重采樣與遺傳算法相結合、并用于解決Amatrice地震的震源參數反演問題是可行的.

由表4中Jackknife獲取的震源參數置信區間的上下限可以看出，頂深、底深及滑動量的Jackknife置信區間的下限均為負值.其中滑動量的Jackknife置信區間為(-4.15551,5.86072)即存在95%的概率相信該區間包含滑動量的真值，因此該結果說明Jackknife方法求取的頂深、底深及滑動量的置信區間包含的均值信息是模糊的或不精確的.通過對比Jackknife方法和Bootstrap方法獲取到的震源參數置信區間的區間間隔可知，在置信水平同為95%的情況下，采用Bootstrap方法獲取的區間間隔與Jackknife方法獲取的結果相差較大，Bootstrap方法獲取的9個震源參數置信區間均比Jackknife方法獲取的區間更窄.這表明與Bootstrap方法相比，Jackknife方法反演得到的震源參數的精確度較低，而Bootstrap方法獲取的參數估值更加可靠.

通過表4中的震源參數中誤差結果可以看出，Jackknife方法獲取的震源參數中誤差與Monte Carlo方法所得結果偏離較大，是因為Jackknife方法的重采樣過程僅僅能夠得到等于總觀測數的少量刀切樣本量所致，因此Jackknife方法會低估Amatrice地震的震源參數精度，無法提供反演估值的有效精度信息.Bootstrap方法獲取的震源參數中誤差小于Jackknife方法的結果，盡管Bootstrap方法獲取的震源參數中誤差比Monte Carlo方法的結果略大，但這兩種方法得到的結果仍在符號和數量級上均保持一致、在數值上也較為接近.造成Bootstrap方法的中誤差結果略大的原因可能是方法采樣的隨機性以及非線性反演算法自身無法避免的不穩定性.實驗結果說明了本文震源參數精度評定的Bootstrap方法的有效性，也表明Bootstrap方法在獲取震源參數精度信息方面較Jackknife方法更加精確.從9個參數的中誤差數值結果可以看出，對于Amatrice地震，走向角、傾角、滑動角的中誤差比其他6個震源參數更大，這可能是由于本文所采用的106個GPS觀測臺站的形變量數據中的隨機噪聲對角度量的反演產生較大影響所致，該現象也從側面反映出本文給出的震源參數精度評定的Bootstrap方法能夠合理反映Amatrice地震各個震源參數的精度信息.

3.3 Visso地震

Amatrice地震發生的兩個月后，震中30km范圍內于2016年10月26日又發生一次6.1級地震.根據USGS消息，此次地震的震中位于13.067°E、42.956°N的維索西部約3km地區.地震觸發后，許多學者和機構采用不同的數據和方法對該地震展開深入研究，分析了該地震的觸發機制和周邊斷層的應力變化(Wang et al., 2018)，反演了該地震的震源參數并獲取了參數精度(Xu et al.，2017).

本文獲取了127個GPS臺站形變數據(https:∥ring.gm.ingv.it/)生成Visso地震同震形變場，GPS水平方向的擬合效果較好，殘差RMS為0.17 cm，垂直方向殘差RMS為0.44 cm，與Wang等(2018)的計算結果較為一致.殘差的主要來源是GPS臺站的觀測誤差，而垂直方向殘差偏大是由于垂直觀測精度較水平方向精度低造成的.各臺站觀測得到的水平方向和垂直方向的實際形變值、模型正演得到的預測值及其對應殘差繪于圖6(d、e、f).

依據Okada彈性半空間位錯模型理論及遺傳算法，利用本文設計的精度評定方法對該地震進行震源參數搜索、置信區間計算及震源參數精度信息獲取；并與Jackknife方法得到的參數置信區間和震源參數中誤差、Monte Carlo方法得到的中誤差進行對比分析，從而驗證Bootstrap方法在實際地震的震源參數精度評定問題中的廣泛適用性.參數搜索區間和部分已有文獻研究結果，Jackknife方法和Bootstrap方法的反演估值、置信上下限、區間間隔和參數中誤差，以及Monte Carlo方法獲取的中誤差均列于表5.

對比表5中的遺傳算法、Jackknife方法及Bootstrap方法獲取的Visso地震矩震級可以發現，Bootstrap方法的結果為MW6.04，略小于USGS提供的結果(MW6.10)及Xu等(2017)的結果(MW6.08)，且比遺傳算法和Jackknife方法更加接近已有的文獻成果.由表5中的震源參數95%置信區間結果可以發現，USGS提供的震源參數及Xu等(2017)的研究結果均在本文Bootstrap方法求取的置信區間內.由Jackknife方法計算得到的傾角置信下限為-31.34522、置信上限為107.78006.而Okada彈性半空間位錯模型中的傾角是斷層面與水平面的夾角，其范圍僅為(0,90].說明Jackknife方法獲取的傾角置信區間內包含了一些無物理意義的信息，因此可信程度較低.由表中數據也可看出，Jackknife方法獲取的滑動量置信區間為(-8.68948, 11.53444)，區間長度達20.22392，說明Jackknife方法求取的滑動量置信區間過寬因此可靠性較弱.對比Jackknife方法和Bootstrap方法獲取的Visso地震震源參數置信區間能夠發現，Jackknife方法得到的9個參數置信區間寬度均比Bootstrap方法的結果偏大較多.以上分析說明了Bootstrap方法獲取的震源參數可信程度更高.

表5 Visso地震震源參數搜索區間、反演估值、95%置信區間及參數中誤差計算結果Table 5 Search ranges, inversion values, and confidence intervals with 95% confidence level of source parameters for the Visso earthquake

表5中的震源參數中誤差顯示，Jackknife方法獲取的震源參數中誤差與Monte Carlo方法所得結果差異較大，甚至出現量級不相等的情況.導致該現象是因為Jackknife重采樣容易受原始樣本影響，對于Visso地震，通過Jackknife方法重采樣僅能夠獲取127個少量刀切樣本，使得Jackknife方法獲取的震源參數中誤差因受樣本量的影響進而偏離參數的真實精度信息，因此Jackknife方法在實際震例的震源參數精度評定問題中可能不適用.Bootstrap方法獲取的震源參數中誤差小于Jackknife方法的結果，并且與Monte Carlo方法所獲結果在量級上一致、在數值上相近，說明Bootstrap方法能夠有效地起到對震源參數精度評定的作用，表現出比Jackknife方法更強的優勢.實驗結果證明了本文方法用于Visso地震震源參數精度評定的可行性和有效性，也證明了本文方法具有廣泛的適用性.從參數的中誤差數值結果可以看出，在Visso地震精度評定研究中，分別使用Bootstrap方法和Monte Carlo方法獲取的走向角、傾角、滑動角的中誤差比其它6個參數的中誤差更大，可能是因為本文采用的127個GPS臺站形變量數據中的隨機噪聲對Visso地震的角度量參數反演產生較大影響所致.該現象也再次從側面反映了Bootstrap方法的確能夠合理反映實際震例的各個震源參數的精度信息.

4 結論

為獲取震源參數在統計意義上更加合理的精度信息，若采用傳統的泰勒級數展開方法在計算過程中無法獲取復雜多維非線性函數的偏導數.本文在Okada彈性半空間矩形位錯模型及非線性反演遺傳算法的基礎上，引入了基于重采樣的Bootstrap方法，定義了一種通過實施自助采樣并解算自助樣本便可獲取震源參數中誤差的精度評定方法，并給出了具體的采樣步驟和迭代流程圖.最后將震源參數精度評定的Bootstrap方法用于模擬地震、Amatrice地震和Visso地震實驗中，通過對震源參數進行非線性反演、置信區間計算以及參數中誤差獲取，并與Jackknife方法計算得到的置信區間和中誤差、Monte Carlo方法獲取的中誤差進行對比，對本文方法進行了驗證.

模擬實驗結果顯示，采用本文方法獲取的震源參數均值與模擬真值更為接近；同時，在不同的隨機噪聲水平、觀測點個數以及GPS站點分布情況下，本文方法均能得到合理的精度信息.這表明將基于重采樣的Bootstrap方法用于震源參數非線性反演及精度評定問題是可行且有效的.在Amatrice地震和Visso地震實驗中，本文對比的已有研究結果均在Bootstrap方法求取的置信區間內，且Bootstrap方法獲取的震源參數95%置信區間均比Jackknife方法得到的區間長度更窄.與Monte Carlo方法獲取的震源參數中誤差相比，Jackknife方法得到的參數中誤差與之差異較大.盡管Bootstrap方法獲取的震源參數中誤差略大于Monte Carlo方法的結果，但這兩種方法得到的結果仍在符號和數量級上均保持一致、在數值上也較為接近.實際地震結果驗證了本文震源參數精度評定的Bootstrap方法的有效性和可靠性.Monte Carlo方法的結果依賴于三個方向形變量的平差值以及模擬隨機誤差時的協方差陣，并且該方法必須對總體模型做出具體的假設，例如需要假設擾動項的概率分布為高斯分布，而所得的結論只適用于該特定的數據生成過程.此時，不需要對未知模型及分布做任何假設的Bootstrap方法則表現出更強的適用性.隨著Bootstrap方法引入震源參數精度評定問題，將拓展該采樣方法在大地測量地震反演領域中的應用，對認知震源參數反演精度也具有一定的實際意義.

致謝感謝審稿專家對本文提出的寶貴意見.文中部分圖件是采用開源軟件GMT繪制.