推定最大余震震级法在南北地震带的应用*

2022-06-23解孟雨孟令媛

地震研究 2022年3期

解孟雨，孟令媛

(中国地震台网中心，北京 100045)

0 引言

我国的国土面积虽仅占全球陆地面积的7%，却有全球33%的大陆强震发生在我国(中国地震信息网，2018)。据中国地震信息网(2018)统计，1901—2001年我国大陆因地震死亡人数达到60余万人，最惨烈的是1920年12月宁夏海原8.6级地震和1976年7月河北唐山7.8级地震，死亡人数分别是23.55万人和24.2万人。20世纪以来，科技得到了迅速的发展，但是从发震频度上看，1950年以来中国大陆地区平均每年发生5级地震24次、6级地震4次、7级地震0.6次，地震仍然是人民生命安全和社会发展的巨大威胁。

南北地震带作为中国大陆中部强震相对密集分布的地震带，发生了2008年四川汶川8.0、2010年青海玉树7.1、2013年四川芦山7.0、2017年四川九寨沟7.0、2014年云南景谷6.6、2014年云南鲁甸6.5等一系列强震，这些强震对震中及其附近城市和乡镇均造成了巨大的伤害。例如截至2008年9月25日，2008年汶川8.0地震共造成直接经济损失8 523.09亿元，致使69 227人死亡，17 923人失踪以及374 643人受伤(陶正如，陶夏新，2021)。研究显示南北地震带仍存在发生强震的可能(M7专项工作组，2012)。

强震发生后，快速研判该地震序列的类型及余震强度具有重要的减灾意义。现阶段我国震后应急及地震序列跟踪工作中，针对地震序列类型及余震强度的判定方法，以震中附近历史地震序列对比分析为主，并综合考虑统计规律(中国地震局监测预报司，2007；吕晓健等，2010；蒋飞蕊等，2015；赵小艳等，2021)、值(刘正荣，孔昭麟，1986；马茹莹等，2016)以及推定最大余震震级方法(吴开统等，1984；Shcherbakov，Turcotte，2004)等的计算结果。其中推定最大余震震级方法可以依据地震序列数据实时评估最大余震震级，且计算简单快速，是震后评估最大余震震级水平的常用方法之一。Shcherbakov等(2013)利用该方法估算了全球70个俯冲带地震序列的最大余震(或最大后续地震)震级，结果显示对于大多数序列可以准确估计出最大余震震级，对于震级差为0.2级的震群序列也可以有效估计出最大后续地震震级。苏有锦等(2014)针对全球184个7级以上浅源主-余型地震序列研究了震后不同时间下推定最大余震震级方法的适用性，结果显示估算出的最大余震震级与实际值的差值多为-0.5～0.5，且可以利用震后短时间尺度的余震数据估算出合理的最大余震震级。解孟雨等(2017)基于九寨沟7.0级地震序列的数据，讨论了震级类型、余震目录时间长度对于推定最大余震震级方法的影响，结果显示该方法对面波震级和里氏震级均适用，且余震目录时间尺度较长时，最大余震震级估算值较好，利用震后短期的计算结果也可以作为最大余震震级的合适估计。

为进一步分析推定最大余震震级方法对我国大陆地区强震序列的有效性，分析不同常用最小完备震级和值计算方法组合下的估算结果，本文针对南北地震带开展相关研究，选取该区2000—2020年≥6.0地震序列数据，讨论分析不同方法组合以及数据时间长度对实际最大余震震级的影响，然后依据分析结果提出了最优的推定最大余震震级方法。

1 计算方法

对于现有的地震序列记录，一般满足G-R定律，即：

lg(≥)=-

(1)

式中：(≥)为震级大于的余震个数；和则是正常数，描述了选定区域地震活动性水平的高低程度，它依赖于选定区域和观测数据时间尺度的大小，一般认为值给出了区域总的地震个数，而值刻画了区域内地震大小的相对分布情况，其范围为0.6～1.4(Marzocchi，Sandri，2003；Shcherbakov，Turcotte，2004；Wang，2015；Hamdache，2017)。=1时，地震序列的推定最大余震震级为：

(2)

从式(2)可以看出，只需在地震序列最小完备性震级的基础上，计算出参数和，便可估算序列最大余震震级，这一方法就称为推定最大余震震级法。

最小完备性震级的估算方法主要分为两类(Mignan，Woessner，2012)：基于台站的方法和基于目录的方法。基于目录的方法相对简单明确，适合实现快速计算，其中最大曲率法(MAXC)和拟合优度测试法(GFT)计算较为简单快捷，满足震后最大余震震级快速评估的要求。

MAXC将地震个数(+d≥≥)最大时对应的震级值作为地震目录的，其中为地震震级记录的最小间隔，一般为0.1。MAXC会低估地震目录的，因此更为接近实际真值，一般会添加修正系数0.2，即=(MAXC)+ 0.2(Woessner，Wiemer，2005)，该计算方法称为修正MAXC。

GFT首先定义了衡量拟合G-R定律曲线与实际数据拟合程度的参数值：

(3)

式中：为选定的震级；、为针对≥地震目录计算出的参数值；和分别为每一震级档下观测到的累积地震个数和按拟合G-R定律计算出的累积地震个数。通过改变大小，计算出相应的值，当值大于设定阈值时，取为地震目录的。通常值阈值设定为90%或95%，相应的一般记为-90和-95。使用GFT估算时，存在无法计算出-90和-95的情况，因此在实际使用过程中常采用GFT和修正MAXC的混合方法，依次确定是否可以计算出-95和-90，如果都无法计算出，则利用修正MAXC计算相应的(Wiemer，2001)。

值的计算方法主要分为两大类：最小二乘法和最大似然法。对于最小二乘法(LS)，本文采用常用的线性回归系数最小二乘估计(贾俊平等，2018)，即设定变量和之间满足：

=++

(4)

式中：和为需要确定的常数；为误差项，反映了除和之间线性关系的其他因素的影响。由实测数据和可以得到和的估计值分别为：

(5)

使用最大似然法(ML)估算值，常用的值及不确定度计算如下：

(6)

(7)

=lg [(≥)]+

(8)

(9)

在采用ML计算值时，≥的地震总数(≥)必须足够大，一般为≥30(Aki，1965；贾俊平等，2018)。

为了寻找最优的组合方法，本文采用上述常用方法的组合计算出推定最大余震震级，然后通过对比分析对实际最大余震震级的估算效果，确定更为有效的方法组合。具体而言，分别采用修正MAXC方法以及GFT和修正MAXC的混合方法(以下分别简称为MAXC和GFT)计算，然后分别使用LS和ML计算值，构成4种最大余震震级估算方法组合，即GFT+LS、GFT+ML、MAXC+LS和MAXC+ML。

(10)

式中：为地震序列总个数；为误差的标准差。

2 数据和结果

南北地震带是一条由不同方向、不同性质断裂和褶皱构成的近南北向的复杂构造带和地震活动带，是中国大陆地震活动最为强烈的地区之一。依据地质构造特征和地震活动的空间分布特征，可以给出南北地震带不同的空间范围定义。为收集足够的地震序列数据，同时也为大震活跃区地震序列最大余震的评估提供可行的手段，本文采用以中国大陆强震的空间分布特征定义的南北地震带范围(邵志刚，张浪平，2013)，即(21°～45°N，95°～110°E)。在该研究范围内，2000—2020年中国地震台网目录中共记录到6级及以上地震35次，其中6.0～6.9地震31次，≥7.0地震4次，去掉余震后共有25个≥6.0地震(表1、图1)。由于2000年1月15日姚安6.5级、2000年9月12日兴海6.6级、2001年10月27日永胜6.0级以及2008年11月10日海西6.3级地震序列的现有目录存在记录时间较短、地震个数较少或缺失小震的情况，无法进行推定最大余震震级的计算，因此仅对其余21个地震序列进行研究。

图1 南北地震带MS≥6.0地震震中分布Fig.1 Locations of MS≥6.0 earthquakes in the North-South Seismic Belt

表1 2000—2020年南北地震带MS≥6.0主震基本信息Tab.1 Mainshocks with MS≥6.0 in the North-South Seismic Belt from 2000 to 2020

蒋海昆等(2007)指出绝大多数5.0级以上地震序列的最大余震发生生在震后1a内，其中93%的最大余震发在震后90 d内。在日常震后趋势研判中，我们更关注震后90 d内的最大余震震级。因此针对南北地震带21个≥6.0地震序列，本文分别使用4种和值的方法组合，分别估算主震后1、3、5、7、15、30、90 d时的最大余震震级，然后与震后90 d内实际发生的最大余震震级进行对比。

表2 4种组合方法得出的MS≥6.0地震序列的MC和MILATab.2 MC and MILA for MS≥6.0 earthquake sequences by using 4 combinations

图2 4种组合方法得出的和CM(b)Fig.2 from 4 combinations

表3进一步给出了使用GFT+ML组合方法得出的21个≥6.0地震序列在震后不同时刻的和最大余震震级的估计值。2013年芦山7.0、2019年长宁6.0等13个地震序列的随时间增加而逐步减小(表3中加粗的数字)，符合通常的预期，而其他序列的则随时间增加基本不变或者呈现小幅波动变化。2010—2020年多数≥6.0地震序列<2.0，这可能是由于南北地震带台网的升级和布设使得出现了整体的下降。震后不同时刻的整体上与余震数据时间长度不存在明显的比例关系。如果定义误差绝对值|Δ|≤0.5为方法估计结果正确，则有约71.4%(15/21)的地震序列可以在震后90 d时获得正确的估计结果，同时有57%(12/21)的地震序列可以在震后1 d内给出正确的估计值，以上表明用GFT+ML组合估算最大余震震级有较好的正确率，利用震后短期的数据进行最大余震震级的估计，也具有一定的正确率。

3 讨论

图3 使用GFT+ML组合方法计算景谷MS6.6地震后1 d(a)和90 d(b)余震序列的G-R定律拟合结果Fig.3 Fitting results of G-R law obtained by GFT and ML 1 in day(a) and 90 days(b)after the Jinggu MS6.6 earthquake

表3 使用GFT+ML得出的21个MS≥6.0地震序列的MC和MILATab.3 Results of MC and MILA for 21 MS≥6.0 earthquake sequences by using GFT and ML

在使用LS方法计算值时，较大震级地震的权重较大(吴果，2018)，因此在计算时，若不限制余震数据的震级上限，则相应的G-R定律拟合曲线会存在偏离实测累积频度-震级分布的情况。如利用LS方法对雅江6.0地震后1 d的余震数据进行拟合时，如不限制震级上限，G-R定律拟合曲线存在明显的偏离(图4a)，此时值以及外推最大余震震级均不可信。而按照本文LS方法的计算方法，即自动选取拟合G-R定律和实测累积频度-震级分布误差最小时的震级上限，则可以一定程度上避免拟合曲线偏离的问题，提高拟合结果的可靠性(图4b)。

图4 采用GFT+LS组合方法计算雅江MS6.0地震后1 d不同震级上限的拟合结果Fig.4 Fitting results of G-R law obtained by GFT and LS with different upper limits of magnitude in 1 day after the Yajiang MS6.0 earthquake

整体上使用LS方法进行最大余震震级估计的结果比使用ML时的大。以景谷地震为例，震后1 d使用LS方法计算出的最大余震震级为8.3，远大于实际最大余震震级6.2(表2和图5)。对比图5和图3a可以看出，尽管采用了寻找最小误差的方法确定震级上限，部分改善了LS方法计算时较大震级地震权重较大的问题，但当小震级地震缺失时，仅考虑较小震级范围余震的拟合G-R定律会偏离考虑全部≥余震时拟合出的G-R定律。Sandri和Marzocchi(2007)通过理论分析和数值计算讨论了使用LS方法估算的值情况，结果显示理论上LS方法得到的值会明显低于给定值。因此为了震后快速准确地估算余震序列最大地震震级，ML方法是更合适计算值的方法。

图5 使用MAXC+LS组合方法得到的景谷MS6.6地震后1 d的G-R定律拟合结果Fig.5 Fitting results of G-R law obtained by MAXC and LS combination in 1 day after the Jinggu MS6.6 earthquake