郭 星 潘 華
1)中國北京100081中國地震局地球物理研究所
2)中國北京100082環(huán)境保護部核與輻射安全中心
在實際地震危險性分析中對于地震發(fā)生過程的描述,最常用且數學上最為簡單的模型即泊松過程模型(Cornell,1968).該模型假定在一個斷層或一個震源區(qū)內,地震的發(fā)生符合泊松過程.泊松過程的特點是時間上的無記憶性,即未來發(fā)生地震的可能性與上一次地震的離逝時間無關,這與已知的強震發(fā)生的彈性回跳理論不一致.依據彈性回跳理論(Reid,1910),一次大地震發(fā)生以后,必須積累足夠的能量才能夠在同一地點發(fā)生下一次大地震.因此,若積累速率保持恒定,則距離上一次地震的時間越長,下一次地震發(fā)生的可能性就越大,此即大地震發(fā)生的時間記憶性.顯然,泊松模型不能反映這一特點.
為了描述大地震的記憶性,Utsu(1972)和Rikitake(1974)以及 Hagiwara(1974)基于彈性回跳理論提出一種更新模型,該模型假定地震的發(fā)生符合更新過程,即大地震的發(fā)生具有時間相關性,一次大地震發(fā)生后,還需要很長時間的彈性應力能量的積累才能在同一斷層段上發(fā)生下一次地震.
為了描述未來一段時期ΔT內發(fā)生潛在地震的可能性大小,國外研究人員提出很多用來計算這種條件概率的更新概率模型,包括雙指數分布(Utsu,1972)、高斯分布(Rikitake,1974)、對數正態(tài)分布(Nishenko,Buland,1987)、威布爾分布、伽馬分布(Utsu,1984)和BPT分布(Ellsworth,1999;Matthews et al,2002)等.
對于任一更新概率模型,若強震發(fā)生的概率密度分布函數f(t)和前一次地震的離逝時間T已知,就可以計算分段斷層源上未來一段時期ΔT內強震發(fā)生的條件概率(Wesnousky,1986):
然而,由于強震重復周期較長,而有強震記載的歷史較短,很多斷層源(或分段斷層源)上并沒有上一次強震發(fā)生的時間,這種情況下很難使用更新模型估計強震的發(fā)生概率.然而很多斷層源上已經有數百年甚至上千年的有確切記載的強震平靜期,若采用泊松過程模型則可能會低估強震發(fā)生的危險性.
針對缺少大震資料情況下活動斷裂上強震發(fā)生概率的計算,Matthews等(2002)對于給定的時間界限(已知或認為地震發(fā)生在這個時間點之前),給出了一個描述該時刻處于不同加載狀態(tài)的穩(wěn)態(tài)(steady-state)概率密度函數,然后利用全概率公式和BPT(Brownian passage time)概率密度函數計算得到未來一段時期ΔT內的發(fā)震概率.但該方法僅適用于基于一種隨機加載過程(穩(wěn)定加載附加布朗運動)的BPT模型.
本文對更新模型的條件概率計算方法進行了改進,提出一種以有確切記載的強震平靜期長度Ts為參數的條件概率計算方法.這里Ts要小于上一次地震的離逝時間T,而T是不確定的.本文提出的方法具有通用性,可以用于各種更新概率模型.
本文以東昆侖斷裂帶東段的塔藏段為研究對象,分別利用泊松分布模型與BPT分布模型計算該段的強震發(fā)生概率,并給出有確切記載的強震平靜期長度Ts與強震發(fā)生概率P的對應關系.
若無法確定斷層源上前一次大地震的離逝時間T,但是如果已知研究區(qū)內強震記載完整的初始年份及其距今的時間Ts,根據貝葉斯公式,則可利用強震發(fā)生的概率密度分布函數f(t)推斷斷層源上未來一段時期ΔT內發(fā)生強震的條件概率P(Ts,ΔT).
首先令至今至少存在長度為Ts的平靜期的事件為BTs,而現在的大震離逝時間為T的事件為AT;再以T年前時刻發(fā)生大震的事件為ET,其下一次地震將在未來時段發(fā)生的事件為X.根據上述定義,則有
若當今之前的歷史地震情況未知(無歷史大震和平靜期資料,即是無記憶的),根據等可能決策法(或拉普拉斯決策法),則可認為這個未知期內每一年發(fā)生地震事件(ET)的概率相同,即P(ET)為常數,這里設該常數為ν;若已知大震發(fā)生的概率密度分布函數f(t),則在T年前發(fā)生大震(ET)條件下,其下一次地震將在未來時段發(fā)生(X)的條件概率為
根據貝葉斯公式,在下一次地震將在未來時段發(fā)生(X)的條件下,該大震在T年前時刻發(fā)生的條件概率P(ET|X)為
式中,P(ET)為先驗概率,P(ET|X)為后驗概率,∫∞0(1-∫T0f(t)dt)dT則相當于一個常數.式(2)中,P(ET|X)即為大震離逝時間為T 的概率P(AT).
令P(BTs|AT)表示在大震離逝時間為T(AT)條件下,當前至少存在長度為Ts的平靜期事件(BTs)的概率.其表達式為
然后利用貝葉斯公式即可計算至今至少存在長度為Ts平靜期(BTs)的條件下,大震離逝時間為T(AT)的條件概率P(AT|BTs)為
式中,P(AT)為先驗概率,P(AT|BTs)為后驗概率.若P(AT|BTs)已知,利用全概率公式,則斷層源上未來一段時期ΔT內發(fā)生強震的條件概率為
為驗證本文方法的有效性,利用該方法計算了東昆侖斷裂帶東部塔藏段的強震發(fā)生概率.塔藏斷裂位于東昆侖斷裂帶東段(圖1),屬于“瑪曲空區(qū)”范圍(付俊東等,2012).該斷裂缺乏歷史強震數據,區(qū)域GPS數據結果顯示該區(qū)域可能處于應變積累階段(任金衛(wèi),王敏,2005).
圖1 東昆侖斷裂帶東段的幾何結構及塔藏段位置示意圖(引自李正芳等,2012)Fig.1 The geometrical structure schematic diagram of the eastern segment of East Kunlun fault and the location of Tazang segment(after Li et al,2012)
塔藏斷裂位于岷山強斷隆與甘南強隆區(qū)之間的邊界帶上,為一條全新世活動斷裂,斷裂長度60km,平均水平左旋滑動速率為2.7—2.8mm/a(李正芳等,2012).
東昆侖斷裂帶東部的塔藏段屬于分段明確的斷層源,本研究假定其符合Aki(1984)以及Schwartz和Coppersmith(1984)提出的特征地震模型.特征地震模型中最關鍵的就是對平均特征震級和平均復發(fā)周期這兩個參數的估計.
對于東昆侖斷裂帶東部塔藏段的平均特征震級,本文采用Wells和Coppersmith(1994)給出的走滑型地震的震級MW-破裂長度L的經驗關系式來估計,即
由式(7)計算得到塔藏段的平均特征震級ˉMW為7.2.
利用更新模型計算強震的發(fā)生概率,需要給出特征地震的平均復發(fā)間隔,也叫重復周期.由于塔藏段缺少古地震數據和大震的同震位錯數據,本文采用地震矩釋放率法估計其特征地震的重復周期(Youngs,Coppersmith,1985).
若已知活動斷裂帶的分段和各段的滑動速率,則其平均復發(fā)間隔為
式中:ˉM0為地震矩,可由特征地震的平均震級ˉMW確定;˙M0為地震矩釋放率,由斷裂段平均滑動速率、斷裂段長度、斷裂面寬度、地震比例因子等確定.其中地震比例因子是指發(fā)生特征地震的滑動占總滑動的比例.參照加州概率工作組和美國地震區(qū)劃圖,本文取地震比例因子為0.9(Frankel et al,2002;Working Group on California Earthquake Probabilities,2003).
震級為MW的地震所釋放的地震矩M0可由Hanks和Kanamori(1979)給出的矩震級-地震矩的經驗關系式計算得到,即
式中地震矩M0的單位為dyn·cm,而本文中地震矩的單位采用的是N·m,在計算過程中需要進行單位換算.對于斷層上年平均地震矩釋放率的估計,一般利用平均滑動速率S與斷層發(fā)震面積A的關系式計算得到,即
式中:μ為地殼巖石的剪切模量,本文中μ取3.3×1010N/m2;L和W 分別為發(fā)震斷層的長度和下傾寬度.
當沒有足夠的資料來確定斷層的下傾寬度時,可以采用Wells和Coppersmith(1994)給出的走滑斷層下傾寬度與震級之間的關系式來估計,即
將塔藏段的平均特征震級7.2代入式(11)中得到走滑斷層下傾寬度W 為15km;再根據式(10)計算得到塔藏段的年平均地震矩釋放率˙M0=8.3×1016N·m/a,去掉10%的蠕滑分量,結果為˙M0=7.5×1016N·m/a;最后將˙M0和ˉM0代入式(8)計算得到塔藏段特征地震的平均重復周期為1 059年.
本文所選取的更新模型為比較常用的BPT模型(Ellsworth,1999;Matthews et al,2002).利用BPT模型計算出塔藏段未來50年的強震發(fā)生概率,并與泊松模型的計算結果進行比較.
泊松模型的概率密度函數為f(t)=λexp-(λt),而未來ΔT年內強震發(fā)生的概率為P=1-exp(-λΔT).在泊松分布下,地震復發(fā)的條件概率與離逝時間無關,其未來年內地震ΔT年發(fā)生率是不變的,等于平均復發(fā)間隔τ的倒數(1/τ).
BPT模型的概率密度函數為
式中,τ為斷裂段的平均復發(fā)間隔,α為平均復發(fā)間隔的變異系數.
Ellsworth等(1999)在分析了37個地震序列后發(fā)現,α=0.5可作為活動構造所有不同震級的地震序列復發(fā)間隔變異系數的估計值,因此建議α可直接取0.5.Working Group on California Earthquake Probabilities(2003)的研究結果也表明,盡管在實際計算時采用了對3個不同α的計算結果進行加權平均的方式,即α(權值)為0.3(0.2),0.5(0.5),0.7(0.3),但實際上α取0.5也能恰當地反映出復發(fā)間隔的變異性.綜合上述分析,本文中α值取0.5.
已知BPT模型的概率密度函數f(t)和特征地震的平均重復周期τ,即可利用式(1)計算得到塔藏段不同離逝時間T所對應的未來50年強震發(fā)生的條件概率(圖2).利用式(6)則可計算得到塔藏段不同記載完整的平靜期長度Ts所對應的未來50年強震發(fā)生的條件概率(圖3).
圖2 塔藏段離逝時間T與未來50年強震發(fā)生概率的關系Fig.2 The relation between the length of time(T)since last earthquake and the occurrence probability(P)of large earthquakes for a 50-year exposure period on the Tazang fault segment
圖3 塔藏段上平靜期長度Ts與未來50年強震發(fā)生概率的關系Fig.3 The relation between the length of time(Ts)with historical record of quiet and the occurrence probability(P)of large earthquakes for a 50-year exposure period on the Tazang fault segment
本文同時還給出了利用泊松模型計算得到的條件概率.通過兩種類型的比較可知,只有當Ts與τ的比值非常小時,由式(6)計算得到的條件概率才接近于泊松分布;而當Ts接近于2倍的τ時,條件概率達到最大,之后開始下降.因此,若研究區(qū)內有較長時間的有確切記載的平靜期,則利用泊松模型估計強震發(fā)生概率可能會低估強震發(fā)生的危險性.
鑒于本文目的僅在于方法性探討,因此僅粗略地估計了研究區(qū)域M≥7.0地震資料基本完整的起始年.參考中國大陸分區(qū)地震資料基本完整的起始年分布圖像(黃瑋瓊等,1994),確定塔藏段M≥7.0強震資料基本完整的起始年為1561年.
已知研究區(qū)內M≥7.0強震資料完整的起始年為1561年,即研究區(qū)內距今已有452年的有確切記載的大震平靜期,因此可利用式(6)計算得到塔藏段未來50年強震發(fā)生的條件概率為P(452,50)=0.064 9,而利用泊松模型計算得到的條件概率則為0.046 1,其小于更新模型計算得到的結果.
本文在更新模型的基礎上,提出一種基于有確切記載的大震平靜期Ts的條件概率計算方法.該方法解決了更新模型無法在沒有歷史大震記載斷層上使用的不足,為強震發(fā)生概率的計算提供了一種新的選擇.
從本文的計算結果可以看出,只有當Ts非常小時,即有確切記載的大震平靜期長度較短的情況下,使用泊松模型才不會低估活動斷裂的強震發(fā)生概率.
我國有著悠久的地震史料記載,針對這些歷史地震資料,特別是強震資料完整性的研究,對于強震發(fā)生概率的估計具有重要意義.
本研究旨在提出一種新的發(fā)震概率計算方法,部分震源參數的計算過程比較簡單.因此,本文給出的塔藏段未來50年強震發(fā)生概率僅供參考,要獲得更精確的計算結果,尚需進行更深入的工作.
付俊東,任金衛(wèi),張軍龍,熊仁偉,楊攀新,陳長云,胡朝忠.2012.東昆侖斷裂帶東段塔藏斷裂晚第四紀古地震研究[J].第四紀研究,32(3):473-483.
Fu J D,Ren J W,Zhang J L,Xiong R W,Yang P X,Chen C Y,Hu C Z.2012.Research on Late Quaternary paleoearthquake on Tazang fault on the eastern section of the Kunlun active fault[J].Quaternary Sciences,32(3):473-483(in Chinese).
黃瑋瓊,李文香,曹學鋒.1994.中國大陸地震資料完整性研究之二:分區(qū)地震資料基本完整的起始年分布圖像[J].地震學報,16(4):423-432.
Huang W Q,Li W X,Cao X F.1994.Research on completeness of earthquake data on the Chinese mainland(Ⅱ):The regional distribution of the beginning years of basically complete earthquake data[J].Acta Seismologica Sinica,16(4):423-432(in Chinese).
李正芳,周本剛,冉洪流.2012.運用古地震數據評價東昆侖斷裂帶東段未來百年的強震危險性[J].地球物理學報,55(9):3051-3065.
Li Z F,Zhou B G,Ran H L.2012.Strong earthquake risk assessment of eastern segment on the East Kunlun fault in the next 100years based on paleo-earthquake data[J].Chinese Journal of Geophysics,55(9):3051-3065(in Chinese).
任金衛(wèi),王敏.2005.GPS觀測的2001年昆侖山口西MS8.1級地震地殼變形[J].第四紀研究,25(1):34-44.
Ren J W,Wang M.2005.GPS measured crustal deformation of the MS8.1Kunlun earthquake on November 14th 2001in Qinghai-Xizang Plateau[J].Quaternary Sciences,25(1):34-44(in Chinese).
Aki K.1984.Asperities,barriers,characteristic earthquakes and strong motion prediction[J].J Geophys Res,89(B7):5867-5892.
Cornell C A.1968.Engineering seismic risk analysis[J].Bull Seismol Soc Am,58:1583-1606.
Ellsworth W L,Matthews M V,Nadeau R M,Nishenko S P,Reasenberg P A,Simpson R W.1999.A Physically-Based Earthquake Recurrence Model for Estimation of Long-Term Earthquake Probabilities[R].Washington:U S Geological Survey:99-522.
Frankel A D,Petersen M D,Mueller C S,Haller K M,Wheeler R L,Leyendecker E V,Wesson R L,Harmsen S C,Cramer C H,Perkins D M,Rukstales K S.2002.Documentation for the 2002 Update of the National Seismic Hazard Maps.Open-File Report 02-420[R].Washington:U S Geological Survey:2-31.
Hagiwara.1974.Probability of earthquake recurrence as obtained from a Weibel distribution analysis of crustal strain[J].Tectonophysics,23:313-336.
Hanks T C,Kanamori H.1979.A moment-magnitude scale[J].J Geophys Res,84(B5):2348-2350.
Matthews M V,Ellsworth W L,Reasenberg P A.2002.A Brownian model for recurrent earthquakes[J].Bull Seismol Soc Am,92(6):2233-2250.
Nishenko S P,Buland R A.1987.A genetic recurrence interval distribution for earthquake forecasting[J].Bull Seismol Soc Am,77(4):1382-1399.
Reid H F.1910.TheMechanicsoftheEarthquake,The California Earthquake of April 18,1906[R].Washington:State Investigation Commission,Carnegie Institution of Washington,(2):43-47.
Rikitake T.1974.Probability of an earthquake occurrence as estimated from crustal strain[J].Tectonophysics,23(3):299-312.
Schwartz D P,Coppersmith K J.1984.Fault behavior and characteristic earthquakes:Examples from the Wasatch and San Andreas fault zones[J].J Geophys Res,89(B7):5681-5698.
Utsu T.1972.Large Earthquakes Near Hokkaido and the Expectancy of the Occurrence of a Large Earthquake of Nemuro[R].Coordinating Committee for Earthquake Prediction,7:7-13.
Utsu T.1984.Estimation of parameters for recurrence models of earthquakes[J].Bull Earthq Res Inst,Univ Tokyo,59:53-66.
Wells D L,Coppersmith K J.1994.New empirical relationships among magnitude,rupture length,rupture width,rupture area,and surface displacement[J].Bull Seismol Soc Am,84:974-1002.
Wesnousky S G.1986.Earthquake quaternary faults,and seismic hazard in California[J].J Geophys Res,91(B12):12587-12631.
Working Group on California Earthquake Probabilities.2003.Earthquake Probabilities in the San Francisco Bay Region:2002 to 2031.Open-File Report 03-214[R].Washington:U S Geological Survey.
Youngs R R,Coppersmith K J.1985.Implications of fault slip rates and earthquake recurrence models to probabilistic seismic hazard estimates[J].Bull Seismol Soc Am,75(4):939-964.