蔣長勝 韓立波 郭路杰

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

引言

據(jù)中國地震臺網(wǎng)測定，2014年2月12日17時19分50.3秒，新疆維吾爾自治區(qū)和田地區(qū)于田縣發(fā)生MS7.3強烈地震.鑒于距此次地震震中約100km范圍內(nèi)近年來曾先后發(fā)生2008年3月21日于田MS7.3和2012年8月12日于田MS6.3等地震，此次地震的發(fā)生引起各方廣泛關(guān)注.

不同的構(gòu)造區(qū)域、地震序列類型、震源區(qū)構(gòu)造應(yīng)力水平和區(qū)域大地熱流值等均可能表現(xiàn)為地震序列參數(shù)的差異（Kagan et al，2010），而地震序列參數(shù)的早期特征對震后序列類型快速判定、強余震預(yù)測、斷層愈合過程等研究具有重要的參考價值（蔣海昆等，2007）；此外，目前國際上在對Reasenberg-Jones方法（Reasenberg，Jones，1989）、短期地震概率STEP模型（Gerstenberger et al，2005）、“傳染型地震序列ETES”模型和“傳染型狀態(tài)－速率ERS”模型（Falcone et al，2010）以及時空“傳染型余震序列”（epidemic type aftershock sequence，簡寫為ETAS）模型（Helmstetter etal，2006；Werner et al，2011；Zhuang，2011）等的研究中，構(gòu)建區(qū)域地震概率短期預(yù)測模型時常需使用統(tǒng)一的平均序列參數(shù).因此，對類似強震頻發(fā)的新疆于田地區(qū)開展地震序列參數(shù)研究具有重要的科學價值和現(xiàn)實意義.

本文針對新疆于田地區(qū)近年來發(fā)生的2008年3月21日MS7.3地震、2012年8月12日MS6.3地震和2014年2月12日MS7.3地震，利用目前國際上廣泛采用的ETAS模型（Ogata，1988）開展地震序列參數(shù)研究，并在此基礎(chǔ)上比較分析這些參數(shù)的早期特征.

1 時間序列ETAS模型及其參數(shù)估計方法

ETAS模型假設(shè)所有余震均可按照大森－宇津公式（Omori，1894；Utsu，1961）激發(fā)自己的余震，且震級的分布是獨立的.進一步假定主震的發(fā)生為初始零時刻，在其后的一個觀測時間段［0，T］內(nèi)地震序列｛（ti，Mi）；i＝1，2，…，N｝的強度函數(shù)可表示為（Ogata，1988）

式中，t為主震發(fā)生后的離逝時間；M0為計算所用序列的截止震級；Mi和ti分別為第i個事件的震級與距主震的時間長度；μ為背景地震發(fā)生率；p表示序列衰減快慢；c為主震后余震頻次達到峰值時的時間長度；K表示余震的活躍程度；α表示觸發(fā)次級余震的能力（Ogata，1989，1992），對于震群型序列α一般小于1，而當?shù)卣鹦蛄兄袩o明顯的被激發(fā)的次級余震時，α一般大于1（Ogata，2001）.

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

將式（1）帶入式（2），即可對待估參數(shù)［μ，K，c，α，p］進行最大似然估計.

為考察ETAS模型的擬合效果，一般使用“殘差分析”（Ogata，1988；Daley，Vere-Jones，2003）方法，將地震序列轉(zhuǎn)換為在“轉(zhuǎn)換時間”（τ）域的分布，并考察實際地震序列與理論值的擬合情況.對條件強度函數(shù)λ（t）可采用下式對時間序列｛ti｝進行轉(zhuǎn)換：

這樣，采用上式即可將｛ti｝轉(zhuǎn)化為服從單位速率的穩(wěn)態(tài)泊松分布（Zhuang et al，2012）的｛τi｝.若ETAS模型對數(shù)據(jù)擬合得較好，在轉(zhuǎn)換時間域｛τi｝的地震累積曲線則表現(xiàn)為線性，接近標準穩(wěn)態(tài)泊松過程的理論直線.

2 研究區(qū)和計算使用的資料

本研究所選區(qū)域為新疆于田地區(qū)，空間范圍為34.0°—38.0°N，78.0°—84.0°E.該區(qū)域地處青藏高原北緣，受青藏高原向北推擠和塔里木地塊阻擋的影響，區(qū)域內(nèi)構(gòu)造活動強烈，主要活動斷裂包括康西瓦斷裂、大紅柳灘—郭扎錯斷裂和左旋走滑的阿爾金斷裂帶及其西南端張剪區(qū)的分支斷裂等.歷史上曾發(fā)生多次7級以上地震，其中1900年以來共發(fā)生5次7級以上地震，分別為1924年7月3日民豐MS7??地震、7月12日民豐MS7.2地震、1996年11月19日和田MS7.1地震、2008年3月21日于田MS7.3地震和2014年2月12日于田MS7.3地震，其分布如圖1所示.

本文中使用了中國地震臺網(wǎng)中心提供的《全國統(tǒng)一地震編目目錄》①全國統(tǒng)一地震編目（日報）查詢系統(tǒng).http：∥www.csndmc.ac.cn／newweb／catalog＿direct＿link.htm#.查閱日期：2014年2月17日.，并選用全國地震監(jiān)測臺網(wǎng)數(shù)字化改造完成后的2002年以來ML2.0以上地震事件.為選定新疆于田地區(qū)發(fā)生的2008年MS7.3、2012年MS6.3和2014年MS7.3地震序列，采用緯度－時間圖、經(jīng)度－時間圖及震中分布圖相結(jié)合的方式，根據(jù)地震時空分布的自然邊界選取地震序列.由于本文討論的是地震序列的參數(shù)早期特征，因此僅選取了時空上叢集性較強且易與周圍離散的非序列地震區(qū)分的時段，選取結(jié)果如圖2所示.其中，2008年MS7.3地震序列選取的長度為811天，2012年MS6.3地震序列長度為153天，受觀測時間所限，對2014年MS7.3地震僅選用了主震后3.5天的地震序列.

圖1 新疆于田地區(qū)活動構(gòu)造與地震分布Fig.1 Tectonic and distribution of earthquakes in the Yutian region of Xinjiang

圖2 2002年以來新疆于田地區(qū)地震活動時空分布.（a）緯度－時間分布圖；（b）經(jīng)度－時間分布圖Fig.2 Spatio-temporal distribution of seismicity in the Yutian region，Xinjiang since 2002（a）Latitude-time plot；（b）Longitude-time plot

研究表明，在一些強震發(fā)生后的短期內(nèi)，由于強震波形振幅較大、面波等波列持續(xù)時間較長，隨后大量發(fā)生的震級較小的余震會被“淹沒”，余震區(qū)甚至更大范圍內(nèi)的地震監(jiān)測能力會顯著降低（Iwata，2008）.為確保上述3個地震序列的目錄完整性，研究中使用“震級-序號”法（Huang，2006；蔣長勝，吳忠良，2011；Zhuang et al，2012）確定完整性震級Mc.“震級－序號”法按地震發(fā)生時間的先后順序排序，可用于識別地震記錄信噪比和人為因素等引起的監(jiān)測能力分段和不連續(xù)變化.該方法一方面通過提高Mc數(shù)值確保地震序列的目錄完整性，另一方面還可通過調(diào)節(jié)ETAS模型擬合的起始時間C0，使得有足夠數(shù)量的地震參與計算.通過對2008年MS7.3、2012年MS6.3和2014年MS7.3等3個地震序列的分析表明，震后短期內(nèi)3個地震序列均出現(xiàn)了監(jiān)測能力的顯著下降，且隨時間逐步回升的現(xiàn)象（圖3）.綜合考慮，對3個地震序列統(tǒng)一選定Mc＝ML3.0，并將C0分別設(shè)定為0.312，0.125和0.042天.

圖3 2008年3月21日（a）、2012年8月12日（b）和2014年2月12日（c）3個于田地震序列震級－序號圖.圖中黃色實心圓標出了主震位置，垂直虛線標出了截止震級為ML3.0時地震目錄完整的初始時間位置（C0表示距離主震的時間，單位：d）Fig.3 Magnitude-rank distribution of three earthquake sequences（a）21March 2008Yutian MS7.3earthquake sequence；（b）12August 2012Yutian MS6.3 earthquake sequence；（c）12February 2014Yutian MS7.3earthquake sequence.Yellow dots show the position of mainshock，and the vertical dashed lines indicate the starting time of ETAS model fitting，C0，which accords to the magnitude completeness of ML3.0

3 地震序列參數(shù)估計結(jié)果

下面對3個地震序列分別進行ETAS模型參數(shù)的最大似然估計.作為示例，圖4給出了2014年MS7.3地震序列的累積地震數(shù)與ETAS擬合曲線在“轉(zhuǎn)換時間”域的比較情況，圖5給出了由ETAS模型擬合得到的地震序列條件強度曲線，即單位時間內(nèi)的地震發(fā)生率.結(jié)果顯示，ETAS模型對2014年MS7.3地震序列擬合得較好，參數(shù)估計結(jié)果分別為μ＝0.0000，K＝0.010 7，c＝0.001 2，α＝1.682 0和p＝0.669 2.此外，還利用最大似然法（Utsu，1965；Aki，1965）對3次地震序列的參數(shù)b值進行了估算：

圖4 利用ETAS模型對2014年MS7.3于田地震序列的擬合（a）累積地震數(shù)在“轉(zhuǎn)換時間”（τ）域的擬合情況，圖中深藍色實線為累積地震數(shù)，虛線為ETAS理論曲線；（b）“轉(zhuǎn)換時間”域的地震序列Fig.4 Residual analysis results from fitting the ETAS model to M≥ML3.0 of the 2014Yutian MS7.3earthquake sequence（a）Cumulative numbers of earthquakes（the deep blue curve）against transformed timeτ with a comparison to the their theoretical values（the thick dashed line）；（b）Plot of magnitude against transformed time calculated according to equation（3）

圖5 ETAS模型給出的2014年于田MS7.3地震序列的條件強度曲線（a）和M-t圖（b）.圖（a）縱坐標上的頻次指每天發(fā)生的地震次數(shù)Fig.5 Conditional intensity of the 2014Yutian MS7.3earthquake sequence fitted by ETAS model（a）and M-t plot（b），where cutoff magnitude Mcof the earthquake sequence is ML3.0.The frequency on the ordindte axis in Fig.（a）means the number of earthquakes per day

式中：ˉM為包含Mc及以上地震的平均震級；ΔMbin為震級分檔寬度，一般取為0.1.采用Shi和Bolt（1982）給出的方法對不確定度δb進行如下估計：

式中n為計算所用地震數(shù)目.

由于地震序列參數(shù)與序列持續(xù)時間的選取有關(guān)（蔣長勝等，2013a，b），為研究新疆于田地區(qū)2008年MS7.3、2012年MS6.3和2014年MS7.3地震序列的參數(shù)特征，尤其是參數(shù)的早期特征，分別考察了不同序列持續(xù)時間下的α值、p值和b值結(jié)果.對2008年MS7.3地震序列分別設(shè)定序列持續(xù)時間t＝｛3，4，5，10，20，…，90，100，200，…，800｝天，2012年MS6.3地震序列分別設(shè)定t＝｛70，100，150｝天，2014年MS7.3地震分別設(shè)定t＝｛1，2，3，3.5｝天，參數(shù)計算結(jié)果如圖6所示.由于2012年MS6.3地震序列數(shù)目較少，無法研究其參數(shù)早期特征，因此僅與其它兩次地震進行序列參數(shù)的整體特征比較.

圖6 3個地震序列ETAS模型參數(shù)α值（a）、p值（b）和b值（c）隨序列持續(xù)時間的變化Fig.6 Temporal variation of theαvalue（a），pvalue（b）and b value（c）after the three mainshocks from ETAS fitting for the three earthquake sequences in this study

由圖6a可見，3個序列的α值存在一定差異，其中：2008年地震序列的α值在主震后20天內(nèi)的早期變化較為劇烈，其后穩(wěn)定在1.013—1.303之間；2012年地震序列的α值在1.300—1.636之間；2014年地震序列的α值在主震后3.5天內(nèi)相對穩(wěn)定，分布在1.464—1.727之間.由此看來，2014年地震序列的α值與2012年地震序列接近，但高于2008年地震序列在主震發(fā)生20天后的穩(wěn)定階段，表明其激發(fā)次級余震的能力相對弱于2008年地震序列.

由圖6b給出的p值比較結(jié)果可見：2008年地震序列p值在主震后20天內(nèi)的早期變化同樣較為劇烈，其后穩(wěn)定在1.149—1.295之間；2012年地震序列p值與2008年地震序列較為接近，在1.159—1.192之間；而2014年MS7.3地震序列p值遠低于前兩個地震序列，分布在0.669—0.800之間，展示了相對緩慢的序列衰減過程.

由圖6c給出的b值比較結(jié)果可見：2008年地震序列b值在主震后20天內(nèi)出現(xiàn)持續(xù)上升變化，其后穩(wěn)定在0.65左右；2014年地震序列b值也存在明顯的上升變化，自震后1.0天的0.58上升至震后3.5天的0.67；2012年地震序列b值明顯高于其它兩次序列，約為0.88—0.91.

4 討論與結(jié)論

本文針對新疆于田地區(qū)2008年以來的3次強震，即2008年3月21日MS7.3地震、2012年8月12日MS6.3地震和2014年2月12日MS7.3地震，利用時間序列ETAS模型進行了序列參數(shù)早期特征研究.選用截止震級Mc＝ML3.0，使用最大似然法分別估算了3個地震序列在不同持續(xù)時間的α，p和b等參數(shù)，并進行了比較研究，考察了序列參數(shù)的早期特征，得到如下結(jié)論：

1）盡管主震震中相距在100km內(nèi)，但3次地震的序列參數(shù)存在明顯差異.其中2014年MS7.3地震序列的α值與2012年MS6.3地震序列接近，但高于2008年MS7.3地震序列在主震發(fā)生20天后的穩(wěn)定階段，表明其激發(fā)次級余震的能力相對弱于2008年MS7.3地震序列.

2）2008年MS7.3地震序列的p值與2012年MS6.3地震序列較為接近，但高于2014年MS7.3地震序列.2014年MS7.3地震序列的較低p值表明當前序列衰減緩慢.

3）2008年MS7.3地震序列和2014年MS7.3地震序列的b值在震后均出現(xiàn)持續(xù)上升變化.其中2014年MS7.3地震序列的b值高于2008年MS7.3地震序列，但兩次地震序列的b值均遠低于2012年MS6.3地震序列.

4）2008年MS7.3地震序列的參數(shù)變化較為顯著，α值、p值和b值在震后20天左右趨于穩(wěn)定.而2014年MS7.3地震序列在震后早期的3.5天內(nèi)，α值和p值變化相對穩(wěn)定，僅b值出現(xiàn)持續(xù)上升變化.

本文所討論的新疆于田地區(qū)3次地震序列中，如不考慮截止震級的選取而僅從直覺比較，2014年MS7.3地震序列早期階段的衰減不應(yīng)慢于2008年MS7.3地震序列.由圖6c可見，2014年MS7.3地震序列的b值高于2008年MS7.3地震序列，表明小地震數(shù)量的比例相對較高.而在統(tǒng)一設(shè)定Mc＝ML3.0情況下，會損失掉更多小地震，結(jié)果可能造成了2014年MS7.3地震序列在早期階段衰減更慢的現(xiàn)象.

本文對較小空間范圍內(nèi)發(fā)生的不同地震序列以及同一地震序列不同持續(xù)時間的序列參數(shù)研究表明，序列參數(shù)均可能有較大變化.因此在利用Reasenberg-Jones、STEP和ETAS等短期地震概率預(yù)測模型進行“向前”的實際地震預(yù)測時，需要預(yù)先設(shè)定的區(qū)域平均序列參數(shù)可能與震后早期階段差異較大，并將顯著影響預(yù)測結(jié)果的可靠性.此外，由于2014年新疆于田MS7.3地震序列仍在持續(xù)中，本研究獲得的ETAS模型參數(shù)僅為序列早期階段不完整的結(jié)果，尚不能反映此次地震序列的整體特征.

本文屬“四川省蘆山‘4·20’7.0級強烈地震科學考察”工作，研究中使用了中國地震臺網(wǎng)中心全國地震編目系統(tǒng)提供的《全國統(tǒng)一地震編目目錄》；日本統(tǒng)計數(shù)理研究所莊建倉教授在訪問地球物理研究所期間編制、提供了ETAS模型計算程序，并予以指導(dǎo).劉杰研究員對本文提出了有益建議，評審專家對本文提出了建設(shè)性的意見.作者在此一并表示感謝.

蔣長勝，吳忠良.2011.2010年玉樹MS7.1地震前的中長期加速矩釋放（AMR）問題［J］.地球物理學報，54（6）：1501－－1510.

Jiang C S，Wu Z L.2011.Intermediate-term medium-range Accelerating Moment Release（AMR）prior to the 2010 Yushu MS7.1earthquake［J］.Chinese Journal of Geophysics，54（6）：1501－－1510（in Chinese）.

蔣長勝，吳忠良，韓立波，郭路杰.2013a.地震序列早期參數(shù)估計和余震概率預(yù)測中截止震級Mc的影響：以2013年甘肅岷縣漳縣6.6級地震為例［J］.地球物理學報，56（12）：4048－－4057.

Jiang C S，Wu Z L，Han L B，Guo L J.2013a.Effect of cutoff magnitude Mcof earthquake catalogues on the early estimation of earthquake sequence parameters with implication for the probabilistic forecast of aftershocks：The 2013 Minxian-Zhangxian，Gansu，MS6.6earthquake sequence［J］.Chinese Journal of Geophysics，56（12）：4048－－4057（in Chinese）.

蔣長勝，莊建倉，龍鋒，韓立波，郭路杰.2013b.2013年蘆山MS7.0地震序列參數(shù)的早期特征：傳染型余震序列模型計算結(jié)果［J］.地震學報，35（5）：661－－669.

Jiang C S，Zhuang J C，Long F，Han L B，Guo L J.2013b.Statistical analysis of ETAS parameters in the early stage of the 2013Lushan MS7.0earthquake sequence［J］.Acta Seismologica Sinica，35（5）：661－－669（in Chinese）.

蔣海昆，鄭建常，吳瓊，曲延軍，李永莉.2007.傳染型余震序列模型震后早期參數(shù)特征及其地震學意義［J］.地球物理學報，50（6）：1778－－1786.

Jiang H K，Zheng J C，Wu Q，Qu Y J，Li Y L.2007.Earlier statistical features of ETAS model parameters and their seismological meanings［J］.Chinese Journal of Geophysics，50（6）：1778－－1786（in Chinese）.

Aki K.1965.Maximum likelihood estimate of bin the formula logn＝a－bmand its confidence limits［J］.Bull Earthquake Res Inst Tokyo Univ，43（2）：237－－239.

Daley D D，Vere-Jones D.2003.An Introduction to the Theory of Point Processes：Volume 1：Elementary Theory and Methods［M］.2nd Edition.New York：Springer：17，33.

Falcone G，Console R，Murru M.2010.Short-term and long-term earthquake occurrence models for Italy：ETES，ERS and LTST［J］.Annals of Geophysics，53（3）：41－－50.doi：10.4401／ag-4760.

Gerstenberger M C，Wiemer S，Jones L M，Reasenberg P A.2005.Real-time forecasts of tomorrow’s earthquakes in California［J］.Nature，435：328－－331.doi：10.1038／nature03622.

Helmstetter A，Kagan Y Y，Jackson D D.2006.Comparison of short-term and long-term earthquake forecast models for southern California［J］.Bull Seismol Soc Am，96（1）：90－－106.

Huang Q.2006.Search for reliable precursors：A case study of the seismic quiescence of the 2000western Tottori prefecture earthquake［J］.J Geophys Res，111（B4）：B04301.doi：10.1029／2005JB003982.

Iwata T.2008.Low detection capability of global earthquakes after the occurrence of large earthquakes：Investigation of the Harvard CMT catalogue［J］.Geophys J Int，174（3）：849－－856.

Kagan Y Y，Bird P，Jackson D D.2010.Earthquake patterns in diverse tectonic zones of the globe［J］.Pure Appl Geo-phys，167（6／7）：721－－741.

Ogata Y.1988.Statistical models for earthquake occurrences and residual analysis for point processes［J］.J Amer Statist Assoc，83（401）：9－－27.

Ogata Y.1989.Statistical model for standard seismicity and detection of anomalies by residual analysis［J］.Tectonophysics，169（1／2／3）：159－－174.

Ogata Y.1992.Detection of precursory relative quiescence before great earthquakes through a statistical model［J］.J Geophys Res，97（B13）：19845－－19871.

Ogata Y.2001.Increased probability of large earthquakes near aftershock regions with relative quiescence［J］.J Geophys Res，106（B5）：8729－－8744.

Omori F.1894.On aftershocks of earthquakes［J］.J Coll Sci Imp Univ Tokyo，7：11－－200.

Reasenberg P A，Jones L M.1989.Earthquake hazard after a mainshock in California［J］.Science，243：1173－－1176.

Shi Y，Bolt B A.1982.The standard error of the magnitude frequency b-value［J］.Bull Seismol Soc Am，72（5）：1677－－1687.

Utsu T.1961.A statistical study on the occurrence of aftershocks［J］.Geophys Mag，30：521－－605.

Utsu T.1965.A method for determining the value of bin a formula logn＝a－bmshowing the magnitude frequency for earthquakes［J］.Geophys Bull Hokkaido Univ，13：99－－103.

Werner M J，Helmstetter A，Jackson D D，Kagan Y Y.2011.High-resolution long-term and short-term earthquake forecasts for Califonia［J］.Bull Seismol Soc Am，101（4）：1630－－1648.

Zhuang J.2011.Next-day earthquake forecasts for the Japan region generated by the ETAS model［J］.Earth Planets Space，63（3）：207－－216.

Zhuang J，Harte D，Werner M J，Hainzl S，Zhou S.2012.Basic models of seismicity：Temporal models［G］∥Community Online Resource for Statistical Seismicity Analysis.［2014-02-13］.http：∥www.corssa.org.doi：10.5078／corssa-79905851.