由于地层的埋藏深度、构成岩石的矿物、孔隙充填的流体不同等因素,声波在不同井深地层介质中传播的时差就会有差异,通过声波时差异常可以识别储层[1,2]。声波时差测井曲线具有较高的纵向分辨率,但径向探测距离较近使得井间情况探测不清,地震记录的横向连续性较好,但因纵向分辨率受限而辨别不出较薄的地层。测井与地震资料联合反演可以充分结合二者的优点、取长补短,获得纵向分辨率高、横向连续性好的速度剖面[3]。传统井震联合反演方法是由模型驱动的,假设地层参数与地震响应之间存在先验的褶积算子,需要准确提取地震子波,褶积模型与实际情况的近似度决定了反演结果的准确度[4,5,6]。当声波时差和密度测井曲线不能反映储层与围岩的差异时,传统井震联合反演波阻抗将难以识别储层,测井曲线重构可以提高岩性识别和储层预测的能力[7,8,9]。由于测井曲线重构在非线性数学建模中遇到困难,因此通常是按照线性反演的。然而,测井曲线与地震记录之间的关系并非都是线性的。后来BP (Back Propagation)神经网络被用于测井曲线与地震记录之间的关系拟合,BP神经网络是非线性的,而且还是由数据驱动的,不用提取地震子波[10]。但是BP神经网络不能直接处理二维矩阵,需将多种信息首尾相连拼接成一维向量,而且BP神经网络采用全连接的方式,权值数量庞大。具有稀疏连接和权值共享两大特点的卷积神经网络(ConvolutionalNeuralNetworks,CNN)权值数量急剧减少,而且作为特征提取器的卷积核还可以从图像中提取多种特征。随着深度学习技术的突破,CNN如火如荼地被用于地球物理反演[11,12,13]。从测井曲线到地震记录的映射反演,不受褶积模型的限制,却仍旧依赖于初始模型[14]。从地震记录到测井曲线的映射反演,无需建立初始模型,但输入数据中只有地震记录一种信息,单由地震记录作为CNN的输入泛化能力较弱[15]。于是提出了一种基于机器学习的测井与地震属性联合反演方法,然后进行了模型测试和实例应用,以期获得纵向分辨率高、横向连续性好的速度剖面。

1、方法原理

基于机器学习的测井与地震属性联合反演方法分为两个阶段:第一阶段为从地震记录到测井曲线的映射反演;第二阶段为从测井曲线到地震记录的映射反演。每个阶段又分为两步:第一步为模型训练;第二步为外推反演。

基于机器学习的测井与地震属性联合反演流程如图1所示。图左是第一阶段从地震记录到测井曲线的映射反演,图左上是第一步的模型训练,图左下是第二步的外推反演;图右是第二阶段从测井曲线到地震记录的映射反演,图右上是第一步的正演映射算子逼近,图右下是第二步的初始模型迭代更新。第一阶段从地震记录到测井曲线的映射反演为第二阶段从测井曲线到地震记录的映射反演提供精度较高的初始模型,第二阶段从测井曲线到地震记录的映射反演通过迭代更新初始模型为第一阶段从地震记录到测井曲线的映射反演进一步提高精度。

1.1 从地震记录到测井曲线的映射反演

第一阶段从地震记录到测井曲线的映射反演是以测井曲线x为期望输出、井旁地震道y为输入,让CNN自动寻找从井旁地震道y到测井曲线x的反演映射算子f-1:y→x,即x=f-1(y)。

从地震记录到测井曲线的映射反演分为模型训练和外推反演两步。

1.1.1 模型训练

在模型训练这一步,输入层的井旁地震道在正向传播过程中与隐藏层和输出层的权值和偏置依次作用后得到模型输出,模型输出与期望输出之间的误差从输出层反向传播至毗邻输入层的隐藏层,每经一层就要更新该层的权值和偏置。反复训练,当输出测井曲线与真实测井曲线之间的误差下降到容忍范围后即可停止训练。这时的权值和偏置,即网络模型,就是要找的反演映射算子。

为了解决单一地震记录作为输入泛化能力较弱的问题,在第一阶段从地震记录到测井曲线的映射反演中,将测井曲线中提取的测井属性与地震记录中提取的地震属性构成的二维矩阵作为输入,以提高泛化能力。

图1 基于机器学习的测井与地震属性联合反演流程图  

在反演映射算子逼近过程中,训练误差随轮数的变化趋势如图2所示,当误差下降到容忍范围(<0.01)后停止训练,第1 414轮的整体训练误差是0.009 9。图3展示的是在反演映射算子逼近这一步WELL-09井的训练结果,训练误差是0.007 6。12口实测井轮流作为验证井,平均验证误差是0.016 5。图4展示的是在反演映射算子逼近这一步WELL-04井的验证结果,验证误差是0.015 5。训练结果和验证结果的均方根误差都较小,说明训练的模型可以用于下一步的外推反演。

1.1.2 外推反演

在外推反演这一步,将远井地震记录中提取的地震属性与井间内插外推的测井曲线中提取的测井属性构成的二维矩阵输入上一步已经训练好的网络模型,模型输出即为外推反演的测井曲线。

图2 在反演映射算子逼近这一步训练误差随轮数的变化趋势   

1.2 从测井曲线到地震记录的映射反演

第二阶段从测井曲线到地震记录的映射反演是以测井曲线x为输入、井旁地震道y为期望输出,让CNN自动搜索从测井曲线x到井旁地震道y的正演映射算子f:x→y,即y=f(x)。

从测井曲线到地震记录的映射反演分为正演映射算子逼近(模型训练)和初始模型迭代更新(外推反演)两步。

1.2.1 正演映射算子逼近

在正演映射算子逼近这一步,输入层的测井曲线在正向传播过程中与隐藏层和输出层的权值和偏置依次作用后得到模型输出,模型输出与期望输出之间的误差从输出层反向传播至毗邻输入层的隐藏层,每经一层就要更新该层的权值和偏置。反复训练,当合成地震记录与实际地震记录之间的误差下降到容忍范围后即可停止训练。此时的权值和偏置,即网络模型,就是要找的正演映射算子。

在正演映射算子逼近过程中,训练误差随轮数的变化趋势如图5所示,当误差下降到容忍范围(<0.08)后停止训练,第436轮的整体训练误差是0.0799。图6展示的是在正演映射算子逼近这一步WELL-02井对应的井旁地震道的训练结果,训练误差是0.0661。12口实测井对应的井旁地震道轮流作为验证道,平均验证误差是0.1257。图7展示的是在正演映射算子逼近这一步WELL-11井对应的井旁地震道的验证结果,验证误差是0.0897。训练结果和验证结果的均方根误差都较小,说明这一步逼近的正演映射算子可以用于下一步的初始模型迭代更新。

1.2.2 初始模型迭代更新

在初始模型迭代更新这一步,将第一阶段从地震记录到测井曲线的映射反演结果作为初始模型输入上一步已经训练好的网络模型,模型输出与期望输出之间的误差从输出层反向传播至输入层,沿途所经各层的权值和偏置都不再更新,更新的只是输入层的初始模型。反复迭代,直到合成地震记录与实际地震记录之间的误差下降到预设范围为止。这时可认为,最后一次迭代的模型即为反演的测井曲线。

图3 在反演映射算子逼近这一步WELL-09井的训练结果   

图4 在反演映射算子逼近这一步WELL-04井的验证结果   

图5 在正演映射算子逼近这一步训练误差随轮数的变化趋势   

在初始模型迭代更新过程中,单纯依靠合成地震记录与实际地震记录之间的误差来调整初始模型的更新方向颇具多解性,于是在代价函数中添加了一项正则化约束条件,以降低多解性。加入正则化约束条件的代价函数为:

式中,di0和yik分别为期望输出和第k次迭代的模型输出,即实际地震记录与第k次迭代的合成地震记录;M是数据点的个数;x0和xk分别为初始模型和第k次迭代的模型,最后一次迭代的模型即为反演的测井曲线;λ是约束权重系数。

图6 在正演映射算子逼近这一步WELL-02井的训练结果   

式中,α为学习率。

图8展示的是在初始模型迭代更新这一步WELL-11井的反演结果。从图8中可以看出,反演测井曲线与真实测井曲线吻合得较好。

图7 在正演映射算子逼近这一步WELL-11井的验证结果   

图8 在初始模型迭代更新这一步WELL-11井的反演结果   

2、模型测试

图9展示的是Marmousi-II速度模型,速度范围为1.5～5.0km·s-1,网格点数为241×241,x和z方向上的空间采样间隔分别为25m和12.5m,尺寸大小为6km×3km。从图9所示的Marmousi-II模型中沿着x方向抽取了61条速度曲线,以模拟声波时差变速度测井曲线。有着较高纵向分辨率的测井曲线,不能反映远井地带井间信息。由于实际地震子波随时空推移而发生畸变,因此用于合成地震记录的Ricker子波频率从20Hz按照给定的步长变化到50Hz。第121道的频率是20Hz,第1道和第241道的频率是50Hz,其余道的频率由中间向两边以0.25Hz的步长从20Hz到50Hz递增。利用Marmousi-II模型的速度和密度计算反射系数序列,然后不同频率的Ricker子波分别与各自对应的反射系数序列进行褶积运算以合成地震剖面,如图10所示。具有较好横向连续性的地震记录因为频带有限,所以从图10展示的地震剖面中只能依稀分辨出厚层及断层的大体轮廓,难以清晰分辨出薄层及断面的确切位置。图11和图12分别为从测井曲线中提取的测井属性和从地震记录中提取的地震属性。测井属性是由井间内插外推的测井曲线经平滑滤波后而得到的低频趋势,地震属性是地震记录的平均频率。图13为基于机器学习的测井与地震属性联合反演的速度剖面。可以看出,反演的速度剖面在纵向分辨率和横向连续性上均有了改善。图14为从x方向上1.475km处抽取的一条反演速度曲线与模型速度曲线进行比较。可以看出,反演速度曲线与模型速度曲线之间的误差较小。

图9 Marmousi-II模型   

图1 0 地震剖面   

3、实例应用

该实例来自于加拿大西部的Blackfoot地区,勘探的主要目标是Mannville群的Glauconitic段,储层深度约为1550m (1060ms)。储层中含海绿石的砂岩和页岩充填了河谷,形成了区域性的Mannville地层,河道砂体大约位于1090ms处。数据包括一个3D地震数据体和12口井的声波时差变速度测井曲线,图15展示的是一幅连井地震剖面。根据Rayleigh准则,通常认为,当地层厚度小于地震波长度的四分之一时则不可分辨,因此从地震剖面上难以分辨出薄层。图16为基于机器学习的测井与地震属性联合反演的速度剖面,比地震剖面的纵向分辨率有了很大的提高,图中的黑色曲线为真实测量的声波时差变速度测井曲线,上下两个矩形所框处分别是薄储层和河道砂体,约1090ms处的河道砂体和约1060ms处的薄储层得到了更好地辨认。图17是将WELL-04井作为验证井,其余井作为训练井的验证结果,蓝实线为实测曲线,红虚线为预测曲线,预测曲线与实测曲线之间的吻合程度较高。

图11 测井属性   

图12 地震属性   

图13 测井与地震属性联合反演的速度剖面  

图14 反演速度曲线与模型速度曲线对比  

图15 Blackfoot地区的地震剖面   

图16 Blackfoot地区的测井与地震属性联合反演的速度剖面   

图17 Blackfoot地区的反演速度曲线与实测速度曲线比较   

4、结论

Marmousi-II模型测试和Blackfoot地区的实例应用结果表明,测井曲线纵向分辨率高但横向连续性差、井间情况不易洞明,地震记录横向连续性好却纵向分辨率低、薄层难以分辨,测井与地震联合反演兼有测井和地震的优点,可以获得纵向分辨率高、横向连续性好的地层参数剖面。基于模型驱动的井震联合反演方法,假设地层参数与地震响应之间存在先验的确定性映射算子,而基于数据驱动的井震联合反演方法,则既不用褶积模型,也无需提取地震子波。从地震记录到测井曲线的映射反演,将测井属性与地震属性构成的二维矩阵替代之前的只有地震记录的一维向量后作为CNN的输入,提高了泛化能力,并把反演结果作为下一步的初始模型。从测井曲线到地震记录的映射反演,通过在代价函数中加入正则化约束条件后,降低了以往损失函数中只有地震记录约束下的迭代反演所具有的多解性。

参考文献:

[1]龚爱华.应用阵列声波测井识别苏里格致密气层的方法研究——以盒8段为例[D].青岛:中国石油大学(华东), 2017.

[2]黄蓉. LWF存储式声波测井数据处理与评价研究[D].西安:西安石油大学, 2016.

[4]张玉芬,罗延钟,凌峰.井震联合地震道多尺度反演[J].地球科学——中国地质大学学报, 2001, 26(5):533-537.

[5]张玉芬,李套山,刘小喜,等.多参数约束高分辨率处理方法在东濮凹陷中的应用[J].石油物探, 1999, 38(1):115-120.

[6]张玉芬,刘春华.多参数约束反演及其应用[J].地球科学——中国地质大学学报, 1997, 22(2):219-222.

[7]王宁宁,赵军龙. D区煤系地层声波时差测井曲线重构方法研究[J].河北地质大学学报, 2021, 44(5):80-83.

[8] 怀杰,乔宝强.基于多元回归模型的拟声波时差测井曲线重构方法研究[J].铀矿地质, 2021, 37(3):500-505.

[9]宋维琪,陈伟.测井声波时差反演重构技术研究及应用[J].地震地质, 2009, 31(1):133-140.

[14]张繁昌,刘汉卿,钮学民,等.褶积神经网络高分辨率地震反演[J].石油地球物理勘探, 2014, 49(6):1165-1169.

基金资助:中国石化胜利油田博士后科研项目(YKB2115);

文章来源:安振芳.一种基于机器学习的测井与地震属性联合反演方法[J].河北地质大学学报,2024,47(03):51-58.