热带雨林区水源涵养服务时空演变及驱动因素分析：以海南热带雨林国家公园为例

杜乐山，刘海鸥，刘文慧，王槐睿，张 颖，全占军*

1.北京林业大学经济管理学院，北京 100083

2.中国环境科学研究院生态研究所，北京 100012

3.海南省辐射环境监测站，海南 海口 571126

生态系统服务是指人类从生态系统中获取的惠益，包括供给服务、调节服务、文化服务和支持服务四大类[1].水源涵养服务作为生态系统调节服务的重要组成部分，直接关系到生态系统的韧性和社会经济可持续发展[2].如何科学识别水源涵养服务重要区域，量化水源涵养服务及其时空演变规律[2-4]，一直是学术界讨论与研究的热点.早期主要基于观测数据定量分析水源涵养服务，研究方法包括水量平衡法、降水贮存法、年径流法、林冠截留剩余法、土壤蓄水能力法等，这类方法往往需要大量的定点观测数据，仅适用于小尺度研究；近年来大量基于模型的水源涵养服务量化方法脱颖而出，如SWAT 模型[5-6]、InVEST 模型[2,7-8]、LPJ 动态模型[9]等，特别是InVEST 模型所需的数据易获取、参数调整灵活、评价结果可空间化表达，得到了学术界的广泛认可[10-11].值得注意的是，不同学者界定的水源涵养服务内涵、选取的评估模型与使用的本底数据不同，往往导致水源涵养服务估算结果出现较大差异[12]，因此，基于不同生态系统类型，采用不同方法与评估模型得出的结果往往不具备可比性，但对同一区域开展水源涵养服务时空动态演变分析，对系统掌握区域水源涵养服务情况具有重要意义[13-14].

近年来，海量数据的出现以及机器学习算法的快速发展，极大提升了机器学习在数据分析中的应用，越来越多的机器学习算法在特定领域甚至超过人类，成为众多领域解决现实问题的有效工具.在生态环境领域，支持向量回归模型(support vector regression，SVR)[15-16]、BP 神经网络[17-18]和LSTM[19-20]等机器学习算法已被广泛应用于大气污染预测，随机森林模型[21-22]、神经网络模型[23]也被应用到生态系统服务驱动因素分析并取得了良好效果.在机器学习算法中，基于Adaboost 算法[24-25]和梯度提升树(Gradient Boosting Decision Tree，GBDT)[26]的集成学习法应用则更为广泛；尤其是GBDT，通过逐步训练组合多个决策树来实现分类或回归任务，具有高度可解释性、鲁棒性强、可扩展性好的优点而被优先选用.LightGBM(grandient boosting machine)是实现GBDT 算法的优化版本[27]，具有训练速度快、内存消耗低、准确率高且可以快速处理海量数据等优点，而被广泛认可并应用于健康管理、灾害防控等方面，但在环境领域应用较少[28].

1 研究区概况

海南热带雨林国家公园位于海南岛中部，地理位置108°44 ′32 ″E～110°04 ′43 ″E、18°33 ′16 ″N～19°14′16″N，涉及五指山、琼中、白沙、东方、陵水、昌江、乐东、保亭、万宁9 个县(市)，总面积4 269 km2，约占海南岛总面积的12.1%.国家公园最高点位于五指山，海拔1 867 m；最低点位于吊罗山区域都总河，海拔仅45 m(见图1).该区光照充足，热量条件优越，年均气温在22.5～26.0 ℃之间；年降水量大，雨水充沛，多年平均降水量为1 759 mm，但时空分布不均匀，干湿季节明显；台风活动频数多、强度大、时间长，台风带来的强风、暴雨和风暴潮等对该区域有较大影响.全区森林覆盖率为95.86%，其中，以热带雨林为主体的天然林占76.56%.

图1 海南热带雨林国家公园数字高程Fig.1 Digital elevation model of Hainan Tropical Rain Forest National Park

2 研究方法与数据来源

2.1 水源涵养服务评估

采用InVEST 模型对2000-2020 年海南热带雨林国家公园的产水深度和水源涵养量进行评估.InVEST 模型是由斯坦福大学、大自然保护协会(TNC)与世界自然基金会(WWF) 于2007 年联合开发的，通过地形系数、土壤饱和导水率以及流速系数对产水量进行修正获得水源涵养量[31].InVEST 模型产水模块基于水量平衡原理，认为栅格单元的降水量减去实际蒸散发后的水量即为产水量，在产水量的基础上再考虑土壤厚度、渗透性、地形等因素的影响，计算水源涵养量.模型主要算法如下：

式中：Yil为 栅格单元i中 土地利用类型l的年产水量，mm；AETil为栅格单元i中 土地利用类型l的年均蒸散发量，mm；Pi为栅格单元i的年均降水量，mm.其中，(实际蒸散量与降水量的比值) 是依据Zhang等[32]在Budyko 曲线基础上发展而来，其表达式如下：

式中：wi为栅格单元i中修正的植被年可利用水量与降水量的比值，是描述自然气候与土壤性质的参数；AWCi为栅格单元i的植物有效含水量，由土壤质地和土壤有效深度决定，是田间持水量和萎蔫点之间的差值，mm；Zh为Zhang 系数，是表征降水特征的常数；Ril为栅格单元i中 土地利用类型l的Budyko 干燥指数，表示潜在蒸发量与降水量的比值；pi为栅格单元i的潜在蒸散量，由气候条件决定，mm；kil为 栅格单元i中土地利用类型l的植被蒸散系数，是不同发育期中作物蒸散量(ET)与 潜在蒸散量(ET0)的比值，由植被类型决定.

2.2 时间演变分析

从时间尺度对水源涵养服务的演变规律进行研究，构建水源涵养服务变化指数(water conservation index，WCI) 对水源涵养服务的波动情况进行刻画，用以指征水源涵养服务的相对增益或损失[33].WCI=0表明水源涵养服务没有变化，即无增益无减损；WCI＜0表示有减损；WCI＞0 表示有增益[34].计算公式：

式中：WCIi为栅格单元i水源涵养服务变化指数；WCcur为待估年份的水源涵养量，m3；WChis为初始年份的水源涵养量，m3.

2.3 空间格局分析

空间自相关分析是研究物体在空间变量分布是否具有聚集性和相互作用的重要方法，包含全局空间自相关和局部空间自相关两种[35]，分别揭示空间单元性质与邻近空间上其他属性值在全局空间和局部空间上的相关性[36].采用全局莫兰指数(Moran′s I)描述水源涵养服务的空间自相关特征，当Moran′s I＞0 时，说明WCI 具有正向空间自相关，值越大表示正向自相关性越强；当Moran′s I＜0 时，说明具有负向空间自相关，值越小表示负向自相关性越强；当Moran′s I=0时，说明没有相关性，即呈随机分布[37].Moran′s I 指数计算公式如下：

一个词的意义发生演变，与整个社会的发展分不开。由于词汇是语言三要素中与客观世界联系最为紧密的部分，因此客观世界的发展变化会影响词汇产生发展变化。在这一变化过程中，词义的演变是必然的。

式中：WCi和 WCj分别为栅格单元i和j的水源涵养量，m3；n为研究区样本总量；为区域水源涵养量的平均值，m3；wij为栅格单元i和j基于距离函数邻接关系建立的空间权重矩阵.

空间聚类分析(local indicators of spatial association，Lisa)用于研究区域内水源涵养服务的空间分异规律.Lisa 计算公式[38]：

式中：wij为空间权重矩阵，xi为栅格单元i的属性值，为所有属性值的平均值，n为区域单元的总数.Lisa＞0 表示空间单元的水源涵养服务是高-高值或低-低值的空间聚集；相反，Lisa＜0 表示空间单元的水源涵养服务为高-低值或低-高值的空间聚集[39].在此基础上，利用ArcGIS 10.2 软件空间分析中的热点分析工具(hotspot analysis，基于Getis-Ord Gi*统计指数)，通过计算各斑块之间的得分Z(标准偏差)和概率P，检测具有统计显著性的热点和冷点空间聚集区域，并绘制局部自相关指标Lisa 分布图[40].

2.4 驱动因素分析

为进一步研究海南热带雨林国家公园水源涵养服务的影响因素，参考已有文献[11,39,41-42]，选取自然因素和人类活动两大类共12 个具体指标(见表1).将各栅格数据的空间投影坐标统一为CGCS2000 高斯投影坐标系，通过ArcGIS 10.2 的“值提取至点”工具，将不同分辨率的数据重采样为30 m 分辨率并对齐栅格[11]，对土地利用类型(LULC)指标做类别变量转换，其他指数进行归一化处理.基于ArcGIS 10.2 软件对2000 年、2005 年、2010 年、2015 年和2020 年每年生成2 000 个随机点，5 期共计10 000 个样本，基于Python 语言对数据进行清洗，去除空值后剩余9 135个样本.采用LightGBM 构建决策树，对样本集采用随机且有放回地抽取，并将样本集的80%划分为训练样本、20%作为验证样本，研究各驱动因子对海南热带雨林国家公园水源涵养服务的贡献，并基于Shape 包对模型进行解释性分析.

表1 海南热带雨林国家公园水源涵养服务驱动因素及数据来源Table 1 Driving factors of water conservation and their sources of Hainan Tropical Rain Forest National Park

3 结果与分析

3.1 时间演变分析

2000-2020 年，以5 年为一个时间段，通过InVEST 模型评估得到海南热带雨林国家公园的产水深度和水源涵养量.总体而言，海南热带雨林国家公园产水深度和水源涵养量波动较大，呈现“先下降、后上升”的趋势(见图2 和图3).其中，2000 年、2005 年、2010 年、2015 年和2020 年的产水深度分别为1 264.85、998.72、927.19、950.39 和1 023.28 mm(见图3).2000-2010 年水源涵养量逐期下降，从5.36×109m3降至3.93×109m3，降幅为26.7%；2010-2020 年水源涵养量逐期上升，从3.93×109m3升至4.34×109m3，升幅为10.37%(见图2).

图2 海南热带雨林国家公园2000—2020 年水源涵养量和产水深度的年际变化Fig.2 Annual changes of water conservation and water yield from 2000 to 2020 in Hainan Tropical Rain Forest National Park

图3 海南热带雨林国家公园2000—2020 产水深度的空间分布Fig.3 Distribution of water yield in Hainan Tropical Rain Forest National Park from 2000 to 2020

采用InVEST 模型评估得到2000-2020 年海南热带雨林国家公园的产水深度在927.19～1 264.85 mm 之间，水源涵养量在4.37×109～5.96×109m3之间，水源涵养量非常大，这与该地区是昌化江、万泉河、太阳河、陵水河等河流的发源地相吻合.本研究结果与21 世纪初的部分研究结果[43-45]存在数量级差异，与近期李昂等[29]和张翠萍等[46]的研究结果相近.从表2 可以看出，采用不同方法与模型得到的水源涵养量往往存在较大差异，即使采用同样的模型，使用参数不同也会导致结果存在差异，这是目前学术界面临的普遍问题[55].但这并不影响笔者所得结果的可靠性，笔者通过相同的方法、模型与参数得出的水源涵养服务，在时间和空间上具有可比性，对指导国家公园保护规划编制及综合管理具有参考价值.另外，刘贤词等[53]研究发现南渡江流域水源涵养量在2005-2010 年呈下降趋势，张翠萍等[46]研究发现海南中部山区在1998-2013 年水源涵养量则呈先上升后下降的趋势，韩念龙等[11]研究发现海南岛水源涵养量在1996-2020 年出现下降，李昂等[29]研究发现海南热带雨林国家公园水源涵养量在2000-2018 年先上升后下降(见表2).上述研究均发现海南岛部分区域水源涵养量的时间波动性，与笔者得出海南热带雨林国家公园2000-2020 年水源涵养量呈现先下降后上升的结论部分吻合，水源涵养量除了与降水量密切相关以外，还与21 世纪前10 年自然生态系统面积持续减少[51]、党的十八大以来森林生态系统面积持续增长有直接关系[56].

表2 海南省部分研究区域水源涵养服务相关研究结果对比Table 2 Comparison of research results of water conservation of some research areas in Hainan Province

3.2 空间格局分析

利用ArcGIS 10.2 软件将海南热带雨林国家公园产水深度进行空间可视化，结果表明，产水深度的空间分布格局差异较大，呈现出自北向南逐渐增加、自西向东逐渐增加的分布特点，并且这种空间分布特点在不同年份之间具有较好的一致性(见图3).

海南热带雨林国家公园水源涵养服务总体呈现“东高西低、南高北低”的分布特点，特别是自西向东明显减少，与张翠萍等[46]对海南中部山区的研究结果一致，也与海南岛降水量的分布基本一致，这是因为海南岛属于热带季风海洋性气候，受海南岛中部高海拔的影响，多年平均降水量由东到西递减[11]；实际蒸散发受降水量及土地利用类型的综合影响，林地是海南热带雨林国家公园的主要生态系统类型，国家公园东部和南部是天然林分布的主要区域[46]，植被和土壤截留的降水相对较多，导致区域水源涵养量较高[11]；西部靠近工业园区，天然林相对较少，橡胶园等人工植被覆盖面积较大[46]，植被的实际蒸散发相对较高，导致被截留的降水较少，使其水源涵养量偏低[11].

3.3 自相关分析

综合空间异质性表达能力和运算工作量，确定1 km×1 km 格网作为研究区水源涵养服务驱动分析的空间单元.由图4 可知，2000-2020 年海南热带雨林国家公园水源涵养服务的空间Moran′s I 指数均大于0.93，Z得分大于2.58(均通过1%显著性检验)，表明海南热带雨林国家公园水源涵养服务有较强的空间正相关性，空间聚集特征明显(见表3 和图4).

表3 海南热带雨林国家公园水源涵养服务Moran′s I 指数与冷热点分析Table 3 Moran′s I and hot spot analysis of water conservation in Hainan Tropical Rain Forest National Park

图4 海南热带雨林国家公园2000—2020 水源涵养服务冷热点分布Fig.4 Distribution of hot spot of water conservation in Hainan Tropical Rain Forest National Park from 2000 to 2020

对海南热带雨林国家公园的水源涵养服务进行局部空间差异性分析，通过空间聚类分析图(Lisa)可以明显看出，西部地区出现空间冷点(“低-低”型)，而东部地区出现空间热点(“高-高”型)，且这种分布特征在不同年份之间呈现较好的一致性(见图4).与以往研究结果[11]类似，“高-高”型主要分布在国家公园东部，与水源涵养服务高值区分布一致；“低-低”型主要聚集在国家公园西部，该区域靠近工业园区，叠加降水较少及实际蒸发较高等因素，水源涵养能力较弱；“高-低”型和“低-高”型分布较为零散.另外，2000-2020 年空间热点区域面积基本表现出“先上升后下降而后上升”的趋势，空间冷点区域面积呈现缩小的趋势(见表3)，这与我国近年来逐步加大以国家公园为主体的自然保护地体系保护力度密不可分.

3.4 驱动因素分析

基于LightGBM 模型，对2000-2020 年海南热带雨林国家公园的水源涵养服务进行决策分析，以确定各驱动因子在水源涵养服务方面的贡献.结果显示，各影响因子对水源涵养服务的解释度达98.98%.为进一步解析各驱动因素对水源涵养服务的贡献，将所有采样点通过蜂窝图进行展示，结果如图5 所示.

图5 海南热带雨林国家公园水源涵养服务驱动因素蜂窝图Fig.5 Shap plot of factors driving water conservation in Hainan Tropical Rain Forest National Park

从蜂窝图可以清晰看出，降水量(Pre)、蒸发量(Pet)、土地利用类型(LULC)对海南热带雨林国家公园水源涵养服务的边际贡献度最高，是主要影响因素；水源涵养能力与降水量呈显著正相关，与蒸发量呈显著负相关，同时受土地利用类型影响较为显著.其他因素的边际贡献度从大到小依次为数字高程(DEM)、坡度(Slope)、净初级生产力(NPP)、归一化植被指数(NDVI)、人类足迹指数(Hii)、人口密度(Pop)、温度(Tem)、人均GDP 和夜间灯光(Dmsp).从图5 可以清晰看出，除土地利用类型外，其他人类活动如人类足迹指数(Hii)、人口密度(Pop)、人均GDP、夜间灯光(Dmsp)对水源涵养服务的影响都不大(见图5).

研究发现，自然因素中降水量、蒸发量，以及人类活动中的土地利用类型是海南热带雨林国家公园水源涵养服务的主要驱动因素，这与以往利用地理探测分析得出的降水量、实际蒸散发和土地利用类型对水源涵养时空变化分布具有密切关系[11]的结论一致，也与我国太湖[57]、三江源[58]等其他地区研究结论类似.LighGBM 模型分析同时发现，其他人类活动，如人类足迹指数(Hii)、人口密度(Pop)、人均GDP、夜间灯光(Dmsp)对水源涵养服务的影响不大，一方面是因为自然因素对生态系统服务影响更大[4]，另一方面也说明近年来国家加大对国家公园保护与管控力度，保护区内人类活动稀少导致对水源涵养服务影响不大.本研究引入LightGBM 模型探究海南热带雨林国家公园水源涵养服务的驱动因素，研究结果对国家公园精细化管理具有重要参考价值；但考虑到影响水源涵养服务时空演变的因素错综复杂，除本研究选取的因素外，生态系统演替阶段、林龄[59]以及其他政策制度、保护规划、价值观念等难以量化的因素也会直接影响水源涵养服务[39]，还需要在未来的研究中进一步探讨.

4 结论与建议

4.1 结论

a) 从时间维度来看，2000-2020 年海南热带雨林国家公园水源涵养量波动较大，呈现“先下降、后上升”的趋势，2000 年、2005 年、2010 年、2015 年和2020 年水源涵养量分别为5.36×109、4.23×109、3.93×109、4.03×109、4.34×109m3.

b) 从空间关系来看，2000-2020 年海南热带雨林国家公园水源涵养量的空间格局差异较大，呈现“东高西低、南高北低”的分布特点，且空间分异在不同年份之间具有较好的一致性.总体来看，2000-2020 年Moran′s I 指数均大于0.93，Z得分大于2.58(均通过1%显著性检验)，表明水源涵养服务呈较强的空间聚集模式；西部地区出现空间冷点，而东部地区出现空间热点，且这种分布特征在不同年份之间呈现较好的一致性.

c) 从驱动力分析来看，本文引入LightGBM 机器学习系统研究海南热带雨林国家公园水源涵养服务驱动因素，发现降水量(Pre)、蒸发量(Pet)和土地利用类型(LULC)对水源涵养服务的边际贡献度最高，其中，水源涵养服务与降水量呈正相关，与蒸发量呈负相关，同时受土地利用类型影响较为显著；其他因素的边际贡献度从大到小依次为数字高程(DEM)、坡度(Slope)、净初级生产力(NPP)、归一化植被指数(NDVI)；其他人类活动，如人类足迹指数(Hii)、人口密度(Pop)、人均GDP、夜间灯光(Dmsp)对水源涵养服务的影响不大.

4.2 建议

水是连结生态系统过程与人类活动的重要纽带，水源涵养服务在生态系统服务中占据中心地位.研究发现，海南热带雨林国家公园水源涵养服务具有时间波动性和空间差异性，且受自然因素(主要为降水量和蒸发量)与人类活动(主要为土地利用类型)影响较大.为提升海南热带雨林国家公园水源涵养服务，提出如下建议：

a) 对国家公园进行整体性保护与区域差异化管理.每年海南热带雨林国家公园水源涵养量为3.93×109～5.36×109m3，对保障海南岛饮用水安全和可持续发展至关重要，必须对国家公园进行统筹规划与整体保护，并以国家公园为中心，协同推进构建全岛水生态安全屏障.同时，基于国家公园水源涵养服务的时空异质性特征，在水源涵养量较高的东部地区继续落实严格保护措施，对西部水源涵养量相对较低且出现“低-低”分布特征的地区，建议在不影响自然演替进程的同时适当开展森林抚育措施，以达到国家公园水源涵养最大化的目的.

b) 深入开展气候变化与人类活动监测及其对水源涵养服务的影响研究.国家公园水源涵养服务的时空变化是自然因素(主要为降水量和蒸发量)和人为活动(主要为土地利用类型)共同作用的结果.基于全球气候变化大背景，建议有效利用大数据及人工智能(机器学习)等先进技术，准确把握自然因素波动对水源涵养服务的影响，提前预判水源涵养服务的时空变化趋势并提出保障措施.同时，对国家公园核心区域进行严格保护才能维持源源不断的水源供给，建议充分利用好“五基”协同生态遥感监测体系，集成可见光、红外、雷达、激光等多种传感器[60-61]，全方位做好国家公园人类活动监督管理工作，切实保障海南岛的饮用水安全.

c) 注重方法的科学性与评估的延续性.从科学研究的角度来讲，水源涵养服务评估模型及研究方法多样化更有利于推动学科发展，基于不同研究目的可选用不同方法；但从精细化管理的角度，必须进一步明确研究的规范性、方法的科学性和模型参数的一致性，从而对不同尺度开展时空演变分析，指导保护管理与决策制定.同时，建议对国家公园等重点区域，开展常态化监测和评估以保证结果的延续性，为国家公园科学保护与精细化管理提供决策支撑.