基于非下采样剪切波变换—参数自适应脉冲耦合神经网络的属性融合裂缝预测方法

汤 韦 李景叶* 王建花 薄 昕 耿伟恒 叶 玮

(①中国石油大学(北京)地球物理学院，北京 102249； ②油气资源与探测国家重点实验室，北京 102249；③海洋石油勘探国家工程实验室，北京 100028； ④华北油田公司勘探开发研究院，河北任丘 062552)

0 引言

裂缝是影响油气藏岩石物性参数的重要因素，精确预测裂缝对油气勘探、开发具有重要意义[1-2]。不同尺度的裂缝对地下油气藏的影响程度不同，基于叠后地震数据预测的大尺度裂缝(断裂)通常改变油气的储存状态[3-5]，因此精确预测大尺度裂缝是油气勘探的重要组成部分[3]。

叠后地震数据与储层地质信息具有相关性[6-8]，大尺度裂缝使叠后地震数据波形明显突变或不连续[4]。提取有效的叠后地震属性，可从几何动力学角度预测储层中的大尺度裂缝。常用的叠后地震属性主要有相干体(描述波形相似性)[9-13]、曲率体(表征构造应力引起的地层弯曲程度)[14-16]、倾角体(刻画地层构造变化特征)[17]等。然而，仅仅依靠单一属性很难准确预测地下裂缝分布情况，因此属性融合技术显得极其重要。通过数学方法综合考虑多种属性特征，将它们组合为更全面、可靠的裂缝表征属性，可以降低单一属性预测的多解性[18-20]，提高裂缝预测精度。

目前，属性融合方法较多，且应用较广泛。在属性色彩技术方面，Guo等[21]将属性分量映射到色彩模型中，通过红、绿、蓝(Red-Green-Blue, RGB)和色调、亮度、色饱和度(Hue-Intensity-Saturation, HIS)模型进行属性融合；为了准确确定融合属性，丁峰等[22]结合主成分分析与红、绿、蓝、透明度(Red-Green-Blue-Alpha, RGBA)颜色融合方法分析属性融合，很好地识别了断裂区域。在神经网络方面，徐丽萍[23]利用人工神经网络方法进行非线性融合，以预测储层的发育程度。曹琳昱等[24]将粒子群优化的BP神经网络方法用于属性融合，但需要优选属性训练网络，且会陷入极小值。李全忠等[25]首次利用脉冲耦合神经网络(Pulse Coupled Neural Networks, PCNN)进行地震属性融合，以指示油气发育情况，然而PCNN模型中的大量参数需要人为干预设置，导致融合结果不准确及不能准确突出单属性的细节特征。因此，参数自适应PCNN(Parameter-Adaptive PCNN, PA-PCNN)模型[26]应运而生并用于图像融合领域，很好地避免了人为设置参数带来的误差。此外，基于多尺度几何分析的PCNN模型在图像融合领域的应用效果很好，通过非下采样小波变换[27]或非下采样剪切波变换[28](Non-Subsampled Shearlet Transform,NSST)分割原图像，可更好地提取细节信息进行融合。基于NSST-PCNN模型的多信息融合技术广泛用于医学[29-30]、遥感[31-32]等领域，且融合效果尤为突出。

本文在前人理论研究的基础上，结合NSST和PA-PCNN的优点，建立了基于NSST-PAPCNN的属性融合裂缝预测方法，以降低属性多解性，消除原属性中的部分冗余成分，可最大限度地保留有用信息，为地质解释提供相对准确的裂缝预测数据。首先，通过NSST将多属性分解为高、低频子带，其中高频子带包含更多的裂缝细节信息，低频子带可更好地刻画裂缝轮廓且具有丰富的能量信息。其次，对高频子带运用PA-PCNN模型进行融合，无需人工设置参数，得到更全面的高频数据；结合八邻域的改进拉普拉斯算子加权和(Weighted Sum of Eight-neighbor-hood-based Modified Laplacian，WSEML)与局部能量加权(Weighted Local Energy，WLE)方法对低频子带进行融合，使低频数据更好地保留细节及能量信息，以得到丰富的低频数据。最后，通过逆NSST方法完成属性融合裂缝预测。应用实例较好地证明了本文方法的有效性，通过对比不同方法的属性融合裂缝预测效果，突出了所提方法的优势。

1 理论方法

1.1 NSST

属性融合过程中，为了获得细节信息更丰富的结果，通过数学变换进行多尺度、多方向分解提取相同位置的不同属性、不同空间细节信息进行融合。NSST具有位移不变性及方向灵敏性等特点[32]，计算相对简单，能够准确提取目标边缘及空间信息，在医学和图像处理领域应用广泛。因此，将其应用于裂缝单属性分解，提取细节信息进行融合，可以更全面、准确地预测裂缝。

图1 NSST分解框架

对高、低频子带分别融合后的数据运用逆NSST进行数据组合重构，主要分两个步骤：首先，通过SFBs对各个方向分解的滤波结果进行累加，生成非下采样金字塔；然后，使用重构滤波器对非下采样金字塔由粗到细地进行数据组合重构，获得融合后的数据结果。

1.2 PA-PCNN高频子带融合

PCNN是一种基于文献[33]的皮质模型，是依据哺乳动物大脑皮层中同步脉冲发放的现象提出的，具有全局耦合和脉冲同步等特点，适用于数据融合领域。

PCNN模型是由若干神经元连接而成的单层反馈网络[34]，输入的单个属性中数据的个数与神经元一一对应，每个神经元都受相邻神经元影响，从而进行信息传递和耦合，且主要由接收、调制和脉冲发生三部分组成。在传统的PCNN模型中，存在大量的自由参数设置且非线性较强，求解过程繁琐、计算量大，因此简化的PCNN模型[26]应运而生。

简化的PCNN单个神经元工作流程(图2)可描述为

图2 简化的PCNN工作流程

Fij(η)=Sij

(1)

(2)

Uij(η)=e-αfUij(η-1)+Fij(η)[1+βLij(η)]

(3)

(4)

Eij(η)=e-αeEij(η-1)+VEYij(η)

(5)

式中：Fij(η)为第η次迭代后(i,j)位置的反馈输入项；Sij为(i,j)位置输入数据；Lij为连接项，将周围神经元依据不同权重值相连；Ynl为周围(n,l)位置神经元的输出脉冲值；Wijnl为(i,j)位置输入数据周围(n,l)位置神经元的权重；Uij为内部活动项，e-αf×Uij(η-1)为上一次内部活动的衰减项，e-αf为Uij(η-1)的指数衰减系数；Fij(η)[1+βLij(η)]为反馈输入和连接项的非线性调制，β为连接系数；e-αe为Eij(η-1)的指数衰减系数；Eij为动态阈值项；VL和VE分别为连接项和动态阈值项的幅值。

在式(1)～式(5)中，主要涉及五个自由参数(VL、αf、β、αe、VE)和一个权重矩阵(Wijnl)。在常规方法中，这些参数往往都由人为设定，会导致不准确的融合预测结果。因此本文运用PA-PCNN方法，无需人为设定，根据输入数据自有特征自适应地求解参数进行融合处理。分析式(1)～式(5)可知，β和VL可以整体视为连接项的权值系数，因此定义κ=βVL为连接权重强度。此外，权重矩阵主要计算周围神经元对目标神经元的影响作用，可通过距离定义权值。根据文献[26]可知，所有自由参数都通过自适应计算获得，即

(6)

(7)

(8)

VE=e-αf+1+6κ

(9)

(10)

为了便于求解计算，需将权重矩阵近似为

由上述PA-PCNN过程得到每个属性迭代后在所有神经元处产生的脉冲点火次数，由此可以获得各个神经元处的融合权重系数

(11)

(12)

1.3 WSEML+WLE低频子带融合

经过NSST分解后的各属性数据低频子带包含大量能量和细节信息，因此低频段的融合策略对最终融合质量具有重要影响。传统方法通常直接对低频数据进行简单加权平均融合，然而不同属性在同一位置的能量不同，易造成低频融合结果能量损失。为了在低频融合数据中尽可能保留原数据能量和细节信息，引入WLE和WSEML算法融合低频子带数据[35-36]。

WLE算法根据区域内各数据之间的相关性，最大限度地保留原数据能量信息

LS(i+u,j+v)2

(13)

式中：W是大小为(2R+1)×(2R+1)、半径为R的权重矩阵，矩阵中每个元素值设为22R-r，其中r是(i,j) 位置到矩阵中心四邻域的距离；LS为NSST分解后各属性的低频子带。以一个3×3阶矩阵为例，设R=1，r=0,1,2，则

由于计算效率问题，对NSST分解层数进行一定限制，因此低频子带中仍然包含一定细节信息。为了尽可能地保存原数据中的细节信息，引入WSEML算法

R+1)×EMLS(i+u,j+v)

(14)

其中

EMLS(i,j)=|2LS(i,j)-LS(i-1,j)-LS(i+

1,j)|+|2LS(i,j)-LS(i,j-1)-LS(i,j+1)|+

(15)

结合WLE和WSEML算法，充分提取了原属性数据中的能量和细节信息，可求取低频子带融合权重，最终获得低频子带数据

(16)

(17)

1.4 算法实现流程

根据上述理论算法，总结了NSST—PAPCNN属性融合裂缝预测结构框架的具体流程(图3)。

图3 NSST-PAPCNN属性融合裂缝预测结构框架

(1)从叠后地震数据中提取表征大尺度裂缝的属性(包括相干A、最大曲率B和倾角C)，并对各属性数据进行归一化处理。由于不同的属性其量纲及数量级不同，不能直接融合，需要将多属性数据归一化到[0，1]。

m∈[1,M(k)]

数据分解过程中需要确定分解的层数K及每层方向数M(k)。

(4)对各属性低频子带进行WLE+WSEML融合，最大限度地保留能量和细节信息。融合前，需确定式(13)中的权重矩阵半径R，最终得到融合的低频子带LFusion。

2 应用实例

为了验证NSST-PAPCNN多属性融合裂缝预测方法的可行性和有效性，对M区地震数据进行测试，数据范围为Crossline1～Crossline751、Inline1～Inline751。该区裂缝较发育，大尺度裂缝(小断裂)分布广泛。为了有效预测大尺度裂缝，从叠后地震数据中准确提取了目的层段的相干、最大曲率和倾角等三种属性，并分别做归一化处理，然后进行属性融合处理，具体过程(图4)如下。

图4 M区多属性融合步骤

(1)将三种属性分别进行NSST多尺度、多方向分解处理。测试并分析该区的数据特征可知，分解的层数超过4时，计算较稳定，但是层数过多会导致计算耗时较长。经权衡计算效果及效率，设定分解层数为5，每层的分解方向由细到粗为{16,16,8,8,4}。

(2)对分解后的高频子带数据进行PA-PCNN融合。测试结果表明，迭代次数过多，融合结果差异不明显。因此通过试验，选取迭代次数为110。融合前，需要将PCNN模型中的各项参数初始化为0，并根据实际数据特征自适应地计算各融合参数。通过点火次数确定各属性高频子带数据权重系数，融合得到最终高频子带。

(3)对分解后的低频子带数据进行WLE+WSEML融合。综合考虑计算效率及对结果的影响程度，设定式(13)中权重矩阵半径R为1，计算相应的WLE及WSEML值，并确定各属性低频子带数据权重系数，融合得到最终低频子带。

(4)对分别融合后的高、低频子带数据进行逆NSST重构，获得最终的属性融合结果。

由M区归一化叠后均方根振幅属性切片(图5)可以大致看出裂缝的发育特征，以此验证融合结果。为了更全面、准确地预测裂缝，需要最大限度保留三种属性的轮廓、能量及细节特征。图6为M区目的层段归一化属性切片，图7为相干、最大曲率和倾角多尺度、多方向分解的部分切片。由图可见，本文的NSST-PAPCNN方法能够对三种属性进行多尺度、多方向分解，充分提取了裂缝细节及能量信息，分别获得了低频(图7a)、浅层分解(粗化分解，图7b)、中层分解(逐渐细化，图7c)及深层分解(细化分解，图7d)等信息，从而得到了全面表征裂缝信息的融合结果。

图5 M区归一化叠后均方根振幅属性切片

图6 M区目的层段归一化属性切片(a)相干；(b)最大曲率；(c)倾角属性数值越大，则裂缝越发育，红框区域展示了三种属性的明显差异

图7 相干(左)、最大曲率(中)和倾角(右)多尺度、多方向分解的部分切片(a)低频分量；(b)浅层分解；(c)中层分解；(d)深层分解

为了验证NSST—PAPCNN方法的属性融合裂缝预测效果，对比了不同方法的属性融合裂缝预测结果(图8)。对比图8与图5可见， NSST-PAPC-NN方法具有明显的优势，消除了部分冗余成分(图8a红框区域)，很好地保留了原属性的能量、轮廓及细节信息(图5)，全面刻画了大尺度裂缝的发育情况，为地质解释人员提供了较精确的裂缝预测数据。实际数据测试结果表明，本文提出的基于NSST-PAPCNN的属性融合裂缝预测方法能够更有效地预测裂缝。

图8 不同方法的属性融合裂缝预测结果(a)NSST-PAPCNN；(b)PCNN；(c)RGB红框区域展示了三种结果的明显差异

3 结论与建议

本文基于NSST及PA-PCNN方法提出了一种属性融合方法用于裂缝预测。该方法基于NSST分解算法，将多种属性数据分解为高、低频子带，对融合后的多尺度、多方向高、低频子带进行数据重构，得到最终的多属性融合结果，可进一步提取裂缝的轮廓及细节信息。

针对包含丰富细节信息的高频子带数据，采用PA-PCNN算法自适应调节自由参数，降低人为设置带来的误差。根据各子带数据各元素的点火脉冲产生次数，确定融合系数进行高频子带融合。针对包含能量、轮廓及残存细节信息的低频子带，结合WLE和WSEML算法，最大限度保留低频子带的能量和细节信息，确定融合系数进行低频子带融合。

将融合后的高、低频子带数据进行逆NSST重构，获得最终的融合结果，降低了原属性的多解性，消除了部分冗余成分，很好地保留了原属性的细节和能量信息，能够综合多属性特征，全面、准确地预测裂缝，从而为地质解释人员提供较为精确的裂缝预测数据。运用本文方法对M区属性数据进行测试，并对比了不同方法的属性融合裂缝预测结果，证明基于NSST-PAPCNN的属性融合裂缝预测方法能够更有效地预测裂缝。

尚需指出，所提方法由于需要经过NSST多尺度、多方向分解，并对各尺度、各方向数据分别进行属性融合，导致计算量较大，耗时较长。因此，经权衡计算效率和融合效果，认为NSST分解的层数不能过多。此外，在尽量保留原属性中的有用细节成分时，往往无法保证融合后的所有细节信息都可靠，因此还需要依靠地质解释人员的专业经验及实际工区的数据特征准确预测裂缝，以进一步降低细节信息的多解性。