1. 引言
地震是板块碰撞引起的地壳形变和振动，同时是一种破坏力极强的自然灾害，通常会给地震区域的人带来巨大的生命财产安全损失。我们国家的云南四川地区，由于位于亚欧板块和印度洋板块边界，板块活动强烈，因而地震频发。2008年8.0级汶川地震造成了10万人左右的死亡 [1] ，近期的地震如2022年5月的四川雅安芦山6.1级地震和2022年9月的四川甘孜泸定6.8级地震，给当地的房屋、道路以及公共设施造成了不小的损失。
大量研究表明，地震之前会有显著的地磁异常：1989年发生在美国加州的7.1级地震当天，出现了增强了20倍的地磁异常信号，地震发生前3小时，信号增强了60倍 [2] 。2000年7月1~15日发生在日本伊豆岛的火山和地震群震前40分钟出现了地磁特大异常 [3] 。2008年汶川地震前5小时，地震重灾区北川县就有指南针指向紊乱的现象，在地震发生前13分钟，北川中学有一个班的同学集体观察到指南针不规则转圈的现象 [4] 。这些现象表明突发地磁异常很可能预示着一场地震的发生。
本文采用在研究区域布设的地磁仪器采集到的地磁数据，综合利用深度学习的相关知识，建设起一个可以同时预测地震方位、震级的地震预测系统。在出现地磁异常时，可以利用该系统进行地震预判和预防，不仅为地震预测提供一个可行的研究方向与研究思路，同时也能实际应用于研究区域附近的防震减灾，具有很高的经济价值。
2. 卷积神经网络
2.1. 卷积神经网络结构
机器学习算法可以分为无监督学习和监督学习。其中无监督学习需要从样本中学习并获得规律；监督学习是在训练样本过程中，要对应确定的结果，然后对未知样本进行结果预测。监督学习可分为分类问题和回归问题，分类问题的标注结果为离散值，回归问题的标注结果为连续值。本文所研究的地震预测就是典型的监督回归问题。常用的有监督回归算法主要包括线性回归、支持向量机、梯度提升算法等。神经网络算法是一种常用的机器学习算法。
卷积神经网络(Convolutional Neural Network, CNN)是一种前馈神经网络，与其他神经网络模型相比，CNN 通过自主学习来学习优化过滤器进行特征筛选，而在传统算法中这些过滤器是人为设计的，因此在CNN中对于输入样本不需要有太多的先验知识或人为干预，这样卷积神经网络可以由浅入深进行全面的特征筛选，从而可以有效消除人工特征筛选的局限性，这是CNN的一个主要优势。近年来已经在图像和视频识别、图像分类、图像分割、自然语言处理 [5] 、脑机接口 [6] 等诸多领域中都有应用。
一个CNN由输入层，卷积层，池化层，全连接层和输出层构成(图1展示了一个经典CNN结构图)，除了输入和输出层外的其他层也叫隐含层。其中卷积层用来进行样本特征提取，不同的卷积核可以提取不同的特征；池化层用于降低卷积后特征图像的数据维数，全连接层利用卷积层和池化层提取到的全部特征进行预测。以CNN模型为基础，将多层卷积和多层池化结合产生新的网络模型，可提高网络结构的准确度，经典的卷积神经网络模型有GoogLeNet、AlexNet、VGGNet等 [7] 。
与循环神经网络、BP神经网络相比，CNN在图像识别领域具有以下三点优势：
1) 能够自动学习图像中的特征，使得它们能够很好地捕捉和描述图像的细节；
2) 可以处理大型图像数据集，并且能够从中快速学习；
3) 有很多开源的CNN模型可以供使用，比如VGG、ResNet和Inception等。
目前CNN已经在图像处理领域取得了广泛成果 [8] [9] ，本研究最终选择采用CNN开展后续异常图像的识别工作。
2.2. 评价指标
对神经网络模型的评价依赖于混淆矩阵的构建。
在分类模型的常用指标中，混淆矩阵是能够比较全面的反映模型的性能，并能够衍生出很多其他的指标。以表1的二分类结果混淆矩阵为例，其中：TP：真正例，实际为正预测为正；FP：假正例，实际为负但预测为正；FN：假反例，实际为正但预测为负；TN：真反例，实际为负预测为负。
根据混淆矩阵可以得到召回率和准确率的计算公式：
(1)
(2)
以上为二分类问题中的混淆矩阵、准确率、召回率的相关概念，在多分类评价问题中，只需要把每个类别单独视为“正”，所有其它类型视为“负”，就可以计算出预测结果中TP、FN、FP、TN值，再代入到公式(1)、(2)中就能计算出多分类的召回率和准确率。
3. 地磁异常数据的采集与分析方法
3.1. 地磁异常等值线图
本研究采用在四川凉山自治州地区某区域密集布设的一组地磁传感器网络的磁强数据，总共由25个磁强传感器节点，分布在一个14,000多平方公里的范围内。数据通过物联网实时采集传输到我们的数据中心。通过对磁强数据进行磁异常分析，可以得到整个监测区域的磁异常等值线图。
图2展示了两幅地磁异常等值线图，其中图2(a)是正常情况下的等值线图，图2(b)是震前出现超阈值地磁异常时的等值线图，图像尺寸均为480 × 480。
可以看到，相较于图2(a)展示出来的正常状态下等值线图，异常发生时的等值线图差异很明显，最主要的差别在于异常出现后，等值线图中会出现若干个极值区域，由于等值线图是基于站点异常绘制的，因此可以推断位于极值区域内及附近的站点出现了较大的地磁异常，因此，对地磁异常等值线图进行特征识别，是可以同时将站点位置、异常强度的影响考虑进来的。
(a) 
(b)
Figure 2. Contour map of geomagnetic anomalies
图2. 地磁异常等值线图
3.2. 地磁异常等值线图特征识别
在得到了地磁异常等值线图后，需要通过观察找到等值线图中与地震相关的地磁异常特征。在绘制和观察了2022年2月1日至2022年12月31日共计334天的等值线图之后，可以总结出四条可能与地震发生有关的特征：
1) 极值区域的面积；
2) 极值区域的异常值大小；
3) 极值区域的持续时间；
4) 极值区域的方位。
基于传统数字图像处理的方法可以分别处理和分析这四个特征，特征与对应变量间的关系见表2。
这四个特征分开统计很简单，但是要综合起来分析却很难，因为如果要像第三章构建模型的方式那样通过为四个特征设定阈值的方式构建本章的模型，那么一方面会割裂这四个特征的内在关联，另一方面需要设置四个阈值，受人为影响因素较大，使得模型移植到别的区域的困难增加。
考虑到以上研究过程中存在的困难以及卷积神经网络在图像特征识别方面的良好表现，本文选择采用卷积神经网络对地磁异常等值线图进行处理、识别和分类。将地震预测转化为一个分类问题，将对(震级，方位)的预测转化为具体的类别。
3.3. 网络模型设计
本研究需要同时输出方向和震级2个标签，根据这两个标签可以划分出5个类别，这是一个多分类问题。
VGGNet是在ImageNet竞赛上针对ImageNet数据集而提出的分类网络，图像大小为224 × 224，而本研究的所用到的等值线图像数据集大小为480 × 480，要比ImageNet数据集大，但是类别数要远远少于ImageNet数据集中的类别，本研究对其采取了一些改进措施以适配本研究的数据集，具体表现在：
1) 大大减少了所有卷积层和全连接层中的滤波器数量；
2) 将VGGNet原本的五个卷积块缩减为四个卷积块，去除最后一个全连接层，并添加了Dropout优化策略；
3) 输出层神经元个数为5，对应5个地震类别。
改进后的VGGNet网络仅有12层，我们将改进后的网络命名为VGG-12。VGG-12的网络结构具体结构见图3。
VGG-12由四个卷积块组成，每个卷积块中的卷积层的滤波器大小都为3 × 3的，每个卷积块中的卷积层中的滤波器数量皆相同，四个卷积块的滤波器数量分别为32，64，128，256。最后，通过两个全连接层至输出层，每个全连接层都添加了Dropout层，防止模型过拟合。输出层是一个完全连通的层，有5个神经元对应5种地震类别。
4. 实验结果
4.1. 实验环境
开展本研究的计算机的硬件环境为：CPU型号为Intel(R) Core(TM) i5-9300H，显卡型号为NVIDIA 1660 GeForce GTX Ti。软件环境为：操作系统为Windows 11操作系统，编程语言为Python，集成开发环境为PyCharm，主要三方库有tensorflow、numpy、pandas、PIL、scipy、scikit-learn。
4.2. 实验数据
地震数据来源于中国地震台网(https://news.ceic.ac.cn/)，按照纬度范围26˚N~31˚N，经度范围100˚E~104˚E，时间范围从2022年2月1日~2023年1月1日，有感地震震级Ms > 3.0为限定条件 [10] 查询地震记录，可以发现在监测区域周围400公里以内发生18次地震，监测区域的地磁数据都有反应。具体可见表3。
4.3. 地震数据集构建
由于期间一共发生了18次地震，次数较少，因此地震类别不宜太多，否则会导致落入每个类别中的地震次数过少，训练后的模型只在训练集上表现效果优秀，却难以推广到测试集中，即出现过拟合现象。
地震震级转化为类别比较简单，只要选取一个震级，保证大于这个震级的地震和小于这个震级的地震数量接近5:5即可，根据表3中的震级分布情况，将用于分类的震级设置为4级，18次地震中，有7次地震小于等于4级，11次地震大于4级。
地震方位转化为类别就要用到聚类算法，本研究所用的选择K均值聚类算法，距离采用的是欧式距离。将方位类别确定为2，利用SPSS软件确定这两个类别中心的经纬度为(29.68˚N，102.37˚E)和(27.19˚N，101.97˚E)，而研究区域中心的经纬度为(28.10˚N，102.16˚E)，从两个类别中心相对于研究区域中心的位置来看，可以将方位类别分为南北两个类别。
根据以上震级和方位的分类规则，可以得到4个地震类别，再加上一个未发生地震的类别，共5个类别，具体可见表4。
Table 4. Classification of contour map
表4. 等值线图像类别
按照这种方式进行分类之后，从数量上看，标记为未发生地震的日期为169，有地震发生的日期为165，二者比例接近5:5。
4.4. 数据增强
观察该数据集可以发现存在一个问题：样本量太少。一共334张图片还分成了5个类别，导致真正要识别的地震类别中的图像数量太少，这样训练出来的神经网络泛化能力一般不会太好。为此，本研究采用了诸如亮度增强、对比度增强、水平翻转等数据增强方法提高样本数量。考虑到地震方位是在南北方向，因此在纵向上旋转图像势必会对图像的南北特征产生影响，因此在图像增强时只是做了水平翻转而没有做任何纵向翻转和图像旋转操作。
这样进行了数据增强操作之后，样本总量从334增加至1649，扩大至原来的4.93倍，各个类别的样本数量也有所增加。
5. 分析与总结
5.1. 实验结果分析
本节将利用VGG-12模型开展训练和测试，并将其与基于传统数字图像处理的预测模型进行对比分析。
在深度学习中，数据集的选取对神经网络模型的训练有着直接影响，本研究所用数据集中共有1649张图片，类别共5类，分为0，1，2，3，4，各有833，437，96，169，114张图片，在1649张图像中，选取1319张用于训练，330张用于测试。
使用VGG-12模型进行训练和测试，并利用混淆矩阵、准确率、召回率来对网络性能进行评价，得到混淆矩阵见表5。
各个类别的准确率与召回率为：
Table 6. Evaluation of classification
表6. 分类评价
5.2. 总结
由表6可知，VGG-12对于除了类别2以外的类别识别准确度都达到了80%以上，模型在等值线图像分类上表现良好，其中类别0的识别度最高，其次是类别1，最低的是类别2，推测可能原因有2个：
分类为0和1的样本数量最多，分类为2的样本数量最少，因此模型对于类别0和1的特征识别度更全面和精确，对类别为2的特征识别度受限于样本数要模糊一些；
类别0代表无地震发生时的等值线图，这种等值线图一般不具备在4.1.3.节所总结的4个特征，因此整个图像的像素值变化比较平缓，因此特征数量更少，也更好识别；而类别1代表等级大于4级的地震，因此图像表现刚好与类别0的图像相反，其像素值变化剧烈，图像中的细节更多，因此在充分训练之后，突出特征被网络识别，因此也能更好识别；而类别2代表震级小于4级的地震，其图像中像素的变化居于无地震和强震之间，因此可能同时具有这两种类别图像的特征，如果其中代表强震或者无震的特征稍微多一些，就会被误分到无震和强震的类别中，从而导致该类的识别效果不好。
对单张图像的识别可以将该图像分类，由于震前地磁异常会连续多天出现，因此该次地震对应的震前地磁异常等值线图也会多天表现出相同的特征并被神经网络分到同一个类别中。因此在地磁异常出现后，可以用异常前后连续几天的分类结果来对这次地震的震级和方位的性质进行分析，进而也就实现了对地震的方位和震级实现预测。
6. 结论
本研究提出了一种以VGGNet为原型的改进后的VGG-12模型，该模型相较于原始VGGNet而言有更少的隐含层和参数，可以有效减少训练过程中的计算量。利用该模型对由分类后的地磁异常等值线图构成的数据集进行训练和预测。实验结果表明，该模型可以很好地提取图像中与地震有关的特征信息，以较高的准确率将图像归入它所对应的地震类别中，从而实现对地震方位和震级的预测。
但是受时间与技术的限制，研究过程还仍然存在一些不足之处：
1) 深度学习时样本数量太少，即使经过数据增强后，也只是勉强达到了神经网络模型训练所需的最低样本数量，而且样本数量太少使得要素类别也只能做简单划分；
2) 研究地震时只研究了研究区域附近约350 km及以内的地震，只考虑了距离因素对于研究区域的影响，实际上如果震级够大，地震影响的范围还能进一步提高。
基金项目
国家自然科学基金项目42074074。