相对传统取样分析确定海底底质方法，声学底质分类具有速度快、费用低等优势。特征提取与选择是声学底质分类过程的重要组成部分，对分类精度有着显著影响。可以提取声纳图像灰阶共生矩阵，使用反向散射强度、入射角和功率谱等特征进行海底分类。其中，分形维因具有自相似性和标度不变性，对尺度要求灵活，更具稳健性^[1]，能反映空间结构信息及其变化规律。在海底底质分类特征提取工作中，Tegowski统计了分形维数、反向散射强度和光谱宽度^[2]；Walreea等基于分形维数和光谱^[3]；阳凡林等基于分形维数和共生矩阵^[4]；马飞虎等^[5]基于分形维数、能量谱、反向散射强度的角度依赖性等特征分别进行分类工作。以上研究都只提取了分形维数进行分类，但是，文献[6]表明，在某些情况下，不同纹理具有相同的分形维数，单一分形维数并不能充分刻画图像的分形信息，从而造成分类可靠性降低。为此，本文提出了基于分形理论的多特征底质分类方法。

1 回波强度特征提取 1.1 分形维数

分形几何中，分形维数d是一个描述分形对空间填充程度的统计量^[7]。数学定义如下：设ε为小立方体的长度，N(ε)为用此小立方体覆盖被测形体所得的数目，分形维数d为：

$ d = \mathop {\lim }\limits_{\varepsilon \to 0} \left[{{\rm{log}}N(\varepsilon )/{\rm{log}}\left( {\frac{1}{\varepsilon }} \right)} \right] $

(1)

分形维定义有很多方法，比如豪斯多夫维数、计盒维数和分配维数等，本文采用差分计盒维数(differential box-counting，DBC)计算分形维数。

对于一副M×M大小的图像，将其看作三维空间。其中，坐标(x, y)表示原图像像素的坐标位置，z坐标表示灰度值，将原二维图像划分为s×s大小的网格，s可取2~M/2之间的任意整数值，即1 < s≤M/2，则r=s/M。设置s×s×h的盒子覆盖该s×s网格，其中，盒子的高度h满足G/h=M/s，G为图像的灰度等级数目。若第(i, j)个网格像素的最大和最小灰度落入覆盖该网格一列盒子的第l和第k个盒子内部，则：

$ {n_r}(i, j) = l-k + 1 $

(2)

式中，n_r(i, j)表示覆盖第(i, j)网格的盒子数。覆盖整个图像的盒子数N_r为所有n_r(i, j)之和，为：

$ {N_r} = \sum\limits_{i, j} {{n_r}(i, j)} $

(3)

取不同r计算相应的N，利用式(1) 线性拟合计算log (N_r)相对于log (1/r)的斜率，为：

$ D = \lim \frac{{\log ({N_r})}}{{\log (1/r)}} $

(4)

3种典型底质如图 1所示，DBC方法分形维数计算结果如下：沙为2.18，泥为2.33，砾石为2.87。

1.2 空隙

分形维数可以刻画图像纹理信息，但是单一的分形维数并不足够刻画该信息。因此，Mandelbrot提出了空隙概念，以更好地描述图像纹理分形信息^[6]。Kelle等(1989) 定义的基于概率分布函数p(m, L)的空隙特征为：

$ C(L) = \frac{{M(L)-N(L)}}{{M(L) + N(L)}} $

(5)

基于DBC的空隙计算方法为^[8]:

$ M\left( L \right) = \frac{P}{{{N_r}}}, N\left( L \right) = \frac{{{N_r}}}{P} $

(6)

式中，P为图像S所含的像素点数；N_r为DBC方法计算的结果，代入式(5) 得：

$ C\left( L \right) = \frac{{{P^2}-{N_r}^2}}{{{P^2} + {N_r}^2}} $

(7)

空隙特征在-1~1之间，盒子尺寸过大过小都不利于底质的分类，实际计算中要进行多次计算，确定最佳区分效果的盒子尺寸。

1.3 多重分形

多重分形测度研究物理量或其他量在几何支集上的分布，多重分形能够从全局和局部两方面充分描述图像的奇异性结构，根据差分盒计维数的多重分形维计算过程为：

$ \mathop {\lim }\limits_{r \to 0} {\mu _r}(x)\sim {r^a} $

(8)

根据DBC的计算方法，令r=s/M，因此：

$ {\mu _r}\left( {i, j} \right) = {n_r}\left( {i, j} \right)/{N_r} $

(9)

式中, n_r(i, j)和N_r为DBC方法计算的结果，令

$ \chi \left( {q, r} \right) = \sum {[{\mu _r}\left( {i, j} \right)]^q} \approx {r^{\tau (q)}} $

(10)

这个函数的作用是通过学习所有的q阶矩去估计所有分布的奇异性。则

$ (q-1)D(q) \equiv \tau (q) = \mathop {\lim }\limits_{r \to 0} \frac{{{\rm{In}}\chi (q, r)}}{{{\rm{In}}r}}, q \ne 1 $

(11)

计算过程中，q的取值可以为(-∞, +∞)中的任意数，在实际预测时，取值不可能遍及整个区域，因此，需要从中选择出对预测有效的数据。实际计算中，当q的取值越过一定范围，再增大对计算结果已无显著影响时，可以确定q的取值范围^[9]。

1.4 灰度基本统计量

基本的统计算法主要包括均值、标准差和分位数。均值为：

$ u = \frac{1}{{mn}}\sum\limits_{i = 1}^m {\sum\limits_{j = 1}^n {{x_{ij}}} } $

(12)

标准差为：

$ \sigma = \sqrt {\frac{1}{{mn}}\sum\limits_{i = 1}^m {\sum\limits_{j = 1}^n {{{\left( {{x_{ij}}-u} \right)}^2}} } } $

(13)

式中，m、n表示每个矩形块有m行n列；变量x_ij表示第i行j列的振幅值。分位数用来描述随机变量的概率分布, 本文用中位数表述反向散射强度分布。

2 特征向量主成分分析及分类

将回波强度图像分成小矩形块，对每个小矩形块计算其特征值，包括分形维数、空隙特征、多重分形维、均值、标准差、分位数，得到原始特征值矩阵为：

$ \boldsymbol{X} = \left[{\begin{array}{*{20}{c}} {{x_{11}}}& \cdots &{{x_{1p}}}\\ \vdots &{}& \vdots \\ {{x_{n1}}}& \cdots &{{x_{np}}} \end{array}} \right] $

(14)

主成分Y_i是p个原始指标向量的线性组合，并且要求Y_i的方差尽可能大，Y的协方差为零。

$ {\boldsymbol{Y}_i} = {a_{i1}}{\boldsymbol{X}_1} + {a_{i2}}{\boldsymbol{X}_2} + \cdots + {a_{ip}}{\boldsymbol{X}_p} $

(15)

$ \mathit{\boldsymbol{Y}} = \mathit{\boldsymbol{X \boldsymbol{\varPhi} }} $

(16)

式中，Φ为系数矩阵，由(a_i1, a_i2, …, a_ip)组成。

计算主成分需要对原始数据进行标准化，即

$ {x_0} = \frac{{x-{\rm{ave}}}}{{{\rm{std}}}} $

(17)

式中，ave为X的均值；std为X的方差。之后，计算X₀的协方差矩阵为：

$ \boldsymbol{R} = \frac{{\boldsymbol{X}_0^{\rm{T}}{\boldsymbol{X}_0}}}{{n-1}} $

(18)

求出协方差矩阵R的特征值λ_i和特征值向量Φ_i。根据λ_i取前k个主成分，得到k个主成分向量后，对其进行K-均值聚类，并分析聚类结果。本文采用轮廓系数作为聚类评价工具，轮廓结合了凝聚度和分离度，可以以此来判断聚类的优良性，单点的轮廓系数计算公式为：

$ {s_i} = \frac{{{b_i}-{a_i}}}{{\max ({a_i}, {b_i})}} $

(19)

式中，a_i表示第i个点到本类其他点的平均距离；b_i是第i个点到其他类点平均距离的最小值；s_i的值范围为(-1，1)。平均轮廓系数指标解释如下：0.71~1.00为分类结果非常好，0.51~0.70为分类结果一般，0.26~0.50为分类结果脆弱，小于0.26为分类结果失败。

3 底质分类实验结果与分析

实验数据来自胶州湾某海域，使用EM3000多波束测深仪获得。多波束回波数据和测深数据经过预处理，获得地理坐标下的回波强度数据，在相邻的波束脚印内作内插，计算出每个栅格的坐标值，栅格的声强值作为地理编码。实验海域钻孔取样采集的样本大致可分为砾石、沙、沙质泥、泥4类，从这4类样本中分别各随机选取10份，以采样点为中心选取4 m×4 m的区域统计特征值：分形维数、空隙、多重分形维、均值、标准差和分位数。特征向量主成分分析结果和K-均值聚类结果如图 2所示。为分析本文方法的优越性，利用传统单分形特征方法开展了分类实验。传统方法只考虑一个分形特征即分形维数，结合均值、标准差、分位数时，特征向量主成分分析结果和K-均值聚类结果如图 3所示，全特征分类方法和传统单分形特征方法分类的可靠性统计结果如表 1所示。

从图 2和图 3对比可以看出，本文的多特征方法样本类间散点较传统单分形特征方法更为集中，而在后续的聚类算法分析时，样本数据的集中度对提高类间区别性有很大影响：数据的集中度越高，类间区别性越大，从而分类结果精度越高。因此可以看到，本文方法在增加空隙、分形维两个特征值后，对样本特征向量的数据结构具有一定的改善作用，而且从样本分类图可以更清楚地看到，多特征分形维方法的样本类间距离更大，同一底质类样本数据更为集中。

由表 1的分类可靠性统计结果对比表明，本文的多特征方法分类精度检验结果总体由传统方法的85%提高到92.5%，其中，淤泥由80%提高到90%，沙质泥由70%提高到80%，沙由90%提高到100%，除砾石外，其他底质分类精度都得到一定的提高，其中，淤泥与沙质泥二者的分类错误有明显的减少。究其原因，主要是由于淤泥和沙质泥的回波强度图像纹理较为接近，而针对纹理接近的回波图像，传统方法统计的单一分形特征值和基本灰度统计值对其分类都有一定的欠缺性，因分形维数、均值、方差和中位数等都未考虑到图像的细节信息，更偏重于对图像整体特征的统计，从而造成图像细节信息的缺失，这些不足导致实际中在对纹理较为接近的底质进行分类时，分类错误更容易发生，增加分类难度。针对图像局域信息和图像不同层次信息的缺失、空隙和多重分形特征值有一定的优越性。空隙用来测量实际质量和期望质量之间的差异，当纹理细密时，空隙值小, 而纹理粗糙时，空隙值大；多重分形则能够从全局和局部两方面充分描述图像的奇异性结构，对于图像表现出来的各向异性和不一致尺度而言, 多重分形含有可供识别的连续谱指数。因此，在传统单分形统计特征向量组基础上增加空隙和多重分形维特征值后，图像纹理有更详尽的表达，新图像特征值向量组更加完备，从而区分精度得到提高。

从表 1总体分类精度也可以看出，本文顾及空隙、多重分形维特征值后，由于兼顾了图像全局和局部信息，不仅对较近纹理区分有一定的优越性，提高了近似底质分类的可靠性，而且总体分类精度也得到了提高。

为进一步验证本文的多特征分类方法的聚类性，将其与传统方法分类结果的轮廓指标进行对比，如图 4所示，二者平均轮廓系数分别为0.602和0.589，多特征方法的平均轮廓系数指标均值较大，根据平均轮廓系数指标解释分析可以认为，本文提取的多特征向量在数据结构方面有所改善，凝聚度和分离度都有不错的效果，聚类具有优良性。

4 结束语

通过对本文实验结果进行分析，可以得到如下结论：① 不同海底底质类别的分形维数不同，而且分形维对图像尺度要求灵活，可以作为底质分类特征；② 分形特征和灰度特征有一定的相关性，需对特征空间作主成分分析以组成新的图像特征向量组，更充分表达回波图像信息，得到可靠的底质分类结果；③ 本文提出的多特征底质分类方法顾及分形维数、空隙和多重分形等特征，考虑了图像局域信息和图像不同层次信息，对于相近分形维数的图像有较好的区分效果。