常云蕾,廖静娟,张 丽
(1. 中国科学院空天信息创新研究院数字地球重点实验室,北京 100094;2. 可持续发展大数据国际研究中心,北京 100094;3. 中国科学院大学,北京 100049)
红树林是生长在热带及亚热带海岸潮间带的湿地木本植物群落,是全球生产力最高的生态系统之一(林鹏,1997),在防风消浪、促淤保滩、维持生物多样性、固碳储碳和海产品供给等方面发挥着重要作用,为沿海生态环境保护、社会发展及经济进步提供着重要支撑(张和钰 等,2013;Roy,2014;郑艺 等,2019)。然而,由于红树林地处海陆交界的特殊地带,是全球最脆弱的生态系统之一,水产养殖、城市发展、环境污染、森林砍伐、海平面上升等问题时刻威胁着红树林的生长(Pham et al., 2016)。在过去的半个多世纪中,全球红树林面积急速下降,损失已经超过1/3(Alongi,2002)。21 世纪以来,全球红树林面积减少的速度虽有所减缓,但每年仍以0.4%左右的速度减少(Hamilton et al., 2016)。按照这种趋势发展,全球红树林将会逐渐消失,其独特的经济、社会、生态价值也将荡然无存,严重威胁到人类的生存与发展(贾明明 等,2021)。因此,开展红树林的长时序变化趋势分析并定量探究红树林损失与潜在驱动因素的关系对红树林的保护、恢复及可持续管理具有重要意义。
遥感技术具有覆盖范围广、重复周期短、信息量大、限制条件少等优势,应用遥感技术进行红树林监测可以有效避免传统实地调查方法中工作量大、效率低的问题(周振超 等,2018)。目前,遥感技术在红树林的监测已被广泛采用,并获得诸多成果(赵玉灵,2010;Kuenzer et al., 2011; Alejandro et al., 2019; Halder et al., 2021)。然而,受限于红树林自身的复杂性,已有的红树林遥感动态监测研究还存在一些不足。首先,只采用少数时相的遥感数据,其时间跨度较大,在刻画红树林动态变化过程中精细程度不足,红树林演变分析的不确定性增加;其次,在处理大尺度长时序遥感影像时需要很大的计算资源,在目前已有的大范围红树林动态监测中,较少开展长时间序列的红树林趋势变化分析。本研究首先基于已有的红树林分布数据集,开展1990—2020年东南亚红树林面积变化的时空分布特征分析;然后基于Google Earth Engine 云计算平台①Google Earth Engine(GEE)是谷歌公司在2012 年建立的用于全球范围地理空间分析的云计算平台(Patel et al., 2015),为大尺度和长时序的遥感资料分析提供高效便捷的处理手段,在农业和林业应用、水体和城镇提取以及潮间带监测等领域都发挥重要作用(Gorelick et al., 2017; Shelestov et al., 2017),采用Landsat长时序影像数据和红树林分布数据集,结合Theil-Sen中值趋势分析和Mann-Kendall检验方法,分析红树林增强植被指数(Enhanced Vegetation Index, EVI)的时空变化趋势;最后,利用地理探测器定量分析自然因子和社会因子对东南亚红树林面积变化的影响。以期为东南亚红树林的保护、恢复和利用提供科学合理的理论依据和决策支持。
东南亚位于亚洲东南部,气候条件主要为热带季风和热带雨林气候,其地理位置和热带型的气候条件为红树林提供了适宜的生长环境。东南亚红树林具有极高的生物多样性,为数百万人提供了赖以生存的生态系统服务(Richards et al., 2016)。其中,有红树林生长的国家共10 个,包括东帝汶、菲律宾、柬埔寨、马来西亚、缅甸、泰国、文莱、新加坡、印度尼西亚和越南。泰国、柬埔寨、缅甸及越南红树林生长区为热带季风气候,一年中有明显的干湿两季,年均降水量在1 000 mm以上,盛行西南季风。东帝汶、菲律宾、马来西亚、文莱、新加坡和印度尼西亚红树林生长区为热带雨林气候,分布在赤道两侧、南北纬10°之间,终年高温多雨,年降水量在2 000 mm以上,没有明显的旱季且年温差较小。
将1990、2000、2010、2015 和2020 年东南亚有红树林生长的所有区域及其周边1 km缓冲区作为红树林分布区(图1)。该分布区内红树林种类丰富、结构复杂,是世界上生物多样性最为丰富的红树林生态系统(Gandhi et al., 2019)。
图1 研究区位置Fig.1 Location of the study area
1.2.1 遥感数据
1)红树林分布数据 采用的红树林分布数据时间范围为1990、2000、2010、2015 和2020 年。其中,1990、2000、2010 和2015 年的红树林分布数据是Guo等(2021)基于深度学习模型提取,2020年的红树林分布数据由廖静娟(2022)课题组以2010年的结果为基准进行目视解译得出。上述数据均经过高精度人工修正,1990、2000、2010、2015和2020 年红树林分布数据的总体精度分别为87%、90%、88%、88%、91%。该数据集均是利用1990和2000年的Landsat5 TM 影像、2010年的Landsat5 TM 和Landsat7 ETM+影 像、2015 和2020 年 的Landsat8 OLI影像生产的。
2)EVI 分布数据 NDVI(归一化植被指数,Normalized Difference Vegetation Index) 和EVI 在监测和评估植被数量和条件的时空变化方面表现出良好的动态范围和敏感性(Liu et al., 1995)。然而在湿润环境下高密度植被监测中,NDVI 会出现过饱和的现象,而EVI改进了大气或者生长地表所造成的误差;此外,与NDVI相比,EVI对针叶/阔叶树冠结构更加敏感,尤其在高生物量区域,NDVI会逐渐饱和,而EVI 对冠层变化仍然保持敏感(Huete et al., 2002)。红树林属于常绿阔叶林,是世界上生产力和碳密度最高的生态系统之一,其结构属性(冠层高度)具有高度空间异质性,因此,选择采用EVI分析红树林生长变化趋势。
使用的EVI分布数据为GEE云平台存储的EVI年际均值数据②https://developers.google.com/earth-engine/datasets/catalog。其中,1990—2012 年由Landsat5 TM(Collection 1 Tier 1)生成、2013年由Landsat7(Collection 1 Tier 1)生成、2014—2020 年由Landsat8(Collection 1 Tier 1)生成。
为了更好地确定不同时期的研究区范围,选取1990、2000、2010、2015 和2020 年红树林出现过的区域作为1990—2020年整个研究时段的研究区范围。利用1990、2000、2010、2015 和2020 年合并后红树林范围对EVI分布数据进行掩膜处理,从而获取1990—2020年红树林EVI逐年均值分布数据。
1.2.2 辅助数据 辅助数据主要为驱动因子数据,包含影响红树林面积变化的年均降水(X1)、年均温(X2)、地形(X3)、养殖池面积(X4)、红树林距道路的最短距离(X5)以及人口(X6)。降水、气温数据来源于英国大气科学中心③https://crudata.uea.ac.uk/cru/data/hrg,分辨率为0.5°;数字高程模型(Digital Elevation Model, DEM)数据采用30 m 分辨率的SRTM 高程数据④http://www.geodata.cn/data/;养殖池数据采用中国科学院南京地理与湖泊研究所提取的1990、2000、2010、2015 和2020 年的30 m 分辨率水产养殖区空间分布范围数据集(罗菊花,2019);矢量道路数据集和人口公里格网数据集分别来源于全球生物多样性模型(Global biodiversity model⑤www.globio.info)和世界人口数据集(WorldPop⑥www.worldpop.org)。
采用地理探测器进行红树林变化驱动因素分析。由于地理探测器的自变量为类型量,故需要对驱动因子数据进行分类。自然断点法可以有效避免人为引入的干扰,最大限度地提高类之间的差异(Chen et al., 2013)。数据离散化的分类效果可以通过地理探测器的q统计量确定,q值越大,分类效果越好(王劲峰 等,2017);因此,借助R 语言,将自然断点法的分类等级数设置为5~12类,筛选q值最大的空间尺度作为地理探测器自变量等级划分的标准。q值最大时不同驱动因子的划分结果为:年均降雨、年均温、高程划分8级;养殖池面积变化10 级;红树林距道路的最短距离9 级;人口7 级(图2)。所有数据均为2020年分布数据,其中红树林距道路的最短距离(d)以点状数据显示,即以蜂窝网为单位统计红树林面积,计算每个蜂窝网距道路最短距离,以点的形式对距离分级显示可以与其他驱动因子的分布具有更相似的视觉效果。
图2 1990—2020年东南亚红树林面积变化驱动因子分级Fig.2 The grading of driving factors of mangrove area change in Southeast Asia from 1990 to 2020
将Theil-Sen 中值趋势分析和Mann-Kendall 检验方法结合,能很好地应用于长时间序列数据的趋势分析中,目前在植被的长时间序列分析中已得到广泛应用(Sen, 1968; Tucker et al., 1991)。
Theil-Sen中值趋势分析是一种稳健的非参数统计趋势分析方法,可以有效规避数据误差,减少数据异常值的影响(袁丽华 等,2013)。其基本原理是计算n(n- 1)/2 个数据组合的斜率中位数,公式为:
式中:EVIi和EVIj分别为红树林第i年和第j年的EVI值。当SEVI> 0时,表明该时间序列数据呈增长趋势,否则呈现退化趋势。
Mann-Kendall是一种非参数统计检验方法,该方法不需要样本数据服从一定的分布,具有抵抗数据误差的能力(Milich et al., 2000)。常被用来判断趋势的显著性,公式为:
式中:EVIi和EVIj分别为红树林第i年和第j年的EVI 值;n为时间序列的长度;sgn 为符号函数;Z取值范围为( -∞, + ∞)。在给定的显著性水平α(一般取值为0.05)下,当|Z|>u1-α2时,表明该序列具有显著性变化。
地理探测器是探测空间分异性,以及揭示其背后驱动力的一组统计学方法(王劲峰 等,2017)。采用地理探测器中的因子探测器和交互作用探测器探究不同驱动因素对红树林变化的影响及交互作用。
因子探测器可以定量描述不同驱动因素的相对影响力大小,公式为:
式中:q代表驱动因子驱动力,值域为 [0, 1],q值越大驱动力越强,反之则越弱;L为因子的分层;Nh和N分别为层h和全区红树林变化分布及因子对应的单元数;σ2h和σ2分别为层h和全区红树林变化方差。
交互作用探测器用于评估不同影响因子共同作用时,对红树林面积变化是否存在交互作用以及交互作用的类型和强度,评估依据参考文献(王劲峰等,2017)。
东南亚红树林总面积在1990—2020年呈显著下降趋势(图3),由1990 年的520.55 万hm2减少至2020年的373.77万hm2,整体上每年以1.1%的平均速度消失。其中,1990—2015年持续下降,2015—2020年面积有所回升,表明相关政策干预在红树林恢复及保护方面的重要作用(Friess et al., 2016)。
图3 1990—2020年东南亚红树林面积变化趋势Fig.3 Change trends of mangrove area in Southeast Asia from 1990 to 2020
分别统计1990、2000、2010、2015 及2020 年东南亚存在红树林的各个国家红树林面积(表1)并计算其变化率(图4)。可知,1990—2020年各个国家红树林面积整体上均呈下降趋势,其中,新加坡、东帝汶、越南下降最为严重,面积下降率达到69%以上,其次为泰国、缅甸、马来西亚、柬埔寨、文莱、印度尼西亚和菲律宾。2015—2020 年,部分国家的红树林面积呈现少量上升趋势,包括菲律宾、文莱、印度尼西亚、马来西亚、缅甸,但新加坡和越南红树林面积仍在急剧下降,柬埔寨和泰国下降趋势有所减缓。
图4 东南亚各国红树林面积变化率趋势Fig.4 Trends of mangrove area change rates in Southeast Asian countries
表1 1990—2020年东南亚不同国家红树林面积Table 1 Mangrove area in different countries of Southeast Asia from 1990 to 2020 万hm2
正六边形具有很低的面积周长比,可以有效减少由于格网形状的边界效应带来的样本偏差,故采用正六边形统计红树林面积变化率,并对其分级显示,分为红树林完全消失、显著下降、轻微下降、显著增长、轻微增长、稳定不变。分别统计边长为10、15、20、25、30、35、40 km 的蜂窝网数据,结果表明边长为30 km 时更能体现红树林空间变化的分异特征,因此选取30 km 边长的蜂窝网形式表示1990—2020 年红树林面积变化率的空间分布(图5)。
从图5可以看出,东南亚各区域红树林面积变化存在明显的空间分布特征,从整体上来看,红树林面积严重下降,甚至大部分红树林消失;而增长的区域呈零星分布,且大多分布在菲律宾以及印度尼西亚西部和东部,越南北部也存在少量红树林增加区域。
图5 1990—2020年东南亚红树林面积变化率空间分布Fig.5 Spatial distribution of mangrove area change rate in Southeast Asia from 1990 to 2020
以蜂窝网为基本单位对30年红树林面积变化率分级统计,结果表明:1990—2020年,红树林完全丧失区域高达40.11%,变化率在-99%~-50%的区域占28.57%,在-49%~-5%的占8.79%,在-4%~-1%占13.95%;面积变化增加的区域(变化率在1%~99%)共占8.18%;面积不变的区域即面积变化率为0%的区域仅占0.4%。
为验证年际EVI均值分布数据的可用性,利用GEE中存储的时间分辨率为8 d的EVI分布数据集,按季节(12-2 月、3-5 月、6-8 月、9-11 月)计算其EVI均值分布。结果表明,按照季节划分的红树林EVI均值变化趋势与年均变化趋势整体上具有一致性。因此,采用Theil-Sen中值趋势法分析红树林年际EVI均值分布的时间序列变化趋势具有一定可行性。
结合Theil-Sen 中值趋势分析和Mann-Kendall检验方法,依托GEE 平台,探求1990-2020 年东南亚红树林EVI变化趋势的空间分布特征。由于研究时间跨度大,基本不存在SEVI值为0的区域,因此根据实际情况,取-0.005和0.005作为临界值,划分红树林面积改善、稳定不变及退化的区域,并判断Mann-Kendall检验在0.05置信水平上红树林趋势变化的显著性。叠加Theil-Sen中值趋势分析和Mann-Kendall检验结果,得到像元尺度上红树林趋势变化数据,并将结果划分为5种类型(表2)。东南亚整体红树林面积退化趋势明显,退化区域占该区域红树林总分布面积的79.25%,稳定不变即没有发生显著变化的仅占0.42%,改善区域仅占20.32%。
表2 红树林EVI变化趋势分级Table 2 Hierarchical statistics of mangrove EVI trend
从图6可以看出,东南亚红树林并非呈单一的变化趋势,而是多种变化趋势共存,但总的趋势是退化。选择东南亚红树林广泛分布的保护区或受多种因素影响发生显著变化的部分区域作为典型研究区,分析其时空变化趋势:
1) 伊洛瓦底江三角洲(Ayeyarwady Delta)(图6-a):该区域是缅甸最大的红树林分布区,红树林变化趋势受孟加拉湾和安达曼海影响较大,沿海红树林主要受飓风影响(Veettil et al., 2018)。退化区域集中分布,尤其在红树林大面积分布区退化趋势更为明显,显著退化面积远远多于轻微退化。改善区域分散于红树林大面积分布区边缘及小面积分布区,其中,边缘主要呈显著改善,轻微改善主要位于红树林小面积分布区。
2)高岗(Koh Kong)(图6-b):该区域为柬埔寨红树林物种最丰富和分布面积最大的地区。从整体上看,该地区1990-2020年红树林退化的区域明显多于改善的区域,显著退化的区域在空间上由中心向四周扩展,尤其在东部的延伸最为明显。改善区域少量分布于红树林分布边缘,且主要呈轻微改善趋势。
3)吉奥(Can Gio)(图6-c):该区域为越南四大红树林保护区之一。该保护区红树林改善和退化伴随分布,退化趋势更为明显,四周主要呈显著退化趋势,并由四周向中心逐渐扩展。红树林改善区域分散于保护区四周和中心,且四周有少量显著改善区域,中心主要为轻微改善。
图6 1990—2020年东南亚红树林EVI变化趋势空间分布(a. 缅甸:伊洛瓦底江三角洲;b. 柬埔寨:高岗;c. 越南:吉奥红树林生物圈保护区;d. 马来西亚:双溪哈吉多拉尼;e. 文莱:斯里巴加湾;f. 印度尼西亚:宾图尼湾红树林自然保护区)Fig.6 Spatial distribution of EVI changes in mangrove in Southeast Asia from 1990 to 2020 (a. Myanmar: Ayeyarwady Delta;b. Cambodia: Koh Kong; c. Vietnam: Can Gio Mangrove Biosphere Reserve; d. Malaysia: Sungai Haji Dorani;e. Brune: Bandar Seri Begawan; f. Indonesia: Bintuni Bay Mangrove Nature Reserve)
4)双溪哈吉多拉尼(Sungai Haji Dorani)(图6-d):该区域为马来西亚半岛西海岸的红树林分布重要栖息地。该地区红树林轻微退化和轻微改善区域分布比较广泛,且呈明显的交叉分布。显著退化的区域零星分布于中部,并在北部有少量集中分布。显著改善的区域大多位于轻微改善的边缘,在南部区域较为明显。
5)斯里巴加湾(Bandar Seri Begawan)(图6-e):该区域为文莱主要的红树林分布区域。1990-2020年该地区红树林整体上呈明显的退化趋势,退化的区域远远大于改善的区域。其中,显著退化和轻微退化占比相当,红树林分布边缘受各种因素的影响退化严重,内部多呈轻微退化,改善的区域仅在该区域北部有零星分布。
6)宾图尼湾(Bintuni Bay)(图6-f):该区域是印度尼西亚红树林自然保护区。区域内红树林变化趋势具有明显的空间分布特征,轻微改善的区域主要分布在该地区中部,退化的区域主要位于东部和西部,其中西部退化更为明显,呈显著退化,而东部红树林呈轻微退化趋势。
采用地理探测器中的因子探测器分析不同驱动因 子 对1990-2000、2000-2010、2010-2015 及2015-2020年4个时段红树林变化的影响程度(表3)。驱动因子的p值均为0,表明所选驱动因子对红树林面积变化的空间变化特征均存在显著影响。总体而言,1990—2020年养殖池面积的增加是东南亚红树林面积下降的主要驱动因素,驱动力均>30%。其次,红树林距道路的最短距离和人口变化也具有较强的影响力;其中,红树林距道路的最短距离在1990—2015 年q值均>20%,人口对红树林影响的驱动力在2010-2015 年达到30%,2015-2020 年为24.5%,说明红树林距道路的最短距离越短,人口越多,对红树林面积下降也存在显著影响。年均温、降水以及地形对红树林变化的影响程度相对较低。1990-2000年养殖池面积快速扩张,对红树林面积下降的驱动作用达到48%;2000年后,随着养殖池面积缓慢下降,红树林受养殖池的影响也逐渐减缓。年均降雨、年均温、地形和红树林距道路最短距离的驱动力相对稳定,而人类因生活及生产活动需求的不定性,对红树林变化的影响程度较大,人口数量对红树林变化的驱动力在2000-2020年先增后减。分别统计各驱动因子在1990-2020年内驱动力的均值,结果表明,红树林变化受驱动因素的影响程度排序主要为:养殖池面积>人口>红树林距道路最短距离>年均温>年均降雨>地形。
表3 1990—2020年东南亚红树林面积变化驱动因子驱动力Table 3 Driving forces of mangrove area change drivers in Southeast Asia from 1990 to 2020
交互作用探测器用于识别不同驱动因子对红树林空间分异的交互作用(图7),不同因子的相互作用驱动力均高于单一因子,其作用方式包括非线性增强和双因子增强2种。这表明1990—2020年红树林面积变化的空间异质性受多种驱动因子的共同作用,单一因子并不能完全解释红树林的空间变化特征。具体地,养殖池面积与其他因子之间的交互作用明显强于其他因子之间的交互作用,其中与人口的交互作用在所有的交互作用中具有最明显的增强效果,且交互作用类型均为双因子增强;红树林距道路最短距离、人口与其他因子的交互作用包含非线性增强和双因子增强2种;此外,单一因子驱动力相对较弱的年均降雨、年均温和地形在与其他因子的交互作用中,影响强度明显增强,且主要呈非线性增强效果。总体而言,不同驱动因子对红树林影响的交互作用均呈现增强的作用,没有出现减弱的情况,增强效果相对明显的有:养殖池面积-人口、养殖池面积-红树林距道路最短距离、人口-红树林距道路最短距离、养殖池面积-年均降雨、养殖池面积-地形、人口-年均温。
图7 1990—2020年东南亚红树林时空变化驱动因子的交互作用探测器结果Fig.7 Interaction of drivers of spatio-temporal changes in mangrove in Southeast Asia from 1990 to 2020
基于GEE云平台EVI数据和东南亚红树林数据集,结合Theil-Sen 中值趋势分析和Mann-Kendall检验方法,探讨东南亚红树林1990-2020年时空变化趋势,并利用地理探测器定量分析其变化的驱动因素,得到的主要结论为:
1)从时空分布特征上,1990-2020 年东南亚红树林面积整体呈下降趋势,在30年间平均每年以1.1%的速度消失。以蜂窝网为基本单位统计的面积变化率表明:红树林完全丧失区域高达40.11%;变化率在-99%~-5%的区域共占37.36%,在-4%~-1%的区域占了13.95%,面积变化增加的区域(变化率在1%~99%)共占8.18%,面积不变的区域即面积变化率为0%的区域仅占0.4%。
2)从时空变化趋势上,东南亚红树林在1990—2020年退化趋势明显,退化的区域明显多于改善的区域。在选择的典型研究区中,高岗(柬埔寨)和斯里巴加湾(文莱)退化趋势最为明显;吉奥(越南)和伊洛瓦底江三角洲(缅甸)地区的红树林退化和改善伴随分布,但整体上退化的区域多于改善的区域;双溪哈吉多拉尼(马来西亚)和宾图尼湾(印度尼西亚)在红树林分布边缘呈退化趋势,但在红树林分布中心有所改善,呈现恢复趋势。
3)从红树林变化驱动因素看:1990—2020 年红树林面积变化的空间异质性受多种驱动因子的共同作用,单一因子对红树林变化的影响程度为养殖池面积>人口>红树林距道路最短距离>年均温>年均降雨>地形,不同因子的相互作用均呈现增强的作用。其中,养殖池面积与其他因子之间的交互作用明显强于其他因子之间的交互作用。
东南亚红树林总面积在1990-2020年呈明显的下降趋势,虽然在相关政策的干预下红树林面积在2015年起有所回升,然而,由于已确定的策略和方法是泛热带的,只针对东南亚红树林最受威胁的部分国家,包括菲律宾、印度尼西亚、越南等(Friess et al., 2016),因此仅在实施红树林保护行动的部分国家,红树林有增加趋势,其余国家红树林面积仍呈严重的下降趋势。因此,应参考红树林面积变化率空间分布,着重对面积严重下降的区域采取相关保护措施。
红树林湿地动态驱动力是一个复杂的系统,既包括自然因素又包括人为因素,其变化是两者共同作用的结果(毛丽君 等,2011)。自然因素包括气温、海温、降水、地形以及红树林自身特性,对红树林面积变化起重要的控制作用(晁碧霄 等,2021);人文因素主要包括滨海地区的城镇化建设及快速增长的海岸带开发活动所导致的红树林退化(Costa et al., 2015)。在本研究中,人为因素(养殖池扩张及人口增长)是导致东南亚红树林面积减少的重要原因,而自然因素(年均温、年均降雨及地形)对红树林的影响相对人类活动所造成的影响较小,这与目前已有的对东南亚各个国家红树林的驱动因素研究(DasGupta et al., 2013; Latiff et al.,2014; Richards et al., 2016; Hai et al., 2020)相符。人为因素在东南亚红树林面积变化中起主导作用,尤其在养殖池大面积建设国家,如印度尼西亚、越南、菲律宾和泰国,养殖池面积的扩张是影响红树林退化的最主要因素(Huitric et al., 2002; Walters et al., 2004; Kusmana et al., 2014; Veettil et al., 2018),此外,人口增长所导致的耕地及建筑用地扩张也是主要的驱动因素。
本研究探讨了东南亚红树林长时间序列的空间变化趋势,并定量分析了其变化的驱动因素,可为生物多样性的保护、沿海经济的发展、生态环境的可持续改善提供重要支撑。但还存在需要深入分析之处:在驱动因素分析中对驱动因子考虑不够全面,在后续研究中将加入地表覆盖数据,从而更加科学完善地探讨东南亚红树林覆盖变化的驱动因素。