基于深度卷积神经网络的2021年5月21—22日云南漾濞地震和青海玛多地震震级估算

朱景宝，宋晋东，李山有*

1 中国地震局工程力学研究所 地震工程与工程振动重点实验室，哈尔滨 150080 2 地震灾害防治应急管理部重点实验室，哈尔滨 150080

0 引言

地震预警是减轻地震造成人员伤亡和财产损失的有效手段之一.地震预警依据P波传播速度比S波传播速度快的原理，利用台站记录到的P波初期信息对地震基本参数和破坏区域进行估算，抢在破坏性地震波到达之前向公众发布预警信息，进而让人们提前几秒或数十秒采取紧急避险措施(马强，2008; 宋晋东，2013).世界上一些地震较为活跃的国家和地区已经建立了地震预警系统，如：美国(Allen et al., 2009; Kuyuk et al., 2014)、日本(Hoshiba et al., 2008; Kamigaichi et al., 2009)、意大利(Satriano et al., 2011; Colombelli et al., 2020)、中国台湾(Wu and Teng, 2002; Hsiao et al., 2009)、墨西哥(Aranda et al., 1995; Suárez et al., 2018)等.同时，中国地震预警系统也在进行验证测试(Peng et al., 2011; Zhang et al., 2016).

在地震预警系统中，地震破坏区域的准确估算以及可靠的预警信息发布依赖于震级估算的结果，因此准确而快速地估算震级，是地震预警系统的重要工作之一，对地震预警也有重要的意义.传统震级估算方法主要根据P波触发后单一特征参数与震级的经验关系，进而估算震级(张红才等, 2012).常用的特征参数包括幅值类参数(Wu and Kanamori, 2005; Wu and Zhao, 2006)、周期类参数(Kanamori, 2005; Böse, 2006; Huang et al., 2015)和能量类参数(Reed and Kassawara, 1990; Nakamura, 2003; Festa et al., 2008).尽管上述特征参数已经被证明与震级存在一定的线性关系，但是某一特征参数包含与震级相关的信息较单一，导致震级估算结果的误差较大，并且震级估算结果也存在一定程度的小震高估和大震低估问题.因此，一些研究人员尝试采用双参数对震级进行估算，相对于单参数震级估算结果，有效提高了震级估算的准确性(Cuéllar et al.，2018).

近年来，人工智能已被引入到地震预警震级估算研究中.一些学者也尝试将多个地震预警参数与机器学习方法相结合，进而对最终震级进行估算(Münchmeyer et al.，2020；朱景宝等，2021a).Zhu等(2021)使用日本K-net台网数据，融合多种类型的12个P波初期特征参数(幅值类参数、周期类参数、能量类参数、衍生类参数)作为卷积神经网络的输入，构建了P波到达后3 s时间窗下的DCNN-M模型，进而用于地震预警震级估算，其结果优于传统的单参数τc方法和Pd方法.Ochoa等(2018)基于支持向量机回归算法，在单台情况下，使用P波触发后5 s的信号得到的25个参数(与方位角相关的参数、与震级相关的参数、与震中距相关的参数)作为输入，进而对里氏震级进行估算，其结果表明：对于ML4.5以下的地震，震级估算有较小误差，且平均绝对误差为0.19个震级单位.此外，还有一些学者尝试将地震波形数据与机器学习方法相结合，进而改进震级估算结果.胡安冬和张海明(2020)使用日本Kik-net和K-net台网强震动数据建立了NN网络用于地震预警震级估算，将波形的频谱信息作为模型输入，并发现在3 s时间窗内，该方法优于传统的单参数τc方法，同时也认为多特征组合输入的使用将会得到更优异的震级估算模型.Mousavi和Beroza(2020)基于30 s的波形输入，构建了一个由卷积神经网络和循环神经网络组成的震级估算模型，有效改善了震级估算结果.然而，已经建立的人工智能地震预警震级估算模型对于中国地震的震级估算是否可行亟待研究.

2021年5月21日和22日，在中国云南漾濞和青海玛多分别发生了MS6.4地震和MS7.4地震，同时还伴随一些前震和余震的发生.据大理市抗震救灾指挥部和青海省应急管理厅发布的消息，地震造成一定程度的人员受伤，以及交通、通信、供水等基础设施的破坏.这次地震得到了国内许多研究者的关注，也为人工智能地震预警震级估算模型测试验证提供了数据支持.

本文依据Zhu等(2021)提出的P波到达后3 s时间窗的DCNN-M模型框架，分别使用P波到达后1～40 s不同时间窗的特征参数作为输入进行训练和验证，构建不同时间窗下的单台DCNN-M模型，并采用多台加权平均方法对2021年5月21日和22日发生在中国云南漾濞和青海玛多的主震、余震和前震进行实时震级估算分析，旨在探索DCNN-M模型对这次漾濞地震和玛多地震震级估算的可行性、以及未来在我国人工智能地震预警系统应用的可能性.

1 方法

卷积神经网络主要用于信号数据的处理以及提取数据层次特征.Zhu等(2021)使用日本K-net台网强震动数据并基于深度卷积神经网络构建了P波到达后3 s时间窗下的单台DCNN-M模型，其目的是为了提高地震预警震级估算的准确性；同时，Zhu等(2021)研究结果表明：对于MJMA3～7.5的地震事件，在P波到达后3 s时间窗下，DCNN-M模型的震级估算结果的准确性高于传统的单参数τc方法和Pd方法，且没有出现明显的小震高估和大震低估现象.

1.1 特征参数

由于深度学习的特征提取是黑箱子，提取特征的物理意义并不明确.因此，DCNN-M模型的输入不同于Mousavi和Beroza(2020)、胡安冬和张海明(2020)等采用的波形输入，而是采用与震级相关的特征参数组合作为输入，再通过DCNN-M模型提取这些特征参数的高阶特征，使得该模型具有物理意义(Zhu et al., 2021).该模型的输入是由4类特征参数，共12个特征参数组成的一维向量，幅值类参数包括峰值位移Pd(Wu and Zhao, 2006)、峰值速度Pv和峰值加速度Pa(Wu and Kanamori, 2005)；周期类参数包括特征周期τc(Kanamori, 2005)、峰值比Tva(Böse, 2006)、构造参数(Huang et al., 2015)；能量类参数包括速度平方积分IV2(Festa et al., 2008)、累积绝对速度CAV(Reed and Kassawara, 1990)、累积能量变化率DI(Nakamura, 2003)；衍生类参数(Zhu et al., 2021)包括竖向累积绝对位移cvad、竖向累积绝对速度cvav、竖向累积绝对加速度cvaa.

1.2 DCNN-M模型

DCNN-M模型的网络结构是由四个卷积层、四个Batch Normalization层(BN层)、四个池化层和三个全连接层组成(Zhu et al.，2021).四个卷积层中卷积核的数目分别是124、150、190和250；卷积核的大小为4，移动步长为2，padding类型使用的是“same”；池化层采用的是最大池化，每个池化核的大小是2，移动步长是2；三个全连接层神经元数量分别是250、125、60；最后一个输出层的神经元数目是1，输出结果为预测震级.

本文在Zhu等(2021)提出的P波到达后3 s时间窗下的DCNN-M模型框架以及超参数的基础上，同样采用Zhu等(2021)所使用的训练数据集和测试数据集的日本K-net强震数据，分别计算了P波到达后1～40 s不同时间窗的特征参数作为输入，并对模型进行训练和验证，进而构建了在P波到达后不同时间窗下的单台DCNN-M模型.

由于30 s之后的震级估算结果趋于稳定状态，所以图1展示了P波到达后1～30 s不同时间窗下的单台DCNN-M模型对验证数据集(这里的验证数据集是指Zhu等(2021)文献中所使用的测试数据集)的震级估算结果.图中红色实线是编目震级与估算震级的1∶1线性关系；红色虚线表示±0.5震级单位误差范围；σ是估算震级与编目震级误差的标准差，用于衡量误差的离散性.从图1中可以发现：随着P波到达后时间窗的增加，震级估算误差逐渐减小，且对于6级以上震级低估问题也逐渐改善.

图1 验证集在1～30 s时间窗下震级估算结果Fig.1 Magnitude estimation at 1～30 second time window for validation dataset

1.3 多台加权平均

图2展示了构建的单台DCNN-M模型对训练集和验证集震级估算误差的标准差与P波到达后时间窗的关系.从图2中可以发现：训练集和验证集的震级估算误差标准差比较接近，且随着时间窗的增加逐渐减小，这也说明我们所构建的P波到达后1～40 s的不同时间窗下DCNN-M模型有较好的连续性和泛化能力.

图2 基于DCNN-M模型的训练数据集和验证数据集的震级估算误差标准差与时间窗的关系Fig.2 Relationship between standard deviation of magnitude estimation error and time window for training dataset and validation dataset based on DCNN-M model

我们将基于P波到达后不同时间窗下建立的单台DCNN-M模型以及多台加权平均方法对2021年5月21日和22日发生在中国云南漾濞和青海玛多的主震、余震和前震进行实时震级估算.

根据图2中验证集的震级估算误差标准差随着时间窗的增加而减小，且验证集反映的是模型的性能，同时，标准差反映的是相应时间窗时模型在单台情况下的震级估算的离散程度.因此，随着时间窗的增加，震级估算的标准差越小，说明单台的震级估算准确性越高，所以在多台加权平均计算的时候，权重与验证数据集在P波到达后相应时间窗下的DCNN-M模型震级估算的标准差成反比(Colombelli et al., 2012; Peng et al., 2017).

实时多台震级加权平均计算公式如下：

(1)

(2)

2 数据及处理

中国地震局工程力学研究所为本研究提供了2021年5月21—22日云南漾濞和青海玛多地震的强震动数据，震级范围MS4.0～7.4，共11次地震事件，129条三分向地震记录.在本研究中，没有对数据的震中距以及信噪比进行筛选.另外，DCNN-M模型是基于日本K-net台网建立的，采用的震级标度是MJMA；这两次地震序列的震级标度采用的是MS.依据Peng等(2017)、Kuyuk和Allen(2013)、Olson和Allen(2005)建立多国家多地区震级统计关系对于震级标注的处理方式，即不区分震级标度，本研究将估算震级标度统一记作为M.

图3展示了本研究中震中以及台站的分布，记录到青海玛多MS7.4主震的最近台站是63DAW,震中距为175.6 km，且UD、EW、NS方向的最大加速度分别为-19.1 cm·s-2、46.0 cm·s-2、40.6 cm·s-2；记录到云南漾濞MS6.4主震的最近台站是53YBX,震中距为7.9 km，且UD、EW、NS方向的最大加速度分别为-448.4 cm·s-2，-379.9 cm·s-2，-720.3 cm·s-2(朱景宝等，2021b).

对于本研究中采用的地震记录，首先采用马强等(2013)提出的P波捡拾方法确定P波到时；然后对加速度记录进行一次积分得到速度记录，再对速度记录进行一次积分得到位移记录，并对积分后的结果使用4阶0.075 Hz的巴特沃斯滤波器进行高通滤波，以消除积分带来的低频飘移的影响(Peng et al., 2017; 宋晋东等, 2018).

3 结果

为了验证DCNN-M模型对于漾濞和玛多地震震级估算的鲁棒性，图4和图5分别展示了云南漾濞MS6.4主震和青海玛多MS7.4主震在首台触发后实时的震级估算结果，以及触发台站数随首台触发后时间增加的变化.图4和图5的横坐标是首台触发时间，即第一个触发的台站后n秒.图中黑色实线表示实际震级；黑色虚线表示实际震级±0.5震级单位误差范围；红色实线表示首台触发后，n秒内DCNN-M模型对已触发台站震级估算结果的多台加权平均值；蓝色虚线表示距离震中最近台站的震级估算结果；其他颜色的实线表示已经触发的单个台站的震级估算结果.

图4 云南漾濞6.4级主震的震级估算结果Fig.4 The magnitude estimation of Yunnan Yangbi MS6.4 mainshock

从图4中可以观察到：对于漾濞主震，大部分台站的单台震级估算结果都分布在±0.5震级单位误差范围内，这表明本文训练后的DCNN-M模型具备单台震级估算准确性；首台触发后2 s时，DCNN-M模型的震级估算结果为M6.0；随着首台触发时间的增加，估算震级逐渐接近实际震级，且在首台触发后9 s时，DCNN-M模型的震级估算结果为M6.2；同时，随着首台触发时间的持续增加以及触发台站数的增加，DCNN-M模型的震级估算结果仍在不断接近实际震级.

从图5中可以观察到：对于青海玛多MS7.4主震，大部分台站的单台震级估算结果都分布在±0.5震级单位误差范围内，只有63DAW、63HEN、63XIH、63GUD台站在触发后1 s内震级估算结果在±0.5震级单位误差范围以外，这表明本文训练后的DCNN-M模型具备单台震级估算准确性；首台触发后1 s时，DCNN-M模型的震级估算结果为M6.4，可以判断青海玛多主震为一次大震；首台触发后2 s时，DCNN-M模型的估算震级是M7.0；随着首台触发后时间的增加，在首台触发后9 s时，DCNN-M模型的估算震级是M7.3；且随着首台触发后时间的持续增加和触发台站数的增加，DCNN-M模型的震级估算结果也逐渐接近实际震级.

图6展示了这11次地震事件基于DCNN-M模型的震级估算误差随首台触发后时间的变化.图中红色实线表示青海玛多MS7.4主震的震级估算误差；蓝色实线表示云南漾濞MS6.4主震的震级估算误差；绿色实线表示漾濞余震以及前震的震级估算误差；黑色虚线表示±0.5震级单位误差范围.从图中可以看到：首台触发后1 s时，对于余震以及前震的震级估算误差主要集中在±0.5震级单位范围内，对于云南漾濞MS6.4主震和青海玛多MS7.4主震的震级估算结果分别是M5.5和M6.4，误差接近1.0震级单位，对于玛多主震和漾濞主震在首台触发后1 s时低估的原因可能是大震在发生后几秒内并没有完全破裂，且P波并未携带完整的与震级相关的信息；但是在首台触发2 s以后，云南漾濞MS6.4主震和青海玛多MS7.4主震的震级估算误差都在±0.5震级单位内，且随着首台触发后时间的增加，误差逐渐减小，该结果也与图4和图5的结果相一致.

图5 青海玛多7.4级主震的震级估算结果Fig.5 The magnitude estimation of Qinghai Madoi MS7.4 mainshock

图6 震级估算误差随首台触发后时间的变化Fig.6 The variation of magnitude estimation error with the trigger time of the first station

4 讨论与结论

本研究在Zhu等(2021)使用日本的K-net强震动数据构建的P波到达后3 s时间窗的DCNN-M模型框架和超参数的基础上，分别使用P波到达后1～40 s不同时间窗下的特征参数作为输入进行训练，构建了在P波到达后不同时间窗下的单台DCNN-M模型，并对2021年5月21—22日云南漾濞和青海玛多的MS4.0到MS7.4共11次地震进行震级估算，分析了基于单台DCNN-M模型和多台加权平均方法得到的震级估算结果随首台触发后时间的变化，探索了DCNN-M模型对于这几次地震在地震预警震级估算中的可行性.

通过对云南漾濞MS6.4主震和青海玛多MS7.4主震在首台触发后实时的震级估算结果的观察：在首台触发后2 s时，DCNN-M模型对云南漾濞地震的主震和青海玛多地震的主震的震级估算结果分别是M6.0和M7.0；在首台触发后9 s时，DCNN-M模型对云南漾濞主震和青海玛多主震的震级估算分别是M6.2和M7.3；随着首台触发后时间的增加和触发台站数量的增加，震级估算结果逐渐接近实际震级.同时，这两次地震事件的前震和余震，在首台触发后1 s时，震级估算的误差主要都集中在±0.5震级单位内；随着首台触发后时间的增加，这11次地震事件的震级估算误差逐渐减小.本文的研究结果表明：DCNN-M模型对于2021年5月21—22日发生在中国云南漾濞和青海玛多的11次地震事件在首台触发后10 s内可以得到较为稳定的震级估算结果.这为我国建立人工智能地震预警系统提供了潜在可能.但是，对于深度卷积神经网络震级估算模型DCNN-M在中国地震预警系统中的应用还需要大量的中国数据进行验证和测试.

需要注意的是，由于作为DCNN-M模型输入的幅值类参数、能量类参数和衍生类参数是需要进行震源距修正，统一校正到参考震源距10 km(Zollo et al., 2006; 彭朝勇等, 2013)，其目的是为了消除幅值类参数、能量类参数和衍生类参数因距离衰减对震级估算的影响，所以在实际的地震预警过程中，DCNN-M模型准确的震级估算结果还依赖于地震预警系统的准确定位.此外，由于MJMA和MS这两个震级标度并非完全相同，依据Peng等(2017)，只是近似的将这两个震级标度相等，并且将估算震级统一记作M，而这两种震级标度的差异性可能会对本研究的结果造成一定的误差影响，对于这两个震级标度的准确转换仍需更深的研究和讨论.

同时，本研究是基于日本K-net强震数据训练得到的DCNN-M模型对漾濞和玛多地震的主震、余震和前震进行震级估算.更进一步的改善对漾濞和玛多地震的主震、余震和前震的震级估算结果，我们还考虑了两种方案：一、在原有DCNN-M模型的基础上加上中国强震数据进行迁移学习训练，而是否可以成功应用迁移学习方法改进结果还存在争议；二、使用中国目前已经公开的地震数据，建立属于中国的DCNN-M模型，而对于目前中国的地震数据量是否可以建立泛化能力强且具有稳定性的模型还有待研究.这些方案也是我们对国内地震震级估算在下一步研究中需要探索的地方.

致谢中国地震局工程力学研究所为本研究提供数据支持.文中图件使用通用制图工具GMT(Genetic Mapping Tools)绘制.