地震资料的分辨率直接影响着地下储层的识别精度, 高分辨率的地震资料可以精细地反映地下的构造和岩性信息, 从而指导储层识别和井位部署等工作。如何提高地震资料的分辨率成为地震资料处理中的关键问题之一, 国内外学者相继提出了不同的处理技术来获得高分辨率的地震资料。

常规的高分辨率地震资料处理技术包括反褶积技术和反Q滤波技术。一方面, 反褶积技术是使用较为广泛的高分辨率处理方法, 包括预测反褶积、脉冲反褶积、同态反褶积等方法^[1-2]。在这些反褶积方法的基础上又衍生出了其它方法, 如陆文凯等^[3]利用神经网络技术求取反子波达到反褶积的目的; 蔡连芳等^[4]采用粒子群算法寻找峭度最大化从而实现地震盲反褶积; 王万里等^[5]建立了一个动态的褶积模型, 利用反射系数脉冲进行随机稀疏脉冲反褶积。但这些反褶积方法都需要提取准确的子波信息来求取反褶积算子, 而在实际应用中, 获得准确的子波信息十分困难, 所以这些方法均存在一定局限性。另一方面, 基于地层吸收衰减原理的反Q滤波技术通过补偿地震数据的振幅衰减和频率损失, 从而提高地震资料的分辨率。HALE^[6]最先提出反Q滤波法, 许多学者在其基础上进行改进并广泛应用于实际资料, 如基于短时傅里叶变换(STFT)、基于小波变换(WT)和基于广义S变换(GST)的时频域反Q滤波方法^[7-9]。但这些方法需要假设反射波各频率分量的能量在时间分布上有较高的相似性, 因此在实际应用中表现不佳。此外, 人们陆续提出了基于希尔伯特-黄变换(HHT)的高分辨率方法^[10]、基于压缩感知的高分辨率方法^[11]、基于经验模态分解(EMD)时频分析高分辨率方法^[12]与基于二阶自适应同步挤压S变换时变子波提取的高分辨率方法^[13]等地震资料高分辨率处理方法。但是这些常规的地震资料高分辨率处理方法存在子波提取困难、抗噪能力差、限制条件多、计算步骤复杂等局限性, 且随着地震勘探技术的迅速发展, 传统方法越来越难以满足日益增长的海量地震勘探数据处理的计算。因此, 提出一种智能化的地震资料高分辨率处理方法有着重大意义。

作为机器学习领域里的一个新方向, 深度学习技术一经提出便成为学术界和工业界的研究热点。目前, 深度学习技术在地球物理领域得到广泛应用, 如岩性和流体识别^[14-15]、地震去噪^[16-18]、断层检测^[19-20]、初至拾取^[21]、裂缝预测^[22]和地震资料反演^[23-24]等。此外, 对于地震资料的高分辨率处理, 国内外也有多位学者进行了一系列的研究, CHOI等^[25]随机生成高低频的地震数据进行神经网络训练, 并应用于楔状模型与实际地震数据来提高地震资料分辨率。孙永壮等^[26]利用U-Net网络训练伪反射系数模型数据, 实现了端到端的地震资料高分辨率处理, 即输入到神经网络的数据为低频的低分辨率地震数据, 输出为高频的高分辨率地震数据。以上端到端的方法训练过程都是在时域中完成的, 因而忽视了地震数据的频域信息, 导致网络的学习效果受限, 实际资料测试结果与真实高分辨率数据存在较大差异, 高频部分的数据无法得到有效补偿, 在性能上尚有较大的提升空间。

针对该问题, 我们提出了一种基于S变换与复值U-Net网络(S transform and complex-valued U-Net, STCVU-Net)的地震资料高分辨率处理方法。首先, 通过褶积模型构建高分辨率与低分辨率的单道地震数据, 使用S变换(ST)得到这些地震数据的时频谱作为训练集; 然后, 利用时频谱对设计的复值U-Net网络进行训练; 最后, 将所提出方法应用于模型资料和实际工区资料的高分辨率处理, 并将结果与传统端到端的U-Net网络方法处理结果进行对比, 证明了本文方法的可行性、有效性与适用性。

1 方法原理 1.1 S变换(ST)

S变换^[27]作为一种高时频分辨率的时频分析方法, 结合了小波变换和短时傅里叶变换的特点, 被广泛应用于地震资料处理中。对于地震数据x(t), 其S变换的定义为:

$ S(\tau, f)=\int_{-\infty}^{+\infty} x(t) g(t-\tau, f) \exp (-\mathrm{i} 2 \pi f t) \mathrm{d} t $

(1)

式中: t和f分别是时间和频率; τ是一个控制高斯窗在时间t上位置的平移因子; g(t－τ, f)是高斯窗函数, 具体公式如(2)式。

$ g(t-\tau, f)=\frac{|f|}{\sqrt{2 \pi}} \exp \left[\frac{-(t-\tau)^2 f^2}{2}\right] $

(2)

在S变换中, 高斯窗的高度与宽度会随着频率的变化而变化。S变换有一个重要的特性, 即S变换得到的时频谱对时间的积分等于其时域信号的频谱:

$ X(f)=\int_{-\infty}^{+\infty} S(\tau, f) \mathrm{d} \tau $

(3)

对(3)式两边同时进行傅里叶逆变换即可实现S逆变换, 即:

$x(t)=\int_{-\infty}^{+\infty}\left[\int_{-\infty}^{+\infty} S(\tau, f) \mathrm{d} \tau\right] \exp (\mathrm{i} 2 \pi f t) \mathrm{d} f $

(4)

利用S变换, 可以得到地震数据的时频谱并应用于后续的复值U-Net网络的训练与测试中。通过S逆变换的计算, 可以将网络学习后的高分辨率结果由时频域转换到时间域, 从而得到高分辨率的地震数据。

1.2 U-Net网络

构建U-Net网络进行地震数据的高分辨率处理。U-Net是一个典型的全卷积网络, 包含编码器与解码器两部分^[28]。编码器中主要利用卷积层和池化层进行特征的提取, 而解码器则利用转置卷积层与卷积层将提取到的特征进行还原, 得到原始尺寸的结果。同时在U-Net中还使用了拼接(concatenate)的方式将同一维度的特征进行拼接, 从而减少下采样带来的信息损失。

输入至U-Net网络中的训练数据大小为101×200的时频谱, 根据输入数据的大小设置了结构如图 1所示的U-Net网络。为了提取更深层次的特征信息, 相较于已提出的U-Net高分辨率处理方法^[25-26]设置了更多的通道数, 通道数从64开始、以2的倍数递增至1024。

在U-Net的编码器部分(图 1左侧)共分4个步骤进行特征提取。其中, 每一个步骤均含有两个卷积层, 卷积核的大小均为3×3, 步长为1。在每一个卷积层后设置一个ReLU激活函数, 在每一步的最后设置一个2×2的最大池化层进行下采样操作实现通道数的提升。解码器部分(图 1右侧)同样分为4个步骤, 在每一步中首先设置一个转置卷积层进行上采样操作来恢复数据的大小, 在每一个转置卷积层后同样设置了两个卷积核大小为3×3的卷积层与两个ReLU激活函数。需要注意的是, 在编码器部分设置了4个2×2的下采样层, 所以当输入数据的大小不为16的整数倍时, 数据无法被整除而需进行取整操作, 进而导致当使用相同卷积核大小的转置卷积层时输出数据与输入数据大小无法匹配, 数据无法恢复至原大小。因此, 为得到相同大小的输出数据, 分别设置4个转置卷积层的卷积核大小为2×3, 3×2, 3×3和3×2。在网络的最后采用一个卷积核大小为1×1的输出层进行数据的输出。

1.3 复值U-Net网络

传统的神经网络只能对实值数据进行训练, 但复值信号在地球物理等领域中广泛存在, 为了直接处理这些复值信号, 复值神经网络应运而生。近些年, 复值神经网络越来越受到国内外学者的关注, ARJOVSKY等^[29]通过多次实验证明使用包含复数的循环神经网络(RNN)相较于传统的RNN拥有更鲁棒、更稳定的性能。同时复值神经网络由于其映射能力强和适应性好的特点已被应用于实值分类、音频处理等。因此, 可以构建复值神经网络对复数形式的时频谱进行训练。

将TRABELSI等^[30]提出的深度复值模块拓展到整个网络结构中, 构建由两个实值的U-Net子网络A和B构成的复值U-Net网络, 表示为H=A+iB。时频谱数据可表示为w=x+iy, 其中, x是时频谱的实部; y是时频谱的虚部。复值U-Net网络的数据运算可表示为:

$ \boldsymbol{H} * \boldsymbol{w}=(\boldsymbol{A} * \boldsymbol{x}-\boldsymbol{B} * \boldsymbol{y})+\mathrm{i}(\boldsymbol{B} * \boldsymbol{x}+\boldsymbol{A} * \boldsymbol{y}) $

(5)

复值U-Net网络的计算过程如图 2所示。由图 2可见, 复值的运算被分解为2种不同的实值运算, 并且在运算过程中实部数据与虚部数据会同时输入到2个U-Net子网络, 2个子网络同时对实部数据与虚部数据进行训练。这种共享U-Net子网络的设计有助于整个网络更好地学习到时频谱的实部与虚部之间的关系, 使得预测的高分辨率结果更加准确。

设置合适的损失函数对提升网络的收敛速度和预测精度有着重要的意义。目前, 基于深度学习的高分辨率处理方法多采用均方误差(mean squared error, MSE)损失函数, 其定义为高分辨率标签与预测高分辨率结果的均方误差, 计算公式为:

$\mathrm{MSE}=\frac{1}{n} \sum\limits_i \operatorname{sum}\left[y_i-f\left(x_i\right)\right]^2 $

(6)

式中: n为训练样本数量; y_i为第i个数据的高分辨率标签; f(x_i)为第i个数据的高分辨率预测结果。

针对复值U-Net网络, 构建了复值MSE损失函数。该损失函数同时计算时频域的实部与虚部的均方误差, 即:

$ \begin{aligned} L= & \frac{1}{n} \sum\nolimits_i \operatorname{sum}\left[y_i^{\text {real }}-f\left(x_i^{\text {real }}\right)\right]^2+ \\ & \frac{1}{n} \sum\nolimits_i \operatorname{sum}\left[y_i^{\mathrm{imag}}-f\left(x_i^{\mathrm{imag}}\right)\right]^2 \end{aligned} $

(7)

式中: y_i^real与f(x_i^real)分别为第i个数据时频谱实部的高分辨率标签与预测结果; y_i^imag与f(x_i^imag)为第i个数据时频谱虚部的高分辨率标签与预测结果。

1.4 基于S变换与复值U-Net网络的高分辨率处理

传统的基于深度学习的高分辨率处理方法通过训练大量数据, 建立低分辨率数据与高分辨率数据之间的映射关系。但这些方法对不同分辨率的地震数据的特征提取与映射建立只在时域中进行, 因而忽视了数据在频域中的特征, 导致网络对高频信息学习不到位, 处理结果中高频信息出现一定缺失、准确性不高。为此, 结合S变换与复值U-Net网络的理论, 构建了一种基于S变换与复值U-Net网络的高分辨率处理方法。STCVU-Net的结构框架如图 3所示。首先对单道地震数据进行S变换得到地震数据的时频谱, 再将时频谱的实部与虚部输入至复值U-Net网络中得到高分辨率数据的实部和虚部, 对得到的结果做S逆变换得到单道高分辨率地震数据。图中红线表示实部的运算, 蓝线表示虚部的运算。

STCVU-Net方法结合深度学习数据驱动能力与时频分析方法的时频定位能力, 基于构建的复值U-Net网络对不同分辨率数据的时频谱进行训练, 使得网络在提取地震数据不同频率成分对应关系的同时获取该频率成分变化的发生时间, 从而更准确地对低分辨率数据到高分辨率数据的映射关系进行学习。

1.5 迁移学习

如果直接用合成数据训练得到的网络对实际工区数据进行测试时, 不同的地震数据及地震子波之间存在的差异往往会导致测试结果不够准确。实际工区的数据通常不足以重新训练网络, 在这种情况下, 引入迁移学习(transfer learning, TL)^[31]的方法可以使已训练的网络更适用于其它工区的数据。

不同的地震数据之间在存在着一定差异的同时仍存在一定的相关性, 迁移学习方法可以利用这种相关性显著减少地震数据之间的差异引起的不利影响^[32]。图 4为对复值U-Net网络进行迁移学习的示意。图 4中, 上半部分为预训练阶段, 下半部分为使用新数据进行迁移的阶段, 相同的颜色代表着网络层被冻结, 参数被固定。可以看到, 在微调过程中, 整个网络除了最后两个卷积层的参数可被调整外, 其余部分的参数均被固定。通过迁移学习的方法, 可以用少量新工区数据对预训练网络的部分参数进行更新, 从而使网络能够更准确地对新工区数据进行高分辨率处理。

2 训练与测试 2.1 网络训练与验证

首先, 利用褶积模型合成单道地震数据, 并对其进行S变换得到相应的时频谱作为训练集对网络进行训练。图 5为部分Marmousi2模型的反射系数, 大小为549×3400, 从中随机抽选50道反射系数, 使用这些反射系数与低频子波进行褶积得到低分辨率的单道地震数据作为输入数据, 而标签为反射系数与高频子波进行褶积得到的高分辨率数据。为了增强网络的泛化性, 在合成低分辨率数据时所用的子波为主频在5~20Hz内随机选取的雷克子波, 合成高分辨率数据所用的子波是主频为50Hz的雷克子波; 同时为了增强网络的抗噪性, 在低分辨率数据中随机添加不同方差的高斯噪声使信噪比(signal noise ratio, SNR)为-10dB至50dB, 信噪比的定义为:

$ \mathrm{SNR}=10 \log \frac{\|\boldsymbol{s}\|_2^2}{\|\boldsymbol{s}-\tilde{\boldsymbol{s}}\|_2^2} $

(8)

式中: $ \boldsymbol{s}$和$\tilde{\boldsymbol{s}} $分别表示真实地震数据与含噪地震数据。

通过不同的主频子波、不同信噪比的设置与滑动时窗的方法, 共构建了1800道、采样数为200的地震数据作为训练集, 另外构建了200道地震数据作为验证集。模型框架用Pytorch进行构建, 优化器为Adam自适应矩估计(adaptive moment estimation), 学习率被设定为0.0005, 批大小为4, 迭代次数为100, 整个训练与测试过程在一个GPU(NVIDIA RTX3090)上进行。

在训练结束后, 引入决定系数(coefficient of determination)来评估高分辨率处理结果的准确率。决定系数是表示回归模型中真实值和预测值的方差比的统计度量, 通常用R表示。R越接近1, 表示预测值越接近真实值, 预测效果越准确。R的公式为:

$ R=1-\frac{\sum\limits_i\left(y_i-\hat{y}_i\right)^2}{\sum\limits_i\left(y_i-\bar{y}_i\right)^2} $

(9)

式中: y_i是真实值; $\hat{y}_i $为真实值的平均值; $\bar{y}_i $是预测值。

为了说明所提出的STCVU-Net方法的优势, 额外构建了一个具有相同结构、超参数和训练集的时域端到端一维U-Net网络来进行对比。在该网络中, 单道地震数据直接输入到网络中进行训练。训练完成后, 在验证集上比较了时域端到端U-Net方法与STCVU-Net方法取得的R分别为0.9551和0.9636。从验证集上随机选取一道, 该道的低分辨地震数据、高分辨地震数据与两种网络的预测结果如图 6所示。图 7展示了该单道数据预测结果的时频谱和时频谱的实部与虚部。由图 7可见, 两种方法均能够对高频成分进行较好的恢复, 但STCVU-Net方法预测结果与高分辨率地震数据更吻合, 时频谱数值更精确。

2.2 模型数据测试

使用两个训练好的网络对整个低分辨率Marmousi2模型进行测试。图 8a是使用主频为20Hz雷克子波与反射系数通过褶积模型形成的低分辨率地震剖面; 图 8b是使用主频为50Hz雷克子波与反射系数进行褶积形成的高分辨率地震剖面; 图 8c和图 8d分别是使用传统U-Net方法和使用STCVU-Net方法得到的预测结果。同时计算了整个剖面上传统U-Net方法与STCVU-Net方法预测结果的R分别为0.9262与0.9522。可以看出, 比起低分辨率数据, STCVU-Net方法的预测结果的同相轴精细程度得到显著提升, 层间信息更加丰富, 局部小构造信息也清晰可辨, 地震资料的分辨率得到了显著提高。图 8e至图 8h分别为图 8a至图 8d中黑框区域的局部放大显示。从红色箭头所指部分可以看出, 相较于传统U-Net方法, STCVU-Net方法的预测结果同相轴更连续, 分辨率更高, 与实际高分辨率数据更吻合。图 9a和图 9b分别为2种方法与真实高分辨率数据的差值剖面。图 9c和图 9d为图 9a和图 9b中黑框区域的局部放大。由图 9可以看出, STCVU-Net方法的预测结果与实际高分辨率数据之间的差值更小, 预测结果更准确。进一步对比原始高低分辨率地震剖面与2种方法预测结果的频谱(图 10)可知, STCVU-Net方法显著提高了模型数据的主频和频带宽度, 其预测结果的频谱相较于U-Net方法更贴合高分辨率标签数据的频谱。

为了更直观地对两种方法的结果进行比较, 在预测结果中随机选取一道进行对比(图 11)。由图 11可以看出, U-Net方法的预测结果在一些地震道上与合成的高分辨率地震数据之间存在着较大误差, 未能完全恢复地震数据的高频信息。而STCVU-Net方法的预测结果误差更小、更加真实可信, 处理效果优于传统端到端U-Net方法。

为测试STCVU-Net方法的抗噪性, 在Marmousi2模型数据中随机加入了不同方差的高斯噪声, 使信噪比分别为2dB, 5dB, 10dB和20dB(图 12a至图 12d)。使用STCVU-Net方法对不同信噪比的数据进行高分辨率处理, 结果如图 12e至图 12h所示。

在不同信噪比条件下, STCVU-Net方法均能在进行高分辨率处理的同时对数据中的随机噪声进行有效地压制, 压制效果会随着信噪比的提高而提升。同时对使用不含噪声的训练数据训练得到的网络进行预测, 结果如图 12i至图 12l所示。由图 12可以看出, 当训练数据中含有一定噪声时, 不同信噪比下的高分辨率预测效果均得到了提高。因此, 在训练数据中加入适当噪声可以有效提高STCVU-Net方法的抗噪能力和预测精度。

2.3 实际数据测试

选取某工区的实际地震资料进行测试。为了评价高分辨率处理方法的效果, 将实际资料的低通滤波数据(图 13a)作为低分辨率数据进行输入, 将实际数据(图 13b)视为高分辨率数据进行对比验证。工区共有718条主测线, 901条联络测线, 此外还有18口井的测井资料。时间范围限定为0.736~1.134s, 采样间隔为2ms。

图 14为2种方法实际资料预测结果对比。由于地震子波之间存在差异, 当直接将Marmousi2合成数据训练所得网络用于实际资料测试时, 预测结果出现较大偏差(图 14a)。因此引入迁移学习方法, 使用测井资料中的反射系数信息与提取的实际数据高低频子波构建了15道高分辨率与低分辨率数据对(图 15), 并利用该数据对两个已训练网络进行微调后, 选取一个未参与微调过程的地震数据剖面进行测试, 结果如图 14b与图 14c所示, 图中黑色实线为测井反射系数与提取的高频子波进行褶积得到的合成记录。进一步计算相应预测结果的决定系数R分别为0.9134和0.9452。

可以看出, 2种方法均能提高地震数据的分辨率, 相比之下, 使用STCVU-Net方法得到的预测结果与真实数据的一致性更好, 不仅同相轴在横向上更连续, 而且构造特征更加清晰。进一步对比每个数据的频谱(图 16)可以看出, 未经微调的网络的预测结果的频谱与真实高分辨率数据频谱的吻合度较差; U-Net方法在一定程度上拓宽了频带, 但未能对高频部分完全进行补偿; 而STCVU-Net方法能对实际数据的高频部分进行最准确的补偿, 有效提升了主频, 拓宽了频带范围。

3 结论

本文结合深度学习的数据驱动能力与时频分析方法的时频定位能力, 提出了一种基于S变换与复值U-Net网络(STCVU-Net)的地震资料高分辨率处理方法, 使网络同时学习到地震数据中不同频率成分的变化与该频率成分变化的发生时间, 从而更准确地学习到低分辨率数据到高分辨率数据的映射关系, 有效且精确地进行高分辨率处理。模型数据与实际数据的测试结果表明, STCVU-Net方法相较于常规端到端的U-Net方法, 能够更准确地进行地震资料高分辨率处理, 并且STCVU-Net的高分辨率处理结果同相轴形态清晰, 高频信息丰富, 具有较大的应用潜力。