CN106896407A - 一种基于近似负熵的微地震信号初至拾取方法 - Google Patents

Abstract

本发明提供一种基于近似负熵的微地震信号初至拾取方法，属于地球物理技术领域。根据微地震信号自身特点设计合适的近似负熵函数，结合负熵谱曲线极值差的加权函数设定阈值，进而实现信号区域与噪声区域的准确划分及微地震信号初至的准确拾取。本发明的有益效果是根据微地震信号自身的特点，引用信息论中近似负熵理论提出了一种基于近似负熵的微地震初至拾取新方法，与现有拾取方法相比，新方法能够实现较低信噪比条件下初至到时的准确拾取，信噪比可低至‑12dB，并且计算过程简单，速度较快，准确性高。

Description

一种基于近似负熵的微地震信号初至拾取方法

技术领域

本发明属于地球物理技术领域，尤其是指微地震信号处理方法。

背景技术

页岩气、煤层气等非常规油气勘探是未来油气产业的重点、是全球都在关注的课题。油气勘探越来越多的使用水力压裂及微地震检测技术。在微震信号识别或震源定位工作中，由于微地震信号具有数据量大、频率分布广、信噪比低等特点，微地震初至到时的检测和拾取是一个重要而困难的问题，很难找到一种精确而快速的微地震初至拾取方法。

近年来，国内外学者先后提出多种基于单分量或三分量的时域或频域初至拾取算法，其中具有代表性的方法有Akaike信息准则法、长短时能量比法等。Akaike信息准则法方法是常与自回归技术一起使用，这种方法计算比较简单。Sleeman(1999)提出了寻找两段平稳时间序列最优分割点的方法，可以直接计算数据AIC并将结果最小值作为P波初至点。2011年，St-Onge等人将该方法应用于微地震初至拾取当中，在弱噪声条件下拾取效果好，但在强噪声条件下无法准确识别微地震信号。长短时能量比法是一种时域的震相识别算法，它利用微地震信号振幅、能量、包络等特性，按照短时帧和长时帧划分特征函数进而计算新时间序列的平均值之比，比值的最大值对应地震记录初至点。后来很多学者陆续对其进行改进，2009年，Wong等人提出了修改能量比法。长短时能量比法受时窗长度影响较大，低信噪比条件下合适窗长的选择是该类方法的一个难点。近几年来，多种初至拾取算法不断被提出，如移动时窗标准差法(Akram,2011)，功率谱包络STA/LTA法(Xiantai,et al,2011)，动态阈值法(Akramand Eaton,2012)，极化坐标方法(Moriya,2008,2009)、高阶统计量(HOS)和AIC(Küperkoch et al.,2010)，分形算法(Liao et al.,2010)、互相关法(Eaton，2016)、聚类法(Zhu et al.，2016)、小波变换和AIC结合方法(Zhang et al.,2003)、多尺度小波分析和STA/LTA联合法(Capilla,2006，Rodriguez，2011，Li et al，2016)等。以上算法都能够实现初至拾取，但在强噪声、弱有效信号条件下，分别存在计算量大，受窗长影响较大，准确度降低等问题。因此，寻找一种计算速度快、拾取精度高的微地震初至拾取算法是目前该研究领域专家们不断探讨和研究的问题。

发明内容

本发明提供一种基于近似负熵的微地震信号初至拾取方法，以解决针对较低信噪比条件下微地震信号初至时刻点拾取困难的问题。该方法根据微地震信号自身特点设计合适的近似负熵函数，结合负熵谱曲线极值差的加权函数设定阈值，进而实现信号区域与噪声区域的准确划分及微地震信号初至的准确拾取。

本发明采取的技术方案是，包括以下步骤：

1)根据微地震信号自身的特性，设计近似负熵函数：

具有n个采样点的随机序列y＝{y_i,i＝1,2,…,n}的负熵可以表示为：

J(y)＝H(y′)-H(y)， (1)

其中y′是与随机序列y具有相同均值和方差的高斯随机序列，H(y)是序列y的熵，同理H(y′)是序列y′的熵，J(y)就表示随机变量y的负熵；

构造基于最大熵准则(最小负熵)的近似负熵函数，其对信号相对敏感，表达式描述如下

其中k_i是大于零的常数，G_i(·)是对比方程，是一个正交的非二次函数，x是公式(1)中随机变量y经过归一化后得到的零均值、单位方差序列，x′是满足零均值，单位方差的标准正态分布序列；E[·]表示数学期望，i表示选取的对比方程的个数，i值越大，对负熵的近似程度越接近，通常情况下n＝2，系数k₁＝k₂＝1；

对于任意随机序列ζ＝{ζ_i,i＝1,2,…,n}，选取的非二次对比方程G_i(·)为公式(3)和(4)中的两种形式：

两个对比方程G₁(ζ)和G₂(ζ)均为分段线性计算对比函数，G₁(ζ)具有良好的稳健性，可以用于检测独立分量的超高斯性，G₂(ζ)适用于估计随机变量的亚高斯性，由于序列x′是已知的标准正态序列，通过公式计算可得于是公式(2)中的近似负熵表达式就变为

式中NJ(x)代表随机序列x的近似负熵值，

2)对地震信号进行归一化和分帧处理，并计算各帧内的负熵值，绘制负熵谱曲线；

根据公式(2)的要求，在计算地震数据近似负熵之前，需要对微地震序列y＝{y_i,i＝1,2,…,n}进行归一化处理，得到零均值、单位方差的随机序列，归一化公式如下

式中y＝{y_i,i＝1,2,…,n}对应n个微地震数据点，和σ²分别代表数据序列y的均值和方差，序列是由n个构成的与序列y同维的数据序列，x为归一化后的零均值、单位方差微地震数据序列；

然后，选取帧长F_s和帧移F_l对归一化后的地震数据进行分帧处理，再根据式(5)中的近似负熵函数分别计算出各个帧内微地震数据的近似负熵值，并根据各帧的计算结果绘制负熵谱曲线；

3)设置阈值函数找到噪声与信号的分界点，得到微地震数据初至时刻点的具体位置；

将近似负熵谱曲线的最大值和最小值的差值进行加权，再通过一个加权因子来调节阈值的敏感度和准确度，具体表达式如下

thr＝α(max(NJ)-min(NJ))， (9)

其中thr为计算得到的阈值，α为加权因子，能够调节插值的高度，根据负熵谱和阈值线的交点检测到初至所在的帧，再根据公式(10)便可得到初至对应的具体数据点A_t；

其中F是初至所在的帧数，F_s和F_l分别代表帧长和帧移，表示对向上取整的结果，从而确定初至到时的所在位置。

本发明的有益效果是根据微地震信号自身的特点，引用信息论中近似负熵理论提出了一种基于近似负熵的微地震初至拾取新方法，与现有拾取方法相比，新方法能够实现较低信噪比条件下初至到时的准确拾取，信噪比可低至-12dB，并且计算过程简单，速度较快，准确性高；该方法根据微地震信号自身特点设计合适的近似负熵函数，结合负熵谱曲线极值差的加权函数设定阈值，进而实现信号区域与噪声区域的准确划分及微地震信号初至的准确拾取。

附图说明

图1(a)单道实际微地震数据及检测结果，含噪的微地震信号信噪比非常低，有效信号几乎无法识别，图中虚线标出位置为有效信号的初至时刻点；

图1(b)近似负熵谱曲线和阈值函数，从负熵谱曲线可以看出噪声段与信号段的负熵值存在明显的差别，图中虚线标注位置为阈值函数值；

图2(a)信噪比为-3dB的含噪记录，该含噪记录信噪比不是很低，有效信号段位置可以看出，图中用三种不同线型分别标注出三种方法的拾取结果，可以看到三种方法的拾取位置大致相同；

图2(b)AIC拾取结果，在AIC曲线中存在一个明显的最小值(最低点)，可以很容易的找到初至位置；

图2(c)STA/LTA拾取结果，在STA/LTA曲线中存在一个明显的最大值(峰值)，可以很容易找到初至位置；

图2(d)近似负熵拾取结果，在负熵曲线中，信号段与噪声段存在明显差别，通过阈值可以找到初至位置，由于此时加入的是高斯白噪声，因此噪声段负熵值为零；

图3(a)信噪比为-12dB的含噪微地震记录，该含噪记录信噪比很低，有效信号段位置几乎不能看出，图中用三种不同线型分别标注出三种方法的拾取结果，可以看到三种方法的拾取位置差别很大；

图3(b)AIC方法拾取结果，在AIC曲线中存在一个最小值(最低点)，可以找到初至位置，但是这个最小值不是准确的初至位置，也就是说在低信噪比条件下，该方法的函数曲线准确性下降；

图3(c)STA/LTA方法拾取结果，在STA/LTA曲线中不存在明显的最大值(峰值)，找到的曲线最大值也不是初至时刻的准确位置，因此在低信噪比条件下，该方法的函数曲线准确性下降；

图3(d)近似负熵方法拾取结果，在负熵曲线中，信号段与噪声段仍然存在明显差别，通过阈值可以找到初至位置，由于此时加入的是高斯白噪声，因此噪声段负熵值为零；

图4 AIC、STA/LTA和近似负熵三种方法对100道不同信噪比的合成微地震信号进行拾取的绝对误差图，从图中可见随着信噪比的不断降低，AIC和STA/LTA方法拾取误差明显增加，AIC方法从-6dB开始误差明显增加，STA/LTA方法在-10dB时误差明显增加，而近似负熵方法在-12dB时误差才开始略有增加；

图5(a)实际的三分量微地震数据的X方向分量，该方向数据信噪比较低；

图5(b)实际的三分量微地震数据的Y方向分量，该方向数据信噪比较高；

图5(c)实际的三分量微地震数据的Z方向分量，该方向数据信噪比较高；

图6(a)实际的三分量微地震数据X方向分量拾取结果，由于该方向数据信噪比较低，三种方法的拾取结果存在很大差异，近似负熵方法拾取位置较另两种方法准确；

图6(b)实际的三分量微地震数据Y方向分量拾取结果，由于该方向数据信噪比较高，三种方法的拾取结果接近；

图6(c)实际的三分量微地震数据Z方向分量拾取结果，该方向数据信噪比较高，但噪声类型有所区别，STA/LTA和近似负熵方法拾取结果接近，AIC方法拾取结果很不准确。

具体实施方式

1)根据微地震信号自身的特性，设计合适的近似负熵函数；

对于具有单位方差随机变量来说，分布越“随机”，熵越大，负熵越小，具有n个采样点的随机序列y＝{y_i,i＝1,2,…,n}的负熵可以表示为：

J(y)＝H(y′)-H(y)， (1)

其中y′是与随机序列y具有相同均值和方差的高斯随机序列，H(y)是序列y的熵，同理H(y′)是序列y′的熵，J(y)就表示随机变量y的负熵；

由于微地震数据具有数据量大、频率分布广、信噪比低等特点，构造基于最大熵准则(最小负熵)的近似负熵函数，其对信号相对敏感，表达式描述如下：

其中k_i是大于零的常数，G_i(·)是对比方程，是一个正交的非二次函数，x是公式(1)中随机变量y经过归一化后得到的零均值、单位方差序列，x′是满足零均值，单位方差的标准正态分布序列。E[·]表示数学期望，i表示选取的对比方程的个数，i值越大，对负熵的近似程度越接近，通常情况下n＝2，系数k₁＝k₂＝1；

在整个式子中对比方程G_i(·)的选择非常重要，它直接影响负熵近似的精确程度，在选取对比方程G_i(·)时，函数值增长不能太快，计算应该简单、快速，并且保证确定分布是可以计算的，为了得到更加稳健的近似并且满足上述选择准则，对于任意随机序列ζ＝{ζ_i,i＝1,2,…,n}，选取的非二次对比方程G_i(·)为公式(3)和(4)中的两种形式：

两个对比方程G₁(ζ)和G₂(ζ)均为分段线性计算对比函数，G₁(ζ)具有良好的稳健性，可以用于检测独立分量的超高斯性，G₂(ζ)适用于估计随机变量的亚高斯性，由于序列x′是已知的标准正态序列，通过公式计算可得于是公式(2)中的近似负熵表达式就变为

式中NJ(x)代表随机序列x的近似负熵值，这种近似负熵的算法不仅简单，计算速度快速，还具有较好的统计特性，很适合应用于微地震信号初至拾取当中。

2)对地震信号进行归一化和分帧处理，并计算各帧内的负熵值，绘制负熵谱曲线；

本发明将近似负熵思想应用于微地震信号初至到时拾取，通过对微地震信号进行归一化、分帧处理后，在各个帧内计算微地震数据的近似负熵值，获得当前帧的统计信息，并根据各帧的计算结果绘制负熵谱曲线，设计一种合适的阈值函数，找到噪声与信号的分界点，从而确定初至到时的所在位置；

首先，根据公式(2)的要求，在计算地震数据近似负熵之前，需要利用公式(6)，(7)和(8)对微地震序列y＝{y_i,i＝1,2,…,n}进行归一化处理，得到零均值、单位方差的随机序列；

式中y＝{y_i,i＝1,2,…,n}对应n个微地震数据点，和σ²分别代表数据序列y的均值和方差，序列是由n个构成的与序列y同维的数据序列，x为归一化后的零均值、单位方差微地震数据序列；

由于微地震信号具有随机性，在计算近似负熵值之前需要对微地震数据进行分帧处理，帧长和帧移的确定直接影响拾取效果，帧长越大，负熵谱曲线越平滑，但是不能太长，以保证拾取的精度。图1(a)为一道512点的实际微地震随机数据，对其进行归一化处理并分帧，为了保证实验数据的一致性，帧移和帧长的大小分别选为38和3；接着根据式(5)中的近似负熵函数计算每帧内地震数据的近似负熵值，并画出各帧对应的近似负熵谱曲线，如图1(b)所示；从负熵谱曲线可以明显看到噪声段和有效信号段的负熵值存在明显差异，从而说明了近似负熵方法的可行性。

要想实现精确地拾取初至，就需要设置合适的阈值函数，从近似负熵谱曲线中找到初至所在的具体位置；阈值选择的方法直接影响拾取的敏感度和准确度，阈值算法有许多种，如何设计一种即简单又可靠的阈值函数也是本文方法的关键点；通过分析近似负熵谱曲线的特点，设计了一种简单的阈值函数，将近似负熵谱曲线的最大值和最小值的差值进行加权，通过一个加权因子来调节阈值的敏感度和准确度，具体表达式如下：

thr＝α(max(NJ)-min(NJ)),0<α<1， (9)

其中thr为计算得到的阈值，α为加权因子，能够调节插值的高度，实际中α应取很小值，这样保证拾取的位置有足够小的延迟，对图1(b)中的负熵谱曲线加权因子α＝0.2，相应的阈值线为图1(b)的点划线所示，根据负熵谱和阈值线的交点就可检测到初至所在的帧，再经过式(10)便可得到初至对应的具体数据点A_t：

其中F是初至所在的帧数，F_s和F_l分别代表帧长和帧移，表示对向上取整的结果，从而确定初至到时的所在位置。对图1(b)中检测到的交点标出的帧数利用式(10)进行转换，得到图1(a)中地震数据初至点的具体位置，如图1(a)中虚线所示。

下面通过具体的实验例来说明进一步说明本发明方法的效果。

实验例1合成记录

模拟生成一张采样频率为1000Hz的100道合成微地震记录，每道512个采样点，有效信号主频为300Hz，对该纯净地震记录加入不同强度的高斯白噪声，得到信噪比从-1dB到-12dB的合成含噪微地震记录。为了验证信噪比对初至拾取方法性能的影响，首先从信噪比为-3dB的含噪记录中选取一道初至到时为160ms单道含噪微地震数据，如图2(a)所示。分别利用AIC方法、STA/LTA方法和近似负熵法对该道含噪地震数据进行初至点拾取，拾取结果分别如图2(b)、图2(c)和图2(d)所示。

图2(b)是AIC方法的拾取结果，其最小值处代表在最小平方意义下噪声信号和微地震信号的拟合度最差，该点即对应了图2(a)中含噪数据的初至点位置。图2(c)是STA/LTA方法的拾取结果，其最大值代表短时窗和长时窗的比值最大，即信号此时的前后段差异最大，所以STA/LTA的最大值对应于待检测含噪数据的初至点位置。图2(d)为近似负熵方法得出的负熵谱曲线，曲线和阈值函数的交点即为图2(a)中含噪数据的初至点所在的帧，根据公式(10)中帧与数据点的对应关系即可得出初至点位置。最后在图2(a)中分别用不同的线型标出了三种方法拾取到的初至点位置。从图2中三种方法的拾取结果可可看出，AIC方法的拾取结果为162ms，STA/LTA方法的拾取结果为157ms，近似负熵方法的拾取结果为159ms，与含噪记录的实际初至位置160ms相比，误差均小于0.5％，因此，对于信噪比为-3dB的含噪地震数据，三种方法均可实现初至点的准确拾取。

接着，从信噪比为-12dB的含噪地震记录中选取一道初至到时仍为160ms的单道含噪地震数据，如图3(a)所示。分别利用AIC方法、STA/LTA方法和近似负熵方法对图3(a)中的含噪地震数据进行初至点拾取。图3(b)是AIC方法的拾取结果，图3(c)是STA/LTA方法的拾取结果，图3(d)为近似负熵方法的拾取结果。从图3(a)中标注的三种方法拾取的初至点位置可知，AIC方法的拾取结果为22ms，STA/LTA方法的拾取结果为212ms，近似负熵方法的拾取结果为159ms，与待检测含噪记录的实际初至位置160ms相比，AIC和STA/LTA方法的拾取误差分别为13.5％和5.1％，而近似负熵算法拾取结果与实际的初至点160ms误差为0.1％。

表1给出了利用AIC、STA/LTA和近似负熵三种方法对信噪比从-1dB到-12dB的含噪记录进行初至点拾取的统计数据。

表1 AIC、STA/LTA和AN三种方法在不同信噪比的条件下的拾取结果

从表中数据的统计结果可知，当信噪比从-1dB到-7dB时，三种方法的平均拾取误差分别为：0.5％,0.3％和0.1％，此时认为三种方法拾取效果都比较好；而当信噪比低于-8dB时，AIC方法的拾取误差高达16.4％，且当信噪比低于-10dB时，STA/LTA方法的拾取误差也高达19.5％，但是近似负熵方法在信噪比为-1dB到-12dB情况下平均拾取误差均小于0.1％，可见新方法的准确性和可靠性。

图4给出了三种方法对100道微地震合成记录在不同信噪比条件下进行初至拾取的绝对平均误差折线。

实验例2实际记录

图5为一组实际的三分量微地震数据，该组地震记录共含15道，每道512个采样点。分别用AIC、STA/LTA和近似负熵方法对其进行初至时刻提取。从含噪数据可以看出，不同分量上的信噪比和信号幅度都不一致，对三个分量方向上的数据分别利用AIC，STA/LTA和近似负熵方法进行初至点拾取。分别从X，Y，Z分量拾取结果中提取了其中一道数据，其单道波形如图6所示。从图6中的单道波形图也可以看出当信噪比较高时，三种方法拾取到的初至点位置接近，拾取误差较小，而当信噪比较低时，AIC和STA/LTA两种方法的拾取位置明显存在较大误差，近似负熵方法具有较好的准确性和可靠性。

Claims (2)

1.一种基于近似负熵的微地震信号初至拾取方法，其特征在于，包括下列步骤：

1)根据微地震信号自身的特性，设计近似负熵函数：

具有n个采样点的随机序列y＝{y_i,i＝1,2,…,n}的负熵可以表示为：

J(y)＝H(y′)-H(y)， (1)

其中y′是与随机序列y具有相同均值和方差的高斯随机序列，H(y)是序列y的熵，同理H(y′)是序列y′的熵，J(y)就表示随机变量y的负熵；

构造基于最大熵准则(最小负熵)的近似负熵函数，其对信号相对敏感，表达式描述如下：

$J\left(y\right)=\underset{i=1}{\overset{n}{\&Sigma;}}{k}_i{\{E\&lsqb;G_i\left(x\right)\&rsqb;-E\&lsqb;G_i\left(x^{\&prime;}\right)\&rsqb;\}}^2,---\left(2\right)$

其中k_i是大于零的常数，G_i(·)是对比方程，是一个正交的非二次函数，x是公式(1)中随机变量y经过归一化后得到的零均值、单位方差序列，x′是满足零均值，单位方差的标准正态分布序列；E[·]表示数学期望，i表示选取的对比方程的个数，i值越大，对负熵的近似程度越接近，通常情况下n＝2，系数k₁＝k₂＝1；

对于任意随机序列ζ＝{ζ_i,i＝1,2,…,n}，选取的非二次对比方程G_i(·)为公式(3)和(4)中的两种形式：

$G_1\left(\mathrm{\&zeta;}\right)=-\exp \left(\frac{-{\mathrm{\&zeta;}}^2}{2}\right),---\left(3\right)$

$G_2\left(\mathrm{\&zeta;}\right)=\frac{1}{4}{\mathrm{\&zeta;}}^4,---\left(4\right)$

两个对比方程G₁(ζ)和G₂(ζ)均为分段线性计算对比函数，G₁(ζ)具有良好的稳健性，可以用于检测独立分量的超高斯性，G₂(ζ)适用于估计随机变量的亚高斯性，由于序列x′是已知的标准正态序列，通过公式计算可得于是公式(2)中的近似负熵表达式就变为

$NJ\left(x\right)={\{E\&lsqb;-\exp (-\frac{x^2}{2})\&rsqb;-\sqrt{\frac{1}{2}}\}}^2+{\{E\&lsqb;\frac{1}{4}{x}^4\&rsqb;-\frac{3}{4}\}}^2,---\left(5\right)$

式中NJ(x)代表随机序列x的近似负熵值，

2)对地震信号进行归一化和分帧处理，并计算各帧内的负熵值，绘制负熵谱曲线；

选取帧长F_s和帧移F_l对归一化后的地震数据进行分帧处理，再根据式(5)中的近似负熵函数分别计算出各个帧内微地震数据的近似负熵值，并根据各帧的计算结果绘制负熵谱曲线；

3)设置阈值函数找到噪声与信号的分界点，得到微地震数据初至时刻点的具体位置；

将近似负熵谱曲线的最大值和最小值的差值进行加权，再通过一个加权因子来调节阈值的敏感度和准确度，具体表达式如下

thr＝α(max(NJ)-min(NJ))， (9)

其中thr为计算得到的阈值，α为加权因子，能够调节插值的高度，根据负熵谱和阈值线的交点检测到初至所在的帧，再根据公式(10)便可得到初至对应的具体数据点A_t；

其中F是初至所在的帧数，F_s和F_l分别代表帧长和帧移，表示对向上取整的结果，从而确定初至到时的所在位置。

2.根据权利要求1所述的一种基于近似负熵的微地震信号初至拾取方法，其特征在于，步骤2)中对地震信号进行归一化和分帧处理如下：

根据公式(2)的要求，在计算地震数据近似负熵之前，需要对微地震序列y＝{y_i,i＝1,2,…,n}进行归一化处理，得到零均值、单位方差的随机序列，归一化公式如下

$\stackrel{\&OverBar;}{y}=\frac{1}{n}\underset{i=1}{\overset{n}{\&Sigma;}}{y}_i,---\left(6\right)$

${\mathrm{\&sigma;}}^2=\frac{1}{n-1}\underset{i=1}{\overset{n}{\&Sigma;}}{(y_i-\stackrel{\&OverBar;}{y})}^2,---\left(7\right)$

$x=\frac{y-\stackrel{\&OverBar;}{y}}{\mathrm{\&sigma;}},---\left(8\right)$

式中y＝{y_i,i＝1,2,…,n}对应n个微地震数据点，和σ²分别代表数据序列y的均值和方差，序列是由n个构成的与序列y同维的数据序列，x为归一化后的零均值、单位方差微地震数据序列。