林利明 郑 颖 史 浩 李忠城 时靖雪 师素珍 郭俊超
1.中联煤层气有限责任公司 2.中海油研究总院有限责任公司 3.“煤炭资源与安全开采”国家重点实验室·中国矿业大学（北京）
随着国内油气勘探程度不断加深，非常规油气勘探目标的储层条件日益复杂，尤其是薄储层含气性预测面临诸多挑战。传统的反演技术作为储层预测的重要手段，由于受地震资料频带的限制，难以满足薄储层预测的精度要求。鄂尔多斯盆地东缘致密气资源丰富，二叠系石盒子组河流相储层，砂岩厚度薄、纵向多期叠置、横向空间变化快，薄储层含气性预测困难，通过高精度反演进行薄储层预测成为致密气高效勘探开发的关键技术之一。
地质统计学反演是随着地质统计学理论发展起来的高精度反演技术，最早由Haas等[1]在1994年提出，该方法结合地质、地震和测井等多种先验信息，纵向分辨率高。早期出现的方法有克里金技术和随机模拟技术，后来发展了以序贯指示和序贯高斯为主导的反演算法，目前应用较为广泛的为Torres-Verdin等[2]提出的马尔科夫链—蒙特卡罗算法（以下简称MCMC法），MCMC法不仅忠实于地震数据，而且具有高效收敛、全局寻优的特点。2001年，Seok-Hoon等[3]将运用地质统计学建立的先验信息与利用地震资料得到的似然函数结合，采用马尔科夫扰动模拟对得到的后验概率分布进行采样，最终得到反演结果。2005年，Contreras等[4]在随机反演过程中采用MCMC法，使随机反演在叠前地震数据中得到较好的应用。2006年，Hansen等[5]为降低反演的复杂性，将确定性反演问题与地质统计学理论相结合，取得了较好的效果。2010年，Cordua等[6]将贝叶斯理论与地质统计学先验信息相融合，实现了AVO随机反演。之后，部分学者研究了在岩石物理约束下进行地质统计学反演，如2010年，Bosch等[7]将岩石物理与地质统计学相结合进行储层预测的方法进行了系统总结。2016年，Connolly等[8]将大量伪井与角度叠加地震进行匹配，基于岩石物理模型进行随机反演，获取了高质量储层特征。2018年，Figueiredo等[9]将线性化的岩石物理模型合并到多元多模态先验分布中，基于贝叶斯反演得到储层的弹性参数、物性参数及岩相。
国内地质统计学的发展时间比较晚，但随着国内学者在随机反演方面的不断研究，地质统计学反演也取得了极大的发展。2005年，印兴耀等[10]介绍了一种贝叶斯—序贯高斯模拟方法，能够整合垂向分辨率与储层厚度相当的地震数据和测井数据，取得了较好的效果。2009年，黄捍东等[11]采用地震非线性随机反演方法，以测井、地质信息作为约束条件，对薄砂岩储层进行了预测。2012年，刘占族等[12]在常规波阻抗反演基础上，应用地质统计学反演提高了反演的纵向分辨率，精细描述了煤层气薄储层的空间分布形态。随着计算机的发展，地质统计学反演从前期的以叠后反演为主逐步发展到了叠前反演。2010年，孙月成等[13]研究了叠前AVA同步反演与随机模拟理论有机结合的反演方法，并将其应用于薄层预测。2020年，李炳颖等[14]针对中深层薄层精细预测难题，通过分析岩相概率和变差函数求取可靠的地质统计参数，应用叠前地质统计学反演技术降低了反演的不确定性。与常规反演相比，地质统计学反演结果在测井和地震数据等先验信息不足时没有明显优势，近年来不少学者在相控反演方面做了相关研究。2013年，张恒山等[15]结合沉积微相认识，开展相控地质统计学反演，并基于MCMC法提供了一种相对较准确的高精度岩性、物性预测结果。2021年，段新意等[16]提出了适用于先验信息不完备地区的相控地质统计学反演方法，取得的高精度反演结果对薄层有较好的反映。徐伟等[17]将含砂率引入到地质统计学反演中，有效地解决了三角洲相薄储层预测的困难。
笔者在地质统计学反演的基础上，结合叠前确定性反演以及地质统计学反演的优势，首先在储层岩石物理特征分析的基础上，优选出能反映储层岩性、含气性的敏感弹性参数，将叠前反演弹性参数体转换为岩相概率密度，从而作为地质统计学反演的空间约束项，提出了一种针对薄储层预测的相控地质统计学叠前反演技术，应用于鄂尔多斯盆地东缘临兴中区的致密砂岩含气性预测，取得了较好的预测效果。
1.1 基于MCMC法的地质统计学反演基本原理
基于MCMC法的地质统计学反演技术融合了随机数值模拟理论和地震反演技术[18]，打破了原来地震数据驱动反演的纵向分辨率极限1/4波长的限制，成为目前最为有效的薄储层预测技术之一[19]。图1为基于MCMC法的地质统计学反演工作流程。
图1 基于MCMC法的地质统计学反演工作流程图
反演问题从概率统计的角度可以看成是一种贝叶斯估计问题，即通过加入先验信息来进行约束，得到后验概率分布的一些性质。MCMC法是在贝叶斯框架下，基于Metropolis-Hastings算法，建立一个平稳分布的马尔科夫链，通过不断迭代使马尔科夫链在高概率密度区收敛，即贝叶斯方法中的最大后验估计，从而获得目标分布的样本，然后基于样本做出各种统计推断[21-22]。
其核心是马尔科夫链的构建，无后效性是马尔科夫链最重要的性质，即t＋1时刻的状态只依赖于t时刻的状态，与t时刻之前的状态无关[23]。马尔科夫链的构建过程[24]：
第一步，定义地震正演过程，即
式中d表示地震记录；
f表示正演算子；
m表示待估参数，无量纲；
e表示噪声，无量纲。
第二步，利用贝叶斯框架表征未知参数m的后验概率密度，即
式中p(m)表示先验概率密度；
p(d|m)表示似然函数。
第三步，假设噪声e是平均值为0、方差为σe的高斯随机噪声，可得到观测地震记录的似然函数，即
式中σe表示方差。
为了得到所估参数的后验概率分布p(m|d)，利用Metropolis-Hastings算法生成马尔科夫链。其建议分布函数取均匀分布，即
式中mt表示马尔可夫链在t时刻的状态值，无量纲；
m*表示t＋1时刻的潜在转移点，无量纲；
δ表示步长，无量纲。
式中π(m)表示马尔科夫链的平稳分布。
为了使马尔科夫链最终收敛于未知参数m的后验概率分布，平稳分布应为。
1.2 相控地质统计学叠前反演思路
常规地质统计学反演对地震数据和测井数据的依赖性较大，在钻井较少的情况下，误差往往较大，当不同岩性地层的弹性参数有严重叠置时，常规叠后地质统计学反演不能有效地反映横向岩性展布规律。因此，如何综合应用地震、地质、测井等多个维度的数据信息，获得兼具预测性和高分辨率的储层参数，是解决致密砂岩薄储层条件下储层预测的关键。地质统计学反演结果的好坏依赖于地质概念模型和随机函数模型的选择[20]，而叠前相控地质统计学反演是在地质统计学反演的基础之上，结合地质先验认识，提取符合沉积储层变化规律的岩性、含气性敏感弹性参数体，把敏感参数体转化为砂岩储层概率密度，并作为储层相约束地质统计学反演。其流程如图2所示，主要包括以下关键步骤：①岩石物理是连接地震和测井的桥梁[25]，开展岩石物理分析，优选出能反映储层岩性、含气性的敏感弹性参数。②建立岩相概率密度体。结合已钻井情况、地质先验信息，利用贝叶斯判别原理，把地震数据驱动的敏感弹性参数反演结果转换为符合地质认识的岩相概率密度，作为储层空间约束项，添加到三维地质统计学叠前反演过程中。③统计学参数分析。变差函数是用来描述储层属性空间变化的一种方法，能否求得理想的变差函数是地质统计学建模及反演的关键[14]。由于变程与储层的非均质性密切相关，是地质统计学反演过程中重要的参数。因此，变差函数拟合图分析的关键问题是得出三维方向上的变程[26]。④协模拟物性参数，通过地质统计学反演得到叠前弹性参数体后，要进一步得到各小层能直观反映气藏特征的物性参数（如孔隙度、渗透率、饱和度等），可以通过协模拟技术来实现，其可以模拟得到两个甚至多个变量之间的空间分布，且模拟结果符合变量之间的相关性信息。
图2 相控地质统计学叠前反演技术流程图
1.3 协模拟技术
地质统计学协模拟是以储层参数的统计规律为基础，通过随机模拟的方式，建立储层参数与弹性参数关系，从而实现储层物性参数预测的同步处理技术。与传统的高斯模拟相比，直接序贯协模拟不需要对原始变量数据进行任何非线性变换，因此可以较好地再现弹性参数和物性参数之间的关系，但仅可以再现两个呈线性变量之间的关系，当变量之间关系较为复杂时，此方法就不再适用[27]。基于此，Horta等[28]在其基础上进行了改进，提出了考虑联合分布的直接序贯协模拟方法，该算法基于贝叶斯法则，首先运用直接序贯模拟方法模拟得到所有次级变量，再将其作为条件数据，协模拟得到主变量，可以有效地表征非线性相关关系变量之间的统计分布特点。
其运用同位简单协克里金估计[29]计算局部条件均值和方差，对条件分布进行抽样得到模拟值。假设Z1表示已知变量，Z2表示待模拟变量，在估计Z2的局部均值和方差时，将Z2的已知原始数据和模拟Z1得到的待估点x0的Z1变量同时作为条件数据。Z2在待估点x0处的协克里金均值Z2(x0)*和协克里金方差σ2(x0)可表示为：
式中x0、xα分别表示待估点和临近点；
Z2(x0)*、Z2(x0)分别表示待模拟变量Z2在x0处的协克里金均值和变量值 ；
Z2(xα)表示Z2在xα处的变量值 ；
λxα表示Z2在xα处的权重系数；
λx0表示已知变量Z1在x0处的权重系数；
Z1(x0)表示Z1在x0处的变量值；
n1、n2分别表示变量Z1、Z2的均值；
σ2(x0)表示Z1在x0处的协克里金方差。
权重值可用式（9）进行计算，即
式中CZ2(xα,xβ)表示待模拟变量Z2在临近点xα处与xβ处之间的协方差；
CZ1Z2(x0,xα)表示待模拟变量Z2与已知变量Z1在临近点xα处与待估点x0处之间的协方差；
CZ2(x0,xβ)表示Z2在x0处与xβ处之间的协方差；
CZ1Z2(xβ,xα)表示Z2与Z1在xα处与xβ处之间的协方差；
CZ1(x0)表示Z1在x0处的方差；
CZ1Z2(x0)表示Z1与Z2在x0处之间的协方差。
采用联合分布直接序贯模拟方法可以利用已钻井的地质统计学特征进行随机扰动，从而更好地呈现非线性相关变量之间映射关系。另外，其抽样时采用的是变量之间的全局联合概率分布函数，可以更大程度地提高模拟的精度[30]。
2.1 研究区地质概况
鄂尔多斯盆地致密砂岩气藏是典型的岩性气藏，是中国非常规天然气的主力产区[31-33]。临兴中区位于鄂尔多斯盆地东缘，三维地震满覆盖面积746.29 km2，已提交致密砂岩气探明储量超过千亿立方米。临兴中区主要目的层是二叠系石盒子组，主要为三角洲平原沉积，其内部河道砂体呈树叉状分布，致密砂岩油气资源分布广泛[34-35]。储层岩性以岩屑石英砂岩为主，单砂体平均厚度小于15 m，含气砂体平均厚度小于5 m（图3），有效储层孔隙度介于3.0%～14.0%，渗透率介于0.1～5.0 mD，属于典型的致密砂岩储层[36]。储层具有明显的“低孔隙度、低渗透率、单砂体厚度薄、砂体纵向叠置、横向变化快”的特点，常规储层预测技术受限于地震分辨率，难以刻画小于5 m的含气薄层砂岩，所以该区致密砂岩薄储层的含气性预测十分困难。
图3 石盒子组砂体厚度分布直方图
2.2 储层岩石物理特性
对临兴中区的录井岩屑资料和测井曲线进行整理，将石盒子组划分为3个主要岩相类别：泥岩、致密干砂岩、含气砂岩。经过多井岩石物理统计分析，砂岩波阻抗与泥岩波阻抗存在部分叠置，单一参数难以区分，需要开展多参数进行砂岩、泥岩有效区分。从众多岩石物理交会组合当中，优选出能通过叠前反演直接获取、又具备明显物理意义的纵波波阻抗和伪泊松比（纵横波速度比，即vp/vs），并绘制石盒子组的纵波波阻抗和伪泊松比统计分布直方图（图4）和交会图（图5）。图5中蓝色散点为测井解释泥岩，绿色散点为测井解释的致密干砂岩，红色散点为测井解释含气砂岩。从交会图来看，可以初步认为伪泊松比具有良好的砂泥岩性区分度，含气砂岩也具备伪泊松比小于1.8的明显低异常特征。因此，可以将伪泊松比作为岩性、含气性敏感因子[37]，进行储层岩性概率密度计算。
图4 石盒子组弹性参数统计直方图
图5 石盒子组纵波波阻抗和伪泊松比交会图
2.3 建立储层岩相概率密度体
概率密度函数描述的是特定岩性对应的岩石物理参数分布特征[38]。将弹性参数反演结果用于岩相概率密度体提取的前提条件是反演得到的弹性参数必须能够反映储层岩性，研究区岩性符合率超过60%，表明砂岩与泥岩具备一定的区分度。然后依据地质沉积背景控制已钻井的权重，最终得到一个既能反映地震空间特征，又符合地质认识的三维岩相概率密度体。
对井点数据进行岩性统计，建立不同岩性对应的岩石物理参数响应范围及定量的概率分布函数，再通过云变换（描述两个参数之间的相关关系）将地震数据驱动的叠前同步反演得到的弹性参数体进行三维数据交会，借助贝叶斯判别原理即可将弹性参数体转换成地震分辨率下的岩相概率密度体。图6为临兴中区石盒子组的盒二段、盒四段、盒六段砂岩概率密度切片，其中虚线为河道沉积。从图6可以看出，暖色的红色与黄色为砂岩发育概率高值区，冷色的蓝色和绿色为砂岩发育概率低值区，研究区物源方向为北偏西，盒二段、盒四段、盒六段砂岩储层呈条带状分布，砂地比特征具有较明显差异，其中盒二段切片显示砂岩概率密度偏高，西北部富砂，储层发育；
盒四段中部砂岩概率密度偏高，中部富砂，砂岩储层发育；
盒六段砂岩整体概率密度偏低，呈现窄条带状分布。岩相概率密度与前期地质沉积认识较为匹配。
图6 石盒子组主力储层岩相概率密度切片图
2.4 测井统计参数分析确定变程、变差
变差函数描述的是横向和纵向地质特征的结构和特征尺度[38]。其中，横向变差函数大小主要影响储层的横向展布规律，纵向变差函数大小主要影响储层的厚度[16]。
通过对石盒子组的已钻井资料和地质信息统计分析，获得各个层段岩性的概率密度函数和变差函数。图7为临兴中区盒四段的测井统计变差函数，由于砂岩薄层、薄互层发育，90%含气单砂岩厚度小于5 m，因此，砂岩的垂向变程设置为1～4 ms。垂向上含气砂岩空间非均质性极强，同一个井台不同定向井砂岩发育差异大，通过统计得出85%单砂岩的宽度介于300～1 500 m。因此，横向变程设置需要根据各个层段的沉积特征而定。富砂层段与富泥层段的变程、变差函数特征不同，可以根据单井分析结果或者是储层预测结果判断选择合适参数[39]。
图7 石盒子组测井统计岩性变差函数图
2.5 反演效果
基于MCMC法的相控地质统计学反演经过了15次随机模拟，得到了15个等概率的伪泊松比、纵波波阻抗、密度以及岩性数据体。图8是一条石盒子组连井反演结果剖面。其中LX-A井、LX-B井、LX-C井是反演井，参与反演约束过程；
LX-D井、LX-E井、LX-F是验证井，是后期新钻井，作为盲井用来检测反演符合率。过井反演剖面图中投影的曲线是自然伽马曲线以及测井岩性解释结果：蓝色代表泥岩，黄色代表砂岩储层。通过与过井点处的测井解释岩性结果对比可以看到，反演结果与实际钻遇的砂泥岩符合程度较高，不但可以准确识别出LX-A井盒四底部厚度为5 m的砂岩，还可以预测无井地区的砂岩展布，验证井LX-D井盒四段两套叠置砂岩也得到了很好的刻画，展示相控反演的高分辨率特征。
图8 石盒子组储层相控地质统计学反演典型剖面图
根据研究区岩石物理分析结果，纵横波速度比和泥质含量具有一定的相关性，纵横波速度比越大，泥质含量越大；
速度比越小，泥质含量越小。两者交会分析相关度为0.88，标准偏差为0.1。在砂岩内部，纵波波阻抗和孔隙度呈负相关关系，纵波波阻抗越大，孔隙度越小，纵波波阻抗越小，孔隙度越大，两者交会分析相关系数为0.91，标准偏差为0.02。同样的，在物性基本一致的情况下，密度可反映含气饱和度的变化，密度与含气饱和度呈负相关性。因此，可利用纵波波阻抗、纵横波速度比和密度3个弹性参数，对泥质含量、孔隙度和饱和度进行约束模拟。同时利用已钻井的地质统计学特征，对孔、渗、饱参数进行空间采样，获得各个变量之间的全局联合概率密度函数，来提高模拟的精度。
在此基础上，采用协模拟技术将储层相控地质统计学反演得到的伪泊松比和纵波波阻抗结合测井解释的孔隙度、渗透率、含水饱和度建立关系，通过协模拟进一步获得能够更加直观反映气层特征的物性参数。图9是协模拟得到的石盒子组连井物性参数剖面，图9-a是孔隙度剖面，图9-b是渗透率剖面，图9-c是含水饱和度剖面。从物性参数剖面分析，LX-A井主力气层为盒二段，发育厚度为5.6 m的气层，孔隙度为12.5%，渗透率为2.4 mD，含水饱和度为35.0%，该层段测试日产气量为17.3×104m3；
LX-B井主力气层为盒四段，发育厚度为4.9 m的气层，孔隙度为10.2%，渗透率为1.6 mD，含水饱和度为39.0%，该层段测试日产气量为3.3×104m3；
LX-C井盒六段，发育厚度为4.9 m的气层，孔隙度为11.2%，渗透率为1.2 mD，含水饱和度为40%，该层段测试日产气量为3.6×104m3。验证井LX-D井仅在盒六段钻遇2.0 m薄气层，该层段测试日产气量为0.7×104m3；
LX-E井在盒四段钻遇7.0 m气层，储层物性条件好，孔隙度为13.2%，渗透率为2.5 mD，含水饱和度为31%，盒四段气层单层压裂，测试日产气量为50.7×104m3；
LX-F井分别在盒二段、盒四段钻遇6.2 m、5.7 m气层，且储层物性条件较好，孔隙度为10.9%，渗透率为0.8 mD，含水饱和度为40%，两层分压合试，日产气量为2.6×104m3。反演井与验证井与实钻气层情况能够基本吻合，而且预测储层物性参数与测试产能情况关系密切，储层物性条件越好，产能越高。
图9 石盒子组协模拟物性参数剖面图
最后将地质统计学反演的岩性结果与储层物性结合，优选出孔隙度大于8%、含水饱和度小于40%的含气有利区，即可得到各个目的层段含气砂岩厚度图。图10为盒四段反演厚度预测平面图，图10-a是砂岩厚度平面图，图10-b是气层厚度平面图。从图10可以看出，在反演井中，LX-A井盒四段气层不发育，LX-B井、LX-C井盒四段发育厚度分别为5.0 m、4.5 m的气层，LX-B井与LX-C井均位于有利区内，反演预测与实际情况完全吻合。验证井LX-D井位于区块西南部，不在盒四段预测含气砂岩有利区内，实钻揭示该井盒四段虽然钻遇砂岩，但是储层物性较差，测井解释气层为0；
验证井LX-E井、LX-F井钻前预测位于含气砂岩有利区内，预测盒四段含气砂体厚度分别为7.0 m、5.7 m，实际钻井证实LX-E井、LX-F井位于河道主体部位，储层发育，钻遇的气层厚度分别为6.8 m、5.0 m，气层预测与实钻吻合较好。
图10 临兴中区盒四段反演厚度预测平面图
1）MCMC法通过不断迭代进行统计分析，具有效率高、精度高、适用范围广的特点，其算法的先进性保证了反演结果误差较小，对储层预测具有重要意义。
2）相控地质统计学叠前反演薄储层含气性预测技术流程特点在于利用岩性敏感弹性参数建立了三维空间储层先验约束，弥补了常规地质统计学反演横向不确定性强的特点。利用协模拟技术可以将相控地质统计学叠前反演得到的弹性参数体与测井解释的孔隙度等物性参数建立关系，能够更加直观地反映储层含气特征。
3）在相控地质统计学叠前反演的基础上，结合协模拟技术，在鄂尔多斯盆地东缘的临兴中区实现了厚度为5 m的薄气层的有效识别，经过新井验证，含气性吻合率情况较好，证明相控地质统计学叠前反演可以应用于致密砂岩薄储层含气性预测，有较高的可靠性，具有推广应用价值。
猜你喜欢概率密度气层岩性一种识别薄岩性气藏的地震反射特征分析云南化工(2020年11期)2021-01-14连续型随机变量函数的概率密度公式数学学习与研究(2020年15期)2020-11-28基于孔、裂隙理论评价致密气层石油地质与工程(2019年4期)2019-09-10K 近邻分类法在岩屑数字图像岩性分析中的应用录井工程(2017年1期)2017-07-31涪陵页岩气层构造变形量化研究获新进展录井工程(2017年1期)2017-07-31船底凹槽气层流动特征及影响因素研究武汉理工大学学报（交通科学与工程版）(2015年5期)2015-12-05Hunt过程在Girsanov变换下的转移概率密度的表示公式数学年刊A辑(中文版)(2015年1期)2015-10-30随机变量线性组合的分布的一个算法河北建筑工程学院学报(2015年2期)2015-04-29随机结构-TMD优化设计与概率密度演化研究振动工程学报(2015年2期)2015-03-01苏里格气田致密砂岩气层识别难点及方法评述天然气勘探与开发(2015年1期)2015-02-28