河套地震带的b值时空变化特征分析

韩晓明 张文韬 王树波 赵星 包金哲 李娟

1）内蒙古自治区地震局，呼和浩特市哲里木路80号 010010

2）安徽省地震局，合肥 230000

0 引言

关于b值的研究源于古登堡和里克特提出的G-R关系（Gutenberg，et al，1954）

式中，参数a和b可根据一定时空范围内地震目录统计估算求取，b值反映了大、小地震的比例关系。更多的研究者通过大量的观测事实和实验证明了G-R关系式的成立，并指出了b值更多的特性和物理意义（Scholz，1968；李全林，1978；黄永祥，1982；耿乃光，1986；张智等，1987）。在此过程中，关于b值物理意义的争论也相当激烈。先是Scholz（1968）指出b值主要代表介质内部应力水平的高低，即随着介质应力水平的提高b值会减小；其后茂木青夫（1981）用声发射实验反驳了Scholz关于b值由应力状态决定的观点，并提出了b值由岩石性质决定的结论。再后Atkinson等（1984）根据双扭实验中的岩石声发射过程再次观测到了b值与应力状态的依存关系。此后的相关实验研究也都从各个角度阐述了影响b值的决定性因素，归结起来主要还是应力状态和介质性质二者为主。在利用b值判断区域应力场状态时，实验室内的不同样品类比于不同区域地质构造体间的差异（这句话比较费解，重新组织一下），因此，b值的区域差异性是应该存在的（李全林，1978；马鸿庆，1982；黄伟琼等，1989；王晓青等，1996）。根据上述观点，对同一区域，可根据b值时序变化判断应力场在不同时段内的变化情形；而在同一时段内，可针对不同区域比较b值的区域性差异。

近些年来，随着对b值研究的不断深入和应用范围的不断扩大，b值在地震预测研究方面的作用逐渐显现出来，研究者通过大量震例讨论了地震前震中附近的b值变化规律（李全林，1978；马鸿庆，1982；黄伟琼等，1989；焦远碧，1998；王晓青等，1996；吴忠良，2001；陈培善等，2003；易桂喜等，2004a、2004b；易桂喜等，2010）。可见，b值不仅是一个简单的统计分析参数，且具有明显的物理意义和广泛的应用领域。

本文以1970年以来发生在河套地震带的ML≥1.5地震事件作为研究对象，在最小完整性震级Mc评估的基础上进行b值求解，并对余震剔除和地震样本数等因素对b值的影响进行了定量分析。文中着重从b值的时序变化过程和空间分布情况进行分析。时序分析主要结合河套地震带6级以上中强地震活动实况，探讨了地震前后b值的变化规律；空间分析则根据河套地震带的断裂构造展布和活动情况，结合b值的空间分布进行区域应力场状态判定。

1 研究区的选取

b值的区域特征是由不同区域的介质构造状态差异引起的（尹祥础等，1987），选取某一地震带或构造单元进行b值分析，是消除b值区域差异性或者保证b值的背景值相对稳定和统一的有效途径，且进而能够较为真实和准确地提取出b值的时空变化信息。本文选取河套地震带作为研究区域，正是出于这方面的考虑。河套地震带既是一个地震事件比较集中、分布密度相对较高的地震带，又是一个地下介质构造状态比较单一和稳定的地质单元，可为b值计算提供比较丰富的地震样本数目和相对统一的应力积累背景水平。

在鄂尔多斯周缘众多断裂带中，地处北缘的河套断陷带规模最大，构造活动较为强烈。它位于长期隆起的内蒙地轴南侧的阴山隆起和鄂尔多斯隆起之间，西起狼山山前断裂，东至和林格尔断裂，南起鄂尔多斯北缘断裂，北至阴山山前断裂，长约440km，宽40～80km，总体走向近东西。经过渐新世、晚第三纪和第四纪3个演化时段复杂的构造运动，河套断陷带地层、地貌和构造标志显著，地表多为全新世沉积物覆盖层，地貌特征明显（国家地震局鄂尔多斯周缘活动断裂系课题组，1988）。河套地震带内的断裂活动控制着山地的抬升和盆地的沉降，其新构造运动以垂直差异运动为主，这与断裂两侧地块的拉伸及旋转运动有关。河套地震带边界断裂及其所控制的凸起、凹陷宏观布局均呈左行斜列，断陷带中最大沉降地段位于狼山山前和包头凸起东南侧，受NE向的断裂控制，断陷带两侧的阴山隆起和鄂尔多斯块体相对呈现为左旋拉张运动。河套地震带1970年以来MS≥4.7地震的震源机制结果显示（曹刚，2001），在NE向压应力和NW向张应力的作用下，除呼-包盆地东侧边缘位置出现的1次正断层活动震例外，其余16次中强地震的震源断层均表现出走滑为主的错动特征，这是河套地震带呈现出统一的水平作用为主的走滑型断层机制的重要证据（图1）。而在b值计算中，保证研究区受到统一的应力场控制是提高计算精度的有效途径（李全林，1978），这也是本文选择河套地震带作为b值研究区域的主要依据。

图1 1970年以来河套地震带M S≥4.7地震震源机制解空间分布

河套地震带是华北地震活动区北部的一条重要的地震活动带。进入20世纪20年代后，河套地震带成为鄂尔多斯块体周缘活动最为强烈的区域，先后发生6.0级以上地震6次，且主要发生在1970～1998年间。此时段内，河套地震带成为华北第4个地震活动期6级以上地震活动的主体区域。伴随着中强地震的持续活跃，河套地震带1970年以来的小震记录比较丰富，空间分布基本均匀。同时，针对该地区强烈的地震活动的地震台站建设亦迅速推进，经过“八五”、“九五”和“十五”，数字化观测系统的逐步完善，台站密度逐渐加大、台网布局日趋科学合理、监测能力逐步提升，这些基础条件都为开展该地区的b值分析提供了良好平台（图 2）。

图2 1970年以来研究区M L≥1.5地震以及断裂和地震台站的分布

2 b值的计算

b值的计算主要有线性最小二乘法、极大似然法、非线性最小二乘法、矩估计方法等4种。其实，矩估计方法最终得出的b值表达式与极大似然法相同，因此，矩估计方法可看做是极大似然法中的1种。就计算误差而言，极大似然法最小，线性二乘法次之，非线性最小二乘法最大（张建中等，1981）。就计算量而言，极大似然法最小，线性二乘法次之，非线性最小二乘法最大。但是，对于每1组计算，随着样本数目的增加，上述方法的计算结果的差异会越来越小，尤其当样本数目超过100时，4种方法基本表现出等同的计算效果（张建中等，1981）。换句话说，只要保证有足够多的地震样本数目参与计算，选择哪种计算方法似乎显得不太重要。本文将采用最小二乘法进行b值时间扫描计算，采用极大似然法进行b值空间扫描计算。

2.1 时间扫描计算

黄伟琼等（1989）分析了b值统计中的影响因素后认为，去除震级频度拟合曲线中小震级端的“掉头”部分能够得到更加合理的结果。显然，这个“掉头”点恰恰对应的是最小完整性震级Mc，因此，进行b值计算以最小完整性震级Mc作为起算震级M0是比较合理的。韩晓明等（2015）采用最大曲率法（Woessner et al，2005）求得河套地震带的最小完整性震级Mc＝ML2.7。能否简单的将该震级作为起算震级？答案是否定的，主要原因有2个：①1970年以来地震事件的震级分布表明（图3），ML2.0～2.5之间包括了2291次地震，占到地震总数的47.52%，如果按照最小完整性震级Mc值来进行地震样本截取，将大大减少每组参与计算b值的样本数，并最终影响计算结果的精度和可信度。②1970年以来的最小完整性震级并不恒定，随着时间的演进，其非均匀性变化特征显著，1990～2012年的最小完整性震级可达ML1.9（刘芳等，2013）。经综合考量，将本文中b值计算的起算震级确定为ML2.0。

图3 余震剔除前后震级分布（a）和时间分布（b）对比统计（1970～2012年，M L≥1.5）

2.2 最小二乘法时间扫描计算

首先考虑地震样本数及各震级档的分布情况。假如震级分布不均匀且样本数较少，则应适当提高震级分档标准，例如采用ΔM＝0.3作为分档步长；如果震级分布均匀，样本数不够多，则采用较小的震级分档标准，一般采用ΔM＝0.1作为分档步长；如果震级分布均匀，且样本数目充足，可采用微分频次进行分档（国家地震局预测预防司，1997）。根据1970年以来河套地震带的地震活动情况，本文采用ΔM＝0.1进行震级分档，以保证有较多的数据参与计算，间接提高结果的可信度。

设m为震级分档总数，Mi为第i档震级，Ni为第i档震级的实际地震数，则有

b值的标准偏差

其中

式中，a、b为线性最小二乘法拟合估计参数。

考虑到b值的估计过程更多的是考察“独立地震事件”，因此，本文在时间和空间扫描计算中都采用余震剔除后的地震目录（图3）。实际计算中，设定每次计算的窗长包含100个地震，滑动步长包含20个地震。结果显示，1970～2012年，河套地震带的b值经历了3次中期时长的趋势性下降过程：1974～1985年、1993～1999年、2007～2012年。由表1可见，第1次下降过程中，研究区分别发生了1976年4月6日和林格尔6.2地震、1976年9月23日巴音木仁6.2地震和1979年8月25日五原6.3地震。第2次下降过程中发生了1996年5月3日包头西6.4地震（图4）。第3次下降过程始于2007年，虽然2010年出现了小幅转折上升，但是整体变化仍处于均值线以下，且继续保持低b值状态至今。根据前2次趋势性下降过程中河套地震带中强以上地震的发生情形，不知目前河套地震带的低b值状态是否具有一定的指示意义。

表1 河套地震带1970年以来M S≥6.0地震

2.3 极大似然法空间扫描计算

相对于最小二乘法，极大似然法的求解过程比较简明，主要基于极大似然法和震级频度关系式给出b值的计算公式（国家地震局预测预防司，1997）

式中，N为地震总数，M0为起算震级，Mi为第i个地震的震级。

空间扫描计算同样采用余震剔除后的地震目录，考虑到河套地震带的区域形状和地震分布密度，用固定网格进行扫描计算会使结果图件出现多处空白，因此这里采用固定半径的圆形扫描窗口进行滑动计算。河套地震带整体呈长条形，东西长约660km，南北宽130～180km，因此圆形扫描半径的选取主要照顾南北宽度的范围。设定扫描半径（窗长）设为55km，滑动步长5km，每个扫描窗口内的地震样本数目不少于20个。

图4 b值（含方差）时序变化和6级以上地震的对应关系

图5 b值空间等值域分布

b值空间扫描图像显示，五原西北部、临河西北地区、包头至西山咀一带低b值异常比较突出（图5）。其中临河以北地区的主要活动断裂为狼山山前断裂，该断裂展布于狼山东南麓，由一系列向盆地倾斜的正断层或阶状正断层组成，从地貌、构造活动特征和第四系厚度来看，该断裂北段活动较为强烈，而低b值分布表明该区目前可能仍处于第四系以来的持续强烈活动过程中，区域地壳介质仍处于高水平应力作用状态。与之关联的是，五原西北部处于狼山山前断裂和色尔腾山山前断裂结合部或转折地段，该地段山势陡峭、第四纪断层和断崖十分发育，是狼山-色尔腾山山前断裂系活动最为强烈地段，较易产生应积累而使得b值出现低值异常。因此，b值空间扫描结果整体表现出与构造相依的分布特征。

重力、航磁和石油勘探等资料显示（国家地震局鄂尔多斯周缘活动断裂系课题组，1988），河套地震带并非单一构造单元，其内部存在2个次一级构造凸起，即包头凸起和西山咀凸起。这两大凸起整体位于河套地震带所辖的临河盆地和呼包盆地的中间，其构造位置十分特殊，且地震活动较为强烈。河套地震带的最大地震（公元849年包头西7.0级地震）和次大地震（1996年包头西6.4级地震，该地震也是距今最近的1次中强地震）都发生在该区域，使其成为全地震带迄今为止活动水平最高的区域。根据1970年以来相对完备的地震记录，包头至西山咀一带的地震频度和地震强度是河套地震带最高的地区，表明该区现今的构造活动依然强烈。震源机制解分析结果表明，包头-西山咀一带是河套地震带主压应力由北东转向北东东的力学节点，是全地震带受力最为复杂的地区；如果临河盆地和呼-包盆地受到外力作用而发生整体运动，那么包头-西山咀一带势必承受来自2个盆地的共同作用而成为应力积累的聚集点，这可能就是使该地区成为河套地震带6～7级强震集中区的动力根源。该区目前出现的低b值异常表明，上述区域正处于高应力集中状态，地壳介质非均匀特征相对较为显著，存在应力大释放的可能性。

3 影响b值计算结果的因素探讨

3.1 余震删除对b值计算的可能影响

在地震活动性研究中，余震剔除的必要性和合理性还在不断的讨论中，但是作者认为，在进行测震学指标计算过程中应该对限定“时-空-强”范围内的地震事件进行余震剔除，以此来确保“独立”的地震事件参与计算和分析。基于这种思路，本文采用G-C方法（陈凌等，1998）进行余震剔除（图3），并对余震剔除前后b值的时序变化制做叠加曲线，并进行了对比分析。

1970年1 月～2012年9月，河套地震带共计发生ML≥1.5地震4821次，进行余震剔除后保留了3520次“独立”地震事件，余震剔除比率27%。从余震剔除前后的对比统计图来看，除去4次6级以上中强地震的余震序列的调制影响，余震剔除前后，震级-频度分布和震级-时间分布皆保持稳定变化，整体构架未发生显著“畸变”（图3）。

从余震剔除前后b值时序变化的叠加曲线可以看出：①2000年以前，余震剔除前后b值时序变化叠加曲线“吻合度”较差，甚至在1975～1980年期间叠加曲线走势出现严重紊乱，这与该时段内研究区连续发生的3次6级以上中强地震的余震序列干扰有关（图4，图6）；②2000年以后，余震剔除前后b值时序曲线“吻合”较好，各个时段内基本保持同步变化，其原因有二。一是2002年以来，河套地震带一直保持有18个地震台，且空间布局未变，台网运行十分稳定（韩晓明等，2015）（图7）；二是2000年以来，研究区地震活动主要以中小地震为主，无中强地震序列及显著震群活动，地震活动未在时间和空间上形成丛集而对b值计算造成干扰。

进一步分析2000年以前的b值时序曲线可以发现，余震剔除前后，低b值阶段（低于均值线）叠加曲线的时序变化差异相对较大，高b值阶段（高于均值线）差异相对较小。根据河套地震带地震活动的起伏变化情况分析认为，2000年以前的高b值阶段对应的时段为1985～1995年，此间地震带未发生中强以上地震，所选取地震事件的震级和时序分布未受到地震丛集等“干扰”现象的影响而发生非均匀性变化；叠加曲线出现较大差异的低b值阶段对应的时段为1975～1980年，此间河套地震带连发1976年和林格尔6.2级地震、巴音木仁6.2级地震和1979年五原6.0级地震（表1），这3次6级以上地震成为影响河套地震带地震活动的“干扰”事件，且其余震序列在时间、震级和空间分布上同时形成“丛集”现象，从而改变了余震剔除前后河套地震带地震活动的整体构架，使得b值在余震剔除前后产生了较大变化（图 6）。

图6 余震剔除前后的b值（含方差）时序变化曲线

图7 1970年以来研究区地震台站个数的时序变化

3.2 计算样本的震级分布和数目对b值的影响

一般来讲，在测震学指标的求解过程中，足够多的地震样本数是保证计算结果可信度的重要一环，b值计算同样如此。为了分析不同震级档包含地震数对b值计算的影响，选取1970年以来的ML≥1.5地震作为考察对象，未进行余震剔除，运用线性最小二乘法计算了b值及其标准方差δb。结果显示，b值随着震级大小出现阶性段起伏变化，ML1.5～3.6之间，b值与震级呈现同步增大，ML3.7～4.4的地震出现反向变化，ML4.5以上地震又表现为同步上升，b值的计算误差δb也表现出了同样的变化过程。这应该与每个震级档所包含的地震事件的样本数目有关，震级优势分布范围（ML2.0～2.5）对应的b值计算误差δb相对较小，包含较少地震样本数目的震级范围对应的δb相对较大，δb与参与计算的地震样本数目呈现出明显的反向变化（图8）。

图8 地震样本数目对b值及δb的影响

4 结论和讨论

b值是一个历久弥新的测震学指标，伴随着数字地震资料的日益完善和丰富，b值在地震活动预测、区域应力场判定中的作用愈加重要，效果愈加明显。河套地震带具有相对稳定而统一的应力场环境和相对均匀的台网布局，是开展b值研究的理想场所。通过分析研究得到以下认识：

（1）河套地震带的b值时间分析结果显示，1970年以来，河套地震带的4次6级以上地震均发生在b值趋势下降过程中，从b值与区域应力场的关系分析，b值的趋势性下降反映了研究区应力水平的增强，也即这4次地震均发生在河套地震带应力水平逐步增强过程中，符合地震发生的动力学过程。从空间扫描结果看，位于临河盆地两端的临河-五原地区、包头-西山咀凸起是河套地震带应力水平较高区域，结合区域地质构造特征和地球动力学环境分析认为，研究时段内河套地震带的b值空间差异可能由地壳介质体性质和应力环境共同作用，相对来讲，b值的构造相依特征更为明显。

（2）b值计算中的起算震级选取过程显示，综合考虑Mc值的分时段变化以及保证足够多的地震样本依然是提高研究结果可信度的重要环节。在b值计算中，不能简单的以研究时段内的最小完整性震级作为起算震级，还应该考虑到研究时段内的地震震级构成、各震级档包含的地震数目等，综合确定地震事件的截取震级。因此，为提高b值的计算精度，应该综合考虑截取的震级下限和参与计算的地震样本数2个因素，并尽可能在两者之间找到平衡。

（3）b值计算中的余震剔除是必要的。余震剔除可以尽可能地避免或消弱中强地震的余震序列、震群活动等时空丛集地震事件的干扰，剥离出相对独立的“主震”事件进行计算，以期更加准确的表现出“时、空、强”限定范围内的b值变化。

致谢：审稿专家的建议对本文的完善和提高帮助很大。文中b值空间分布图由甘肃省地震局冯建刚提供的程序完成，其余图件由Matlab程序设计平台、GMT和Zmap程序包绘制完成。作者在此一并致谢！