研究区位于南渡江、昌化江、万泉河3大河流的上游，地理位置18°42′35′′—19°14′16′′N、109°12′48′′—109°48′42′′ E，涉及五指山、琼中、白沙、昌江4个县（市），总面积1 780 km²，约占国家公园总面积的41.7%。其中，包含的国家公园五指山分区和黎母山分区面积分别约为537和541km²，包含的鹦哥岭分区中西部区域面积为701 km²。

该区山地、丘陵地形广布，山体陡峭，地形起伏较大，海拔为180～1 867 m，其中，五指山为海南最高峰，主峰海拔达1 867 m。该区年均气温在17～24℃，为热带海洋性季风气候，年降水量大，雨水充沛，年均降水量1 800 mm以上，但时空分布不均，干湿季节明显。土壤以热带红壤和砖红壤为主，局部高海拔山地分布有山地赤红壤和山地黄壤。植被类型以热带天然林为主体，是海南热带雨林国家公园天然林保存最完整、连片分布程度最高的区域之一。

本研究参照邓雅文等^[9]提出的基于多特征指标和层次聚类分析的河源区自动划分方法，结合海南3大河流（南渡江、昌化江、万泉河）的水系特征和海南热带雨林国家公园管理边界，首次明确界定海南“三江源”区域。首先，基于高程数据，采用均值变点法确定最佳汇流累积量阈值，提取子流域单元作为河源区划分的基本单元。其次，从水文条件（汇流累积量、河网密度）、地形地貌（高程、坡度）和植被条件（植被覆盖度、土地利用类型）3个维度构建多特征指标体系，对各子流域单元的指标均值进行极差标准化处理。然后，采用ward聚类分析法进行层次聚类分析，聚类数设置为2～9，将8次聚类结果中被判定为河源区的频次≥50%的子流域单元进行组合，得到初步划定的海南“三江源”区域。最后，将初步划定结果与海南热带雨林国家公园管理边界（五指山分区、黎母山分区、鹦哥岭分区）进行叠加分析和人工目视修正，确保3大河流源头均包含在内且保持分水岭边界的完整性。最终确定的海南“三江源”区域总面积1 780 km²，包含黎母山全境、五指山全境和鹦哥岭中西部（图1）。

图  1  海南热带雨林国家公园范围及“三江源”区域图

Figure 1.  National Park of Hainan Tropical Rainforest and the “Sanjiangyuan” of Hainan Island

本研究采用2013、2018、2023年3年的多源数据开展研究。构建InVEST模型所需输入数据包括年降水量、年潜在蒸散发、植物根系限制层深度、植物可利用水分含量、土地利用类型等。根据模型计算得到逐年产水量，再利用地形指数、流速系数和土壤饱和导水率数据，对逐年产水量进行修正得到逐年水源涵养量。模型计算中用到的所有栅格数据利用Arcgis10.8软件统一重采样至空间分辨率为30 m，使用投影坐标系Krasovsky_1940_Albers。此外，本研究使用的数据还包括归一化植被指数（NDVI）、温度、高程等，主要用于水源涵养量驱动因素分析（表1）。

本研究所采用的土地利用数据来源于地球资源数据云平台（www.gis5g.com），空间分辨率为30 m，根据LUCC分类体系将土地利用类型分为6个一级类，25个二级类。采用二级分类体系，一级分类（6类）基于土地利用属性划分，包括耕地、林地、草地、水域、建设用地和未利用土地；二级分类（25类）依据土地覆被自然属性细分（耕地分为水田/旱地，林地分为有林地/灌木林等）。该体系由中国科学院地理科学与资源研究所刘纪远^[12]提出，作为中国土地变化研究的基准框架，其数据已应用于国土规划、生态评估及全球变化研究，分类精度通过人机交互验证确保可靠性，在适用性方面具有重要的现实意义。

根据海南水资源公报（https://swt.hainan.gov.cn）及相关文献研究成果^[3,4,8,13]，当Z取值为12时，InVEST模型产水模块的模拟效果好。

根据InVEST模型用户指南、联合国粮农组织（FAO）的《作物蒸散量作物需水量计算指南》，以及参考相关文献的研究成果^[3,4,8,13]，InVEST模型输入所需的生物物理系数见表2。

图  2  水源涵养模型参数

Figure 2.  Water conservation model parameters

采用InVEST模型中的产水模块对海南“三江源”区域产水量进行评估。InVEST模型产水模块基于水量平衡原理，认为栅格单元的降水量减去实际蒸散发后的水量即为产水量^[3]。模型相关算法如下：

式中，$ {Y}_{xj} $为j类土地利用类型栅格单元x的年产水量，单位mm；$ {AET}_{xj} $为j类土地利用类型栅格单元x的年实际蒸散量，单位mm；$ {P}_{x} $为栅格单元x的年降水量，单位mm。

基于Budyko的水热耦合平衡假设^[14]计算水量平衡的蒸散部分。

式中，$ {w}_{x} $为描述自然气候与土壤性质的非物理参数；$ {AWC}_{x} $为栅格单元x的植物有效利用水分含量，单位mm；Z为Zhang系数，又称季节常数，表示研究区降水特征的常数，取值范围1～30；$ {R}_{xj} $为j类土地利用类型栅格单元x的Budyko干燥指数，表示潜在蒸发量与降水量的比值；$ {k}_{xj} $为j类土地利用类型栅格单元x的植被蒸散系数；$ {ET}_{0x} $为栅格单元x的潜在蒸散量，单位mm；$ MaxsoilDept\text{h} $为最大土壤深度，单位mm；$ RootDept\text{h} $为根系深度，单位mm；$ {PAWC}_{x} $为栅格单元x的植物可利用水分含量，单位mm，取值范围0～1。

基于世界土壤数据库，根据周文佐^[15]提出的经验估算模型计算PAWC。

式中，$ sand $为土壤砂粒含量（%）；$ silt $为土壤粉粒含量（%）；$ clay $为土壤黏粒含量（%）；$ OM $为土壤有机质含量（%）。

在产水量的基础上，通过地形指数、地表流速系数、土壤饱和导水率对产水量进行修正得到水源涵养深度（water retention）^[4]。相关算法如下。

式中，$ {Retention}_{xj} $为j类土地利用类型栅格单元x的水源涵养深度，单位mm；$ {Vel}_{xj} $为j类土地利用类型栅格单元x的流速系数；$ {{{K}_{sat}}}_{x} $为栅格单元x的土壤饱和导水率，单位m·d^-1；$ {TI}_{x} $为栅格单元x的地形指数；$ drainagearea $为集水区栅格数量；$ {soil}_{{{dept\text{h}}_{x}}} $为栅格单元x的土壤深度，单位mm；$ {slope}_{x} $为栅格单元x的百分比坡度。

土壤饱和导水率（K_sat）采用Cosby土壤传递公式^[16]计算。

在ArcGIS软件中，根据每个栅格的水源涵养深度（mm），使用分区统计工具按海南“三江源”区域及其分区（五指山、黎母山、鹦哥岭）面积分别计算水源涵养量。水源涵养量（m³）计算公式：

式中，$ {V}_{x} $为栅格单元x的水源涵养量，单位m³；$ {Retention}_{x} $为栅格单元x的水源涵养深度，单位mm；A为单个栅格单元的面积（本研究采用空间分辨率为30 m的数据，栅格面积为900 m²）。

本研究以InVEST模型模拟得到的水源涵养量为因变量，采用最优参数地理探测器方法探究水源涵养量的空间分异特征及其主要驱动因素，最优参数地理探测器主要基于R语言软件里的GD包来实现。

地理探测器是探测空间分异性，以及定量揭示其背后驱动因素的一种统计学方法，包含4个探测器：分异及因子探测、交互作用探测、风险区探测以及生态探测^[10]。最优参数地理探测器是基于传统地理探测器改进得到的模型，它优化了空间数据离散化过程和空间分析尺度，确定了地理探测器模型的最佳参数组合，显著提升了空间分析的准确性和有效性^[11]。本研究采用最优参数地理探测器计算所有连续型数据在不同分级方式下以及不同间断数下的值，使用自然断点分类、分位数分类、相等间隔法、几何间隔分类和标准差分类5种分类方式，选择q值最高的参数组合（分类方法与分类数）进行离散化，q值最高的组合即为驱动因素的最优空间离散化，从而探测研究区水源涵养量的空间分异性和影响因子对水源涵养量变化的解释力。

在参数优化的基础上，本研究主要使用模型中的因子探测器和交互作用探测器。通过因子探测分析单一影响因素对研究区水源涵养量的解释程度，通过交互作用探测识别双因素共同作用时对研究区水源涵养量的解释程度，解释力的强弱通过比较q值的大小来反映。双因子交互作用的强弱及类型详见表3。具体计算公式：

式中，$ q $为某影响因子对水源涵养量空间异质性的影响力，取值[0,1]，$ q $值越大，说明解释程度越强；$ \text{h=1,\cdots,}L $为变量Y或因子X的分层（Strata），即分类或分区；$ N $和$ {N}_{\text{h}} $分别为全区和层h的单元数；$ \sigma _{\text{h}}^{2} $和$ {\sigma }^{2} $分别是层h和全区的Y值的方差。

2013—2023年，以5年为一个时间段，通过InVEST模型评估得到海南“三江源”区域的水源涵养量。综合降水量和潜在蒸散量变化来看，总体而言，海南“三江源”区域2013、2018和2023年3期的降水量和水源涵养量都呈现“先下降，后小幅回升”的变化特征，潜在蒸散量整体变化幅度较小，降水量和潜在蒸散量整体具有相反的空间分布特征（图3）。降水量在2013年最高（1 860 mm），2018年下降（1 250 mm），2023年（1 366 mm）有所回升。潜在蒸散量整体变化幅度较小，2013年最低（1 233 mm），2018年略升（1 282 mm），2023年最高（1 304 mm）。

图  3  海南“三江源”降水量、潜在蒸散量和水源涵养深度时空分布特征

Figure 3.  Spatiotemporal distribution characteristics of precipitation, potential evapotranspiration and water retention of the“Sanjiangyuan” of Hainan Island

海南“三江源”区域2013、2018和2023年的年均水源涵养量分别为1 005、463和551 mm，水源涵养量分别为19.4×10⁸、8.97×10⁸、10.66×10⁸ m³。2013到2018年水源涵养量降幅为53.8%,2018到2023年升幅为18.8%，但仍明显低于2013年水平（图3、表4）。五指山、鹦哥岭和黎母山各区的水源涵养量及单位面积水源涵养量变化情况亦呈现类似趋势。

空间分布上，研究区水源涵养量和降水量高的地方都集中分布于南渡江、昌化江和万泉河等主要河流源头的中高海拔山区，即海拔高、蒸散水平低、降水条件良好的鹦哥岭的西南部、黎母山东北部及五指山。相比之下，鹦哥岭的西北部较低（图3）。

基于最优参数地理探测器对研究区水源涵养驱动因素进行探测，定量评估各驱动因子对研究区水源涵养量变化的解释力。选取降水量（Pre）、潜在蒸散发（Pet）、温度（Tem）、归一化植被指数（NDVI）、高程（DEM）、坡度（Slope）、土地利用类型（LULC）这7个因子，分别计算2013、2018和2023年这些因子对水源涵养量空间分异的决定值q，从而揭示不同年份各因子的影响强度差异。显著性检验结果表明，大部分驱动因子的P值均低于0.05，通过显著性检验，表明所选取的驱动因子对研究区水源涵养量的空间分异具有一定的解释力。但是，在2013年，NDVI、坡度和土地利用类型的P值均大于0.05，未通过显著性检验，说明该年份植被覆盖、地形起伏及土地利用类型对水源涵养空间分异的影响尚不显著。

各因子在不同年份对水源涵养量空间分异的解释力不同，整体来看，从大到小排序为降水量>高程>温度>潜在蒸散发>土地利用类型>NDVI>坡度（表5、图4）。降水量的解释力最强，对水源涵养量的空间分异起主导作用，q值为0.525～0.654；其次为高程，q值为0.498～0.621；温度和潜在蒸散发两个因子的q值分别为0.440～0.593和0.416～0.541。土地利用类型、NDVI和坡度的解释力较弱，q值分别为0.088～0.107、0.030～0.076和0.030～0.035（表5）。

图  4  海南“三江源”水源涵养量单因子探测结果

Figure 4.  Single factor detection of water conservation of the “Sanjiangyuan” area of Hainan Island

交互作用探测结果显示，所有驱动因子之间的交互作用均表现为双因子增强或非线性增强，表明因子的交互作用均可增强对水源涵养量空间分异的解释力（图5）。降水量与高程交互作用的解释力最强（0.779<q<0.804），其次为降水量与温度的交互作用（0.768<q<0.770）以及降水量与潜在蒸散发的交互作用（0.758<q<0.770）。土地利用类型、NDVI和坡度单因子对水源涵养量空间分异的解释力较小（0.030<q<0.107），但分别与降水量、潜在蒸散发、温度和高程交互时解释力显著增强（0.431<q<0.786），其中，最强的是土地利用类型分别与降水量、潜在蒸散发、温度和高程的交互（0.553<q<0.786）。

图  5  海南“三江源”水源涵养量影响因子交互作用探测结果

Figure 5.  Detection of factors influencing water conservation of the “Sanjiangyuan” area of Hainan Island

通过OPGD模型揭示了土地利用类型（LULC）是影响水源涵养空间分异的重要交互因子。为进一步理解LULC的具体变化过程及其对水源涵养量的潜在影响，本节对研究期内（2013—2018年，2018—2023年）的土地利用类型转移特征进行分析。

整体来看，2013—2018年土地利用各地类之间存在转换，但变化量较小，大部分面积保持不变，林地为研究区主体土地利用类型（表6）。林地转出面积达到1 567.68 km²，其中除1 552.28 km²保持林地外，仍有13.18 km²转为草地、1.44 km²转为耕地，同时存在向水域与建设用地的小范围转移，分别为0.34 和0.44 km²，说明林地在维持主体地位的同时仍存在向周边地类的小幅退化趋势（表6）。草地主要向林地转移，转移面积为13.69 km²，同时仅有少量转向耕地、水域和建设用地（均小于0.2 km²）（表6）。耕地、水域、建设用地的转移量都较小，大部分面积保持不变（表6）。值得注意的是，2018年建设用地虽在区域内占比极低，但其面积在2013—2018年期间呈现增长趋势，由2013年的0.47 km²增加至2018年的0.95 km²，增长幅度超过一倍（表6）。

2018—2023年，总体而言，土地利用类型变化以林地和草地之间的微幅转换为主，向建设用地的转移规模较小，土地利用类型变化整体趋于稳定（表7）。

综上所述，2013～2018年间，研究区土地利用变化主要表现为林地的少量转出和建设用地的扩张，这与同期水源涵养量的下降趋势在时间上吻合，暗示了人类活动增强可能加剧了水源涵养功能的退化；而2018～2023年土地利用趋于稳定，为国家公园建设后水源涵养功能的恢复提供了有利的下垫面条件。

海南“三江源”区域是整个海南的水源涵养核心区域，维系着整个海南水资源的生态平衡与安全，定量评估其水源涵养功能并探明水源涵养功能驱动因素，对掌握海南饮用水安全、推动全省绿色可持续发展具有关键作用。

本研究结果表明，在时间尺度上，2013—2023年间海南“三江源”区域水源涵养量整体呈现出“先下降、后小幅回升”的阶段性变化特征，该变化特征与海南热带雨林国家公园尺度和海南主要流域尺度的相关研究结果相似^[3,7-8]。主要原因是2013年的降水量明显高于2018年和2023年，这种降水量的时段差异是驱动研究区水源涵养量变化的主导因素。作为典型的热带季风气候区，海南降水总量大但不同时期变化不稳定，降水的变化直接影响区域水源涵养量的变化起伏。自2019年海南热带雨林国家公园试点建设以来，植被覆盖度与生态系统完整性逐步提升，海南“三江源”区域水源涵养量呈现小幅回升趋势，体现了国家公园建设在促进水源涵养功能恢复方面的积极效应。

在空间尺度上，水源涵养功能表现出高值集聚特征。水源涵养高值区主要集中分布于南渡江、昌化江和万泉河等主要河流源头的中高海拔山区，整体呈现出海拔较高、降水条件优越且蒸散发水平较低的典型特征。类似的空间分异规律已在海南热带雨林国家公园及其他山地河流源头区的水源涵养评估研究中得到验证^[3,17]。从水文过程机理的角度分析，高海拔与复杂地形对区域内潮湿气流的抬升有积极作用，使得该地区降水量极为丰富，为水源涵养功能的形成提供了充足的水量基础。同时，高海拔地区气温相对较低，潜在蒸散发受到明显抑制，在高植被覆盖和稳定森林结构的共同作用下，降水更易通过下渗和土壤蓄存方式被有效保留，从而增强了区域水分调蓄与水源涵养能力。相关研究表明，高海拔地区较低的潜在蒸散发和有利的植被结构是其水源涵养功能高值形成的重要水文机制^[18]。

最优参数地理探测器结果在统计意义上为前述基于时空格局和水文过程的成因解释提供了定量支持。单因子探测结果表明，降水量在研究区水源涵养空间分异中具有最高的解释力，是其主要控制因子。高程次之，潜在蒸散发和温度的解释力亦保持在较高水平。同时，双因子探测结果表明，降水量与高程交互作用的解释力最强，其次为降水量与温度的交互作用以及降水量与潜在蒸散发的交互作用。即在热带季风气候和复杂山地地形共同作用的背景下，研究区水源涵养空间分异总体表现为以降水条件为主导、地形因子起强化作用，同时受潜在蒸散发和温度等热量条件共同调节的控制特征。类似的水量服务驱动机制已在季风气候控制下的山地地区相关研究中得到验证^[19-20]。

尽管土地利用类型、NDVI 和坡度在单因子探测中的解释力较弱，但其与降水量、潜在蒸散发、温度及高程交互后解释力显著增强。在热带季风气候控制下的山地河流源区，水源涵养首先受降水输入及能量条件约束，并在地形梯度作用下形成基本空间格局，而土地利用类型和植被状态主要通过调节蒸散过程、入渗能力和水分再分配路径，对气候—地形效应起放大或削弱作用。类似的“弱单因子—强交互”特征在相关研究中得到证实^[21]，强调了土地利用类型与植被因子的调节作用。

需说明的是，本研究未直接纳入人类活动强度指标。一方面，土地利用类型在一定程度上可作为人类活动的间接表征，其交互作用结果表明，土地利用方式变化可能通过影响植被结构和下垫面条件，对水源涵养功能产生间接调节作用。另一方面，研究区尺度相对较小，样本数量有限，进一步增加解释变量可能削弱地理探测器结果的稳健性与统计可靠性。因此，本研究重点关注气候、地形等自然因子及土地利用因子对水源涵养功能的影响。

研究发现，海南“三江源”区域水源涵养能力具有时段差异性和空间差异性，为提升海南“三江源”区域水源涵养服务，提出如下对策建议：

1）加强河流源头区生态空间的保护与优化布局，划定并严守水源涵养生态保护红线。研究发现，海南“三江源”区域水源涵养高值区主要集中分布于中高海拔山地源头区，其空间格局受降水输入与地形条件共同控制，因此，水源涵养能力提升应以源头区生态空间为优先保护对象，实施严格管控。

2）增强水源涵养功能对气候变化的适应与调节能力。研究发现，海南“三江源”区域水源涵养功能主要受降水输入及潜在蒸散发、温度等气候因子控制。在全球气候变化大背景下，气候变暖驱动水循环加速，极端降水事件增多且蒸散发显著增强，打破了降水、蒸发、径流与地下水补给之间的动态平衡，从而削弱了水源涵养功能的稳定性^[22]。因此，未来对于海南“三江源”区域水源涵养功能的优化提升必须充分考虑气候变化的影响。

基于自然的解决方案（Nature-based Solutions，NbS）强调通过保护和恢复生态系统自身的调节功能，以自然过程应对气候变化和水安全等社会挑战^[23]。在该理念的指导下，建议加强源头区森林生态系统的连续性与完整性，通过实施天然林保护，维持连片分布的森林覆盖和完整的林分结构，增强森林对降水的自然调蓄作用，从而提升水源涵养功能在气候变化背景下的稳定性与适应能力。此外，可结合地形条件布局小尺度、分散式的生态蓄水单元（如适当移动山间石块和倒木，增加沟谷水潭的数量及其容水量，降低流速和增加局域水量维持能力），提高径流滞留和洪峰削减能力，既留住雨季水量，也补给枯水期用水，从而提升局域水源涵养与水文系统稳定性。

3）健全水资源生态补偿机制，构建支撑水源涵养功能稳定提升的长效激励体系。海南“三江源”区域作为典型的热带山地河流源头区，承担着全岛重要的水源涵养与水安全保障功能，其水源涵养服务具有显著的空间外溢性。受生态保护红线管控和天然林保护等措施约束，源头区在土地利用、产业发展和资源开发等方面，面临着较高的机会成本。而水源涵养服务带来的供水保障、防洪减灾和水质改善等综合收益则主要由下游城镇与社会经济系统共享。在气候变化背景下，极端降水事件增多与蒸散发增强进一步放大了源头区生态系统调蓄功能的重要性，也对持续投入生态保护提出了更高要求。因此，有必要通过制度化的水资源生态补偿机制，将水源涵养服务的受益关系转化为稳定的经济激励，以保障以源头区生态保护和天然林保护为核心的水源涵养提升措施得以长期有效实施。

本研究以海南“三江源”区域为研究对象，基于2013、2018和2023年3期数据，运用InVEST模型与最优参数地理探测器（OPGD），系统评估了水源涵养量的时空演变特征及其主要驱动因素。

1）2013、2018和2023年，海南“三江源”区域水源涵养量由19.4×10⁸ m³降至8.97×10⁸ m³，后回升至10.66×10⁸ m³。

2）水源涵养高值区主要集中在鹦哥岭西南部、黎母山及五指山降水条件优越且蒸散发水平较低的中高海拔山区。

3）降水量是水源涵养空间分异的主导因子（q值为0.525～0.654），其次为高程（0.498～0.621）、温度（0.440～0.593）和潜在蒸散发（0.416～0.541）。交互作用探测表明，降水量与高程交互作用的解释力最强（q值为0.779～0.804）。

4）水源涵养量的提升需聚焦中高海拔源头区生态保护与气候适应性管理，加强天然林保护与植被结构和沟谷水潭优化，加强对极端气候的缓冲能力；同时完善水资源生态补偿机制，为水源涵养能力持续提升提供制度保障。