裂缝型油气藏广泛分布于陆相盆地, 油气储层面积大, 实际地质储量丰富, 是增产上储的重要领域^[1]。但裂缝型油气储层的形成过程复杂多样, 裂缝类型和尺度不尽相同, 非均质性强, 勘探难度极大^[2]。传统岩心定向分析^[3]和井壁成像技术分辨率高, 能够精准确定裂缝发育程度与类型^[4], 但仅代表单点信息, 对三维空间裂缝带的识别能力有限。面对上述问题, 许多学者利用纵横波和转换波之间的参数差异来预测裂缝, 但其识别目标一般为区域构造和与成岩构造有关的大裂缝带。地震属性分析通过提取地震波的反射特征(频率、振幅、相位等)进行裂缝带预测。越来越多的基于地震属性进行裂缝带预测的方法得到发展。相干属性^[5-6]、曲率分析技术^[7-8]、蚂蚁体追踪^[9-10]、反射强度属性^[11-12]以及基于这些属性发展的多尺度属性技术也用于不同尺度的裂缝综合预测^[13]。这些地震属性预测方法各有优势, 但受区域地质差异和测井资料等影响, 仅利用单一地震属性进行裂缝预测往往会产生多解性与假异常情况。因此, 基于测井信息约束的多属性融合技术是当前关注的重点, 但如何融合则是需要研究的关键问题。

近年来, 机器学习发展迅速, 在地球物理领域引起广泛关注。其中很多有代表性的算法, 例如支持向量机、神经网络、线性回归、朴素贝叶斯、随机森林等算法已被引入地震学领域, 在流体识别^[14]、速度谱自动拾取^[15]、地震数据噪声压制^[16-17]方面均有不错的应用效果。2006年, HUANG等^[18]提出了一种新型学习算法——极限学习机(extreme learning machine, ELM), 它基于单隐层前馈神经网络构建的算法, 仅有一个隐藏层的简单结构却具有不错的学习效果。与传统的单隐层前馈神经网络不同, 极限学习机的输入权重和偏置均随机设定, 再利用广义逆矩阵理论来计算输出权重。相比于其它一些机器学习算法, ELM在不失一定的学习精度的前提下, 具有快速学习、强泛化能力、人为干预较少的优点^[19]。潘华贤等^[20]研究了ELM进行储层渗透率预测的能力并与支持向量机(support vector machine, SVM)算法进行了对比, 结果表明ELM算法较SVM算法泛化能力强, 运算效率高。目前, 已有一些学者将机器学习应用于裂缝预测, 但大多都局限于测井数据的预测^[21]。李军等^[22]将测井数据与地震数据相结合, 利用BP神经网络综合预测裂缝发育带。受制于BP神经网络的固有缺陷^[23], 该方法学习效率较低、易陷入局部拟合且当预测区域深度过大时, 效果不理想。何健等^[24]将近似支持向量机(proximal support vector machine, PSVM)算法用于裂缝分类。作为支持向量机的改进算法, PSVM有着较高的运算效率和准确度。ELM因其较强的泛化能力, 运算效率优于SVM, 但与PSVM相比, 两者的分类效果还需要进一步研究。为此, 本文引入ELM算法, 比较ELM算法与PSVM算法基于测井数据和井旁道地震数据的分类效果, 再对实际研究区域进行裂缝综合预测。

1 方法原理 1.1 极限学习机

ELM是架构在单隐层前馈神经网络基础上的算法, 神经网络的输入权值和偏置均采取随机赋值的方式, 并在最小二乘准则的框架下, 利用Moore-Penrose广义逆计算输出权值。因此, 较传统的基于梯度下降学习理论的神经网络, ELM有快速收敛、不易陷入局部极值等优点^[25], 适合对区域尺度大, 数据繁杂的裂缝带进行预测。ELM的原理如下。

给定包含N个任意样本(x_i, t_i)的数据集, 输入层节点数为n, 输出层节点数为m, 其中x_i=[x_i1, x_i2, …, x_in]∈Rⁿ, t_i=[t_i1, t_i2, …, t_im]∈R^m。对于一个激励函数为g(x), 且有K个隐藏节点的单隐藏层的神经网络:

$ \sum\limits_{i = 1}^K {{\mathit{\boldsymbol{\beta }}_i}} g\left( {\left\langle {{\mathit{\boldsymbol{w}}_i},{\mathit{\boldsymbol{x}}_j}} \right\rangle + {b_i}} \right) = {\mathit{\boldsymbol{o}}_j}\quad j = 1,2, \cdots ,N $

(1)

式中:g(x)可选用Sigmoid函数、Gaussian函数等; w_i=[w_{i, 1}, w_{i, 2}, …, w_{i, n}]^T为第i个隐藏节点与输入节点间的权值向量; β_i=[β_{i, 1}, β_{i, 2}, …, β_{i, m}]^T为第i个隐藏节点与输出节点间的权值向量; b_i是第i个隐藏节点的偏置; 〈w_i, x_j〉表示w_i和x_j的内积; o_j为输出值。

ELM的网络结构如图 1所示。

已知单隐层神经网络的学习目标是使输出的误差最小, 即存在β_i, w_i和b_i使得:

$ \sum\limits_{i = 1}^K {{\mathit{\boldsymbol{\beta }}_i}} g\left( {\left\langle {{\mathit{\boldsymbol{w}}_i},{\mathit{\boldsymbol{x}}_j}} \right\rangle + {b_i}} \right) = {\mathit{\boldsymbol{t}}_j}\quad j = 1,2, \cdots ,N $

(2)

公式(2)用矩阵表示为:

$ \mathit{\boldsymbol{H\beta }} = \mathit{\boldsymbol{T}} $

(3)

式中:H表示为神经网络隐层的输出矩阵; β为输出权重矩阵; T为期望输出矩阵。

$ \begin{array}{l} \mathit{\boldsymbol{H}}\left( {{w_1},{w_2}, \cdots ,} \right.\left. {{w_K},{b_1},{b_2}, \cdots ,{b_K},{x_1},{x_2}, \cdots ,{x_N}} \right)\\ \;\;\;\;\;\;\;\;\;\;\;\; = {\left[ {\begin{array}{*{20}{c}} {g\left( {{w_1} \cdot {x_1} + {b_1}} \right)}&{g\left( {{w_2} \cdot {x_1} + {b_2}} \right)}& \cdots &{g\left( {{w_K} \cdot {x_1} + {b_K}} \right)}\\ {g\left( {{w_1} \cdot {x_2} + {b_1}} \right)}&{g\left( {{w_2} \cdot {x_2} + {b_2}} \right)}& \cdots &{g\left( {{w_K} \cdot {x_2} + {b_K}} \right)}\\ \vdots & \vdots & \vdots & \vdots \\ {g\left( {{w_1} \cdot {x_N} + {b_1}} \right)}&{g\left( {{w_1} \cdot {x_N} + {b_2}} \right)}& \cdots &{g\left( {{w_K} \cdot {x_N} + {b_K}} \right)} \end{array}} \right]_{N \times K}} \end{array} $

(4)

通常, 期望找到${\mathit{\boldsymbol{\widehat w}}_i}, {{\hat b}_i}$和${\mathit{\boldsymbol{\widehat \beta }}_i}$, 从而实现训练单隐层神经网络的目的, 使得:

$ \begin{array}{*{20}{c}} {\left\| {\mathit{\boldsymbol{H}}\left( {{{\mathit{\boldsymbol{\widehat w}}}_i},{{\hat b}_i}} \right){{\mathit{\boldsymbol{\widehat \beta }}}_i} - \mathit{\boldsymbol{T}}} \right\| = \mathop {\min }\limits_{\mathit{\boldsymbol{w}},b,\mathit{\boldsymbol{\beta }}} \left\| {\mathit{\boldsymbol{H}}\left( {{\mathit{\boldsymbol{w}}_i},{b_i}} \right){\mathit{\boldsymbol{\beta }}_i} - \mathit{\boldsymbol{T}}} \right\|}\\ {i = 1,2, \cdots ,K} \end{array} $

(5)

这等价于求解最小化损失函数:

$ E = \sum\limits_{j = 1}^N {{{\left[ {\sum\limits_{i = 1}^K {{\mathit{\boldsymbol{\beta }}_i}} g\left( {{\mathit{\boldsymbol{w}}_i} \cdot {\mathit{\boldsymbol{x}}_j} + {b_i}} \right) - {\mathit{\boldsymbol{t}}_j}} \right]}^2}} $

(6)

传统的基于梯度下降算法的神经网络可以解决这类问题, 但需要在迭代过程中不断调整参数^[26]。根据ELM理论, 输入权重值w_i和隐藏层偏置b_i被随机确定后, 隐层的输出矩阵H就不会再变化, 恒为常数矩阵。这时, ELM的训练过程等效为求解Hβ=T的最小二乘解${\mathit{\boldsymbol{\widehat \beta }}_i}$。如果隐层节点数K等于训练样本数N, 则矩阵H是方阵而且可逆, 当输入权值和隐藏层偏置随机赋值时, ELM可以以零误差逼近训练样本。则Hβ=T的最小范数二乘解为:

$ \mathit{\boldsymbol{\widehat \beta }} = {\mathit{\boldsymbol{H}}^{†} }\mathit{\boldsymbol{T}} $

(7)

式中:H^†为隐层输出矩阵H的Moore-Penrose广义逆矩阵。

1.2 几种常见的裂缝带预测属性

当利用机器学习算法解决分类问题时, 选取的输入特征需要相对独立。本文选用曲率、相干和反射强度属性作为输入特征, 它们具有不同的物理意义。

曲率属性属于地震几何属性, 能够描述地质体的几何变化, 由于其对断裂和裂缝等反应敏感, 可以用于预测裂缝的发育程度。曲率属性包括面曲率和体曲率属性, 面曲率属性仅利用了层位数据, 而体曲率属性充分利用了地震体的倾角、振幅等信息。考虑到实际三维工区的需要, 我们选取体曲率属性作为本文的输入特征之一。

相干属性由地震数据中相邻地震道之间的相似性来确定, 能识别并表征断层及裂缝发育带的分布规律和延展形态。相干体技术有3类, 一是基于互相关的相干算法, 二是基于相似性分析的相干算法, 三是基于特征值分析的相干算法。相较于第3类算法, 前两类算法的抗噪性相对较差, 不利于压制噪声, 计算时窗选取不当易产生人为假像。本文采用基于特征值分析的相干体技术, 具有较高分辨率。

反射强度属性即瞬时振幅是振幅类属性, 为时间域内所选样点各道在不同时刻对应的振幅大小。反射强度属性可以作为断层解释的依据, 有助于对大断层的检测。

我们利用ELM进行裂缝带预测时, 首先根据测井解释报告等资料将工区大致划分成裂缝欠发育区、裂缝较发育区和裂缝发育区3部分, 再按区域将3种属性数据按划分区域分别制作成训练样本, 处于裂缝欠发育区域样本视为1类样本, 裂缝较发育区域样本视为2类样本, 裂缝发育区域样本视为3类样本, 以此构建一个三分类问题, 利用ELM学习, 得出网络模型, 判别由目的层全部区域内曲率、相干和反射强度属性制成的预测集, 最终实现裂缝带预测。

由于测井或地震数据输入单位不一, 有些数据变化范围较大, 有些数据变化范围较小, 导致神经网络训练时间长, 收敛慢, 还有可能使得输入属性的作用权重不同, 影响训练结果。因此, 要将网络训练的目标数据映射到激活函数值域。本文做法是将属性值归一化到[-1, 1]。归一化公式为:

$ Y = 2 \times \frac{{X - {x_{\min }}}}{{{x_{\max }} - {x_{\min }}}} - 1 $

(8)

式中:Y为归一化后属性值; X为归一化前属性值; x_min为该类属性最小值; x_max为该类属性最大值。

我们建立了基于ELM的裂缝预测技术流程(图 2)。

2 模型数据实验 2.1 模型设计

本文选用川东北某工区2组测井数据(w1井和w2井)以及对应井旁道地震属性数据来验证ELM算法的分类效果。采用图 2所示的技术流程。在w1井和w2井油气储层裂缝发育、欠发育区域分别选取自然伽马(GR)、声波时差(AC)、补偿中子(CNL)、浅侧向电阻率(RS)、深侧向电阻率(RD)5种测井数据制作训练数据集和预测数据集(表 1)。这5种测井参数的数值在裂缝发育区和裂缝欠发育区有较为明显的差异, 可以较好地区分裂缝发育区和欠发育区。w1井数据制作成为训练数据集, w2井数据则为预测数据集。另外, 与测井数据集制作方法相同, 在w1井和w2井旁地震道按油气储层裂缝发育、欠发育区域选取常用于揭示裂缝特征的曲率、相干和反射强度等地震属性来制作训练数据集和预测数据集(表 2)。w1井数据制作成为训练数据集, w2井数据则为预测数据集。

基于测井解释结果及测井曲线, 挑选w1井裂缝区域测井数据, 并将裂缝发育状况按裂缝欠发育和裂缝发育分别附上1类和2类标签作为训练数据集。再依据w2井测井解释结果及数据制作预测数据集并附上标签。同理, 挑选w1井旁单地震道地震属性数据, 按裂缝发育状况分别用1类和2类标签标注作为训练数据集。挑选w2井旁单地震道地震属性数据制作预测数据集并附上标签。将训练数据集与预测数据集作为输入, 运用ELM算法进行分类计算, 得到分类结果。

2.2 测试及结果 2.2.1 基于测井数据

基于实验模型的分类效果主要是由分类正确率和计算用时来衡量。分类正确率即此类样本数据划分正确的数量与此类样本总数的比值。

制作的训练数据集由509个裂缝欠发育带样本数据和491个裂缝发育带样本数据组成; 预测数据集包含490个裂缝欠发育带样本数据和501个裂缝发育带样本数据。将训练数据集和预测数据集作为输入使用ELM算法进行分类。预测的结果如表 3。

预测结果表明ELM算法耗时与隐层节点设置数量有直接联系, 隐层节点数设置越多, 算法耗时越长。总体上, 随着隐层节点数量增加, 预测正确率越高。但当隐层节点数大于训练数据集数量后, 正确率增加不明显, 且有波动。实验表明, 隐层节点数设置1000时, 预测数据集正确率能达到94.95%, 已经满足了分类预测的要求。

为了研究ELM算法与PSVM算法的分类效果, 我们利用PSVM算法对相同的训练数据集和预测数据集进行分类。理论上, 分类效果直接受到映射方式的影响, ELM因其随机输入权重和隐层偏置的特性使其有若干种方式投影到高维, 且训练速度极快; 而PSVM映射方式受核函数影响, 加之训练速度相对较慢, 总体性能相对不及ELM。分类效果如表 4。可以看出ELM算法正确率略高于PSVM算法并且耗时较短。随机输入权重和隐层偏置使得ELM算法在处理大数据训练集时有更快的运算效率。

由于揭示裂缝欠发育带和裂缝发育带的测井数据局部有重叠, 因此采用多属性融合的方式划分复杂数据类型能有效提高分类正确率。图 3为PSVM算法与ELM预测结果的声波时差和自然伽马的交会图, 图中橙色圆点和红色叉号对应的声波时差和自然伽马数据均来自测井裂缝发育带, 浅蓝圆点和深蓝叉号数据则来自裂缝欠发育区域, 叉号是预测失败部分。从黑色虚线框所在的声波时差和自然伽马数值交集部分分类效果来看, ELM算法分类效果明显优于PSVM算法, 表明PSVM算法对数据重叠部分的分类效果不如ELM算法。ELM算法在划分重叠部分的数据时能够保持较高的分类正确率。

2.2.2 基于地震属性数据

利用w1井旁道地震属性数据制作的训练数据集由479个裂缝欠发育带样本数据和521个裂缝发育带样本数据组成。利用w2井旁道地震数据制作的预测数据集包含218个裂缝欠发育带样本数据和237个裂缝发育带样本数据。将训练数据集和预测数据集作为输入分别使用ELM算法和PSVM算法进行分类。ELM算法隐层节点数设置为1000。二种算法的分类效果见表 5。

预测结果表明, ELM算法的正确率略高于PSVM算法并且有计算耗时更短的优势。基于地震数据和测井数据的ELM的裂缝分类都能保持较高的正确率。选择适当的样本集, 地震数据3属性的分类效果就能够同测井数据5属性的分类效果相当, 能够满足实际应用的要求。

3 实际应用

测井数据的ELM分类常常能够达到较高的精度, 对岩样中的一条或几条裂缝有一定的识别能力。在尺度大、数据多的地震勘探中, 要做到识别单条裂缝, 常规技术达不到这种精度。但对于发育有较多数量裂缝且具有一定规模的裂缝带区域, 地震数据可以满足裂缝预测的识别精度。针对以上研究内容, 选取川东北某工区须家河组须四段地震资料验证ELM算法用于裂缝带的预测效果。

3.1 研究区概况

川东北地区区域内断层主要发育在九龙山构造转折端, 须家河组须四段主要发育北西和北东两组方向的裂缝。成像测井等资料表明北西向中高倾角裂缝开启性好。研究区内致密砂岩气藏发育, 有4口产气井(w1, w2, w3和w4井), 其中w1井、w2井、w3井为高产气井。受北西向对冲、背冲断层组合影响, 有效裂缝发育规模大, 高产单井产量均达到12×10⁴m³/d。单井日产量及测试段裂缝的地震相类型见表 6。

裂缝的地震相是不同尺度的断裂在地震剖面上所反映的特征。在该研究区域地震相分析, 解释为5类地震相:第1类为大断裂; 第2类为断裂; 第3类为微断裂; 第4类为裂缝; 第5类为基质。第1、2、3、4类裂缝的地震相为研究区域须四段裂缝发育有利相带, 占总数据比例为65%。

3.2 裂缝带预测训练数据集的选取

体曲率属性可以有效表征裂缝发育程度, 相干属性能反映断层空间展布特征, 反射强度属性反映了地震反射波能量强度的变化。尽管这些地震属性能够从不同角度刻画裂缝带的各类特征, 但仅使用单一属性不能全面客观地评价裂缝发育的实际情况, 容易造成多解等问题。因此, 综合利用多种属性, 将体曲率、相干属性和反射强度属性融合, 利用ELM算法进行裂缝发育带的预测。

训练数据集由研究区内4口井的井旁道体曲率、相干属性和反射强度属性数据构成。再基于完钻报告等测井解释资料将裂缝的发育程度分成3类:裂缝欠发育, 裂缝较发育和裂缝发育, 分别对应附上标签1、2、3。

在制作训练数据集时需要考虑以下两点:

1) 训练数据的选取要有代表性。本文是在测井解释结果及测井曲线的基础上筛选数据, 考虑到数据的均衡性, 选取的数据应兼顾主裂缝带和次级裂缝带。分析裂缝发育带、较发育带和欠发育带的数据特征, 根据其差异挑选典型的数据样本, 尽量使得局部重叠的数据最小化。

2) 训练数据中3类样本的比例要合适。在针对整个研究区进行大范围的裂缝带预测时, 需要考虑训练数据中各类样本之间数量的比例。根据地质背景等资料, 研究区域裂缝欠发育带范围远大于裂缝发育带和较发育带且裂缝发育带相比于较发育带规模更小。因此, 选取的训练数据中3类样本的占比需要根据实际情况调整。

3.3 预测效果

不同地震属性对裂缝带识别各有优势, 但同时也存在一定缺陷。图 4为研究区须四段3种地震属性的等时切片。可见, 体曲率属性(图 4a)虽能够反映裂缝发育情况, 但易受噪声干扰, 对地层起伏形态的刻画也受制于人工解释的主观性影响, 所以位于欠发育带的w4井的裂缝发育情况难以分辨。反射强度属性(图 4b)表现出横向变化趋势, 却在一定程度上损失了垂直分辨率, 在大尺度裂缝预测应用中, 存在一定误差, 导致研究区高产气w3井是否在裂缝发育带上难以准确判别。相干属性(图 4c)对数据信噪比有一定要求, 在低信噪比情形下预测会产生假相干, 因此, 预测结果中w4井所在区域的裂缝发育状况不太准确。

利用ELM算法基于研究区地震三属性数据进行裂缝预测。图 4d为ELM算法对研究区裂缝发育带的预测结果。图中白色部分代表裂缝带欠发育区域, 绿色部分代表裂缝带较发育区域, 红色部分代表裂缝带发育区域。预测结果较好地反映研究区内大断裂的基本形态, 准确地将w1井、w2井、w3井及其附近裂缝发育情况预测出来, 这3口高产气井均位于裂缝较发育带上, 且井位附近也有对裂缝发育带分布形态的表征。同时, 将单地震属性难以判断的w4井所在位置裂缝发育情况准确地进行了预测。位于裂缝欠发育带的w4井预测结果符合其低产井的实际开采情况。以上预测结果与测井资料、单井产量等钻探信息一致, 说明ELM算法针对该研究区裂缝带类型预测具有较高的可信度。

4 结论

1) 将ELM算法用于裂缝带综合预测, 相较于单地震属性裂缝带预测, 综合预测精度高, 同时能有效避免单一属性在裂缝时的多解性。研究区裂缝带预测结果表明, ELM算法能有效应用于裂缝带预测, 这也为ELM算法推广到其它地区提供了思路。

2) 因为ELM算法能随机生成输入与输出层间连接权重以及隐层神经元阈值, 极大简化了样本的学习过程, 训练效率显著提高。同时相比PSVM算法, ELM算法预测结果准确率高、用时短。针对数据量大、范围广的工区进行裂缝带分类预测具有较大潜力。

在一定范围内, ELM算法隐层节点个数越多分类准确率越高。当隐层节点个数设置过多时(一般远大于样本数量), 分类效果受节点增多影响就会变小, 有时分类正确率反而会有所降低。此外, 计算时间随节点个数增多而增加。因此, 在实际应用中, 应根据工区需求合理选择节点数量。