中国西部地区利用烈度数据估计地震动参数的方法

陈 鲲 俞言祥 高孟潭 亢川川

1)中国地震局地球物理研究所，北京 100081

2)四川省地震局，成都 610041

0 引言

地震烈度与地震动参数的关系一直是国内外学者的研究热点。而大多数学者只是将地震动参数与烈度这2个变量通过简单回归的函数形式相联系(Gutenberg et al.，1956;Murphy et al.，1977;刘恢先，1978;刘贞荣，1982;李大华，1991;李大华等，1991;Wald et al.，1999a;Khosrow et al.，2001;Kazi et al.，2002;李山有等，2002;Yih-Min Wu，et al.，2003;Atkinson et al.，2007a;王玉石等，2008，2010;Faenza et al.，2010;李敏，2010)。上述研究基本上是用地震动参数估计台站位置处的烈度，而该过程往往是不可逆的，其反变换并不适用于用烈度来估计地震动参数。并且估计场点处地面运动参数如PGA，仅仅使用MMI和PGA的回归，而不考虑震源的任何信息以及从震源到场点处地面运动衰减特性，将存在很大的不确定性。

烈度与地震动参数的对应关系同震级和距离的相关性较强，特别是地震动的低频成分与烈度的关系同震级显著相关，反之高频成分与烈度的关系同距离相关(Atkinson et al.，2000)。如果烈度与地震动参数的对应关系不扣除其震级和距离的相关性，将产生明显的区域性(Atkinson et al.，2007a)。国外学者对烈度与地震动关系进行了大量研究，由于该关系存在显著的区域性，所以有些成果并不适用于中国西部地区。

许多重大的破坏性地震发生在现代地震仪器出现之前或者没有强震仪记录的地区。在重大历史地震强烈影响的地区，人们希望建造的建筑物能抵御这些复发地震的影响，因此地震学家和建筑设计师必须基于宏观烈度数据，利用模型或者间接证据来估计地面运动进行抗震设防。另一方面对于震后的地震动图速报工作，我们也希望大量灾害现场反馈回来的烈度数据可以对烈度速报的地震动参数分布图提供较好的约束。特别是特大地震发生后，发震断层的破裂尺度不能确定时，这些烈度数据对震动图提供约束就更为重要。美国地质调查局Did You Feel It系统通过互联网收集烈度信息，利用经验关系与地震动参数相联系，进而生成震动烈度图为地震应急服务(Wald et al.，1999b;Wald et al.，2005;Atkinson et al.，2007b)。

近年来，随着中国强震记录数据的增多，也逐步开展了中国西部地区烈度数据与地震动参数关系的研究。本文的主要目的是希望利用烈度点数据对地震动参数分布图提供一定的约束，进一步减小快速产出的震动图的不确定性。

1 数据

本文统计所采用的强震加速度记录，必须能够同时获得相应地震的烈度信息。我们收集了6次地震加速度记录和烈度数据共计186条加速度记录，其中包括宁洱2007年6月3日MS6.4地震8条，汶川2008年5月12日MS8.0地震151条，盈江2008年8月21日MS5.9地震4条，攀枝花2008年8月30日MS6.1地震14条，乌恰2008年10月5日MS6.8地震4条，姚安2009年7月9日MS6.0地震5条。各强震台站的烈度数据来自中国地震局应急救援司烈度评定报告的烈度等值线资料，其中地震烈度V度区为推测范围，根据Ⅵ度区宽度向远处拓展得到。186条强地震动记录资料所对应的强震观测台站中，1个台站位于烈度Ⅹ度区，1个台站位于烈度Ⅸ度区，13个台站位于烈度Ⅷ度区，35个台站位于烈度Ⅶ度区，72个台站位于烈度Ⅵ度区。其中4个台站位于烈度Ⅸ度区与Ⅷ度区分界线附近稍稍偏于Ⅷ度区，按照其位于Ⅷ度区处理。另有64个台站位于对应地震的Ⅵ度区以外，但与Ⅵ度区的距离小于Ⅵ度区的宽度，按照其位于V度区处理。各加速度记录的烈度分布见表1和图1。Ⅹ度及以上烈度主要是以地表震害现象作为评定标准，并且本文烈度Ⅸ和Ⅹ度的台站记录只有1个，因此本文不考虑根据Ⅸ和Ⅹ度及以上烈度数据来估计地震动参数。

为了增加研究数据和样本量，每条地震记录2个水平向加速度分量分别作为独立的样本。每条加速度记录经过基线校正，截止频率为0.02Hz的高通滤波后，求得各烈度档峰值加速度(PGA)、峰值速度(PGV)、5%阻尼比1s和3s反应谱谱值的均值和标准差。此外，希望能够反映场点处地面运动参数可能跨2个不同的烈度档的情况，因此合并邻近的2个烈度档(如烈度Ⅵ和Ⅶ)的地面运动参数数据，一起进行处理，并计算这些合并数据的均值和标准差，其结果示于表2。

表1 研究选用加速度记录的烈度分布Table 1 Numbers of strong ground motion records for each intensity

2 烈度估计地震动参数的贝叶斯方法

表2 各烈度档地震动参数的均值和标准差Table 2 Mean and standard deviation of ground motion parameters for each intensity level

区别于仅仅将一种参数与另一种参数的分布信息联系起来的方法，贝叶斯方法是基于贝叶斯概率理论，允许将一些独立的背景信息考虑到计算地面运动档的概率过程中(Ebel et al.，2003;Cua et al.，2008)。该方法认为围绕地震震中的各场点不同烈度档的地震动参数分布反映了地震的强度及地面运动局部传播的特性。用烈度估计地震动参数的贝叶斯方法可用以下公式表述

式(1)中:P(GM|I)是烈度I引起的地面运动GM的条件概率，而P(I|GM)是地面运动GM引起烈度I的条件概率。各烈度档的地面运动参数的对数值围绕其均值呈正态分布(Murphy et al.，1977;Faenza et al.，2010)。这样各烈度档的P(I|GM)可由烈度与地面运动参数的数据集经验性地确定，示于表2，确定方法详见下一节。P(GM)是可能在场点处被观测到的地面运动参数值的“先验”概率。P(GM)的先验概率是从当地衰减关系的标准差中得到，并且衰减关系可以用来计算该场点处任一地震动参数值的概率。分母是在场点可能经历的地面运动参数范围内进行累加。实际上，贝叶斯方法的计算过程是先选定一系列可能的地面运动参数值GMi，然后利用已知的烈度观测值I和式(1)来计算每个地面运动参数GMi的Pi(GMi|I)的概率。最大可能的概率值Pi(GMi|I)所对应的GMi值就是场点处最有可能经历的地面运动参数值。

衰减关系是地震震级、场点到震中或者发震断层的距离及场点处近地表地震波速的函数，其中明确包含了震级大小、从震源到场点距离及场点处场地条件的信息。由于贝叶斯方法将震源信息和地面运动局部传播的特性考虑进来，因此，能够有效地去除地震动参数与烈度对应关系中震级与距离的相关性。

3 P(GM|I)确定

图1表示各烈度档地震动参数的分布，沿常数烈度Ⅵ(图1中水平线)的分布代表P(GM|MMI)概率分布，沿常数PGA(图1中垂直线)表示P(I|GM)的分布。图2是图1峰值加速度PGA数据的3D直方图。由图2可以看出，受本文收集到的各烈度档的地震记录数量的限制，相应烈度档的PGA分布并不光滑，但总体的趋势仍近似满足对数正态分布。实际上，求解P(I|GM)的方法是认为每一烈度档的地面运动围绕某一均值呈对数正态分布，这样P(GM|I)能被解析计算。图3表示图1，2数据每个烈度档的PGAi的对数正态分布。图4表示几个PGAi档P(I|PGAi)概率分布，其中每一个PGAi的P(I|PGAi)分布都经过正规化处理。

图2 图1峰值加速度(PGA)3D频度分布图Fig.2 3D frequency distribution map for peak ground acceleration(PGA)in Fig.1.

4 贝叶斯方法的应用

图3 图2数据的对数正态分布图Fig.3 Logarithmic normal distribution of data in Fig.2.

图4 峰值加速度与烈度对应的P(I|PGA)分布图Fig.4 P(I|PGA)maps corresponding to peak acceleration and intensity.

2011年4月10日和8月11日分别发生了四川炉霍5.3级地震和新疆伽师5.8级地震。这2次地震均获得了台站观测值和相应的烈度值，因此可以利用这些数据来验证贝叶斯方法从烈度观测值中估算地震动参数的精度。从四川炉霍地震资料中获取了能确定烈度的6组加速度记录，结果见表3(国家强震动台网中心，2011)。其中最大地震动峰值加速度为炉霍地办台观测到的－359.863cm·s－2;从伽师地震资料中获取了能确定烈度的4组加速度记录，结果见表4(国家强震动台网中心，2011)。2次地震的烈度等值线资料来自中国地震局震后实际烈度调查数据。本次研究使用的衰减关系为汪素云等(2000)回归的中国西部PGA衰减关系以及肖亮(2011)回归的中强地震活动区(伽师地震)和川藏地区(炉霍地震)的PGA长轴衰减关系(以下称汪素云2000年回归的衰减为“衰减2000”，肖亮2011年的衰减关系为“衰减2011”)。

表3 根据M S5.3炉霍地震烈度数据用3种方法估算的地震动参数Table 3 Ground motion parameters estimated using intensity data of the M S5.3 Luhuo earthquake

5 比较与讨论

表3，4中所说的3种方法，第1种是利用衰减关系直接估计的结果;第2种是利用表3各烈度档地震动参数均值估计的结果;第3种是利用贝叶斯方法估计的结果，其中使用了衰减2000和衰减2011两套衰减关系，分别简称为Bayes1和Bayes2。

表4 根据M S5.8伽师地震烈度数据用3种方法估算的地震动参数Table 4 Ground motion parameters estimated using intensity data of the M S5.8 Jiashi earthquake

为了比较根据烈度数据用各种方法估算的地震动参数的精度，计算了2次地震用5种方法得到的估计值与观测值之间的均方根(rms)，所用公式如下:

式(2)中:Ypred，i表示第 i个台站的地震动参数估计值(本文使用 PGA)，Yobs，i为第i个台站的观测值，n为台站数量。2次地震的均方根结果示于表5。

从表3，4可以看出，用2种衰减关系直接估计的结果，衰减2011的PGA估计值基本上都大于衰减2000的结果。总体上说，衰减2011的估计值更接近台站观测值。对于炉霍地震震中距＞50km的2种衰减的PGA估计值基本一致，从图6，7可以更清楚地看出这种趋势。从表5可以看出，衰减2011的均方根值0.418 8小于衰减2000的0.464 3，因此衰减2011的估计效果优于衰减2000。并且用表2中各烈度档的PGA均值估计的均方根为0.267 3，远小于直接使用衰减关系估计的均方根。其原因可能是2种衰减关系均是基岩场地上的衰减关系，如果再考虑局部场地效应则会使估计值更接近真实的台站观测值，从而进一步提高估计的精度。而对于基于2种衰减关系用贝叶斯方法进行估计的PGA结果而言，从均方根大小可以看出，2种贝叶斯方法估计的均方根是5种方法中最小的。其中，基于衰减2011比基于衰减2000的贝叶斯方法的均方根小，为0.249 9。因此，基于衰减2011的贝叶斯方法的估计结果优于基于衰减2000的贝叶斯结果，说明贝叶斯方法对于其先验概率P(GM)的依赖性较强。如果贝叶斯方法中的先验概率再考虑局部场地效应，那么会进一步减小估计值与观测值之间的残差。

图5 5种方法的PGA估计值对比图Fig.5 Comparison of PGA estimated from five methods.

表5 用5种方法估算的2次地震总体参数的均方根Table 5 The rms of the population parameters of two earthquakes estimated by five methods

图6 四川炉霍地震3种峰值加速度估算方法结果比较Fig.6 Comparison of PGA estimated from three methods for the Luhuo earthquake.

图7 新疆伽师地震3种峰值加速度估算方法结果比较Fig.7 Comparison of PGA estimated from three methods for the Jiashi earthquake.

对于各个台站，估计的PGA值与观测到的PGA值差别最大的是四川炉霍地震炉霍地办的强震台。主要是因为该台的PGA位于衰减2000衰减关系3倍标准差之外，并且炉霍地办观测到的PGA与周围的观测到的PGA极不协调，可能与中强地震近场峰值加速度的高频脉冲有关(Ebel et al.，2003;Yi-Min Wu et al.，2003;陈鲲等，2012)，该PGA值并不能代表真实的烈度值。

图5表示2次地震PGA观测值与5种方法估计值的对比。从图5可以看出，衰减2011的贝叶斯方法的结果总体上效果最好，PGA只是在360cm·s－2左右差别较大。如前所述，这可能是因为四川炉霍地震炉霍地办的强震台中强地震近场峰值加速度的高频脉冲所致。

图6，7表示基于衰减2011衰减关系的贝叶斯方法、表3各烈度档对应的PGA均值以及PGA观测值的对比。总体上说，2次地震基于衰减2011衰减关系的贝叶斯方法的估计结果更接近实际的观测值。

6 结论

中国目前的地震台站布设密度，一次地震只能获取有限数量的台站记录。仪器记录仅仅是地震动参数值，所表示的是地面的振动程度，而烈度是一种宏观的震害现象。虽然烈度与地震动参数之间存在不小差别，但随着中国地震烈度速报工作的发展，震后较短时间内可能获得大量的烈度数据，也可以对烈度速报的震动图(如PGA，PGV等)提供一定的约束。特别是对于特大地震，在震源破裂的尺度不能确定的情况下，场点到断层破裂面的距离不能确定，这样宏观的烈度数据能更有效地控制地震的震动程度，是对台站观测数据的一种补充。综合以上结果得出以下结论:

(1)根据烈度数据用贝叶斯方法估算峰值加速度是可行性。虽然本文以峰值加速度为例说明贝叶斯方法的估计精度，但该方法同样适合于其他地震动参数的估计。

(2)用贝叶斯方法估算的参数精度与先验概率(衰减关系)显著相关，如果选用考虑局部场地效应的衰减关系会进一步提高地震动参数的估算精度。