长江流域产水功能对土地利用变化的响应及其驱动因素

胡文敏,杨睿瀚,贾冠宇,殷梓强,李 毅, 沈守云, 李 果,*

1 中南林业科技大学, 长沙 410004 2 湖南省自然保护地风景资源大数据工程技术研究中心, 长沙 410004 3 湖南农业大学, 长沙 410152

随着全球气候变化加剧和人口城市化的迅速增长,全球性水资源矛盾及水安全问题日益突出,已成为限制国家和地区可持续发展的关键问题,其中土地利用/覆被变化对水文过程的影响成为该领域的研究重点[1—2]。由于LULC会影响地表雨水的截留、下渗、蒸发等水文要素及其产汇流过程,进而造成区域产水量在时间过程和空间分布上呈现多样性特征,土地利用变化被认为是影响区域生态系统功能的重要驱动[3—4]。2020年国际生物多样性和生态系统服务政府间科学政策平台发布的《全球生物多样性和生态系统服务评估报告》中指出水资源功能退化是造成生态环境损失和生态系统服务衰退的重要因素,未来随着全球气候变暖,水资源功能退化风险将很难管控,但土地利用作为人类的可控行为,有助于减缓水资源的退化[5—6]。因此了解水资源对土地利用变化的响应机制有重要意义。目前,揭示地表产水功能如何响应土地利用方式,以及两者之间的变化规律尤其是量变过程还存在不足,本研究重点在于揭示产水功能对土地利用变化的响应及其驱动因素。

流域是陆地表层生态系统重要的组成单元,其产水功能作为流域生态系统的关键要素,在整个流域生态系统发育中具有决定作用[7—8]。长江流域是土地利用变化以及生态功能变化研究的热点区域,中国高度重视长江流域的生态环境保护,生态环境部印发的《长江流域水环境质量监测预警办法(试行)》(2018),提出推进长江流域水环境质量持续改善,此外长江流域水生态环境保护作为核心内容被列入“十四五”规划[9]。因此,探讨长江流域土地利用变化对生态系统产水功能的影响,对于深化土地利用/覆被变化的水环境功能效应,健全长江流域水资源的生态安全格局具有重要意义。

从地理学上看,长江流域表现为从河口至高原的跨尺度特征,具有大尺度、跨区域的特点,以长江流域全域为研究对象,解决的是流域的跨尺度研究问题。目前针对长江流域的产水功能研究主要集中在长江流域的局部典型区域,如长江源头[10]、蓄滞洪区[11]、三角洲河口[12]、中下游平原[13]等,虽然这些地域的研究具有典型性,但主要从区域角度考虑,在全域或跨尺度的地理探讨上还存在不足。并且,当前对长江全流域尺度的水资源研究主要涉及径流量变化[14]、水资源量分布特征[15]、水资源综合管理[16]等,对长江流域全尺度的产水功能研究鲜见。国际上研究土地利用/覆被变化的水功能效应主要使用模型法,如InVEST模型[17]、SWAT 模型[18—19]、DHSVM模型[20]等。其中InVEST模型产水模块,主要基于水量平衡原理,根据土地利用、气候及土壤等要素实现对生态系统产水功能的定量、动态及可视化评估,适用于不同尺度、不同区域类型的生态系统产水服务功能的评估,模型在英国典型集水区[21]、洞庭湖湿地[22]、横断山区[23]、太湖流域[24]、三江源[25]、海南岛[26]和武汉市[27]等区域均取得较好的模拟效果。因此,本研究拟采用InVEST模型来探讨长江流域产水量对土地利用变化的响应模式,并结合多元统计方法以期揭示影响区域产水功能的主要驱动因素。

1 数据来源与研究方法

1.1 研究对象

长江是中国第一大河,干流全长约6300km,流域面积约180万km2,占中国国土面积的18.8%。长江流域作为中国重要的水源涵养区,承担了中国约39.2%的水源涵养功能[28],流域地形地貌复杂、土地开发利用程度高,一直以来都是土地利用和环境变化研究的热点,也是重要生态经济发展区,其土地利用变化对流域水环境的影响深刻关系到国家的水生态安全。长江流域聚集了超过全国总量40%的人口和国内生产总值(GDP)[29],近年来,气候变化和社会经济改变了流域的水文条件和水资源时空分配格局,导致部分区域水生态系统受损,水生态安全面临较大威胁。为探讨长江流域的子流域功能变化,研究基于长江水利委对水系的划分,将其划分为12个子流域[30],以宜昌市划分上游和中游,以湖口县划分中游和下游,如图1所示。由于1990年至2015年是长江流域经济的高速增长时期,经历了从“八五”至“十三五”宏观区域土地利用、生态环境、社会发展的剧烈转变,因此,以该时期为研究时段对于揭示产水功能与土地利用变化的关系具有代表性。

图1 研究区域地理分布示意图Fig.1 Geographical distribution map of the study areaHS:长江源头;JS:金沙江水系;MTS:岷沱水系;JLS:嘉陵江水系;WS:乌江水系;UMS:上游干流水系;LDS:洞庭湖水系;HJS:汉江水系;LPS:鄱阳湖水系;MMS:中游干流水系;LMS:下游干流水系; LTS: 太湖水系

1.2 数据源与数据处理

1.2.1数据来源

研究数据源主要包括两个方面,一是计算产水量所需模型参数,二是计算影响产水量变化的主要驱动指标。具体数据包括：年平均降水量(PRE)、年平均潜在蒸散量(PET)、土地利用现状数据、土壤根系最大埋藏深度、植物可利用水量、蒸散发相关系数、流域和子流域矢量分区、季节参数Z、年平均温度(TEMP)、植被覆盖度(NDVI)、人口密度(TPOP)、GDP公里网格、实际蒸散量(AET)。其中,年平均降水量和年平均温度来源于国家气象科学数据中心(http://data.cma.cn/),单位分别为mm和℃,采用了1990年至2015年的月值数据,得到年平均降雨量和年平均温度。年平均潜在蒸散量(PET)同样来源于国家气象科学数据中心,单位为mm,通过月日照时间和月平均气温,结合Ham2公式计算流域的潜在蒸散量[31],公式为：

(1)

式中SD为日照时间(h),Ta为平均气温(℃)。土地利用现状数据来源于中国科学院资源环境科学与数据中心,采用土地利用类型一级划分标准：耕地(AL)、林地(WL)、草地(GL)、水域(WA)、建设用地(CL)和未利用土地(UL)六个类型,见图2。土壤根系最大埋藏深度(mm)、植物可利用水量(mm)来源于国际土壤数据库(HWSD)的土壤参考深度(mm)。季节参数Z借鉴文献中关于长江流域的参考经验值0.5取值[28]。归一化植被指数、人口密度(万人/km2)和GDP公里网格(万人/km2)数据来源于中国科学院资源环境科学与数据中心。实际蒸散量(mm)通过InVEST模型产水量模块计算。

图2 1990—2015年长江流域土地利用类型分布图Fig.2 Distribution map of land use types in the Yangtze River Basin from 1990 to 2015

1.2.2产水量模型

本研究采用InVEST模型的产水量(Water yield)模块计算长江流域的产水量。InVEST模型产水量模块基于Budyko水热耦合平衡假设,由各个栅格单元的降水量减去实际的蒸散发量从而得到该栅格单元的产水量[32],计算公式如下：

(2)

(3)

(4)

(5)

式中,Yxj为j地类的栅格单元x的产水量深度(mm)。AETxj为j地类栅格单元x的年实际增散发量(mm)。Px为栅格单元x的年降水量(mm)。Rxj为潜在蒸散发与降雨量的比值。ωx为气候与土壤属性关系的非物理参数。AWCx为植物可利用水量(mm)。Z即Zhang系数,表示研究区降水量特征。ET0x为栅格单元x的潜在蒸散发。Kxj为某植被类型的蒸散系数。研究计算的平均产水量为平均产水体积(m3),国际上较为通用,主要通过栅格单元的产水量深度乘以像元面积计算,本研究中采用1km×1km像元大小,总产水体积为所有像元平均产水体积的总和。

1.2.3主成分与聚类分析

为分析长江流域产水量驱动因子,研究选取了13个指标进行了探讨分析,利用主成分分析提取了前3个主成分,主要因子包括社会经济、气候变化、土地利用3种因素类型,具体为子流域的耕地、林地、草地、水域、建设用地和未利用地六个地类面积、PRE、PET、AET、TEMP、NDVI、TPOP、GDP,计算采用Origin 9.1软件。为分析长江流域产水量驱动因子的性质和贡献类别,研究结合系统聚类分析方法,将所选取的13个指标进行了指标聚类[33]。系统聚类方法中相似度方法选择类平均法,聚类数为3。此外,为探讨所选因子对产水量的效应,采用了Pearson相关分析和边际效应分析方法。

1.2.4边际效应分析

边际效应分析可以定性定量描述自变量变化对因变量的影响程度,适合在同样协变量因子下揭示自变量因子的梯度水平[34—35]。本研究中对平均产水量的边际效应进行了分析,主要通过各指标对产水量的边际效应体现。分析基于部分依赖(Partial dependence)函数进行计算[34],其结果可以解释某个特征和目标值y的边际效应关系f,具体实现主要通过在R语言中构建平均产水量因子与各指标构建随机森林模型(Random Forest,RF),然后利用部分依赖函数求出特征指标(驱动力指标)的边际效应。随机森林边际效应函数公式如下：

(6)

2 结果与分析

2.1 长江流域产水量时空变化特征

图3显示了长江流域1990年至2015年各时期的平均产水量空间分布特征,按照标准差划分为低、较低、中等、较高、高5个等级,单位：103m3/km2。可以看出：(1)从空间分布上看长江流域平均产水量在空间分布上呈现出明显的西北低、东南高的特点；(2)低等级区呈团状分布,其它等级区域都呈现东北至西南向的带状分布；(3)高等级产水区域主要集中在流域中下游,中等级产水区域主要集中在流域中游,低等级产水区域主要集中在流域上游；(4)从产水量最高值分布上看,最高值出现在2000年,最低值出现在1990年。

图3 1990—2015年长江流域平均产水量等级空间分布特征/(103m3) Fig.3 Spatial distribution characteristics of average water yield grade in the Yangtze River Basin from 1990 to 2015

图4 长江流域各子流域1990至2015年间各期平均产水量特征 Fig.4 Characteristics of average water yield of each sub basin in the Yangtze River Basin from 1990 to 2015

从时间变化上来看(图4和表1)：(1)各期平均产水量排序为：2000年(933.0×103m3/km2)>1995年(892.7×103m3/km2)>2015年(889.9×103m3/km2)>2005年(881.9×103m3/km2)>1990年(867.0×103m3/km2)>2010年(860.4×103m3/km2)；(2)平均产水量变化趋势上,1990至2000年间平均产水量呈升高趋势,2000至2015年间产水量呈先降低再升高趋势；(3)长江流域各期总产水量排序特征为：2000年(1657.0×109m3)>1995年(1585.4×109m3)>2015年(1580.5×109m3)>2005年(1566.2×109m3)>1990年(1539.7×109m3)>2010年(1528.1×109m3),总产水量变化趋势与平均产水量一致。

表1显示了长江流域各子流域1990至2015年间各期平均产水量特征,其中趋势采用各期平均产水量的拟合曲线,可以看出：(1)长江流域上、中、下游的平均产水量均呈现出前10年(1990—2000年)上升,后15年(2000—2015年)先下降再上升的趋势,与流域总体产水量变化趋势一致；(2)从各子流域平均产水量变化上看,上游的长江源头与流域总体的产水量变化趋势差异较大,前5年(1990—1995年)呈下降趋势,而后20年(1995—2015年)呈上升趋势；(3)上游子流域中,乌江水系和上游干流水系平均产水量较高,长江源头最低；(4)中游子流域中,鄱阳湖水系和洞庭湖水系为最高,汉江水系最小,下游子流域差异较小；(3)中下游大部分子流域平均产水量超过1000×103m3/km2,鄱阳湖水系和洞庭湖水系最高,均超过1300×103m3/km2和1100×103m3/km2。因此,不同子流域之间的空间差异性较为明显,且洞庭湖与鄱阳湖是平均产水量的主要贡献流域。

表1 长江流域平均产水量变化特征/(103m3/km2)

表2显示了长江流域各子流域1990—2015年的总产水量变化特征,可以看出(1)总产水量变化趋势与平均产水变化趋势基本一致,但在变化程度上有所差异,如总产水量上金沙江水系、中游干流水系产水量后15年(2000—2015年)下降趋势更明显;洞庭湖水系、鄱阳湖水系前10年(1990—2010年)上升趋势更明显。(2)此外,洞庭湖、鄱阳湖和金沙江三个水系对流域总产水的贡献超过了50%,太湖、乌江和下游干流水系不足15%。

表2 长江流域总产水量变化特征/(109m3)

2.2 土地利用变化对产水量的影响

2.2.1土地利用变化与产水量变化联动特征

土地利用变化与产水量变化关系较为复杂,研究主要从土地利用类型水平上探讨两者的联动关系,如图5所示。可以看出：(1)总体上,长江流域农用地(AL)面积逐年降低,建设用地(CL)和水域(WA)面积逐年升高,林地(WL)、草地(GL)、未利用地(UL)面积变化起伏较大,且农用地和建设用地对平均产水的贡献较大,水域、林地、草地较小；(2)总体体积上,农用地与林地总产水体积最大,为主要贡献地类,建设用地由于面积占比小,虽然平均产水量高,但对总产水体积影响较弱；(3)各土地利用类型上的平均产水量与类型面积变化关系存在较大差异,图中灰色背景区域代表两者变化趋势一致的时期,结果显示农用地上两者的变化趋势在2000—2010年一致,均呈下降趋势和正相关关系,其它时期呈相反变化和负相关关系。建设用地与水域类型上两者变化趋势在1990—2000年和2010—2015年期间一致,呈正相关关系,林地类型上两者变化趋势仅在1990—1995年一致,呈正相关关系,其它时期呈相反变化差异。草地类型上两者变化趋势在1990—2010年一致,呈正相关关系,仅在2010—2015年变化相反,呈负相关系。未利用地类型上,两者变化趋势在1990—2005年一致,其它时期呈负相关。这种平均产水量与类型面积变化的联动关系的差异一方面由于土地利用覆被变化的导致,另一方面与气候变化和人类活动对产水量的影响有密切关系,尤其是产水量在不同流域的时空变化差异。

图5 土地利用类型水平上产水量与类型面积变化的联动特征Fig.5 Linkage characteristics of water yield and type area change at the level of land use type

2.2.2产水量对土地利用变化的响应

图6显示了长江流域不同时期各子流域产水量变化的差异。总体上看：(1)农用地、建设用地、林地和水域类型上的平均产水量变化过程较为一致,在1990—2010年和2010—2015年期间,其产水量的变化呈增长趋势,而在2000—2010年期间,产水量变化呈降低趋势。草地与未利用地类型上与以上土地利用类型差异主要体现在1990—1995年期间,这一时期产水量变化在草地与未利用地类型上呈略微降低状态;(2)从各土地利用类型上的子流域产水量变化与长江流域总体产水量变化对比来看,各时期的变化具有显著差异,土地类型与产水量具有阶段性特征,其中农用地、建设用地的平均产水量变化差异最为显著,其次是林地与水域类型,说明产水量对农用地与建设用地类型的变化响应较为敏感,而林地与水域类型的响应次之,草地与未利用地类型的响应敏感性较弱。

2.2.3土地利用变化与产水量变化的协同关系

图6 产水量对土地利用变化的响应Fig.6 Response of water yield to land use change

图7揭示了研究期间各土地利用类型上产水量变化与类型面积变化的协同关系,可以看出：(1)农用地、草地、未利用地类型上的平均产水量与其面积变化呈正相关,相关系数分别为：7.28、70.04和354.2,说明农用地、草地、未利用地每增加1km2面积,产水量将增加7.28m3、70.04 m3、354.2 m3,具有正向协同效应；(2)而建设用地、林地和水域面积上的平均产水量与其面积变化呈负相关,相关系数分别为：-1.86、-59.00、-265.64,说明建设用地、林地和水域每增加1km2面积,产水量将减少1.86m3、59.00 m3、265.64 m3,具有负向协同效应；(3)在影响程度上,水域与未利用地呈现较强的协同效应,建设用地虽然呈负相关,但相关性较弱,呈较弱的协同效应。

图7 1990—2015年不同土地利用类型上产水量变化与类型面积变化的协同关系Fig.7 Synergistic relationship between water yield change and area change of different land use types from 1990 to 2015

2.3 产水量变化的影响因子综合分析

2.3.1因子变化对产水量的敏感性分析

考虑土地利用变化、气候变化以及人类活动对产水量的影响,进一步探讨不同影响因子对产水量的敏感性,研究统计了长江流域各子流域的平均产水量(WY)、AL、WL、GL、WA、CL、UL的地类面积、平均降水量、平均潜在蒸散量、平均实际蒸散量、平均温度、年平均植被覆盖度、年平均人口密度、年平均生产总值,分析了因子变化对产水量的敏感性。

图8反映了从总体阶段1990—2015年上各因子与产水量变化的关系,可以看出：(1)PRE、PET对产水量在α=0.01水平上具有显著相关性,相关性系数分别达到0.998、-0.481。AET、WL、GL、UL在α=0.05水平上具有较显著相关性,相关性系数分别为-0.262、0.267、-0.289、0.294,说明平均降水量、平均潜在蒸散量两个因子对产水量具有高敏感性,平均实际蒸散量、林地、草地、未利用地对产水量具有较强敏感性；(2)AL、CL、NDVI、TPOP相关性系数均小于0.2,与产水量呈弱相关性,但与其它产水因子有高或较高相关性系数,说明该类因子有可能通过其它产水因子间接作用于产水量,通过合理调控这些变化有可能对产水量产生积极影响。由于因子变化量采用的是各子流域的变化量,其相关性变化反映了长江流域各子流域间的空间差异性,可以看出上述显著性因子中,PRE和PET在0.01水平上表现出较强的显著性,说明不同子流域上PRE和PET的空间差异性较为明显。

图8 1990—2015年长江流域产水量驱动因子敏感性Fig.8 Sensitivity of water yield drivers in the Yangtze River Basin from 1990 to 2015* 表示α=0.05 水平上相关性显著 ** 表示α=0.01 水平上相关性显著;WY：年均产水量；PRE：年均降水量；PET：年均潜在蒸散量；AET：年均实际蒸散量；TEMP：年均温度；NDVI：年均植被覆盖度；TPOP：年均人口密度；GDP：年均生产总值；AL：农用地面积；CL：建设用地； WL：林地；WA：水域；GL：草地；UL：未利用地

图9为产水量与驱动因子之间的边际效应,图中横坐标为各驱动因子的变化范围,纵坐标为效益值f,分析时把其它因子变量作为协变量,效应反映了单一驱动因子变化时与产水量变化的梯度关系,有助于理解驱动因子对产水功能的影响。可以看出随着驱动因子的变化,产水量变化呈现不同的效应模式,土地利用类型上AL、GL、UL呈现正的效应关系,WL、WA、CL呈现负的效应关系。其它驱动因子成分上,PRE、NDVI、GDP和TPOP呈现正的效应关系,PET、AET、TEMP呈现负的效应关系。此外,平均产水量与其驱动因子之间的线性关系表现出不同程度的稳定性,从线性变化结果上看,UL、PET、PRE对产水量的影响呈较为稳定,说明长江流域城市化、蒸散和降雨对产水量的影响较为稳定。而其它因子如AL、WL等呈现复杂线性关系,说明农用地和林地对产水的影响关系较为复杂,这也进一步反映出通过人为作用对农用地、林地等进行结构调整,可以改善区域产水功能,但气候因子如温度升高,会降低地表的产水功能,当人口和经济提升达到一定水平,产水功能有可能提升,反映出通过人为策略对区域水资源进行保护有可能改善区域产水功能。

2.3.2影响产水量主要驱动因素分析

图10显示了利用主成份分析和聚类分析方法探讨影响长江流域产水量的主要驱动因素,主成份分析结果将1990—2015年各期的产水量驱动因子分为三个主成分：PC1、PC2、PC3,聚类分析方法采用平均距离聚类方法。

从图10结果可以看出,PC1的贡献率范围为(28.78%,48.46%),PC2的贡献率范围为(19.41%,28.63%),PC3的贡献率范围为(10.73%,17.13%)。PC1贡献率较大可认为是产水量的主要驱动因素,PC2和PC3贡献率较低可认为产水量的次要驱动因素或间接驱动因素。研究期六个阶段的累积贡献率分别为79.90%(1990—2015年)、70.46%(1990—1995年)、70.98%(1995—2000年),76.44%(2000—2005年)、81.57%(2005—2010年)、72.29%(2010—2015年),综合说明各阶段主成分因子的贡献对产水影响的解释力超过70%。主要驱动因素第一主成分PC1在各个阶段呈现不同贡献强度,1990—1995年PC1主要贡献因子为AL、TPOP、WL和PRE,1995—2000年主要贡献因子发生了变化,WL、CL、GDP、TPOP、EMP占主要因子,2000—2005年主要贡献因子变为WA、CL、GDP、NDVI、PET,2005—2010年主要贡献因子变为AL、PET、EMP、NDVI、WL、UL、AET,2010—2015年主要贡献因子变为GL、PRE、TPOP、GDP、CL、TEMP,各阶段呈现的因子贡献变化说了不同阶段产水量主要驱动力的差异,总的来说1990—2015年期间,PC1主成分主要贡献因子为WA、CL、TPOP、GDP、PRE。同样在次要驱动因素中PC2主成分中,主要贡献因子为PET、AET、GL、UL、WL。PC3主成分对产水的贡献较小。

图10 基于主成分与聚类方法的产水量主要驱动因素分析图Fig.10 Analysis chart of main driving factors of water yield based on principal component and clustering method

3 讨论与结论

3.1 长江流域产水与土地利用时空变化响应模式

通过对长江流域产水与土地利用变化的关系分析,有助于更清楚地识别流域的土地利用变化模式,及其产水响应模式。长江流域土地利用变化与区域城市化、人类活动以及政策有重要联系,研究揭示了流域中建设用地面积逐年增加,呈1990—2000年缓慢增加, 2000—2015年加剧增加态势,且农用地面积变化与建设用地呈相反变化趋势,这种变化与杨阳等[36]、李威等[37]土地利用演变研究规律相似,但较为特殊的是流域林地、草地、水域面积和产水变化具有多样性,林地与草地作为典型植被类型,在长江流域实施重大生态恢复过程中起到重要作用。就林地与草地而言,植树造林、退耕还林、生态林保护、退化草地恢复、退田还湖等对林地和草地的转移产生重要影响,虽然植被恢复提升了流域水资源利用效率,但结果表明林地与草地面积依然呈现轻微降低趋势,因此土地生态恢复政策总体上的作用还不明显,且从土地利用方式上看长江流域依然面临由于土地利用配置的差异引起的产水功能空间分布不均,且流域间异质性程度加剧问题,这就需要未来流域的水资源调控和土地利用配置更为科学合理,尤其对于产水贡献最大的洞庭湖、鄱阳湖、金沙江流域。

流域产水体现了区域水源涵养与人类活动的协调关系,不同土地类型上产水上呈现较大差异,从平均产水量排序上看,建设用地>农用地>林地>水域>草地,其结果与胡砚霞等[38]对汉江流域的研究结果一致,但不同时期的变化幅度具有差异,产水变化特征较为一致的是后15年(2000—2015年),所有土地类型上的产水量呈现逐步增加态势,这期间一方面受全球气候变化影响,气候变暖、降雨、融雪等加速了流域水循环,另一方面土地利用政策导向从注重大力开发到生态保护的变化也产生了积极作用[39],尤其在“十三五”与“十四五”期间推动长江流域水环境保护与经济带环境保护上作用明显,但是由于不同子流域的土地覆被与产水的变化的空间异质性,在具体政策实施中导致土地利用和类型转移方向发生改变,进而引起土地下垫面产水功能的变化[40],通过产水与土地利用类型变化的关系可以看出,农用地变化在长江流域产水功能中起到重要作用,虽然其相关性系数较低,但由于农用地面积占比大,其总体产水量变化较其它用地明显。总体来说,长江流域产水变化对其土地利用变化具有土地类型上的多样性响应模式,农用地、草地、未利用地具有正向协同效应,林地和水域具有负向协同效应,建设用地相关性小,对产水的响应不明显,这些特征在对于长江流域生态系统服务功能维持、生态功能区划、以及不同地理区域的土地利用范式上具有一定参考意义。

3.2 影响长江流域产水量的主要与潜在驱动因素

通过对产水量主要驱动因素的分析为进一步揭示影响流域主要驱动力提供了参考。从研究方法上来说,InVEST模型的产水量模块基于Budyko水热耦合平衡假设,将蒸散量视为水量损失,一般来说降雨量与潜在蒸散量、实际蒸散量具有较高相关性[11, 24],可以视为影响产水的直接驱动因素,但由于长江流域自身的空间跨度大,地形、经纬度、海拔、区位的变化高,导致各地区的气候条件差异性也较大,所以虽然作为直接驱动因素,但在各子流域的变化中呈现差异性,其中潜在蒸散量、实际蒸散量与降雨量都表现为东北低、东南高的趋势,其特征在主成分分析中亦体现,其中降雨量作用贡献较大。主成分PC1可以视为主导因素,其中除降雨量外,土地利用类型的水域、建设用地对产水的作用较大,这也反映了虽然水域与建设用地与总体产水的响应上中并不明显,但由于子流域之间的相互作用起到了总体平衡的作用,总体上间接抵消了类型变化对产水变化的效应,这一结果与王保盛等[41]对闽三角产水量模拟研究结果相似。因此,重视不同子流域之间的作用,合理规划和权衡子流域的功能定位对于调控长江流域产水功能有重要意义,未来可以通过基于不同子流域产水的差异特点,权衡合理的水资源功能区划,进一步提高区域与宏观政策的实施效率。

此外,人口密度、GDP作为PC1主要成分其作用亦较大,说明社会经济与城市化的发展亦是主要驱动因素,这与夏庆云[42]、凌子燕[43]等研究一致,虽然该因素不是直接驱动因子,但其对产水的影响呈现间接式的作用。PC2的贡献虽然较PC1低,但对产水仍具有一定影响,不同的是PC3成分的贡献最低,虽然对产水的作用较小,但通过不同阶段的主成分分析可以看出,农用地、植被覆盖度、温度发生了转变,对产水的影响作用变大,这也说明该类因子可以作为一种潜力因素影响产水功能,虽然它们与平均产水量相关性弱,敏感性小,但在调节产水功能上亦具有潜在功能,如温度在调节降雨和蒸散的同时,对于地表植被生物量以及土壤持水力的变化具有约束作用,从而影响植物[44]和土壤的水分利用[45]。另外从边际效应上看,受温度突变影响,温度上升有可能降低产水的时空差异,且有温度与降雨、蒸散具有负相关关联,因此降雨量和蒸散的减少亦具有相同作用。需要指出的是总体流域与子流域之间亦存在驱动力的差异,如上游、中游、下游之间存在因子贡献的差异,本研究主要是从流域总体上考虑。总之,从产水的主导驱动因子分析来看,长江流域水域、建设用地、人口密度、GDP、降雨可以作为主要驱动力,农用地、植被覆盖、温度可以作为潜在驱动力,主要驱动力和潜在驱动力的区分对于长江流域整体与子流域的宏观和微观政策引导具有指导意义。

3.3 结论

研究基于InVEST模型产水模块探讨了长江流域产水量对土地利用变化的响应模式,探明了流域产水的时空变化特征,揭示了影响流域产水的主要驱动因素,研究结果与讨论表明：(1)长江流域产水量变化存在一定规律,首先在空间分布上呈现出明显的西北低、东南高的特点,且不同子流域之间的空间差异性较为明显,呈现中游>下游>上游的分布,农用地与林地为平均产水主要贡献类型,洞庭湖、鄱阳湖和金沙江为总产水的主要贡献流域,这种差异性对于流域产水功能合理区划有重要意义；(2)流域土地利用变化也存在自身特点,表现为建设用地面积变化呈上升趋势,而农用地呈下降趋势,其它类型呈阶段性和多样性变化,这种变化导致长江流域产水变化对其土地利用变化的响应呈类型上的多样性,从产水量的边际效应上看农用地、草地、未利用地具有正向协同效应,林地和水域具有负向协同效应,建设用地对产水的效应不明显；(3)通过主成分与聚类分析方法探讨的产水的主导驱动因子上看,长江流域建设用地、水域、人口密度、生产总值、降雨为主要驱动力,农用地、植被覆盖、温度为潜在驱动力。平衡不同子流域之间的产水,合理区划子流域的功能定位有助于科学调控长江流域产水功能。研究采用的InVEST模型是基于水热平衡原理,在产水功能中反映的是水资源的纵向功能,不足之处在于未考虑水资源的横向功能,比如区域的径流输送、地下水的传输等。因而未来研究还可以结合水资源横纵方向综合探析其产水功能的效应及其驱动机制。