融合Landsat8与Sentinel-2数据的红树林物候信息提取与分类

薛朝辉，钱思羽

河海大学地球科学与工程学院,南京 211100

1 引言

红树林生长在热带、亚热带区域，是海岸潮间带滩涂上的木本植物群落，在维持滨海湿地生产力、生物多样性以及保护海岸海堤、防风御浪、改善海湾环境、固碳等方面发挥着巨大的作用（范航清和王文卿，2017；杨盛昌等，2017；生农等，2021；林珲和张鸿生，2021）。然而，红树林生态系统由于沿海开发利用遭到严重破坏，因其抗干扰能力较差，导致修复困难。因此，开展红树林的合理利用与保护成为全球关注的热点问题之一。

科学准确地获取红树林时空分布信息是保护和修复红树林生态系统的基础和前提。红树林通常分布在富含淤泥的海岸潮间带，难以进行大规模野外调查和采样，而遥感技术具有观测范围广、信息量大、获取信息快、更新周期短、可比性强等优点，已成为红树林监测的重要手段（周振超等，2018；Wang 等，2019）。其中，Landsat 8 与Sentinel-2 等中高分辨率多波段光学遥感数据已被广泛应用到红树林资源探测中（Tamura 和Kikushima，2008；Giri 等，2011；Vo 等，2013）。然而，红树林遥感探测存在以下难点：红树林与陆生植被高度相似，地理空间上无自然隔离，容易造成混分、漏分现象；红树林会受到周期性潮水的浸没，导致长期受到海水淹没的红树林较难提取。针对以上难点，国内外学者已开展了大量的研究。

首先，针对红树林与陆生植被难以区分的问题，张雪红（2016）利用归一化差值湿度指数和修正的归一化差值池塘指数对广西山口国家红树林生态自然保护区的红树林进行提取，成功将红树林与陆生植物区分。为了更加充分地利用红树林与陆生植被的光谱差异信息，Baloloy 等（2020）利用Sentinel-2 数据设计了专门用于红树林提取的植被指数。由于红树林与陆生植被叶片结构存在较大区别，有研究融合光学和雷达数据提取红树林（Zhen 等，2018；Bunting 等，2018；Thomas等，2018；周振超，2019；Zhao 和Qin，2021）。此外，利用红树林为常绿植被的生物特性也能将其与陆生植被进行区分，如张威等（2015）选择Landsat 秋冬季的影像进行分类，成功将红树林和互花米草区分。进一步，有研究利用红树林与陆生植被在不同季相的差异提高了分类精度，如路春燕等（2019）利用红树林、米草以及滩涂在不同季相的光谱特征差异进行分类；Li 等（2019）利用Sentinel-2 构建归一化植被指数时间序列对福建漳江口保护区进行红树林提取。然而，以上研究仅利用红树林常绿植被的生物特性，Chen 等（2017）不仅使用了红树林的常绿特性，同时结合了潮汐淹没的地学特性，成功提取中国2015年红树林分布图。综上所述，已有研究根据红树林的生物特征，利用多源和多时相数据将红树林与陆生植被进行有效区分。

其次，针对红树林受到周期性潮水浸没难以识别的问题，贾明明（2014）通过确定高潮和低潮线建立淹没红树林指数成功提取被海水淹没的红树林。此外，也有学者考虑红树林的高度进行区分，如Xia 等（2018）利用高潮和低潮时期的遥感影像中更容易反应光谱差异的近红外波段，并结合归一化植被指数，构建了淹没红树林识别指数，分别提取了低矮和高挺的红树林；徐芳等（2020）对比分析了出露的红树林、高潮水位淹没的红树林、海水水体的光谱特征，提出了归一化潮间红树林指数，成功提取福建省龙海九龙江口湿地的红树林。综上可见，学者们通过研究和对比被海水淹没红树林的特征，较好地解决了红树林受到潮水浸没难以识别的问题。

此外，为了能更加高效地提取红树林，也有学者利用红树林的地学特性，如建立海岸线缓冲区（Long 和Giri，2011；刘大召等，2019）和筛选低海拔或坡度平缓地区（Peng 等，2020；梁浩等，2016；董迪等，2020）等方法，建立红树林可能存在的感兴趣区域来进一步提高分类精度。

以上研究虽然能有效提取红树林，但仍存在一定的不足：（1）红树林生长在热带、亚热带地区，光学传感器易受云、雨、雾等不利天气因素的影响，对红树林提取带来困难。微波传感器虽然能够弥补光学影像受天气影响的限制，但雷达数据本身受噪声影响严重，信息较弱，重访周期长，数据成本较高，难以广泛使用。（2）红树林属常绿阔叶林，提取时最为容易与其他常绿植被混淆，而以往研究仅针对红树林与陆生植被的区分进行探讨，并未对红树林和常绿植被的区分进行详细研究。已有研究表明，将Landsat 8 OLI 和Sentinel-2 MSI 进行融合重建可以有效解决遥感数据存在的时空数据缺失问题，能够为红树林监测提供更高质量的数据源（Gao和Zhang，2021）。红树林与陆生植被的物候信息存在差异，特别地，考虑物候信息能够有效区分红树林与常绿植被（Pastor-Guzman 等，2018；Berlanga-Robles 和Ruiz-Luna，2020；Mandal 等，2020），但目前尚未有利用物候特征进行红树林提取的研究。

据此，本文基于GEE 云平台，以孙德尔本斯为研究区，获取2016年全年的Landsat 8 OLI 和Sentinel-2 MSI 数据（共计395 景影像），进行红树林提取研究。首先，基于最小二乘回归构建两个传感器在相同指数之间的关系，并对Sentinel-2 MSI数据计算得到的指数进行线性变换，获得融合后的数据集，以此重建时间序列数据，之后根据J-M距离选择可分性较高的EVI和LSWI指数时间序列；其次，对两个指数的时间序列数据进行S-G滤波处理，并分别提取生长期始期等13种物候信息；最后，将两个指数的物候信息进行特征级联，采用随机森林方法进行分类，提取研究区红树林范围。

2 研究区与数据资料

2.1 研究区概况

孙德尔本斯拥有全球最大的单块潮汐嗜盐红树林，该研究区位于恒河三角洲，地处印度和孟加拉国交界，地理范围为21°30′N—22°30′N，88°15′E—89°53′E （Dutta 等，2015；Ghosh 等，2016），研究区区位如图1所示。孙德尔本斯河道纵横交错，由许多泥滩和小岛组成（Ghosh 等，2015）。该区域属于热带季风气候，降雨量为1600—2000 mm，温度为11—37 ℃。区域内植被资源丰富，通常可识别的植被类型有常绿灌木、常绿林、落叶林、耕地和草本植被等。红树林通常沿海岸呈块状分布，能够对该地区起到防风御浪的作用，但红树林与沿海植被生长较为紧密，给提取研究带来一定难度。

图1 研究区地理位置Fig.1 Location of the study area

2.2 数据资料

2.2.1 遥感数据

研究基于GEE平台获取2016年全年的Landsat 8 OLI 数据和Sentinel-2 MSI 数据，共计395 景影像，如表1所示。所获取的Landsat、Sentinel-2 数据已经过辐射校正，并且没有几何畸变的问题，无需进行辐射定标、大气校正和几何校正。根据数据集自带的质量评价信息，采用波段逻辑运算进行去云处理。

表1 本文使用的遥感数据Table 1 Remote sensing data used in this research/景

2.2.2 样本数据

为了更加客观准确地选择地面样本点，选取以下3类数据集作为参考。

（1）2016年全球红树林地图数据（Global Mangrove Watch 2016，以下简称“GMW2016”）：数据来自全球红树林监测项目，该项目主要目的是向缺乏红树林监测系统的国家提供红树林范围和变化图，为防止红树林的进一步损失和退化提供帮助（Thomas 等，2017；Bunting 等，2018）。项目使用ALOS PALSAR和Landsat数据生成全球红树林地图数据，空间分辨率为25 m，用户精度为97.5%。本文选择2016年全球红树林地图数据。

（2）2016年哥白尼全球地表覆盖类型数据（Copernicus Global Land Cover Map， 以下简称“CGLS-LC100”）：数据由哥白尼土地监测服务团队提供，包括2015年—2019年全球土地覆盖分类图，空间分辨率为100 m，总体精度为83.3%。该产品将植被分为常绿针叶林、常绿阔叶林、落叶针叶林、落叶阔叶林、裸露/稀疏植被等类别。虽然该产品空间分辨率较低，但对植被类型的区分较为详细，故在样本选择时仅作为区分植被类型的参考。本文选择2016年地表覆盖分类图，孙德尔本斯地区分类图如图2所示。

图2 孙德尔本斯地表覆盖类型图Fig.2 Land cover map of the Sundarbans area

（3）Google earth 历史影像：本文通过Google Earth Pro软件查看2016年孙德尔本斯地区的历史影像。在对每个地类进行矩形感兴趣区选择时，利用Google earth历史影像对每个感兴趣区进行目视解译。

本文根据研究区植被覆盖类型进行分类，包括红树林、常绿植被、非常绿植被、耕地、无植被地区和水体，如表2所示。

表2 研究区地表覆盖类型Table 2 Land cover types in study area

为简化人工选择的过程、保证样本均匀合理地分布，先选择每类纯净的矩形感兴趣区，再生成随机样本点。首先，利用GMW2016 选取研究区内红树林矩形感兴趣区，对于其他地表覆盖类型则参考CGLS-LC100数据，结合Google earth历史影像进行矩形感兴趣区选择；之后，利用Google earth历史影像结合目视解译再次确保样本的准确性；最后，在每类矩形感兴趣区内随机选择100个样本点，共600个样本，得到分类样本数据，如图3所示。

图3 研究区样本点分布图Fig.3 Distribution map of sample points in study area

2.3 计算平台

谷歌地球引擎是2010年由谷歌公司、卡内基梅隆大学和美国地质调查局联合推出的云计算地理空间处理计算平台，该平台存储了海量遥感数据，可以使用基于网络的代码编辑器或开源的Python 界面进行数据处理，并拥有强大的计算能力，能够大大减少地理空间数据检索、处理和分析的时间，可为大规模遥感数据处理提供支持（Chen 等，2017；付东杰等，2021；Gorelick 等，2017；刘凯等，2019）。目前，已有较多研究利用GEE 平台对红树林进行遥感监测研究（Tieng 等，2019；王子予等，2020；Liu 等，2021）。本文在Windows 10、Intel Core i7 CPU、64 GB RAM 系统环境下运行GEE平台。

3 研究方法

首先，采用最小二乘回归融合Landsat 8 OLI和Sentinel-2 MSI数据，分别重建EVI和LSWI时间序列；其次，采用S-G 滤波对时间序列进行平滑处理，提取物候信息后将特征进行级联；最后，采用随机森林进行分类。本文采用的技术路线如图4所示。

图4 技术路线图Fig.4 Technical workflow of this study

3.1 指数计算与数据融合

红树林属常绿阔叶林，植被指数可有效区分红树林和非常绿植被，同时也可在一定程度上区分其他常绿植被。此外，红树林生长环境受潮汐变化影响，具有特殊的水体淹没信息，采用水体指数可进一步区分红树林与其他常绿植被。为此，本文同时利用植被指数和水体指数对红树林进行提取研究，如表3所示。

表3 本文采用的指数Table 3 The indexes used in this study

由于Landsat 8 OLI 和Sentinel-2 MSI 传感器具有相似的光谱特性，将这两种数据融合具有全球覆盖以及中等空间分辨率、高时间分辨率地面观测的潜力（Zhang 等，2018）。本文使用最小二乘回归法将Landsat 8 OLI与Sentinel-2 MSI进行融合，从而获得稠密时间序列数据。首先，筛选2016年全年的有效数据，并对数据进行去云处理，同时将Sentinel-2 MSI 数据重采样至30 m 并计算4 种指数；其次，筛选出每张Landsat 8 OLI 影像前后8 d内日期最接近的Sentinel-2 MSI 影像，并获取对应的指数数据；最后，利用匹配的指数数据拟合线性最小二乘回归模型，以此对Sentinel-2 MSI 数据进行转换，重构时间序列数据。最小二乘回归数据融合原理如下：

式中，xi为Sentinel MSI 植被指数集；yi为Landsat OLI植被指数集，每个像素相应数据集的植被指数按照时间顺序排列。

式中，ri为残差，分别对β1，β2求偏导，求得在残差最小的情况下每个像素最佳的β1，β2。

3.2 时间序列选择与滤波处理

为统一数据的时间间隔并在一定程度上削弱大气影响，采用最大合成法构建8 d 时间序列。为有效评估不同时间序列对后续分类的效果，采用J-M距离计算类间可分性，对比得到更有利于区分不同地物的时间序列。采用S-G 滤波（Savitzky 和Golay，1964）对合成后得到的时间序列曲线进行平滑处理，消除曲线中的异常值。

在滤波过程中，需要设定两个参数：分别是窗口长度（N）和多项式拟合阶数（k）。为达到最优的平滑效果，需对S-G 滤波的参数进行寻优。通过保真性和降噪能力两方面评价滤波效果，其中，保真性利用均方根误差评定，该值表示平滑前后时间序列的差异程度，差异越小，保真性越好；降噪能力通过最值统计进行评定，若滤波后最大值和最小值被去除，则说明达到了去除噪声的目的。

3.3 物候信息提取

本文选取了生长期始期、生长期末期、生长期长度、基值等13 种物候信息，如表4所示。其中，对于生长期始末期，使用最大分离法（Descals等，2021）来确定，其余物候信息根据TIMESAT 3.3 中的定义（Jönsson 和Eklundh，2004）进行提取；此外，对于生长期的积分值利用逐段计算面积的方法来近似代替。

表4 物候信息描述Table 4 Phenological metrics description

3.4 分类方法及精度评价

本文采用随机森林算法进行分类，该方法分类精度高，模型训练时间少，且对训练样本数量、质量以及样本不均衡问题的敏感性较低，在遥感影像分类中应用广泛（Talukdar 等，2020；张露洋等，2021）。研究区包含6 种地类，共计600 个样本点，其中红树林样本数为100个。采用随机抽样方法，将每类样本打乱后抽取50 个，共计得到训练样本和验证样本各300个。通过混淆矩阵计算用户精度、总体分类精度和Kappa系数，作为分类结果的评价指标。

4 实验与分析

4.1 指数计算与分析

基于Landsat 8 OLI 和Sentinel-2 MSI 遥感影像分别计算4 个指数，如图5所示。从图5 中可以看出，对于植被指数，EVI 和NDVI 均可有效区分红树林、水体和无植被地区，一定程度上能够区分红树林和耕地。对比图5（a）和（b）、（e）和（f），红树林和其他常绿植被的指数差异不太明显。对于水体指数，对比图5（c）和（d）、（g）和（h），红树林、非常绿植被和其他常绿植被LSWI 差异较NDWI 大，但LSWI 较难区分红树林和水体，然而二者能被NDWI有效区分。本实验主要利用水体指数区分红树林和其他常绿植被，故选择LSWI 更为有效。综上所述，初步分析水体指数中选择LSWI指数的区分度更高，但具体指数选择以及区分能力的判断还要通过指数构建的时间序列来分析。

图5 不同类别样本指数年平均值Fig.5 Index annual mean of different samples

4.2 融合有效性分析

通过计算融合前后时间序列的相关性分析融合可靠性。从表5 可以看出，相关系数均在0.7 以上，融合前后时间序列具有较强的相关性。时间序列融合效果较好，具有较高的相关性。

表5 数据融合相关性分析Table 5 Data fusion correlation analysis

通过对比数据融合前后有效观测值的个数分析融合有效性，如图6所示。从图6（a）和（b）可以看出，数据融合后研究区内全年可获得的有效观测值个数、明显提高。融合前，较多区域全年有效观测值个数仅在0—20；通过数据融合，研究区内全年有效观测值个数基本在20 以上。为验证全年有效观测值的时间分布能够满足构建时间序列的要求，对每个季度融合前后有效观测值个数进行统计对比。从图6（c）—（f）可以看出，通过数据融合，能够保证每个像素在每个时间段内均存在有效观测值，可以构建较为稠密的时间序列。

图6 数据融合前后有效观测值对比Fig.6 Comparison of the number of good observations before and after data fusion

4.3 基于J-M距离的指数可分性评价

通过J-M 距离计算不同指数类别间的可分性，如表6所示。可以看出，EVI 和LSWI 对红树林、常绿植被和非常绿植被这3 个类别的可分性相较NDVI 和NDWI 要高；不同指数对红树林、耕地、裸地和水体这4个类别的可分性均较高，本文选择EVI和LSWI指数构建时间序列。

表6 指数可分性评价结果Table 6 Results of index separability evaluation

4.4 S-G滤波参数选择与结果分析

本文对红树林、常绿植被、非常绿植被和耕地4 个类别的EVI 和LSWI 指数的时间序列进行滤波参数寻优，如图7所示。首先，根据RMSE 寻优曲面，选取RMSE较低即保真性较好的参数；再对经过该参数平滑的时间序列进行最值统计，分析去噪能力；最后综合保真性和去噪能力选择最优参数。

图7 红树林S-G滤波参数选择Fig.7 Selection of S-G filtering parameters in mangroves

表7给出了EVI、LSWI时间序列S-G参数寻优结果，综合判断选取S-G 滤波窗口长度为9，阶数为3。图8 给出了红树林平滑后的时间序列，可以看出平滑效果较好，能够在保留原始时间序列特征的同时达到去噪效果。

图8 红树林时间序列Fig.8 Time series of mangrove

表7 S-G滤波参数统计Table 7 S-G filter parameter statistics

4.5 物候信息提取与结果分析

4.5.1 EVI时间序列物候信息

利用EVI 时间序列提取物候信息，如图9—10所示。从图9 可以看出，生长期始期、生长期长度、基值、生长期内最大值、生长期振幅、整体生长期积分值、生长期始末期的值区分度较高，红树林和其他植被的差别较为明显。首先，对比红树林和其他常绿植被，红树林的生长期始期要明显晚于其他常绿植被，故整个生长期的长度也相对较短，生长期始末期的值要略高于其他常绿植被；在生长期内，红树林的最大值要略高于其他常绿植被，但生长期内的振幅要比其他常绿植被略小。其次，对比红树林、非常绿植被和耕地，物候信息有明显不同。红树林生长期较长，而非常绿植被和耕地都具有明显的季节性，生长期均较短，进入生长期之后较短时间就达到最大值，故两者生长期内的物候尺度，如积分值，增长下降率都有明显差异。

图9 EVI物候信息图Fig.9 The map of EVI phenology metrics

从各类植被中选取具有代表性的样本点描述具体物候信息，如图10所示。首先，对比红树林和其他常绿植被，红树林生长期始期约在第150天左右，常绿植被在第50—100天之间，红树林生长期的开始时间要比常绿植被晚两个月左右，而两者的生长期末期基本都在第250天之后，时间较为一致，红树林生长期末期的EVI值要比常绿植被高大约0.1；在生长期内，红树林和常绿植被基值均在0.35—0.4，红树林生长期内的面积明显小于常绿植被。其次，对比红树林、非常绿植被和耕地，非常绿植被的生长期长度在100 d 左右，比红树林短的多。非常绿植被生长期始末期的值在0.45 左右，与红树林较为相近，而耕地生长期始末期的值较低，旱地在0.3左右，水田则大约为0。

图10 样本点EVI物候信息（a—m具体含义参见表4）Fig.10 EVI phenology metrics of sample points

综上所述，虽然红树林属于常绿植物，在物候信息上与常绿植被有明显的相似性，但也仍有不同之处。不同植被类型EVI时间序列的物候信息有明显区别，可有效区分红树林和其他植被类型。

4.5.2 LSWI时间序列物候信息

通常物候信息都是通过植被指数构建时间序列进行提取，但潮汐变化会导致红树林淹没状态的不同，月内、年内的潮汐变化会影响到红树林的水分信息，通过LSWI 大致记录潮汐变化构建时间序列是区分红树林和其他常绿植被的关键。本实验通过最大合成法来选择每个时间段LSWI 的最大值，这在一定程度上能够减少日内高潮、低潮的影响。避免某些时间段受潮汐影响本该有较强的水分信息，但由于选取了该时段内某天低潮时刻的值而削弱了水分信息。利用LSWI 时间序列提取物候信息，如图11—12所示。从图11 可以看出生长期始末期、生长期长度、基值、生长期中期、生长期结束时的下降速度、生长期始末期的值区分度较高，红树林和其他植被的差别较为明显。首先，红树林生长期始中末期均早于其他植被，生长期始末期的值明显高于其他植被；其次，在生长期内，红树林基值高于除水田外的其他植被，振幅略小于其他植被。

图11 LSWI物候信息图Fig.11 The map of LSWI phenology metrics

从各类植被中选取具有代表性的样本点描述具体物候信息，如图12所示。首先，红树林生长期始期在第50—100 天，其他类型植被均在第100 天之后，较红树林始期早一个月左右；红树林生长期末期在第150—200 天，明显早于其他植被150 天左右；红树林生长期始末期的值在0.4左右，明显高于其他内陆植被，水田生长期始末期的值在0.55 左右；其次，在生长期内，红树林基值约0.35，要明显高于除水田外的其他植被，最大值在0.5左右，比常绿植被和非常绿植被高0.1左右。

图12 样本点LSWI物候信息（a—m具体含义参见表4）Fig.12 LSWI phenology metrics of sample points

综上所述，由于红树林受到海水周期性侵蚀，具有较强的潮汐信息，不同植被类型LSWI 时间序列的物候信息有明显区别，可有效区分红树林和其他植被类型。

4.6 分类结果与分析

4.6.1 随机森林参数选择

随机森林算法中决策树的个数对分类结果有较大影响，需对该参数进行最优化选择。本文设置决策树3，5，10，20，50，100，150，200，分别运行10 次进行统计分析。图13 显示了Kappa 系数随决策树个数的变化情况，当决策树个数由3—150 变化时，精度有明显地提升；当决策树个数达到150 时，分类精度最高；当决策树个数增加至200 时精度有明显地下降。为此，将决策树的个数设置为150进行分类。

图13 随机森林参数选择Fig.13 Random forest parameter selection

4.6.2 数据融合效果对比

表8给出了融合前后的精度对比结果，可以看出，通过数据融合总体精度由89.44% 提升到91.02%，Kappa 系数由0.873 提升到0.892。图14显示了仅使用Sentinel-2 MSI 进行分类的结果，与GMW2016 数据集对比可以看出，由于时间序列较为稀疏导致的红树林和其他植被的时间序列存在相似性，在内陆地区存在较多的红树林错分现象。由此可知，通过数据融合构建了更加真实的时间序列，从而能够提取更为准确的物候信息，提高了不同植被之间的区分度，分类精度显著提高。

表8 数据融合精度对比Table 8 Data fusion accuracy comparison

图14 仅使用Sentinel-2 MSI数据分类结果Fig.14 The classify results of only using Sentinel-2 MSI dataset

4.6.3 不同输入特征的红树林提取结果对比

为了分析EVI 和LSWI 以及物候信息对红树林提取结果的影响，实验设置不同的特征组合进行对比：（1）分别输入［EVI，LSWI］的时间序列和［EVI，LSWI］时间序列的物候信息进行对比；（2）分别将EVI、LSWI、［EVI，LSWI］时间序列提取的物候信息作为输入进行对比。

根据第1组实验对比，如表9所示，利用EVI、LSWI时间序列共96个特征进行分类得到总体精度为89.10%，Kappa 系数为0.869，而利用EVI、LSWI时间序列的物候信息共26个特征进行分类得到的总体精度为91.02%，Kappa 系数为0.892，精度有明显地提升。结果与GMW2016 数据集进行对比，从图15（a）可以看出，利用时间序列直接进行分类，不仅在红树林与其他植被过渡区域内有大量的错分漏分现象，且研究区北部也有大量的植被被错分为红树林。由此可见，物候信息能够反应植被随季节的变化特征，对时间序列进行物候信息的提取可以更好地显示红树林的季节特征，提高不同地物类别的区分度。

图15 不同特征分类结果对比Fig.15 Comparison of different features classification results

表9 不同输入特征的分类精度对比Table 9 Comparison of classification accuracy of different input features

根据第2 组实验对比，如表9所示，分别利用EVI 和LSWI 时间序列的物候信息作为特征输入进行分类总体精度为76.91%和81.33%，Kappa 系数分别为0.723 和0.776，而同时使用EVI 和LSWI 时间序列的物候信息使得总体精度分别提高14.11%和9.69%，分类精度显著提升。结果与GMW2016数据集进行对比，从图15（c）可以看出，仅利用EVI时间序列的物候信息进行分类时，虽然能够较好地区分水体和红树林，但红树林与其他类别具有相似的植被信息，导致产生大量错分现象；根据图15（d），仅利用LSWI时间序列的物候信息进行分类时，一些无植被沿海沼泽地区被错分为红树林，但红树林叶片含水量以及土壤湿度与其他植被的显著差异大大减少了陆地植被的错分现象。由此可见，仅利用EVI时间序列的物候信息无法较好地区分红树林和其他植被，但未添加植被特征信息也无法将红树林与无植被沿海沼泽地区较好区分。因此，同时使用EVI 和LSWI 时间序列的物候信息可以更加充分地利用红树林的物候特性，从而能够更为精确的提取红树林。

从图15（b）可以看出，利用EVI、LSWI时间序列的物候信息分类效果较好，仅在红树林与陆地植被的过渡区域存在少量错分现象。

综上所述，物候信息能够更好地表征时间序列随季节变化的规律性，利用物候信息进行分类有效增强了不同植被的可分性，从而提高红树林分类精度；此外，同时使用EVI、LSWI 时间序列的物候信息进行分类较好地弥补了单独使用一个特征指数的不足，有利于红树林的精确提取。

4.6.4 红树林提取结果分析

从图16 可以看出，本文方法提取的红树林大多呈斑块状分布在沿海滩涂及附近地区，内陆地区无明显的红树林分布，但有较多耕地、常绿植被和非常绿植被。分类图中主要的错分区域为红树林和陆地植被过渡区域，有部分陆地植被被错分为红树林，但整体结果较可信。

图16 2016年Sundarbans地区分类图Fig.16 Map of the Sundarbans area in 2016

4.6.5 相关研究结果对比

近年来有较多学者针对孙德尔本斯地区进行红树林提取研究，对孙德尔本斯地区红树林相关研究进行整理，如表10所示。

表10 相关研究结果对比Table 10 Comparison of related studies

由表10 可知，近年来对孙德尔本斯进行红树林的提取研究中，都未考虑利用多源遥感数据和红树林的物候信息，但也有研究中取得了较为理想的精度。

在小范围提取研究中，Gupta 等（2018）考虑到红树林受潮汐淹没的特性，利用NDVI 和NDWI指数构建红树林综合识别指数；Jia 等（2019）也利用此特性构建指数，得到了较为理想的精度；除了利用多光谱数据，Kumar 等（2019）利用EO-1高光谱数据，在提取亨利岛及周边红树林的基础上对红树林种群进行了分类。

但以上研究的提取范围都较小，在较大范围的提取研究中，较多研究会预先对研究区进行圈定。如Dan 等（2016）去除高于30 m 地区和非植被地区并根据海岸线建立缓冲区；Islam等（2019）选择的研究区内不存在除红树林外的植被区域，这大大提高了红树林的可分性。

本研究对整个孙德尔本斯地区进行红树林提取，未对研究区进行预先圈定，排除了人工选择的干扰，且研究区内包含各种常绿、非常绿等易于红树林混淆的植被。与已有研究不同，本文通过：（1）融合Landsat 8 OLI 和Sentinel-2 MSI 数据构建更加可靠的时间序列；（2）从时间序列中提取物候信息增大红树林与其他植被的可分性；（3）联合使用EVI 和LSWI 时间序列的物候信息弥补单一指数的不足，客观准确地对红树林进行了提取，也得到较为理想的精度。

5 结论

本文基于GEE 平台获取Landsat 8 OLI 和Sentinel-2 MSI 影像数据集进行红树林提取研究。首先，基于最小二乘回归构建两个传感器在相同指数之间的关系，重建时间序列数据，之后选择可分性较高的EVI 和LSWI 指数时间序列；其次，对两个指数的时间序列数据进行S-G 滤波处理，并分别提取13 种物候信息；最后，将两个指数的物候信息进行特征级联，利用随机森林对孙德尔本斯地区2016年的红树林进行提取，总体精度为91.02%，Kappa 系数为0.892，红树林提取效果较好。本文主要结论如下：

（1）融合Landsat 8 OLI 和Sentinel-2 MSI 数据可以构建更加可靠的时间序列，有利于得到更准确的物候信息，融合后总体精度提高1.58%，Kappa系数提高0.02。

（2）物候信息可以增大红树林与其他植被的可分性，利用物候信息进行分类可使总体精度提高1.92%，Kappa系数提高0.02。

（3）联合EVI 和LSWI 时间序列的物候信息较好地弥补了单一指数的不足，有利于红树林的精确提取，与单独使用这两个指数相比，总体精度分别提高14.11%和9.69%，Kappa 系数分别提高0.17和0.12。

本文充分顾及单一来源遥感数据监测红树林的不足、不同指数对红树林的区分特点以及红树林的物候特性，融合Landsat 8 OLI 与Sentinel-2 MSI 时间序列数据，联合EVI 和LSWI 时间序列的物候信息，对孙德尔本斯地区红树林进行提取研究。与已有的研究相比，本研究的主要特色和不同体现在：一方面，利用EVI 和LSWI 时间序列描述红树林的物候信息，解决了红树林与其他常绿植被难以区分的问题；另一方面，首次将物候信息引入到红树林遥感提取中，能够更客观、准确地提取红树林，验证了物候信息对红树林提取的有效性。此外，本方法利用红树林特性对其进行提取，理论上来说能够普遍使用于各个地区，未来将会对该方法进行推广并验证其有效性。论文不足之处在于仅利用了一年的观测数据，物候信息可能会有一定的偏差，后续研究将利用更长周期的时间序列数据进行红树林提取研究。

志 谢 感谢全球红树林监测组织和哥白尼土地监测服务团队发布的GMW2016和CGLS-LC100数据。