統(tǒng)計地震學的兩個基本定律對應(yīng)力的依賴關(guān)系＊

C lém entNarteau Svetlana Byrdina Peter Shebalin Danijel Schorlemm er

余震衰減速率Λ由大森-宇津定律描繪,即Λ(t)=k(t+c)-p,這里 t是距離主震的時間,k是余震序列的發(fā)生率,p是冪律指數(shù),c是冪律余震衰減速率起始前的延時。對于 c值的爭議很大[4-10],研究人員[11-12]斷定 c接近于零,大于零值就是 (至少部分是)由于主震后 (特別是主震尾波期間)數(shù)據(jù)不完整造成的人為結(jié)果。2004年加州帕克菲爾德M6.0地震[13]和 2004年日本新瀉縣中越M6.6地震[8]后的研究表明,主震后立即對微震進行仔細識別,估算出的 c值較小,但仍非零。我們對南加州和日本的研究結(jié)果表明,將分析范圍限定在較小主震震級和較大余震震級會得出代表早期余震系統(tǒng)性特征的非零 c值,而不是一個完全的人為結(jié)果。

對于南加州,我們使用 1984—2003年有震源機制解[14]并經(jīng)重新定位的地震目錄和南加州地震臺網(wǎng) (SCSN)自 2003年以來的地震目錄。我們只選用那些解的誤差 ＜0.2、臺站分布系數(shù)≥0.5的高質(zhì)量地震事件。這一子目錄包含 12000個事件。為識別主震,我們采用基于冪律關(guān)系的時間和空間窗算法[15-16],此關(guān)系反映了地震現(xiàn)象的分級特性。沒有震源機制解的相應(yīng)余震按相同規(guī)則取自全部 SCSN目錄 (1984—2007年, ＞385000次事件)。

和南加州一樣,我們對于日本也使用兩個不同的地震目錄。我們將 F-Net地震臺網(wǎng)的目錄與 JMA(日本氣象廳)的目錄相結(jié)合,分離出自 1997年以來主震 (F-Net)及其各自余震 (JMA)的震源機制。為剔除日本不同火山區(qū)的群發(fā)性地震活動,我們還去除了距主震時間的幾何平均值超過 4小時的所有余震序列。因此也就去除了那些余震衰減速率不明顯的群發(fā)性空間地震活動。

為了避免由重疊記錄引起的人為結(jié)果,我們重點考慮由中等震級主震產(chǎn)生的余震序列,而不考慮大地震。同樣,我們只考慮較大余震,并根據(jù)主震時間對其進行疊加處理,以彌補每一序列中地震數(shù)量的不足。尤其是必須分別確定主震和余震的兩個震級范圍:我們將主震的最小震級閾值設(shè)定為因為這是重新定位的目錄的總體完整性水平[3]。然后,我們確定為重疊記錄開始對南加州目錄中M≥1.8余震的地震頻度-震級分布產(chǎn)生系統(tǒng)性影響的主震震級 (見附 2:方法)。余震震級范圍是一個可以改變從而產(chǎn)生不同 c值的自由參數(shù)。在所有情況下,我們都只選擇完整性地方震級小于(余震的最小震級閾值[17])的余震序列。對于每一次主震,利用擬合優(yōu)度檢驗法[18]對過去 3年間半徑為 15 km的圓內(nèi)這一完整性震級進行計算。由此,值不僅可以確定我們在單個序列中考慮的事件數(shù)量,而且可以確定被考慮的序列數(shù)量 (例如,南加州取則會在最初的 12小時期間產(chǎn)生近 2800次主震和 4000多次余震)。使用 A ki規(guī)則中的滑動角λ對主震的斷裂模式進行分級,即走滑型事件λ≈0°或λ≈180°,逆沖型事件λ≈90°,正斷層事件λ≈ -90°[3]。在此,我們既用滑動角也用滑動窗,其中滑動角的寬度為 60°,滑動窗的步長為 5°。這就產(chǎn)生了 73個主震等級。對于每一個等級和不同的震級區(qū)間,我們都根據(jù)余震距離各自主震的時間對余震進行疊加處理,并利用由模擬退火[4]實現(xiàn)的連續(xù)極小化和最大似然法對 c值進行估算。此外,我們還利用蒙特卡羅 (Monte Carlo)方法評估最大似然法估算 c值的不確定性,證實對第一個和第二個滑動角分別進行處理得出的觀測結(jié)果都是一致的。我們也用優(yōu)化程序的不同起始時間進行了一系列系統(tǒng)性的計算,以測試c值的穩(wěn)定性,因此去除了與目錄不完整相關(guān)的潛在的人為因素 (見附 2:方法)。

對于南加州,圖 1(a、b)顯示了 c值相對于余震震級分別為 1.8～2.8和 2.4～3.4的滑動角的演化。對于這些震級區(qū)間,圖 1(c～h)顯示了每一滑動角所選余震的數(shù)量以及所選主震和余震的平均震級。c值的不確定性對疊加中所用事件的數(shù)量有很強的依賴關(guān)系 (圖 1c、d),但即使是在震級閾值很高的情況下也會出現(xiàn)信號。走滑型地震的 c值比逆沖斷層型的高,但比正斷層型的低。這意味著逆沖斷層型地震的冪律余震衰減速率起始前的延時較短,正斷層型的較長,走滑型的居中。兩個震級范圍的特性也是如此,只不過會出現(xiàn)一般性的偏移值,這可由 c值對余震平均震級的依賴關(guān)系來解釋 (見附 2:方法)。逆沖型、走滑型和正斷層型地震之間 c值的差別不能用余震平均震級的變化(即相應(yīng)的 b值)來解釋,因為當按照滑動角繪制時,它呈現(xiàn)的不是一個系統(tǒng)性模式 (圖1e、f)。然而,兩個震級范圍的主震的平均震級都表現(xiàn)出類似特性 (圖 1g、h)。與走滑型地震相比,逆沖斷層型地震的主震較大,而正斷層型地震的主震較小。這是南加州的普遍特點,可由地震頻度-震級分布的差異解釋。

圖 1 南加州地區(qū)斷裂模式對冪律余震衰減速率起始前延時的影響。(a、b)分別與余震震級范圍(1.8,2.8)、(2.4,3.4)的滑動角相對應(yīng)的 c值的對數(shù)。紅、綠、藍分別表示走滑型、正斷層型和逆沖斷層型事件的滑動角區(qū)間。誤差棒由Mon te Ca r lo方法估算,相當于每個點上 500個合成余震序列的最大似然法估算 c值的 16%～84%的分位數(shù)。(c～h)主震和兩個余震范圍內(nèi)余震的統(tǒng)計資料。(c、d)所選余震的數(shù)量(N)。(e、f)所選余震的平均震級。(g、h)所選主震的平均震級

對于日本,為了證明缺失的余震不會影響 c值的估算,我們也對涉及震級為 3.0～4.5的 77次主震的 JMA目錄的完整性與Peng目錄[7]的完整性進行比較。然后,正如在南加州一樣,我們也發(fā)現(xiàn)相同的余震和主震震級范圍內(nèi) c值對滑動角有著相同的依賴關(guān)系 (圖 2a、b)。因為 F-Net目錄只包含少數(shù)幾次M＜3.5事件[3],所以所選主震的平均震級比南加州的要大。因此,估算的 c值總體較小。

圖 2 在日本地區(qū)斷裂模式對冪律余震衰減速率起始前延時的影響。(a、b)分別與余震震級范圍(1.8,2.8)、(2.4,3.4)的滑動角相對應(yīng)的 c值的對數(shù)。對 Peng[7]和 JMA目錄進行比較可以發(fā)現(xiàn),對于 1.8≤MA≤2.8的余震,優(yōu)化程序在主震發(fā)生后 5×10-4天開始。誤差棒由Mon te Car lo方法估算,相當于每個點上 500個合成余震序列的最大似然法估算 c值的 16%～84%的分位數(shù)

幾乎在所有情況下,地震頻度-震級分布都遵循一個以 b值 (Gutenberg-Richter頻度-震級關(guān)系指數(shù))為特性的冪律關(guān)系。對于用于識別日本和南加州主震的相同目錄,b值會因斷裂模式的不同而發(fā)生系統(tǒng)性變化。正如文獻 [3]中所描述,走滑型事件的 b值(～0.9±0.05)確實比逆沖型事件的 (～0.7±0.05)高,而比正斷層事件的 (～1.1±0.05)低。將兩個觀測結(jié)果結(jié)合起來看,c值和 b值都依賴于可由差異剪應(yīng)力解釋的震源機制,具體結(jié)果如下。

根據(jù)摩爾-庫倫 (Mohr-Cou lom b)理論,均勻介質(zhì)中受均勻三軸應(yīng)力 (σ1,σ2,σ3)作用的剪切破壞應(yīng)發(fā)生在最優(yōu)取向平面上。這些平面與中等剪應(yīng)力σ2平行,且與最大壓應(yīng)力σ1形成角Θ。因此,斷裂運動狀態(tài)受控于阿蒙頓 (Amonton)定律:

這里μ是靜摩擦系數(shù),σn和τ是在斷層面上確定的正應(yīng)力和剪應(yīng)力。假設(shè)其中一個主應(yīng)力是垂直的,另外兩個位于一個水平面內(nèi),那么就有可能表示出相對于μ而言斷層作用出現(xiàn)時的角Θ和差異剪應(yīng)力 (σ1-σ3)。在此情況下,主要剪應(yīng)力的取向不僅決定著斷裂機制,而且也決定著差異剪應(yīng)力的大小。例如,垂直應(yīng)力對應(yīng)的是正斷層作用的最大壓應(yīng)力和逆沖斷層作用的最小壓應(yīng)力。因此,對于相同的靜摩擦系數(shù),逆沖斷層作用的差異剪應(yīng)力要比正斷層作用的高得多。在這些擠壓和伸展狀態(tài)之間,走滑斷層作用會取 (σ1-σ3)的中間值。因此,c值與滑動角之間的關(guān)系表明余震衰減速率起始前的延時與差異剪應(yīng)力呈顯著的正相關(guān)。

考慮到余震源自主震附近應(yīng)力的階梯狀擾動,有關(guān)余震產(chǎn)生的各種模型預示這種擾動的幅度控制著短時期內(nèi)余震衰減速率的非冪律狀態(tài)的持續(xù)時間。例如,我們可以看到一個增加的應(yīng)力步長在速率-狀態(tài)摩擦模型[20]、靜疲勞模型[4,21]和破壞機理模型[22-23]中是如何加速冪律衰減速率的起始的。此外,經(jīng)典裂紋擴展實驗表明應(yīng)力擾動的幅度是遠場外加應(yīng)力的遞增函數(shù)[24]。圖1a、b和圖 2a、b中所表現(xiàn)出的 c值的特性印證了差異應(yīng)力、應(yīng)力重分布強度以及余震序列初始階段之短暫性之間的這種關(guān)系。然而,個別余震序列卻依然有高度分散,這極有可能是由應(yīng)力不均勻性效應(yīng)和測量不確定性引起。因此,c值對差異剪應(yīng)力的依賴關(guān)系只能通過余震的疊加來描繪。

b值決定著小地震和大地震數(shù)量之間的比率,其變化可與成核狀態(tài)下地震破裂能(低 b值 )否 (高 b值)傳播相關(guān)。由此,不同應(yīng)力狀態(tài)下 c值與 b值演變之間的正相關(guān)在地震傳播和余震產(chǎn)生期間都表現(xiàn)出相同的對應(yīng)力的依賴關(guān)系。在需要較高的差異剪應(yīng)力才能引起滑動的區(qū)域,孔隙度會降低,斷裂帶不同區(qū)段之間的相互作用增強,裂隙在不同方向上的大范圍內(nèi)蔓延速度更快[25-26]。于是,地震傳播期間,在較高圍限應(yīng)力下的過程區(qū)內(nèi)會呈現(xiàn)較短的破裂、滑移或破壞的特征時間尺度。因此,逆沖斷層型地震有可能比正斷層型地震傳播得更遠,這一結(jié)果與擠壓狀態(tài)下的 b值可能比伸展狀態(tài)下的 b值高的觀測結(jié)果一致。如果我們假設(shè)由于結(jié)構(gòu)和組分的不均勻性而導致破裂終止,那么同樣的推理依然有效,在高應(yīng)力環(huán)境下使余震成核的時間也就更短。由此可更快地實現(xiàn)穩(wěn)態(tài)松弛,逆沖斷層型地震的 c值也就比正斷層型的 c值更小。

地震頻度-震級分布和余震衰減速率對應(yīng)力的依賴關(guān)系為地震危險性分析和地震物理學提供了新的視點,因為它表明地震活動性的兩種獨立情況可用來約束脆性孕震地殼內(nèi)的應(yīng)力狀態(tài)。例如,地震活動背景和余震發(fā)生的準系統(tǒng)性可能不僅有助于量化沿主要斷層的應(yīng)力積累[27],而且還可以量化一次潛在地震的特征震級[28]。

附 1:方法概述

為識別主震,我們沒有選取那些震級小于M的地震,即在一次震級為M的地震事件發(fā)生后的 0.125×100.55M天、一個半徑為0.02×100.50Mkm的圓內(nèi)發(fā)生的地震[15-16]。利用同樣的空間尺度標準,在較大事件前不到 12小時內(nèi)發(fā)生的地震被視為可能的前震,并被從主震目錄中刪除。這樣一來,在最初的 12小時內(nèi)就不會出現(xiàn)重疊的余震序列。與所選主震相對應(yīng)的沒有滑動角的余震取自全部 SCSN目錄 (1984—2007年, ＞385000次事件)和 JMA目錄 (1997—2007年, ＞1100000次事件)。特別值得注意的是,我們選取了震級為M的主震后 12小時內(nèi)在一個半徑為0.02×100.50Mkm的圓內(nèi)發(fā)生的震級較小的地震。最后,只有 2%的所選余震被記錄在有震源機制解的地震活動目錄中。因此,這里所分析的數(shù)據(jù)集不同于確定 b值和滑動角之間的關(guān)系所使用的地震活動目錄[3]。

我們研究的只是來自特定震級范圍的地震,即使是在余震序列的早期,這些地震的資料也非常齊全。有如下兩個主要約束條件:(1)主震應(yīng)該足夠小,小到可產(chǎn)生較短尾波,這樣才能對隨后的余震進行可靠的探測;(2)余震應(yīng)該足夠大,以確保從余震序列的早期開始該特定震級水平上資料的完整性。這些震級范圍由地震學和統(tǒng)計學約束確定 (見附 2:方法)。

附 2:方法

不用區(qū)分斷裂模式,我們分別求出余震不同震級范圍的 c值。由此,我們觀測到冪率衰減之前的特征時間隨著平均余震震級的增大而持續(xù)增加。在加州和日本的單個地震事件中也有類似相關(guān)性報道[5,9]。這些觀測結(jié)果可以解釋圖 1a、b和圖 2a、b中不同余震震級范圍的平均 c值的變化。那么,平均 c值的絕對變化可以完全脫離 c值對滑動角的依賴關(guān)系。

在用最大似然法估算 c值時,優(yōu)化程序取決于兩個值,一個是 Tcomp,即目錄完整的時間 (所有余震都有記錄),另一個是 Tstart,即開始擬合的時間。很顯然,如果 c＜Tcomp,估算的 c值就不會被認為是可靠的。如果c＜Tstart,c值誤差就會很大,而且還會隨著Tstart以指數(shù)方式增大。在 Tstart＜Tcomp＜c的情況下,c值可能會因為目錄的不完整而被高估。因此,估算 c值的首選條件是Tcomp＜Tstart＜c。

對于日本,我們利用一個由人工挑選的目錄[7]來確定不受余震缺失影響的 c值和Tcomp值。對于不同余震震級范圍,我們證實c＞Tcomp。然后,采用 JMA目錄中列出的地震事件,我們用增大的 Tstart值向 c值 (即由人工挑選出來的 c值)的合理估算靠攏。這種靠攏取決于記錄事件的比例,也取決于它們的時間分布。對于 2.4≤MA≤3.4的余震,我們發(fā)現(xiàn) JMA目錄可被視為完整,我們所取的 Tstart值可以非常小。對于 1.8≤MA≤2.8的余震,我們推薦取 Tstart=5×10-4天,因為此時記錄事件的比例從零迅速轉(zhuǎn)變?yōu)?1。那么,還是這些地震事件,正斷層型事件的 c值較大,逆沖型事件的較小,走滑型事件的居中(圖2a)。

最重要的是,對于任何其完整性水平尚不確定的余震目錄,可以使用優(yōu)化程序的起始時間來估算 c值的穩(wěn)定性,并證明缺失的余震不會對該參數(shù)的估算有任何干擾。本研究中,我們對南加州的情況進行了類似分析。我們發(fā)現(xiàn)在 Tstart=10-4天 (9 s)、Tstart=5 ×10-4天 (45 s)、Tstart=10 ×10-4天 (90 s)和 Tstart=14×10-4天 (2 m in)情況下,所有等級的滑動角的 c值估算都沒有明顯變化。由于MA≥1.8的余震目錄在一次MM＜4.5的主震后 2m in時是完整的[10],所以 c值的這種穩(wěn)定性清楚表明缺失的余震在這一特例中 (即 Tcomp＜Tstart＜c,包括不同滑動角)不會對該參數(shù)的估算產(chǎn)生影響。

[1]GutenbergB,R ichterC F.Frequency of earthquakes in California.Bu ll.Seismo l.Soc.Am.,1944,34:185-188

[2]Utsu T.A ftershocks and earthquake statistics.J.Fac.Sci.Hokkaido Univ.,1965,Ser.V II3:379-441

[3]Schorlemm erD, Wiem er S, WyssM.Variations in earthquake-size distribution across different stress regim es.Natu re,2005,437:539-542

[4]Narteau C,Shebalin P,HolschneiderM.Temporal lim itsof the power law aftershock decay rate.J.Geophys.Res.,2002,107,B2359,doi:10.1029/2002JB001868

[5]Shcherbakov R,TurcotteD L,Rund le JB.A generalized Omori’s law for earthquake aftershock decay.Geophys.Res.Lett.,2004,31,L11613,doi:10.1029/2004GL019808

[6]V idale J E,Peng Z,IshiiM.Anom alous aftershock decay rates in the firsthundred seconds revealed from the H i-netborehole data.Eos Trans.AGU 85(FallMeet.Supp l.),2004,abstr.S23C-07

[7]Peng Z G,V idale J E,IshiiM,et al.Seism icity rate imm ediately befo re and afterm ain shock rup tu re from high frequencywaveform s in Japan.J.Geophys.Res.,2007,112,B03306,doi:10.1029/2006JB004386

[8]Enescu B,Mori J,MiyasawaM.Quantifying early aftershock activity of the2004m id-N iigata Prefecture earthquake(MW6.6).J.Geophys.Res.112,B04310,do i:10.1029/2006JB004629(2007).

[9]Nanjo K Z,Enescu B,Shcherbakov R,etal.Decay of aftershock activity for Japanese earthquakes.J.Geophys.Res.,2007,112,B08309,doi:10.1029/2006JB004754

[10]Kilb D,Martynov V,Vernon F.A ftershock detection asa function of tim e:results from theANZA seismic network following the 31October2001ML5.1 Anza,Califo rnia,earthquake.Bull.Seismo l.Soc.Am.,2007,97:780-792

[11]Kagan Y Y.Short-term p roperties of earthquake catalogs andmodels of earthquake source.Bull.Seismo l.Soc.Am.,2004,94:1207-1228

[12]Lo lliB,GasperiniP.Com paring differen tmodelsof aftershock rate decay:the role of catalog incompleteness in the first tim es afterm ain shock.Tectonophysics,2006,423:43-59

[13]Peng ZG,V idale JE,Houston H.Anomalousearly aftershock decay rate of the 2004MW6.0 Parkfield,California,earthquake.Geophys.Res.Lett.,2006,33,doi:10.1029/2006GL026744

[14]Hauksson E.Crustal structure and seismicity distribution ad jacent to the Pacific and North America plate boundary in sou thern Califo rnia.J.Geophys.Res.,2000,105:13875-13903

[15]Gardner J,Knop off L.Is the sequence of earthquakes in southern Californiaw ith aftershocks removed Poissonian?Bull.Seismo l.Soc.Am.,1974,64:1363-1367

[16]Reasenberg P.Second-order moment of central California seismicity,1969—1982.J.Geophys.Res.,1985,90:5479-5495

[17]Schorlemmer D, Woessner J.Probability of detecting an earthquake.Bull.Seismo l.Soc.Am.,2008,98,2103-2117

[18] Wiemer S, WyssM.Minim um magnitude of completeness in earthquake catalogs:examples from Alaska,the western United States,and Japan.Bu ll.Seismol.Soc.Am.,2000,90:859-869

[19]Sibson R H.Frictional constraintson thrusts,wrench and no rm al fau lts.Natu re,1974,249:542-544

[20]Dieterich J.A constitutive law for rate of earthquake production and itsapp lication to earthquake clustering.J.Geophys.Res.,1994,99:2601-2618

[21]Scholz C.Microfractures,aftershocks,and seism icity.Bu ll.Seismo l.Soc.Am.,1968,58:1117-1130

[22]Shcherbakov R,TurcotteD L.A dam agemechanicsmodel foraftershocks.Pure App l.Geophys.,2004,161:2379-2391

[23]Ben-Zion Y,Lyakhovsky V.Analysis of aftershocks in a litho spheric model with seismogenic zone governed by dam age rheology.Geophys.J.Int.,2006,165:197-210

[24]Atkinson B K.Subcritical crack grow th in geological materials.J.Geophys.Res.,1984,89:4077-4114

[25]Amitrano D.Brittle-ductile transition and associated seismicity:experimental and numerical studies and relationship with the b-value.J.Geophys.Res.,2003,B108,2044,doi:10.1029/2001JB000680

[26]O’ConnellD R H,Ma S,A rchuleta R J.Influence of dip and velocity heterogeneity on reverse-and no rm alfau lting rup ture dynam ics and near-fault ground motions.Bu ll.Seismo l.Soc.Am.,2007,97:1970-1989

[27]Narteau C,Shebalin P,Ho lschneiderM.Loading rates in California inferred from aftershocks.Non lin.Process.Geophys.,2008,15:245-263

[28]Schorlemm erD, Wiem er S.Microseism icity data forecast rup ture area.Nature,2005,434:1086

[29] WellsD L,Coppersm ith K J.New empirical relationships among magnitude,rupture length,rupture width,rupture area and surface displacement.Bull.Seismol.Soc.Am.,1994,84:974-1002