随着生态系统服务功能研究的迅速发展, 生态系统水源涵养功能受到越来越多的关注, 并成为国内外研究的热点之一^[1-3]。水源涵养是陆地生态系统重要生态服务功能之一, 是生态系统在一定的时空范围和条件下, 将水分保持在系统内的过程和能力, 在多种因素的作用下(如生态系统类型、地形、海拔、土壤、气象等)具有复杂性和动态性特征^[4]。生态系统水源涵养功能的内涵随着人们对生态系统与水关系认识的不断深入而变化^[5]。从早期对水源涵养功能的研究单指生态系统对河流水量的影响, 到现在全面的理解水源涵养功能的各种表现形式, 包括生态系统的拦蓄降水、调节径流、影响降雨量、净化水质等, 人们越来越认识到生态系统水源涵养功能的重要性。

国内外学者均对流域产水量及水源涵养功能进行了广泛的研究, 但是研究的侧重点各不相同。国内学者主要是根据研究需求, 通过不同的计算方法对流域生态系统水源涵养功能及其时空变化进行研究, 如丁程锋等^[6]利用InVEST模型, 结合气候、土壤、地形等因素, 通过假设云杉分布的不同情景, 对乌鲁木齐河流域水源涵养功进行了定量评估和分析；刘璐璐等^[7]根据综合蓄水能力法, 结合森林资源二类调查数据, 分析了南北盘江流域森林生态系统的水源涵养功能及其时空分布格局；而国外学者则更注重研究影响流域产水量和水源涵养功能的水文、气候要素, 如Pessacg等^[8]研究了降水因子的不确定性对位于干旱地区的阿根廷丘布特河流域产水量的影响；Marquèsa等^[9]运用InVEST模型对西班牙东北部的Francolí流域产水量进行了评估, 并分析了气候及降水模式的变化对该区域水生态系统服务所产生的影响。

石羊河流域位于我国西北干旱区生态环境脆弱地带, 而水资源是制约该区域经济发展、生态安全的关键因素^[10], 其亦是退耕还林还草较早实施的区域。因此, 研究该流域在土地利用变化情境下水源涵养功能和空间差异, 对区域水土保持治理、合理配置水资源、促进生态系统管理和保护、科学实施退耕还林还草工程具有重要的现实意义。

许多专家和学者对石羊河流域的水资源问题已进行了一些的研究, 但大多数都是进行定性分析^[11-13]。InVEST(Integrated Valuation of Ecosystem Services and Tradeoffs)模型可以评价当前状态和未来情景下生态系统服务的量和价值^[14], 能够基于不同情景(如在不同土地覆被情景下)定量评估流域尺度产水量以及水源涵养功能^[15-16], 近年来, 已经被越来越多的应用于流域水源涵养评估中, 但该模型还需要有针对性的改进和完善, 才能使其更适用于特定研究区域^[17]。流域尺度水源涵养功能研究目前多集中于森林植被^[18-20], 但流域生态系统植被类型并不仅仅是单一的森林植被, 草地、耕地等地类都是流域的重要组成部分, 对流域产水量及水源涵养功能都产生着重要影响^[21]。故本文运用InVEST模型, 定量评估石羊河流域产水量、水源涵养量及各地类水源涵养能力, 分析其空间差异及原因, 通过对模型相关参数的本地化修正来增强模型的区域适用性, 并进行水源涵养重要性空间分级, 模拟不同退耕情景下石羊河流域产水量及水源涵养量, 同时通过情景模拟后结果的差异, 基于最小模糊度法, 确定合理退耕还林还草空间分布模式, 以期为石羊河流域退耕还林还草工程、水资源保护和合理利用提供理论依据。

1 研究区概况

石羊河流域位于西北干旱区东端, 河西走廊东部, 祁连山北麓, 介于101°22′—104°16′E、36°29′—39°27′N间, 流域面积4.16万km²(图 1)。石羊河起源于南部祁连山, 消失于巴丹吉林和腾格里沙漠之间的民勤盆地北部, 主要由大靖河、古浪河、黄羊河、杂木河、金塔河、西营河、东大河、西大河这8条河流组成, 多年平均径流量15.75×10⁸ m³, 产流面积1.1×10⁴ km²。多年平均降水量222 mm, 人均水资源量约700 m³, 耕地亩均水资源量约250 m³, 属典型的资源型缺水地区。石羊河流域内地貌类型复杂, 水热资源的分布不均。石羊河流域是我国内陆河流域中人口最密集、水土资源开发利用程度最高的流域之一, 也是构建国家生态保护屏障的重要地段^[22]。由于特殊的自然地理条件和不合理的人类活动致使石羊河流域出现了地下水天然补给量少、区域性地下水位下降、土地沙化、盐渍化加重、植被退化等一系列生态环境问题^[23-25]。

2 数据与研究方法 2.1 数据来源及预处理

本研究所需的主要数据有：研究区2015年土地利用/覆被类型数据、土壤数据、2005—2015年气象数据、地形数据等, 具体数据来源及预处理见表 1。

2.2 研究方法 2.2.1 产水量模拟

InVEST产水量模块是基于水量平衡原理, 用各栅格的降水量减去实际蒸散发后的水量, 进而得到该栅格产水量。具体计算公式如下：

(1)

式中, Y_{x, j}为j类土地利用/覆被类型、栅格x的产水量；AET_{x, j}为j类土地利用/覆被类型、栅格x的年实际蒸散量；P_x为栅格x中的年降水量。是布德科曲线(Budyko curve)的近似值, 其计算公式如下^[26]：

(2)

式中, R_{x, j}是土地利用/覆被类型j、栅格x处的布德科干燥度指数, 它是潜在蒸散与降水量的比值。ω_x是改进的、无量纲的植被可利用水量与年预期降水量；用于描述自然的气候-土壤属性的非物理参数。其计算方法如下：

(3)

式中, AWC_x是植被可利用的体积含水量(单位mm), 其值由土壤质地和有效土壤深度决定。Zhang系数是表征降雨分布和深度的季节性参数, 其取值范围在1—10之间, 对于夏季降雨为主的季风区流域Zhang系数接近于1^[16]。

布德科干燥度指数R_{x, j}的计算公式如下：

(4)

式中, ETo_x是栅格x内的潜在蒸散量, k_{x, j}表示植被的蒸散系数。

2.2.2 水源涵养计算

用InVEST模型计算产水量之后, 根据DEM计算地形指数D, 考虑不同土壤的渗透性、不同土地利用类型地表径流流速系数的影响, 再计算水源涵养量, 计算公式如下：

(5)

式中, WR为多年平均水源涵养量(mm)；D为地形指数, 无量纲；K_soil为土壤饱和导水率(cm/d)；V为流速系数, Y为产水量。

(6)

式中, Watershed pixel count为集水区栅格数量, 无量纲；Soil depth为土壤深度(mm)；Percent slope为百分比坡度。

2.2.3 输入参数及本地化处理

模型需要输入的参数(图 2)及本地化处理过程如(表 2)所示。

3 结果与分析 3.1 产水量空间格局分析

由InVEST模型产水量公式(1), 计算出石羊河流域及上游8条子流域的产水深度和产水量总值(图 4)。结果表明：研究区单元平均产水深度值为60.90 mm/a, 产水总量为24.71×10⁸ m³/a。

研究区单元产水深度分布具有明显的空间差异性, 与降水量南多北少的分布情况一致(图 2), 与蒸散量南低北高的分布情况相反, 总体呈现出南高北低的趋势, 出现这种空间分布格局主要是由于垂直地带分异因素造成的, 南部是石羊河流域上游区域, 位于祁连山北部, 主要是山地地貌, 属山区河流, 平均海拔3000 m左右, 植被以森林、草地为主, 单元平均降水量达251.50 mm, 比全流域平均值高44.18 mm, 单元平均实际蒸散量为119 mm, 比全流域平均实值低26 mm；上游8条子流域单元平均产水深度分别为黄羊河127.56 mm、杂木河126.51 mm、大靖河123.58 mm、金塔河121.42 mm、古浪河109.08 mm、西营河105.03 mm、东大河99.94 mm、西大河92.16 mm。位于中部和北部的石羊河流域中下游地区, 海拔较低, 植被以耕地为主, 平均降雨量为175.63 mm, 比全流域平均值低31.76 mm, 单元平均实际蒸散量为176 mm, 比全流域平均实值高31 mm；与研究区南部相比, 中部及北部降水量偏小、蒸发量偏高, 故产水量相对较小。

3.2 水源涵养功能总体分析 3.2.1 水源涵养功能空间格局分析

基于产水量模拟结果, 再根据公式(5)计算得到研究区水源涵养量。得出, 研究区单元平均水源涵养量为23.97 mm, 平均水源涵养总量为1.35×10⁸ m³。由图 4可知, 研究区单元平均水源涵养量空间分布格局与单元产水深度空间分布格局基本一致, 也呈现出南高北低的趋势。

水源涵养高值区主要集中在石羊河流域南部上游区域, 该区域产水量大, 主要植被森林、草地的郁密度高, 形成良好的土壤结构及通风状况, 其土壤下渗、持水能力较强^[33-34], 水土保持功能较好, 水源涵养功能强。研究区上游8条子流域单元水源涵养量由高到低分别为杂木河62.11 mm、黄羊河56.11 mm、大靖河54.52 mm、金塔河50.74 mm、古浪河48.90 mm、西营河48.71 mm、西大河42.59 mm、东大河41.79 mm。水源涵养总量由高到低分别为古浪河2.22×10⁷ m³、西营河1.87×10⁷ m³、大靖河1.69×10⁷ m³、东大河1.64×10⁷ m³、西大河1.61×10⁷ m³、杂木河1.23×10⁷ m³、黄羊河0.89×10⁷ m³、金塔河0.84×10⁷ m³。石羊河流域上游8条子流域水源涵养功能相差不大, 与整个流域相比, 水源涵养功能较强。

水源涵养功能低值区位于研究区中部及北部区域, 该区域地势相对平坦, 是人类活动集中区, 耕地分布较多, 东北部地区受腾格里和巴丹吉林沙漠影响, 大部分土壤沙化严重, 水源涵养量较低；且该区域自然植被面积少, 未利用地面积大、分布广, 降水截留能力差, 故总体水源涵养能力较差；但该区域散布着少量水源涵养较高的区域, 结合土地覆被类型图(图 2)可以看出, 这些区域属于林地分布区。

3.2.2 不同地类的水源涵养功能分析

利用ArcGIS空间统计分析工具, 分别计算出不同地类单元水源涵养量和水源涵养总量, 结果见表 3。

各地类的水源涵养能力由其单元水源涵养量来体现。研究区不同地类单元水源涵养量依次为：林地＞草地＞耕地＞建设用地＞未利用地＞水域。林地的土壤孔隙度较大, 能够明显改善土壤结构, 水源涵养能力最强；建设用地、未利用地、水域水源涵养能力均很低。水源涵养总量依次为：草地＞林地＞耕地＞未利用地＞建设用地＞水域, 草地的水源涵养能力较强, 且所占面积较大, 因此水源涵养总量最大；虽然未利用地的单元水源涵养量与水域、建设用地相差不大, 但由于其所占面积很大, 因此水源涵养总量远大于水域和建设用地, 可见水源涵养总量的大小与各种地类面积有很大关系^[35]。

3.3 水源涵养功能重要性空间分级

对区域水源涵养功能重要性进行分级是合理分配、利用和保护水资源的前提。而定量准确的评估区域水源涵养功能又是科学的对水源涵养功能进行重要性空间分级的前提。

利用InVEST模型对研究区水源涵养功能进行定量评估, 再根据《国家生态保护红线—生态功能红线划定技术指南(试行)》中涵养水源重要生态功能区的分级方法, 使用ArcGIS软件中的分位数分类方法, 对石羊河流域2015年水源涵养功能重要性进行5级分类。水源涵养重要性评价能够做到定量化和空间化, 识别的重要区更为精细。

根据水源涵养功能的重要性分为5个重要性等级(表 4)分别为极重要、高度重要、中等重要、较重要及一般重要。

结合图 1、图 5和表 4可知, 研究区水源涵养功能重要性分级结果总体呈现出由西南向东北递减的阶梯分布趋势, 其中一般重要区域所占面积最大, 比例为53.94%。主要分布在中部和东北部的石羊河中下游区域, 包括民勤县、金昌市、永昌县和武威市中北部地区；极重要和高度重要区域所占比例分别为10.41%和11.99%, 主要分布于西南部的石羊河上游区域, 包括肃南县和天祝县；水源涵养功能较重要、中等重要区域, 所占比例为11.77%和11.89%, 主要分布于研究区下游和中上游的过度地带, 包括永昌县南部与肃南县接壤的地区、武威市南部与天祝县接壤地区以及古浪县南部大部分地区。

研究区水源涵养重要性分级与水源涵养功能空间分布格局基本一致, 水源涵养极重要和高度重要地区降水较多、蒸发较小, 分布着8条子流域、地表水系相对发达, 并且植被覆盖度高, 是石羊河流域水源涵养的重要生态区。根据分级结果, 可将极重要和高度重要区域作为重点和优先保护的地区, 并依此来制定生态保护、水资源开发利用等政策。

3.4 基于水源涵养功能的退耕还林还草模式

退耕还林还草是推进生态文明建设、实现可持续发展的重要举措, 其具体规划实施时需要考虑社会、经济、生态等多种因素, 是一项长期复杂的工程。为了更好的研究不同退耕还林还草模式的植被效应对区域水源涵养功能产生的不同影响, 本研究利用InVEST模型分别模拟不同土地利用变化情境下水源涵养量的差异, 通过最小模糊度法得出合理的退耕还林还草模式。

3.4.1 情景模拟

本研究设定了3种变化情景：

情景1：研究区2015年现状实际情况(图 5)；

情景2：退耕还林情景(将研究区内所有耕地纳为林地)；

情景3：退耕还草情景(将研究区内所有耕地纳为草地)。

在以上3种情景设置中, 仅改变土地利用覆被类型, 其他数据均与现状保持一致。

为了在统一标准下, 对3种情景水源涵养量进行对比分析, 采用水源涵养功能重要性空间分级的分类方法, 进行分级, 结果如图 6所示。

分别计算3种情景下水源涵养总量及单元水源涵养量, 结果如表 5所示。

不同情景下的单元产水量顺序为：情景1＞情景3＞情景2。可以看出, 随着林地、草地面积的增加, 产水量有所减少, 且退耕还林情景产水量最少。因为森林植被对土壤水分需求量大, 蒸散发量也大于草地和农作物, 进而导致土壤水分亏缺, 故退耕还林情景产水量最低。不同情景下的单元水源涵养量和水源涵养总量顺序均为：情景2＞情景3＞情景1。在退耕还林情景下, 单元水源涵养量增加到34.17 mm, 此时林地面积最大, 水源涵养总量增加到1.93×10⁸ m³；在退耕还草情景下, 单元水源涵养量增加到26.88 mm, 此时草地面积最大, 水源涵养总量增加到1.52×10⁸ m³。与前文不同地类单元水源涵养功能的分析结果：林地＞草地＞耕地相吻合。

3.4.2 退耕还林还草模式

将情景1、2、3的结果进行空间叠置, 可以发现, 在研究区不同退耕方式对水源涵养的影响不尽相同。有些区域退耕还林对其原水源涵养量的提升很大, 且明显高于退耕还草的水源涵养量, 这种情况即可视为“适宜退耕还林区”；而有些区域退耕还林与退耕还草对其原水源涵养量的提升相差不大, 考虑到退耕还林的成本高于退耕还草, 故把这些区域视为“适宜退耕还草区”；还有些区域无论是退耕还林、还是退耕还草与其原本耕地的水源涵养量相比均无太大差别, 故将这些区域视为“保留耕地区”。

由于“退耕还林”、“退耕还草”、“保留耕地”3种类型的概念本质上是模糊的, 其水源涵养量差异值是有交叠的, 很难从中明确判断出其所属类型, 因此, 本研究通过最小模糊度法来计算其隶属函数, 从而确定3种类型的分类条件。用情景2与情景1进行栅格计算, 在结果中抽取30个样点, 并将情景1、情景2、情景3及情景3与情景1的栅格计算结果等作为判断其所属类别的先验知识, 样本点的选取均匀的分布在流域的上中下游区域, 结果如表 6所示。

按照统计学的观点, 评价指标统计数据的均值最能反映这个指标的平均特征, 偏离均值越大的值, 其所反映的特征偏离平均特征越大, 于是以情景2与情景1水源涵养的差值作为变量来描述耕地像元适合“退耕还林”“退耕还草”“保留耕地”这3个模糊概念时, 图 7所示的是像元所属类型的三角形隶属函数, 需要进一步确定待定参数t₁和t₂。

然而, t₁和t₂都是可变的, 其中t₁∈(14.0, 36.4), t₂∈(36.4, 66.3)(14.0、36.4、66.3分别是3种类型情景2与情景1水源涵养差值的平均值), 因此, 3个模糊集的模糊度达到最小值时, 参数t₁和t₂取值最佳。

设论域X={x₁, x₂, …, x₃₀}={3.8, 4.6, …, 99.6}, 表示“保留耕地”、“退耕还草”、“退耕还林”这3个模糊概念的模糊集分别为A、B、C, 选择模糊熵作为模糊度的度量, 则可以建立如下模型：

求解上述模型：t₁=32.8, t₂=65.0, 得到描述模糊概念“保留耕地”、“退耕还林”和“退耕还草”的3个模糊集如图 8所示。

通过最小模糊度的分类条件, 可以确定研究区内适宜保留耕地的区域、适宜退耕还草的区域以及适宜退耕还林的区域。分布模式如图 9所示, 可以看出, 适宜保留耕地的区域主要分布于石羊河中下游地区, 该区域地势平坦, 降水量小、蒸发量大, 单元水源涵养量小, 无论是退耕还林还是退耕还草对其水源涵养量提升均不大, 且就实际情况来说, 该区域人口密集, 是主要农业区, 将其耕地保留也较为合理；适宜退耕还草的区域分布在石羊河上游和中下游的过度地带, 该区域水源涵养能力较差, 虽然有一定的提升空间, 但退耕还林与退耕还草对其水源涵养量的提升差异不大, 综合考虑到退耕还林、还草的经济成本及政策因素(退耕还林国家要支付每公顷补助22500元, 而退耕还草为每公顷15000元), 故将其划为适宜退耕还草的区域；适宜退耕还林的地区位于石羊河流域南部, 该区域降水较多, 海拔较高, 适宜林地生长, 退耕还林后的水源涵养量明显高于现状及退耕还草后的结果, 故这些区域属于适宜退耕还林区。

4 结论与讨论 4.1 结论

(1) 基于InVEST模型定量评估了石羊河流域产水量及水源涵养功能, 并分析了其空间差异, 得出：模型中的Zhang系数为2.1时, 研究区产水量模拟效果最佳；研究区2015年单元平均产水深度值为60.90 mm, 产水总量为24.71×10⁸ m³；单元平均水源涵养量为23.97 mm, 水源涵养总量为1.35×10⁸ m³；单元产水深度与单元平均水源涵养量空间分布格局基本一致, 呈现出南高北低的趋势；研究区不同地类单元水源涵养能力以林地最强, 草地的水源涵养总量最大。

(2) 在对研究区水源涵养功能进行定量评估的基础上, 进行了水源涵养功能重要性空间分级, 得出：一般重要区域面积最大, 占总面积的53.94%, 主要位于研究区北部的石羊河中下游地区, 极重要和高度重要区域共占总面积的22.40%, 主要位于研究区南部, 即石羊河流域上游地区。

(3) 对不同退耕模式进行情景模拟, 得出：研究区内耕地全部纳为林地时, 单元水源涵养量及水源涵养总量均最高；研究区内耕地全部纳为林地或草地时, 产水量均有所减少, 且全部纳为林地时, 产水量最少；根据最小模糊度法, 基于水源涵养功能角度, 确定了退耕还林还草的合理模式。

4.2 讨论

(1) 虽然本文基于InVEST模型定量评价了石羊河流域水源涵养功能, 并对其进行重要性分级, 但是由于流域生态系统的复杂性, 在未来的研究中还需要进一步定量考虑流域上、中、下游的水资源关系, 定量分析中下游对上游区域水资源的需求, 分层次、分等级的研究流域水资源分配和生态补偿等问题。

(2) 通过本文研究可以发现, 随着林地、草地面积的增加, 位于干旱半干旱地区的研究区产水量有所减少, 会对下游绿洲用水带来问题, 并且可能会影响植被生长发育, 甚至导致群落衰败和生态系统退化, 因此, 应遵循自然演替规律, 重点关注区域土壤水管理, 科学合理的制定退耕还林还草政策。