青海地區(qū)地震序列參數(shù)的早期特征研究

余娜 張曉清 袁伏全

青海省地震局，西寧市興海路1號 810001

0 引言

地震序列是地震活動的主要表現(xiàn)形式之一，是統(tǒng)計地震學(xué)的重要研究對象。不同的構(gòu)造區(qū)域、地震序列類型、震源區(qū)構(gòu)造應(yīng)力水平和區(qū)域大地熱流值等均可能表現(xiàn)為地震序列參數(shù)的差異（Kagan et al，2010），且余震序列的衰減特征和激發(fā)余震的能力可通過序列的統(tǒng)計參數(shù)來表征，因此，獲得準確可靠的地震序列參數(shù)對快速判定地震序列類型、研究震源區(qū)特征和評估后續(xù)地震危險性等均具有重要的參考價值（蔣海昆等，2007；Guo et al，1997）。目前國際上對地震序列參數(shù)的計算主要采用時間序列的“傳染型余震序列”（epidemic type aftershock sequence，簡稱 ETAS）模型（Ogata，1988、1989；Zhuang，2011）。蔣海昆等（2006）系統(tǒng)研究了中國大陸5級以上地震序列的統(tǒng)計特征，并從物理上討論了不同條件下的序列衰減和余震激發(fā)問題。蔣長勝等（2013a、2013b、2013c）對1976年唐山地震序列、2013年岷縣漳縣MS6.6地震序列和蘆山 MS7.0地震序列利用ETAS模型進行了研究。余娜等（2016a、2016b）采用時間序列的ETAS模型計算了2010年玉樹MS7.1、2016年門源MS6.4地震序列。

青海地區(qū)是西北地震活動最為活躍的地區(qū)之一，是地震學(xué)家高度關(guān)注的地震危險區(qū)域，歷史上多次強震的發(fā)生，積累該區(qū)域內(nèi)更多、更可靠的地震序列參數(shù)，有利于震后地震趨勢估計和余震危險性研究。

1 研究資料

1.1 時間序列ETAS模型

用于擬合地震序列的 ETAS模型（Ogata，1988、1989；Zhuang，2011），其假設(shè)所有的余震均可按照大森-宇津公式（Omori，1894；Utsu，1961）激發(fā)自己的余震

式中，f（t）為t時刻對應(yīng)的地震序列發(fā)生率；t為主震發(fā)生后的離逝時間；μ為背景地震發(fā)生率；p表示余震序列衰減的快慢；c為主震后余震頻次達到峰值時所對應(yīng)的時間；K為常數(shù)，用于描述余震的活躍程度。ETAS模型假定主震的發(fā)生為初始零時刻，在其后的一個觀測時間段［0，T］內(nèi)的地震序列｛（ti，Mi）；i＝1，2，…，N｝的強度函數(shù)可表示為（Ogata Y，1988）

其中，Mi和ti分別為第i個事件的震級和發(fā)生時間；M0為參考震級，一般可取截止震級。式（2）中的參數(shù)p表示地震序列衰減的快慢，p越大衰減越快，反之則越慢；參數(shù)α表示觸發(fā)次級余震的能力（Ogata，1989），α越大表示觸發(fā)次級余震的能力越弱，反之亦然。對于震群型序列，一般情況下 α＜1，而當?shù)卣鹦蛄兄袩o明顯的被激發(fā)的次級余震時，α一般大于 1（Ogata，2001）；參數(shù)c為數(shù)值極小的正常數(shù)，一般用于調(diào)節(jié)當?shù)卣鹦蛄性谥髡鸢l(fā)生后極短時間內(nèi)余震記錄不完整，以及避免式（2）中分母為0的情況。為確保參數(shù)擬合的穩(wěn)定性，當區(qū)域內(nèi)背景地震發(fā)生率不高時，可通過設(shè)置背景地震的發(fā)生率μ＝0來降低擬合難度。

1.2 研究區(qū)域和序列目錄完整性分析

本文選取的“青海地區(qū)”的空間范圍為 32.0°～39.0°N，89.5°～103.0°E。該研究區(qū)域位于青藏高原中北部及東北緣，地質(zhì)構(gòu)造規(guī)模大、活動性強，主要有阿爾金斷裂帶、祁連山-河西走廊活動斷裂帶、昆侖山斷裂帶、巴顏喀拉斷裂帶和唐古拉斷裂帶等，以NNW、NW-NWW向斷裂帶最為發(fā)育。地震具有頻次高、強度大、分布廣的特點（馬玉虎等，2006；曾秋生，1999）。2009年以來研究區(qū)主震震級MS≥5.0地震有12個，其中，5.0≤MS＜6.0的有8個，6.0≤MS＜7.0的有3個，MS≥7.0的有1個。受ETAS模型計算對地震序列質(zhì)量的要求，本文使用了記錄較好的5級以上地震序列共9個（圖1），分別為2009年8月28日大柴旦MS6.4、2010年3月24日西藏聶榮MS5.8、2010年4月14日玉樹MS7.1、2011年6月26日囊謙MS5.2、2013年1月30日雜多MS5.2、2013年9月20日門源MS5.1、2015年11月23日祁連MS5.0、2016年1月21日門源MS6.4和2016年10月17日雜多MS6.2地震。

本研究所采用的地震目錄為中國地震臺網(wǎng)中心提供的《全國統(tǒng)一正式目錄》。青海地區(qū)地震目錄存在著明顯的空間不均勻現(xiàn)象，說明不同時期、不同地區(qū)地震記錄水平存在明顯的差異。在實際的地震序列選取中，采用了蔣長勝等（2014）區(qū)分地震序列與背景地震的時-空分布界限的“自然邊界法”，通過緯度-時間圖、經(jīng)度-時間圖及震中分布圖相結(jié)合的方式選取地震序列。圖2給出了2009年大柴旦MS6.4地震序列的選取結(jié)果。依照上面的選取方法，選取了其余8個地震序列。此外，由于在一些強震發(fā)生后的短期內(nèi)，主震的波形振幅較大，面波等波列持續(xù)時間較長，可出現(xiàn)將隨后大量發(fā)生的震級較小的余震“淹沒”而導(dǎo)致余震區(qū)甚至更大范圍內(nèi)的地震監(jiān)測能力顯著降低（Akaike，1974）的現(xiàn)象。為確保地震序列的完整性，采用“震級-序號”法（Huang，2006；蔣長勝等，2014；Zhuang et al，2012）確定了每個地震序列完整性震級Mc（圖3）?！罢鸺?序號”法按地震發(fā)生時間的先后順序排序，地震密度較大的區(qū)域的連線大致為截止震級Mc的時序變化。

圖1 研究區(qū)9個地震序列主震的空間分布

圖2 2009年大柴旦MS 6.4地震序列的選取結(jié)果

2 計算結(jié)果

2.1 地震序列ETAS模型擬合情況

ETAS模型參數(shù)一般使用最大似然法進行估計。在擬合時間［S，T］范圍內(nèi)（0＜S＜T），似然函數(shù)L的形式為

將式（2）帶入式（3），即可對序列參數(shù)［μ，K，c，α，p］進行最大似然估計。作為示例，圖4給出了2009年海西MS6.4地震序列在Mc＝ML1.9情況下的ETAS模型參數(shù)結(jié)果，獲得ETAS模型參數(shù)為 μ＝0.0000，K＝0.0280，c＝0.0035，α＝1.2411和 p＝1.2650。為進一步考察 ETAS模型的擬合效果，常使用“轉(zhuǎn)換時間”域的“殘差分析”法（Ogata，1988）?！皻埐罘治觥狈椒蓪?fù)雜的地震序列｛ti｝轉(zhuǎn)換為齊次泊松過程，并在“轉(zhuǎn)換時間”域｛τ｝中進行擬合效果分析，相應(yīng)的轉(zhuǎn)換過程為

這種經(jīng)過轉(zhuǎn)化的｛τ｝域的數(shù)據(jù)被稱為“殘差點過程”（RPP）。若在殘差分析中，轉(zhuǎn)換時間域的RPP在｛τ｝上的累計頻次近似直線，并與理論直線較為接近（Zhuang et al，2012），則ETAS模型擬合效果較好。以2009年大柴旦MS6.4地震序列為例，圖4（a）給出了主震之后0～91天和Mc＝ML1.9時的“殘差分析”結(jié)果，可以看出，序列累計地震數(shù)與理論ETAS擬合曲線較為接近，ETAS模型對此地震序列的擬合較為理想；圖4（b）為相應(yīng)的｛τ｝域的 M-τ分布。研究區(qū)9個地震序列參數(shù)的估算結(jié)果見表1。

2.2 序列參數(shù)穩(wěn)定性分析

ETAS模型參數(shù)與截止震級的選取有一定的關(guān)系（蔣長勝等，2013a、2013b）。本文以2009年8月28日大柴旦MS6.4、2010年4月14日玉樹MS7.1地震為例，研究了不同截止震級對地震序列參數(shù)的影響。2009年大柴旦MS6.4地震選取截止震級Mc＝ML1.8、1.9、2.1，發(fā)現(xiàn)地震序列參數(shù)α值和p值在Mc＝ML2.1時變化較小，誤差也?。▓D5）。2010年玉樹MS7.1地震選取截止震級Mc＝ML1.5、1.6、…、2.1，當截止震級 Mc從 ML1.5到2.1逐漸增大時，α值和p值隨截止震級的變化均有較大變化，參數(shù)k值逐漸減小，α值總體上有增加的趨勢，p值隨截止震級的增大而減小，但是變化幅度不大。因此，在考察地震序列參數(shù)的早期特征時，須考慮在何種截止震級條件下進行討論。

由于板內(nèi)地震所處的構(gòu)造運動速率不同，故地震序列持續(xù)時間存在明顯差異。因此，所選取的地震序列的擬合截止時間不同，獲得的ETAS模型參數(shù)可能不同，擬合參數(shù)的標準差也將存在差異（Iwata，2008）。以 2010年4月14日青海玉樹 7.1級地震為例，設(shè)定 tend＝［1.0，2.0，3.0，…，42.0］，當 Mc＝ML1.5時，ETAS模型參數(shù) α值和 p值均在震后 14天出現(xiàn)突跳變化，標準差幅值較大，尤其是激發(fā)次級余震能力的α值變化最為劇烈，p值在震后18天內(nèi)有一定幅度波動。b值變化相對穩(wěn)定，但在震后14天內(nèi)有逐漸增加的變化，隨著序列持續(xù)時間的增加α值減小。地震早期ETAS模型參數(shù)的α值、p值變化幅度較大，其后，趨于相對穩(wěn)定。

圖3 研究區(qū)9個地震序列的“震級-序號”法處理結(jié)果

圖4 2009年8月28日大柴旦MS 6.4地震序列擬合情況的“殘差分析”

表1 研究區(qū)9個地震序列主震參數(shù)和地震序列ETAS模型參數(shù)值

圖5 2009年8月28日大柴旦MS 6.4地震 ETAS模型計算結(jié)果

3 青海地區(qū)地震序列參數(shù)的時變特征

地震序列參數(shù)的早期變化比較劇烈，反映了地震序列主震發(fā)生后應(yīng)力快速調(diào)整的過程，余震也主要集中發(fā)生在地震早期階段。為了研究地震序列參數(shù)在時間上的變化特征，本文將地震序列不同截止條件下的ETAS模型參數(shù)穩(wěn)定時間和震級進行統(tǒng)計（表2），可以看出，地震序列以不同的截止震級作為研究對象時，地震序列的穩(wěn)定時間 Tα和Tp是不同的，也就是說，在一個地震序列中，地震序列參數(shù)α值和p值的穩(wěn)定時間Tα和 Tp與截止震級沒有明顯的正相關(guān)關(guān)系，當 Mc＝ML1.8時，Tα＝19，Tp＝19；當 Mc＝ML1.9時，Tα＝17，Tp＝21；當 Mc＝ML2.1時，Tα＝23，Tp＝19。地震序列參數(shù) α值和 p值的穩(wěn)定時間Tα和 Tp與主震震級間沒有明顯的相依關(guān)系，主震震級較大的其穩(wěn)定時間也可能較短，震級較小的地震穩(wěn)定時間反而較長。比如玉樹MS7.1地震，地震序列參數(shù)α值和p值的穩(wěn)定時間為：當Mc＝ML2.0時，Tα＝15，Tp＝15；雜多 MS5.2地震，當 Mc＝ML2.0時，Tα＝18，Tp＝13；反映了玉樹 7.1地震震后應(yīng)力迅速調(diào)整的特征。

研究區(qū)不同震級檔地震序列參數(shù)的最長穩(wěn)定時間和最短穩(wěn)定時間分別為：地震序列主震震級 MS5.0～6.0，Tα，min＝9，Tα，max＝47，Tp，min＝13，Tp，max＝45；主震震級 MS6.0～7.0，Tα，min＝17，Tα，max＝57，Tp，min＝11，Tp，max＝26。

表2 地震序列參數(shù)穩(wěn)定時間統(tǒng)計表

4 討論與結(jié)論

為系統(tǒng)考察青海地區(qū)地震序列的ETAS模型，并對地震序列的震后早期特征進行分析，本文運用中國地震臺網(wǎng)中心提供的《全國統(tǒng)一正式目錄》，對青海地區(qū)的地震序列進行了篩選，發(fā)現(xiàn)青海地區(qū)地震目錄存在著明顯的空間不均勻現(xiàn)象，不同時期、不同地區(qū)地震記錄水平存在明顯的差異。在實際的地震序列選取中采用“自然邊界法”選取了研究區(qū)2009年以來的9個地震序列，利用“震級-序號”法確定了每個地震序列的完整性震級Mc，使用最大似然法對每個地震序列進行參數(shù)估算，獲得了如下認識：

（1）為檢測地震序列參數(shù)的穩(wěn)定性，考察了截止震級和擬合截止時間的選取對ETAS模型參數(shù)影響。研究發(fā)現(xiàn)，截止震級Mc對α、k和p值有一定的影響。在考察地震序列參數(shù)的早期特征時，須考慮在何種截止震級下進行討論；地震早期ETAS模型參數(shù)的α值、p值變化幅度較大，反映了地震主震發(fā)生后震源區(qū)應(yīng)力的快速調(diào)整過程，其后，地震參數(shù)的α值和p值趨于相對穩(wěn)定。

（2）青海地區(qū)地震序列的ETAS模型參數(shù)α值、p值的穩(wěn)定時間Tα、Tp與截止震級間沒有明顯的正相關(guān)關(guān)系，與主震震級間也無明顯的相依關(guān)系，主震震級較大的地震事件穩(wěn)定時間可能較短，震級較小的地震事件其穩(wěn)定時間也有可能較長。

（3）青海地區(qū)不同震級檔地震序列參數(shù)的最長穩(wěn)定時間和最短穩(wěn)定時間分別為：地震序列主震震級 MS5.0～6.0，Tα，min＝9，Tα，max＝47，Tp，min＝13，Tp，max＝45；主震震級 MS6.0～7.0，Tα，min＝17，Tα，max＝57，Tp，min＝11，Tp，max＝26。

由于研究區(qū)域強震數(shù)目和ETAS模型計算條件所限，僅對9個地震序列進行了估算，并對其時間特征進行了討論，且所得出的結(jié)論僅反映部分統(tǒng)計特征，今后還需要繼續(xù)積累更多的震例進行更深入研究。

致謝：研究中使用了中國地震臺網(wǎng)中心“全國地震編目系統(tǒng)”提供的“統(tǒng)一正式目錄”，日本統(tǒng)計數(shù)理研究所莊建倉教授編制、提供了本文所使用的ETAS模型和余震短期概率預(yù)測程序，在此一并表示感謝。