基于多特征融合的遥感图像河流提取

付宝晶，李自立

（广西师范大学电子工程学院，广西桂林 541000）

0 引言

河流是自然环境中最常见的地形特征之一，也是人类生产生活的重要依靠。从遥感图像中进行河流区域提取，对掌握某地区的水文特性、环境保护和发展建设等都有着十分重要的理论和实际意义。

近年来，河流提取是图像分割的热点，研究算法大致可分为像素阈值法［1］、纹理特征法［2］、形状轮廓法、机器学习法等。阈值法是图像分割中的经典思想，刘炜［3］等通过创建特征空间进行全局阈值分割粗提水域分布信息。此类方法通常能进行有效分割，但在图像复杂或者目标灰度范围接近的区域，很难分割目标。基于纹理特征分割在光学图像中有很好的应用，Moslem Ouled Sghaier［4］等基于局部纹理特征和全局形态学算子的方法对河流和湖泊进行提取和区分。此类方法能够对图像进行全局的纹理分析，但是也忽略了局部纹理细节，也就使得它不能在纹理和像素之间建立联系。形状轮廓提取法在河流提取中也有良好的表现，张祝鸿［5］等利用河流的边缘形状特征运用笔画宽度变换与几何特征集结合的方法提取河流。此类算法依赖于检测到的河流边缘，在受到桥梁、山体阴影的干扰下会影响河流结构的准确性。机器学习近年来在河流提取方面有很大的优势和突飞猛进的发展，李宁［6］等将精致Lee 滤波与卷积结合，基于河道几何优化的权值卷积核，提出一种新型的河道提取网络模型。章斯腾［7］等提出增强河流剖面特征自动提取河流的神经网络，能够观测河流网络的动态变化。其中非监督分类凭自然聚类的特性，在图像尺寸较大、纹理较为丰富的区域容易检错。监督学习效果较好，但效果和输入的训练样本数量和训练次数有关。

综合现阶段的有代表性的研究结果，基于某单一特征或方法进行特定河流区域提取时，存在无法避免的局限性，致使提取的结果不完整或者干扰较大。针对以上的方法局限，本文提出一种基于纹理、光谱和几何形状多特征融合的方法。综合考虑图像河流区域所蕴含的多种信息，获取全面有效的特征描述，大幅提升河流提取的准确率和完整度。

1 遥感图像数据分析

本文采用高分一号拍摄的台山市遥感影像进行实验，该数据的空间分辨率为2 m，实验选取的图像大小为1 500×1 500 像素。图像包含蓝光通道、绿光通道、红光通道、近红外通道4 个波段，如图1所示。蓝光波段的水体衰减系数最小，对水体穿透能力最强，常用于研究河流、海洋水下环境状况等；绿光波段对水体有一定的穿透力，常用来反映水下地形、沿岸沙洲等特征；红光波段位于叶绿素的主要吸收带，常用来区分植被生长状况、覆盖率等；近红外波段位于水体的强吸收区，吸收了水体绝大部分的能量［8］。对比4 个通道在河流流域上的成像特点，红外波段中的水体成像部分与其他地物差异较明显，并且在计算角二阶矩时也能使水体能量保持最大化。因此，本文选取近红外波段图作为河流流域提取的处理对象。

图1 遥感影像四通道图Fig.1 Four bands images of remote sensing image

2 研究方法

首先，基于河流本身具有一定宽度且一致连续的特点，在一段相对较短的距离中，该段河流的任意两个像素可以认为具有相似的灰度统计特性。本文通过计算灰度共生矩阵［9］（Gray Level Co-occurrence Matrix，GLCM）来反映图像中的河流纹理信息，选取河流与其他地物对比效果最明显的角二阶矩图作为下一步的输入。然后进行OTSU［10］双阈值分割滤除绝大部分的非河流区域，但还存在一些小斑点噪声，通过构造几何特征滤波器（Geometric feature filter，GFF）进行筛选去噪以得到完整的河流流域图。本文提出的算法流程图如图2所示。

图2 多特征融合算法流程图Fig.2 Flow chart of multi-feature fusion algorithm

为了将提取过程描述的更清晰，图3 为每个过程的具体示意图。（a）是输入的第四通道红外波段图；（b）是通过计算灰度共生矩阵得到的角二阶矩特征图；（c）是OTSU 算法进行双阈值分割得到的二值图像；（d）是经过连通域框选分组后的结果图；（e）为通过几何特征滤除后的得到的河流效果图。

图3 算法过程示意图Fig.3 Illustration of the algorithm process

2.1 纹理特征获取

基于输入的红外波段图，按照以下步骤进行灰度共生矩阵的计算。

（1）灰度级量化。其函数表达式为：

（2）特征值计算。以目标像素点为中心，向邻域的0°、45°、90°、135°四个方向都进行纹理特征计算，最后取得4 个特征值的均值作为最终的结果以消除旋转噪声因子。通过计算可得到以下矩阵：

（3）特征选取。计算得到的灰度共生矩阵有多个统计量，包括均值、方差、对比度、角二阶矩等。其中选取的角二阶矩谱计算公式为：

图4为各统计量结果示意图。

图4 灰度共生矩阵各统计量示意图Fig.4 Illustration of GLCM’s statistics

2.2 光谱特征获取

基于得到的角二阶矩图，用OTSU 算法基于下式进行图像的光谱特征分析。

式中：σ是类间方差；w0前景像素占图像总像素的比例；u0是前景平均灰度；w1背景像素占图像总像素的比例；u1是背景平均灰度。

2.3 几何特征获取

基于得到的二值图像，先进行连通域标记［11］。每一个单独的连通域为一个标识的目标块。进一步可以获取这些标识块的质心、内接圆、外接矩形等几何参数。示意图如图3（d）所示，每个红色矩形框选的区域为一个连通域。

河流在遥感图像中呈现线性细条状，且具有一定的宽度，根据此特征构建几何特征滤波器，滤除图像中的非河流区域噪声。几何特征滤波器包括连通域长度、连通域矩形度、连通域面积与其外接矩形面积比。

（1）连通域长度l。与其他噪声相比，河流区域的l值较大。因此通过设定滤波器阈值Tl以滤除小块非河流区域噪声。

（2）连通域矩形度η。用来描述河流的矩形特征，设第i个连通域的矩形度为ηi

式中：Wi，Ll分别为第i个连通域的宽和长；η的取值范围为（0，1］，越接近于0，矩形度越高。认定ηi≤Tη的区域为河流区域。

（3）连通域面积与外接连通矩形标记框面积比β。河流具有弯曲细长的特征，连通域矩形框选某一连通域河流时，面积会远远大于该矩形框中河流的面积。第i个连通域面积与其外接连通矩形标记框面积比为：

式中：εi为连通域面积；Si为连通域标记框面积，可由Wi和Li相乘得到。

通过设定Tβ这一阈值来滤除掉不符合该特征的区域。示意图如图5所示。

图5 连通域与外界矩形面积示意图Fig.5 Examples of external rectangle area

基于上述特征定义几何特征滤波器g={l，η，β}，对连通域标记后的图像进行筛选滤除，以得到无噪声干扰的河流流域图。

3 提取结果

基于高分一号卫星获取的台山市遥感影像，该影像共有4个波段，空间分辨率为2 m，运用本文设计的方法进行河流区域提取并与其他算法进行性能比较和评估。

3.1 阈值选取

在几何特征分析中涉及到的阈值选取对本算法的提取效果至关重要。这里对几何特征滤波器各参数进行讨论，寻求最优阈值，以获取最好的提取效果。

河流在连通域中体现为细条状，由图6（c）～（f）可知，合理设置Tl数值可以有效滤除非河流的小斑点噪声。考虑到河流受到横穿公路和桥梁以及周围面积较大湖泊的影响，因此设置Tl=100，即连通域长度小于100的将会被当作噪声滤除。矩形度会滤除因仅考虑连通域长度而忽略掉的方形水田、建筑物等地物。当ηi≥0.6 时能有效滤除方形噪声且保持河流完整流域，如图6（g）～（j）所示。连通域面积与外接连通矩形标记框面积比是刻画连通域在矩形框中的填充度，可以有效滤除弯曲度较小的噪声，河流连通域在其矩形标记框内的所占比一般是一个较小的数值，βi≥0.35时可以有效滤除山体、矩形水田等地物的干扰，如图6（k）～（n）所示。

图6 几何特征滤波器参数分析Fig.6 Parameter analysis of Geometric Feature Filter

3.2 处理结果分析

为了突出本文算法的特点，选取了以下常见方法进行提取结果对比。包括在遥感影像分割中常用的OTSU 算法、近年来有学者提出的SWT 算法以及全卷积神经网络分割方法中的UNet 算法［12］。为了验证本文算法的普遍适用性，选取了流经不同地貌的河流河段。截取的影像大小都为1 500×1 500像素，其中图1是流经林地的河段，河流上方有公路横穿，公路与河流交接地带可能对提取产生干扰。图2 是流经城市区域的河段，此处的楼房依河而建，离河流较近的楼房对河流提取会有一定的干扰。图3 是流经山地区域的河段，山体阴影对河流提取存在影响。图4 是流经农田区域的河段，水田以及灌溉沟渠对河流提取的影响较大。分别使用的OTSU+GFF、SWT+GFF、U-Net以及本文方法提取对比图分别如表1的第2～6行所示。

表1 算法提取河流效果对比图Tab.1 Algorithm to extract river effect comparison

OTSU 结合GFF 算法能够提取到完整的河流形状，但也包含大量的非河流区域噪声。将其与几何滤波器结合，能够滤除绝大数的非河流区域。但河流在林地区域的公路、城市区域的楼房、山地区域的山体阴影处产生了粘连噪声，是因为该方法仅考虑了图像的灰度特征，在有与其灰度相当地物的干扰时，很难将其分开。

SWT 结合GFF 算法是依据图像中平行的边缘线进行河流提取的，当截取的河岸一边有丢失时，会导致该部分的河流无法提取，如表1 第三行第一列流经林地区域的提取图所示。此外该算法主要依据的是图像的边缘信息，在依河而建的城市区域以及农田区域，楼房水田、水渠等与河流距离很近，容易产生边界不明显的现象，因此很容易出现噪声。

U-Net 网络使用的训练数据集是Sentinel II 卫星拍摄的位于不同气候带的河流影像，包含2982 对400×400 像素的彩色遥感图像和相应的河流标签，图像的空间分辨率为10 m。初始学习率设为0.000 01，动量设为0.9，衰减系数设为0.000 000 01，梯度下降法的最小批次设为2，设置训练迭代次数分别为10、50。输入上述影像对应的彩色图像用得到的训练模型进行分割，结果显示当epoch=10 的模型对流经林地、城市区域的河流分割效果较好，河流提取较为完整，但存在少数粘连噪声以及少量河岸缺失的现象。耕地区域的提取结果中有大量非河流噪声，而在山地区域提取失败。当epoch=50时，流经林地、城市区域的河流出现了过拟合现象，与河流相似地物如湖泊被当作河流提取出来，而此时流经山地、耕地区域的河流得到了近乎完美的提取效果。

本文算法综合兼顾了遥感影像中河流的纹理、灰度、几何特征，能够最大程度的保持河流流域的完整性并能有效抑制其它算法出现的噪声粘连现象。

3.3 结果定量评价

在提取到的河流流域图中自动生成10 000个随机点，设N0为河流像元个数，Nr为正确提取的像元个数，Ne为错误提取的像元个数，Nl为漏提取像元个数。基于以上4个参数，得到两个定量评价公式。

（1）准确率。

（2）完整度。

根据上述公式对表1中每一幅图像分别计算其准确率和完整度，基于随机点数的不确定性，进行多次实验取平均值，结果如表2所示。

表2 定量分析结果统计Tab.2 Statistics of quantitative analysis results

从定量分析结果中可以看出，将本文算法应用到四种地形时，无论从准确率还是完整度方面都有明显的优势。在山地与耕地区域得到了近乎完美的河流流域，在流经林地和城市区域由于公路和桥梁的影响，导致该方法在断连处有很小的局部河流信息丢失。OTSU 算法受到灰度近似的桥梁、公路和建筑物的影响，会出现几个特征滤波器无法滤除的粘连噪声，因此在林地、城市和山地区域的准确率不高。而在耕地区域，OTSU 算法有着不错的效果，精度和完整度都可达到99%以上。SWT算法受到岸边水田、灌溉沟渠和楼房的影响，在进行笔画宽度变换时，与河道构建成一个连通域，几何特征滤波器对该部分噪声失效，所以该算法在林地、城市和耕地区域的准确率不高。其中SWT算法在林地区域的完整度最低，是因为所给图像缺失了河岸的一边，而该算法非常依赖边缘的提取，导致该河段有较多信息丢失。U-Net算法的输出效果与迭代训练次数有着密切的联系，其中当epoch=10时，流经林地、城市的河流都有着不错的提取效果，但是在山地和耕地区域的河流提取误差较大甚至提取失败。当加大迭代训练次数后，山地和耕地区域河流的精确度得到极大的提升，准确率和完整度达到了98%以上。而此时流经林地和城市区域的河流提取结果中包含了大量湖泊噪声，致使河流提取的准确率下降。

4 结论

高分辨率遥感卫星图像中河流有多种特征，本文综合考虑河流包括纹理、光谱以及几何3 个方面特征进行验证分析。选取流经不同地貌区域的河段实验并进行结果的定性和定量分析，结果表明，本文方法是一种适应性很强的方法，能够在不同地貌环境下准确提取河流。通过与对照组的对比分析，该方法能最大程度减少灰度近似地物以及河岸缺失的影响，有着较好的性能表现。

但是该算法还存在一定的缺陷：一是由于灰度共生矩阵的计算量较大，因此在实验中的计算速度会比较慢。如何进一步进行图像预处理缩减灰度共生矩阵的计算量是下一步要解决的关键问题。二是在几何特征滤波器的参数阈值设定需要人为的干涉与判断，怎样自适应确定滤波器的参数是下一步工作的重点。