中国是淡水养殖大国, 其淡水渔业产量约占世界淡水渔业产量的60%^[1]。养殖水体面积是衡量水产养殖能力的一个重要指标, 传统的人工养殖水体信息调查存在工作量大、数据更新较慢等缺点, 难以准确、客观地反映养殖现状。利用遥感技术可实时、快速、准确地获取养殖水体的时空分布及动态变化, 有助于水产养殖产业的科学规划与管理, 对快速提升渔业信息化和科学化管理具有重要的意义^[2]。

目前, 常用于提取水体信息的遥感数据主要有MODIS、TM和ETM⁺、Sentinel系列卫星以及国产的HJ-1A/B等。Rao等^[3]基于MODIS数据(MOD09A1和MOD13Q1)并采用随机森林模型, 对2000—2016年长江流域进行水体信息提取, 并分析其变化趋势, 其结果对水资源利用和洪水监测具有非常重要的意义。Pekel等^[4]采用300万景TM和ETM⁺遥感影像, 研究了1984—2015年全球陆表水体变化, 为水资源利用管理提供了参考。Xing等^[5]采用Sentinel-1卫星数据对2016年洞庭湖水面面积进行月度变化监测, 结果表明洞庭湖水面面积变化明显, 7月水面面积是12月的近2.8倍。Lu等^[6]采用2009年5月HJ-1A/B影像, 结合水体指数法对北京及周边地区进行水体提取, 取得了较好的效果。但是以上研究中采用的遥感影像空间分辨率均不足, 在提取面积较大的湖泊和水库等水体时满足要求, 但在提取面积较小的养殖水体时存在较大误差。随着高分辨率遥感技术的进步, 可以获得的高分辨率遥感数据日益丰富。许多研究者正在开展基于高分辨率遥感数据的水体提取研究, 如SPOT、WorldView-2以及国产ZY-3和GF-1卫星数据等。陈文倩等^[7]以2013年12月GF-1卫星数据为基础, 提出了单波段阈值法与构建阴影水体指数相结合的决策树水体信息提取方法, 对新疆特克斯河上游支流进行水体提取, 精度得到较大的提高。Yao等^[8]基于2012—2013年ZY-3卫星数据, 提出了一种城市水体提取方法, 在青岛、阿克苏、福州等城市的实验结果表明, 该方法相比NDWI方法精度明显提高。这些高分辨率遥感数据在一定程度上克服了混合像元对水体信息提取精度的影响, 已经成为水产养殖水体遥感动态监测的重要数据源。

国内外学者曾提出很多利用遥感技术提取水体信息的方法, 如单波段阈值法、多波段谱间关系法、多波段组合指数法、遥感影像分类方法以及深度学习技术等^[9-11]。单波段阈值法对遥感影像中单个波段设置阈值, 运算速度较快, 但是阈值大小直接影响水体提取精度, 由于不同区域水体性质(如含沙量、覆冰状况等)差别较大, 该方法存在一定的局限性。多波段谱间关系法利用遥感数据的多波段光谱信息, 通过建立波段关系函数来判别水体和非水体, 其水体提取精度比单波段阈值法高^[12]。多波段组合指数法通过构建波段运算函数来判别水体和非水体。其中最有代表性的方法为归一化差异水体指数(normalized difference water index, NDWI)。Xu^[13]提出MNDWI (modified NDWI)方法, 利用中红外代替NDWI中的近红外波段, 解决了NDWI在抑制建筑物信息方面的不足。王小标等^[10]提出一种MBWI (multi-band water index)方法, 相比NDWI具有更高的精度。但是对于IKONOS、QuickBird、GF-1等高分辨率遥感影像, 由于都缺少中红外波段, 也限制了这些方法的应用。现阶段, 许多研究者提出把遥感影像分类方法应用到水体信息的提取中^[8-9], 如随机森林(random forest, RF)、支持向量机(support vector machine, SVM)以及BP神经网络(back propagation neural network, BP)模型等。这类方法的优点是不需要人为设定分类阈值, 消除了阈值设定对水体信息提取的影响, 而且也更容易结合遥感影像的光谱特征和纹理特征, 只需要少量样本进行模型训练即可实现对水体高精度提取。深度学习技术在计算机视觉领域的目标检测和分类应用中表现出强大的性能, 其在遥感影像分类和目标检测中也得到广泛的应用^{[11, 14]}, 但是该方法需要大量已标注样本进行模型训练, 而大量已标注的样本较难获得, 该类模型的水体提取精度还有待验证。

基于上述可行性分析, 本研究选择安徽省巢湖市为实验区, 采用国产GF-1和ZY-3高分辨率遥感影像为数据源, 利用NDWI和纹理特征, 发展一种集成随机森林、支持向量机和BP神经网络的分类模型, 再利用谱间关系法进行阴影检测和剔除, 以数学形态学滤波进行分类后处理, 来实现水产养殖水体信息的提取与动态监测应用, 为科学合理开发和利用水产养殖用地提供依据。

1 实验区和遥感数据源 1.1 实验区概况

实验区选择安徽省巢湖市(图 1), 经纬度范围为31.27°~32.00°N, 117.42°~117.96°E。区域内水体类型丰富多样, 包含了中国五大淡水湖之一的巢湖部分水域。沿湖河流众多, 其中较大的河流有杭埠河、丰乐河、白石天河、派河、南淝河等, 且分布不均, 除此还有大量的水库、坑塘水面等, 是长江下游的重要生态湿地。另外实验区内地物覆盖类型也十分多样, 不仅包含大量水体, 同时包含城区以及沈家山、银屏山等山区, 存在大量城市建筑阴影和山区阴影对比区。利用遥感技术获得该区域内水产养殖水体的面积, 可以全面掌握巢湖市的养殖规模及其动态变化, 监测结果对指导当地渔业生产管理、合理规划产业布局具有现实指导作用。

1.2 遥感数据源

本研究以GF-1卫星数据和ZY-3卫星数据为数据源。GF-1是中国高分辨率对地观测系统的首颗卫星, 选择全色波段(2 m)和多光谱波段(8 m)数据进行融合, 获得空间分辨率为2 m的多光谱影像。ZY-3卫星是中国首颗高精度民用立体测绘卫星, 选择全色波段(2.1 m)和多光谱数据(5.8 m)进行融合, 获得空间分辨率为2.1 m的多光谱影像。

研究选取2016年1月26日2景、2016年4月21日1景GF-1影像数据, 经过拼接得到实验区。选取2013年3月28日2景ZY-3影像数据, 拼接得到实验区。此外, 在水体动态变化监测应用时, 为了保证GF-1和ZY-3的空间分辨率一致, 将ZY-3融合后数据重采样为2 m空间分辨率。

2 研究方法

本研究中的淡水养殖水体主要是指规则型池塘、山塘水库以及大水面等封闭型水体, 不涉及河流、进排水渠等流动性水体及公园内观赏水体^[1]。

由于NDWI能够较好地反映水体的光谱特征; 同时水体的纹理比非水体的纹理更加细腻, 也可用于判别水体和非水体; 加之集成分类模型是由两个或两个以上分类器组合而成, 可以有效地对单个分类器的分类信息进行互补, 提高分类精度, 因此基于上述思路, 本研究设计了如下水产养殖水体信息提取及其动态变化监测的技术流程, 如图 2所示。

2.1 影像预处理

影像的预处理选用ENVI 5.3中针对GF-1和ZY-3数据的处理模块, 具体内容包括辐射校正、几何校正、正射校正、影像融合、影像拼接和裁剪等。辐射校正包括辐射定标和大气校正。几何校正选取遥感影像与具有准确地理坐标卫星影像间的同名点, 进行重采样纠正几何偏差。正射校正对影像同时进行倾斜改正和投影差改正, 将影像重采样成正射影像。影像融合处理高分辨率全色影像和低分辨率多光谱影像, 获得高分辨率多光谱影像。在此基础上, 通过影像拼接和裁剪获得实验区影像。

本研究精度验证方法, 是在Arcgis10.3中通过人工目视解译, 把实验区内的地物分为水体和非水体两种类型, 其结果作为模型的训练样本标签与精度验证标准。

2.2 NDWI和纹理特征生成

本研究采用NDWI和纹理特征组成的特征向量作为分类依据。在水体中可见光范围内反射率最高的是蓝绿波段, 最低的是近红外波段, 依据这一特性建立的NDWI可以最大程度地突出水体信息, 抑制非水体信息, 表达式如下^[15]:

$ \mathrm{NDWI}=\frac{B_{2}-B_{4}}{B_{2}+B_{4}} $

(1)

式中, B₂为绿波段, B₄为近红外波段, 分别对应GF-1卫星影像第2波段和第4波段。

NDWI仅考虑单个像元, 而纹理特征可以表达其邻近像元之间的相互关系。本研究采用ENVI 5.3中二阶概率统计纹理分析工具, 选择均值(mean)、方差(variance)两种纹理特征, 计算窗口选择5×5。由于多光谱遥感影像波段数目多, 且各波段之间具有较强的相关性, 因此本研究在提取遥感影像纹理特征时, 首先对影像进行主成分分析, 并选取第一主成分进行影像纹理分析。

2.3 集成分类模型构建

本试验以MATLAB为平台, 编程实现随机森林、支持向量机、BP神经网络以及集成分类模型4种方法的运算。

随机森林算法是Breiman^[16]在2001年提出的以多个决策树为基础的集成式分类器。该算法通过从原始训练样本集中有放回地重复随机抽取一定数量的训练样本, 生成新的训练样本集合, 每个训练样本集合对应一个决策树, 待分类结果由所有决策树采用众数投票法决定。该方法的关键就是设置合理的决策树个数, 该值越大模型的收敛性越好, 但训练时间也会显著增加, 且当树的个数过多时, 模型会出现过拟合现象。

支持向量机是20世纪90年代由Vapnik提出的基于统计学理论的机器学习方法, 采用非线性变换将输入空间转换到高维空间, 并在高维空间构建最优分类超平面, 使得不同种类样本之间能够最大化分离。该方法已经广泛应用在遥感影像分类中并取得较高的精度^[17]。

BP神经网络是一种按照误差逆传播算法训练的多层前馈人工神经网络^[18]。整个训练过程由信息的正向传播和误差的反向传播两个子过程交替进行来完成, 在此过程中根据梯度下降算法, 不断调整网络中各层神经元的权值, 当输出层结果与期望结果之间的差值很小或者达到训练设定阈值时, 模型训练结束并保存模型参数。

集成分类模型可以获得比单个分类模型更好的分类结果以及鲁棒性。本研究选择现有研究中分类精度较高的随机森林、支持向量机和BP神经网络作为子分类器, 采用众数投票法决定待分类像元的属性。3种分类方法中两种或两种以上的模型都判断为水体时, 集成模型才最终判断为水体。

2.4 水体提取信息后处理

利用集成分类模型算法能够准确、快速地提取出养殖水体信息。但是由于山区阴影、城市建筑阴影以及信号噪声的干扰, 提取出的养殖水体往往仍可能存在一定的误差, 需要对养殖水体提取结果进行后处理, 具体包括阴影检测及剔除和数学形态学处理。

在遥感影像中, 由于水体与阴影的光谱特征和纹理特征比较接近, 对水体进行提取时, 部分阴影也会被误认为是水体。现有研究发现, 阴影在ZY-3卫星各波段的光谱均值表现为B₁+B₄ > B₂+B₃, 而水体在各波段的均值特征为B₁+B₄ < B₂+B₃^[19], 考虑到GF-1和ZY-3的波段设置相同, 也可以根据此特性建立基于GF-1的增强型阴影水体指数(enhanced shadow water index, ESWI)进行阴影剔除, 表示如下:

$ \mathrm{ESWI}=\frac{B_{1}+B_{4}}{B_{2}+B_{3}} $

(2)

在水体提取结果的基础上, 通过计算ESWI进行阴影剔除, 当ESWI大于1时为阴影, ESWI小于1时为真实水体。

利用数学形态学方法可以去除养殖水体提取后产生的异质小斑块及椒盐噪声。数学形态学是一种非线性滤波方法, 运算时是将结构元素在图像的范围内平移, 同时施加交、并等基本几何运算, 以达到对二值图像的处理。膨胀和腐蚀是数学形态学的基础, 可以用集合的方式定义如下:

$ 膨胀:X \oplus S=\{x | S+x \cup x \neq \phi\} $

(3)

$ 腐蚀:X \ominus S=\{x | S+x \subseteq X\} $

(4)

式中, X为待处理图像; S为结构元素。其他形态学运算可由以上两种基本运算组合产生, 主要有开运算和闭运算。先腐蚀后膨胀为开运算, 开运算可以去除孤立的小点和毛刺; 先膨胀后腐蚀为闭运算, 闭运算可以填平小孔, 弥合小裂缝。

2.5 精度评价指标

本研究精度评价指标由分类误差混淆矩阵计算而得, 其基本统计量包括总体精度(p_t)、Kappa系数、漏分误差(p_o)和错分误差(p_c)。表达式分别如下:

$ {p_t} = \frac{{\sum\limits_{i = 1}^n {{X_{ii}}} }}{N} $

(5)

$ {\rm{Kappa}} = \frac{{N\sum\limits_{i = 1}^n {{X_{ii}} - \sum\limits_{i = 1}^n {\left( {\sum\limits_{j = 1}^n {{X_{ij}}\sum\limits_{j = 1}^n {{X_{ji}}} } } \right)} } }}{{{N^2} - \sum\limits_{i = 1}^n {\left( {\sum\limits_{j = 1}^n {{X_{ij}}\sum\limits_{j = 1}^n {{X_{ji}}} } } \right)} }} $

(6)

$ {p_{\rm{o}}} = \frac{{{N_{{\rm{wn}}}}}}{{{N_{{\rm{wt}}}}}} $

(7)

$ {p_{\rm{c}}} = \frac{{{N_{{\rm{nw}}}}}}{{{N_{{\rm{wl}}}}}} $

(8)

式中, N为样本数; X_ij为类别i被分为类别j的样本数; N_wn为水体错分为非水体的像元数; N_wt为真实水体像元个数; N_nw为非水体错分为水体的像元数; N_wl为分类得到的水体像元总数。

3 结果与分析 3.1 水体提取结果

水体提取结果如图 3所示, 在面积较大的水体上, 随机森林、支持向量机、BP神经网络以及集成分类模型水体提取结果和人工提取结果保持一致。例如, 巢湖水面, 模型提取的水体轮廓清晰, 边界完整, 和人工目视解译提取的水体结果相差不大。而对小面积水体, 如水库和坑塘, 4种方法表现出的差异性较大, 尤其在研究区右边部分, 随机森林、支持向量机和BP神经网络法均存在误分现象, 但集成分类模型则和人工提取结果更为接近。另外, 人工提取结果和模型法提取结果的最大差异主要表现在对河流的处理结果上。由于文中定义的养殖水体不包含河流等开放性水体, 而模型法提取水体是基于像素特征而得, 因此无法判断水体像元是属于养殖水体还是非养殖水体, 这是模型法提取养殖水体的最大误差源。

选取人工目视解译结果为验证数据, 分别对4种方法的水体提取结果建立混淆矩阵并进行精度评价。结果如表 1所示。集成分类模型的总体精度为97.4%, Kappa系数为0.87, 均高于其他3种模型。集成分类模型的漏分误差和错分误差分别为3.7%和6.4%, 该方法在漏分误差和其他3种方法相差不大的情况下, 错分误差得到了较大的减小, 说明集成分类模型可以很大程度地提高水体信息的提取精度。集成分类模型法的水体提取面积为423.9 km², 而人工目视解译法提取水体面积为421.3 km²。4种方法中, 集成分类模型提取的水体面积误差最小。

3.2 集成分类模型提取细节

集成分类模型提取水体结果细节如图 4所示。在水库坑塘方面, 自动提取结果和人工提取结果差别主要在水体的边界上以及干涸半干涸的水体部分。造成这种差异的原因主要有两方面: (1)人工提取水体时, 提取原则是按照水体的最大外边界提取水体面积, 在水陆交界处存在混合像元的影响, 而自动化提取水体时是按照纯水面提取水体。另外, 一般选择的训练样本为纯水面像元, 导致水陆混合像元会被分类为非水体。(2)在人工提取水体时, 基于丰富的经验判断, 把部分半干涸的库塘也归为养殖水体(图 4中红框位置), 而自动提取水体时, 由于该类水体特征较弱, 会把这类库塘归为非水体, 导致模型提取水体和人工提取水体之间的不同。

3.3 阴影剔除结果

由于水体和阴影在光谱和纹理上存在相似性, 很多分类模型在提取水体时常常把阴影也提取出来。选择实验区中存在建筑阴影和山体阴影的区域, 利用ESWI方法进行阴影检测和剔除, 其结果如图 5所示, 其中红框中为主要的阴影区。结果表明, 利用ESWI方法可以在保持水体边界的情况下去除绝大部分建筑阴影和山体阴影。

3.4 养殖水体面积动态变化监测

掌握养殖水体面积的动态变化, 是水产养殖水体资源遥感监测的重要任务之一。本研究分别采用实验区2013年ZY-3卫星影像和2016年GF-1卫星影像进行了养殖水体提取, 并对比分析两个年度间养殖水体面积的变化, 结果如图 6所示。图 6(a)中巢湖水体面积基本无变化, 而图 6d中面积较小的坑塘水面变化较大, 通过对比图 6bZY-3影像和图 6cGF-1影像, 发现其原因是巢湖面积较大, 水资源丰富, 受人类活动和季节变化的影响较小, 而较小的坑塘水面较容易受到人类活动和季节变化的影响而产生或者消失。另外, 由于监测区域的GF-1和ZY-3数据是由不同幅影像拼接而成, 影像的几何校正精度以及配准精度也会影响变化监测的结果。2013年研究区内养殖水体面积为396.7 km², 2016年养殖水体面积为423.9 km², 面积增长了6.9%。

4 讨论 4.1 模型的选择与构建

现有研究表明单一的遥感影像分类方法存在明显的不足, 如分类精度较低、错分和漏分现象严重等^[20]。多分类器集成分类模型可以实现不同分类模型的优势互补, 实现整体分类精度的提高, 而子分类器之间的差异性是集成分类模型的关键^[21]。Feng等^[9]的研究表明随机森林、支持向量机和BP神经网络是现有遥感影像分类模型中精度较高的3种模型, 其精度优于常用的最大似然分类模型, 且3种模型的分类原理具有较大的差异性。本研究采用的集成分类模型综合了随机森林、支持向量机和BP神经网络3种精度较高的模型, 3种模型的分类结果相互独立, 最后采用众数投票法获得最终分类结果。表 1中的结果表明集成模型的总体分类精度和Kappa系数均优于单个模型。

集成分类模型不仅考虑了影像的光谱信息, 也考虑了影像的纹理信息。纹理信息是对遥感影像中地物空间分布规律的一种特征表达, 是一种比光谱信息更高层的特征^[22]。水体在遥感影像中较为平滑, 区域灰度变化较小, 其纹理特征和非水体有较大差别。周亚男等^[22]和Feng等^[9]的研究结果表明加入纹理信息可以有效提高水体提取精度。同时, 纹理特征作为一种区域统计信息, 可以减弱“椒盐噪声”的影响, 具有较强的稳定性。

4.2 水体提取误差消除

现有研究表明混合像元的存在会导致小水体提取误差增加^[23], 尤其对于MODIS、Landsat等低空间分辨率遥感影像, 混合像元的影响更加显著^[3]。本研究基于GF-1高分辨率遥感影像, 可以大大减弱混合像元对水体提取精度的影响。同时, 阴影也是水体提取误差的主要来源, 这是由于阴影和水体的光谱曲线存在高度相似性^[24]。目前针对GF-1遥感影像的研究中普遍存在对阴影的剔除效果不理想的情况^[25], 而Feyisa等^[26]提出的利用Landsat影像中的5个波段数据建立的阴影剔除指数并不适合仅有4个波段的GF-1遥感影像, 本研究采用谱间关系法建立的增强阴影水体指数虽然无法完全消除阴影的影响, 但是可以剔除绝大多数城市建筑阴影和山体阴影。最后采用数学形态学滤波方法对提取结果进行分类后处理, 消除细碎的异质斑块, 从而更进一步提高了精度。

虽然本研究采用的方法利用较少样本即可取得较高的精度, 但是模型是基于像素的水体提取方法, 只能把遥感影像分类为水体和非水体, 无法对养殖水体与河流、水渠以及城市水体进行区分。另外, 深度学习技术在遥感影像水体提取方面表现出巨大的应用潜力, 现有研究表明深度神经网络模型在大量训练样本的条件下, 其分类结果优于现有的传统的浅层机器学习模型^{[14, 27]}。在人工标注样本数量有限的情况下, 综合使用传统机器学习算法和深度学习算法, 进一步构建精度更高、鲁棒性更强的水体信息提取模型, 将是今后十分重要的研究方向。

5 结论

本研究通过高分辨率遥感数据的全色波段和多光谱波段进行融合, 得到分辨率更高的多光谱影像, 在此基础上采用集成分类模型对水体进行了提取, 并利用人工目视解译结果对其进行了精度验证。实验区的研究结果表明, 采用集成分类模型提取水体, 其总体精度为97.4%, Kappa系数为0.87, 漏分误差为3.7%, 错分误差为6.4%。对监测区域的养殖水体变化监测结果表明, 2016年相比于2013年水产养殖面积增加了6.9%。研究结果误差已达到全国水产养殖水体资源遥感监测项目制订的≤±10%技术要求。目前, 该方法已成为水产养殖水体资源遥感动态监测的分析工具之一, 并在全国范围内的监测工作中加以应用。