Temporal and spatial variation characteristics of annual water yield and its influencing factors in Hangzhou based on InVEST model

浙江省杭州市位于中国东南沿海北部、长江三角洲南沿(29°11′~30°33′ N，118°21′~120°30′ E)，属于亚热带季风气候，2019年全市平均气温约18.1 ℃，年降水量为1 647.7 mm。夏季高温多雨，冬季低温干燥。杭州市地形复杂多样，东部为浙北平原地区，西部为浙西丘陵地带，京杭大运河和钱塘江穿城而过。截至2019年，全市下辖10区2县1市，分别为上城区、下城区、拱墅区、西湖区、江干区、滨江区、萧山区、余杭区、临安区、富阳区、桐庐县、淳安县和建德市。

本研究所采用的数据包括Landsat遥感影像数据、气象观测数据、土壤数据和地形高程数据。采用中国科学院地理科学与资源研究所提供的长三角城市群2000、2005、2010、2015年共4期土地利用数据为基础数据，将地类划分合并为耕地、林地、草地、水体、建设用地和未利用土地等6类。气象数据由国家气象科学数据中心提供的浙江省基准地面气象观测站及自动站共计18个站点的2000—2015年观测数据集，去除异常值后在ArcGIS平台下进行投影，使用空间分析工具中的克里金(Kringing)插值分析得到2000、2005、2010和2015年的降水量、潜在蒸散量分布图，其中潜在蒸散量按Penman-Monteith (PM)公式计算得到，实际蒸散量由InVEST模型产水模块运行计算得到。利用中科院寒区旱区数据中心提供的中国土壤类型图提取杭州市土壤根系深度、物理组成、有机质等信息，在ArcGIS平台栅格数据空间叠加运算下得到土壤有效含水量(AWC)^[24]。采用中国科学院地理所提供的数字高程模型(DEM)计算得到流域矢量数据。所有数据进行栅格化处理，转为同一投影的栅格图像，空间坐标系为Krasovsky_1940_Albers，像素单元大小为30 m×30 m。

利用InVEST模型的产水模块，考虑实际蒸散量，建立年降水量与Budyko热液耦合平衡的关系，估算杭州市不同土地利用类型的产水功能。其核心算法是运用水量平衡法结合气候、地形和土地利用类型计算得出每个栅格的产水量。产水量为区域内各栅格的降水量减去实际蒸散量。产水量计算公式如下：

式(1)中：Y_x为栅格x的年产水量(mm)；T_AEx为栅格x的年实际蒸散量(mm)；P_x为栅格x的降水量(mm)。T_AEx根据Budyko曲线近似得到，公式^[25]为：

式(2)~(5)中：R_x为栅格x的Budyko干燥度指数，是潜在蒸散量(T_PEx)与降水量(P_x)的比值。C_Kx为植被蒸散系数，在不同的植被类型下该值是不同的，表示不同发育期植物潜在蒸散量与参考蒸散量(T_E0x)的比值，具体取值参考模型默认范例。T_AEx通过T_E0x直接计算，并由T_E0x和P_x共同决定其取值。w_x是土壤可利用水量与降水量的比值。根据DONOHUE等^[26]的定义：

式(6)中：C_AWx为栅格x的土壤有效含水量(mm)，Z为季节常数，基于杭州市自然地理特征及相关参考文献将其取值为6.5^[27]。本研究年产水总量为杭州市所有栅格单元产水量的总和。

以城市原点(紫薇园，30°15′00″ N，120°10′12″ E)为城市中心、以5 km为间距，建立覆盖杭州市主城区的10个同心圆梯度带^[28]。为消除因面积不同产生的结果差异，采用ArcGIS空间分析工具Zonal对每个梯度带要素进行空间分布均值统计分析。

选取年降水量(Pre)、年实际蒸散量(AET)和年均气温(Tem)作为影响年产水量的气象因子指标，根据主成分分析(PCA)从密度、面积、边缘、形状、多样性和聚散性6个方面选取斑块密度(PD)、最大斑块指数(LPI)、边缘密度(ED)、景观形状指数(LSI)、香农多样性指数(SHDI)和蔓延度(CONTAG)等6个景观格局指数作为影响年产水量的景观格局指标。对研究区进行1 km×1 km网格化处理并去除边缘值，分别构建气象因子和景观格局指数网格数据集，其中景观格局指数由土地利用类型数据通过Fragstats 4.2软件计算得到。利用SPSS 25中Pearson分析计算得到年产水量与气象因子、景观格局指数的相关性。

由表1可知：2000—2015年间杭州市年产水总量与平均值均呈先下降后上升趋势。就总量来看，2005年最小，为9.81×10⁹ m³，比2000年下降了23.24%，2015年达到最大值(24.48×10⁹ m³)；15 a平均增长率为6.10%。就栅格单元平均值来看，最小值为582.42 mm (2005年)，最大值为1 457.24 mm(2015年)，4 a相对误差为5.90%，2000—2015年平均增长率为6.17%。

从不同土地利用类型来看，2000—2015年间各类型土地栅格单元年产水量均呈先下降后上升趋势(表1)。建设用地高于其他土地利用类型，最大值为1 796.25 mm (2015年)，水体产水量低于其他土地利用类型，最小值为214.30 mm (2005年)。2000—2015年间，各类型土地年产水量增长率从大到小依次为水体(181.72%)、未利用土地(95.98%)、林地(90.10%)、耕地(86.97%)、建设用地(79.42%)、草地(77.78%)。

从空间区域上来看，2000—2015杭州市年产水量呈现东北高、西南低的分布特征 (图1)。高值区主要分布在东北部的主城区，且呈聚集扩张趋势，栅格单元多年平均产水量为1 042.68~1 254.85 mm；其他区县市(富阳区、临安区、桐庐县、淳安县、建德市)相对较小，为933.57~1 009.06 mm (图2)。

Figure 1.  Spatial distribution of water yield in Hangzhou in 2000－2015

Figure 2.  Changes of water yield of counties and districts in Hangzhou

从梯度变化来看，年产水量呈现随着与城市中心距离的增加先增大后减小趋势(表2)，最高值均出现在距离城市中心10 km处。对比距离城市中心5~15 km和40~50 km的年产水量发现：2000、2005、2010和2015年的差异率分别为4.50%、16.50%、9.25%和10.31%，这表明年产水量的城乡变化在2005年前后出现拐点。

杭州市各土地利用类型转移矩阵表明(表3)：2000—2015年间土地利用变化最剧烈的为耕地与建设用地之间的转化，耕地面积减少了2.75%，建设用地面积增加了3.67%，15 a间有3.24%耕地转化为建设用地。结合表1可知：杭州市2000—2015年间建设用地面积比例增加显著，导致杭州市2015年产水量显著大于2000年产水量。

从景观格局来看，2000—2015年间杭州市存在景观明显破碎化和异质性(表4)。15 a间，PD先下降后上升，LPI下降1.23、ED上升1.07，表明景观破碎化程度明显；LSI上升3.69，表明斑块形状较不规则、复杂度不断增加；CONTAG下降1.77，表明景观连通性处于下降趋势，反映景观破碎化增加；SHDI由0.95上升至1.00，表明景观异质性增加。从梯度变化来看，随着与城市中心距离的增加，LPI、CONTAG总体呈现先减小后逐步增加的趋势，最大斑块面积和景观连通性先减小后增大；ED、LSI、PD、SHDI总体呈现先增加后逐步减小的趋势，表明随着与城市中心距离的增加，景观破碎化和异质性先增大后减小(图3)。其中，CONTAG和SHDI的拐点出现在10 km处，表明在距离城市中心10 km处，景观连通性最小、异质性最大；PD、ED、LSI的拐点出现在20 km处，表明在距离城市中心20 km处，景观破碎化最为严重；LPI的拐点则在15和30 km处，表明在距离城市中心15和30 km处，最大斑块面积减小导致一定程度的景观破碎化。因此，在距离城市中心10~20 km处存在最为显著的景观破碎化和异质性。结合年产水量发现：单位面积产水量拐点出现在距离城市中心10 km处，与景观格局的梯度变化趋势较为一致。

Figure 3.  Changes of landscape pattern indices along gradient in Hangzhou　　　　

Figure 4.  Changes of climate factors along gradient in Hangzhou

从时间上看(表5)：2000—2015年间，杭州市年降水量最小的是2005年，最大的是2015年。与2000年相比，2005年降水量减少183.54 mm (13.90%)，2015年增加696.21 mm (52.75%)。年均气温最低的是2000年，最高的是2015年，相差0.31 ℃。年均实际蒸散量最小值和最大值分别出现在2005和2010年；与2000年相比，2005年实际蒸散量减少9.21 mm (1.64%)，2010年增加2.31 mm (0.41%)。结合表1可知：降水量变化与年产水量趋势一致。从空间上看，各县区气象因子存在差异(表6)。主城区与富阳区年降水量为1 540~1 545 mm，其他4区县市均高于1545 mm。主城区年均气温为17.35~17.39 ℃，其他5区县市为17.34~17.38 ℃。主城区年实际蒸散量明显低于其他5区县市，其中最低值出现在下城区(289.10 mm)，最高值出现在淳安县(617.13 mm)。从城乡梯度来看，在距城市中心50 km内，年降水量变化不大；随着距离增加，年均气温逐渐下降，年实际蒸散量先减小后增大，最小值出现在距离10 km处(图4)。结合2000—2015年产水量的空间分布特征和梯度变化趋势可知：年产水量空间变化与年实际蒸散量表现相似。

由表7可知：2000—2015年，杭州市年产水量与降水量相关性最强，相关系数为0.959，与气温呈极显著正相关(P＜0.01)，与实际蒸散量呈极显著负相关(P＜0.01)。在景观格局指数方面，年产水量与PD相关性最强，相关系数为0.153，与ED、LSI、SHDI呈极显著正相关(P＜0.01)，与LPI呈显著负相关(P＜0.05)，与CONTAG呈极显著负相关(P＜0.01)。

InVEST模型水量平衡模型原理认为：降水量和实际蒸散量是估算产水量的2个重要参数，其中降水量作为气象因子的主要变量，由自然条件决定，受人为因素影响较小。也有研究^[15]认为：景观格局比气象因子对年产水量的影响更大。本研究发现：2000—2015年杭州市年产水量与降水量的相关性大于其与指数，与顾晋饴等^[13]、窦攀烽等^[14]结果一致。气温与年产水量具有一定的相关性，而实际蒸散量由于计算方法受到气象因子、土地利用和景观格局的综合影响，相关系数较小。年产水量与各景观格局指数相关系数均较小，但显著性较强，说明虽然景观格局不是影响杭州市年产水量变化最重要的影响因素，但仍值得关注。提示不同地区或生态系统间影响产水量的决定性因素应该是不同的。

研究发现：2000—2015年间，杭州市年产水量总体上表现为先降低后增加，建设用地的年产水量高于其他各土地利用类型。建设用地是研究区变化规模和速度最快的土地利用类型，由不透水性的下垫面组成，植被覆盖少，不利于有效拦截大量降水，同时地表蒸散发能力较弱，人类活动更为剧烈；建设用地的增加一定程度上引起产水量的增加^[11]。相比之下，水体单位面积产水量最小。研究区水域广布(千岛湖、西湖、西溪湿地、钱塘江)，水体下渗程度小于其他土地利用类型，但由于水面蒸散速度与强度高^[29]，因此水体的年产水量也最小。这也可能是研究区生态系统质量整体提升所带来的效应^[30]。林地通过林冠层、枯枝落叶层和土壤层拦截降雨，并通过蒸腾作用散失水分。快速城市化导致研究区部分高质量城市绿地损失，景观破碎化和异质性显著增大，因此林地提供的产水量不稳定^[15]；而为建设绿色城市、生态城市新增的城市林地大部分处于幼龄阶段，这可能是林地产水量在2000—2015年间增长率大于耕地和草地的原因，也说明了城市林地恢复及提供优质服务需要较长时间，具有较大的提升空间。

由主城区梯度变化可知：景观格局在距离城市中心10~20 km处破碎化程度尤其严重，表明2000—2015年间此区域是人类活动最为密切、景观破坏最为严重的区域。距离城市中心10 km处的单位面积产水量最大，表明这是年产水量在梯度上变化的拐点；气象因子在该区域内变化不大，表明景观格局的变化将重新分配城市产水服务，对城市生态平衡和可持续发展产生积极或消极影响^[31]。该结果提示：快速发展的城市化进程中，合理布局规划城市用地^[32]，可避免因大规模集中建设导致的生态系统服务功能的降低。

目前InVEST模型产水模块已被广泛应用，与邱问心等^[33]研究相比，本研究利用Penman-Monteith 公式计算年蒸散量具有一定的可行性。《杭州市水资源公报》数据显示：杭州市2000、2005、2010和2015年实际产水深度均值分别为0.63、0.54、1.13和1.43 m。本研究结果对比验证发现：4个年份相对误差为5.90%，说明模拟结果相对可信。但其中部分输入参数的选择仍会影响结果的准确性，选择不同的经验公式、空间插值方式和土地利用类型分类均会得到不同的结果。本研究采用研究年份的当年气象数据进行模型估算，未考虑到气象因子作用于环境过程的时间滞后性^[34]；对景观格局的分析仅考虑了景观水平，对斑块和类型水平的研究有待进一步研究。未来可通过评估其他指标与景观格局的空间相关性来进一步研究生态系统服务与城乡景观格局的关系，为实现区域可持续规划提供科学依据。