基于小流域的汾河上游生态服务价值演变

王怡帆, 师学义

(中国地质大学(北京) 土地科学技术学院, 北京 100083)

生态系统服务是指由生态系统提供的，有助于人类生存和生活质量提升的各类产品和服务功能的总和。在不同的研究尺度下，人们对生态系统服务的理解不尽相同[1]，综合而看，生态系统服务类型主要分为供给服务、支持服务、调节服务和文化服务4种[2]。生态系统服务为人类和地球上各类生物提供稳定和舒适的生活环境，其功能多样，包括食物生产、水土保持、气候调节、防风固沙等。生态服务价值评估是以货币的形式显化自然资源的经济价值，从短期看，有利于帮助人们正确利用生态系统，从而保障自身宜居的生态环境；从长远看，可以作为政府决策者制定发展战略和衡量发展效益的工具，有助于人类福祉和可持续发展。

20世纪末，Costanza等[3]13位科学家将生态系统服务功能价值进行估算求和。2001年，为了促进生态系统保护和可持续发展，千年生态系统评估(MA)对全球生态系统进行首次多层次综合性评估[4]。国内最具有影响力的是谢高地等[5-7]基于Costanza研究提出的模型，制定出我国生态系统当量因子表，并且于2008年重新修订。2015年谢高地首次选择除NPP外的其他生态指标(降雨量和土壤保持量)对服务功能进行了修正，将单位面积价值当量因子静态评估修订为动态评估，并将其运用到时间和空间上的生态系统服务评估。

本文的研究区域是山西省汾河上游的宁武县和静乐县，区域资源丰富，但在过分追求经济发展的同时，常常伴随着生态环境越发脆弱的问题，地形破碎、水土流失严重、自然灾害现象时有发生，严峻的生态问题不容小觑。2018年以来，山西省为贯彻落实习近平总书记“一定要高度重视汾河的生态环境保护，让这条山西的母亲河水量丰起来、水质好起来、风光美起来”的重要指示精神，着力开始实施以汾河中上游为治理核心的山水林田湖草生态保护修复工程。目前，学界对于生态服务价值的研究，大多基于大尺度区域，如以县、乡镇为评价单元，而本文将汾河上游地区按分水岭和下游河道出口断面划分为数个小流域，以小流域为评价单元，开展汾河上游地区2006年、2011年、2016年生态系统服务价值评估演变与空间动态分析，以此建立起汾河上游重点保护区生态账户，为生态保护修复工程的布局提供依据，有利于最终形成“山青、水绿、林郁、田沃、和美”的生态格局。

1 研究区概况

研究区域位于山西省忻州市中部，地处黄土高原东沿的汾河上游地区，东依原平市和忻府区，西隔管涔山与五寨、岢岚县为邻，南与太原市娄烦县、古交市相连，北靠神池县、朔州市朔城区。地理坐标介于东经111°42′43″—112°36′55″，北纬38°7′57″—39°8′41″，两县境内汾河流域总面积3 441 km2，其中宁武境内1 395 km2，静乐境内2 046 km2，汾河自北向南流经芦芽山和云中山，纵贯宁武、静乐两县。区域山峰高耸、岩石裸露、沟壑纵横、土少石多。宁武、静乐两县涉及28个乡镇845个行政村，区域总土地面积达3 981.29 km2，境内总人口为32.29万人，人口密度达81.10人/km2。

2 试验材料与方法

2.1 数据来源及处理

本文所用数据包括卫星遥感数据、野外实地调查数据以及其他数据来源。其中卫星遥感影像数据来自于“地理空间数据云”平台，2006年、2011年所使用的遥感数据影像由Landset5卫星拍摄，2016年所使用的遥感数据影像由Landset8所拍摄，扫描宽度与Landset5相仿。本文利用ENVI软件对3期卫星遥感影像进行预处理，以汾河上游1∶5万的地形图为基准实现图像几何校正，之后再对其镶嵌和裁剪。利用ArcGIS等相关软件，通过人工目视解译与监督分类相结合的方法，对研究区的遥感影像进行信息分类。根据国家相关分类标准和已有的调查资料，并结合之后研究的需要，拟将汾河上游地区的土地利用类型划分为6大类，分别为耕地、草地、林地、建设用地、水域和未利用地。本文将针对这6种土地类型的生态服务价值演变进行深入的研究与评价。

2.2 研究方法

2.2.1 小流域的划定 国内外研究中基于小流域的研究较多，但普遍只是选取某一个固定小流域开展相关研究，如郑国权等[8]以广东省清远市瑶安小流域为例，对小流域综合治理效益进行定量评价。大多没有系统地将研究区域划分成众多小流域进行研究。人依水而居，本次研究的评价单元之所以选取小流域，是因为小流域已经成为水土流失治理、生态系统修复等研究工作的重要地理单元。在定义上，小流域通常是指集水面积在小于50 km2以下的相对独立和封闭的自然汇水区域，水利上通常指是在一个县属范围内的流域。

本次以小流域为单元进行综合评价和管控，有利于因地制宜改善生态环境，保护水土资源，进而提升人居质量。本次小流域的划分遵循以下标准：(1) 每个小流域的面积控制在3～50 km2；(2) 小流域的边界须经过山顶部、鞍部，须与等高线处处垂直；(3) 小流域可分为完整型、坡面型和区间型3种类型[9]。

2.2.2 系数修正与ESV评估 本文以谢高地等[7]2015年修订的中国生态系统单位面积生态系统服务价值当量表为基础，结合研究区实地情况修正当量表系数，其中修正指标包括社会经济条件、自然资源稀缺度等，由此获得汾河上游生态系统服务价值当量因子表(表1)。

表1 汾河上游生态系统服务价值当量因子

经过测算，粮食单产扣除成本后的价值即为农田的自然产量经济价值，在没有人力投入的自然状态下，农田自然生态系统提供的经济价值是现有粮食平均单价的1/7，公式如下：

En=P×Qn×1/7

(1)

式中：En为1个标准当量因子的生态系统服务价值；P为粮食作物的平均价格；Qn为所有粮食作物的平均产量。

生态服务价值计算公式为

Vij=Ak×VCij

(2)

(3)

式中：Vij为研究区地类j的第i种生态服务功能价值；Ak为地类k的面积；VCij为校正后地类j的第i种生态服务功能价值系数；ESV为总生态服务价值。

2.2.3 敏感性指数 为了衡量生态系统服务价值(ESV)随时间变化对于生态服务功能价值系数(VC)的依赖程度，本文特引入敏感性指数(CS)来验证评估结果的可信程度[10]。如果CS>1，说明价值系数每增加或减少1，生态服务价值将对应变化的幅度大于1，表明ESV具有弹性，评估结果不具有可信度；如果CS<1，说明价值系数的变化不会很大程度地影响生态服务价值的变化，ESV缺乏弹性，评价结果具有可信度。敏感性指数越大，说明其价值系数越准确[11]。本研究通过各地类的价值系数分别增加或减少50%来计算敏感性指数，以此来判断ESV的可信度。敏感性指数(CS)计算公式如下:

(4)

式中：CS为生态敏感性指数；ESVj，ESVi分别为变化后和变化前的生态服务价值总量；VCjk，VCik分别为变化后和变化前的生态服务价值系数；k为各个土地利用类型。

2.2.4 空间自相关与冷热点演化

(1) Global-Moran′sI分析。Moran′sI用于确定某个范围内的空间对象之间是否存在相关性，其数值分布在[-1，1]。若Moran′sI>0表示空间正相关性，其值表明空间正相关性的明显程度，若Moran′sI<0表示空间负相关性，其值越小，则空间差异越大，若Moran′sI=0，则表明空间呈随机性，不具有相关性[12]。空间自相关的Moran′sI统计用如下公式计算：

(5)

(6)

(7)

统计的zi计算公式如下：

(8)

E[I]=-1/(n-1)

(9)

V[I]=E[I2]-E[I]2

(10)

(2) Getis-Ord GeneralG分析。高/低聚类工具与Global-Moran′sI不同，是对于给定的研究区域，可以测量其聚集密度，进而得出高/低聚类空间分布模式[13]，这一工具是一种推论统计。计算公式如下：

(11)

式中：xi，xj为要素i和j的属性值；∀j≠i代表要素i和j不是相同值。

统计的ZG计算公式如下:

(12)

(13)

V[G]=E[G2]-E[G]2

(14)

如果Z得分为正，观察GeneralG指数大于期望GeneralG指数，则表示数据在高值区域聚类；如果Z得分为负，期望GeneralG指数大于观察GeneralG指数则表示数据在低值区域聚类。

(3) 冷热点Getis-Ord Gi*统计。这一工具能够对点或面数据进行分析，提取出具有一定统计显著性的高值(热点)和低值(冷点)在空间上的聚类，从而对相邻区域间的要素或属性的空间关联模式进行测度，在时间变化的冷热程度上反映区域生态服务价值的变化特征[11]。其中z得分越大，表明热点的聚类程度越大，即越紧密。对于统计学上的显著性负z得分，z得分越小，则表明冷点的聚类程度越大。其计算公式如下：

(15)

(16)

(17)

3 结果与分析

3.1 生态服务价值评估结果

3.1.1 生态服务价值系数 通过计算汾河上游2006—2016年粮食作物种植面积与总产量，得出研究区在这一区间内粮食平均单产量为3 150.24 kg/hm2，同时，为消除不同年份价格变化差异，以山西省2016年粮食作物平均价格1.93元/kg为标准，最终求得汾河上游单位面积农田生态系统食物生产服务的经济价值为868.56元/hm2。由此得到的汾河上游生态服务价值系数见表2。

表2 2006-2016年汾河上游生态系统服务价值系数 元/hm2

3.1.2 生态服务价值时空演变 (1) 总体生态系统服务价值变化。根据生态系统服务价值计算公式和2006—2016年3个时相各土地利用类型面积，结合对应的生态服务价值系数，对汾河上游3期不同土地利用类型的生态系统服务价值分别进行计算。由计算可得，2006年、2011年、2016年的生态服务价值总量分别为45.99亿元、40.72亿元、39.67亿元，呈现出持续减少的变化特点。以5年时间尺度为跨度，减少量分别为5.27亿元，1.04亿元，由表3可知，减少率分别为11.47%，2.6%，其中2006—2011年生态服务价值减少程度相当剧烈。2011—2016年减少程度相对缓和(表3)。

表3 2006-2016年汾河上游地区各地类ESV变化

从研究时段不同土地利用类型生态服务价值来看，虽然各地类的价值占比在数值上发生了变化，但整体结构上保持稳定，不同地类生态价值的大小关系保持不变：草地>林地>水域>耕地>未利用地>建设用地。其中，草地生态服务价值连续多年保持第一。在3个时期所占比例分别为37.83%，39.59%，38.1%，而建设用地生态服务价值最少。草地作为汾河上游的优势景观，贡献了大量的生态服务价值，其次是林地。虽然耕地面积也很大，但因其生态服务价值系数较低，综合而言对汾河上游生态服务价值的贡献不是很大。

从不同阶段生态服务价值变化来看，2006—2011年生态服务价值总量减少程度高达11.47%。从各地类角度分别观察，耕地、林地、建设用地、未利用地的生态服务价值均有不同程度的增加,而水域和草地的服务价值均略有减少。所以这一阶段生态服务价值的减少主要是受到水域面积缩减的影响。2011—2016年生态服务价值总量变化幅度较小，共减少了2.75%。其中耕地、林地的生态服务价值增大，其余地类的价值量减小。林地价值的增量最大，为4 951.09万元，草地价值的减量最大，高达10 378.47万元，减少速率最大的是未利用地，其生态服务价值累计减少了66.2%。这一阶段服务价值的增量主要来源于林地，得益于研究区内国家天然林保护、三北防护林、京津冀风沙源综合治理等工程的扎实推进，以及地区的“两山”造林、“两网”绿化、“两区”增绿、“两林”富民、“双百”精品、“双保”管护等省级林业6大工程的有效顺利实施。

从总体空间格局来看(图1)，研究区内各地区生态服务价值的分布和变化呈现出较为显著的空间差异现象，研究区以汾河、东碾河、鸣河为主的8大河流网状分布，东部为管涔山，西部为芦芽山，中部为汾河川，山川地区生态环境好，总体生态服务价值较高，而各乡镇中心周边区域生态服务价值较低，与农田和建设用地占比较大有密切关系。从时间尺度来看ESV分布，研究区西南部生态服务价值减少程度极为显著，主要原因是随着开发建设力度的加大，导致原来ESV高和较高地区的水域面积减少，具有生态服务价值负效应的建设用地面积加大。应加大保护力度，采取生态修复工程来不断提高其生态服务价值。

图1 2006-2016年汾河上游生态系统服务价值等级空间分布

(2) 单项生态系统服务价值变化。由不同生态系统服务功能的价值变化来看(表4)，2006年、2011年、2016年各单项生态服务价值贡献率整体结构分布相对比较稳定。2006年等级次序为：水文调节>气候调节>土壤保持>生物多样性>气体调节>净化环境>美学景观>原料生产>食物生产>水资源供给>维持养分循环，2011年与2016年等级次序相同，相比2006年只是后几项功能顺序发生些许变化。10年来水文调节功能一直是汾河上游的主导生态系统服务功能，价值量均大于90 000万元，而气候调节和土壤保持的贡献率分别位于第二位和第三位，研究时段上述3种服务功能的综合贡献率加和已超过总体服务价值的一半。原因是草地和林地作为汾河上游的优势景观，具有较高的水文调节、气候调节和土壤保持的功能。水资源供给和维持养分循环的贡献率始终处于最后两位，合计贡献率不超过3%，这说明汾河上游调节服务功能远大于其他服务功能。持续增加的功能有食物生产、原料生产和维持养分循环。其中食物生产价值量在单项功能中增加最多，净增607.88万元，研究后期(2011—2016年)增速明显加快，主要因为具有较高食物生产功能的耕地快速增加所致。单项功能水文调节价值量减少最多，研究期内共减少了44 969万元，减少幅度最大的是水资源供给功能，主要原因是具有较高水文调节功能的水域面积的减少以及有水源负效应的建设用地增加。

表4 汾河上游2006-2016年不同生态系统服务功能的价值变化

将各项单项生态系统服务功能聚类成为4大服务，分别是供给服务、调节服务、文化服务和支持服务[14]，得到生态服务类型空间分布图(图2)，由图可得汾河上游水体分布较为广泛，对于内部的小气候改善及生态系统调节服务具有重要影响。供给服务整体价值量增加最多，研究后期(2011—2016年)增速明显加快，主要因为中部地区具有较高食物生产功能的耕地快速增加所致。研究区西南地区文化服务和支持服务等级分布明显降低，在一定程度上反映出人类的社会经济活动极大程度地负向影响着汾河上游生态服务价值。但是由于土地利用方式和面积的分布，单项生态服务功能价值结构分布没有较大变化，保持相对稳定。

3.2 生态敏感性分析

通过公式计算得出各地类的敏感性指数如图3所示，2006—2016年汾河上游敏感性指数变化较小，其中草地的敏感性指数最高，为0.483～0.633，未利用的敏感型指数最低，为0.001。所有地类ESV的敏感性指数在研究时段内都小于1,表明ESV缺乏弹性,研究区的生态系统服务相对于价值系数变化缺乏弹性,本文采用的价值系数符合汾河上游实际情况，计算结果可信度高。

3.3 空间自相关与冷热点演化

3.3.1 空间自相关分析 研究区域2006年、2011年、2016年ESV的Moran′sI值都大于0(Z>0)，表明该区域生态服务价值的空间分布存在一定程度的空间自相关性，但值并不大，表明整体空间自相关性程度不是很显著。随着年份的增长，Moran′sI值先减小后增大，且2016年Moran′sI值达到最大，表明在研究时段中，该区域生态服务价值的空间集聚程度先减小后增大，空间自相关性也存在先减小后增大的发展趋势。

3.3.2 空间分布模式 分析汾河上游生态服务价值的高/低聚类指数可得(表5)，G值与E(G)值均较小且相当接近，且研究时段中，Z得分均为负，E(G)值大于G值，表示数据在低值区域聚类。随着年份的增长，Z(G)值呈增长趋势，说明ESV空间分布及其动态变化的低值聚集程度有所减弱，表明随着时间的变化，ESV整体空间分布形态也发生了一定程度上的变化。

2006—2011年，研究区ESV增值热点区和损失冷热点区域变化显著，整体分布有集中范围扩大的趋势，表明此阶段ESV发生了剧烈的变化，其中位于研究区中南部的石家庄镇、娘子神乡、丰润镇等地由热点区转化为冷点区，说明此区域ESV急剧下降，主要原因是城市扩张、城镇化的推进占用了大量农用地。而ESV热点区集中分布于研究区西北部和东南部地区，主要原因是区域生态环境良好，分布有芦芽山国家级自然保护区、冰洞国家地质公园、天池等，该区域群山苍茫，植被茂盛，人烟稀少，人为活动对该区域干扰较小，土地较少被开发破坏。

图2 汾河上游2006-2016年单项生态系统服务空间分布

图3 生态系统服务价值敏感性指数变化

表5 汾河上游空间自相关与高/低聚类指数

2011—2016年，ESV增值热点区和损失冷热点区的分布在这一研究区间变化不大。ESV热点区主要管涔山地区以及万辉沟，这与小流域尺度的ESV分布趋势保持一致，研究区西南部地区ESV冷点区面积有扩大趋势，大多位于丰润镇、神峪沟乡区域内，这与当地工业化、城镇化的推进密不可分，因此需要重视对冷点区乡镇生态系统的监测工作，采取一系列保护措施，最大程度上控制冷点区范围，遏制其范围扩散。在冷点区内，通过改变地类的方式不断促进生态价值的提高，进而促进汾河上游地区环境与经济统一发展。

4 结 论

(1) 汾河上游地区2006年、2011年、2016年的生态服务价值总量分别为45.99亿元、40.72亿元、39.67亿元，呈现逐年递减的变化特点，2006—2016年共减少了6.31亿元，减少了13.74%，ESV减少的主要原因是水域和草地面积的减小，说明人类活动一定程度上对生态服务价值产生了负向影响。研究时间区域各地类生态服务价值之间的关系保持稳定，为草地>林地>水域>耕地>未利用地>建设用地。各地类生态敏感性指数变化较小，且在研究时段内都小于1，说明生态系统服务价值对价值系数变化缺乏弹性，研究结果具有强可信度。

图4 汾河上游ESV冷热点分布格局

(2) 研究区内各地区ESV空间格局存在较为明显的差异现象，研究区以汾河、东碾河、鸣河为主的8大河流网状分布，东部为管涔山，西部为芦芽山，中部为汾河川，山川总体生态服务价值较高，而各乡镇中心及周边区域因农地和建设用地所占比例较大形成较大范围范围的ESV低值区。因此，为了提高汾河上游整体的生态系统服务价值，应该以小流域为基础，针对各地生境现状，实施生态修复工程来改善生态环境条件。

(3) 单项生态系统服务功能中，水文调节功能的贡献率一直居于首位，价值量均大于9亿元，气候调节和土壤保持的贡献率分别位于第二位和第三位。2006—2016年上述3种服务功能的综合贡献率分别为62.64%，60.81%，60.53%，超过总体服务价值的一半。而生态系统服务价值最低的水资源供给和维持养分循环功能，合计贡献率不超过3%。

(4) 研究时段中区域的ESV存在一定程度的空间自相关性，但相关性程度不是很显著，整体数据在低值区域聚类，随着年份的增长，ESV空间分布及其动态变化的低值聚集程度有所减弱，空间分布在一定程度上发生显著变化。国家天然林保护等工程和“两山”造林、“两网”绿化等省级林业工程的实施，在汾河上游ESV增值热点区中发挥了重要作用，但是研究时段后期ESV损失冷点区面积有扩大趋势，这与当地工业化、城镇化的发展密切有关，这些区域往往是城乡建设土地开发的重要区域，在促进当地经济发展的同时，应当注重生态服务价值的损失，严格控制冷点区范围。

(5) 本文基于遥感影像解译，囿于数据的获取，仅考虑了两个时间区段内3个特定年份生态系统服务价值的变化，对于详细的每年生态系统服务价值变化情况考虑欠缺，无法准确体现渐变型区间段变化情况，这在一定程度上会影响测算结果的精确性。且在对生态服务价值系数进行修正的过程中，结合实地选取的修正指标所需数据不尽充分，在后期研究工作中，应当进行更多的实地调研，收集更加完善的数据。深入分析生态系统服务价值的年际变化特征，结合专家问卷法等获取更加准确的当量因子，并将次研究与生态补偿研究紧密结合起来，健全生态保护的生态补偿机制，实现汾河上游生态效益价值、社会效益价值和经济效益价值的有效统一。