基于分形理论的水合物回波特征提取及分类识别

徐海良，蒋昌健，杨放琼，孙思聪

(1.中南大学机电工程学院，湖南长沙，410083；2.中南大学高性能复杂制造国家重点实验室，湖南长沙，410083)

天然气水合物(natural gas hydrate)是目前最具开采前景的替代能源[1-3]。目前普遍认为储存在天然气水合物中的碳质量约是目前已探明的所有化石燃料中碳质量总和的2倍[4-5]，其中海底蕴藏着占目前全球探明总量90%以上的天然气水合物，海底将成为未来天然气水合物开采的主要场地[6]。前期勘探结果表明[7-8]，海底天然气水合物主要有以下4种赋存形式：1)伴随少量沉积物的大块天然气水合物(我国采集的天然气水合物实物样品即为此种形式)；2)以球粒状分散于细粒沉积物中的天然气水合物；3)分布在较大沙石孔隙间的天然气水合物；4)填充在岩层裂缝中的天然气水合物。而在海底浅层沉积物中，以第一种和第二种水合物为主要赋存形式。本文主要针对这种形态的天然气水合物识别问题进行研究，通过对海底天然气水合物成藏区域进行实时探测和底质识别，以区分天然气水合物和其他海底沉积物，提高天然气水合物开采效率。国内外针对天然气水合物探测、识别方法的研究主要集中在地球物理方面，特别是对如何采用地震法识别天然气水合物进行了大量研究。在勘探方面，地震资料似海底反射层(bottom simulating reflector,BSR)是用来证明海底存在水合物的最常见的证据[10]。然而，针对断裂和泥火山等构造有关的水合物矿藏，其BSR 特征往往不明显或不存在。沙志彬等[11]认为在难以利用常规地震剖面判别水合物的异常特征时，利用AVO(amplitude versus offset)反演、波阻抗反演、地震瞬时属性和能量半衰时剖面，能较好地解释水合物的地球物理异常特征，为识别水合物提供了有力的证据。吴其林等[12]发现利用宽频地震数据可以很好地刻画水合物薄层，同时也能较好地分辨厚层，更容易分析水合物矿层纵向腹肌规律变化。然而，以上地震勘探法均无法对海底底质进行实时识别，且BSR 的存在并不一定标志着天然气水合物存在。在测井法识别天然气水合物方面，孟大江等[13]利用测井资料较高的纵向分辨率和地震资料较高的横向分辨率，通过井震联合反演获得了天然气水合物准确的空间展布形态；靳佳澎等[14]通过对测井数据分析识别了水合物层,再通过合成地震记录将测井上识别的水合物层与地震反射特征及属性异常进行对比研究，利用高分辨率三维地震资料描述了水合物层分布范围，并分析了水合物层形成的主要控制因素。获取以上测井数据较复杂，需要对岩层进行钻井操作，利用高分辨率三维地震资料不适宜用于海底浅层天然气水合物的实时探测与识别。近年来，海洋可控源电磁法(controlled source electromagnetic method,CSEM)作为探测海底天然气水合物海洋油气藏的新方法得到了快速发展，探测时，拖曳的水平电偶极子作为发射源，位于海底的接收站接收包含海底以下水合物储层信息的电磁响应，由此可实现对水合物储层和储级范围的识别[15]。然而，将拖曳的电偶极子作为发射源不利于同时进行开采作业，无法使探测、识别与开采同时进行。为了实现天然气水合物的高效采集，天然气水合物识别需在海底采矿车上实时完成，上述方法多为后处理技术，因此，不适用于本文研究。而声波作为一种弹性波，具有指向性和传播性好以及传播速度稳定等特点，可应用于海底底质识别的研究。孟金生等[16]指出海底沉积物结构和物理性质的不均匀性导致声学回波信号幅值以及相位发生变化，可从回波形状中提取有效的分类特征。王正垠等[17]采用20～140 kHz的线性调频脉冲信号提取回波子带能量特征，对5类湖底沉积物进行了识别研究；KIM等[18]将相似性指数(SI)作为回波信号的特征量，分别对黏土、泥、细砂和粗砂这4种沉积物进行了分类识别，均取得了很好的识别效果。若要从声波信号中准确提取出底质信息，则还需提取回波特征。针对回波特征提取的多种处理方法存在着诸多差异性。任新敏等[19]通过相空间重构的方法，从海底回波信号中提取了3个特征量并有效地将2类海底底质分开，其结果与实际海底底质样品分析结果较吻合。PACE 等[20]提出了利用宽带回波信号识别海底介质类型的谱特征和倒谱特征法；王正垠等[17]从宽带回波信号中提取了频谱特征及宽带回波包络特征，实现了对3 种类型湖底介质的分类识别，平均准确率接近75%；栾经德等[21]研究了水下目标宽带信号散射回波频域特性，认为宽带信号可用于水下目标特征的提取和目标识别。在以上方法中，相空间重构法需提取相空间中的相似序列重复度、熵以及相平面一阶矩，后期处理过程较繁琐，效率较低；而通过频域提取特征的方法会损失部分时域信息，不利于提高识别的准确率。分形理论则广泛应用于信号处理技术中，一个信号的非线性特征和不规则程度可用分形维数来度量。回波信号的波形图像包含丰富的海底底质物理特性，可通过分形理论获得不同海底底质的声学回波特征，将天然气水合物与其他底质成分区分开来，以达到识别的目的。郝秋实等[22]利用分形布朗运动的小波系数方差估算了噪声信号的分形维数，认为分形维数可以作为一种固有特征对信号进行描述；周木荣[23]通过提取海底钴结壳调查区域内23 种底质回波信号的分形特征，实现了深海钴结壳的有效分类和识别，正确率可达85%。综上可知，相比于相空间重构法和频域特征提取法，采用提取分形维数的方法可直接提供波形结构特征，避免了因信号转换而产生的信息损失，其提取方式更直接，准确性和效率更高，可实时、便捷地对回波信号进行处理，有利于集矿机采矿作业时对海底矿区天然气水合物进行实时探测与识别，以达到提高采矿效率的目的。

1 分形理论

分形维数是表征一个分形标度不变性的基本不变量，常用的分形维数计算方法主要分为计盒子法和小波变换法。而在实际应用中，只通过单一分形维数来刻画一个复杂分形信号的分形特征往往不够全面，而需要采用多重分形来描述对象的分形特征。

设在n维欧氏空间Rn中，有1 个n维集合F，用μ表示集合F在Rn中的测度。将分形对象(F，μ)划分成N个尺度为δ的单元，设μi为每个单元的测度，其与单元尺度δ之间的关系可用标度指数αi来表征，即

具有相同标度指数的单元构成子集合Fα，Fα的s维Hausdorff测度可以表示为

Hs(Fα，μ(α))满足以下关系：若s＞f(α)，则Hs(Fα，μ(α))=0；若s=f(α)，则Hs(Fα，μ(α))为有限值；若s＜f(α)，则Hs(Fα，μ(α))=∞；其中f(α)是子集合Fα的Hausdorff维数。

若测度μ(α)在区间[α，α+Δα]内有N(α)(即覆盖盒子数)个单元，则

因此，只有当满足N(α)∝δ-f(α)时，Fα的s维Hausdorff 测度Hs(Fα，μ(α))才能取到有限值，从而，Fα的Hausdorff 维度f(α)可以通过下式计算得到：

f(α)是对分形集进行定量描述的重要指标，用于描述多重分形局部维数的连续谱。若研究对象具备单分形特征，测度为均匀分布，则f(α)是α的δ函数；若研究对象具备多分形特征，测度呈非均匀分布，则在f(α)是α的单峰凸函数[24]。

f(α)在具体应用中可能无法直接通过计算求得。而用多重分形方法描述信号局部特征的另一个理论是基于信息论角度提出的广义维数D(q)。设第i个单元的平均概率为pi，设定1个参数q，当q＞＞1，且pi＞pj时，piq＞＞pjq。由于pi对测度μ和维数D(q)的贡献度不同，因此，较大的概率被赋予较高的权值。定义概率测度的q阶矩如下：

定义广义r维测度如下：

假设存在一临界指数τ(q)，其大小依赖于pi的q阶矩，τ(q)可以表示为

其中：临界指数τ(q)为质量指数。广义Reny 维数D(q)的定义如下：

在计算时域信号广义维数的过程中，δ一般取为信号的采样间隔。通过改变q，按式(10)可计算并绘制出不同q对应的D(q)-q曲线，用最小二乘法计算曲线中无标度区的斜率，斜率的绝对值即为该q值对应的广义维度D(q)，将不同的D(q)组合最终得到1个分形信号的广义维数谱。

2 回波信号采集

为尽可能地模拟海底天然气水合物的物理环境，本文在实验研究前，制作天然气水合物和不同类型沉积物实物的模拟样品，通过采集一系列模拟样品的水下超声回波信号，对海底天然气水合物进行分类识别，以验证本文所采用的研究方法对水下天然气水合物识别的有效性与可行性。

2.1 实验系统

中南大学机电工程学院深海技术研究与开发实验室搭建了1 个超声波试验平台(实验平台示意图如图1所示)，主体是1 个长×宽×高为2.0 m×1.6 m×1.8 m的矩形水池。水池上边沿平行安装有2根平行的轨道，上面的行走平台可以前后移动以改变探测位置。平台中央安装有超声波换能器的固定装置，可用于调节换能器距水池底部的距离。本文实验中，调节超声波换能器与模拟样品的垂直距离为1 m。

图1 实验装置结构示意图Fig.1 Schematic diagram of experimental device structure

实验中所用的水下超声回波信号采集系统主要由超声波换能器、发射控制系统、回波接收系统、主控系统、工控计算机等部分组成，如图2所示。

图2 超声回波信号采集系统结构示意图Fig.2 Schematic diagram of structure of ultrasonic echo signal acquisition system

选用谐振频率为150 kHz的水下收发一体式超声波换能器，当换能器与水池底部距离1 m时，超声波在底面的辐射直径约为0.05 m，只有将实验模拟样品放置在此区域中间，换能器才能接收到有效的回波信号。发射控制系统由信号发生器(RIGOL DG1022Z)、功率放大器(西安Agitek ATA-2021)和阻抗匹配电路组成。回波接收系统由前置放大电路、时间增益控制(TGC)电路、带通滤波电路、高速模拟信号采集卡(PCI-1714U)和工控计算机组成。主控系统的核心是AT89C52 单片机，控制超声波信号的产生、接收增益的调整以及信号采集卡的同步采集。

2.2 实验样品及实验方法

天然气水合物的物理性质见表1[25]。天然气水合物结构中85%是水分子，具有多孔隙的性质，其物理性能和冰的物理性能接近。根据多孔介质冻融过程的相变规律[26]，将冰在4 ℃的环境中放置3 h，使其孔隙结构被液态水饱和，此时冰为双相多孔介质，其结构与实际天然气水合物结构类似。因此，在实验中，天然气水合物的样品可以用冰代替，实际测得样品密度为0.89 g/cm3，与天然纯水合物的密度大致相等。

表1 冰与天然气水合物的物理性质[25]Table 1 Physical properties of ice and natural gas hydrate[25]

海底沉积物一般按照粒度划分为泥、砂和砾3种主要组分，根据沉积物中各组分所占质量分数的不同，天然气水合物成藏区域的沉积物可以使用Folk分类图解进行描述，如图3所示。结合已知成藏区域内沉积物的情况，在实验中需要建立的9类模拟样本依次为泥、砂、砾、泥质砂、泥质砾、砂质泥、砂质砾、含砾泥和含砾砂。这9类模拟沉积物样本都是通过泥、砂、砾这3种基础沉积物按一定质量比例混合而成。根据国际通用粒度划分标准，要求泥的粒径小于0.062 5 mm。实验中，选用硅酸盐水泥代替沉积物中的泥；砂的粒径要求为0.062 5～2.000 mm，实验中选用粒径为0.08～0.10 mm 的细砂和1.00～2.00 mm 的粗砂混合而成；要求砾的粒径大于2.00 mm，实验中选用不同粒径(3.00～12.00 mm)的砾岩混合而成。根据这9类模拟样本中各组分的质量比，确定各样本中泥、砂、砾所需的质量如表2所示。

图3 Folk沉积物粒度三角图解法Fig.3 Triangular diagram method of Folk sediment grain size

表2 模拟沉积物样本的各组分质量Table 2 Component content of simulated sediment samples kg

制作方法如下：将各组分充分搅拌混合均匀后，浇注在统一的模具中，同时对样本进行适当夯实处理，并人为地对样本表面进行处理，使其粗糙程度不同；之后，将样本在室内环境下养护1周，直到样本胶结硬化。

回波信号采集实验的具体过程如下。

1)将模拟样本放在水池底部的中间位置，调整行走机构的位置以及超声波换能器的高度，使超声波能够垂直入射，待水池水面平静后开始实验。

2)主控系统启动信号发生器产生150 kHz的方波信号，经过功率放大器后，驱动换能器发出超声波，同时启动采集卡以8 MHz 的采样频率开始采样。

3)调整行走机构的位置，使入射声波落在每种模拟样本表面的不同位置，重复步骤2)，采集各模拟样本的超声回波信号。

4)对每种模拟样本回波信号分别截取2 000个采样点，以获得模拟样本回波信号的信号特征。各模拟样本采样信号的波形结构如图4所示。

图4 回波信号波形图Fig.4 Echo signal waveforms

3 回波信号特征提取及分类识别结果

3.1 广义维数特征提取

根据式(1)～(9)，采用MATLAB 编写广义维数D(q)的计算程序，对上述10 个样本的回波波形提取广义维数谱。在计算广义维数的过程中，阶数q的取值范围和间隔对计算精度和计算量都会产生影响。只有q的取值范围足够大才能充分反映分形信号的细节特征，但随着q增加，计算量也会成倍增加，且当q超过一定值时，继续增加q对计算结果不会产生显著影响。另外，q的取值间隔也不能太大，否则会造成对应D(q)的个数较少，这样便失去了用多重分形的方法提取信号分形特征的意义。本文中，为了保证计算的精确度，q的取值范围选为[-40，40]，间隔为2。

冰与单一介质沉积物以及混合介质沉积物回波信号的广义维数谱如图5所示。

图5 冰和沉积物样品回波信号的广义维数谱Fig.5 Generalized dimensional spectrum of echo signals from natural gas hydrate and sediment samples

由图5可知冰的广义维数谱与其他9类沉积物的广义维数谱有较大差异，这表明通过提取回波信号的广义维数特征可达到区分冰(天然气水合物模拟物)和其他沉积物的目的。

3.2 神经网络分类识别

概率神经网络是一种基于核密度函数估计法和贝叶斯分类准则的神经网络模型，属于径向基函数神经网络的一个分支。它以核密度函数作为隐含层激活函数，具有传统前馈神经网络所没有的优点，其在故障分析、模式识别和特征分类领域已得到广泛应用。概率神经网络包含有4层基本网络，分别是输入层、隐含层、求和层和输出层，其结构模型如图6所示。

图6 概率神经网络结构图Fig.6 Structure diagram of probabilistic neural network

1)输入层的作用是将样本接收后输入到样本层的各个节点，其节点个数与特征矢量x的维数相同。

2)对隐含层将输入样本xi乘以加权系数ωi后求和。每个节点对应1个中心，其个数与输入样本个数相同。

3)求和层由m个节点组成，将同一类样本经隐含层输出之后求和，从而得到样本概率密度的核密度函数估计结果。

4)输出层由m个竞争神经元组成。比较所有输出层神经元后验概率密度，概率密度最大的那一类输出1，其余输出0。

隐含层中高斯核函数宽度系数σ的取值会影响概率神经网络的性能，这也是神经网络设计过程中唯一需要确定的自由参数，它表征不同模式概率密度函数之间的重叠程度。σ要足够大才能包含所有输入样本所在的区间，但σ过大时会影响神经网络的计算精度，因此，需选用不同的参数σ，分别进行神经网络识别实验。

概率神经网络进行模式分类的过程主要分为训练和测试两步。训练的具体流程如下：

1)确定训练样本的总类别数m，获取各种类的特征向量样本集，经归一化处理后输入到神经网络。

2)根据样本参数确定神经网络的总体结构。

3)确定不同的高斯核函数宽度系数σ。

4)选用一定数量的样本数据构成训练样本集，只需要将样本集输入训练1次，便可建立最终的概率神经网络模型。

概率神经网络训练完成之后，便可以将待识别的样本输入到网络中进行识别测试，实现贝叶斯最优化分类。

在回波信号采集的实验中，共提取包括天然气水合物模拟物(冰)在内的10类模拟样本的回波数据，每一类样本各采集100组数据。对以上数据提取广义维数特征，可以得到41 维广义维数特征向量集，从向量集中随机抽取60 组作为概率神经网络的训练样本，另外40 组作为测试样本。下面分别对未经降维优化和经过核Fisher 判别分析法(kernel Fisher discriminant analysis,KFDA)降维之后的分类识别结果进行比较分析。

3.2.1 基于广义维数特征的PNN分类识别

直接将实验中提取的回波样本的41 维广义维数谱D(q)作为概率神经网络的样本集，将10 类试样共600个训练样本输到网络中进行训练，确定各个节点的权值系数。将剩下的400个测试样本作为输入，训练完成后神经网络的识别结果如表3所示。

表3 基于广义维数特征的分类识别Table 3 Classification and recognition based on generalized dimension features

3.2.2 基于KFDA特征优化的PNN分类识别

从向量集中随机抽取60 组作为概率神经网络(PNN)的训练样本，另外40组作为测试样本。利用核Fisher判别分析法对训练样本进行降维优化，得到9 维最优特征向量集和1 个最优判别矢量w。选择不同宽度系数σ，将60组训练样本作为输入数据对神经网络进行训练，建立识别分类模型。

通过对40 组测试样本的判别可以评价识别分类模型的性能。首先将40 组测试数据按照最优判别矢量w投影，得到9 维测试样本最优特征向量集。然后，将该向量集作为神经网络的输入，训练神经网络完成后输出分类结果，据此可判断该测试样本属性。识别结果如表4所示。

表4 基于KFDA特征优化的分类识别Table 4 Classification and recognition based on KFDA feature optimization

从表3和表4可知：从原始回波数据中提取的广义维数特征的维数偏高，使得神经网络模型的结构过于庞大(输入层节点单元数为41)，含有大量对分类识别作用不大的特征，造成最终水合物的平均识别准确率为68.75%，非水合物的平均识别准确率为67.93%，识别效率和准确率较低。但经过KFDA 特征优化处理过后，神经网络模型网的结构得到大大简化(输入层节点单元数为3)，计算时间缩短到未经降维的1/4，水合物的平均识别准确率达83.13%，非水合物的平均识别准确率达81.45%。

4 结论

1)利用分形维数可以有效提取出回波信号特征。引入的多重分形理论将冰及其余沉积物的复杂回波信号进行降维优化，提高了特征的分类精度，可以在保证特征分类精度的前提下，提取出对分类有效的信息，剔除冗余的特征维度。另外，从广义维数谱看，冰的广义维数谱相比于其他类沉积物有较大区别，说明本文所采用的方法对天然气水合物及沉积物的分类识别具有可行性。

2)采用概率神经网络分别对冰与9类沉积物试样回波信号的原始广义维数谱以及经过核Fisher判别分析法优化之后的特征进行分类识别实验，得到对未经优化的原始广义维数谱的平均识别准确率为68.8%，降维优化处理过后的平均识别准确率为82.3%，证明了核Fisher 判别分析法适用于水下沉积物回波信号的特征降维优化，能有效提高分类识别的准确率。

3)不同宽度系数对应的冰和沉积物的识别准确率不同。实际上，两者的识别准确率往往不能同时达到最高。在后续研究中，建议通过设计优化相关算法，改善实验条件，增加样本数量，以提高对天然气水合物(模拟样品)实时探测并提高识别的准确率。