水产养殖池是利用天然水面或人造池塘,在人为控制下繁殖、培育和收获水生动植物生产活动的场所^[1]。《2022中国渔业统计年鉴》显示,2021年我国内陆养殖面积为4.98万km²,养殖产量为3 183.27万t,其中池塘养殖产量占内陆养殖产量的52.26%。池塘养殖已成为中国水产养殖的主要形式和水产品供应的主要来源,在保障优质动物蛋白供给、促进农业增效和农民增收、扩大农村就业等方面发挥着重要的作用^[2]。然而,水产养殖池的快速扩张也对内陆湿地生态系统造成了显著影响,如造成生物多样性下降、水体污染、土壤酸化和农田盐碱化等。因此,迫切需要精准的、及时的内陆水产养殖池时空动态变化监测信息,以实现内陆水产养殖的可持续发展与科学管理。

近年来,面向遥感大数据的云存储和云计算技术在过去几年间得到了迅速发展,特别是谷歌地球引擎(Google Earth Engine,GEE)云计算平台,大大促进了各类遥感研究的发展。应用云平台进行遥感大数据的处理和分析,已经成为养殖池分布信息快速提取的关键技术途径^[3]。受到数据可获取性、覆盖能力等因素的影响,Landsat和Sentinel-1/2是目前研究中使用最多的数据源^[1]。其中,研究者们基于Landsat的时间跨度大、覆盖范围广和可免费获取等优势,在中国沿海、越南等地区做了各种面向对象的分类研究^[3⇓-5]。Sentinel-1具有可以穿透云层、植被和土壤等优点,在绘制水产养殖分布时能够更为精细地检测到养殖池对象的线性结构细节。因此,越来越多的研究者以Sentinel-1作为数据源在越南、印度、东南亚等地区做了不同尺度的研究^[6⇓⇓-9]。和Landsat相比,Sentinel-2能够提供更具优势的影像(10 m空间分辨率,2~6 d重访周期),能够有效区分水产养殖池的范围,监测其年内的变化; 和Sentinel-1相比,Sentinel-2具有更丰富、有用的影像信息与信噪比,能够增强水产养殖池及其周围水坝的信息。综上所述,Sentinel-2在湿地分类制图上具有巨大潜力。

在内陆地区,养殖池与其他地表水体具有相似的光谱特征,且在水色、形状、面积等方面表现出较强的空间异质性,导致在内陆区域对其进行准确的遥感制图与监测十分困难。目前,已有众多的监督分类方法,如决策树、随机森林、深度学习等,在绘制不同国家或地区养殖池的研究中^[3⇓-5]被广泛使用。然而,养殖池与其他地表水体具有极其相似的光谱和纹理特征,因此仅仅依靠监督分类的方法无法有效区分养殖池与其他地表水体,需要复杂的数据处理才能获得准确的水产养殖池产品。一些研究在监督分类的基础上进一步应用了基于知识的决策树,用于区分其他地表水体的干扰^[6]。这些研究通常只使用了基于单幅影像获取的空间特征,忽略了时间变化特征,难以区分养殖池与其他在大小、形状、特征上极为相似的地表水体及其他地物,如水库、坑塘、水稻田等。此外,以往的研究主要针对滨海的养殖池,而内陆的养殖池相比于滨海区域,其土地覆被构成更为复杂,干扰更多(如水田等),提取难度更大。因此,需要一种适用于内陆的,能够考虑到时空变化特征,准确区分养殖池与其他地表水体及地物的方法。

针对上述问题,本研究将以江汉平原为研究区,提出了一种基于GEE与时间序列Sentinel-2数据的江汉平原养殖池智能提取的方法,并进行变化监测,实现了2016—2022年逐年的江汉平原养殖池精确提取及时空格局分析。

江汉平原(29°26'~31°37'N,111°14'~114°36'E)位于湖北省的中南部,是中国海拔最低的平原之一^[10]。江汉平原属于亚热带季风气候,河网稠密、湖泊众多,水域面积广大,水域面积占其总面积的18%,其中湖泊面积达1 605.4 km²。《2022中国渔业统计年鉴》显示,2021年湖北省水产养殖池面积为5 163.42 km²,常年位居全省第一。鉴于该平原的城市与其养殖池的分布情况,本文提取湖北省境内的长江、汉江流域周围50 m等高线所包罗的范围,以此为基础再外扩至该范围所属的县级行政区,包括22个区、15个县与县级市、4个省直管市,总面积约57 286.6 km²(如图1所示)。由于地理和养殖方式的差异,江汉平原养殖池通常具有不同的形状和大小。

研究基于GEE平台获取江汉平原2016—2022年的Sentinel-2数据,共计3 058景影像(表1)。为了减轻云污染的影响,在合成无云影像时对云量百分比进行了限制(<20%)。然后使用Sentinel-2云掩模算法对指定时间和空间范围内影像进行计算,以中值合成方法重构最小云量合成高质量影像^[11],从而生成了2016—2022年每年的高质量合成影像。与原始影像相比,使用中值影像不仅减小了数据集的大小,还可以获取给定时间段内较为稳定的像元,从而合成更好质量的影像,这使得分析更加容易和快速。此外,本研究还使用上述方法,依次生成了每个月的高质量合成影像,来创建2016—2022年的月度时间序列图像。

为了增强研究区内水产养殖池及其周围水坝的信息,本文从Sentinel-2影像中提取了4个光谱指数: 归一化差异植被指数(normalized difference vegetation index,NDVI)^[12]、归一化差异水指数(normalized difference water index,NDWI)^[13]、修正的归一化差异水指数(modified normalized vegetation index,mNDWI)^[14]和归一化建筑物指数(normalized difference built-up index,NDBI)^[15]。其中,NDWI和mNDWI可以大大增强养殖池的信息,而NDVI和NDBI可以增强养殖池周围水坝的特征,这些水坝通常被植被、砾石或不透水表面覆盖。指数定义公式分别为:

式中: ρ_nir为近红外波段(B8,中心波长842 nm)的反射率; ρ_red为红波段(B4,中心波长665 nm)的反射率; ρ_swir为短波红外波段(B11,中心波长1 610 nm)的反射率; ρ_green为绿波段(B3,中心波长560 nm)的反射率。

本文通过实地调查的方式,结合谷歌地球高分辨率图像,利用目视解译的方法获取了2022年共计1 000个样本数据,然后基于2016—2022年的图像对采集的样本数据进行了检查和修改,最终获得了2016—2022年的样本数据。在本研究中,影响养殖池正确分类的因素主要有一些与其具有相似的光谱特征,且在水色、形状、面积等方面表现出较强的空间异质性的其他水体的干扰。因此,本研究将研究区地表水体分为几个类别,包括水产养殖池、水田和其他没有明显水坝、形状不规则的地表水体,如河流和水库等(表2)。这些地表水参考样本是通过一系列实地调查所获取的,最终样本数据进一步分为养殖池和非养殖池,样本数量如表3所示。此外,本研究还基于最新的谷歌地球手动绘制了研究区内300个养殖池面状样本数据,这些样本数据最终用于确定养殖池分割尺度和生成分类规则。

为了提高分类精度,本文引入数字高程模型( digital elevation model, DEM )、 高分辨率谷歌地球影像数据和3个常用的全球10 m空间分辨率土地覆盖产品, 包括Dynamic World ( https://dynamicworld.app ), ESA_Worldcover ( https://esa-worldcover.org/en ) 和FROM-GLC10( http://data.ess.tsinghua.edu.cn ), 作为辅助数据来评估本文提取结果的准确性。

本文技术流程见图2。在GEE云计算平台的支持下,基于Sentinel-2时间序列影像,结合K均值聚类算法(K-means)与层次决策树分类算法,本研究开发了一种面向对象的江汉平原养殖池智能提取与变化监测的方法,构建了10 m空间分辨率的江汉平原养殖池数据集(2016—2022年)。首先基于GEE获得所有可使用的时间序列影像,之后基于K均值聚类算法(K-means)对影像进行分割,再提取对象的光谱、纹理、地形和时间等特征构造特征空间,后计算变化指数,最后计算混淆矩阵并验证其精度结果。

基于GEE云计算平台,本研究采用K-means算法^[16]对研究区影像的光谱波段以及特征波段进行图像聚类分割,为了增强其聚类结果的精度,使用了1.2.1节中的4种光谱指数来增强研究区养殖池及其周围的水坝信息。K-means算法是经典的无监督学习算法,具有时间复杂度低、算法伸缩性强等优点,在图像分割领域有非常广泛的应用^[15,17],K-means算法的关键参数是“K”,它与目标地物类别个数有关。不同的中心类簇个数会带来不同的聚类分割结果。本文主要采用实地调查的方式,结合高空间分辨率谷歌地球影像和GEE,通过不同的聚类分析来确定研究区主要土地覆被类型及种类。如图3所示,分别分析了研究区3种典型养殖池,在不同的K值变化下的分割结果。发现在天门市的养殖池,当K=9时分割效果明显,K取到最优解,当K=10时,数据集中的数据对象所属类簇不再发生显著变化;在荆州市的养殖池,当K=8时分割效果明显,当K大于8(取到9)时分割效果变化不大;在仙桃市的养殖池,当K=9时分割效果明显,K取到最优解。因此,最终本研究选用9作为初始类簇中心个数。

本研究基于层次决策树算法进行江汉平原养殖池的提取。决策树是一种用于分类和回归任务的非参数监督学习算法。它是一种分层树形结构,由根节点、分支、内部节点和叶节点组成^[18]。决策树能够直接体现数据的特点,具有运行速度快、准确性高的优点,可以在相对短的时间内对大型数据源做出可行且效果良好的结果。本文结合研究区影像的光谱指标、形状、地形(坡度、高程)、纹理等条件设计层次分类树来提取江汉平原水产养殖池。决策树的详细描述如下:

1)利用NDWI和OTSU算法^[19]将处理生成的景观对象分为水对象和非水对象2个部分。OTSU算法又名最大类间方差法,是一种常用且高效的算法,在土地覆盖分类中得到了广泛应用^[20⇓-22]。其原理是根据影像的灰度特性,将图像分割成背景和前景2类,选择使其类间方差最大、类内方差最小的分割阈值为最优阈值,进行图像的自动二值分类^[21]。

2)分析水体对象的NDWI时间序列特征,排除水田等季节性水体。虽然水田在大小、形状和水体特征上与养殖池十分相似,但是在一年中的某些时候,水田具有明显的季节性特征,定期有水面覆盖。基于此特性,分析了江汉平原水体对象的NDWI时间序列特征,如图4所示,可以很明显地看出水田与养殖池的区别,最后在结果中剔除水田等一些季节性水体。

3)根据水体对象的形状特征,进一步提取江汉平原水产养殖场对象。水产养殖池通常位于地势低洼的地区,并有伴有规则形状的水坝。这些独特的地形和形状特征可以有效地将其与其他永久水体区分开。因此计算了1.2.2节构建的300个养殖池面状样本数据的地形和形状特征,包括高程(elevation,E)、坡度(slope,S)、面积(area,A)、紧凑度(compactness,C)、景观形状指数(landscape shape index,LSI)和矩形度(rectangularity,R)(表4),确定最佳分割尺度,生成分类规则集并构建了决策树分类集指标。表中: P为周长; A为面积; A_mbr为外部矩形的面积。利用这些特征指标有效地将养殖池与其他水体的区别开来。其中,计算得到的最佳分割尺度阈值分别为坡度(T1=0.82)、紧凑度(T2=0.58)、高程(T3=38.2)、LSI(T4=1.3)、矩形度(T5=0.82)和面积(T6=0.01)。层次决策树分类模型见图5,养殖池与其他水体的形状特征分析见图6。

4)从上一步的结果中剔除水库。水库和水产养殖池具有相似的大小、形状和长期特征,水库与水产养殖池的根本区别在于:水库通常是孤立的,而水产养殖池则呈聚集分布。基于此特性,本文根据其空间关系的孤立性来识别水库。具体步骤如下:首先,对每个对象按其边界向外延伸50 m,确定其与周围物体的空间关系;然后判断其周围是否有物体与它相交,如果有,就确定它是一个水库。根据判断结果进一步对水库进行识别和剔除,最终得到准确的江汉平原水产养殖池对象。

本研究基于1.2.2节获得的2016—2022年所有样本数据对本文提取结果进行准确性评估。评估指标包括总体准确性(overall accuracy,OA)、Kappa系数精度、用户精度(user’s accuracy,UA)和生产者精度(producer accuracy,PA),得到的混淆矩阵及指标如表5所示。结果表明K-means聚类算法得出的分类结果与地面验证样本高度一致。2016—2022年,每年分类的OA高于91%。江汉平原养殖池的PA和UA均大于90%。为了进一步评估K-means聚类算法的详细形状和轮廓的准确性,将本研究中生成的江汉平原水产养殖池与基于0.54 m空间分辨率的谷歌地球图像的养殖池面状样本数据进行了比较,如图7所示,发现江汉平原养殖池的轮廓和形状与通过人工数字化获得的高精度结果非常吻合。图7(e)—(g)所示,发现了在面积、LSI和紧凑度方面,本研究所提取的结果与实际结果具有很强的相关性,相关系数均大于0.912。

研究区水产养殖池分布广泛,总面积为2 126.43 km²。主要集中在荆州市、仙桃市、武汉市、孝感市、潜江市等5市(图8)。这些地区的水产养殖池占整个江汉平原研究区总数的93.87%以上。如图8(b)和(c)所示,荆州市水产养殖池面积最大(919.48 km²,43.24%),其次是仙桃市(316.29 km²,14.87%); 宜昌市水产养殖池面积最小(16.13 km²,0.76%),其次是荆门市(35.74 km²,1.68%)。

江汉平原养殖池在2016—2022年间的动态变化呈现出明显的空间异质性,整体呈现增加的趋势,总面积从1 947.43 km²增加到2 126.43 km²,增长率为9.19%(179.00 km²)(图9和表6)。其中,养殖池面积增加最显著的地方是荆门市,增长率为41.48%(平均每年增加约1.50 km²),其次是荆州市(38.89%,36.78 km²)、宜昌市(28.13%,0.51 km²)、孝感市(15.23%,5.37 km²)。相比之下仙桃市、天门市、武汉市和潜江市养殖池面积出现了减少的趋势,平均每年损失分别约8.90 km²,1.90 km²,5.28 km²和2.50 km²,减少率分别为16.46%,14.47%,11.41%和8.47%。

本文基于时间序列Sentinel-2密集时间序列影像和GEE云计算平台,提出了一种结合K均值聚类(K-means)与层次决策树分类算法的养殖池提取与变化监测方法。与传统的基于机器学习的分类方法相比,本文的方法充分考虑了水产养殖池的空间异质性,结合光谱、纹理、地形、形状等特征可以很好地识别养殖池与其他干扰的水体;与其他传统的基于像素的方法相比,本文使用的是一种面向对象的方法,在提取养殖池时,可以利用到更多的形状、大小等特征信息,可以对养殖池进行准确的分割,提高了养殖池分类精度。此外,本研究还结合最佳尺度分割分析和时间特征,有效区分了养殖池与其他容易混淆的水体,最终获得了江汉平原地区空间分辨率为10 m的水产养殖池数据集(2016—2022年),总体分类精度达到91.90%以上。因此,本研究的方法可以很容易地适用于在全国开展养殖池提取与动态监测,所得数据集对支持江汉平原地区生态保护和可持续发展目标的评估具有重要的借鉴价值和现实意义。

如表7所示,本研究数据集与3种常用的全球10 m分辨率土地覆盖产品(Dynamic World,ESA_Worldcover,FROM-GLC10)进行比较。为确保可比性,本文选择了2016—2022年间变化较小的区域进行比较。3种数据集都是采用基于像素的监督分类算法,该方法依赖于大量的训练样本,并且容易受到地物复杂性的影响,产生“椒盐现象”。这种方法很难准确识别养殖池周围的水坝和狭窄的排水渠。相比之下,本文数据集是使用一种基于面向对象的智能提取的方法,该方法结合了K-means和层次决策树分类算法,从而避免了“椒盐”问题的影响,并能够准确检测大多数水产养殖池周围的水坝等信息。不足的是:在获取样本数据时,仅采集了不同形状、大小的养殖池样本数据,未考虑到不同功能的养殖池,未做到更细致的分类,在以后的研究中,划分样本数据需要考虑到这一因素; 在提取养殖池时,Sentinel-2影像容易受到研究区的云雾干扰,导致在识别养殖池时造成数据信息的缺失,而雷达影像作为主动式成像,能够很好地克服这一缺陷^[23]。因此,在今后的研究中,将联合Sentinel-1雷达影像与Sentinel-2多光谱影像作为数据源,依靠优势互补提高水产养殖池制图精度,实现更高频次、长时间特征的变化检测^[24]。

本研究基于GEE与Sentinel-2密集时间序列影像,提出了一种结合K均值聚类(K-means)与层次决策树分类算法的养殖池提取与变化监测方法,实现了2016—2022年逐年的江汉平原养殖池精确提取及其时空格局分析。主要结论有:

1)结合K-means与融入时间变化特征的层次决策树算法能够实现精准的养殖池分类,2016—2022年的总体分类精度均达到91.90%以上,Kappa系数达到0.84以上。

2)2022年江汉平原的水产养殖池面积为2 126.43 km²。其中,43.24%的养殖池集中分布于荆州市,宜昌市养殖池面积最小,仅占0.76%。

3)江汉平原养殖池在2016—2022年期间的动态变化呈现出明显的空间异质性,整体呈现增加的趋势,总面积从1 947.43 km²增加到2 126.43 km²,增加了9.19%。其中荆门市和荆州市的养殖池面积增加最为明显,分别增加了41.48%和38.89%;而仙桃市和天门市的面积出现了减少的趋势,分别减少了16.46%和14.47%。

本研究的方法为养殖池的精准提取研究提供了参考,所得数据集对支持江汉平原地区生态保护提供依据,对支持可持续发展目标的评估具有重要的借鉴价值和现实意义。