0 引言

地震勘探在地下地质构造研究、油气资源勘查等领域起着重要的作用。在实际地震资料采集过程中，有效信息通常与复杂的背景干扰混合在一起，因此需要处理之后地震资料才可用于地质解释。噪声压制是地震资料处理中的重要一环，目的是提高地震数据的信噪比，改善成像质量。而随机噪声是一种不规则噪声，频率、传播方向不固定，很难直接去除。因此，随机噪声压制的研究具有重要的意义。

目前，地震数据随机噪声压制方法大致可分为两类，即基于模型驱动和基于数据驱动。基于模型驱动的去噪方法有f-x域预测去噪^[1-3]、S变换去噪^[4-5]、曲波变换去噪^[6-7]、小波阈值去噪^[8-14]等。该类方法提出时间较早，以较为常用的小波阈值去噪方法为例，其发展至今已有二十余年。1995年，Donoho^[8]首次提出小波阈值去噪方法，奠定了小波变换在信号去噪领域的地位。之后，张宇等^[9]、王勇等^[10]利用小波阈值去噪方法完成了地震数据随机噪声的压制。在此基础上，不同学者对小波阈值去噪方法进行了不同程度的改进。付燕^[11]在小波阈值去噪基础上，对随机噪声控制的尺度一上的小波系数进行第二次多尺度小波变换，并将二次小波变换后的尺度一上的小波系数置零，完成了对随机噪声的压制。刘军等^[12]将广义交叉验证函数作为一种新的阈值函数，优化了小波阈值的选取过程。除了对小波阈值去噪方法进行改进外，还有学者将小波阈值去噪方法与其他理论相结合，对地震资料进行联合去噪。夏洪瑞等^[13]将小波阈值去噪与加权叠加理论相结合，提出了小波时空域变阈值去噪方法。赵迎等^[14]将小波阈值去噪与完备总体经验模态分解理论相结合，提出了一种相对保幅、有效的去噪方法。

不论是小波阈值去噪方法，还是其他基于模型驱动的去噪方法，其发展时间相对较长，技术逐渐趋于成熟。但是，压制噪声的初始条件及数学原理均是从频率差异、区域统计规律、振幅差异等某一特定角度出发，因此，方法本身存在一定的局限性。

近年来，随着计算机硬件水平的提高，基于数据驱动的去噪方法得到了较快发展。主要代表为不同深度学习框架下的噪声去除方法，如韩卫雪等^[15]提出了基于深度学习卷积神经网络的地震数据随机噪声去除算法，并通过对不同类型实际地震数据去噪测试，验证了算法的可行性和有效性；王钰清等^[16]改进了卷积神经网络的训练策略，提出了基于数据增广和卷积神经网络的地震随机噪声压制方法；李海山等^[17]将残差卷积神经网络应用于叠前随机噪声压制，取得了较好的效果。卷积神经网络主要依靠卷积层提取地震数据的纹理特征压制噪声，但当网络层深度增加时，会导致非凸目标函数产生局部最优解。另外，由于卷积神经网络自身结构相对简单，存在模型训练速度慢、对训练样本要求高等问题^[18]。生成对抗网络(Generative Adversarial Networks, GAN)^[19]在一定程度上解决了这些问题。生成对抗网络作为一种新的无监督学习算法，具有较强的自主学习和生成样本能力，降低了对训练样本的要求^[20-21]。Wang等^[22]、Dong等^[23]将GAN用于地震数据的噪声压制，取得了一定的效果。Radford等^[24]将卷积神经网络与GAN相结合，提出了深度卷积GAN，此网络具有更强的特征学习能力及生成能力，但网络训练不稳定。针对此问题，Zhu等^[25]在GAN的基础上提出循环一致性生成对抗网络(CycleGAN)，通过引入循环一致性损失强化训练过程的稳定性，从而更好地实现网络对抗学习的过程。

因此，本文提出一种基于CycleGAN的地震数据随机噪声压制方法。通过对简单模型、复杂模型的地震剖面及实际地震数据的测试，并与经典的小波阈值去噪方法对比，验证本文方法的可行性和有效性。

1 基本原理 1.1 生成对抗网络

GAN由Goodfellow等^[19]提出，是一个通过对抗过程估计深度学习中生成模型的新框架。GAN主框架由一个生成器G(Generator)和一个判别器D(Discriminator)组成(图 1)。

GAN的基本思想来源于博弈论中的一个经典问题——“零和博弈”。在GAN中，参与“零和博弈”的双方为生成器G和判别器D，主要通过对抗学习的方式训练、优化模型，达到纳什平衡^[26]，此时模型参数达到最优。GAN目标函数为

$ \min\limits_{G} \max\limits_{D} V(D, G)=E_{\boldsymbol{y} \sim P_{\text {data }}(\boldsymbol{y})}[\lg D(\boldsymbol{y})]+ \\ E_{\boldsymbol{x} \sim P_{\boldsymbol{x}}(\boldsymbol{x})}\{\lg [1-D(G(\boldsymbol{x}))]\} $

(1)

式中：V(D, G)表示关于G和D的价值函数；y为输入的真实数据；x为服从高斯分布的随机变量；P_x(x)为x的概率分布；P_data(y)为y的概率分布；x~P_x(x)表示随机变量服从P_x(x)的采样；y~P_data(y)表示真实数据服从P_data(y)的采样；D(y)、D[G(x)]分别表示将y、G(x)判别为真实数据的概率，其中G(x)为x经生成器G后生成的数据；E_x(·)表示在x采样下“·”的期望值。

目标函数可以是最大化判别器判定正确的可能性，同时训练生成器使lg{1-D[G(x)]}最小。

1.2 循环一致性生成对抗网络

CycleGAN的设计初衷是在不需要其他额外信息的前提下，将图像从源领域映射到目标领域。相较于传统GAN，CycleGAN主要有两点改进^[25]：

(1) 在CycleGAN里，输入数据可以为不成对的训练数据，无须建立训练数据间一对一的映射，即可实现输入数据到目标数据之间的转换；

(2) 引入循环一致性损失函数强化训练过程。

CycleGAN的基本结构是由两个GAN组成，并且构成一个环形网络(图 2a)。在CycleGAN内的两个GAN(图 2b、图 2c)共享两个生成器G_{X, Y}、G_{Y, X}(G_{X, Y}表示该生成器作用是将X域内的数据转换至Y域，下标“_{X, Y}”不具有实际物理意义，仅为一种命名方式，G_{Y, X}同理)，并且各自带一个判别器D_x、D_y。其中，为避免生成器G_{X, Y}、G_{Y, X}将域内的所有数据都转换为另一个域内的某一数据，CycleGAN使用循环一致性损失(Cycle Consistency Loss)作为约束，即x和y从源域转换到目标域后，同样也可以从目标域返回源域^[27]。

由于CycleGAN的结构发生变化，目标函数也随之发生变化。CycleGAN的目标函数由对抗损失和循环一致性损失两部分组成。

(1) 对抗损失。在CycleGAN中，对于GAN1，目标函数为

$ \operatorname{los}_{\mathrm{GAN} 1}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, D_{\boldsymbol{y}}, \boldsymbol{x}, \boldsymbol{y}\right)=E_{\boldsymbol{y} \sim P_{\text {data }}(\boldsymbol{y})}\left[\lg D_{\boldsymbol{y}}(\boldsymbol{y})\right]+ \\ E_{\boldsymbol{x} \sim P_{\boldsymbol{x}}(\boldsymbol{x})}\left\{\lg \left[1-D_{\boldsymbol{y}}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}(\boldsymbol{x})\right)\right]\right\} $

(2)

同理，对于GAN2，目标函数为

$ \operatorname{loss}_{\mathrm{GAN} 2}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, D_{\boldsymbol{x}}, \boldsymbol{y}, \boldsymbol{x}\right)=E_{\boldsymbol{x }\sim P_{\text {data }}(\boldsymbol{x})}\left[\lg D_{\boldsymbol{x}}(\boldsymbol{x})\right]+ \\ E_{\boldsymbol{y} \sim P_{\boldsymbol{y}}(\boldsymbol{y})}\left\{\lg \left[1-D_{\boldsymbol{x}}\left(G_{\boldsymbol{Y}, \boldsymbol{X}}(\boldsymbol{y})\right)\right]\right\} $

(3)

(2) 循环一致性损失。在对GAN1、GAN2进行对抗训练时，仅使用对抗损失不能保证学习得到的映射会将一个输入x生成期望的输出y。为了进一步提高对抗训练的准确性，减少可能的映射函数空间，学习得到的映射函数G_{X, Y}(·)、G_{Y, X}(·)还应满足循环一致性条件，如图 2d所示，将X域(含噪声数据)内的数据x转换至Y域(去噪后数据)后，再转换为x′，即x→G_{X, Y}(x)→G_{Y, X}[G_{X, Y}(x)]→x′，使得x与x′接近，甚至相同，这样可以避免模型网络将所有X域内数据转换为Y域内某一固定数据；同理如图 2e，将Y域内的数据y转换至X域后，再转换为y′，即y→G_{Y, X}(y)→G_{X, Y}[G_{Y, X}(y)]→y′。循环一致性损失为

$ \operatorname{loss}_{\text {cycle }}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, G_{\boldsymbol{Y}, \boldsymbol{X}}\right)= \\ \mathrm{E}_{\boldsymbol{x} \sim P_{\mathrm{data}}(\boldsymbol{x})}\left\{\left\|G_{\boldsymbol{Y}, \boldsymbol{X}}\left[G_{\boldsymbol{X}, \boldsymbol{Y}}(\boldsymbol{x})\right]-\boldsymbol{x}\right\|_{1}\right\}+ \\ E_{\boldsymbol{y} \sim P_{\text {data }}(\boldsymbol{y})}\left\{\left\|G_{\boldsymbol{X}, \boldsymbol{Y}}\left[G_{\boldsymbol{Y}, \boldsymbol{X}}(\boldsymbol{y})\right]-\boldsymbol{y}\right\|_{1}\right\} $

(4)

式中‖·‖₁表示L₁范数。

循环一致性损失使输入的x、y与输出的x′、y′之间保持高度一致性，从而使生成器的学习更具有目标性，可防止网络对抗性学习的退化^[28]。

最终的损失函数由上述三部分组成，即

$ \operatorname{loss}_{\mathrm{all}}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, G_{\boldsymbol{Y}, \boldsymbol{X}}, D_{\boldsymbol{x}}, D_{\boldsymbol{y}}\right)=\operatorname{loss}_{\mathrm{GAN} 1}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, D_{\boldsymbol{y}}, \boldsymbol{x}, \boldsymbol{y}\right)+ \\ \operatorname{loss}_{\mathrm{GAN} 2}\left(G_{\boldsymbol{Y}, \boldsymbol{X}}, D_{\boldsymbol{x}}, \boldsymbol{y}, \boldsymbol{x}\right)+\lambda \operatorname{loss}_{\mathrm{cycle}}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, \mathrm{G}_{\boldsymbol{Y}, \boldsymbol{X}}\right) $

(5)

式中λ为控制系数(通常为10)，用来控制对抗损失和循环一致性损失在最终损失中的比重。最终的目标函数为

$ \min\limits_{G_{\boldsymbol{X}, \boldsymbol{Y}}, G_{\boldsymbol{Y}, \boldsymbol{X}}} \max\limits_{D_{\boldsymbol{x}}, D_{\boldsymbol{y}}} \operatorname{loss}_\text{all}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, G_{\boldsymbol{Y}, \boldsymbol{X}}, D_{\boldsymbol{x}}, D_{\boldsymbol{y}}\right)= \\ \operatorname{loss}_{\mathrm{GAN} 1}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, \mathrm{D}_{\boldsymbol{y}}, \boldsymbol{x}, \boldsymbol{y}\right)+\operatorname{loss}_{\mathrm{GAN} 2}\left(G_{\boldsymbol{Y}, \boldsymbol{X}}, D_{\boldsymbol{x}}, \boldsymbol{y}, \boldsymbol{x}\right)+ \\ \lambda \text { loss }_{\text {cycle }}\left(G_{\boldsymbol{X}, \boldsymbol{Y}}, \mathrm{G}_{\boldsymbol{Y}, \boldsymbol{X}}\right) $

(6)

2 CycleGAN构建

在CycleGAN中，本文用深度残差网络^[29](Deep Residual Network，ResNet)构建生成器网络，用马尔柯夫判别器^[30](PatchGAN)构建判别器网络。

2.1 生成器网络

在深度学习过程中，卷积神经网络通过不断堆叠卷积层进行特征提取，并且随着层数增大，特征也变得更丰富。但是，层数与准确率并不是单纯的正比关系，即不是神经网络中卷积层数越大越好，因为当达到一定层数时会出现网络退化现象。而ResNet可有效解决此问题。

本文利用ResNet网络构建CycleGAN的生成器(图 3)，与图像处理领域中的RGB三通道数据不同，地震数据类似于单通道灰度图。因此，网络的输入数据为256×256的矩阵，使用64个7×7大小的卷积核(Conv)进行卷积处理，经批量标准化(BN)、ReLU激活函数、下采样，进入由6个残差块组成的ResNet网络；然后经过上采样、64个7×7大小的卷积核进行卷积处理及Tanh激活函数，输出数据为256×256的矩阵。

2.2 判别器网络

普通GAN判别器输出一个实数，是对输入数据整体进行判别的结果；而PatchGAN输出一个矩阵，矩阵中的每个元素代表输入数据中某一部分的判别结果，因而具有更高的分辨率及准确性。

本文判别器(图 4)使用PatchGAN，与生成器网络输入数据相同，判别器网络的输入数据大小同样为256×256的矩阵，使用不同层深、卷积核大小均为4×4的卷积层处理，每次卷积处理后，需经批量标准化、斜率为0.2的Leaky ReLU激活函数处理，最终输出30×30的矩阵，该矩阵的每个元素值代表某部分输入数据为真实数据的概率。

3 数据测试

本文研究使用Dell台式工作站，Linux64位操作系统，版本为Ubuntu.18.04，GPU为NVIDIA. GTX1070。为验证CycleGAN对地震数据随机噪声的去除效果，分别利用理论数据(简单模型、复杂模型)和实际数据进行测试，同时对上述地震数据进行小波阈值去噪(综合对比应用不同小波及不同阈值的去噪结果，本文最终选用了sym3小波、固定式阈值、软阈值函数，小波分解层数为3层，见附录A)，通过结合不同的信号去噪效果评估指标(本文使用信噪比SNR、均方根误差RMSE两种评估方法，见附录B)，对比、评估地震剖面噪声压制效果；并抽取单道数据进行频谱分析，通过剖面细节进一步明确CycleGAN方法的噪声压制效果。

3.1 理论数据测试 3.1.1 简单模型数据测试

简单模型无噪声剖面如图 5a所示，添加一定程度的高斯白噪声得到含噪声剖面(图 5b，SNR为0.31)。剖面共250道，每道时长2000ms，包含水平、倾斜、不连续和孤立等多种同相轴。对于Cycle-GAN算法，需建立合适的训练集进行训练，在对模型进行去噪试验时，训练集样本分别选取2750幅无噪声和含随机噪声剖面(图 6)。其中，含噪声剖面与无噪声剖面无需一一对应，每幅训练集剖面中含有水平、倾斜、不连续和孤立等不同类型同相轴。

利用制作好的训练样本集，经历100Epoch，用时约40h，可得到训练好的网络结构，通过该网络对含噪剖面进行随机噪声压制测试。图 7为该网络的损失函数曲线。由图也可看出，本文网络中的各部分达到了一种相对平衡的状态，网络可以完成对剖面内随机噪声的压制。

从图 8可以看出，小波阈值去噪方法与Cycle-GAN去噪方法均对噪声进行了一定程度的压制，但后者去噪效果明显优于前者。由表 1可知，Cycle-GAN去噪后数据SNR为10.38，明显大于小波阈值，RMSE也小于小波阈值。

对比剖面细节，图 8b红框内同相轴在去噪后形态存在部分畸变，对比剖面内其他同相轴及训练集样本，可知是由训练集内与其类似剖面训练样本较少，导致网络对于该形态的同相轴去噪能力不够所致。抽取每幅剖面第60道数据进行频谱分析，结果如图 9所示。由图可见，原始无噪声地震数据的频谱较为光滑，而加入随机噪声之后的含噪声地震数据频谱发生了较大变化；通过小波阈值去噪后，频谱虽有所改善，但仍与原始数据频谱差别较大；而CycleGAN去噪后的数据频谱毛刺状干扰基本消除，并且频谱特征与原始地震剖面基本一致。

3.1.2 复杂模型数据测试

使用Overthrust模型(图 10a)的反射系数与主频为40Hz的Ricker子波进行褶积得到无噪声地震剖面(图 10b)，剖面共500道，每道时长640ms。对无噪声数据添加高斯白噪声得到含噪声数据(图 10c，SNR为0.65)。为了制作足够数量的训练样本，可将数据进行一系列的处理，如上下倒转、左右倒转、插值等。由于此噪声模型相对复杂，需要选取比简单模型更多的训练样本，本文选取10000幅含随机噪声剖面和10000幅无噪声剖面作为训练样本，图 11为部分训练样本示例。

利用训练好的网络，对含噪声数据进行去噪测试，结果如图 12所示。对比图 12a与图 12b可知，在图 10c中受随机噪声影响较为严重的同相轴(图 10c中红色箭头处)，经过两种去噪方法处理后情况均有所改善，但CycleGAN去噪效果明显优于小波阈值；对比图 12c与图 12d可见，图 12c中除随机噪声外，还存在部分同相轴，可知小波阈值去噪算法将部分有效信息去除，而CycleGAN去噪结果中无此现象。同样，由表 2也可以看出，本文去噪算法性能更佳。

对于Overthrust模型地震剖面的细节部分，如图 12a和图 12b中红框内区域，抽取第240道数据进行频谱(图 13)分析，在90~180Hz范围内，可发现CycleGAN去噪方法几乎完全消除了随机噪声；在0~90Hz范围内，CycleGAN去噪方法也能够很好地保留有效信号，最大限度地消除了随机噪声的影响。

3.2 实际数据测试

选用中国东部A工区实际地震数据进行测试。对于受噪声影响较小的数据进行常规噪声压制处理，即可获得成像较为清晰的地震剖面。从中选取10幅作为训练集中无噪声剖面，图 14为其中的两幅剖面，每幅剖面有250道，每道时长为1024ms。通过对无噪声剖面进行上下反转、左右反转及插值处理，扩充训练数据集数量，得到5000幅无噪声剖面。对无噪声剖面添加一定程度高斯白噪声，得到5000幅含噪声剖面。利用上述数据作为训练集进行网络训练。取该工区内另一随机噪声严重的剖面(图 15a，SNR为0.55)为测试数据，测试数据共250道，每道时长为1024ms，去噪结果如图 15b、图 15c所示。

对比图 15b、15c可以发现，由于原始数据在进行小波阈值去噪时，小波分解得到的小波系数无法对有效信号进行完美的表达，因此，虽然“雪花状”随机噪声得到压制，但其中部分连续同相轴变为多个不连续块状同相轴(如图 15a、15b内红框内部分)，有效信号损失较为严重。图 15d、图 15e分别为两种方法去除的噪声，其中，小波阈值去噪方法去除的噪声(图 15d)存在类似条状同相轴，相比之下，CycleGAN在压制随机噪声的同时，有效信息也得到了最大程度的保留。表 3为实际地震数据去噪效果对比，由表 3也可以看出，本文去噪方法性能更佳。

4 结论

本文提出一种适用于地震数据随机噪声压制的循环一致性生成对抗网络，该网络生成器由ResNet网络构成，有效避免了梯度爆炸(消失)及网络退化问题；判别器由PatchGAN构成，可以更好地对生成结果进行识别；同时，该网络在传统对抗损失的基础上引入循环一致性损失，提高了训练过程的稳定性。通过训练好的网络分别对简单模型数据、复杂模型数据及实际地震数据进行随机噪声压制测试，并以常规去噪算法(本文应用小波阈值去噪方法)作为对照组，利用去噪前后数据信噪比、均方根误差等参数，证明了本文构建的网络去噪方法优于小波阈值去噪方法；通过对比单道频谱分析结果进一步说明CycleGAN去噪效果优于常规去噪算法，验证了本文方法在地震数据随机噪声压制方面的可行性。

目前，深度学习算法在地震资料处理领域的应用才刚刚起步，还存在许多问题亟待加以改进。本研究下一步工作是：如何在维持当前去噪效果的同时，缩短网络训练时间，改善CycleGAN去噪的泛化能力，提高本文方法的实用性。

附录A 固定式阈值和软阈值函数计算公式

文中小波阈值去噪的阈值计算规则为固定式阈值，阈值函数为软阈值函数。

去噪过程中首先需要估计阈值，然后保留大于阈值的系数，舍弃小于阈值的系数。阈值选择过大则可能去除有效信号，过小则不能完全消除噪声，因此需要选择一个合适的阈值规则。文中固定式阈值为

$ \lambda=\sigma \sqrt{2 \ln N} $

(A-1)

式中：σ为噪声的标准方差；N为信号的长度。

在完成阈值计算后，需要阈值函数对小波系数进行处理，文中使用的软阈值函数计算公式为

$ \hat{\omega}(j, k)= \begin{cases}0 & |\omega(j, k)|<\lambda \\ |\omega(j, k)-\lambda| \operatorname{sgn}[\omega(j, k)] & |\omega(j, k)| \geqslant \lambda\end{cases} $

(A-2)

式中：$\hat{\omega}(j, k)$为阈值函数处理过的小波系数；ω(j, k)为含噪信号分解后的小波系数；sgn(·)为符号函数。

附录B 信噪比和均方根误差计算公式

信噪比(Signal to Noise Ratio，SNR)是有效信号能量与噪声能量之比，为振幅比的平方。假设有效信号为s(t)，去噪后信号为g(t)，则

$ \mathrm{SNR}=\frac{E_{\mathrm{s}}}{E_{\mathrm{n}}}=\frac{\sum\limits_{i=1}^{N} s^{2}(t)}{\sum\limits_{i=1}^{N}[g(t)-s(t)]^{2}} $

(B-1)

式中E_s、E_n分别为有效信号能量、噪声能量。通常认为，SNR越高，去噪后的信号残留噪声能量越小，去噪效果越好。

均方根误差(Root Mean Squared Error，RMSE)是均方误差的算术平方根，用来表示有效信号与去噪后信号之间的差异

$ \mathrm{RMSE}=\sqrt{\frac{1}{N} \sum\limits_{i=1}^{N}[g(t)-s(t)]^{2}} $

(B-2)

通常认为，RMSE越小，去噪效果越好。