基于机器学习的地表破裂自动识别和填图
——以2021年青海玛多MW7.4级地震为例
2023-03-15曾宪阳刘静王伟姚文倩吴静刘小利韩龙飞王文鑫邢宇堃杜瑞林杨绪前
曾宪阳, 刘静,*, 王伟, 姚文倩, 吴静, 刘小利,韩龙飞, 王文鑫, 邢宇堃, 杜瑞林, 杨绪前
1 中国地震局地质研究所地震动力学国家重点实验室, 北京 100029 2 天津大学表层地球系统科学研究院,地球系统科学学院, 天津 300072 3 杭州五代通讯与大数据研究院, 杭州 310012 4 中国地震局地震研究所地震大地测量重点实验室, 武汉 430071
0 引言
地表破裂带是地震破裂最直观的现象.大地震的地表破裂带精细填图是详细记录地表破裂带的几何形态、同震位移及其空间分布规律的重要工作(比如,Ambraseys and Tchalenko, 1970; Tchalenko and Ambraseys, 1970;Sharp et al., 1982; King and Nábělek, 1985; Spotila and Sieh, 1995; Treiman et al., 2002; Haeussler et al., 2004; Liu-Zeng et al., 2009; Quigley et al., 2012; Choi et al., 2018; 韩龙飞等, 2022),为探索地震发生机理、断层破裂过程及运动方式提供重要基础(Sieh et al., 1993; Haeussler et al., 2004; Oskin et al., 2012;Choi et al., 2018;Olson et al., 2019; DuRoss et al., 2020; Scholz, 2002, 2019).比如,有助于构建和检验大地震破裂扩展模型,断裂滑动模型以及动态破裂模拟和强地面运动建模等(Zachariasen and Sieh, 1995; Day and Ely,2002; Olsen et al., 2008; Xu et al., 2016; Lo Yi-Ching et al., 2018; Klinger et al., 2018; Hough et al., 2019; Pollitz et al., 2019; Goldberg et al., 2019; Lozos and Harris, 2020).此外,通过地表破裂观测到的同震位错也是构建地震位移与震级和破裂长度的经验关系(邓起东等,1992;Wells and Coppersmith,1994; Wesnousky, 2008)及研究古地震复发习性的重要参数(Schwartz and Coppersmith, 1984; Sieh, 1996; Liu-Zeng et al., 2006),对预测未来地震的破裂尺度提供重要的理论基础,为震后分析地震灾情、评估断层的地震危险性等工作提供重要依据(Mignan et al., 2015).
然而,对地震地表破裂的精细填图并不容易.在以往对青藏高原内部及周边大地震的震后地表破裂调查中,研究人员一般通过野外实地考察的方式近距离对地表破裂进行详细观测,虽然能够在局部范围内获取精度较高的地形地貌数据,但是会遗漏断层外更细微的变形,同时需要大量的人力物力,也受到地区的限制,很难快速对整个地表断裂进行精细填图.近年来,借助米级分辨率卫星光学影像或无人机影像对单次地震的地表破裂或多次实践累计地貌的特征进行较详细描述(Klinger et al., 2005; Xu et al., 2006; 刘静等,2013;Li et al., 2016; Bi et al., 2017; Ren and Zhang, 2019).由于地表破裂往往比较细碎,米级分辨率卫星影像不能对微小破裂成像,且对垂直位错很难进行分辨.所以获取的数据在精度和准确度上仍然不足再现地表破裂的复杂度和丰富程度.得益于无人机在地球科学领域的广泛普及和应用,利用无人机摄影测量技术对强震造成的同震地表破裂进行迅速和全面的影像采集是大势所趋(Pierece et al., 2020; Ponti et al., 2020; Liu-Zeng et al., 2022; 王文鑫等, 2022),此项技术极大推动了地表破裂带填图的完整性和精细程度,尤其是针对一些无法涉足的区域.地震发生后第一时间对同震地表破裂进行高精度无人机影像采集对于确定同震地表破裂的精确位置、范围、最大同震位移等都起到了前所未有的推动作用(Ponti et al., 2020; Liu-Zeng et al., 2022).
传统处理高分辨率无人机影像通常采取人工绘制地表破裂的方法,非常耗时,同时高度依赖个人经验,往往会造成漏画、过度解译等情况(Rodriguez Padilla et al., 2022; Mattéo et al., 2021).随着计算机技术的快速发展,机器学习(Machine learning)为处理这类遥感影像数据提供了更多的可能 (Bergen et al., 2019).近些年来,机器学习已经被逐步广泛地应用在地球科学领域的多个方面(Tingdahl and De Rooij, 2005; Meier et al., 2007; Adeli and Panakkat, 2009; Lary et al., 2016; 周永章等,2018;龚健雅和季顺平, 2018),如地震震级预测(Adeli and Panakkat, 2009)、大洋海底地形中海山的自动识别和填图(Valentine et al., 2013)、地形特征识别(Li et al., 2017, 2020; Du et al., 2019)、河道形态 (Beechie and Imaki, 2014)、滑坡陡坎结构(Dunning et al., 2009)、滑坡易发性研究 (Yilmaz, 2010;Brenning, 2005)、冰川流向识别(Smith et al., 2016)、冰缘地貌(Wainwright et al., 2015)、海岸地貌和海洋沉积物特性(Li et al., 2011; Martin et al., 2015)、地表变形和表面过程的复杂性和相互作用(Hu et al., 2021)等.目前机器学习在震后地表破裂自动识别上的研究相对较少,仅有个别学者对陆地断层和裂隙等构造要素的识别进行了尝试(Mattéo et al.,2021).因此,我们基于2021年5月22日玛多MW7.4地震发生后无人机航拍采集的地表破裂的高精度影像数据,采用机器学习的方法——基于卷积神经网络(CNN)(LeCun et al., 1999)的有监督深度学习算法,构建了地表破裂识别模型;通过迁移学习训练模型,使用模型进行运算,完成了对研究区地表破裂的快速识别;通过对比自动识别的结果与人工填图的结果,分析了当前模型的优点和不足,为后续利用该方法进行地表破裂的快速自动识别提供了新的优化策略.
1 研究背景
2021年5月22日,在青海省果洛藏族自治州玛多县发生了MW7.4地震,造成大量房屋倒塌、部分桥梁等基础设施受损.玛多地震发生在此前未被识别出来的江错断裂上,位于东昆仑断裂以南70 km处.震中位于(34.59°N,98.34°E),地震震源机制解为一条近EW走向、具有显著拉张分量的左旋走滑错动事件,发震断层走向N102°E,倾向S,倾角81°,滑动角-11°(中国地震台网中心,2021)(图1).震后我们对玛多地震进行野外调查,同时于5月24日—6月15日对地表破裂全段进行了高精度无人机影像(UAV)数据采集,发现玛多地震震后形成了长达158 km的地表破裂带 (王文鑫等,2022;姚文倩等,2022).地表破裂带西起鄂陵湖南(97.6135°E,34.7459°N),向E终止于昌麻河乡(99.2698°E,34.5102°N),依据其走向变化自西往东依次可划分为鄂陵湖段、野马滩北段、野马滩南段、黄河乡段、东哦段、昌玛河段和江错分支段等8段(姚文倩等,2022;刘小利等,2022;Liu-Zeng et al., 2022),主要分布在海拔4200~4600 m的山区、草原、沼泽、冲积扇和部分沙丘地区,主要沿N105°E走向展布,由一系列张裂缝、张剪裂缝、剪切裂缝以及挤压鼓包、挤压脊和裂陷等雁行状组合而成(李智敏等,2021;潘家伟等,2021;邵延秀等,2022),近地表植被主要以低矮的高寒沼泽灌丛—嵩草为主(李宁云,2018),因此植被对遥感图像上裂缝的识别影响较少.在影像上,这些地表破裂大都表现为深灰色细条带,其宽度和长度分别为几十厘米至 1~10 m和20~500 m(邵延秀等,2022;姚文倩等,2022;刘小利等, 2022).
图1 玛多地震构造背景图Fig.1 Seismogenic tectonic background map of the Madoi earthquake
此次玛多MW7.4地震是继汶川地震之后我国发生的震级最高的一次地震事件,同时也是获得地表破裂带厘米级无人机影像资料最为详细的地震.玛多地震地表破裂的高精度填图和特征研究为震后评估地震断层危害性等工作提供重要依据,对于地震破裂过程的精细研究具有重要的科学意义.
2 方法
2.1 无人机图像采集与正射影像(DOM)的生成
我们于5月28日—6月8日用纵横大鹏CW-15垂直起降固定翼无人机沿着地表破裂带进行了UAV图像采集(图2),无人机可以垂直起降,续航时间3 h,续航里程180 km,最高起降海拔>4500 m,抗风能力强(王文鑫等,2022).无人机在CW Commander平台上操作,配备JOUAV CA103全画幅正射相机(4200万像素),操作中使用CW Commander航空成像应用软件进行飞行规划和监测,以及采集航空照片,其中航空照片的航向、旁向重叠率均为80%,地面比例尺大多设置为1∶300(在地形复杂的破裂东端设置为1∶600).同时在研究区内设置临时基站并使用RTK (Real-time Kinematic)测量其地理坐标,定位精度为厘米级.数据采集后,使用Photoscan软件对UAV数据进行整理,生成数字正射影像(DOM).影像的质量取决于航拍区域地形、植被、作业时的风向、飞行高度、相邻照片的重叠度、特定区域照片数量等外界条件.尤其是在高原地区,地势起伏大,地形复杂,天气多变,不同飞行区域的无人机飞行高度不同,以上诸多因素导致影像分辨率有所差异(Liu-Zeng et al., 2022).最终生成的DOM(图1)空间分辨率为3~7 cm/pix,包含红、绿、蓝三个光谱带,其中研究区东端山地地形较为复杂,起伏度大,分辨率在5~7 cm/pix.之后使用ArcGIS10.2提取相应的山影图、坡度图等,通过对各要素图层进行拉伸、滤波、变换等处理,以便更好的识别地表破裂,提升数据解译能力.
图2 玛多地震地表破裂带无人机图像采集区Fig.2 UAV image acquisition area of the surface rupture zone of the Madoi earthquake
2.2 训练数据集准备
深度学习中的监督学习算法需要大量的标签数据来训练网络,这些数据泛化性质量的好坏决定着模型鲁棒性的强弱.因此我们使用ArcGIS在DOM上对地表裂缝采用勾画折线的形式进行了大量的手工标注,获得鄂陵湖段西端A区和野马滩北段B区(图2)约40 km2的无人机图像数据与之对应的约800条裂缝标注.为了减少遗漏、避免重复、方便复查,对整个测区采取统一的解译标志,将影像分为等大小格网,依次对每个格网采取1∶1000、1∶500、1∶250等比例尺显示,保障最小裂缝识别长度<0.1 m,按照移动窗口法拉网式、高密度人工解译最终获得沿总体走向长度约160 km、累积破裂长度超310 km的精细地表破裂数据.解译期间,各解译人员采取统一的解译标志,并及时沟通.并于9月在野外对标注数据进行现场确认以及调整,减少错标误标,减轻了这些因素带来的不利影响.数据现场确认后,将获取到的玛多震后遥感影像与人工标注信息清洗制作为PASCAL VOC2007格式的数据集.
PASCAL VOC2007数据集在检测、分割等任务上催生了很多优秀的计算机视觉模型.玛多震后遥感影像的地表裂缝识别可以划分为二分类的分割任务,一类为背景,另一类为地表裂缝.在计算机编程语言中,可以利用地理空间数据抽象库(GDAL)读取矢量数据中的标注信息,根据标注信息生成标签图片(图3).由于人工标注是在ArcGIS软件中采用勾画折线段的形式展示,而非多边形.因此利用标注数据直接生成的标签是由连续的单点像素组成,较为狭窄,无法完全覆盖地表裂缝.我们选择利用形态学中的膨胀算法处理,对标注的线膨胀40个像素元,把点序列形成的线转成裂缝区域,这样的处理方式可以使当前不能识别的数据变成可以利用的数据集,同时可以使标签更好的拟合裂缝的延展方向与空间展布.但是利用膨胀算法处理可能会导致图像裂缝中一些局部纹理细节丢失严重,从而影响某些像素点准确定位.最后切割为(512×512)大小的标签图片,形成2135份数据集,并把数据集按7∶2∶1的比例划分为训练集、测试集、验证集.
图3 地表裂缝数据集(a) 原图; (b) 标签图,为8位彩图.Fig.3 Surface fracture data(a) Original drawing; (b) Label map, 8-bit color map.
2.3 使用卷积神经网络识别裂缝区域
2.3.1 卷积神经网络模型结构
训练神经网络往往需要大量的训练集样本,而我们的训练集仅仅只有约1495份,在训练集较少的情况下易产生过拟合的问题.因此我们使用卷积神经网络VGG-Unet作为识别裂缝区域的深度学习网络,U-net采用了跳跃连接(skip-connection)结构对编码器和解码器中的各级分辨率特征进行连接,这种结构能使模型更好的利用图像的低级纹理特征,该方法曾在数据集样本较少的医疗影像结构分割任务中有很好的表现.本文使用VGG-16图像分类网络的图像特征提取模块作为U-net的主干特征提取器.VGG-16是一个经典的图像分类网络,在使用中舍弃VGG中的分类模块,在保留了高效的特征提取能力的同时,避免了VGG庞大的分类模块对整体性能的影响.同时考虑到裂缝提取是一个纹理特征提取任务,我们使用U-net对提取到的特征进行解码分类,U-net的优势在于能很好兼顾纹理特征和语义特征,对于纹理类的目标识别有很好的效果(图4).
图4 基于VGG-16模型的U-net Encoder模块Fig.4 U-net Encoder module based on VGG-16 model
首先,我们将切分后的RGB图像输入模型(图5),图像经由编码器模块(Encoder Module)的卷积层(Conv)滤波和池化层(Pooling)下采样得到不同尺寸的特征图,这些特征图被保存在Features Map中,小尺寸的特征图会在Decoder Module中经过上采样层转换成与较大尺寸的特征图相同的大小,然后与较大尺寸的特征图进行特征融合.这种做法可以让网络在小尺寸的特征图中通过高级语义特征,以较大的感受野搜索裂缝区域,同时融合大尺寸的特征图中的低级纹理特征,补全裂缝区域的细节.最终,融合后的特征通过一个Conv 2d层进行特征压缩,将特征压缩为两个通道,用于分别表示每个像素点分别为裂缝或非裂缝的权重.
图5 基于上采样和skip-connection的U-net解码器模块Fig.5 U-net decoder module based on upsampling and skip-connection
2.3.2 评价标准
(1)
2.3.3 训练
预训练已经在许多实验中被证明是一种有效提高训练结果稳定性和加快网络参数拟合的方法.对于Encoder部分,我们使用经由ImageNet数据集训练的VGG网络,在舍弃全连接层参数后作为迁移学习预训练网络,而对于Decoder部分,我们选择随机初始化参数重新训练.为了防止随机参数Decoder产生的较大LCE值使网络反向传播时对预训练网络参数造成破坏,在训练的前50个Epoch,我们锁定Encoder模块的参数,训练Decoder使其获得基本的识别能力以抑制LCE值.在该过程中,由于Adam优化器具有计算高效、占内存少、收敛快等优点,我们使用Adam优化器、批数据量8条、学习率0.0001和学习率下降率0.92作为训练参数.在第50个Epoch后,我们解锁Encoder模块的参数,对网络整体参数进行微调.此时由于Encoder需要参与网络参数更新,导致网络参数数据增加,同时因显存容量有限,我们将批数据量改为4条,其余训练参数不变.
2.4 后处理
在通过深度学习网络识别到裂缝区域后,为了将预测结果放进ArcGIS中进一步观察裂缝分布与走向,可在识别到的裂缝区域内对裂缝进行近一步精细刻画.注意到裂缝所在位置像素灰度值快速变化的特征,我们使用Canny算法在裂缝区域中提取裂缝实体,该方法主要分为以下四个步骤:(1) 首先利用高斯滤波,将图像转化为灰度图并使用高斯算子对灰度图像进行滤波,以消除高斯噪声,防止高斯噪声在后续过程中被误识别为裂缝. (2) 使用一阶离散微分算子在滤波后的灰度图像x和y方向分别计算近似梯度,用该值表示图像在该点x和y方向上灰度值变化的剧烈程度.为了获得清晰、明确的裂缝,我们需要找到图像中灰度变化最快的方向,即近似梯度最大的方向.为此,需要通过x方向的近似梯度Gx,y方向上的近似梯度Gy,来计算梯度幅值G和梯度方向θ,计算公式如下:
(2)
(3)
在获得幅值G和梯度方向θ后,将θ离散为8个方向,对应每个像素周围8个像素方向.并根据离散后θ指向的方向,找到该方向幅值G的最大值,最大值所在的像素点即判断为候选裂缝.(3) 最后为了排除一些灰度变化不明显的纹理噪声,同时保证裂缝的连续和裂缝末端的完整,我们采用双阈值的滞后阈值处理进行较长边缘裂缝的链接.通过对梯度幅值G设置高门限G1为120、低门限G2为80.将候选裂缝进行检测,其中检测值大于G1的定为强裂缝;检测值在G2与G1之间,且与确定为强裂缝的像素邻接的定为弱裂缝;检测值小于G1的定为非裂缝.最终将强裂缝、与强裂缝邻接的弱裂缝作为识别结果.
为了能够在ArcGIS软件中观察到裂缝识别结果的几何位置信息,需要对裂缝识别结果的图片数据进行计算并转化为文档数据,最后结合投影文件信息,形成用点来描述裂缝几何位置的ESRI Shapefile文件.
2.5 评价指标
因为裂缝数据是以矢量的方式保存,直接判断裂缝识别的矢量结果缺乏客观性,因此我们将矢量数据栅格化后膨胀的裂缝标签来验证预测裂缝的准确度.我们将DOM上的裂缝标注和裂缝预测结果分别进行膨胀,两者的膨胀结果分别视作真实的裂缝区域和预测的裂缝区域.我们根据预测结果和实际结果的一致性将识别结果以像素级分为以下四类.TP(真正例):正确预测裂缝区域的像素数;FP(假正例):预测为裂缝区域但非真实裂缝的像素数;TN(真反例):正确预测非裂缝区域的像素数;FN(假反例):预测为非裂缝但实际是裂缝区域的像素数.在此基础上,我们计算查准率(Precision)和召回率(Recall),前者为预测为裂缝区域中预测正确的比重,后者为真实裂缝区域中预测正确的比重.查准率体现了模型捕获裂缝区域的能力,召回率体现了模型捕获裂缝区域的质量.查准率和召回率的计算方式如下:
(4)
(5)
为了减少非裂缝区域误识别的同时尽可能获取所有的裂缝区域,我们希望模型同时拥有较好的查准率指标和查全率指标,因此我们使用Dice系数指标作为最终全面评估识别效果的指标.Dice系数的计算方式如公式(6)所示:
(6)
在Dice系数计算公式中,由于准确率和召回率取值为0~1,因此无论其中一个指标有多好的表现,如果另一个指标表现过差都会使Dice系数指标变差.同时Dice系数指标在分子乘上系数2,保证Dice系数指标取值范围在0~1之间.
这些评价指标为我们的工作提供了重要的帮助,其中查准率和召回率为我们调整训练方法和修改超参数提供了重要的方向指导,Dice评分为评估我们的方法结果的总体准确性提供了一个定量的标定.
3 结果
3.1 实验结果
本实验中需要对113 GB的数据进行计算,因此要求实验设备具有强大的计算能力.在实验硬件搭配中,采用了1 TB固态硬盘储存数据与RTX 3060 12 GB显卡做计算力支持.我们采用随机梯度下降的方式训练网络,设置总Epoch为300.训练过程中平均每个Epoch耗时2分52秒,总共耗时12小时36分.为了衡量模型的性能,在每个Epoch分别计算出训练集小批次与验证集小批次的平均Dice系数指标及损失值(图6).
图6 图像分割评测指标Dice、损失值变化曲线Fig.6 Image segmentation evaluation index Dice,LCE value curve
从图6可以看出模型参数拟合良好,最优性能模型对验证集进行预测的Dice约为87.5%.选取在验证集上表现最优性能的模型参数(第263 Epoch),并在验证集上进行识别将标注区域与机器识别区域结果区域勾绘出来(图7c),紧接着对模型的识别结果进行预测区域和地面实况区域之间相似性数据统计(图7d),其中TP的像素点个数为79430,FN的像素点个数为3292,FP的像素点个数为13051,由此计算出Precision=0.85、Recall=0.96、Dice=0.9.
图7 模型识别结果(a) 数字正射影像; (b) 人工识别地表破裂,其中红线为人工标注的地表破裂; (c) 标注区域与机器识别区域对比,其中红线为标签标注区域, 蓝线为机器识别区域; (d) 机器识别结果,其中绿色区域为真正例(TP)、红色区域为假反例(FN)、蓝色区域为假正例(FP).Fig.7 Model recognition results(a) DOM image; (b) Artificially identified surface rupture map,where the red line is the artificially marked surface rupture; (c) Comparison of marked area and machine recognition area, where the red line is the labeling area, and the blue line is the machine recognition area; (d) Machine recognition result map,where the green area is TP, the red area is FN, the blue area is FP.
按照上述识别率统计分析方法,随机在不同区域的验证集中抽取10张图进行识别与数据统计,结果如表1所示,各验证集的Dice系数在0.7~0.91之间,这也说明了我们的结果鲁棒性较好.与DOM上的裂缝进行对比(图7)表明,大多数裂缝能够很好地被识别.
表1 验证集样例识别数据统计Table 1 Verification set sample identification data statistics
3.2 测试结果
选取B-2(19839×16384)的测试集区域(即未在训练集与验证集中出现的区域)进行地震裂缝识别,为了观察模型效果,对B-2区域进行了人工标注,并将识别结果与人工标注进行对比计算,得出Precision、Recall、Dice计算结果如图8.
图8 B-2区域识别结果数据Fig.8 B-2 test area identification map
根据图8中的数据,该方法对B-2区域的识别结果为Precision=0.417、Recall=0.663、Dice=0.512.为了更好地观察识别效果与地震裂缝的识别位置,将识别结果按不同颜色进行展现,可以更直观地观察混肴矩阵中TP、TN、FP、FN各个样例的分布(图9).由于误差不可避免,在这种分辨率高的图片进行人工标注工作时,会产生漏标、错标等错误因素,从而影响数据统计结果,使Precision、Recall、Dice的计算值变小.从图9(a、b)中可以看出蓝色FP区域为人工漏标区域,而模型预测结果为裂缝,说明模型具有一定的泛化性能.图9(c、d)中红色FN区域的裂缝纹理杂乱呈斑块状,模型错误预测为背景.因此,性能优良的模型依赖于良好数据集的建立,需要尽可能多、地质纹理种类多、标签正确的数据集.尽管我们期望模型在其他数据中的表现如同训练集上一样良好,但是显然,模型只能尽可能地趋近于人类期望结果,我们依旧可以根据检测结果,观察裂缝分布.
图9 B-2测试区域标签图与识别结果对比图(a,b,c,d) 分别为4个不同的区域,红线为标签标注区域,绿色区域为真正例(TP),红色区域为假反例(FN),蓝色区域为假正例(FP) .Fig.9 B-2 test area label map and identification map comparison chart(a,b,c,d) are four different areas , the red line is the labeling area, the green area is TP, the red area is FN, and the blue area is FP.
4 讨论
4.1 假反例分析
召回率(Recall)的大小受真正例(TP)、假反例(FN)的影响,其主要因素为真实标签的质量,此外高分辨率无人机影像每张约1 GB的储存容量,图片较大使得标注工作量也十分巨大,这需要我们在进行人工标注地表破裂时更加细致,并对标注进行清洗.
我们对A-4区选取的数据集识别(图10),从图10a中可看出FN(假反例)区域主要覆盖在一些纹理细小、长度短的裂缝上,图10b中FN(假反例)区域并没有覆盖在真实裂缝上,这与人工标注时将非裂缝标记为裂缝所导致的误差有关,图10c中FN(假反例)区域主要覆盖在一些纹理线路杂乱、呈斑块状的裂缝上.因此,裂缝纹理特征不明显或聚集成斑块状、真实标签标注不细致等因素,模型将裂缝预测为背景的概率会变大,从而导致FN数量增多,召回率偏低.
图10 A-4测试区域标签图与识别结果对比图(a,b,c) 分别为3个不同的区域,红线为标签标注区域,绿色区域为真正例(TP),红色区域为假反例(FN),蓝色区域为假正例(FP) .Fig.10 A-4 test area label map and identification map comparison chart(a,b,c) are three different areas, the red line is the labeling area, the green area is TP, the red area is FN, and the blue area is FP.
4.2 假正例分析
地表破裂的精确率(Precision)与假反例数量主要受到训练集样本影响,为了使模型具有更好的准确性与泛化性,往往需要准备尽量多、且包含不同纹理地形的训练数据集.相反,当训练数据集样例较少或样本不均衡时,模型不能学习到最佳的参数并获得更好的泛化能力,会使得假正例预测结果增多.
从图11a中可看到FP(假正例)区域覆盖在标签的周围,模型预测区域比标签区域广,该种情况可视为模型预测正确.而图11b中可看出FP(假正例)区域覆盖在真正裂缝上,这是由于人工漏标导致的误差,可通过细致勾绘标签减少误差.图11c中可看出FP(假正例)区域覆盖在类似于对比图中的裂缝纹理上且未在真实地裂缝周围,这种情况视为模型错误预测,可在训练数据集中放入一批类似于裂缝纹理的图片,标签设为全背景进行训练.显然,上述3种情况都是影响精确率数值变低与假正例样例增多的主要因素.
图11 B-3测试区域标签图与识别结果对比图(a,b,c) 分别为3个不同的区域,红线为标签标注区域,绿色区域为真正例(TP),红色区域为假反例(FN),蓝色区域为假正例(FP) .Fig.11 B-3 test area label map and identification map comparison chart(a,b,c) are three different areas, the red line is the labeling area, the green area is TP, the red area is FN, and the blue area is FP.
4.3 深度学习在地表破裂识别中的优势与局限
同震地表破裂带的精细填图可以为理解断裂带复杂几何结构、动态破裂过程和古地震复发习性等提供重要的基础信息(Oskin et al., 2012;Klinger et al., 2018; Choi et al., 2018; 刘小利等,2022;姚文倩等,2022).而随着海量高精度地形数据的便捷化日益获取,如何采取更高效可靠的(半)自动化定量处理和分析也逐渐引发了众多科研人员的关注和投入(Palamara et al., 2007;Walker et al., 2016;姚文倩等,2019;韩龙飞等,2019;Mattéo et al.,2021).相较于传统的人工识别和专家识别系统而言,本次研究基于深度学习的算法对玛多地震同震破裂实现了自动识别和提取,拓展了利用深度学习在同震地表破裂研究领域上的广泛应用.总体而言,深度学习在构造地貌自动识别研究中具备以下优势:(1)相较于传统的人工识别方法,基于深度学习的自动识别技术可以提升数十倍效率,且可以降低标注过程中漏检误检的概率.(2)深度学习的自动识别结果是基于严格统一的判别标准,相对于人工识别更多依赖于判读者的主观认知来说,后者难以标准一致,因而深度学习的结果具有较为稳定的识别标准和和更为优异的质量.(3)传统手动识别中要求研究人员具有一定的专业背景知识,相对而言,深度学习网络可以具有较高的自主学习能力,因此不需要过多强调研究人员的图像知识、地质知识等来设计优秀的特征模式.
然而需要注意的是,目前深度学习在地学领域的运用仍存在一些问题,主要体现在:(1)深度学习网络的推广性尚且有限.深度学习网络的识别能力是基于对已有数据的学习,当构造地貌类型发生其他变化时,需要重新收集相关的高精度地形数据对网络进行二次训练,需要追加一定的成本.(2)尽管深度学习网络总体上可以实现较高准确度的识别,但在较小的粒度上仍有不同程度的误识别和漏识别现象.因此,目前深度学习更适用于宏观构造地貌的精细刻画.为了克服以上问题,未来我们需要在以下几个方面进行改进和创新:(1)尽可能采用最先进的算法,尝试突破已经训练好的模型的应用局限性;(2)继续改进算法,利用大型图像处理器更加快速地处理覆盖较大区域的数据集,减少运算时间.(3)结合递归神经网络和CNN,利用多波段信息,达到更准确识别的效果.
4.4 展望
深度学习方法也可拓展到其他地貌和地表过程等定量研究中.近年来,移动摄影测量技术(SfM)以其方便快捷、低成本获得数公里到数十公里范围的厘米级三维地形数据和光学影像(魏占玉等, 2015; 毕海芸等, 2017),激光雷达扫描(LiDAR)能快速、高精度地获得更大范围的三维地形数据,这些都极大地推进了活动构造和构造地貌的量化研究(Hudnut et al., 2002; Gold et al., 2009; Cowgill et al., 2012; Oskin et al., 2012; 刘静等, 2013; 任治坤等, 2014).随着大规模LiDAR和SfM技术的应用,海量密集地貌点云和影像数据获得越来越便捷,数据呈现指数级增长.对于海量高精度地形描述数据的有效挖掘,一方面在三维数据可视化方面尝试应用虚拟现实技术研究目标区域的构造和地貌特征(Cowgill et al., 2012; Fowler et al., 2012; 王鹏等,2020),另一方面基于高分辨率地形数据对地表地貌和地质岩性的识别也逐渐从定性描述向定量分析转变.
以往研究人员通常以感性差异的裸眼经验进行人工地貌解译,或以简单的计算方法和指标(如地形起伏度和粗糙度)反映不同期次的地貌单元发育特征 (韩龙飞等,2019; Frankel and Dolan, 2007; Regmi et al., 2014),进行不同地貌面特征的参数表征和自动化提取后,开展半人工填图.机器学习为基于高精度海量数据的大范围地貌定量研究提供新的机会.深度学习网络更容易更快速地将遥感影像和DEM、坡度等其他多类型的数据集进行有效地挖掘,提高地貌识别性能,以帮助更好地了解地表过程(Du et al., 2019).
总之,随着算法和模型不断的改进和创新,机器学习会在地貌和地表过程的研究中有着更广泛的运用.机器学习技术的应用近年来已经涉及地学的各个领域,为地球深浅部过程研究提供了直接的手段,开辟了新视野.未来在算法、精度、应用场景和批处理便利性等方面还有更大的发展空间.
5 结论
我们在MW7.4 玛多地震发生后第一时间沿整个地表破裂带进行了高分辨率无人机影像数据采集,获取的数据集能够很好的识别地表破裂,为探索地震发生机理、断层破裂过程及运动方式提供基础信息.本文提出了一种基于机器学习算法,从高分辨率无人机影像数据中自动检测同震地表破裂的能力,采用卷积神经网络(CNN)的有监督深度学习算法构建了地表破裂识别模型,并应用于鄂陵湖段和野马滩北段约40 km2的无人机航片数据集,快速、高效地绘制了该段的地表破裂.尽管存在一些误报和漏识别,但与人工解译结果对比表明,深度学习算法与自动化处理系统相结合能够有效地绘制地震所产生的地表破裂,为震后分析地震灾情、评估地震断层危害性等工作提供重要依据,也为未来研究大地震地表破裂提供新思路.本次采用的卷积神经网络模型为后续利用深度学习进行地表破裂的快速自动识别提供了一种可能的优化策略.同时展示了机器学习在地震地质、地表过程和地貌等定量研究中的巨大优势和广阔前景.
致谢感谢评审专家以及编委对于本文提出的问题和建议;感谢青海省地震局、成都纵横自动化技术有限公司在数据采集过程中给予的支持和帮助.在此一并表示感谢!