水下地磁辅助导航系统采用地磁标准差、粗糙度、相关系数、坡度标准差和地磁信息熵等参数作为评价候选地磁匹配区适配性能的度量。不同参数反映了候选地磁匹配区不同方面的性能，能够反映地磁图的离散程度、局部起伏状态、光滑程度等固有属性^[16]。地磁图通常是用一系列的离散点来描述地磁场在空间中的变化特征。

设某地磁场区域的经纬度跨度为$M\times N$网格（$M$为纬度跨度，$N$为经度跨度），$f(i,j)$为网格点坐标$(i,j)$处的地磁异常值。从宏观起伏特征（地磁费歇信息量、地磁标准差、累加梯度均值）、微观破碎特征（地磁信息熵、粗糙度、坡度标准差）及自相似特征（相关系数）等7种地磁图特征参数的定义如下^[17]。

1）地磁标准差（$\sigma$）。$\sigma$描述的是候选匹配区地磁场的整体起伏变化，其数值越大，相应地磁特征就越明显。计算公式为：

式中，$\overline{f}=\frac{\mathrm{1}}{MN}{\displaystyle \sum _{i=\mathrm{1}}^M}{\displaystyle \sum _{j=\mathrm{1}}^N}f\left(i,j\right)$表示候选匹配区地磁异常的均值。

2）地磁粗糙度（$r)$。$r$描述的是候选匹配区地磁场的平均光滑程度，其数值越大，表示该候选匹配区的地磁信息量越丰富。计算公式为：

式中，$r_x$表示该区域$x$方向粗糙度；$r_y$表示该区域$y$方向粗糙度。。

3）地磁相关系数（$\rho$）。$\rho$表示地磁数据之间的相关程度，其数值越小，表明该区域地磁信息越丰富。计算公式为：

式中，${\rho }_x$表示该区域$x$方向的相关系数；${\rho }_y$表示该区域$y$方向的相关系数。

4）地磁坡度标准差（$S_{\sigma })$。$S_{\sigma }$描述的是候选匹配区地磁坡度的整体情况。计算公式为：

式中，$S(i,j)$表示地磁坡度值；$\overline{S}$表示地磁坡度均值。

5）地磁信息熵（$H$）。$H$描述的是候选匹配区的平均信息量，$H$值越小，代表地磁场特征越独特。计算公式为：

式中，$p(i,j)=\frac{f(i,j)}{{\displaystyle \sum _{i=\mathrm{1}}^M}{\displaystyle \sum _{j=\mathrm{1}}^N}f(i,j)}$。当对数的底数取2时，熵的单位为bit。

6）地磁费歇信息量（Fisher information content，$\mathrm{F}\mathrm{I}\mathrm{C})$。$\mathrm{F}\mathrm{I}\mathrm{C}$可度量地磁场所包含的信息量，其值越大，表明地磁信息越丰富。计算公式为：

式中，$\nabla h(i,j)$表示在地磁场区域网格点$(i,j)$处的梯度向量。

7）累加梯度均值（$T$）。累加梯度代表地磁场的复杂度，其均值描述的是地磁图单位面积的地磁变化。计算公式为：

式中，$T_x(i,j)$表示网格点$(i,j)$处东向的梯度值；$T_y(i,j)$表示网格点$(i,j)$处北向的梯度值。

§1.1介绍的7个地磁图特征参数分别从不同角度刻画了地磁图的特征，但参数间会存在共线性问题，进而使所得结果具有较高的预报方差和病态。同时，相关特征参数会增加输入维数，计算量也会大幅提高^[18]。因此，需要借助降维的思想进行地磁图特征参数的筛选，将原有的多个特征参数进行线性变换，得到少数几个新的无关联变量来近似描述原有特征信息，称为PCA方法。

设候选地磁匹配区的个数为$m$，特征参数个数为$n$，具体流程为：

1）特征参数的同趋势化处理。地磁图特征参数既包含如粗糙度、地磁标准差等高优指标，又包含如相关系数、信息熵等低优指标。为确保主成分具有相同的趋势性，需要采用取负法进行同趋势化处理，即$x{\text{'}}_{ij}=-x_{ij}$，$x_{ij}$表示低优指标。处理后的评价矩阵表示为$\mathit{X}\text{'}={\left[x{\text{'}}_{ij}\right]}_{m\times n}$$(i=\mathrm{1,\; 2}\cdots m;j=\mathrm{1,\; 2}\cdots n)$。

2）特征参数的标准化处理。各特征参数的量纲存在不一致现象，容易导致奇异数据的出现，并且标准化处理后仍能保持原始数据的基本特征。计算公式为：

式中，${\overline{x}}_j=\frac{\mathrm{1}}{m}{\displaystyle \sum _{i=\mathrm{1}}^m}x{\text{'}}_{ij}$表示第$j$个特征参数的均值；${\sigma }_j=\sqrt[]{\frac{\mathrm{1}}{m-\mathrm{1}}{\displaystyle \sum _{i=\mathrm{1}}^m}(x{\text{'}}_{ij}-{\overline{x}}_j{)}^{\mathrm{2}}}$表示第$j$个特征参数的标准差。处理后的评价矩阵表示为$\mathit{Z}={\left[z_{ij}\right]}_{m\times n}(i=\mathrm{1,\; 2}\cdots m;j=\mathrm{1,\; 2}\cdots n)$。

3）标准化评价值相关矩阵$\mathit{R}$。计算公式为：

式中，$r_{kj}=r_{jk}$，$r_{jj}=\mathrm{1}$，且$r_{kj}=\frac{\mathrm{1}}{m}{\displaystyle \sum _{i=\mathrm{1}}^m}{z}_{ik}{z}_{ij}$。

4）$\mathit{R}$的特征值、特征向量和方差贡献率的求解。求解得到的特征根${\lambda }_g$按从大到小的顺序排序，其表征各主成分在综合评价中的作用大小，同时求得特征向量${\mathit{L}}_g$。方差贡献率为每个主成分表征原评价指标变量的信息量，表示为：

通过累计贡献率$\alpha (k)$可以确定$k$的值，这里判定准则为$\alpha (k)={\displaystyle \sum _{g=\mathrm{1}}^k}{\alpha }_g\ge \mathrm{99}\mathrm{\%}$。

5）主成分个数确定。计算结果矩阵$\mathit{T}=\mathit{Z}\mathit{L}$，$\mathit{L}$表示特征矩阵，$T_i$表示第$i$个主成分，这里取前$k$ ‍个主成分。

BP神经网络是一种非线性不确定性数学模型，通过不断修改网络的权值和阈值使得均方误差最小，最终实现数据的高精度拟合或分类。但BP算法采用的是梯度下降法，初始权值的选取不同会导致网络发散或陷入局部极值点。遗传算法（genetic algorithm，GA）作为一种随机搜索算法，在搜索过程中并不需要借助梯度信息，而是将网络的权值和阈值作为种群个体仿照基因编码，经过选择、交叉和变异等遗传操作迭代进化得到种群最优个体，解码最优个体获得BP神经网络的近似最优解。这种利用GA对BP神经网络进行优化的模型称为GA-BP神经网络模型^[19-21]。

GA的参数设置为：最大进化代数为100，种群规模为50，优良子辈个数为2，交叉概率为0.6，变异概率为0.001。

选取分类正确率（C_i）、平均正确率（A）及全局分类正确率（O）作为评估指标^[22]，分别表示为：

式中，$N\left(R={\omega }_i\right)$表示真实类别为${\omega }_i$的样本数；$N(P={\omega }_i\bigcap R={\omega }_i)$表示真实类别为${\omega }_i$的样本中被实际正确筛选出的样本个数$\left(i=\mathrm{1,\; 2}\right)$；$C_i$表示单类分类正确率，表征分类模型对单一类别的分类性能；A表示平均正确率，表征单一分类性能的平均度量；O表示全局分类正确率，表征分类模型的整体分类性能。

借助特定地磁图特征参数和适配性标签，通过对GA-BP神经网络训练，得到两者的网络映射。将候选地磁匹配区的地磁图特征参数作为输入量，经过神经网络分类获取对应的适配性标签，实现地磁适配区的智能选取。

神经网络构建和训练过程如下：

1）将某一特定水下区域作为候选地磁匹配区进行研究，计算得到该区域的若干个地磁图特征参数，借助PCA提取出新的独立特征参数。

2）利用地磁滤波匹配算法，本文选用无迹卡尔曼滤波（unscent‍ed Kalman filtering，UKF）算法在任务区域开展地磁匹配定位实验，选取匹配概率$P_r$^[22]作为衡量适配性指标。将1或0设置为适配标签，结合精度要求，当$P_r\ge \mathrm{0.7}$时，类别${\omega }_{\mathrm{1}}$标定为1；当$P_r<\mathrm{0.7}$时，类别${\omega }_{\mathrm{2}}$标定为0。

3）将前两步中的独立特征参数及其对应的适配标签作为BP神经网络训练样本的输入向量，通过GA优化BP神经网络参数，经过多次迭代计算和样本训练，获取可靠稳定的神经网络。

4）将待研究区域的地磁图划分成若干子区域，重复步骤1）获得若干独立特征参数，并将其输入到训练好的神经网络中，实现对整个研究区域地磁适配性的识别和分类。

2017年5月，美国国家地球物理数据中心利用卫星、船载、航空和地面磁测数据发布的全球地磁异常格网数据（Earth magnetic anomaly grid at 2 arc-minute resolution, version 3, EMAG2 v3）模型，相较于EMAG2（2009）有重大更新。该模型分辨率仍为2′，重点针对船载和航空磁测数据进行更新，新推出了海洋区域的海平面格网数据产品，并利用克里格算法进行海洋稀疏航线地磁异常格网的插值改进，极大提高了全球海洋区域的地磁异常格网模型精度^[23]。因此，EMAG2 v3模型能够为AUV水下高精度长航时自主导航提供良好的背景场条件，更有利于仿真验证与分析。本文选取EMAG2 v3模型区域范围为125ºE~135ºE、25ºN~35ºN的中国东海部分海域，利用该部分海域的地磁异常信息展开研究和仿真实验分析。

将上述海域等分成400块格网作为样本，每个样本的格子尺寸为$\mathrm{15}\times \mathrm{15}$，如图 1所示。分别计算各子区地磁图特征参数，利用PCA将若干原始地磁图特征参数转化为新的主成分特征参量；通过对GA-BP神经网络训练，实现对各子区的类别划分，给出地磁适配区/非适配区标签。从400 ‍块子区中随机选取120块区域作为训练样本，剩下280块区域作为测试样本，进行自组织优化分类方法性能评估实验。

图  1  候选地磁匹配区划分

Figure 1.  Division of Candidate Geomagnetic Matching Area

借助PCA依次提取各子区的7个地磁图特征参数，按照特征属性的同趋势化处理、标准化处理、相关矩阵及其方差贡献率等的求解和综合评价等步骤进行主成分分析，求得每个主成分对应的贡献率和累积贡献率如图 2所示。

图  2  PCA主成分方差信息分布图

Figure 2.  Distribution Map of PCA

由图 2可得，前4维主成分包含原始数据99%的方差信息，后3维主成分为冗余信息。稳健起见，选取前5维主成分作为独立输入向量开展后续神经网络训练。

由于输入向量维数降低，使得BP神经网络训练中的神经元个数减少，训练效率也相应得到提高。神经元个数的确定与网络的性能息息相关，需要确定最佳网络神经元的个数。目前该问题并没有统一的解决方案，通常由经验公式或者试算来确定。通过BP神经网络训练来确定最佳神经元个数即隐层节点数，如图 3所示。

图  3  全局分类正确率随隐层节点数变化曲线

Figure 3.  Variation Curve of O with the Number of Hidden Layer Nodes

由图 3可得，隐层节点数为7时，全局分类正确率最高，因此这里神经元个数确定为7。

利用本文提出的基于PCA和GA-BP神经网络结合的候选地磁匹配区自组织优化分类方法，从特征选取（feature selection，FS）和参数优化（parameter optimization，PO）两个角度来提高候选地磁匹配区的分类正确率。划分成4种分类方法：（1）无特征选取，无参数优化（Nono）；（2）无特征选取，利用GA优化BP神经网络（PO）；（3）利用PCA特征选取，无参数优化（FS）；（4）利用PCA特征选取，利用GA优化BP神经网络（FSPO）。不同分类方法的性能如表 1所示。

由表 1可得，Nono方法的O相较于其他3种方法是最低的，根本原因在于选取的7个地磁特征参数包含较多的相关信息，造成地磁特征信息的冗余，并且BP神经网络参数没有根据实际情况设置，导致分类结果并不理想。同Nono方法对比，PO和FS方法分别将O提高了6.43%和5%，表明无论从特征选取角度还是参数优化角度，两种方法均能有效改善分类性能。本文提出的FSPO方法分类性能表现最好，其在完成独立特征选取的前提下，同时通过GA来优化神经网络参数，使得O达到96.79%，将分类错误率维持在较低水平，提高了候选地磁匹配区的分类精度。

从测试样本中随机选取3个候选地磁匹配区进行匹配验证分析，并计算相应区域的匹配概率，按照自组织优化分类方法的一般步骤进行适配区分类，结果如表 2所示。

由表 2可知，基于GA-BP神经网络得到的分类结果与相应区域的匹配概率具有良好的一致性，即分类标签为地磁适配区的区域2实际也具有较好的适配性能，分类标签为地磁非适配区的区域1和区域3的实际适配性能则一般。上述方法给出的结果可信，可将区域2选作地磁适配区。

为进一步验证该自组织优化分类方法所得结果的正确性和可靠性，采用基于地磁异常的UKF滤波算法对上述随机选取的3个区域进行仿真实验，惯导航迹仿真参数设置如表 3所示。航迹仿真中AUV的运动状态可以描述为：匀加速运动→匀速运动→转弯→匀速状态。起始点分别选定在相应子区域的左下角点位置。航迹匹配结果如图 4所示，候选地磁匹配区沿经度和纬度方向的绝对匹配误差统计如表 4所示。

图  4  候选地磁匹配区的航迹匹配结果

Figure 4.  Track Matching Results of Candidate Geomagnetic Matching Areas

综合图 4航迹匹配结果及表 4经纬度方向绝对匹配误差统计，区域2作为地磁适配区具有良好的匹配效果，匹配航迹最接近真实航迹，其经纬度方向的绝对匹配误差也是最小的；而区域‍1和区域‍3作为地磁非适配区出现误匹配现象，相较于区域1，区域3的匹配航迹漂移更为严重。这进一步验证了匹配概率、适配性能及分类标签具有很好的一致性，同时验证了上述方法分类结果的可靠性。

将各子区域经过PCA获取的独立特征参数输入到训练好的神经网络中，输出得到整个区域的适配性标签。图 5为中国东海部分海域的地磁适配区栅格化选取结果。

图  5  地磁适配区分布

Figure 5.  Distribution of Geomagnetic Matching Areas

图 5中，空白区域表示地磁适配区，阴影区域表示地磁非适配区。统计可得，该区域地磁适配区的数量为283个，地磁非适配区为117个。图 5给出的分类结果可为AUV航迹规划问题提供支持，即在航迹规划前，基于所赋予航行任务的先验知识，可以在AUV允许航行范围内优先选择若干地磁适配区，这为AUV的高精度长航时自主导航提供重要保障。

本文围绕地磁场适配性分析问题，首先，利用PCA选出独立的地磁特征参量；然后，在此基础上利用GA从参数优化角度对BP神经网络的分类性能进行优化；最后，设计了一种基于自组织优化分类的地磁适配区选取方法，实现整个研究区域候选地磁匹配区的自动识别和分类。在上述仿真条件下，经验证，该自组织优化分类方法具有较高的分类精度和可靠性，约束AUV尽可能经过地磁适配区能够有效提高其导航定位精度，能够为AUV完成特定的水下任务提供重要保障。

本文主要贡献如下：（1）利用PCA进行候选地磁匹配区的特征选择，保证了神经网络训练输入参量的独立性，避免了神经网络结果存在较高的预报方差和病态，同时降维也减小了网络计算量。（2）将GA应用到BP神经网络初始的权阈值优化，解决了BP神经网络极易陷入局部极小值的问题，加快了网络的收敛速度，最终得到BP神经网络的近似最优解。