应用统计物理方法分析川滇地区地震危险性

王 凡 沈正康 王 敏 王阎昭 陶 玮

1)中国地震局地质研究所，地震动力学国家重点实验室，北京 100029

2)北京大学地球与空间科学学院，北京 100871

0 引言

川滇地区位于青藏高原的东缘，构造运动强烈。受印度板块NNE向推挤和青藏高原E向挤出的动力学作用控制，川滇地区构造形变场主要表现为绕东喜马拉雅造山结的顺时针大回转及与之相伴的大小地块间的断层错动，其主要活动断裂包括构成青藏高原与四川盆地边界的龙门山断裂带，高原内部的龙日坝和岷江断裂带，构成川滇“菱形地块”东边界与北边界的鲜水河－小江断裂带及南边界的红河、南华－楚雄－建水、无量山断裂带，菱形地块内部的丽江－小金河断裂带，南部的龙陵－瑞丽断裂带，以及西南的怒江、澜沧江断裂带等。伴随着强烈的构造运动，本地区地震活动频繁，据中国地震台网中心记录(图1)，自1970年以来该地区发生5.5级以上地震达280次之多。

图1 川滇地区活动构造及1500年以来5.0级以上地震Fig.1 Tectonic map of Sichuan-Yunnan region and M≥5.0 earthquakes of 1500—2010.

根据历史记载，大地震通常会造成巨大的人员伤亡和财产损失。因此，川滇地区的地震危险性分析一直是地震预测研究的热点(闻学泽，2001;易桂喜等，2002，2004，2006，2008;蒋长胜等，2008，2010)。易桂喜等(2002)和徐锡伟等(2005)利用简单的点过程统计方法研究地震复发行为，结果表明，这一区域活动断裂带的强震复发表现出趋于随机的丛集行为，其过程不存在良好的准周期性;同时，基于地震复发间隔与震级的统计研究表明，该地区所发生的地震也不具有“强度可预报”或“时间可预报”的性质。因此，基于断裂带整体强震复发性质的中长期地震危险性评估具有一定困难。

自20世纪80年代末，国际上逐渐将现代统计理论和方法应用于地震预测问题的研究并取得多方面的进展，已成为物理学和地震学之间一个十分活跃的交叉领域。Kagan等(1994)提出了基于地震空间相关性的方法，根据历史地震目录预测未来地震发生的概率，并将此方法应用于日本(Jackson et al.，1999;Kagan et al.，2000)、中国(Rong et al.，2002)以及全球(Kagan et al.，2011)。Rundle等将图像信息学(Pattern Informatics)算法应用于地震活动性图像分析，并尝试对加州的地震做出预测(Rundle et al.，2000，2003;Tiampo et al.，2002)。蒋长胜(2008)将此算法应用于川滇地区，其回溯性检验表明该方法对于川滇地区的中长期地震活动具有一定的预测能力。Keilis-Borok等(1990)利用TIPs(Time of Increased Probability)方法对意大利中部地区进行了长期强震预测。Huang等(2001)利用RTL(Regional-Time-Length)方法研究了M 7.2阪神地震造成的区域地震危险性变化。除利用地震目录进行统计分析的方法外，Shen等(2007)提出了基于地壳形变场建立地震危险性预测模型的方法，研究发现，过去50a以来南加州地区的强震分布与应变率分布有很高的相关性。以上这些研究证实了基于统计理论的方法应用于不同观测数据进行地震危险性分析的有效性。

本文分别尝试应用基于地震空间相关性和地壳形变场的方法分析川滇地区的长期地震危险性，并通过对预测结果的回溯性检验分析该方法的有效性。

1 基于地震空间相关性的方法

1.1 基于地震空间相关性的地震危险性模型

基于地震空间相关性的地震危险性预测方法由Kagan等(1994)最先提出，其基本思想是根据历史地震目录建立统计学模型，从而得到未来发生地震的概率。假定地震发生概率与历史地震发生频度成正比，地震震级－数目分布遵从修正的古登堡－里克特关系(Bird et al.，2004)，地震发生为泊松过程，由此可以根据区域地震目录给出强震发生的统计预测。地震概率密度函数可表示为

其中，Λ为单位时间、单位面积内发生地震的概率密度，θ和φ分别表示经度和纬度，m是震级，G为地震震级-数目分布函数，F为地震空间分布概率函数，是由根据历史地震目录得到的地震空间分布概率函数与背景场地震概率C的加权和。之所以采用加权和是为了避免某些没有历史地震的地区未来发生地震的概率为零。地震空间分布概率函数是历史地震目录的平滑函数，本文采用的平滑函数如下:

其中φi和θi分别是地震i震中的纬度和经度，rs是尺度参数，它与很多因素有关，如地质构造、地壳圈层厚度以及地震定位的准确度等，因此，不同区域的rs是不同的。

利用1500—2010年5.0级以上地震目录提供约束，给出了川滇地区未来发生5.0级以上地震的概率(图2)。在基于地震空间相关性的预测方法中，每一次地震都是独立的，结果大震的前震和余震都作为独立事件参与统计。因此，只选取5.0级以上的地震，既能保证有充分的数据约束预测模型，又可以在一定程度上去除部分前震以及余震以获得尽可能合理的结果。此外，我们还在一定的范围内遍历尺度参数rs，结果表明，当rs=15时预测模型能更好地符合检验阶段的地震目录。由图2可以看出，地震危险性较高的地区包括龙门山断裂中2008年汶川地震破裂段及邻域、虎牙断裂、鲜水河断裂、龙陵地区、马边断裂及丽江－小金河断裂中－南段。此外，地震危险性较高地区还包括小江断裂、楚雄－建水断裂南端、无量山断裂中段及邻域地区等。

图2 由1500—2010年地震目录约束得到的5级以上地震发生的概率，背景色表示地震概率Fig.2 Earthquake long-term potential based on earthquake catalog from 1500 to 2010.Colors show the long-term probability of earthquake occurrence.

1.2 模型结果的回溯检验

基于地震空间相关性的预测方法能否准确估计川滇地区长期地震危险性?为验证其有效性，下面采用2种方法对模型的预测结果进行了回溯性检验。方法一为基于似然函数的检验方法，即用似然函数比来检验预测模型与实际发生地震的符合程度(Kagan et al.，1977)。其中，j是网格序号，N是网格数，vj是归一化的网格内发生地震的概率，τj为归一化的网格面积，n为检验阶段地震目录中的地震数，λi为预测模型得到的地震i的震中所在网格的单位面积单位时间内发生地震的概率，ξ为单位面积、单位时间内发生地震的统计值，I0为基于预测模型的地震目录，I1为基于检验阶段的地震目录，两者之差反映了预测模型的准确性(详细介绍见Kagan，2009)。

检验方法二为区域地震概率与地震发生数的累积函数检验(Kagan，2009)。将研究区分成小的地域单元，按照地域内地震发生概率由高到低排列，排列后的概率函数按面积积分，同时将单元内发生的地震数累积。将2个函数在终点归一化，得到概率积分函数与地震数累积函数。如果地震发生率符合地震概率预测模型，则2条函数曲线应有很好的相似性;如果实际地震发生曲线在概率曲线之上，表明地震的发生概率高于概率模型预期，即在概率模型预期地震发生率高的地区实际地震发生数大于概率预期值，而在概率模型预期值低的地方则相反。

图3 利用历史目录的地震预测模型Fig.3 Earthquake forecast test using historical catalog data.

将1500年以来5.0级以上地震分为1500—1999年和2000—2010年2个部分，第1部分用于约束地震概率估计，所得结果见图3a;第2部分用于回溯检验。第1种检验方法的结果列于表1，结果中I0、I1差异较大，表明预测模型不够准确。第2种检验方法的结果见图3b，显示实际地震发生曲线位于概率模型曲线之下，同样表明概率模型(图3a)的预测结果并不好。考察过去10a地震发生的空间分布，发现确实有不少地震发生在预测概率相对高的地区，但仍有相当多的地震发生在预测概率较低的地区，其中特别引人注目的是2008年的汶川地震。综其原因可能是预测方法只是在大样本空间统计意义上有效，当统计样本数较少时，如本例，小概率事件出现的可能性大大增加，因此对统计结果的影响极大。普遍认为发生汶川地震的龙门山断裂带的特征地震复发周期在千a以上(Densmore et al.，2007;Burchfiel et al.，2008;Shen et al.，2009;Zhang et al.，2010)，有历史记录的500a以来极少有强震发生，故本方法得到的汶川地震发生地预测概率值很低，预测并不成功。

表1 基于似然函数的检验结果Table 1 Forecast testing catalog and information scores

针对上述预测模型存在的样本数较少的问题，进行了另一个回溯检验，将1500—1972年的地震目录用于约束地震概率估计，得到的预测概率如图4a所示，1973—2010年的地震目录用于模型检验。检验方法一的结果见表一，I0和I1二者差异约为12%，表明模型预测效果较好。检验方法二得到的结果如图4b所示，从中可以看出预测函数曲线与实际地震发生函数曲线符合得较好，虽然二者在10%高概率地区的范围仍有一定差距，但在20%高概率地区的范围已相当接近。结果进一步显示50%强震发生在模型给定20%高概率地区，而80%强震发生在50%高概率地区。2种检验都表明在较为充分的样本空间内基于地震空间相关性的地震预测方法可以很好地反映地震发生率，虽然在中、短期地震预测中有一定局限性，但对于长期地震预测不失为一种有效方法。

图4 a背景为利用1500—1972年地震目录约束得到的5.0级以上地震发生的概率，黑色圆圈为1973—2010年发生的5.0级以上的地震;b区域地震概率与地震发生数的累积函数检验，红线为地震预测概率积分函数，蓝线为实际发生地震的累计数Fig.4 The same as Fig.3 except that the earthquake catalog data are of 1500－1972 and the M ＞5.0 events are of 1973 －2010 in(a)，and the test curves in(b)are plotted accordingly.

2 基于地壳形变场的方法

2.1 基于地壳形变场的地震危险性模型

基于地壳形变场的地震危险性模型(Shen et al.，2007)假定区域长期地震危险性与区域构造应变率成正比，地震发生数遵从修正的古登堡－里克特震级－地震数统计关系，则地震发生概率可表示为

其中A(x)为大地测量得到的地表应变率场，G(m)为修正的古登堡－里克特震级－地震数统计关系。本文收集了研究区内1999—2011年所有可利用的GPS观测资料，采用非线性模型模拟台站位置变化，尽可能地消除2001年MW7.8可可西里地震、2004年MW9.2苏门答腊地震以及2008年MW7.9汶川地震的同震位移和震后形变的影响，得到了反映区域构造运动的速度场(图5)，具体方法见王敏等(2003)和Shen等(2011)。

由区域构造运动速度场求解应变率场，采用了一种连续滑移回归方法(Shen et al.，1996;Jackson et al.，1997)。其主要原理是对于每个计算点，集合周边台站速度场数据，用最小二乘法估算平均应变率场，所用速度场数据依据相邻台站的方位角大小和与计算点的距离加权，其中距离加权函数表示为exp(－Δ2/σ2)的形式，Δ为台站至计算点距离，σ为平滑因子。不同σ所对应的总权重和应变率估值的形式误差之间存在负相关关系，这一负相关关系中的最佳平衡点即为最优的σ值。

图5 由GPS观测获得的区域水平运动速度场(相对稳定的欧亚板块)Fig.5 GPS horizontal velocity field with respect to the stable Eurasian reference frame.

图6 背景为应变率第二不变量，白色菱形为1500—1972年间M＞6.0地震的震中，白色圆圈为1973年以来M＞5.5地震的震中Fig.6 Geodetic strain rate field and spatial distribution of strong earthquakes.The color background is the second invariant of strain rates，and white diamonds and white circles are epicenters of M ＞ 6.0 events occurring in 1500－1972 and M ＞5.5 events occurring in 1973－2010，respectively.

2.2 模型结果的回溯检验

为了对地震的空间分布与应变率场的相关性有更客观的认识，应用1.2节中区域地震概率与地震发生数的累积函数检验方法对上述预测模型进行回溯性检验。图7a和图7b分别为1973年以来5.5级以上地震和1500年以来6.0级以上地震与应变率场的关系。在有地震仪器记录之前(即1500—1973年)地震发生地点与震级量度都存在一定不确定性，并且可能还有未被记录的地震。考虑川滇地区历史上居住人口较为稠密，过去500a中6.0级以上地震被遗漏的可能性不大，而导出应变率场的平滑因子也在几十km范围，因此地点不确定性的影响不会很大。从图7a可以看出，1973年以来发生在低应变率地区的5.5级以上地震的数目明显多于基于应变率场的预期，说明应变率场不能作为过去37年强震发生地点的指标，而图7b则显示1500年以来6.0级以上地震的空间分布与应变率场存在很好的相关性，即在数百a时间尺度上，应变率的高低可以作为强震发生地点的指示因子。该地区10a尺度的地壳应变率场实际上反映了本地区长期的构造形变场，而这一形变场的主体为弹性形变，只能通过地震、特别是强震来释放，但其释放时间需要数百a。

图7 应变率与地震发生数的相关关系Fig.7 Geodetic strain rate vs.earthquake count.

3 讨论与结论

Rong等(2002)应用基于地震空间相关性的预测方法研究了整个中国大陆的地震危险性，采用的地震目录为1978—1999年，截止震级为M＞5.4，检验时间区间仅为2a(2000—2002年)。相比之下，我们研究的川滇地区范围小且地震频度大，约束模型和检验模型的数据区间都远远大于前者。因为使用的参数不同，2个模型的预测结果不能严格比较，但从回溯检验的结果不难看出(本文图4b与Rong等(2002)图3)，本文模型预测结果与实际发生地震情况符合程度明显更好。这也再一次说明基于地震空间相关性的预测方法对历史地震目录有很强的依赖性，历史地震目录越长，记录越完整，模型的预测结果越可靠。在较为充分的样本空间内，基于地震空间相关性的预测方法对于长期地震预测具有一定的有效性。

地震发生的过程很短，也就是几十s，但是地震孕育的过程很长，通常要几百a，甚至上千a。在几十a时间尺度内，地震的发生有丛集性的特征，而基于地震空间相关性的方法所根据的地震目录只有近几十a比较完整，所以对这一时段有较高的有效性;而数十a的GPS观测获得的区域应变率场更能代表区域长期的应变率场，在这个时间尺度内，地震的发生主要表现为泊松过程。不难理解，相比过去30a发生的地震，过去500a间发生的地震与基于地壳形变场的地震预测结果有更好的一致性。因此，区域应变率场可以作为长期地震危险性预测的指标，但需要特别关注2个方面:其一，由数十a尺度的大地测量资料获得的区域应变率场应该能够代表区域长期的应变率场;其二，大地测量观测得到的应变通常包含了弹性应变和非弹性应变，要注意区分这2种应变，因为只有弹性应变能是与地震相关的。

本文的研究表明，基于地震空间相关性和地壳形变场的地震预测方法对于川滇地区的长期地震活动具有一定的预测能力，这为地震危险性分析和地震预测研究提供了新的线索。但无论如何，本文仅是概率性的地震预测方法在川滇地区的初步应用探索，未来应当考虑开展更深入的研究。