基于单一图像生成对抗神经网络方法在沉积相建模中的应用

李少华，史敬华，，于金彪，王 军，周传友，喻思羽

（1.长江大学地球科学学院，湖北 武汉 430100；2.中国石化胜利油田分公司勘探开发研究院，山东 东营 257015）

储层建模包括多个步骤，不同学者有不同的划分方法，如二步建模［1］、三步建模以及多步建模［2］，沉积相建模是其中的一个重要环节。通常会先建立沉积相模型，然后在相模型基础上建立储层物性参数模型。为了评估建模参数不确定性的影响，通常还需要建立概率储量［3-5］。一般来说，不同相带物性参数的分布特征是不一样的，在相模型基础之上建立物性参数模型能够更加合理地刻画物性参数的空间分布，中国学者称之为相控建模［6］，已得到了广泛的应用［7-8］，并且这个概念及应用领域也被进一步扩大，如多级相控建模、碳酸盐岩相控建模、地质模式约束建模等［9-11］。有多种建立相模型的技术，包括基于两点地质统计学的方法（序贯指示模拟、截断高斯模拟），基于目标的建模方法、多点地质统计学方法（MPS，Multiple-points Statistics Methods）和基于沉积事件（过程）的建模方法［12］。自1993 年提出MPS 方法，在相建模方面展示出良好的应用前景，该方法综合了两点地质统计学建模方法和基于目标建模方法的优点，既能较好地刻画复杂地质体的几何形态，又容易满足已有的条件数据。MPS 方法逐渐成为研究的热点，中外学者对经典的MPS 方法提出了很多改进［13-17］。MPS方法基于训练图像进行模拟，例如经典的基于概率的单一正规方程模拟（Snesim）方法，其核心技术包括：一是对训练图像建立搜索树，二是从搜索树中获取条件概率［18］；基于样式的Simpat 方法，通过计算数据事件与训练图像中样式的差异进行模拟［19］，能够更好地再现地质模式。针对一个工区通常基于一个训练图像进行模拟，除非是对该工区地质模式认识上存在较大不确定性，可以采取几张训练图像代表不同的地质模式［20］。基于训练图像MPS 方法的重点是提取训练图像的模式，但是这些算法不能生成真实的非线性模式，也不能展现地质推断的可变性和不确定性，DUPONT认为基于GAN模拟可以解决这些问题［21］。

1968年，IVAKHNENKO 公开了感知器（Perception）的概念，被认为是深度学习的起源［22］。2014年，IAN 等提出了GAN 方法［23］，主要包含生成器和判别器，两者相互博弈，该方法受到广大学者的关注，并在其基础上做了大量的改进研究。近年来，人工智能方法在油气领域的应用越来越广泛，例如以支持向量机回归和XGBoost 回归算法为基础的储层裂缝开度预测，结合卡尔曼滤波和人工神经网络进行油藏井间连通性研究，基于测井的岩相自动识别，应用多次数据吸收集合平滑算法进行自动油藏历史拟合［24-28］。ERIC 等提出基于训练图像的空间生成对抗神经网络地质统计学反演方法，该方法能够快速生成二维和三维的非条件模拟结果［29］。COIFFIER 等使用生成对抗神经网络将二维样本合成为三维地质图像，该方法实现简单，在标准GAN的生成器和判别器之间引入随机截取采样步骤［30］。ZHANG 等使用基于沉积结构的生成对抗神经网络方法生成逼真的三维相模型，该方法可以再现各种概念地质模型，同时具有满足约束条件（如井资料）所必需的灵活性［31］。目前利用GAN 方法建模需要大量的训练图像，基于15 000 张训练图像学习开展模拟预测，尽管可以得到与训练图像非常逼真的模拟结果，但是大量训练图像的生成与学习训练增加了该方法应用的难度。SinGAN 方法能够从单个自然景观图像生成新的模型，克服了传统的生成对抗神经网络需要大量的图像进行训练的限制［32］。为此，笔者尝试将SinGAN 方法应用到相建模中，以辫状河储层为例，并与经典的基于样式的Simpat 方法进行对比，SinGAN方法展示出良好的应用前景。

1 SinGAN方法基本原理

GAN 方法可以生成与期望地质模式高度匹配并且满足条件数据的模型，但是训练数据必须基于大量的特定训练集［24］。SinGAN 方法能够从单个自然图像中学习并生成模型［26］。

SinGAN 方法是一种金字塔式的全卷积轻量权重的GAN 方法，与常见的GAN 方法类似，其重要改进之处是训练样本采用的是单个训练图像，核心在于通过不同尺度来捕获复杂图像结构的统计特征。模型包含金字塔式的生成器和训练对应的金字塔式图像。降采样是一种图像处理技术，降采样的结果是原始图像的缩略图。训练过程是从最粗尺度开始，多次迭代到最细尺度。所有的生成器和判别器具有相同的感受野，因此随着尺度的逐步增加能够捕获到更小尺寸的图像结构。在最粗尺度上由生成器根据高斯分布的噪声生成的图像，其表达式为：

在最粗尺度下，有效感受野大约是训练图像尺寸的一半，因此生成器生成的是图像的布局和全局结构。每一个更细尺度的生成器会增添上一个生成器没有生成的细节部分。因此，除了空间噪声，每一个生成器还要接受上一个尺度的上采样图像数据，其表达式为：

所有的生成器结构类似，特别之处在于在噪声进入卷积层之前，噪声会被添加到上一个尺度生成的图像数据中，这既保证GAN 方法未忽视噪声，又避免完全随机。卷积层目的是生成上一个尺度所缺乏的细节，其表达式为：

卷积核大小为3×3，在最粗尺度上每个卷积块使用32 个卷积核，每增加一个尺度，卷积核数量增加16个。因为生成器是全连接的卷积网络，通过修改噪声的维度，可以生成任意大小和长宽比例的图像。

训练过程采用从最粗尺度到最细尺度的顺序模型。一个GAN 方法训练完毕，就保持固定状态。训练损失包括对抗损失和重建损失，其表达式为：

对抗损失用来减少训练图像采样后的图像中块的分布和生成样本中块的分布之间的差异，重建损失确保存在可以产生训练图像降采样后的图像的特定噪声图集。每一个生成器与一个马尔科夫判别器耦合，判别器将输入的每个重叠块分为真和假两类。SinGAN 方法定义的损失是整个图像的损失，这样网络可以学习边界条件。马尔科夫判别器的结构和生成器中的由BatchNorm-LeakyReLU 的5个卷积块组成的全连接卷积网络相同，因此块的大小也相同。要生成原始的图像x，必须保证存在一组特定的输入噪声集。重建损失的表达式为：

当尺度数量为N时，重建的图像在训练的过程中还有一个重要的作用，重建的图像决定了每一个尺度上噪声的标准差。

SinGAN方法的结构尺度数量会影响模拟结果。当尺度数量较少时，粗尺度的有效感受野更小，只能捕获细微的特征。随着尺度数量增加，更高质量的结构和整体目标布局会更好地保存下来，图像结构随着尺度数量增加而变化。

2 实例应用与分析

N 气田位于中国东海陆架盆地某凹陷的中央反转构造带中北部，为中国东海大型低渗透气田，已经勘探的3套气藏储层的厚度均大于100 m，总气藏的厚度大于400 m。其主力层段花港组上段主要发育灰色泥岩、粉砂质泥岩和泥质粉砂岩等，且夹有薄层细砂岩；花港组下段主要发育厚层细砂岩，局部发育中粗砂岩，并夹杂少量粉砂岩及泥岩。

为了检验提出方法的合理性，选取了N 气田2个有代表性的小层沉积微相图作为训练图像，包括河道、心滩和泥岩，训练图像的大小为256×159个网格单元。其中图1a 平稳性较好，图1b 较差，在训练图像上部与下部具有较大差异，上部的河道间泥岩更为发育。为对比该方法的适用性，采用2 组不同的训练图像。测试硬件环境包括：CPU 采用第八代英特尔酷睿i7 处理器，频率为3.7 GHz，运行内存为16 G，显卡采用英伟达GeForce GTX1 080 Ti。

图1 2个小层的沉积微相图Fig.1 Microfacies maps of two layers

2.1 利用MPS方法建立沉积微相模型

Simpat方法由ARPAT 提出［19］，与基于概率统计的Snesim 方法不同，该方法将图像重建的思想引入到储层建模，主要模拟步骤包括：①以数据样板扫描训练图像，构建样式库。②在模拟过程中，计算数据事件与样式库中所有样式的相似性，Simpat 方法采用Manhattan 距离计算数据事件与样式的相似度。③将最相似的样式整体粘贴到模拟网格中。④继续模拟，直至所有网格都被模拟。

利用Simpat 方法分别采用3×3，5×5，7×7，9×9，11×11，13×13，15×15，17×17 共8 组数据样板进行测试。由基于相对平稳训练图像Simpat 方法模拟结果（图2）和基于非平稳训练图像的模拟结果（图3）可以直观地看出，数据样板大小对于模拟结果的影响很大，数据样板过小模拟结果很离散，过大容易失真，且计算量增加。对于MPS 方法来说，合理确定数据样板大小对于建模质量的好坏具有重要的影响。另外，相对平稳的训练图像模拟效果比非平稳的训练图像好。

图2 基于相对平稳训练图像的Simpat方法模拟结果Fig.2 Simulation results of Simpat based on stationary training image

图3 基于非平稳训练图像的Simpat方法模拟结果Fig.3 Simulation results of Simpat based on non-stationary training image

2.2 利用SinGAN方法建立沉积微相模型

基于SinGAN 方法的相关理论，搭建了相应的PyTorch 环境。初始化参数设置包括噪声、尺度因子、输出路径、种子数等，本次设置尺度因子为0.75，种子数使用程序随机产生。训练的过程会根据尺度数量以金字塔式的顺序依次从最粗尺度模拟至最细尺度，每个尺度迭代次数可以人为设置，本次设置为2 000 次，目标函数重建损失前的系数设置为10。SinGAN 方法除了具有GAN 方法共有的核心生成器和判别器，还需要进行尺度数量的选择。尺度数量的选取一种是自适应的方式，主要是根据图像的像素大小计算选取，另一种方式是人工指定尺度数量。

首先利用SinGAN 方法对相对平稳的训练图像进行模拟，选择图1a相对平稳的辫状河训练图像进行测试。分别采用尺度数量为3，5 和7 进行模拟计算，设置成图的数量为50 幅，部分模拟结果如图4所示。

从图4 中可以看出，N值对模拟结果的影响很大，N值较小的时候，会出现大量的过渡色，因为SinGAN方法是针对自然景观图像，是连续变量。当N为7 时，效果良好，能够刻画不同沉积微相的结构和空间展布规律，直观对比发现，优于经典Simpat方法模拟结果。

图4 基于相对平稳训练图像的SinGAN方法模拟结果Fig.4 Simulation results of SinGAN based on stationary training image

为进一步验证SinGAN 方法对非平稳训练图像的模拟效果，选择图1b非平稳的辫状河训练图像进行测试。同样采取N为3，5和7进行模拟计算，模拟结果如图5所示。

由图5 可以看出，N值对模拟结果的影响更加明显，N为3时，明显失真，N为5变形也比较严重，N为7 效果较好。进一步增加N值容易导致系统崩溃，一般不要超过7［26］。对比Simpat 方法模拟结果可以看出，SinGAN 方法优于经典的MPS 方法，与训练图像刻画的沉积微相空间结构更相似。因为Sin⁃GAN 方法是基于自然景观图像，利用该方法建立沉积相模型会出现一些过渡色，即使当采用N为7 时也会在个别地方出现，且非平稳的训练图像这种现象更明显一些。这个问题可以采用聚类的方法进行后处理，也就是对带过渡色的图像根据RGB 值进行聚类，判断每个像素的颜色与3 种颜色的差异，将最接近的那种颜色赋值给过渡色的像素，从而得到只有3种沉积相类型的模拟结果。非平稳训练图像基于N为6 时得到模拟结果3 的聚类后处理结果（图6），在不同沉积相类型的接触边界，原始模拟结果（图6a）呈现过渡色，通过聚类后处理的结果（图6b）有效地消除这种过渡色，模型只保留了3种沉积相类型。

图5 基于非平稳训练图像的SinGAN方法模拟结果Fig.5 Simulation results of SinGAN based on non-stationary training image

图6 聚类前后结果对比Fig.6 Comparison of results before and after clustering

2.3 模拟结果定量对比

为了进一步比较SinGAN 方法与经典MPS 方法的样式再现能力，选取了基于Simpat 方法模拟效果较好的2 个模拟结果（分别对应数据样板为13×13和15×15），N为7的SinGAN 方法模拟结果和原始训练图像（相对平稳）进行样式相似性对比分析。将模拟结果中的样式全部提取出来，与训练图像所包含的样式进行比对，采用逐网格属性比对的方式刻画样式的相似性。其表达式为：

通过计算模拟结果的每个样式的相似性，即可统计比较模拟结果与训练图像的样式相似性差异。

由对比结果（图7）可知，SinGAN 方法模拟结果与原始训练图像的样式相似性最好，其次是基于13×13 样板模拟的结果。通过对比分析表明，Sin GAN 方法在再现训练图像内部结构方面具有强大的能力，优于经典的MPS 方法，展示了良好的应用前景。本次仅开展了二维非条件的沉积相模拟，进一步实现条件化的三维模拟将极大推进该方法在储层建模中的实际应用。

图7 数据样式相似性对比Fig.7 Comparison of similarity of data patterns

3 结论

基于单幅自然景观图像的SinGAN 方法应用于沉积相建模，取得了良好的效果，通过直观分析及数据样式再现的定量化对比，SinGAN方法模拟的结果优于Simpat 方法。SinGAN 方法模拟结果不论是相对平稳的沉积微相图还是非平稳的沉积微相图受模拟尺度数量的影响明显。相对而言，对于非平稳的沉积微相图影响更大，需要选取更大的尺度数量才能够达到较好的模拟效果。对于沉积相建模一般选取最粗尺度数量为7，过大可能会导致系统崩溃。当SinGAN 方法用于建立离散的沉积相模型时，即使采用最粗尺度数量为7 也会在个别地方出现一些过渡色，非平稳情况下更为明显。针对这种情况可以采取聚类的方式处理，消除过渡类型，达到模拟离散模型的目的。SinGAN 方法是针对连续变量的，可以直接应用于连续变量的模拟。目前主要进行了二维非条件的沉积相模拟测试，三维空间的模拟可以借鉴已有做法，需要考虑不同尺度相互匹配的问题。如何将条件数据、地震属性等更多的约束信息加入到SinGAN 方法仍是一个难点，但该方法在储层建模领域已经展示出良好的应用前景。

符号解释

Dist——计算模拟结果的第i个样式与PatDBT中的每个样式的距离值，无量纲；

Dn——马尔科夫判别器；

Gn——生成器；

GN——金字塔式的生成器；

i——编号；

Ladv——对抗损失，无量纲；

Lrec——重建损失，无量纲；

n——当前尺度数量；

N——最粗尺度数量；

Pat(i)——模拟结果的第i个样式；

PatDBT——训练图像里所有的样式；

S(i)——训练图像中与模拟结果的第i个样式最相似样式的相似度；

xn——当前尺度数量为n时的图像；

——对应xn生成器生成的图像；

x͂rnec——当前尺度数量为n时噪声集生成的图像；

x͂N——最粗尺度数量为N时生成的图像；

zn——当前尺度数量为n时的空间噪声图像；

zN——最粗尺度数量为N时的空间噪声图像；

α——目标函数重建损失前的系数，无量纲；

↑r——重采样；

ψn——BatchNorm-LeakyReLU 的5个卷积块组成的全连接卷积网络。