青海瑪多MS7.4 地震前b 值時空變化特征1

杜 航 楊 云 鄭江蓉 王 俊 張 揚 宮 杰

（江蘇省地震局, 南京 210014）

引言

北京時間2021 年5 月22 日02：04，青海果洛州瑪多縣發生MS7.4 地震（34.59°N，98.34°E），震源深度17 km。震中位于瑪多黃河鄉，距縣城38 km，距西寧380 km，震中海拔高度4 258 m。地震表現出無明顯前震、余震較多的特征，主震后1 小時內發生MS3.0 以上余震10 次，最大余震為MS4.3，一天內發生MS3.0 以上余震36 次，最大余震為MS5.1，為主震-余震型地震。中國地震臺網中心對震源機制解進行反演得出本次地震為走滑型事件，地震破裂持續45 s，震后30 s 內釋放大量能量，最大滑移量4.5 m，地表破裂超過160 km（鄧文澤等，2021）。本次地震為汶川8.0 級地震以來，中國境內發生的震級最高、地表破裂最長的地震事件。本文將對此次地震震前b值進行研究，尋求該區域b值與應力的變化關系，為該區域后續地震預報提供科學依據和參考。

Gutenberg 等（1944）發現震級M與地震頻次N之間存在統計學關系logN=a?bM，a表示研究區域地震活動水平，b表示研究區域不同震級的相對分布。Kun 等（2013）通過數值模擬研究驗證了震級、頻次的數學關系。研究人員在進行巖石實驗時發現介質的應力水平不斷增加時，b值會呈現逐漸減小的趨勢（Mogi，1962；Scholz，1968；Rivière 等，2018）。Mousavi 等（2017）對整個板塊進行研究，同樣發現b值的變化與介質的應力變化存在一定關系。因此，研究學者使用b值研究區域應力的大小（Nuannin 等，2012；Scholz，2015）。易桂喜等（2006）對龍門山斷裂帶b值進行分析，得到斷裂帶區域的強震危險性；劉艷輝等（2015）對青藏高原東南緣b值進行時空掃描，發現在大地震發生前后b值會出現較大變化，地震發生前b值長時間處于低值；Tormann 等（2015）研究得到全球平均b值約為1，中強震多發地區b值小于1，小地震多發地區b 值大于1； Nanjo 等（2018）利用b值與應力關系，得到海槽地震的初始破裂部位；謝卓娟等（2019）利用b值與應力關系研究首都圈地震活動性；曾憲偉等（2021）通過b值分析得到寧夏吳忠?靈武地區的強震危險性。前人已開展大量區域應力變化與b值之間的關系研究，在此基礎上Nievas 等（2020）提出可以利用b值變化研究區域應力變化，尋求b值變化與地震發生的關系。本文對瑪多MS7.4 地震震前地震序列進行分析，研究地震前b值的變化特征，為研究區域地震預測和地震發生風險提供參考。

1 區域地質背景

新生代以來，青藏高原不斷受到歐亞板塊和印度板塊的擠壓，使青藏高原呈現擴張趨勢，表現出板塊運動活躍、地震多發的特點，21 世紀以來發生7 級以上地震10 次（Molnar 等，2009）。巴顏喀拉塊體是青藏高原的次一級構造單元，該塊體被東昆侖、甘孜-玉樹-鮮水河、龍門山3 條大型斷裂帶包圍，成為中國目前地震活動最為強烈的地區之一（張培震等，2003；Rui 等，2016；聞學澤，2018）。其中，巴顏喀拉塊體東緣最為活躍，十年內相繼發生了2008 年汶川8.0 級地震、2013 年蘆山7.0 級地震、2017 年九寨溝7.0 級地震（圖1）。

圖1 青藏高原主要活動斷裂及歷史地震分布Fig. 1 Simplified distribution of main active faults and historic earthquakes in the Qinghai-Tibet plateau

東昆侖斷裂位于巴顏喀拉塊體北邊界，EW 走向，長約2 000 km，是一條巨型左旋走滑斷裂帶，傾角55°～85°。該斷裂自西向東走滑速率不斷降低，瑪沁段為12.5 mm/a，至瑪曲段降為5 mm/a，至塔藏段已不足3 mm/a（Kirby 等，2007；Ren 等，2013）。有記錄以來，該斷裂發生7 級以上地震6 次，最大震級為2001 年MS8.1 地震，地表破裂達426 km（徐錫偉等，2002；鄧起東等，2003；Xu 等，2006）。

甘孜-玉樹-鮮水河斷裂帶位于巴顏喀拉塊體南邊界，為左旋走滑斷裂帶。學者對該斷裂帶滑動速率有不同見解，聞學澤、石峰等認為甘孜-玉樹段滑動速率為8～14 mm/a，但 Chevalier 等認為甘孜-玉樹段滑動速率僅有5～8 mm/a（聞學澤等，2003；石峰等，2013；Chevalier 等，2018）。學者們普遍認為鮮水河段滑動速率為10～13 mm/a（Zhang，2013；Chen 等，2016）。有記錄以來，該斷裂發生7 級以上地震17 次，最大震級為1997 年MS7.9 地震。

龍門山斷裂位于巴顏喀拉塊體東邊界，為逆沖右旋走滑斷裂帶，由多個近平行斷裂組成。21 世紀以來，該斷裂活動加劇，2008 年汶川MS8.0 地震、2013 年蘆山MS7.0 地震都發生在該斷裂上。

巴顏喀拉塊體除四周存在3 條大型斷裂帶外，其內部也發育著一系列活動斷裂。瑪多-甘德斷裂帶位于東昆侖斷裂帶以南70 km，為全新世活動的左旋走滑斷裂，整體呈NW-SE 走向，長約650 km，水平滑動速率小于6～8 mm/a，活動速率明顯小于位于其南北兩邊的甘孜-玉樹-鮮水河斷裂帶和東昆侖斷裂帶（熊仁偉等，2010）。

2 數據及方法

2.1 地震目錄數據

本次研究地震目錄由中國地震臺網提供，目錄選取時間為2009 年1 月1 日至2021 年4 月30 日。青海瑪多MS7.4 地震研究范圍選取震中所在的青藏高原東北部巴顏喀拉塊體內部，南北范圍選取東昆侖斷裂帶、甘孜-玉樹-鮮水河斷裂帶之間區域（32°N～37°N），東西范圍覆蓋瑪多-甘德斷裂帶且包含本次地震所有余震（96°E～100°E）。

由于巴顏喀拉塊體內部發育著眾多活動斷裂，導致區域內5 級以上地震較為分散。其中瑪多MS7.4 地震及其余震沿瑪多-甘德斷裂帶分布，呈NW-SE 走向，空間分布長約150 km（詹艷等，2021；張喆等，2021；潘家偉等，2021）（圖2）。

圖2 研究區域中強震（MS≥5.0）及瑪多地震余震（ML≥3.0）分布圖Fig. 2 Simplified distribution of major earthquake（MS≥5.0）and Maduo earthquake aftershocks（ML≥3.0）in the study area

在b值計算時，MS震級適用于大范圍、長時間、高震級（一般指一級板塊以上范圍、100 年以上時間、6 級以上震級）范圍內b值計算（陳培善等，2003）。本次研究所選區域較小，時間尺度為12 年，地震目錄99.8%震級小于5 級，因此本次研究采取ML震級進行計算。

地震活動性研究中，余震是否應該剔除尚未有統一答案。眾多學者認為中強地震引發的余震會造成研究區域地震頻次變化，因此在b值變化研究時，應對地震目錄內的余震進行剔除（Keilis-Borok 等，1980；Reasenberg，1985；陳凌等，1998；韓曉明等，2016）。也有學者認為余震去除會影響地震目錄完整性，在研究時應保持地震目錄的原始性，選取最小完整性震級以上地震進行計算即可（史海霞等，2018）。為研究余震序列對b值計算結果的影響，本文采用完整地震目錄和去余震目錄 2 種地震目錄研究b值的時空變化，為準確對比2 種地震目錄得到的b值差異，除地震目錄外其它參數均相同，其中余震剔除時采用G-C 方法，即對5 級以上地震余震時間窗915 天、余震空間窗半徑R0內的余震進行剔除（陳凌等，1998）：

2.2 地震b 值計算方法

目前b值計算方法主要有2 種，最小二乘法和最大似然法。最小二乘法是直接通過震級-頻次的關系確定b值的大小，該方法簡單易行，但其穩定性受地震頻次影響，當地震較少時結果容易出現偏差（Weichert，1980；Sandri 等，2007）。最大似然法作為概率方法具有更好的穩定性，可降低少數地震對區域整體b值結果的影響。由于研究區域臺站間距較大，導致區域內最小完備震級較大，地震頻次較少，因此，采用最大似然法計算本區域b值，該方法最早由Aki（1965）提出：

式中，M為震級；Mc為最小完備震級； ?M為震級分檔，本次地震目錄震級分檔為0.1。b標準誤差可由以下公式進行估算（Aki，1965）：

式中，N為樣本量，即地震目錄中大于等于Mc的 所有地震的個數。可以看出N與δb成反比，因此提高樣本數量可降低b值誤差。該方法先將研究區域網格化，然后計算每一個網格的b值。計算b值時需要震級為連續隨機變量，而實際震級計算時，會對震級進行歸檔，使震級只保留一位小數。因此震級分檔?M越小，b值的計算結果越準確，本次地震目錄?M為0.1。

2.3 最小完備震級

最小完備震級通過最大曲率法確定，該方法認為震級-頻次曲線的一階導數最大值對應的震級即最小完備震級（Wiemer 等，2000），研究區域最小完備震級為ML2.1（圖3）。最小完備震級反映該地區的地震監測水平，青海地區雖地震頻發，但人口稀少、臺站密度較小，因此最小完備震級較大。為避免因地震監測能力差異導致地震目錄不完整，影響b值的時空變化研究，本次研究將選取ML2.1 以上地震對瑪多7.4 級地震震源區進行b值研究。

圖3 研究區域2009 年1 月1 日?2021 年4 月30 日震級-頻次圖Fig. 3 The frequency-magnitude distribution of earthquakes in the study area during January, 2009?April, 2021

3 地震b 值時空掃描特征

3.1 b 值時間掃描特征

為研究時間段區域內b值隨時間的變化規律，選用ML2.1 以上地震目錄進行計算，為避免地震數量較少對結果準確性造成影響，以4 個月為時間段進行b值時間掃描，分別計算余震剔除前后b值隨時間的變化（圖4）。由圖4 可知，2 種地震目錄得到的b值整體變化趨勢相同，在5 級地震發生前兩者吻合的較好，但在5 級以上地震發生后的一段時間內存在明顯差異。圖4 中2 條水平直線分別表示2 種地震目錄b值均值，其中完整目錄b值均值為0.86、去余震目錄b值均值為0.85，均值相差0.01。

圖4 b 值時間曲線Fig. 4 The temporal variation of b value

從b值時間掃描特征來看，是否剔除余震對b值時間變化影響不明顯，僅在中強震發生后一段時間內存在差異。這種現象是由于研究區域中強震多表現為主震-余震型地震序列類型，余震主要集中在主震發生后，因此是否剔除余震對中強震發生后一段時間較為明顯，對其他時間段影響較小。同時，雖然中強震余震持續時間較長，但隨時間的推移余震不斷減小，而該區域最小完備震級較大，不管是否剔除余震，震級小于ML2.1 的地震都被刪除，后續較小余震對2 種地震目錄未產生影響，中強震發生后一段時間后2 種地震目錄的b值變化曲線又互相重合。因此，是否剔除余震對中強震余震窗內b值時間變化存在一定影響。

瑪多MS7.4 地震前1 年b值持續低于均值，并處于下降趨勢，b值有上升趨勢時地震發生，擴大到區域內其它5 級以上地震，也發現此規律。同時，2 種地震目錄的b值曲線均存在3 段較長的低b值時間段，且起止時間差異不大，分別為2009 年5 月?2011 年12 月、2014 年5 月?2016 年4 月、2019 年9 月至今，基本為中強震分割點，地震發生后b值明顯上升，短時間內又下降到較低位置，并一直處于較低位置直至下次地震發生。這是由于地震發生前，該區域處于應力積累階段，b值處于較低位置，地震發生后應力得以釋放，區域內b值恢復到正常水平。

3.2 地震b 值空間掃描特征

b值時間掃描結果表明，b值變化與中強震發生存在一定關系，可以作為地震預測時間上的參考。但研究區域較大，時間掃描特征僅能反映研究區域時間尺度上發生地震的可能，無法準確給出地震發生的空間區域。要實現通過b值分析給出研究區內中強震發生的空間預測范圍，就要對b值進行空間掃描，分析研究區域內b值處于異常狀態的地區。由b值時間掃描特征知，是否剔除余震對中強震余震窗內b值存在一定影響。為驗證余震剔除是否對b值空間掃描特征同樣具有影響，在b值空間掃描特征研究中同樣計算ML2.1 及以上震級去除余震地震目錄及完整地震目錄，分別計算2009 年1 月?2021 年4 月、2009 年1 月?2012 年12 月、2013 年1 月?2016 年12 月、2017 年1 月?2021 年4 月的b值空間特征，研究瑪多地震前該區域b值空間變化特征及b值空間分布與中強地震的關系。瑪多MS7.4 地震前去除余震目錄及完整目錄b值空間掃描結果分別如圖5、圖6 所示。

根據我國西部地震多發的特點，本文將研究區域經緯度掃描步長設為0.1°，并以此進行網格化（易桂喜等，2007）。由于研究區域地震較少，網格內地震掃描半徑設為70 km，為保證結果的可信度，僅選取地震數量超過30 個的網格進行計算，以每個網格點為單位計算b值，網格內的b值誤差最大為0.18b，小于地震b值的變化范圍，可信度較高。再根據每個網格點的b值做出研究區域b值空間分布，若網格點地震數量達不到要求則進行空白處理。

圖5（a）、圖6（a）分別為2009 年1 月1 日?2021 年4 月30 日去除余震目錄及完整目錄b值空間掃描結果，可以看出兩者在瑪多震中東南部低b值區域存在差別，圖5（a）中該區域b值約為0.7，圖6（a）中該區域b值約為0.8。區域內整體b值范圍為0.5～2.0，瑪多MS7.4 震中處于低b值區域內部，低b值區域以瑪多MS7.4 震中所在的瑪多-甘德斷裂為中心向兩側擴展，該區域b值整體低于0.6，瑪多MS7.4 震中所在區域為b值最低區域，b值約為0.5，表明該斷裂震前構造應力處于較高狀態。由圖5（a）、圖6（a）可知，b值的空間分布與中強地震的發生有較好的對應，為研究瑪多MS7.4 地震前震中區域b值變化過程，本文將地震目錄時間進行分割。

圖5（b）、圖6（b）分別為2009 年1 月1 日?2012 年12 月31 日去除余震目錄及完整目錄b值空間掃描結果，兩者在稱多西南部低b值區域存在差別，圖5（b）中該區域b值約為0.6，圖6（b）中該區域b值約為0.7。2010 年玉樹MS7.4 地震和2011 年囊謙MS5.2 地震均發生于該區域，兩者差異可能由是否剔除這2 次地震余震引起。圖5（b）、圖6（b）中瑪多MS7.4 地震震中位置空白，震中東部存在低b值區域，區域面積較小，b值在0.6 附近。

圖5（c）、圖6（c）為2013 年1 月1 日?2016 年12 月31 日去除余震目錄及完整目錄b值空間掃描結果，兩者沒有差異。

圖5 去余震目錄b 值空間分布圖像Fig. 5 Spatial distribution of b value with deleting aftershock

圖5（d）、圖6（d）為2017 年1 月1 日?2021 年4 月30 日去除余震目錄及完整目錄b值空間掃描結果，兩者在瑪多震中東南部存在差別，差異區域與圖5（a）、圖6（a）一致。造成該差異的主要原因為2020 年石渠MS5.9 地震發生于此區域，此次地震余震較多，去除余震會引起差異區域b值發生變化，因此該區域是否去除余震b值空間掃描結果呈現較為明顯的差異性。

圖6 完整目錄b 值空間分布圖像Fig. 6 Spatial distribution of b value without deleting aftershock

圖5（b）、（c）、（d）顯示隨著時間推移，瑪多震中附近低b值區域范圍逐漸擴大，由正東不斷向震中擴展。圖5（b）中低b值區域位于東昆侖斷裂帶與瑪多-甘德斷裂帶之間，圖5（d）中低b值區域遷移至瑪多-甘德斷裂帶，呈現南移現象。圖5（d）中低b值區域貫通震中東南區域，該區域b值逐年降低，因此存在發生強震的危險性。

對比2 種地震目錄在不同時間段的b值空間掃描結果，整體差異不大，但在低b值區域存在差異。該區域去余震目錄的b值相對較低，而完整地震目錄b值較高，這將掩蓋中強震區域的低b值特性。對比同一種地震目錄在不同時間段的b值空間掃描結果，可發現瑪多MS7.4 地震發生前低b值區域向震中的遷移現象，表明地震發生前震中附近應力不斷集中。

4 結論及討論

本文根據中國地震臺網提供的2009 年1 月1 日?2021 年4 月30 日地震目錄，對瑪多MS7.4 地震震源區震前地震活動性進行研究，分析區域內中強震發生前b值的變化特征。通過對完整地震目錄和去余震地震目錄分別進行b值時空掃描，發現是否剔除余震將影響中強震余震窗時間段內b值結果，空間掃描結果表明，在低b值區域完整地震目錄b值高于剔除余震地震目錄b值，從而掩蓋了中強震區域低b值特性，而去余震地震目錄可以較好凸顯低b值區域，對中強震預測的指示考意義更強。因此，在進行b值時空變化研究時，尤其是通過空間掃描結果研究中強震的位置時，應剔除余震。對于文中結果整體差異較小的現象，主要是由于該地區最小完備震級較大，大量余震都處于最小完備震級以下，導致2 種地震目錄差異較小，因此所得結果差異也較小。雖然2 種地震目錄差異較小，但在較為關心的低b值區域仍存在一定差異，這將掩蓋中強震區域的低b值特性。因此，后續研究可選則地震監測能力強、最小完備震級小的地區進行去余震地震目錄與完整地震目錄對比的研究，驗證是否隨著余震的增多，2 種地震目錄所得結果差異也隨之增大。具體結論如下：

（1）瑪多MS7.4 地震前1 年b值開始低于均值，且不斷下降，至b值有上升趨勢時地震發生，擴大至區域內其它5 級以上地震亦符合此規律，地震發生后b值明顯上升，短時間內又下降至較低位置，并一直處于較低位置直至下次地震發生。這是由于地震發生前，該區域處于應力積累狀態，b值處于較低位置，地震發生后應力被釋放，區域內的b值即恢復到正常水平。

（2）從b值空間掃描結果看，瑪多MS7.4 地震前震中位于低b值區域中心，該位置為b值最低處。

（3）不同時間段的b值空間掃描結果表明，瑪多MS7.4 地震發生前低b值區域向震中逐漸遷移，表明地震發生前應力向震中附近不斷集中，因此低b值區域的遷移現象可為地震預報提供參考。

（4）余震剔除可以有效避免中強地震余震序列對b值結果的干擾，b值空間掃描時完整地震目錄掩蓋了中強震區域的低b值特性，而去余震地震目錄可以較好凸顯低b值區域，對中強震預測的指示意義更強。因此，在b值時空掃描研究時，應剔除余震，力求更準確的表現b值時空變化特征。