CN106707335A - 一种叠后地震信号波形分类方法 - Google Patents

Abstract

本发明提供一种叠后地震信号波形分类方法，包括如下步骤：(1)提取待分析的地震数据，并进行预处理；(2)使用稀疏编码方法提取地震数据中的特征数据x，并对特征数据x进行降维；(3)将步骤(2)中得到的特征数据x利用聚类算法进行分类。本发明通过稀疏编码降维较好地解决因数据维度变大而影响分类结果本发明引入的新聚类方法，相对于传统自组织神经网络方法能自动确定合适的分类数目，而且对于任意类型和特征的数据，都能获得较好的分类结果。

Description

一种叠后地震信号波形分类方法

技术领域

本发明属于油气储层预测和勘探技术领域，具体涉及一种叠后地震信号波形分类方法。

背景技术

地震波形分类技术在油气储层预测和勘探中发挥了举足轻重的作用，成为石油勘探领域中的重要技术之一。波形分类在地质研究中也发挥着很重要的作用，特别在识别沉积相结构方面。波形分类算法的最终效果取决于在进行分类的时候是否已经拥有当地的实际地质信号以及该地质信号含有的准确地震相个数。在波形分类未引入到石油勘探之前，人们对储层发育有利相带、油气分布的预测主要是通过有经验的地质工作者按照已有的经验及知识来进行分析，采取这样的分析手段需要大量的人力、物力，并且由于人为的主观干预过多，得出的结果不够有说服力，并且实际应用也往往出现问题。

在地震学引入到石油勘探的过程中，人们慢慢发现自采集的地震数据中获取的特征信息能够在石油勘探和开采过程中起到很大的作用，并且在降低人力、物力、成本投入的同时加快了勘探开发速度。三维地震信号波形分类作为地震学的一个重要分支，在沉积相描述、储层预测、油气识别等方面起到了关键作用。三维波形分类通过分析地下信号(这些信号是通过放炮，然后通过传感器采集返回的信号振幅)，然后通过对这些信号进行特征提取，用提取出来的特征来表征原来的信号，再根据这些特征的相似度，利用分类方法把信号分成预先设定类中的某一类，根据它们在平面上的位置信息生成地震相图，地球物理学家和地质学家可以根据地震相的分布规律，揭示出地下的沉积相分布、储层和油气分布等。

传统的信号分类方法在数据维度变大或数据类型复杂时，往往不能得到正确的分类结果。

目前在无监督的波形分类算法中自组织神经网络能够获得较好的结果，基于自组织神经网络的无监督波形分类算法已经成为地震相分析的一种主流算法，其主要是通过不断训练来调节神经元的大小，调节的方向是使得获胜神经元在下一次训练获胜的概率更大，且使得获胜神经元周围的神经元也能在下次训练中优势突出，该算法分类效果较好，且收敛速度快。

稀疏编码的研究起源于人类视觉感知的研究，研究者发现人类视觉神经元对于外界的颜色、尺度、方向、纹理等不同的刺激有着不同的反应，且某些刺激只会导致少数的神经元活动，其他大部分则没有响应，这表现为人脑的稀疏性与简约性。在对语音信号的研究实践中，科学家发现一段信号总是可以分解为一系列子信号的加权，同时这些子信号的系数满足稀疏的特性，即大部分系数为零，信号只由少数几个子信号构成，这与人类视觉机理是一种的，促进了人们对稀疏编码的进一步研究。如果能从信号集合上提取到一组过完备的字典，利用字典中的少数原子信号来表述一个复杂信号，则能实现对信号的稀疏表达，基于稀疏表达的方法实质上省去了传统方法的复杂信号的特征提取过程，却能获得较好的效果。

从数学的角度，稀疏编码用来寻找一组“超完备”的基向量来高效的表示样本数据，假定输入数据X＝(x₁,x₂,…，x_n)^T为n维向量，通过找到这样一组基向量{φ_i}使得：

因{φ_i}是一组超完备的基向量，因而系数a_i可以满足一定的稀疏性要求，即其中大部分为零。超完备基的好处是能够方便地找到数据集中隐含的内部结构与模式，但是带来的问题是系数并不唯一。稀疏性要求的另一个解释是在实际的信号中，图像可以是少量“基本”元素的叠加。m个输入向量的稀疏编码代价函数定义为：

其中，S()为惩罚函数，用来对稀疏性进行约束；式中第一项为重构误差项，这一项使得稀疏编码算法为输入向量X提供一个高拟合度的线性表达式；第二项迫使在X的表达式变得稀疏；λ为平衡因子，用它来调节两项的相对重要性。对于惩罚函数S()而言，最直接的测度应该是L₀范式，即：

S(a_i)＝1(|a_i|＞0) 式(1-3)，

但是由于L₀范式不可导，最常用的仍然选择L₁范式：

S(a_i)＝|a_i|₁ 式(1-4)，

或对数惩罚函数：

S(a_i)＝log(1+a_i ²) 式(1-5)。

但是这样还不够，如果不限制φ_i的范围，那么通过增大φ_i就会有可能使得对a_i的约束变得无意义，故还需要增加对基向量φ_i的约束，得到最终的损失代价函数：

关于稀疏编码的解释除了以上基于超完备基的理论之外，还有从概率角度出发，将图像或信号看做一个生产模型，其最终的结果仍然能回到以上的代价函数上，在此不做更多介绍。

KSVD是稀疏编码理论的具体实现，最早由Michal Aharon提出，算法通过对误差不断降低，使其达到最小。具体的方式是对重构误差项进行SVD(Singular ValueDecomposition)分解，找到其最大的特征向量作为字典的列，是一种典型的迭代优化算法，KSVD算法的具体目标函数为：

其中，Y为训练数据集合，D为超完备的字典，λ为系数，F为常见范数的定义，T₀为给定参数；

为了求解最优化的D和λ采用迭代逼近寻优方法，具体的迭代步骤为：

输入：训练数据集合Y{y₁,...y_i...y_m}，其中y_i为第i个训练信号。

输出：超完备的字典D{d₁,...d_i...d_n}，其中d_i为字典D的第i个原子向量。

Step1：固定字典D，优化系数λ，则问题可以转化为一个稀疏重构问题，典型的算法有MP(Matching Pursuit)或者OMP(Orthogonal Matching Pursuit)算法。

Step2：利用上一步得到的λ，重新更新字典，将原来的目标方程改写为如下形式：

其中，d_k代表字典的第k列(这里单独把字典的第k列列出进行优化，问题转化为选取什么样的d_k可以最大化的降低误差)。

Step3：对E_k进行SVD分解，求取其最大特征值对应的特征向量作为d_k。

Step4：重复以上步骤，直到更新完字典所有列。

上述算法中使用的匹配追踪，或正交匹配追踪算法，被用来寻找已知字典的情况，找到最有稀疏解。首先在原子库中寻找与原信号内积最大的原子，即该原子与原信号最相似，找到该原子后再将原信号减去该原子，再重复迭代，直到满足一定条件。通过匹配追踪算法可以找到基于原子库重构的最佳稀疏解，但是匹配追踪算法仅能保证信号的残余与最近选择的原子正交，却不能保证与前面选择的所有原子均正交，这限制了收敛的速度。在正交匹配追踪中，改进了残差更新公式，使得残差与前面已经选取的所有重构原子正交。

基于密度峰值搜寻的聚类方法由Alex Rodriguez和Alessandro Laio提出，该算法核心思想对于聚类中心有两个要求：1、密度大于邻居点密度；2、与密度更大点间的距离相对更大，下面详细描述该聚类算法：

数据集标签集I_S＝{1,2,...,N}，两数据点x_i与x_j之间距离d_ij＝dist(x_i,x_j)，可以取欧式距离，曼哈顿距离等。对任意数据点x_i定义两个属性，局部密度ρ_i和距离δ_i。

1、局部密度ρ_i：

2、距离δ_i：表示的一个降序排列下标，即

定义：

每个点都能计算出两个属性值，我们定义：

γ_i＝ρ_iδ_i 式(1-12)，

上式表明，γ_i越大说明该点越有可能是聚类中心点。

发明内容

本发明要解决的技术问题是提供一种叠后地震信号波形分类方法，通过稀疏编码降维较好地解决因数据维度变大而影响分类结果本发明引入的新聚类方法，相对于传统自组织神经网络方法能自动确定合适的分类数目，而且对于任意类型和特征的数据，都能获得较好的分类结果。

为解决上述技术问题，本发明的实施例提供一种叠后地震信号波形分类方法，包括如下步骤：

(1)提取待分析的地震数据，并进行预处理；

(2)使用稀疏编码方法提取地震数据中的特征数据x，并对特征数据x进行降维；

(3)将步骤(2)中得到的特征数据x利用聚类算法进行分类。

其中，所述步骤(1)的具体步骤包括：

(1-1)对算法中输入的宽方位角地震数据进行降噪；

(1-2)提取待分析层位段地震数据，每个层位段地震数据的样本点数t，组成训练样本：

X＝{x₁,...x_n}，

其中，x_i为维度t的向量。

其中，所述步骤(2)中：

输入：字典D，目标信号y，稀疏度K，初始化残差r₀＝y，索引集Λ＝φ，t＝1；

输出：y在字典D中的K稀疏逼近所对应的稀疏系数；

具体步骤包括：

(2-1)找出残差与字典矩阵中列内积最大的项，其列号记为λ，即：

(2-2)更新索引集Λ＝Λ∪{λ}，并根据更新后的索引集构建新的重构向量集合

(2-3)对向量Φ_t和x_t进行处理，通过最小二乘得到：

x_t＝argmin||y-Φ_tx_t||₂ 式(2-2)；

(2-4)更新残差r_t＝y-Φ_tx_t，t＝t+1；

(2-5)重复步骤(2-1)-(2-4)，直到满足停止条件：t＞k。

其中，所述步骤(3)包括下述具体步骤：

(3-1)初始化

(3-1-1)给定用于确定截断距离d_c的参数t∈(0,1)；

(3-1-2)计算任意两样本数据间的距离

(3-1-3)将升序排列，计算截断距离

(3-1-4)计算样本数据的局部密度使局部密度降序排列的下标x_i排列生成即

(3-1-5)计算其中，表示所有局部密度比x_i大的数据点x_j中与x_i距离最近的数据点的编号，

(3-2)确定聚类中心并初始化数据点归类属性标记

(3-3)按照ρ从大到小遍历数据对非聚类中心数据点归类；

(3-4)若簇的个数N_c>1，则将每个簇中数据点进一步分为core和halo，core的局部密度大，对应簇核心部分，halo的局部密度小，对应簇边缘部分，设halo部分为若h_i＝0，则x_i属于core部分，反之x_i属于halo部分。

进一步，所述步骤(3-4)的步骤包括：

(3-4-1)初始化h_i＝0；

(3-4-2)为每一个簇生成一个平均局部密度上界

(3-4-3)标记簇halo部分。

本发明的上述技术方案的有益效果如下：

1、对于某些高维叠后地震波形数据，只利用传统方法难以取得较好的效果，本发明采用稀疏编码既能很好地保留原始数据中蕴藏的特征，又能降噪以免对后续算法影响。

2、利用由Alex和Alessandr提出的基于密度峰值搜寻聚类算法，能够实现分类数的自动确定，并能保障得到良好的分类效果。所谓的分类数的自动确定是指：一般的聚类算法需要指定分类数如k均值，而本算法不需要指定分类数。本算法会指定聚类中心，当一个簇中两个点都符合作为聚类中心的条件，该簇当然也就会分成两个簇。这样不断迭代直到没有新的合适的聚类中心，最终聚类结果也就确定下来。

附图说明

图1为本发明中稀疏编码方法的流程图；

图2为本发明实施例一中获得的石炭系顶层的波形分类结果图。

具体实施方式

为使本发明要解决的技术问题、技术方案和优点更加清楚，下面将结合附图及具体实施例进行详细描述。

一种叠后地震信号波形分类方法，包括如下步骤：

(1)提取待分析的地震数据，并进行预处理，具体步骤包括：

(1-1)对算法中输入的宽方位角地震数据进行降噪，降噪方法不限定；

(1-2)提取待分析层位段地震数据，每个层位段地震数据的样本点数t，组成训练样本：

X＝{x₁,...x_n}，

其中，x_i为维度t的向量。

(2)使用流程如图1所示的稀疏编码方法提取地震数据中的特征数据x，并对特征数据x进行降维，其中，

输入：字典D，目标信号y，稀疏度K，初始化残差r₀＝y，索引集Λ＝φ，t＝1；

输出：y在字典D中的K稀疏逼近所对应的稀疏系数；

具体步骤包括：

(2-1)找出残差与字典矩阵中列内积最大的项，其列号记为λ，即：

(2-2)更新索引集Λ＝Λ∪{λ}，并根据更新后的索引集构建新的重构向量集合

(2-3)对向量Φ_t和x_t进行处理，通过最小二乘得到：

x_t＝argmin||y-Φ_tx_t||₂ 式(2-2)；

(2-4)更新残差r_t＝y-Φ_tx_t，t＝t+1；

(2-5)重复步骤(2-1)-(2-4)，直到满足停止条件：t＞k。

步骤(2)用文字描述为：利用KSVD算法对抽取到的信号进行字典学习，将学习到的字典可以作为原输入的稀疏表达。

利用正交匹配追踪算法以及上一步获得的字典对所有数据进行处理，使得每个叠后地震信号有一个稀疏表达，并将获得的稀疏系数作为叠后地震信号的特征编码。将叠后信号的每一道进行稀疏编码，利用步骤(2-5)中条件停止时输出的x_t作为叠后地震信号的特征。

(3)将步骤(2)中得到的特征数据x利用聚类算法进行分类，包括下述具体步骤：

(3-1)初始化

(3-1-1)给定用于确定截断距离d_c的参数t∈(0,1)；

(3-1-2)计算任意两样本数据间的距离

(3-1-3)将升序排列，计算截断距离

(3-1-4)计算样本数据的局部密度使局部密度降序排列的下标x_i排列生成即

(3-1-5)计算其中，表示所有局部密度比x_i大的数据点x_j中与x_i距离最近的数据点的编号，

(3-2)确定聚类中心并初始化数据点归类属性标记

(3-3)按照ρ从大到小遍历数据对非聚类中心数据点归类；

(3-4)若簇的个数N_c>1，则将每个簇中数据点进一步分为core和halo，core的局部密度大，对应簇核心部分，halo的局部密度小，对应簇边缘部分，设halo部分为若h_i＝0，则x_i属于core部分，反之x_i属于halo部分，具体步骤包括：

(3-4-1)初始化h_i＝0；

(3-4-2)为每一个簇生成一个平均局部密度上界

(3-4-3)标记簇halo部分。

实施例1：按照研制的叠后地震波形分类算法，对车排子66井区石炭系地层顶界面叠后地震波形进行处理。地震资料分析表明实验区石炭系地层顶部层段的地震波主频约为30Hz，则本次波形分类的时窗长度确定为35ms(约为地震波记录的一个主周期)，时窗的构成为石炭系顶界上10ms，下25ms。

图2为由本发明的波形分类算法获得石炭系顶层的波形分类结果图，其中分类数最终确定为10类，不同类数对应的地震道波形如图2中左上角所示

根据钻井和实际地震数据获得的不同岩性发育区的波形特征，以及钻井资料、录井、地质资料综合分析，图中不同的颜色对应的地质解释为：

(1)图中的深黑色对应的波形特征为弱振幅、杂乱反射，与钻井、地震实际剖面分析的泥岩发育区和断层发育相吻合。与测井及地质规律对比分析表明，具有线性结构的深黑色揭示了该区石炭系顶部断裂带的平面分布，而相对面积较大的区域与泥岩发育区有较好的吻合。

(2)图中除了深黑色以外的其他颜色指示火成岩的发育区，颜色的差异与火成岩的岩性、以及上覆地层的结构(尖灭)密切相关，需要综合多种因素联合分析。

(3)实验区块石炭系顶部火山岩的发育规律可以大致分为两个区域，以近南北向大断裂为界，东边区域火成岩的岩性相横向展布相对较复杂，且具有南北向的条带分布特征。

以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明所述原理的前提下，还可以作出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (5)

1.一种叠后地震信号波形分类方法，其特征在于，包括如下步骤：

(1)提取待分析的地震数据，并进行预处理；

(2)使用稀疏编码方法提取地震数据中的特征数据x，并对特征数据x进行降维；

(3)将步骤(2)中得到的特征数据x利用聚类算法进行分类。

2.根据权利要求1所述的叠后地震信号波形分类方法，其特征在于，所述步骤(1)的具体步骤包括：

(1-1)对算法中输入的宽方位角地震数据进行降噪；

(1-2)提取待分析层位段地震数据，每个层位段地震数据的样本点数t，组成训练样本：

X＝{x₁,...x_n}，

其中，x_i为维度t的向量。

3.根据权利要求1所述的叠后地震信号波形分类方法，其特征在于，所述步骤(2)中：

输入：字典D，目标信号y，稀疏度K，初始化残差r₀＝y，索引集Λ＝φ，t＝1；

输出：y在字典D中的K稀疏逼近所对应的稀疏系数；

具体步骤包括：

(2-1)找出残差与字典矩阵中列内积最大的项，其列号记为λ，即：

(2-2)更新索引集Λ＝Λ∪{λ}，并根据更新后的索引集构建新的重构向量集合

(2-3)对向量Φ_t和x_t进行处理，通过最小二乘得到：

x_t＝argmin||y-Φ_tx_t||₂ 式(2-2)；

(2-4)更新残差r_t＝y-Φ_tx_t，t＝t+1；

(2-5)重复步骤(2-1)-(2-4)，直到满足停止条件：t＞k。

4.根据权利要求1所述的叠后地震信号波形分类方法，其特征在于，所述步骤(3)包括下述具体步骤：

(3-1)初始化

(3-1-1)给定用于确定截断距离d_c的参数t∈(0,1)；

(3-1-2)计算任意两样本数据间的距离

(3-1-3)将升序排列，计算截断距离

(3-1-4)计算样本数据的局部密度使局部密度降序排列的下标x_i排列生成即

(3-1-5)计算其中，表示所有局部密度比x_i大的数据点x_j中与x_i距离最近的数据点的编号，

(3-2)确定聚类中心并初始化数据点归类属性标记

(3-3)按照ρ从大到小遍历数据对非聚类中心数据点归类；

(3-4)若簇的个数N_c>1，则将每个簇中数据点进一步分为core和halo，core的局部密度大，对应簇核心部分，halo的局部密度小，对应簇边缘部分，设halo部分为若h_i＝0，则x_i属于core部分，反之x_i属于halo部分。

5.根据权利要求4所述的叠后地震信号波形分类方法，其特征在于，所述步骤(3-4)的步骤包括：

(3-4-1)初始化h_i＝0；(3-4-2)为每一个簇生成一个平均局部密度上界

(3-4-3)标记簇halo部分。