吳 果 冉洪流 周 慶 謝卓娟

1)中國地震局地質(zhì)研究所， 活動構(gòu)造與火山重點實驗室， 北京 100029 2)應(yīng)急管理部國家自然災(zāi)害防治研究院， 北京 100085

0 引言

中國領(lǐng)海幅員遼闊， 沿海地區(qū)是中國的經(jīng)濟(jì)重心。隨著近年來“21世紀(jì)海上絲綢之路”戰(zhàn)略構(gòu)想的提出， 濱海及海域的工程建設(shè)和資源開發(fā)日益增多。然而， 中國海域及鄰區(qū)受到歐亞板塊、 印度-澳大利亞板塊、 菲律賓海板塊和西太平洋板塊的相互作用(圖 1)， 構(gòu)造背景復(fù)雜， 地震頻發(fā)(劉光鼎， 1992； 臧紹先等， 1996； 彭艷菊等， 2008)。海域地震嚴(yán)重影響沿海地區(qū)的經(jīng)濟(jì)建設(shè)和社會穩(wěn)定， 也威脅到海洋資源開發(fā)的安全(丁海平等， 2011； 謝卓娟等， 2020)。因此， 分析中國海域及鄰區(qū)的地震活動特征， 建立相應(yīng)的地震活動模型成為日益迫切的需求， 也是中國下一代地震區(qū)劃圖的工作重點之一(李小軍等， 2020)。

圖 1 中國海域及鄰區(qū)震中分布圖(1767BC—2018年10月)Fig. 1 The distribution map of earthquakes in China seas and adjacent areas(1767BC—2018-10).高程數(shù)據(jù)引自文獻(xiàn)(1)https： ∥www.ngdc.noaa.gov/mgg/global/。

地震目錄是地震活動分析的重要基礎(chǔ)資料之一。與大陸地區(qū)不同， 海域的地震監(jiān)測臺網(wǎng)密度有限、 地震資料積累不足(林趾祥等， 1999)。長期以來中國的地震目錄編錄工作都集中于中國大陸及外延的近海地區(qū)(李善邦， 1960； 中央地震工作小組辦公室， 1971； 顧功敘， 1983a， b； 國家地震局震害防御司， 1995； 中國地震局震害防御司， 1999； Chengetal.， 2017)。這些成果在中國的抗震設(shè)防、 地震預(yù)測預(yù)報等工作中發(fā)揮了重要作用， 但是對于遠(yuǎn)海、 俯沖帶及周邊鄰近海域鮮有涉及， 只有吳戈等(2001)對黃海及其沿岸的歷史地震目錄進(jìn)行了編目和研究。依托于國家重點研發(fā)計劃項目“海域地震區(qū)劃關(guān)鍵技術(shù)研究”， 謝卓娟等(2020)通過廣泛收集全球地震臺網(wǎng)的地震資料， 首次編錄了中國海域及鄰區(qū)統(tǒng)一震級標(biāo)度為MS的地震目錄(圖 1)。該目錄為進(jìn)一步分析中國海域及鄰區(qū)的地震活動特征提供了基礎(chǔ)。

空間平滑模型(Smoothed Seismicity Model， SSM)是目前國際上應(yīng)用最為廣泛的基于地震目錄的地震活動模型之一， 該模型通過空間平滑函數(shù)/核函數(shù)(Kernel function)將記錄到的中小地震發(fā)生率向震中附近平滑分配以預(yù)測未來的發(fā)震概率分布(Kafkaetal.， 2000)。SSM自Frankel(1995)提出至今， 已經(jīng)在國際上得到廣泛應(yīng)用， 并成為美國國家地震區(qū)劃圖的一種基礎(chǔ)模型(Frankeletal.， 2000； Petersenetal.， 2008， 2015； Fieldetal.， 2014； Khodaverdianetal.， 2016； Akincietal.， 2018； Assatouriansetal.， 2019)。由于傳統(tǒng)的SSM對同一震級檔的地震的平滑半徑取相同值， 因此又被稱為固定平滑半徑的SSM(Petersenetal.， 2015)。

Helmstetter等(2007)提出了自適應(yīng)空間平滑地震活動模型(Adaptively Smoothed Seismicity Model， ASSM)， 對傳統(tǒng)的SSM做了改進(jìn)： 1)根據(jù)預(yù)測效果自動確定每個輸入地震的平滑半徑； 2)自動調(diào)節(jié)平滑半徑適應(yīng)地震的空間分布密度。ASSM模型連續(xù)2次在“地震可預(yù)測性合作研究”計劃(Co1laboratory for the Study of Earthquake Predictability， CSEP)的5a期中長期地震預(yù)測競賽中取得最佳成績(Zecharetal.， 2013； Straderetal.， 2017； Schorlemmeretal.， 2018)。因此， ASSM迅速地被廣泛應(yīng)用于全球的地震預(yù)測和地震危險性分析中(Kaganetal.， 2010； Werneretal.， 2010， 2011； Zecharetal.， 2010； Wangetal.， 2011； Helmstetteretal.， 2014； Moschetti， 2015； Talebietal.， 2017； Karaca， 2018)， 并于2014年被正式引入美國國家地震區(qū)劃圖(Petersenetal.， 2015， 2020)和美國加州地震概率工作組(Working Group on California Earthquake Probabilities， WGCEP)的UCERF3模型中(Fieldetal.， 2014)。

國內(nèi)和國際一些學(xué)者和機(jī)構(gòu)已經(jīng)多次將SSM用于中國研究區(qū)。Rong等(2002)將SSM用于評估中國大陸的潛在發(fā)震能力； 多位學(xué)者以中國北部或者中國大陸為研究區(qū)， 用SSM建立地震活動模型并以此為基礎(chǔ)開展地震危險性分析(Yangetal.， 2008； Xu， 2019； Fengetal.， 2020)； 中國第5代全國地震區(qū)劃圖中將SSM引入并作為確定潛源空間分布函數(shù)的因子之一(潘華等， 2013)； 最近的研究用SSM來描述中小地震的發(fā)震概率分布， 以配合斷層模型建立川滇地區(qū)或中國大陸的地震活動模型(程佳等， 2020； Rongetal.， 2020)。

國內(nèi)學(xué)者對于ASSM的研究較少， 只有Wu等(2019)將其用于建立青藏高原東部地區(qū)的中長期地震預(yù)測模型。在該工作中， Wu等(2019)針對中國歷史地震記錄時間長、 地震記錄水平隨時空變化的特點， 改進(jìn)了Helmstetter等(2007)的算法， 使其可以充分利用完整性隨時間和空間同時變化的地震記錄， 同時采用空間掃描的方法求取b值， 以更好地體現(xiàn)地應(yīng)力的分布(吳果等， 2019； 王少坡等， 2020)， 從而實現(xiàn)預(yù)測效果的優(yōu)化。

通過上述分析可知， 由于此前缺乏海域統(tǒng)一的地震目錄， 前人對于SSM的研究都是以中國大陸及鄰近地區(qū)為研究區(qū)， 尚未以整個中國海域及鄰區(qū)作為研究對象。只有張力方等(2013)在福州—天津沿線近海地區(qū)做過相關(guān)研究， 而其他海域的相關(guān)工作尚屬空白。此外， Helmstetter等(2007)提出的ASSM在地震中長期預(yù)測中表現(xiàn)優(yōu)異， 但國內(nèi)相關(guān)研究很少， 尚未將其應(yīng)用于中國廣大的海域地區(qū)。

本文基于謝卓娟等(2020)編錄的中國海域及鄰區(qū)的統(tǒng)一震級標(biāo)度的地震目錄， 采用分震級檔累積地震數(shù)-時間曲線分析各個地震帶的地震目錄完整性， 進(jìn)一步采用極大似然法擬合地震活動參數(shù)。在此基礎(chǔ)上， 采用Wu等(2019)改進(jìn)后的自適應(yīng)空間平滑算法， 充分利用完整記錄水平隨時間和空間變化的地震數(shù)據(jù)， 建立中國海域及鄰區(qū)M≥5.0地震年發(fā)生率的分布模型， 并采用概率增益函數(shù)對不同參數(shù)設(shè)置下的模型表現(xiàn)予以評價。最后分析了模型的優(yōu)點和局限性， 并給出使用建議。

1 地震目錄和完整性分析

1.1 地震目錄來源和前處理

本文的地震數(shù)據(jù)來自謝卓娟等(2020)新編錄的中國海域及鄰區(qū)的地震目錄， 該目錄匯集了中國大陸、 中國臺灣、 日本、 韓國、 菲律賓等國家和地區(qū)的地震數(shù)據(jù)， 時間跨度從公元前1767年—2018年10月， 是截至目前數(shù)據(jù)收集最為完整的中國海域及鄰區(qū)地震目錄， 共包含M2.0以上地震事件61285條。為了滿足后續(xù)用于編制海域地震區(qū)劃圖的需要， 該目錄的統(tǒng)計范圍除了渤海、 黃海、 東海和南海等主要海域外， 還向鄰區(qū)有所延伸， 乃至包括中國大陸沿海的地震區(qū)帶(圖 1 中紅色虛線框)。此外， 考慮到面波震級MS是國內(nèi)工程地震界通用的震級標(biāo)度， 為了與前人的工作保持延續(xù)性， 謝卓娟等(2020)采用正交回歸法擬合了不同震級標(biāo)度之間的轉(zhuǎn)換關(guān)系， 將地震目錄的震級標(biāo)度統(tǒng)一為面波震級MS(Xieetal.， 2021)。

類似前人在地震預(yù)測和地震危險性分析時的工作(Helmstetteretal.， 2007； Wangetal.， 2011； Petersenetal.， 2015； Talebietal.， 2017； 程佳等， 2020)， 本文采用國際上廣泛應(yīng)用的Gardner等(1974)提出的算法刪除前、 余震(Wiemer， 2001)。共發(fā)現(xiàn)地震叢集6353個， 包括31835個地震事件。除叢后的地震目錄包含29450條地震記錄， 其中M5.0以上主震分布見圖2a。

1.2 完整性分析

采用記錄時間長但不完整的地震目錄會給地震活動參數(shù)的評估帶來偏差， 無益于提高地震活動模型的質(zhì)量(徐偉進(jìn)等， 2014； Talebietal.， 2017)， 故在建立模型之前需要分析地震目錄的完整性。區(qū)別于Wu等(2019)的做法， 本研究并沒有采用空間掃描的方法求取每個網(wǎng)格點的地震目錄完整性， 這是基于2方面考慮： 1)海域的地震監(jiān)測臺網(wǎng)密度有限， 有地震完整記錄的時間較短(林趾祥等， 1999)； 2)如圖 1 所示， 研究區(qū)內(nèi)的地震分布極不均勻， 在東海、 南海等大片海域地震記錄十分稀少。以上原因?qū)е氯舨捎每臻g掃描法， 很容易遇到地震記錄達(dá)不到數(shù)據(jù)分析所需的樣本量下限的情形。

因此， 本研究延續(xù)中國第5代全國地震區(qū)劃圖的做法(潘華等， 2013)， 即以地震區(qū)帶為統(tǒng)計單元進(jìn)行完整性分析?！昂Ｓ虻卣饏^(qū)劃關(guān)鍵技術(shù)研究”項目組依據(jù)新編制的海域活動構(gòu)造框架圖中主要構(gòu)造單元的邊界對第5代區(qū)劃圖的地震區(qū)帶劃分方案做了部分調(diào)整， 使得地震區(qū)帶的邊界基本與二、 三級構(gòu)造單元的邊界重合(高戰(zhàn)武等， 2021)。主要調(diào)整了東海地震帶、 臺灣東部地震帶和南海地震區(qū)的邊界， 其中對臺灣東部地震帶進(jìn)行了調(diào)整并擴(kuò)充了范圍， 改稱臺灣南-馬尼拉海溝地震帶。同時， 沿琉球海溝劃分出琉球海溝地震帶。調(diào)整后的地震帶分布見圖2b。

圖 2 a 中國海域及鄰區(qū)M5.0以上主震的震中分布圖； b 研究區(qū)內(nèi)的地震帶分布圖(引自高戰(zhàn)武等， 2021)Fig. 2 The distribution map of M≥5.0 main shocks in China seas and adjacent areas(a)； the locations of seismic zones in the study area(b)(Cited from GAO Zhan-wu et al.， 2020).A 華北平原地震帶； B 郯廬地震帶； C 長江下游-南黃海地震帶； D 朝鮮地震帶； E 東海地震帶； F 琉球海溝地震帶； G 華南沿海地震帶； H 臺灣西部地震帶； I 臺灣南-馬尼拉海溝地震帶； J 南海地震帶

由于難以理清廣大海域的臺站分布和記錄能力隨時間的變化， 只能通過分析地震目錄判斷其完整性。常用的方法有G-R關(guān)系曲線法(Gutenbergetal.， 1944)和累積地震數(shù)-時間曲線法， 二者通過觀察曲線的線性程度是否優(yōu)良來判斷某起始震級的地震目錄是否完整(黃亦磊等， 2016； 史翔宇等， 2020)。由于G-R關(guān)系本身是對數(shù)線性關(guān)系， 對于小地震的缺失并不敏感(任雪梅， 2011； 徐偉進(jìn)等， 2014)。同時G-R關(guān)系曲線法需要預(yù)先假定高質(zhì)量地震目錄的時間起點， 不便于區(qū)分不同大小地震完整記錄的起始時間(吳果等， 2014)。通過累積地震數(shù)-時間曲線法可以直觀地觀察累積地震數(shù)隨著時間的變化趨勢， 便于判斷地震完整記錄的起始時間。但累積地震數(shù)-時間曲線法也存在缺陷， 當(dāng)活動水平隨時間波動的不同震級檔的地震數(shù)據(jù)混在一起時， 該方法難以加以區(qū)分(吳果等， 2014)。因此， 本研究采用分震級檔的累積地震數(shù)-時間曲線， 即以0.5個震級單位為間隔分檔繪出累積地震數(shù)-時間曲線， 再判斷該震級檔地震記錄完整記錄的起止時間。

各地震帶分震級檔的完整性分析結(jié)果見表1。由于地震帶和震級檔較多， 受篇幅限制， 本文僅就完整性分析過程中的關(guān)鍵問題進(jìn)行論述：

表 1 中國海域及鄰區(qū)地震帶各震級檔對應(yīng)完整記錄的起始年份Table1 The start years from which the earthquake catalogs within different magnitude intervals are completely recorded for the seismic zones in China seas and adjacent regions

(1)總體而言， 大陸內(nèi)部或近海的地震帶歷史地震記錄悠久， 儀器地震記錄大致始于1970年后， 如華北平原地震帶、 郯廬地震帶、 華南沿海地震帶等； 遠(yuǎn)海的地震帶缺乏歷史地震記錄， 且儀器地震數(shù)據(jù)完整記錄的時間很短， 如東海地震帶和南海地震帶， 這與臺網(wǎng)的布設(shè)時間以及密度、 精度直接相關(guān)。

(2)地震目錄的截止時間為2018年10月， 故2018年的地震記錄尚不完整， 因此表1 中默認(rèn)的截止時間是2017年。華北平原地震帶、 郯廬地震帶、 華南沿海地震帶以及臺灣西部地震帶的小震級檔存在2010年以后數(shù)據(jù)不全的現(xiàn)象。本文以華北地震帶M2.5～2.9的數(shù)據(jù)為例予以展示(圖3a)。該現(xiàn)象與中國地震臺網(wǎng)中心不再公布部分地區(qū)2010年后的小地震數(shù)據(jù)有關(guān)(Chenetal.， 2006)。

圖 3 分震級檔的累積地震數(shù)-時間曲線Fig. 3 The curve of cumulative number of earthquakes versus time for a magnitude interval.a 華北平原地震帶M2.5～2.9； b 朝鮮地震帶M4.5～4.9

(3)華北平原地震帶、 郯廬地震帶、 長江下游-南黃海地震帶和華南沿海地震帶M≥5.0的地震記錄從1500年以來記錄完整(表1)， 該結(jié)果與黃瑋瓊等(1994)、 潘華等(2006)和徐偉進(jìn)等(2014)的分析結(jié)果基本一致。

(4)如圖3b 所示， 朝鮮地震帶的M≥4.5地震記錄的完整記錄時間最早可以追溯到1400年， 這一判斷與翟文杰等(2004)的工作一致。此外， 朝鮮地震帶的地震目錄在1737—1963年約200a間缺少M4.5～5.9的地震記錄。假設(shè)地震活動服從泊松分布， 則可以估算出該現(xiàn)象屬于正常地震活動水平波動的概率約為10-31， 屬于極低概率事件。因此， 可初步判定該現(xiàn)象是由于地震記錄不完整或數(shù)據(jù)收集不充分導(dǎo)致的。然而， 由于本文采用的算法具有很大的靈活性， 同一震級檔對應(yīng)2個乃至更多個完整記錄時段的數(shù)據(jù)仍然可被有效利用。

(5)隨著震級增大， 高震級檔的地震數(shù)逐漸減少， 給完整性分析帶來困難。因此， 對于同一地震帶， 依據(jù)高震級檔的完整記錄起始時間不晚于相對低的震級檔的原則予以約束。

觀察表1 可知， 不同地震帶或不同大小的震級檔對應(yīng)的完整記錄時長差異巨大， 若想要模型充分地利用有限的地震記錄數(shù)據(jù)、 增加樣本量， 對上述分震級檔的地震目錄進(jìn)行完整性分析是很有必要的。

1.3 地震活動參數(shù)擬合

基于上述地震目錄完整性分析結(jié)果， 采用Weichert(1980)提出的方法計算各個地震帶的a值和b值(式(1))及b值的標(biāo)準(zhǔn)差(式(2))。該方法可利用不同震級檔對應(yīng)不同完整時段的地震數(shù)據(jù)(吳果等， 2019)。

(1)

(2)

各地震帶的地震活動參數(shù)擬合結(jié)果見表2， 表中Mmin為擬合參數(shù)時使用的地震數(shù)據(jù)的震級下限； 各地震帶的震級上限Mmax主要參考中國第5代地震動參數(shù)區(qū)劃圖(潘華等， 2013)和高戰(zhàn)武等(2021)的工作， 同時結(jié)合地震目錄中記錄到的最大歷史地震對其進(jìn)行調(diào)整。

表 2 中國海域及鄰區(qū)各地震帶的a值、 b值及b值的標(biāo)準(zhǔn)差Table2 Parameters including a-value， b-value and the standard deviation of b-value for each seismic zone

2 自適應(yīng)空間平滑模型原理

ASSM與傳統(tǒng)的SSM都屬于空間平滑模型， 二者的基本原理相同， 都是利用已經(jīng)發(fā)生的中小地震的分布預(yù)測未來的地震分布， 在一定程度上類似于國內(nèi)常用的地震“原地重現(xiàn)”理論。但空間平滑模型主要針對沒有明確發(fā)震構(gòu)造背景的中小—中強(qiáng)地震， 與“原地重現(xiàn)”理論存在一定差異。此外， 地震不可能絕對精確地“原地重現(xiàn)”， 未來地震的發(fā)生位置具有很大的不確定性。已有的地震密集或者稀疏的地區(qū)只能相應(yīng)地代表該地區(qū)未來發(fā)生地震的可能性更大或者更小， 具有概率含義。因此， ASSM和SSM都需要利用空間平滑函數(shù)將記錄到的中小地震發(fā)生率向震中附近平滑分配， 使地震發(fā)生率的高低變化在研究區(qū)內(nèi)呈現(xiàn)平滑過渡的趨勢， 而不是簡單地“落”在已有地震記錄的位置上。前人已對上述理論的可靠性進(jìn)行了驗證(Kafkaetal.， 2000)， 同時ASSM在CSEP計劃中的優(yōu)異表現(xiàn)也可看作是對該理論可靠性的重要檢驗(Straderetal.， 2017； Schorlemmeretal.， 2018)。

廣大海域地區(qū)缺乏大陸地區(qū)豐富的地震地質(zhì)、 地球物理探測等多學(xué)科研究數(shù)據(jù)， 地震目錄是開展其地震活動性研究最重要的基礎(chǔ)資料。同時， 如圖 1 所示， 廣大海域地區(qū)地震分布十分不均， 例如兩大海溝俯沖帶地震分布的密集程度遠(yuǎn)超過東海和南海地區(qū)。因此， 考慮到ASSM模型是目前國際上應(yīng)用十分成功的基于地震目錄的地震活動模型， 且能充分體現(xiàn)地震分布密集程度的空間變化， 本研究在中國海域及鄰區(qū)建立ASSM模型具有一定的合理性和應(yīng)用價值。

基于前文的地震目錄完整性分析和地震活動參數(shù)擬合的結(jié)果， 進(jìn)一步采用Wu等(2019)改進(jìn)的自適應(yīng)空間平滑算法建立ASSM模型， 下文將簡述其原理。

2.1 學(xué)習(xí)目錄和檢測目錄

ASSM與SSM的最大區(qū)別在于ASSM將整個地震目錄劃分為2部分(圖 4)， 第1部分為學(xué)習(xí)目錄(Learning catalog)， 第2部分為檢測目錄(Testcatalog)。首先用學(xué)習(xí)目錄建立模型， 再用檢測目錄中的地震評價預(yù)測效果， 從而實現(xiàn)參數(shù)優(yōu)化。學(xué)習(xí)目錄中的地震并非都能被利用， 在使用前要采用各個震級檔對應(yīng)的完整記錄時段對其進(jìn)行篩選(表1)， 篩選后的地震可稱為輸入地震(Input earthquakes)。國際上最新的地震危險性分析研究一般將起算震級定為M5.0(Petersenetal.， 2015)， 因此本研究將檢測目錄的目標(biāo)震級Mt定為M5.0， 即用檢測目錄中M≥5.0的地震測試模型的預(yù)測效果。

圖 4 自適應(yīng)空間平滑地震活動模型中的學(xué)習(xí)目錄和檢測目錄劃分示意圖Fig. 4 Schematic diagram of the division of learning catalog and test catalog in adaptively smoothed seismicity model.

2.2 空間平滑函數(shù)

將研究區(qū)網(wǎng)格化， 網(wǎng)格大小設(shè)為0.1°×0.1°， 然后采用空間平滑函數(shù)將學(xué)習(xí)目錄中的輸入地震發(fā)生率向四周的網(wǎng)格平滑分配以預(yù)測未來的發(fā)震概率分布。本文采用應(yīng)用最為廣泛的高斯平滑函數(shù)：

(3)

2.3 輸入地震對地震年發(fā)生率的貢獻(xiàn)

式(3)中默認(rèn)每個輸入地震對研究區(qū)內(nèi)地震年發(fā)生率的貢獻(xiàn)相等。為了改進(jìn)模型的地震預(yù)測效果， 需要盡可能充分地利用有限的地震數(shù)據(jù)， 擴(kuò)大樣本量。然而如表1所示， 不同地震帶、 不同震級檔的地震記錄來自不同時長的完整記錄時段， 如果不對式(3)加以調(diào)整， 則輸入地震不再符合原有的震級-頻度關(guān)系， 得到的模型也不合理。

為了解決上述問題， 需要求出每個輸入地震對研究區(qū)內(nèi)地震年發(fā)生率的貢獻(xiàn)v(i)。假設(shè)地震i屬于震級檔j， 那么v(i)為

(4)

(5)

利用G-R關(guān)系可以將最小輸入地震震級Mmin以上地震的年發(fā)生率轉(zhuǎn)化為任意目標(biāo)震級Mt以上地震的年發(fā)生率：

(6)

(7)

(8)

式中，Nyear為檢測目錄時長， 一般為5a或10a。

2.4 平滑半徑

ASSM與SSM的另外一大區(qū)別在于對平滑半徑的取值方式不同。傳統(tǒng)方法中， 同一震級檔的地震記錄都采用相同的平滑半徑(Frankel， 1995)。而ASSM中每個輸入地震i的平滑半徑di均取決于其周邊地震的密度， 衡量標(biāo)準(zhǔn)是地震i到第n個最近地震的水平距離(Helmstetteretal.， 2007； Werneretal.， 2011)。因此， 當(dāng)選定任意n時， 地震密集區(qū)的平滑半徑較小， 地震稀疏地區(qū)的平滑半徑較大。

然而， 如前文所述， 本研究中不同地震帶和震級檔的地震都可能對應(yīng)不同的完整記錄時長， 因此需要對求取平滑半徑的算法做相應(yīng)調(diào)整： 令地震i的平滑半徑bi等于以地震i的震中為圓心的圓的半徑， 要求該圓內(nèi)的所有地震的年發(fā)生率的和等于設(shè)定值RS， 本文將RS稱為周邊地震年發(fā)生率。計算時由小到大嘗試不同的RS值， 依次計算出模型結(jié)果并評價其優(yōu)劣， 最終選取令模型表現(xiàn)最優(yōu)時的RS作為最終值。確定RS后， 即可相應(yīng)地計算出每個輸入模型的地震i的平滑半徑bi， 其平均值為

(9)

式中，Ni為輸入地震總數(shù)。

2.5 模型優(yōu)化和概率增益函數(shù)

采用前人推薦的概率增益函數(shù)對模型的優(yōu)劣予以評價(Kaganetal.， 1977； Helmstetteretal.， 2007； Werneretal.， 2010)。假設(shè)每個網(wǎng)格內(nèi)地震的發(fā)生過程服從泊松分布， 均值為μ′(ix，ii)(式(8))， 那么在檢測目錄相同的時間段內(nèi)網(wǎng)格(ix，iy)中發(fā)生n個地震的概率為

(10)

假設(shè)研究區(qū)內(nèi)所有網(wǎng)格都恰好發(fā)生檢測目錄中該網(wǎng)格記錄到的地震數(shù)， 對其發(fā)生概率取對數(shù)似然函數(shù)：

(11)

引入一個空間均一模型(Lunif)作為參考， 在該模型中發(fā)震概率在研究區(qū)內(nèi)均勻分布， 類似于美國國家地震區(qū)劃圖中的背景地震模型(Frankel， 1995)。將檢測目錄中震級大于目標(biāo)震級Mt的地震數(shù)Nt除以研究區(qū)的網(wǎng)格數(shù)Ntg即可得到空間均一模型網(wǎng)格的地震發(fā)生率。

最終， 模型的優(yōu)劣可通過概率增益函數(shù)G體現(xiàn)，G越大則對應(yīng)模型表現(xiàn)越優(yōu)， 選取使G達(dá)到最大時的周邊地震年發(fā)生率RS作為最終的參數(shù)設(shè)置。G的形式為

(12)

3 模型計算結(jié)果及分析

為了檢測不同參數(shù)設(shè)置對模型的影響， 基于上述原理建立了不同的模型， 模型參數(shù)和檢測結(jié)果見表3， 表中展示了該參數(shù)設(shè)置下概率增益函數(shù)G最大時的結(jié)果。選取2013—2017年共5a的目錄作為檢測目錄， 目標(biāo)震級Mt定為M5.0， 目標(biāo)震級以上地震數(shù)Nt共87條。學(xué)習(xí)目錄中的地震最早可追溯到1400年， 對應(yīng)朝鮮地震帶歷史地震記錄的完整記錄起始時間(表2)。學(xué)習(xí)目錄的截止時間為2012年， 與檢測目錄之間互不重疊。

表 3 不同參數(shù)設(shè)置下模型預(yù)測結(jié)果的對比Table3 Forecast results for models with different parameter settings

3.1 輸入地震最小震級

模型Model-1、 Model-2、 Model-3和Model-4對應(yīng)的學(xué)習(xí)目錄輸入地震最小震級Mmin分別為M2.0、M3.0、M4.0和M5.0(表3)。Mmin越小， 輸入地震數(shù)Ni越多， 隨著Mmin的減小Ni從915條增長為19378條。其中表現(xiàn)最優(yōu)的是Model-3， 對應(yīng)的G值為4.48。事實上，G值與周邊地震年發(fā)生率RS直接相關(guān)。上述模型的G值隨RS變化的曲線如圖 5 所示。4條曲線的總體趨勢均為： 隨著RS增大，G值較為快速地增大， 隨后G值增加速率逐漸減小直至達(dá)到最大值， 然后G值開始緩慢下降。

圖 5 不同模型的概率增益函數(shù)G隨周邊地震年發(fā)生率RS的變化曲線Fig. 5 The curves of the probability gain per earthquake(G)versus the annual occurrence rate of surrounding earthquakes(RS)for different models.

從圖 5 可以看出， Model-3的G值始終占優(yōu)， Model-1和Model-2居中， Model-4的最差。這說明當(dāng)模型的Mmin從5.0降低到4.0， 模型的表現(xiàn)顯著優(yōu)化， 這與輸入地震數(shù)從915提高到2264直接相關(guān)。

然而Model-1和Model-2的表現(xiàn)不如Model-3， 這說明對于本研究區(qū)而言Mmin并不是越小越好。分析其原因可知： 琉球海溝地震帶和臺灣南-馬尼拉海溝地震帶的最小完整震級均為M4.0(表1)， 因此當(dāng)Mmin為2.0或者3.0時， 這2個地震帶并沒有增加輸入地震數(shù)， 反而受到其他地震帶小地震過度平滑的影響。同時， 這2個地震帶的地震最為密集(圖 1， 2a)， 對結(jié)果造成影響的權(quán)重最大。因此， 建立模型時需要綜合考慮研究區(qū)內(nèi)的地震分布和各個地震帶最小完整震級的情況。

3.2 ASSM與SSM的效果對比

為了對比自適應(yīng)空間平滑模型(ASSM)和傳統(tǒng)的固定平滑半徑的空間平滑模型(SSM)的表現(xiàn)， 建立了4個SSM模型Model-5、 Model-6、 Model-7和Model-8(表3)， 分別與Model-1、 Model-2、 Model-3和Model-4相對應(yīng)， 區(qū)別在于選用了固定的平滑半徑。參考美國國家地震區(qū)劃圖的做法(Petersenetal.， 2015)， SSM中對M4.0以下輸入地震的平滑半徑取50km，M4.0以上輸入地震的平滑半徑取75km。2類模型的概率增益函數(shù)G的對比見圖 6， 圖中輸入地震最小震級Mmin逐步從2.0上升到5.0， ASSM的G值始終大于SSM， 這說明ASSM相比SSM具有穩(wěn)定的優(yōu)勢。

圖 6 ASSM和SSM的概率增益函數(shù)G隨輸入地震最小震級Mmin 變化的曲線Fig. 6 The curves of the probability gain per earthquake(G)of ASSM or SSM versus the minimum input magnitude(Mmin)of the learning catalog.

3.3 最終模型

經(jīng)過上述對比分析， 選取Model-3作為本研究的最終模型， 該模型采用M4.0以上的地震數(shù)據(jù)建立模型， 對應(yīng)的周邊地震年發(fā)生率RS為0.057， 所有輸入地震的平均平滑半徑bw為92.3km。Model-3預(yù)測的M≥5.0地震的年發(fā)生率分布見圖 7。

圖 7 M5.0以上地震的年發(fā)生率預(yù)測結(jié)果(表3中的Model-3)Fig. 7 Forecasted seismicity rates of M≥5.0 annually in each cell from Model-3 in Table 3.

觀察圖 7 和圖2a可發(fā)現(xiàn)， 圖 7 中地震年發(fā)生率的高低變化與圖2a中地震分布的密集程度一致， 同時地震年發(fā)生率由高到低平滑過渡。年發(fā)生率最高的地區(qū)位于琉球海溝、 中國臺灣、 馬尼拉海溝一帶， 這源于該地區(qū)的俯沖帶背景； 其次是華北平原地震帶的北段， 其中一處顯著高于周邊地區(qū)(圖 7 中的A區(qū))， 這是由于1976年唐山地震后有大量余震出現(xiàn)(萬永革等， 2008； 蔣長勝等， 2013)， Gardner等(1974)提出的除叢算法無法完全有效識別； 渤海、 黃海、 朝鮮半島的地震活動水平與廣東、 福建沿海一帶接近， 浙江沿海、 東海和南海的地震活動水平最低。

3.4 模型的優(yōu)點和局限性

本文的ASSM模型的主要優(yōu)點有：

(1)本模型完全基于歷史和儀器地震記錄數(shù)據(jù)， 可以非常方便快捷地建立模型， 并且ASSM的地震預(yù)測效果已得到了前人研究的證實(Straderetal.， 2017； Schorlemmeretal.， 2018)。

(2)本模型在利用地震數(shù)據(jù)時具有較強(qiáng)的適用性和靈活性： 通過對算法進(jìn)行調(diào)整， 完整記錄水平隨時間和空間變化的地震數(shù)據(jù)都可以被本文的模型充分利用。當(dāng)遇到類似于朝鮮地震帶的情況， 即同一個震級檔有2個以上完整記錄時段(圖2b， 表1)， 在式(4)中扣除記錄水平較差的時段的時長即可， 這在程序中很容易實現(xiàn)。

(3)本文的ASSM與傳統(tǒng)的SSM相比， 預(yù)測效果具有穩(wěn)定的優(yōu)勢。最新的美國國家地震區(qū)劃圖也從單獨使用SSM改為同時使用ASSM和SSM(Petersenetal.， 2015， 2020)。

本文的ASSM模型的不足和局限性有：

(1)本模型只能給出目標(biāo)震級Mt以上的地震年發(fā)生率， 不能給出震級上限， 而震級上限是地震活動和地震危險性分析中的重要參數(shù)。

(2)本模型基于歷史和儀器地震記錄數(shù)據(jù)， 不能體現(xiàn)活動構(gòu)造、 地球物理、 大地測量等學(xué)科的研究成果， 例如斷層的空間分布、 滑動速率、 閉鎖程度等(程佳等， 2020)。近年來CSEP的研究發(fā)現(xiàn)， 同時考慮地震目錄和活動斷層等多種資料和手段的混合模型的預(yù)測效果要優(yōu)于單一的預(yù)測模型(Rhoadesetal.， 2018； Taronietal.， 2018)， 因此建議在進(jìn)行地震危險性分析時將本模型與其他模型相結(jié)合使用。

(3)大地震的能量積累和釋放具有一定的時間相關(guān)性(Matthewsetal.， 2002； Lietal.， 2017)， 而本文的ASSM模型假設(shè)網(wǎng)格內(nèi)地震的發(fā)生服從時間泊松分布， 屬于時間獨立的地震活動模型。

4 結(jié)論

通過建立中國海域及鄰區(qū)自適應(yīng)空間平滑地震活動模型， 并比較不同參數(shù)設(shè)置時模型的表現(xiàn)， 得到以下幾點結(jié)論：

(1)通過算法的調(diào)整， 本文的模型可以充分利用完整記錄水平隨時間和空間變化的地震數(shù)據(jù)， 顯示出較強(qiáng)的適用性和靈活性。

(2)當(dāng)輸入地震最小震級設(shè)為M4.0時模型的表現(xiàn)最優(yōu)。這提示建立模型時并非輸入地震最小震級越小、 地震數(shù)越多， 模型的表現(xiàn)就越好， 而是需要綜合考慮研究區(qū)內(nèi)地震的分布情況和各個地震帶的完整記錄水平。

(3)本模型與傳統(tǒng)的固定平滑半徑的平滑模型相比具有更大的概率增益函數(shù)值， 這種優(yōu)勢不受模型的輸入地震震級的影響。

(4)本文的模型存在不能計算震級上限和反映斷層分布、 未能考慮大震復(fù)發(fā)的時間相關(guān)性等不足和局限性。因此， 本模型可以單獨用于描述中小—中強(qiáng)地震的發(fā)震概率， 也可作為確定潛在震源區(qū)的空間分布函數(shù)的重要因子之一。但在進(jìn)行地震危險性分析時， 建議結(jié)合地震地質(zhì)、 GPS等其他學(xué)科的成果形成混合模型， 從而進(jìn)一步提高模型的適用性和有效性。本文的工作將為中國海域地震危險性分析和建模提供依據(jù)。