基于深层神经网络压制多次波

宋欢， 毛伟建* ， 唐欢欢

1 中国科学院精密测量科学与技术创新研究院计算与勘探地球物理研究中心， 武汉 430077 2 大地测量与地球动力学国家重点实验室， 武汉 430077

0 引言

在地震勘探中，多次波既可以看成是一种能量很强的干扰噪声，也可以看成是与一次波类似的有用信号.将多次波视为干扰噪声是因为多次波的存在容易混淆一次波的能量，甚至出现多次波同相轴掩盖一次波同相轴的情况，从而影响地震成像，使反射波形态发生畸变，导致地震资料的偏移成像效果降低，甚至导致对地震资料解释和地震介质构造的错误认识.因此，有效压制多次波已成为地震勘探中的一个重要问题.另一方面，从地震波场传播过程来看，多次波同一次波一样也是地下反射层的反射，是在地下界面反射或海水面反射一次以上的地震反射波，也蕴含着地层结构信息.而且，多次波与一次波相比，在地下行走的射线路径更长，覆盖区域更广，能照射到一次波无法照明到的陡倾角地层，因此能携带更加丰富的地下结构信息.已有研究表明，将多次波视作有用信号进行偏移成像时，其对应的成像范围较一次波的成像范围明显加宽，成像质量也明显提高，特别是能提高盐丘下部成像质量和照明度(Liu et al., 2011，2016；刘伊克等，2015，2018；刘学建和刘伊克，2016)，这对认清受多次波影响而成像效果不佳区域的地质构造特征具有重要的工程应用价值.但是，对多次波进行成像之前，通常需要对多次波进行有效预测，因此也涉及部分多次波压制内容.

目前，已有的多次波压制方法根据原理的不同大致可分为两类：滤波法和波动理论法.滤波法主要是利用多次波和一次波的差异性完成分离(Weglein et al., 2011).差异性主要体现在两方面.一是多次波具有周期性，利用这一特性可对多次波进行预测进而达到压制多次波的目的，如预测反褶积方法(Treitel, 1969; Morley and Claerbout, 1983).但值得注意的是，只有在零偏移距或近偏移距处，多次波才呈现周期性，所以这类方法会随着偏移距的增大而降低压制效果.二是多次波和一次波的动校正速度差异，利用这一速度差异可以将动校正(NMO)后的多次波和一次波变换到不同域进行分离，在变换域中多次波和一次波比在时空域里更容易分离.常见的方法有F-K滤波(Sengbush, 1983; 吴战培, 1995)、聚束滤波(胡天跃等, 2000, 2002)和Radon变换(Hampson, 1986; Foster and Mosher, 1992; Sacchi and Ulrych, 1995; Hargreaves et al., 2001; Trad et al., 2003；熊登等，2009; Lu, 2013; Li and Yue, 2017)等.滤波法的原理较为简单，在一次波和多次波差异明显的情况下可以达到较好的压制效果，而当差异性较小时，该方法的多次波压制效果就会受到很大的限制.波动理论法则是从波动方程出发，根据地震波的传播特点建立多次波模型实现对多次波的预测，然后从原始数据中将预测的多次波自适应相减得到衰减后的地震数据，所以波动理论法也称为预测相减法.当前发展较成熟的波动理论方法有波场外推法(Loewenthal et al., 1976; Wiggins, 1988; Jiao et al., 2002)、反馈迭代法(Berkhout, 2006; Berkhout and Verschuur, 2006; Verschuur and Berkhout, 1997)和逆散射级数法(Weglein et al., 1997; 陈小宏和刘华锋, 2012; 刘伊克等, 2014).波动理论法不需要对地震波场作先验假设，特别是反馈迭代法和逆散射级数法不需要任何附加信息就可以压制多次波，因此适用性广泛，但方法本身计算量大，过程复杂，且震源子波估计的准确性和近偏移距数据的缺失会影响该方法的压制效果.所以，滤波法和波动理论法都有其各自的优缺点(张广利等, 2016).此外，无论使用滤波法还是波动理论法压制多次波都涉及大规模科学计算.在日常生产中，随着宽方位、高密度地震数据采集方式的出现，往往需要对海量地震勘探数据进行多次波压制，导致无论使用滤波法还是波动理论法都面临着计算耗时长的问题.因此，需要发展一种新的技术来提高多次波压制效率.

近年来，得益于深度学习技术的发展、海量数据资源的积累和计算能力的增强，人工智能技术引发了学术界、工业界的广泛关注和高度重视(赵改善, 2019).人工智能就其本质而言，是对人的意识、思维的一种模拟.在人工智能领域最活跃的是机器学习，在机器学习中最活跃的是深度学习.深度学习因能学习样本数据的内在规律和特征，使机器能够像人一样具有分析学习能力，能够识别文字、图像和声音等数据，已成为人工智能或机器学习领域研究与应用的热点，被越来越多的学者应用在勘探地球物理中，涉及领域主要包括地震构造解释、地震地层解释、地震相识别、地震反演、数字化岩石物理分析等.在多次波压制领域，Siahkoohi等(2019)采用生成式对抗网络对自由表面多次波进行了压制.为了取得较好的多次波压制效果，在测试阶段，该方法需要Surface-Related Multiple Elimination (SRME)算法提供初步的多次波预测结果作为输入.

为了继续探索和挖掘深度学习在多次波压制研究领域的潜能，本文提出了一种基于深层神经网络的多次波压制技术，采用的深层神经网络是一种改进的具有卷积编码器和卷积解码器的U-net网络.U-net网络被Hinton和Salakhutdinov(2006)提出, 当时主要用于压缩图像和去噪.Ronneberger等(2015)首次将其用于分割医疗图像后，该网络从此广泛应用于图像分割领域.然后，本文在大量带标签训练集上采用监督学习的方式训练定义的深层神经网络.训练成功的网络模型具备分离多次波和一次波的能力，即将含多次波的地震数据输入训练成功的神经网络，可直接快速输出多次波压制后的地震数据，从而避免常规多次波压制方法涉及的大规模计算.最后，利用工业界模型数据多角度验证了本文方法的有效性.

1 深层神经网络

1.1 定义网络模型

基于深层神经网络实现多次波压制的过程中，网络模型的定义及训练是关键.本文在TensorFlow深度学习框架下定义网络模型，使用的模型是一种改进的具有卷积编码器和卷积解码器的U-net网络，其结构如图1所示.相对于全连接层编码器而言，卷积编码器使用卷积层来学习数据特征，能够避免训练参数过多、难以训练和过拟合的问题.而且卷积编码器可以学习数据特征的稀疏表达，将高维地震数据通过非线性映射得到低维特征，通过设置合适的损失函数，使得提取的低维特征仅来自一次波信号，从而压制地震数据中的多次波信号.卷积解码器可以将包含低维特征的信号映射回高维空间，重构多次波压制后的地震数据.

从图1可知，本文使用的深层神经网络由编码和解码两部分组成.其中，编码部分包含8个下采样模块，分别为sequential、sequential_1、…、sequential_7，用于提取地震数据中的一次波信号特征.解码部分包含8个上采样模块，分别为sequential_8、sequential_9、…、sequential_14 和conv2d_transpose_7，用于将编码部分提取的一次波低维特征还原解码到原图尺寸，最终得到多次波压制后的地震数据.为了更好地恢复原始地震数据中的一次波信息，在编码器和解码器之间采用跳跃式连接，即在每一级的上采样过程中，利用拼接层(Concatenate)将编码和解码过程中对应位置上的特征图按通道进行融合，使得解码器在进行上采样时能够保留更多高层特征图蕴含的高分辨率信息，从而提高图像分辨率.

图1 深层神经网络结构图Fig.1 The structure diagram of deep neural network

每个下采样模块由一个序贯模型(Sequential)代表，该Sequential模型由三个网络层线性堆叠而成，分别是卷积层、批归一化层(Batch Normalization, BN)和激活函数.卷积层采用的卷积方式是“same”，使用“same”方式需要对原数据进行零填充，以确保卷积前的数据与卷积后的数据具有相同大小.为了达到下采样的目的，卷积层使用的划动步长取值为2.激活函数能实现去线性化，增加模型对复杂问题的表达能力.常见的激活函数有tanh函数、sigmoid函数、ReLU函数.在下采样模块中，本文从两方面考虑选用ReLU函数作为激活函数.一，ReLU函数能解决神经网络训练中一个非常致命的问题：梯度消失.ReLU函数将小于0的部分直接置0，大于0的部分即为输入，如式(1)所示，它既实现了非线性变换，同时输入大于0的部分梯度为1.所以从本质上来说，ReLU解决梯度消失的原则是靠梯度等于1，从而避免连续相乘的结果衰减.二，根据ReLU的定义，信息只能在输入大于0的区域进行传播，这带来了另一个优点就是稀疏性，稀疏性对网络性能的提高有帮助.在卷积层和激活函数之间还安插了一个BN层，目的是为了防止输入或者网络中间某一层的输出处于激活函数梯度很小的区域，这会导致很慢的学习速率甚至陷入长时间的停滞.对训练的每一批数据进行归一化处理，即减去均值再除以方差之后，可将数据移到激活函数梯度变化敏感的中心区域，这不仅是一种对抗梯度消失的手段，同时能提高神经网络的训练效率.

(1)

8个上采样模块中的前7个同样由Sequential模型代表，不过，该Sequential模型由四个网络层线性堆叠而成，分别是反卷积层(Conv2DTranspose)、BN层、Dropout层和ReLU激活函数.最后一个上采样模块是一个自带tanh激活函数的反卷积层，它也是输出层.在编码部分，原始输入地震数据经历8次下采样后，分辨率依次缩小了1/2、1/4、1/8、1/16、1/32、1/64、1/128、1/256.对于最后一层的输出数据，需要进行256倍的上采样，以得到原图一样的大小.本文的上采样是通过反卷积实现的，采用的反卷积方式同样是“same”，使用的划动步长同样为2.除了BN层和ReLU函数之外，解码部分使用的Sequential模型还包含一个Dropout层.Dropout的直接作用是减少中间特征的数量，降低每个特征之间的关联性，从而增强神经网络模型的泛化能力.

1.2 训练网络模型

本文采用监督学习的方式训练神经网络.监督学习最重要的思想是，在标注好的训练集上，通过不断优化神经网络参数，使得网络模型输出的预测结果尽量接近真实结果.优化神经网络参数的过程就是神经网络的训练过程.为了取得一个好的训练效果，需要准备合适的数据集、定义合适的损失函数、选取合适的优化算法和超参数等，它们与模型是否能训练成功息息相关，并直接决定了模型质量.下面将一一对这些影响因素进行详细介绍.

1.2.1 准备数据集

为了训练一个能成功压制多次波的模型，首先需要提供大量具有代表性的带标签数据.带标签数据是指输入数据对应的真实输出是已知的.本文制作的带标签数据集，输入数据是原始的含多次波的地震数据，标签是与原始数据一一对应的不含多次波的地震数据.从带标签数据集中选取部分数据作为训练集，输入定义好的神经网络，可训练出最优的模型参数，剩余数据可作为测试集，用于评估模型性能.不过，在训练神经网络之前，需要对标注好的数据集做进一步处理以用于深度学习算法，包括统一地震数据尺寸和归一化.

在设计神经网络时，会根据将要采用的神经网络和原始地震数据大小，定义Input层输入数据尺寸和Output层输出数据尺寸.一旦确定神经网络结构，用于训练和测试该网络的所有输入数据尺寸必须与该网络Input层的数据尺寸一致.如果不一致，需要做相应的调整.具体来说，如果输入的地震数据尺寸大于Input层数据尺寸，剪切地震数据中多出的行和列；如果输入的地震数据尺寸小于Input层数据尺寸，通过补零的方式增加行和列.

对数据集进行归一化处理最根本的目的是让模型更好更快的收敛.由于链式法则的累积效应，过大/过小的输入值会导致梯度爆炸或消失，所以将输入数据维护在一个大小合适的区间有助于加速训练.常见的数据区间有[0,1]和(-1,1)两种，具体采用哪种主要取决于网络结构.因为本文定义的神经网络结构对应的输出层使用了一个tanh激活函数，其对应的输出结果分布在(-1,1)之间，所以本文采用的归一化方法将输入数据归一化到(-1,1)之间，计算公式如下所示：

(2)

其中，x为样本数据，max(abs(x))为样本数据绝对值的最大值，xscale为样本数据归一化后的值.

俗话说：“苍蝇不叮无缝的蛋”，因为有机可乘，才有漏洞可钻。因此完善各种制度、机制是加强会计人员职业道德建设的保障。

1.2.2 定义损失函数

损失函数度量的是网络模型输出的预测值和真实值之间的误差，其中模型在训练数据集上的误差被称为训练误差，在测试数据集上的误差称为测试误差.一般来说训练误差是否变小代表着模型的训练过程是否收敛，而测试误差是否足够小直接与模型对未知样本预测能力的好坏相关.所以，损失函数不仅定义了优化问题，还对模型是否能训练成功起着至关重要的作用.

分类问题和回归问题是监督学习的两大种类.由于本文需要解决的是对具体数值的预测，属于回归问题.对于回归问题，最常用的损失函数是均方误差(Mean Squared Error，MSE).本文直接采用MSE作为损失函数，它的定义如下：

(3)

其中，yi为原始地震数据矩阵y中第i个元素的标签值，y′i为神经网络给出的第i个元素的预测值，n为该原始地震数据矩阵中所有的元素个数.

1.2.3 选取优化算法

神经网络的优化可分为两个阶段：第一阶段先通过神经网络计算得到预测值，并根据损失函数的定义计算误差；第二阶段通过后向传播算法计算损失函数对每个参数的梯度，然后根据梯度和学习率使用优化算法更新每个参数.最常见的优化算法是随机梯度下降法(Stochastic Gradient Descent, SGD)，及后来在其基础上发展起来的小批量梯度下降法、冲量梯度下降法等.SGD及其变种都以相同的学习率更新每个参数，这会导致一个问题：如果学习率太小，梯度很大的参数会收敛很慢；如果学习率太大，接近最优值的参数可能会不稳定.为了避免这个问题，本文使用Adam优化算法，该算法采用的是自适应学习率，即每个参数在优化的过程中能根据自身梯度大小的变化采用不同的学习率.不过，我们还是需要手动设置初始学习率.除了初始学习率，在训练模型的过程中还需要手动设置的参数有迭代次数(iterations)、代数(epochs)等.为了选取合适的超参数，需要在小数据集上有耐心地进行反复实验.

2 数值实验

为了证明利用深层神经网络模型压制多次波的想法是可行的、结果是可靠的，本文开展了一个数值实验，即利用含多次波和不含多次波的带标签训练集训练深层神经网络，然后利用训练的网络模型直接压制地震数据中的多次波.

2.1 实验数据介绍

使用的实验数据是Chevron Gulf of Mexico 2D Elastic Synthetic Seismic Benchmark(Chevron GOM 二维弹性合成地震数据)，该数据是在非均匀分层且包含多个盐丘的复杂地质条件下合成的，如图2 所示，复杂地质地震数据更有利于验证训练的网络模型压制多次波的能力.图2中出现的多种标识将在结果部分说明.

图2 合成Chevron GOM 二维弹性地震数据使用的速度模型Fig.2 The velocity model used by the Chevron GOM two-dimensional elastic-seismic-synthetic benchmark data

2.2 实验步骤及参数设置

为了证明训练的深层神经网络能代替常规方法压制多次波，并取得好的多次波压制效果.本文首先对大小为6400的带标签数据集进行预处理，包括按照定义的输入层数据尺寸统一所有的地震数据大小，和按照式(2)对地震数据进行归一化.然后在预处理后的包含320个带标签数据的训练集上通过采用Adam优化算法最小化目标函数式(3)来训练深层神经网络，获取的损失值曲线如图3所示.在训练过程中使用的初始学习率为2×10-4、迭代次数为320、代数为200.从图3可知，前100代，损失值以较快的速度从-9.5左右减小到-10.6左右，然后随着代数的增加损失值逐渐趋于平稳，在200代时损失值为-10.7左右，说明训练的神经网络收敛了.最后，利用测试集评估训练的神经网络压制多次波的能力.

图3 模型训练时对应的损失值曲线Fig.3 The training loss curve

3 结果与讨论

图4—7分别展示了利用训练的神经网络对测试集中四个具有代表性、且分布在不同地理位置的CMP道集进行多次波压制的结果.这四个CMP道集分布在图2中的四个红色圆点处，分别位于38750 m、54750 m、74000 m和95000 m.图2中四条虚线经过的区域分别代表了四个CMP道集对应的地下地质情况，其中第一条、第三条和第四条虚线都经过盐丘本身，第二条位于两个盐丘之间，说明选取的这四个CMP道集对应的地下地质情况都很复杂，通常它们会比简单地下地质条件下获取的CMP道集面临更大的多次波压制难度.此外，由于四个CMP道集具有不同的局部地下地质环境、且间隔较远，导致不同CMP道集中的一次波(或多次波)具有不同的形态特征，因此也会增加神经网络压制多次波的难度.

图4 第一个CMP道集对应的多次波压制结果(a) 含多次波的输入数据； (b) 不含多次波的真实结果； (c) 深层神经网络的预测结果； (d) 图4b与图4c之间的误差.Fig.4 The multiples-attenuation results of the first CMP gather(a) The input data with multiples； (b) The target data without multiples； (c) The predicted results of DNN； (d) The bias between Fig.4b and Fig.4c.

图5 第二个CMP道集对应的多次波压制结果(a) 含多次波的输入数据； (b) 不含多次波的真实结果； (c) 深层神经网络的预测结果； (d) 图5b与图5c之间的误差.Fig.5 The multiples-attenuation results of the second CMP gather(a) The input data with multiples； (b) The target data without multiples； (c) The predicted results of DNN； (d) The bias between Fig.5b and Fig.5c.

图7 第四个CMP道集对应的多次波压制结果(a) 含多次波的输入数据； (b) 不含多次波的真实结果； (c) 深层神经网络的预测结果； (d) 图7b与图7c之间的误差.Fig.7 The multiples-attenuation results of the fourth CMP gather(a) The input data with multiples； (b) The target data without multiples； (c) The predicted results of DNN； (d) The bias between Fig.7b and Fig.7c.

图4—7中的子图(a)是输入数据(含多次波的CMP道集)，子图(b)是真实结果(不含多次波的CMP道集)，子图(c)是利用训练的神经网络获取的预测结果(多次波压制后的CMP道集)，子图(d)是图(b)减去图(c)的误差图.比较图4—7中的图(c—d)可知，四组预测结果与其对应的真实结果都基本吻合，说明训练的神经网络能在保留一次波信号的同时很好地压制复杂地质地震数据中的多次波，而且对多次波形态特征不同的CMP道集也能取得好的多次波压制效果，说明训练的神经网络还具有很好的泛化能力.泛化是指：在数据同分布的假设下，从训练数据中学习到的特性，能有效预测没有见过但和训练数据来源于同一分布的数据.

为了更加充分地了解训练的神经网络压制多次波和重构一次波的能力，本文将以上四个CMP道集对应的真实结果与预测结果的零偏移距道进行比较，结果如图8所示.图中的蓝色和红色曲线分别是真实的和预测的零偏移距道的幅值随时间采样点的变化曲线.选取的时间采样点数从800到2000，对应的时间范围从3.2 s到8 s.在这个时间范围内，信号的幅度较强.从图8可知，整个时间段的比较结果可分为两部分，从采样点数800到黑色箭头处，一次波能量较强，且预测值与真实值高度吻合，说明训练的神经网络在压制多次波时能完全重构该时间段的一次波信号；从黑色箭头处到采样点数2000，预测值与真实值吻合得较差些，但我们也能看出吻合趋势，特别是在一次波能量相对强的位置，说明训练的神经网络在压制多次波时能部分重构该时间段的一次波信号.图中的黑色箭头处是多次波首次出现的位置，从含多次波的CMP道集中可以看出，该位置以下多次波的能量较强，而一次波的能量大幅衰减且同相轴的连续性较差.本文提出的基于深层神经网络压制多次波的方法通过学习有别于多次波的一次波信号特征，从原始地震数据中重构一次波.但是，当能量较强的多次波和能量较弱的一次波重叠时，一次波的特征不明显，使得一次波的特征学习不到位，最终导致训练的网络难以完全重构这部分一次波，所以图8中多次波较强的地方一次波有一些损伤，导致预测的一次波与真实的一次波难以高度吻合.

图8 零偏移距道的比较结果(a) 第一个CMP道集； (b) 第二个CMP道集； (c) 第三个CMP道集； (d) 第四个CMP道集.Fig.8 The zero-offset trace comparison(a) The first CMP gather； (b) The second CMP gather； (c) The third CMP gather； (d) The fourth CMP gather.

更进一步，本文比较了多次波压制前和压制后的零偏移距剖面，该剖面由1680个CMP道集的零偏移距道组成，结果如图9所示.该零偏移距剖面的时间范围从2.4 s到8 s，覆盖距离是2.1 km，从2.9 km到5 km，具体位置如图2中的红色方框所示，该方框包含一个完整的盐丘，有利于验证训练的神经网络压制复杂地质地震数据多次波的能力.图9a是多次波压制前的零偏移距剖面，与不含多次波的真实的零偏移距剖面(图9c)进行对比可知，图9a中的白色箭头处存在明显的幅度较强的多次波.利用训练的神经网络对零偏移距剖面中的多次波进行压制，结果如图9b所示.从图9b可知，训练的神经网络能很好地压制白色箭头处的多次波，且保幅重构大部分一次波能量，只有小部分能量较弱、连续性较差的深层一次波能量损失了，损失原因在图8说明中有介绍.

图9 多次波压制前和压制后的零偏移距剖面(a) 多次波压制前的零偏移距剖面； (b) 利用训练的深层神经网络压制多次波后的零偏移距剖面； (c) 不含多次波的真实的零偏移距剖面.Fig.9 The zero-offset profiles before and after multiples attenuation(a) The zero-offset profile before multiples attenuation； (b) The predicted zero-offset profile obtained by DNN after multiples attenuation； (c) The target zero-offset profile without multiples.

此外，在内存为8 GB、CPU为E5-1603的普通台式计算机上可完成整个深层神经网络的训练，耗时24 h 35 min.为了生成多次波压制后的零偏移距剖面(图9b)，在同样的计算环境下利用训练的神经网络对1680个CMP道集进行多次波压制处理，耗时仅22 min，完成一个CMP道集的多次波压制平均仅需要0.8 s，说明深层神经网络具有很快的多次波压制速度，而且数据量越大，利用深层神经网络压制多次波的时间优势越明显，适用于三维地震数据的多次波压制.不过，利用深层神经网络压制多次波的时间是三种时间的总和：准备标签数据、训练网络模型和压制多次波.一旦模型训练成功，压制多次波的时间很短，但准备标签数据和训练网络模型需要的时间相对较长.为了缩短时间，后续研究可利用CPU/GPU异构协同并行计算完成标签数据的准备和网络模型的训练.

本文方法可以用于实际数据，但不能将基于Chevron模型数据训练的深层神经网络直接用来处理实际数据，因为实际数据与Chevron数据不属于同一分布.在处理实际数据之前，需要利用实际数据制作带标签数据集，然后采用与Chevron数据相同的步骤重新训练网络.对于Chevron模型数据而言，含多次波和不含多次波的地震数据是已知的，容易制作带标签训练集.但对于实际数据而言，只有含多次波的地震数据是已知的.为了制作带标签数据集，需要采用常规方法压制实际数据中的多次波，将多次波压制后的地震数据视为标签数据，然后与原始含多次波的地震数据对应起来完成带标签训练集的制作.在处理实际数据时，为了节约时间，还可以利用迁移学习，在基于Chevron模型数据训练的深层神经网络上做适当微调，避免重头开始训练网络.利用迁移学习处理实际数据也是我们后续研究的重要方向.

4 结论

本文提出了一种基于深层神经网络的多次波压制技术，将含多次波的地震数据输入训练成功的神经网络，可直接输出多次波压制后的地震数据.通过比较多次波压制前、后的CMP道集、零偏移距道振幅和零偏移距剖面可知，本文定义的深层神经网络能较好地压制复杂地质地震数据中的多次波，且对测试集中具有不同多次波形态特征的地震数据也能进行多次波压制处理，具有较好的泛化能力.同时，训练成功的深层神经网络能在很短的时间内完成一个CMP道集的多次波压制，具有较高的多次波压制效率.文中的结论是根据工业界认可的模型数据得出，由于缺乏实测数据的验证，后续需要收集大量不同工区采集的实测地震数据进一步验证文中的结论.

致谢感谢Chevron提供的实验数据(Chevron Gulf of Mexico 2D Elastic Synthetic Seismic Benchmark).