湿地被誉为“地球之肾”、“生物超市”和“基因库”, 与森林、海洋并称为全球三大生态系统, 在气候调节、促淤造陆、蓄洪防旱等方面发挥着重要作用^[1]。黄河三角洲湿地是由黄河冲淤和海水动力共同作用的结果, 是我国暖温带保存最完整、最广阔和最年轻的湿地, 是典型的河口湿地^[2]。目前, 黄河三角洲湿地受入海水沙量减少、人类活动和全球气候变化等多重因素的影响, 正面临湿地萎缩、景观破碎化及生物多样性降低等问题 ^[3]。因此, 开展黄河三角洲湿地地物分类研究, 对湿地资源的保护、开发和利用具有重要意义。

高光谱数据具有光谱分辨率高、光谱波段多以及“图谱合一”的特点, 针对不同的地物类型可以表现出不同的光谱信息, 基于上述特点, 许多学者利用高光谱数据开展了湿地地物分类研究^[4-5]。王建步等^[6]将线性光谱混合分析模型与归一化植被指数(NDVI)和归一化水体指数(NDWI)相结合, 实现了高光谱图像的湿地分类。初佳兰等^[7]利用ISODATA非监督分类方法与最大似然法(ML)和支持向量机(SVM)法监督分类法, 提出了一种基于众数赋值的高光谱图像地物分类方法, 有效提高了湿地分类精度。李世波等^[8]对高光谱数据进行主成分分析(PCA)降维后, 利用随机森林和径向基内核支持向量机分类, 得到了较稳定的湿地分类精度。然而, 上述研究多是利用人工提取的特征进行分类, 提取特征质量依赖人工经验, 自动化与智能化程度不高。

近年来, 深度学习技术在计算机视觉领域的巨大成功, 使其在遥感图像处理领域也得到了广泛应用^[9-11], 卷积神经网络CNN(Convolutional Neural Network)作为其代表, 已被广泛用于湿地分类研究^[12-13], 提高了分类算法的自动化与智能化。Hu等^[14]提出了基于多目标CNN和决策融合模型的高光谱分类算法, 使用卷积神经网络进行特征提取, 引入联合投票法和模糊隶属度的决策融合模型进行湿地分类。Chen等^[15]利用深度置信网络, 引入共轭梯度更新, 使用以中心像素的光谱块特征进行分类, 在GF-5黄河三角洲影像数据上取得了较好效果。需要注意的是, 基于CNN的算法关键是对高光谱数据进行特征提取, 然而多数算法的特征提取依赖固定尺寸、固定感受野的卷积核, 这使得在特征提取过程中, 信息只局限在小范围内, 更多的全局信息得不到利用, 限制了分类精度的提高。

图卷积神经网络GCN(Graph Convolutional Neural Network)因其灵活的感受野为解决上述问题提供了新思路^[16-17]。Hong等^[18]设计了MiniGCN网络用于高光谱图像分类, 提高了分类精度。Ding等^[19]设计了联合上下文的多尺度图卷积网络, 同样在高光谱分类中取得不错效果。然而, 上述方法的图结构都是通过计算欧几里得距离来构造的, 欧几里得距离只能反映光谱值之间的绝对差异, 忽略了光谱波形间的相对差异。在湿地中, 水体、潮滩、植被等界限并不明显, 仅仅利用光谱之间的绝对差异无法完成分类, 单一的图构造方法很难完成湿地地物分类任务。

为解决上述问题, 本文提出了基于双路图卷积的黄河三角洲湿地分类模型, 称为DGCNWC算法, 构建欧式图和余弦图挖掘光谱间的绝对差异与相对差异; 利用图卷积提取特征, 获取全局信息; 在CHRIS和GF-5数据集上取得了较好结果。

1 双路图卷积网络

尽管当前基于GCN的算法在湿地分类中已经取得了不错效果, 但其对湿地复杂地物的边界区分不明显。因此, DGCNWC算法利用全局信息, 聚焦光谱值与光谱波形特征提取, 之后进行特征融合, 开展湿地分类研究。算法整体框架图如图 1所示。DGCNWC算法包括两个模块: 图结构数据构建模块、特征提取与融合模块。

1.1 欧式图与余弦图的构建

与基于CNN的方法相比, 基于GCN的方法具有更灵活的感受野, 能够获取全局信息, 这得益于数据的图结构表示。初始图结构直接决定了能够提取的特征的质量。因此, 本文设计了欧式图和余弦图两种图结构用于表达湿地高光谱影像的图结构数据。

欧式图表示光谱之间的绝对差异, 利用欧几里得距离来度量, 欧式图可以直接反映每个像素之间光谱值的差异, 衡量光谱的相似性。余弦图则表示的是光谱曲线之间的相对差异, 利用余弦相似度度量, 余弦图体现的是每个像素光谱波形之间的差异, 衡量光谱波形之间的相似性。光谱值差异有利于判断不同的地物类别, 光谱波形差异可以较好的区分潮滩、植被等有明显差异的地物界限。

基于上述原理, 假定高光谱影像为$H = \left( {{h_1},} \right.$ $ \left. {{h_2}, \cdot \cdot \cdot ,{h_n}} \right) \in {{\text{R}}^{{\text{1}} \times {\text{1}} \times {\text{D}}}} $, 其中h_i表示高光谱影像的每一个像素, D表示影像的光谱维度, n表示影像的像素个数。对高光谱影像的每个像素为中心, 取矩形邻域可以得到$ H'{\text{ = (}}{{{h'}}_1},{{{h'}}_2}, \cdot \cdot \cdot ,{{{h'}}_{{n}}}{\text{)}} \in {{{R}}^{{{r}} \times {{r}} \times {{d}}}} $, 其中h′_i表示以h_i像素为中心, 取矩形邻域后的像素块, r表示邻域的边长, d表示经过主成分分析(PCA)降维后的特征维度。对于每个像素块h′_i, 都有两个无向图表示其欧式图与余弦图。

将任一像素块h′_i的欧式图定义为G_ed=(V_ed, E_ed)。其中, V_ed代表像素块h′_i中的所有像素, 表示图G_ed的顶点, E_ed是像素块内每个像素间光谱值的绝对差异, 由欧几里得距离度量, 采用高斯核函数^[20], 其计算公式如下:

$ E_{ij}^{ed} = {e^{ - \lambda {{\left\| {{{h'}_{ii}} - {{h'}_{ij}}} \right\|}^2}}}, $

(1)

式中, h′_ii, h′_ij表示像素块h′_i中的像素, λ为尺度系数, 一般设为0.2^[20], $E_{ij}^{ed}$表示像素块h′_i中的像素h′_ii与其他像素的欧几里得距离。

对于任一像素块h′_i的欧式图G_ed, 其初始特征矩阵$X = \left[ {{x_1},{x_2}, \cdot \cdot \cdot ,{x_n}} \right] \in {R^{d \times n}}$, 其中, x_i表示像素块内每一像素的光谱值, d表示每个像素x_i的光谱维度, n表示像素块内的像素个数。其邻接矩阵$ {A_{ed}} \in {R^{n \times n}} $, 邻接矩阵表示的是像素块内的每一个像素与其他像素的邻接关系, 其计算公式如下:

$ a_{ij}^{ed} = \left\{ {\begin{array}{*{20}{c}} {E_{ij}^{ed}rank\left( {E_{ij}^{ed}} \right) \geqslant k} \\ {1\;\;\;\;\;\;\;\;\;\;\;\;\;\;i = j} \\ {0\;\;\;\;\;\;\;\;\;\;\;\;others} \end{array}} \right. , $

(2)

式中, $ a_{i j}^{e d}$表示像素块h′_i中的像素h′_ii欧式图的邻接系数, $ {rank}\left(E_{i j}^{e d}\right) \geqslant k$, 表示$ E_{i j}^{e d}$降序中最大的k个值, k为经验参数, 一般为10~30且小于r^{2 [21]}, 当$ a_{i j}^{e d}=0$时, 认定两个像素无连接关系。

同理, 任一像素块h′_i的余弦图定义为G_cos=(V_cos, E_cos)。V_cos代表像素块h′_i中的所有像素, 表示图G_cos的顶点, E_cos是像素块内每个像素间光谱值的相对差异, 由余弦距离度量, 其计算公式如下:

$ E_{i j}^{\cos }=e^{-\lambda\left(\frac{h_{i i}^{\prime} \cdot h_{i j}^{\prime}}{\left\|h_{i i}^{\prime}\right\|\left\|h_{i j}^{\prime}\right\|}\right)}, $

(3)

式中, $ h_{i i}^{\prime} \cdot h_{i j}^{\prime}$表示像素块h′_i内两个像素光谱值组成的向量的点积, $ \left\|h_{i i}^{\prime}\right\|\left\|h_{i i}^{\prime}\right\|$表示两个上述两个向量的模相乘, λ为尺度系数, 同公式(1), 一般设为0.2。

对于任一像素块h′_i的余弦图G_cos, 其初始特征图仍为X, 其邻接矩阵为$ A_{\cos } \in R^{n \times n}$, 参数选择与欧式图相似, 计算公式如下:

$ a_{i j}^{c o s}=\left\{\begin{array}{cc} E_{i j}^{c o s} & \operatorname{rank}\left(E_{i j}^{c o s}\right) \geqslant k \\ 1 & i=j \\ 0 & \text { others } \end{array} .\right. $

(4)

1.2 特征提取与融合

欧式图与余弦图构造完成后, 将利用其特征矩阵和邻接矩阵输入到特征提取与融合模块。特征提取模块由两路网络组成, 一路为两层图卷积网络用于提取欧式图的绝对差异光谱特征; 另一路为图U-Net网络, 用于提取余弦图的相对差异光谱特征。融合模块由一层卷积神经网络和两层全连接网络组成, 用于融合绝对差异光谱特征和相对差异光谱特征, 利用融合的特征进行最终的湿地地物分类。

对于欧式、余弦图的邻接矩阵, 其对角线的值表示各像素的自相关性。为了提高模型训练稳定性, 采用自连接的邻接矩阵$\tilde{A}$扩大其自相关性, 表达式如下:

$ \tilde{A}=A+I_n, $

(5)

式中, I_n代表对角线为1, 其余值为0的单位矩阵。

对欧式图, 受限于图卷积信息传播特性, 三层及以上图卷积会出现“过渡平滑”现象, 导致信息不具有可分性, 一层图卷积信息提取不充分, 因此采用双层图卷积进行绝对差异光谱特征的提取。图卷积层通过不断聚合相邻节点的信息实现邻接信息的传递, 在不断传递过程中, 可以遍历所有节点, 从而获取光谱特征的全局信息。图卷积方程定义如下:

$ X^{l+1}=\sigma\left(\tilde{D}-\frac{1}{2} \tilde{A} \tilde{D} \frac{1}{2} X^l W^l\right), $

(6)

式中, X^l为当前l层的特征矩阵, X^l+1表示经过图卷积后l+1层的特征矩阵, W^l表示可学习的权重参数, σ(ּ)表示激活函数, 本算法采用的是ReLU(ּ)。$ \tilde{D}$为对角矩阵, 表示邻接矩阵的度矩阵, 其对角元素计算公式为:

$ \tilde{D}_{i, i}=\sum\nolimits_j \tilde{A}_{i, j}, $

(7)

对余弦图, 因其蕴含了光谱波形的相似性信息, 只通过图卷积层进行信息传递不能充分挖掘到不同地物的边界信息, 因此, 欧式图双层图卷积的基础上, 引入了四组图池化与图反池化操作, 最后利用残差结构加速模型的收敛。其中, 图池化操作如图 2所示, 其方程定义如下:

$ y=X^l P^l/\left\|P^l\right\| , $

(8)

$ i d x={rank}(y, k), $

(9)

$ \tilde{y}=\operatorname{sigmod}(y(i d x)), $

(10)

$ \tilde{X}^l=X^l(i d x)\text{，} $

(11)

$ A^{l+1}=A^l(i d x, i d x), $

(12)

$ X^{l+1}=\tilde{X}^l \odot(\tilde{y})\text{，} $

(13)

式中, P^l是第l层的一个可训练的投影向量, X^l是第l层的特征矩阵, $ \|\cdot\|$用于获得P^l的欧几里得范数, y是将图节点特征矩阵映射到了一维空间的结果。rank(ּ)是一个排序操作, 对y进行降序排序, 获取最大的k个值的坐标, 同时也为池化后的图节点坐标。利用公式(10)~(13)获取池化后的特征矩阵与邻接矩阵。

图反池化利用“零”填充, 对于反池化的图, 构建一个与其反池化后同样尺寸的节点值全为0的图, 随后将反池化图的节点值填入, 其填入位置对应其池化层的位置, 保证每层反池化的图尺寸都能与相应池化层的图尺寸相同。反池化的方程定义如下:

$ X^{l+1}={distribute}\left(0_{N \times C}, X^l, i d x\right) \text {, }, $

(14)

式中, 0_N×C为初始0矩阵, 与X^l+1尺寸相同; distribute(ּ)操作为利用X_l, 在与idx相对应的位置更新0矩阵。

同时, 为了更好地利用特征, 使用了Nonlocal注意力机制关注有效信息, 获得最终的相对差异光谱特征。余弦图特征提取的图U-Net流程如图 3所示。

最终, 为了融合欧式图提取到的绝对差异光谱特征F_ed与余弦图提取到的相对差异特征F_cos, 特征融合模块定义如下:

$ F_{{fusion }}=F C\left({conv}\left(F_{e d} ; F_{c o s}\right)\right) , $

(15)

式中, conv(ּ)表示卷积操作, FC(ּ)表示全连接。

利用卷积层学习融合参数, 利用全连接层进行降维, 得到融合特征F_fusion。最后将其输入到Softmax分类层, 预测每一地物分类概率, 得到最终分类结果。

2 实验分析与讨论 2.1 数据集

为了验证所提的双路图卷积网络在黄河三角洲湿地地物分类的有效性, 本文分别选取了CHRIS卫星和GF-5卫星的两景高光谱遥感影像进行实验。

CHRIS数据集区域位于黄河入海口新老入海径流清八汊和清水沟交界处, 采集于2012年6月, 像素尺寸为510×511, 共有18个波段, 其光谱范围406~1 036 nm, 光谱分辨率从5.9~44.1 nm不等, 空间分辨率17 m, 包含芦苇、互花米草、滩涂、水体、裸地以及柽柳碱蓬混生区6个地物类型, 地物真值图由当年夏季实地踏勘所得, 以站点或断面的形式开展滨海湿地地物类型种类、覆盖度和边界范围等信息调查。

GF-5数据集覆盖了黄河三角洲自然保护区滨海湿地, 采集于2018年11月, 像素尺寸为462×617, 共有150个波段, 光谱范围为390~1 029 nm, 光谱分辨率为3.67~4.81 nm不等, 空间分辨率为30 m, 包含互花米草、芦苇、柽柳、盐地碱蓬、潮滩芦苇、裸潮滩、盐碱滩以及水体8个地物类型, 由当年秋季实地踏勘所得, 以站点或断面的形式开展滨海湿地地物类型种类、覆盖度和边界范围等信息调查。

每个数据集, 针对每个地物类型随机抽取200个作为训练数据, 其余地物类型用于测试。因此CHRIS数据总量为260 610个, 训练样本为1 200个, 测试样本为259 410个; GF-5数据集总量为285 054个, 训练样本为1 600个, 测试样本为283 454个。

2.2 实验环境及评价指标

本文所有实验开展的计算机环境为: Windows 10操作系统、i9-10900K CPU、NVIDIA GeForce RTX 3090 GPU以及128G内存。本文所提双路图卷积算法使用Pytorch实现。

为了有效评价算法性能, 本文采用总体精度(Overall Accuracy, OA)、平均精度(Average Accuracy, AA)以及Kappa系数(K)作为算法评价指标。

2.3 算法性能对比与讨论

为了充分验证所提算法, DGCNWC算法在黄河三角洲湿地分类的有效性, 本小节对DGCNWC算法开展欧式支路与余弦支路的单支路消融实验, 同时与SVM、3DCNN^[22]、DGU^[21]算法进行对比实验。SVM算法是机器学习的经典算法; 3DCNN算法是利用3D卷积神经网络, 实现光谱和空间联合的算法, 在多个公共数据集上取得了优异性能; DGU算法是利用欧几里得距离构建光谱图和空间图的图卷积算法。两组数据集每种地物类型分类精度、测试训练样本比以及总体实验结果如表 1、表 2所示, 分类结果的可视化图像见图 4、图 5。

(1) PCA维度分析

为了验证PCA降维后的数据对分类性能的影响, 本小节探究了PCA降维后不同的维度对分类性能的影响, 因为CHRIS只有18个波段, 数据冗余度较小, 不对其进行降维, 只对GF5数据集进行PCA降维, 实验结果如表 3所示。

由表 3结果可以看出, PCA对GF5数据集分类结果的影响, 在60~110维度间呈现先上升后下降的趋势。这表明在80维度, 是GF5数据降维的最佳维度, 分类精度取得了80.61%; 在60或70维度, 由于降维过多, 导致光谱信息产生了丢失, 精度下降; 同时在90~110维度, 因为保留成为过多, 光谱信息中仍存在冗余, 也导致精度下降。因此后面实验中, GF5数据经PCA降维后维度采用80。

(2) DGCNWC算法消融实验

在消融实验中, 分别研究了欧式支路、余弦支路对分类精度所产生的影响。实验结果表明, 在CHRIS数据集中, 欧式支路的总体分类精度为83.16%, 平均分类精度为86.79%; 余弦支路的总体分类精度为77.62%, 平均分类精度为82.98%。均低于DGCNWC算法的总体精度88.18%, 平均精度90.28%。这验证了所提算法利用光谱值绝对差异与光谱波形差异的作用。同时, 在各地物类型分类精度上, 欧式支路绝大部分优于余弦支路, 而在柽柳碱蓬混生区、互花米草的分类上, 余弦支路分别以81.29%、96.87%优于欧式支路。因其生长习性原因, 柽柳碱蓬混生区与互花米草是CHRIS数据集地物边界较为复杂的两类, 这表明余弦支路可以利用波形之间的差异性识别不同地物类别边界, 验证了余弦图是对欧式图的有效补充。

针对GF5数据集, DGCNWC算法以总体精度80.61%, 平均精度78.30%的结果优于欧式支路和余弦支路, 同样揭示了双路图卷积对黄河三角洲湿地分类的有效性。需要注意的是, DGCNWC算法在芦苇、潮滩芦苇、盐地碱蓬等地物上的分类精度较单分支网络提升了7%~34%不等, 如此大的分类差距在CHRIS数据集上没有发生。考虑到GF5原始数据含有150个波段, 波段多, 数据复杂, 在对其进行大气校正、辐射校正时容易引入噪声, 其数据质量与CHRIS数据相比稍差, 因此单分支网络受到影响导致精度偏低。而双路网络受此影响较小, 各地物精度都比较高, 表明融合双支路后, DGCNWC算法有着不错的鲁棒性。

(3) 不同算法对比实验

在对比实验中, 开展了机器学习算法SVM、基于CNN的算法3DCNN、基于GCN的算法DGU在分类精度上与DGCNWC算法的对比研究。实验结果显示, 在CHRIS数据集上, DGCNWC算法的总体分类精度达88.18%, 高于SVM的83.98%, 3DCNN的83.33%, DGU的84.93%, 同时在GF5数据集上, DGCNWC的总体分类精度都分别高于SVM、3DCNN、DGU有7%、4%、10%左右。这些现象表明, DGCNWC算法面对各种类型的分类算法在黄河三角洲湿地分类上都有较为明显的优势, 具有精度更高的分类性能。

同样的结果, 在Kappa系数的结果中也有体现, 在CHRIS和GF5数据集中, DGCNWC算法分别以0.844 7、0.749 2均高于SVM、3DCNN、DGU算法, 表明DGCNWC算法的分类结果与实际地物结果具有更高的一致性。

另一方面, SVM、3DCNN、DGU算法在两个数据集上的分类性能各有优劣, 在CHRIS数据集中, DGU总体精度分别高于SVM、3DCNN有0.95%、1.6%, 为三个算法中分类性能最好的; 在GF5数据集中, 3DCNN算法以76.29%的总体精度高于SVM、DGU的73.03%、70.85%。这表明, 深度学习算法在分类性能上确实优于传统机器学习算法, 在黄河三角洲湿地分类领域, 基于CNN的算法与基于GCN的算法具备一定的竞争性。然而, 考虑模型的稳定性发现, 在CHRIS数据集, SVM、3DCNN、DGU的总体精度与平均精度的插值分别在2.04%、3.24%、4.59%, DGCNWC为2.1%。在GF5数据集中, 这些数据分别为1.34%、3.42%、2.56%、2.31%。可以看出SVM、DGCNWC对不同地物类型的鲁棒性更好, 某一地物分类精度不易出现极端值。同样表明DGCNWC在兼顾高分类精度的同时, 又有较好的鲁棒性。

图 6表示了DGCNWC对地物边界区分的效果, 如图 6所示, 以左上角红框标出区域为例, SVM、3DCNN、DGU的地物边界不清晰, 地物混杂。DGCNWC的地物边界较上述3种方法有明显提高, 地物边界清晰, 没有斑点分布。这表明DGCNWC在一定精度范围内解决了地物边界不易区分的问题。

3 结论

针对当前黄河三角洲湿地地物分类算法中存在的全局信息利用不足, 地物边界不易区分等问题, 本文提出了一种基于双路图卷积的黄河三角洲湿地分类算法, DGCNWC算法。DGCNWC算法包含图结构数据构建模块、特征提取与融合模块两部分。图结构数据构建中, 设计欧式图和余弦图, 欧式图的图节点的邻接关系利用欧几里得距离度量, 表示光谱值之间的绝对差异, 余弦图支路利用余弦相似度表示图节点的邻接关系, 体现光谱波形之间的差异, 用以区分地物类型边界。特征提取与融合中, 对欧式图利用两层图卷积层做特征提取, 对余弦图使用图U-Net网络进行特征提取, 最后DGCNWC算法将两路特征进行融合。相较于现有算法, 在黄河三角洲湿地分类的任务中取得了较为满意的结果。实验表明, 在CHRIS、GF5数据集中分类的总体精度分别达到了88.18%、80.61%。

值得注意的是, 虽然DGCNWC算法学习了光谱值差异与光谱波形差异特征, 但在大尺度遥感影像下只考虑光谱特征是不够的, 应当引入空间相关性。因此未来的工作将在此基础上提取空间相关性特征, 进一步提高特征的鉴别能力, 从而提高黄河三角洲湿地地物分类的精度, 更好的服务黄河三角洲的保护、开发和利用。