页岩气作为典型的非常规油气, 已成为我国油气勘探开发的重要领域。而页岩总有机碳(TOC)含量是评价该类气藏的关键指标之一, 该指数的高低决定了产气能力的大小^[1-2]。岩心总有机碳含量分析数据离散且数量有限, 而利用地震资料直接预测页岩储层总有机碳含量的研究文献也比较少, 因此为了有效预测页岩储层的总有机碳含量, 需要测井和地震资料相结合, 提高总有机碳含量预测精度。

目前, 利用测井资料预测总有机碳含量曲线的研究有: SCHMOKER^[3]利用体积密度法; AUTRI等^[4]利用声波时差、电阻率和放射性资料; PASSEY等^[5]利用ΔlogR法; ALTOWAIRQI等^[6]利用总有机碳含量与弹性参数的线性关系, 采用多元回归分析方法; 以及孟召平等^[7]利用BP神经网络算法计算总有机碳含量曲线。尽管这些方法对总有机碳含量曲线的预测有一定实际意义, 但都只是根据岩心总有机碳含量数据与测井曲线的相关性, 在二维平面中开展页岩总有机碳含量测井曲线预测, 无法定量预测总有机碳含量的三维空间分布。利用地震资料获取总有机碳含量的研究目前有: 金吉能等^[8]、尹俊等^[9]利用地震多属性反演技术预测总有机碳含量平面分布, 虽然总有机碳含量预测值与实测值较为吻合, 但该方法存在较强的多解性; 许杰等^[10]基于地质统计学反演的地震储层预测技术对页岩储层总有机碳含量进行预测, 但该技术没有考虑页岩储层的特殊性, 只适用于常规气藏储层预测; 陈祖庆^[11]、邓宇等^[12]和王健等^[13]通过总有机碳含量与密度交会分析, 建立总有机碳含量与密度的线性拟合方程, 得到计算总有机碳含量的经验公式, 再结合叠前地震反演获得的密度数据体计算出总有机碳含量数据体, 但该方法只考虑总有机碳含量与密度的关系, 存在相关性不够、预测精度低的问题, 无法满足页岩储层精细预测的需求^[14]。

针对上述研究存在的不足, 以四川盆地渝西地区Z井区为例, 基于总有机碳含量与弹性参数存在多因素、多变量的相关性, 提出了一种基于粒子群优化(particle swarm optimization, PSO)的支持向量机(support vector machine, SVM)计算页岩气储层总有机碳含量预测算法。虽然国内一些学者如张成龙等^[15]已经用支持向量机算法进行总有机碳含量预测, 但是他们的研究只是针对特定井段上的总有机碳含量测井曲线预测, 并没有利用地震反演数据体获得总有机碳含量数据体, 进行总有机碳含量平面和三维空间分布的预测。因此, 本文尝试将粒子群优化支持向量机算法运用于三维空间中总有机碳含量数据体的求取, 该算法是借助于地质、测井和地震等资料建立总有机碳含量与敏感参数复杂的非线性关系, 在非均质性较强的页岩气储层中, 利用该算法进行总有机碳含量的预测, 有效提高页岩气储层总有机碳含量的预测精度。

1 工区概况

渝西地区中部铜梁、大足、璧山一带是页岩勘探有利地区, 其中弥陀场向斜、蒲吕场向斜及西山背斜页岩气资源储量丰富, 显示渝西地区龙马溪组—五峰组具有丰富的资源。Z井区位于西山构造西侧、蒲吕场向斜中部, 区内地腹构造平缓, 断裂不发育, 是页岩气勘探有利区^[16-17]。

沉积相带是控制页岩气富集的关键因素, 据研究^[18], 将龙马溪组划分为浅水泥质陆棚亚相、半深水陆棚亚相和深水陆棚亚相3种沉积亚相类型。其中, 深水陆棚亚相页岩具有有机质含量高、岩石脆性高、有机质孔发育和孔隙保存条件好等特征, 是龙马溪组页岩气富集的最有利相带。渝西地区Z井区主要发育深水陆棚-浅水泥质陆棚亚相, 表明该区具有广阔的页岩气勘探开发前景。

由钻井和测试资料可知, 研究区钻达五峰组的井有5口, 即Z1井、Z2井、Z3井、Z5井和Z6井。主要目的层为龙马溪组龙一₁—五峰组, 根据总有机碳含量值将研究区页岩储层划分为3类: 总有机碳含量>3%为Ⅰ类优质页岩, 2% < 总有机碳含量 < 3%为Ⅱ类优质页岩, 总有机碳含量 < 2%为Ⅲ类非优质页岩。由Z井区龙马溪组+五峰组优质页岩储层参数表(表 1)可知, Z1井储层段Ⅰ类优质页岩厚度为9.4m, Z2井储层段Ⅰ类优质页岩厚度为8.3m, Z3井储层段Ⅰ类优质页岩厚度为13.5m, Z5井储层段Ⅰ类优质页岩厚度为12.2m, Z6井储层段Ⅰ类优质页岩厚度为10.9m。可以看出, Z井区主要目的层储层段Ⅰ类优质页岩薄, 储层预测难度大, 需要精细的储层预测支撑该地区页岩气整体勘探开发。

由Z6井测井曲线(图 1)可以看出, 优质页岩储层发育段与上、下围岩在纵波速度、横波速度、自然伽马、密度和纵横波速度比等方面有明显差异, 优质储层(龙一₁—五峰组段)具有明显的高伽马、低速度、低密度和低纵横波速度比等特征^[19-21]。

对研究区5口井龙马溪组实测岩心总有机碳含量数据和测井资料进行敏感参数分析(图 2), 发现总有机碳含量与密度、自然伽马和纵横波速度比等参数相关性较好, 表现为密度、纵横波速度比与总有机碳含量呈负相关性(图 2a和图 2c), 自然伽马与总有机碳含量呈正相关性(图 2b), 而总有机碳含量与杨氏模量的相关性较差(图 2d)。原因在于: ①泥页岩中有机质的密度为1.03~1.10g/cm³, 围岩基质中黏土骨架的密度为2.30~3.10g/cm³, 富含有机质的页岩密度小于泥页岩密度, 有机质含量的增加导致页岩密度降低, 因此, 总有机碳含量随着密度的增加而减小, 与密度曲线呈负相关; ②通常泥质含量越高, 自然伽马放射性也越高。页岩中泥质含量较高, 高放射性物质的含量也更多。总有机碳含量高的烃源岩具有粒度小、表面积大的性质, 可以吸附较多的放射性元素, 导致总有机碳含量随着自然伽马的增加而上升, 因此总有机碳含量与自然伽马曲线呈正相关; ③页岩比砂泥岩致密, 孔隙度小, 页岩的纵波速度介于泥岩和砂岩之间。声波在有机质中传播速度较慢, 有机质含量的增加导致页岩纵波速度降低, 因此总有机碳含量随着纵横波速度比的增加而减小, 与纵横波速度比曲线呈正相关。

综上所述, 选择自然伽马、密度和纵横波速度比作为敏感参数, 进行总有机碳含量预测^[22]。

常规方法中, 页岩总有机碳含量预测是通过建立总有机碳含量与密度的交会关系, 由叠前反演得到的密度参数数据体计算得到^[23-25]。研究区总有机碳含量与密度、自然伽马、纵横波速度比和杨氏模量的拟合公式如表 2所示, 总有机碳含量受多种地质因素影响, 与测井数据存在复杂的非线性关系, 仅用单一的密度通过线性交汇法求总有机碳含量时准确性较低。为了得到更准确的总有机碳含量, 本文提出利用支持向量机算法计算总有机碳含量。

2 方法原理 2.1 支持向量机原理

支持向量机(SVM)是VAPNIK等^[26]在1996年提出的以统计学习理论和结构风险最小化为基础的新型机器学习方法。基本思想是利用核函数, 将低维空间的自变量x通过非线性关系φ(·)映射到高维空间, 然后在高维空间使用线性判别函数解决低维空间的非线性问题, 通过建立决策函数输出预测值y。

对于总有机碳含量预测, 给定样本集{(x_i, y_i), i=1, 2, 3, …, m}, x_i是二维空间的属性向量, 为{自然伽马、密度、纵横波速度比}曲线的随机输入值; y_i是对应的预测目标值, m是样本数。

对于低维非线性空间的页岩总有机碳含量预测模型, 通过非线性函数φ(x)将样本集数据x映射到高维线性空间, 在高维线性空间中建立决策函数, 即:

$ \mathit{\Phi}(x)=\boldsymbol{\omega} \cdot \varphi(x)+b $

(1)

式中: ω是权值向量, b是偏置参数。

基于统计学习理论的结构风险最小化思想, 算法可以利用极小化优化模型来建立决策函数, 即

$ \min \frac{1}{2}\|\boldsymbol{\omega}\|^{2}+C \sum\limits_{i=1}^{m}\left(\zeta_{i}+\zeta_{i}^{*}\right) $

(2)

约束条件为: $\left\{\begin{array}{l}y_{i}-\boldsymbol{\omega} \cdot \varphi(x)-b \leqslant \varepsilon+\zeta_{i} \\ {[\boldsymbol{\omega} \cdot \varphi(x)]+b-y_{i} \leqslant \varepsilon+\zeta_{i}^{*}} \\ \zeta_{i}, \zeta_{i}^{*} \geqslant 0 \quad(i=1, 2, 3, \cdots, m)\end{array}\right.$式中: C为惩罚因子, 代表对超出误差ε的样本的惩罚程度; ζ_i, ζ_i^*为松弛因子; ε为误差上限。

引入拉格朗日乘子α_i和α_i^*, 将最小化优化模型转化为对偶优化问题进行求解, 即对Q(α)求极大值:

$ \begin{aligned} &\max : Q(\alpha)=-\frac{1}{2} \sum\limits_{i, j=1}^{m}\left(\alpha_{i}-\alpha_{i}^{*}\right)\left(\alpha j-\alpha_{j}^{*}\right) \cdot \\ &K\left(x_{i}, x_{j}\right)+\sum\limits_{i=1}^{m} y i\left(\alpha_{i}-\alpha_{i}^{*}\right)-\varepsilon \sum\limits_{i=1}^{m}\left(\alpha_{i}+\alpha_{i}^{*}\right) \end{aligned} $

(3)

约束条件为:

$ \left\{\begin{array}{l} \sum_{i=1}^{m} \alpha_{i}=\sum_{i=1}^{m} \alpha_{i}^{*} \\ 0 \leqslant \alpha_{i} \quad \alpha_{i}^{*} \leqslant C \quad(i=1,2,3, \cdots, m) \end{array}\right. $

式中: Q(α)是拉格朗日函数的对偶目标函数, K(x_i, x_j)=φ(x_i)·φ(x_j)为核函数, 当α_i×α_i^*=0且α_i, α_i^*非同时为0时, 所对应的x_i为支持向量(SV)。

本文选择径向基函数作为核函数:

$ K\left(x, x_{i}\right)=\exp \left(-\frac{\left\|x-x_{i}\right\|}{2 \sigma^{2}}\right) $

(4)

因此, 算法计算中涉及惩罚因子C和核函数参数σ两个主要参数^[27]。

2.2 粒子群算法原理

粒子群算法是通过模拟鸟群觅食行为发展起来的一种基于群体协作的寻优处理算法, 其中, 每一只“鸟”表示一个粒子, 鸟群寻找的“食物”就是所求的最优解^[28]。在此次总有机碳含量预测中, 采用粒子群算法对“惩罚因子C和核函数参数σ”这两种粒子进行初始化处理, 第i个粒子的坐标为x_i^t=(x)_i1^t, x_i2^t, …, x_iD^t), x_id^t∈[L_D, U_D], D个目标解的坐标为未优化前的目标解。第i个粒子运行速度为v_i^t=(v_i1^t, v_i2^t, …, v_iD^t), v_id^t∈[v_{min, D}, v_{max, D}], 第i个粒子最优化位置为p_i^t=(p_i1^t, p_i2^t, …, p_iD^t), 种群最优化位置为p_g^t=(p_g1^t, p_g2^t, …, p_gD^t)。则第s维度的粒子运动速度和粒子坐标的迭代结果可表示为:

$ \begin{gathered} v_{i D, s}^{t}=\beta v_{i D, s}^{t}+c_{1} r_{1}\left(p_{i D, s}^{t}-v_{i D, s}^{t}\right)+ \\ c_{2} r_{2}\left(p_{g D, s}^{t}-v_{g D, s}^{t}\right) \\ x_{i D, s}^{t+1}=x_{i D, s}^{t}+v_{i D, s}^{t+1} \end{gathered} $

(5)

式中: v为粒子速度, β为惯性权重, c₁、c₂为学习因子, r₁和r₂为[0, 1]的随机数。

β具有平衡全局和局部搜索的能力, β的值越大越具有较好的全局收敛能力, 反之则有较强的局部收敛能力。学习因子c₁和c₂分别控制着粒子寻找个体最优位置和全局最优位置的能力。

参数优选中, 将自然伽马、密度和纵横波速度比作为输入量回归预测总有机碳含量, 若预测的总有机碳含量与实测总有机碳含量之间的误差在误差范围内, 则结束参数选优程序, 此时粒子群的最优位置即为惩罚因子C和核函数参数σ的最优解; 若预测的总有机碳含量与实测总有机碳含量超出了设定的误差范围, 则继续执行程序, 直至两者的误差在设定的误差范围内。通过粒子群算法寻找最优支持向量机参数, 利用寻找的最优参数构建总有机碳含量预测的支持向量机模型。因此粒子群算法作为优化领域中的有效算法, 在支持向量机参数寻优过程中可以取得良好的效果^[29]。

2.3 粒子群优化支持向量机算法实现步骤

本文提出的粒子群优化支持向量机算法实现步骤如下。

1) 将模型的适应度函数定义为训练样本的预测误差, 即最小预测误差所对应的粒子位置便是函数所对应的最优解;

2) 随机初始化粒子群位置和粒子速度;

3) 计算每个粒子的适应值;

4) 将当前适应值和个体历史最佳位置适应值进行对照。若前者值较高, 则用当前位置代替个体历史最佳位置;

5) 将当前适应值和全局最佳位置适应值进行对照。若前者值较高, 则用当前位置代替全局最佳位置;

6) 利用粒子群算法进行粒子速度和位置的更新;

7) 如未满足结束条件, 则返回步骤3)。反之则执行下一步。最后得到最优参数C和σ值;

8) 将最优参数C和σ值进行SVM算法训练, 再运用验证集验证算法的预测精准度, 得到储层总有机碳含量的预测结果。

相应的算法流程如图 3所示。

3 实际应用及效果分析 3.1 数据实验与分析

由图 2分析可知, 对总有机碳含量反应较为敏感的测井曲线主要为自然伽马、密度和纵横波速度比测井曲线。因此, 选取上述3条曲线数据作为总有机碳含量识别的原始样本输入数据, 输出结果为粒子群-支持向量机算法计算的总有机碳含量曲线。从已知5口井中, 选择Z2、Z3、Z5和Z6井作为训练井, 参与粒子群优化支持向量机算法计算; Z1井作为验证井不参与计算, 用来检验算法的精准性。随机在训练井选取400个龙一₁—五峰组段样本数据作为粒子群优化支持向量机算法的训练集, 用于预测总有机碳含量曲线; 验证井选取100个样本作为验证集, 训练集和验证集数据是分别来自5口井的自然伽马、密度和纵横波速度比测井解释数据。

粒子群算法的参数设置依据为: 核函数采用径向基函数; 种群规模影响着算法的搜索能力和计算量, 依据模型大小和算法性能, 通过多次经验尝试, 把粒子群算法初始规模设置为20可以达到很好的求解效果; 惩罚因子和核函数参数的范围影响算法的精度和预测率, 参数过大, 降低泛化能力, 增加训练时间, 反之则预测精度较低。为了保证算法的稳定性和高效性, 惩罚因子C的粒子群搜索区间设置为[0.1, 100.0], 核函数参数σ的粒子群搜索区间设置为[0.01, 1000.00], 粒子群的最大进化代数为200;合适的c₁, c₂可以加快收敛且不易陷入局部最优, 根据研究区运算条件, 设置学习因子c₁=1.5, c₂=1.7。图 4为粒子群算法对支持向量机进行参数优化的适应度曲线, 计算得到最优参数C=1024和σ=0.25。

将搜索到的最优参数C和σ代入支持向量机算法, 利用粒子群优化支持向量机算法对训练集进行训练, 总有机碳含量曲线的预测正确率达到了95%以上, 如图 5所示, 表明算法参数设置合理, 进行总有机碳含量曲线预测可以获得较高的精确度。

3.2 总有机碳含量曲线预测

利用粒子群优化支持向量机算法进行总有机碳含量曲线预测, 相关的预测结果如图 6所示, 并与表 2中总有机碳含量与密度的关系式拟合的总有机碳含量曲线进行对比。图 6中第6道红色实线为粒子群优化支持向量机算法预测的总有机碳含量曲线, 第7道黑色实线为总有机碳含量与密度关系式拟合的总有机碳含量曲线, 黑色圆点为岩心总有机碳含量值, 可以看到粒子群优化支持向量机算法具有极强的非线性逼近能力, 能真实反映总有机碳含量与测井参数的关系, 由粒子群优化支持向量机算法预测的总有机碳含量曲线与岩心总有机碳含量吻合度较高。表 3给出了粒子群优化支持向量机算法在5口井中总有机碳含量预测正确率统计, 计算正确率较高, 均超过了95%, 说明该算法计算结果精度较高, 可用于页岩气储层总有机碳含量的地震预测。

3.3 总有机碳含量地震预测

基于总有机碳含量敏感参数统计分析结果, 密度、自然伽马和纵横波速度比反演数据体是总有机碳含量数据体地震预测的基础。密度反演数据体和纵横波速度比反演数据体可以利用叠前同时反演直接获得, 自然伽马反演数据体可以通过测井参数反演得到。图 7a、图 7b和图 7c分别是过Z1和Z6井连井密度、纵横波速度比和自然伽马反演剖面, 井旁反演结果与测井曲线相比, 吻合度均较高。图 7a中龙马溪组龙一₁—五峰组底储层段密度较低, 主要在2.50~2.62g/cm³, 呈红黄色条带。图 7b中龙马溪组龙一₁—五峰组底储层段纵横波速度比较低, 主要分布于1.40~1.64, 呈红黄色条带。图 7c龙马溪组龙一₁—五峰组底储层段自然伽马值较高, 主要分布于150~190API, 呈红色条带。

在密度、自然伽马和纵横波速度比反演数据体的基础上, 利用粒子群优化支持向量机算法就可以计算出总有机碳含量反演数据体。图 8为Z1和Z6连井总有机碳含量反演剖面, 红色部分为优质页岩, 总有机碳含量相对较高, 龙马溪组龙一₁—五峰组底储层段优质页岩总有机碳含量较高, 主要分布于2.0%~3.5%, 剖面从上到下总有机碳含量逐渐增大, 五峰组储层特征明显, 横向上展布稳定。Z1和Z6井井旁总有机碳含量预测结果与井的分析结果较为一致, 表明反演结果真实可靠。图 9为利用密度与总有机碳含量拟合公式计算的Z1和Z6连井总有机碳含量反演剖面, 图中, 五峰组储层特征不明显, 表明利用粒子群优化支持向量机算法预测总有机碳含量的效果比利用密度法预测更好。

根据总有机碳含量预测数据体, 绘制出龙一₁段—五峰组底优质页岩总有机碳含量平面分布(图 10), 从图中可以看出, 工区内龙一₁段—五峰组底优质页岩总有机碳含量整体较高, 大多在2.0%~2.7%之间变化。Z1, Z2, Z3, Z5和Z6井的预测评价总有机碳含量值分别为2.27%, 2.38%, 2.49%, 2.57%和2.42%, 与实测结果基本一致, 相对误差较小。

页岩气开发的靶体主要是龙一₁各小层, 为了更好地展示各小层的总有机碳含量, 再次绘制龙一₁³段顶—五峰组底平面分布(图 11)。由图 11可见, 工区内龙一₁小层优质页岩总有机碳含量同样整体较高, 分布在2.0%~2.7%。

表 4对比了龙马溪组龙一₁段—五峰组底优质页岩总有机碳含量预测结果与实际测井解释结果, 可以看出, 本区龙马溪组优质页岩各项特征均比较稳定, 优质页岩总有机碳含量预测结果与测井解释结果误差较小, 绝对误差在0.03%以内, 相对误差在1.4%以内, 优质页岩总有机碳含量预测结果可靠。

4 结论

利用粒子群优化算法得到支持向量机算法的最优参数, 利用最优参数创建总有机碳含量预测模型, 将算法应用于四川盆地渝西地区Z井区优质页岩总有机碳含量预测, 具有较好的应用效果。

1) 粒子群优化支持向量机算法通过在高维空间进行总有机碳含量与敏感参数的线性拟合, 解决了低维空间中非线性关系, 为页岩气储层总有机碳含量的预测提供了更加准确的数学算法;

2) 粒子群优化支持向量机算法在测井曲线预测中, 预测的总有机碳含量曲线和岩心总有机碳含量吻合度较高, 达到了95%, 表明算法参数设置合理;

3) 通过粒子群优化支持向量机算法对储层总有机碳含量进行地震反演预测, 井旁反演结果、平面预测结果与井位处总有机碳含量曲线较为吻合; 表明算法具有较好的实用性, 能够提高储层预测精度, 对四川盆地渝西地区页岩气勘探开发具有指导意义。