生态系统服务是指人类直接或间接从生态系统中获得的所有惠益^[1-2], 包括供给服务、调节服务、支持服务和文化服务4个方面^[3], 不同生态系统服务间的关系可具体表现为相互增益协同和此消彼长的权衡关系^[4-6]。国内外大量研究表明, 生态系统服务间权衡与协同关系普遍存在^[7-10], 如供给服务与固碳释氧、水源涵养和土壤保持等支持或调节服务间存在权衡关系, 而支持服务与调节服务间往往表现为协同关系^[11]。值得注意的是, 生态系统服务间的关系并不是一成不变的^[12], 景观结构与人类活动等区域差异往往造成生态系统服务之间关系具有复杂的时空变异特征^[10]。

目前, 学界已经普遍认识到生态系统服务间的关系具有区域异质性^{[10, 13]}。如Chisholm等^[14]在南非的研究发现, 植树造林产生的固碳释氧服务有利于缓解全球气候变化, 但是对产水和生物多样性产生负面的影响, 说明碳固定与淡水供给和生物多样性间存在权衡关系。但钱彩云等^[15]研究表明, 甘肃白龙江流域碳固定与淡水供给间为协同关系。武文欢等^[16]对鄂尔多斯市研究发现产水量与土壤保持存在协同关系, 但Jia等^[17]研究则表明陕西北部产水量与土壤保持间为权衡关系。在当前的研究中, 仍缺少采用统一计量标准, 对不同区域间生态系统服务间关系的对比分析。生态系统服务间关系出现时空变异的内在机理有待进一步明晰, 如何识别生态系统服务权衡关系对土地利用、社会发展等驱动因素的响应及服务间关系变化的关键驱动指标, 已成为生态系统服务研究的前沿议题。

在气候变化与城市化进程的推动下, “一带一路”沿线国家土地利用快速变化, 土地资源亟待优化, “绿色丝绸之路”倡议的提出, 为传播我国生态文明理念并协助“一带一路”区域国家应对气候变化提供生态环保合作的基础框架。然而“一带一路”沿线国家处于不同气候区和不同的社会发展阶段下, 不同类型的土地利用方式、土地利用空间格局和不同程度的土地利用强度对生态系统服务产生了不同的作用^[18], 生态系统服务及其相互关系存在区域异质性, 制订服务于洲际的生态政策和实现土地利用方式的优化, 必须充分把握区域差异规律。因此, 本文评估了“一带一路”区域土壤保持、水产量、碳固定和生境质量四种生态系统服务时空分布规律, 探究在不同生态区、不同国家、不同人口规模下生态系统服务间关系的地域分异规律, 旨在为“一带一路”区域多国土地资源的优化配置和区域生态系统服务的提升提供科学依据。

1 研究方法 1.1 研究区与数据来源 1.1.1 研究区概况

“一带一路”是“丝绸之路经济带”和“21世纪海上丝绸之路”的简称, 空间范围以亚欧大陆为核心, 延伸到非洲和大洋洲。区域气候差异明显, 东部季风气候显著, 大陆性气候分布广泛, 多数国家处于干旱和半干旱地区, 南北跨度大, 跨越了热带、温带和寒带。地形复杂多样, 主要为高原和山地, 平原面积狭小分布在大陆周围, 地势起伏大。在气候与地形的驱动下, 主要的生态区有热带和亚热带湿润阔叶林、温带落叶和混交林、沙漠和沙漠灌木、山地草原和灌木、温带草原和灌木和泰加林等。本文研究区包括中国在内的65个国家(图 1), 多数国家为发展中国家, 其总人口超过44亿, 占全世界人口的63%。

1.1.2 数据来源

本文使用的基础数据包括：土地利用数据, 高程数据, 土壤数据, NPP数据, 气象数据, 流域及人口数据(表 1)。利用ArcGIS将所有数据统一为Winkel_Tripel投影坐标, 分辨率为1 km。

1.2 生态系统服务估算 1.2.1 土壤保持

土壤保持量(Soil Retention, SR)为潜在侵蚀量与实际侵蚀量之差^{[5, 15-16]}, 本文采用InVEST模型中的SDR模块, 估算潜在土壤侵蚀量和实际侵蚀量。

(1)

(2)

(3)

式中, SR为土壤保持量；RKLS为在研究区特定条件下裸地的土壤侵蚀量；USLE为考虑植被覆盖管理和水土保持措施下, 现状土地利用的土壤侵蚀量；R为降水侵蚀性因子, 使用了Fournier指数的基本公式^[19]；K为土壤可侵蚀性因子, 计算采用EPIC模型^[20]；LS为坡度坡长因子。C为植被覆盖和作物管理因子, 参考Borrelli等^[21]估算全球土壤侵蚀量的赋值方法；P为水土保持措施因子, 假设所有土地类型都受到水土保持的管理, 将P赋值为1。

(4)

(5)

(6)

式中, P为年降水量；P_i为i月的降水量；SAN为沙粒含量(%)；SIL为粉砂含量(%)；CLA为黏粒含量(%)；C₀为有机碳含量(%)。

1.2.2 产水量

本文采用InVEST模型中的产水模块计算一带一路区域的产水量(Water Yield, WY)^[22]。产水量评估模块基于Budyko水热耦合平衡假设(1974)和年平均降水量数据, 公式如下：

(7)

式中, Y(x)为研究区每个栅格单元x的年产水量(mm)；AET(x)为栅格单元x的年实际蒸散量(mm)；P(x)为栅格单元x的年降水量(mm)。

(8)

式中, PET(x)为潜在蒸散量；ω(x)为自然气候-土壤性质的非物理参数, 是一个经验参数。

(9)

(10)

式中, ET₀(x)为栅格单元x的参考作物蒸散；K_c(lx)为栅格单元x中特定土地利用类型的植物(植被)蒸散系数；AWC(x)为土壤有效含水量(mm)；Z为经验常数。

1.2.3 碳固定

NPP是指绿色植物在单位面积、单位时间内所积累的有机物数量, 是地表碳循环的主要组成部分, 也是判定生态系统碳源/汇和调节生态过程的主要因子^[23]。本研究碳固定(Carbon Storage)服务用植被净初级生产力(NPP)予以指代^[24]。由于NPP产品缺少2015年数据, 本研究使用2014年数据替代。

1.2.4 生境质量

生境质量(Habitat Quality, HQ)地图是由InVEST生境质量模型结合土地覆被和生物多样性威胁因素的信息而生成, 以耕地、道路、城镇和河流作为生境胁迫因子, 打分形成敏感性参数(表 2)^[25-26]。计算公式如下：

(11)

(12)

式中, Q_xj为LULC类型j中的栅格x的生境质量；D_xj为LULC或生境类型j中栅格x的总威胁水平；K和Z为比例因子；H_j为LULC类型j的生境适宜度；R为胁迫因子；Y_r为胁迫因子r所占的栅格数；W_r为威胁因子权重, 值在0—1之间；r_y为栅格y的胁迫因子值；i_rxy为栅格y的胁迫因子值r_y对栅格x的胁迫水平；β_x为栅格x的可达性水平；S_jr为生境类型j对胁迫因子r的敏感性(表 2)。

1.3 生态系统服务间关系分析

为了表征生态系统服务间的权衡/协同关系, 本研究利用SPSS软件, 对2000和2015年研究区土壤保持、产水量、碳固定和生境质量两两服务间的相关性系数进行计算^[27], 正相关表明生态系统服务间为相互增益的协同关系, 负相关则生态系统服务间表现为此消彼长的权衡关系^[5]。

为了研究生态系统服务间的空间权衡与协同关系, 本文以次一级小流域为基本单位, 对生态系统服务进行分区统计并将统计结果赋值到流域矢量图上, 基于GeoDA软件, 利用Weight模块构建空间权重矩阵, 在Space模块下的双变量局部莫兰指数(Bivariate Local Moran′s I)对生态系统服务进行双变量空间自相关分析^[15]。高高集聚和低低集聚表示空间协同, 高低集聚和低高集聚表示空间权衡。Moran′s I指数计算如下：

(13)

式中, z_i是要素i的属性与其平均值(x_i-X)的偏差, w_{i, j}是要素i和j之间的空间权重, n等于要素总数, S₀是所有空间权重的聚合。Moran′s I >0表示空间正相关性, 其值越大空间相关性越明显；Moran′s I < 0表示空间负相关性, 其值越小空间差异越大；Moran′s I = 0, 空间呈随机性。

1.4 生态系统服务关系地域分异规律分析

为了探究“一带一路”生态系统服务间权衡与协同关系的区域差异, 首先采集19123个随机样点(20000个点删除空值), 提取各样点上2000年、2005年、2010年和2015年4种生态系统服务值和人口数量, 然后将生态区和行政区关联到样点并输出为表格。在区域对比中, 本研究选取6个主要的生态区(TSF：热带和亚热带湿润阔叶林；TMF：温带落叶和混交林；TG：泰加林；TGS：温带草原和灌木；MGS：山地草原和灌木；DXS：沙漠和沙漠灌木), 按照地区选取了9个面积较大的国家(俄罗斯、中国、印度尼西亚、印度、哈萨克斯坦、沙特阿拉伯、埃及、乌克兰和希腊), 并将人口数量分为4个等级 < 10、10—100、100—1000和>1000人。在SPSS软件中计算不同生态区、不同行政区和不同人口规模下生态系统服务间的相关性系数, 得到不同地域中生态系统服务间关系的区域响应规律。

2 结果 2.1 生态系统服务的时空变化

在研究中, 将2000年数据作为生态系统服务评估起点, 2015年数据作为生态系统服务评估终点。图 2表示了2000和2015年的4种生态系统服务的空间格局和2000年到2015年4种生态系统服务的增益和减损情况。2000和2015年4种生态系统服务的空间格局整体变化不大, 局部发生了较为明显的变化。土壤保持量为0 t、碳固定值小于0.5 t/hm²和生境质量得分在0—0.2的区域分布大体保持一致, 分布于陆地内部和西南部, 主要位于中国的西北部、蒙古、中亚和西亚地区, 属于干旱地区, 植被为沙漠灌木, 沙漠广布。产水量的低值区(0—200 mm)除和其他3种生态系统服务分布一致的区域, 还包括西西伯利亚地区。土壤保持高值区与山区分布大体保持一致。产水量大于1000 mm的区域分布在南亚、东南亚和东亚, 集中于印度北部、尼泊尔、孟加拉国、不丹, 东南亚国家和中国的南部, 并在中纬度地带大致呈现又沿海向内陆产水量逐渐减少的趋势, 与降水量的分布规律吻合。在东南亚、青藏高原南部边缘和东欧地区碳固定大于10 t/hm², 拥有较高的固碳量主要是由于它们的植被为热带或温带森林, 植被覆盖率高^[26]。生境质量得分在0.90—1的区域分布在亚欧大陆的北部、东北部和东南部, 主要分布在俄罗斯和东南亚地区, 俄罗斯中西部生态威胁因子少, 自然生境质量得分高。

印度河流域上游和恒河流域中下游土壤保持量有明显的增加。2000—2015年印度河流域产水量增加明显。区域大部分地区碳固定有所增益, 主要增益区域分布于南欧和东欧的部分国家, 如希腊和保加利亚。2000—2015年间, 生境质量有改善的区域主要集中在亚欧大陆中部、中国的西北部、阿富汗和巴基斯坦。东亚、东南亚和东欧生境质量降低明显, 尤以中国的东部和中部最为显著, 且表现为自西向东生境质量减损趋于严重的空间态势, 人类活动和土地利用类型的变更可能是导致生境质量发生变化的主要原因。

2.2 生态系统服务的关联特征 2.2.1 生态系统服务关系的数值体现

由表 1可知, 2000、2005、2010和2015年土壤保持、碳固定、产水量和生境质量四种生态系统服务之间均存在着正相关性, 通过了0.01水平上的显著性检验, 呈相互增益的协同关系。其中碳固定与产水量的正相关性最强, 相关性系数大于0.55, 表现为强的协同关系, 体现了水是维系植物生长的基本条件。生境质量与碳固定的相关性在0.25到0.3之间, 呈较强的协同关系, 碳固定量多, 表明高植被覆盖率为生物提供了良好的生存条件。2000年到2015年, 生境质量与产水量的相关性系数呈下降趋势, 协同关系逐渐减弱。生境质量与土壤保持的相关性性最弱, 系数在0.07—0.1之间。

2.2.2 生态系统服务关系的空间表达

基于流域尺度, 对“一带一路”区域生态系统服务间关系进行了空间制图。统计了2000—2015年两两生态系统服务关系呈现为空间权衡或协同的流域数量(表 4)。由表 4可知, 2000—2015年, 大部分流域(>1400)生态系统服务间的空间关系不显著。生境质量与碳固定的空间协同流域与协同流域较多(>160), 且高低权衡流域比低高权衡流域多将近一倍, 但2000—2015年间流域数量没有明显的变化。生境质量与土壤保持的空间协同区有略微减少, 主要集中于高高集聚的流域, 空间权衡区有细微的增加。2000—2015年间, 生境质量与产水量的整体协同区由116个流域增加到150个流域, 西亚和东欧低低协同区增加明显。与此同时, 空间权衡区也有明显增加(2000年：128；2015年：175), 其中高低权衡区在2000—2015年间增加了30流域, 主要位于亚洲中部。碳固定与土壤保持呈现空间权衡和协同的流域数量总体基本不变。碳固定与水产量的协同区较多(2000年：147；2015年：196)且多数为低低权衡, 并且2000—2015年低低权衡增加了49个流域。水产量与土壤保持的空间协同区较少, 低高权衡区明显多雨高低权衡区。

2.3 生态系统服务间权衡/协同关系的地域分布规律 2.3.1 不同生态区下生态系统服务间关系的差异

总体来看, 2000年到2015年TSF(热带和亚热带湿润阔叶林)、TMS(温带落叶和混交林)、MGS(山地草原和灌木)、DXS(沙漠和沙漠灌木)4个生态区内生态系统服务间都为正相关, 表现为协同关系。但生态系统服务间的协同度存在地区差异, TSF(热带和亚热带湿润阔叶林)生态区内生境质量与土壤保持、产水量和碳固定之间呈较强的协同关系, 而土壤保持与产水量之间的协同关系较弱。TMF(温带落叶和混交林)生态区相较于其他生态区, 生境质量与土壤保持(0.14—0.16)、土壤保持与产水量(0.35—0.5)呈较强的协同关系, 但生境质量与产水量和碳固定之间的相关性系数在0.1左右, 协同度较低。MGS(山地草原和灌木)生态区土壤保持与碳固定的协同度比其他生态区要高, 而生境质量与产水量表现为较弱的协同关系。DXS(沙漠和沙漠灌木)生态区生境质量与土壤保持、碳固定与生境质量和产水量之间的协同度较高。TG(泰加林)和TGS(温带草原和灌木)生态区内, 生态系统服务间的关系较为复杂, 具体表现为2000、2005、2010和2015年间, 生境质量与土壤保持, 水产量与土壤保持, 均呈协同关系, 但生境质量与土壤保持的相关性(< 0.1)较水产量与土壤保持的相关性(0.15—0.3)低, TG(泰加林)生态区内生境质量与水产量、土壤保持间的关系四个时期几乎全部为负相关, 绝对值在0.1以内, 表现为较弱的权衡关系。TGS(温带草原和灌木)生态区生境质量与产水量相关性系数小于-0.1, 权衡关系较为显著。不同生态区气候等自然因素存在差异, 导致生态系统服务的空间差异, 从而使服务间的权衡与协同关系在不同生态区出现异质性。

2.3.2 不同国家生态系统服务权衡协同关系的差异

以行政区为单位, 按照地区分别选取9个国家, 统计了2000年、2005年、2010年和2015年生态系统服务的相关性(图 4)。结果发现, 俄罗斯4个生态系统服务间的相关性都较低, 生态系统服务间权衡与协同关系弱。中国相较于其他国家, 碳固定与土壤保持(>0.2)和产水量(>0.6)之间保持着较强的协同关系。印度尼西亚土壤保持与产水量和碳固定间的相关性系数的绝对值小于0.1, 且土壤保持与产水量之间总体保持权衡关系, 2000和2005年土壤保持与碳固定之间也为权衡关系。哈萨克斯坦土壤保持与产水量和碳固定间的协同关系较强(>0.35)。沙特阿拉伯生境质量与土壤保持间协同关系强(>0.4), 而产水量与碳固定间2000、2005和2015年均为负相关。埃及生境质量与碳固定的协同关系强, 但产水量与碳固定为权衡关系。埃及和希腊在4个时期, 产水量与碳固定均为权衡关系。

2.3.3 生态系统服务权衡协同关系对人口规模的响应

在不同人口规模下, 生态系统服务的权衡与协同关系呈现不同的特点。人口数量小于1000时, 生境质量与产水量对人口规模的响应与生境质量与土壤保持对人口规模的响应相似, 表现为人口极少(< 10)区域, 生境质量与土壤保持和产水量的相关性系数较低, 呈较弱的协同关系, 随着人口数量增加到10—100间, 协同关系随之增强, 人口数在100—1000间, 协同关系又有减弱。但人口数量达到1000后, 生境质量与土壤保持在2010和2015年, 协同关系达到最强, 而2000和2005年则表现为权衡关系。生境质量与产水量的相关性系数逐渐降低, 协同关系逐渐减弱, 2000、2005和2010年甚至转为权衡关系。土壤保持与水产量的关系总体表现为随着人口的增加, 协同关系逐渐减弱, 2000和2005年略有不同, 在人口超过1000时, 协同关系最强。随着人口的增加, 水产量与碳固定间的协同关系显著降低。

3 讨论 3.1 生态系统服务的区域响应规律

生态系统服务作为链接自然和人类社会的桥梁和纽带, 自然生态系统和社会生态系统均是影响生态系统服务的重要因素^[28]。生态系统服务的时空分布主要受到自然因素的影响, 如产水量的空间分布格局与降水量分布规律吻合, 亚欧大陆内陆、西亚和北非气候干旱, 降水少, 因此成为产水量的低值区。而水是生命之源, 产水量低的区域, 往往植被覆盖率低, 碳固定低, 进而导致生物生存条件差, 生境质量低。所以, 由于受到降水的影响, 产水量、碳固定和生境质量的低值区分布大体一致。此外, 土地利用变化是导致生态系统系统服务变化的主要原因^[28-29], 2000—2015年间, 中国等东南亚国家经济发展迅速, 对土地利用强度大, 因此这些区域生境质量显著降低。

生态系统服务供给和需求往往在空间上不匹配^[10], 生态系统服务间相互关系出现区域异质性。如热带和亚热带湿润阔叶林生境质量与碳固定为协同关系, 而温带草原和灌木生境质量与碳固定间则为权衡关系, 这可能是由于温带草原区畜牧业发达, 生境质量高时, 对草地的利用率高, 进而导致碳固定降低。必须明确生态系统服务供给与需求的动态变化和驱动机制, 才能对比不同区域间生态系统服务权衡协同的关键环节和主导因素, 明晰不同区域所承担的社会、经济和生态的主体功能^[4]。

权衡研究的最终目的是针对特定区域复杂的生态系统服务进行权衡关系科学管理, 通过对比不同区域生态系统服务间的相互关系, 呈强协同关系的区域, 将对协同关系弱或者呈权衡关系的区域具有空间指引作用。如同为草原生态区, 山地草原和灌木土壤保持与碳固定间为协同关系, 而温带草原与灌木土壤保持与碳固定间则表现为权衡关系, 山地草原处于生态脆弱区, 畜牧等农业活动受到限制, 而温带草原区畜牧业发达, 这可为温带草原区的土地管理决策提供一定指引, 促进生态系统服务间关系朝协同方向发展。当人口规模在10—100时, 各个生态系统服务间的关系都呈较强的协同关系, 因此研究生态系统服务间关系对人口规模的响应对区域制定人口政策有一定指导意义。

3.2 研究优势与不足

目前, 对生态系统服务间关系的研究主要集中在各类生态系统服务间的相互关系、人地耦合和表现类型上, 并且大多集中在小尺度研究生态系统服务间关系的空间格局。本文基于大尺度区域, 利用相关分析和空间自相关分析, 定量表达了生态系统服务间的关系及其关系的空间表达, 并从不同生态区, 不同国家和不同人口规模下, 分析了生态系统服务间关系的地域分异规律。然而, 本文仅分析了生态系统服务权衡协同的空间关系, 和区域间生态系统服务的权衡协同关系存在异质性的这些现象, 未对这种权衡协同的空间关系和区域异质性的驱动因素进行深入探讨。由于数据获取困难, 只评估了生境质量、土壤保持、碳固定和产水量这4种生态系统服务。对生境类型进行敏感性参数赋值时, 不同的植被生态区对应不同的林地类型, 实际赋值时对所有林地统一赋值为1, 运行模型时未考虑不同植被生态区生境得分的差异, 如热带和亚热带湿润阔叶林生境得分应比泰加林得分高, 但为了不影响不同生态区的生态系统服务权衡关系的区域响应结论, 因此在生境质量计算中不直接使用植被生态区作为赋值规则。两种生态系统服务间相互关系可以界定为6—10种类型^[30], 两种生态系统服务间关系可能受共同因素驱动, 或者生态系统服务间的相互作用^[31], 然而这种作用与驱动因素间并非是一成不变的。在未来研究中, 析取驱动因素中影响生态系统服务关系变化的关键阈值, 准确预测“一带一路”区域生态系统服务及其关系的未来情景十分必要。

4 结论

为明晰“一带一路”区域生态系统服务间相互关系的区域差异性, 本文分析了生境质量、土壤保持、产水量和碳固定4种生态系统服务的时空变化及它们间的权衡协同关系的地域分异特点。结果表明：

2000、2005、2010和2015年, “一带一路”区域土壤保持、产水量、碳固定和生境质量这4种生态系统服务的空间分布有所差异, 低值区的分布大致重合, 高值区分布各具特点。2000—2015年间, 生态系统服务的增益与减损区的分布格局不同。4种生态系统服务间的关系在数值上都变现为相互增益的协同关系, 其中碳固定与产水量的协同关系强, 相关性系数大于0.5, 生境质量与土壤保持的协同关系较弱, 相关系系数小于0.1。从生态系统服务权衡协同关系的空间表达来看, 区域异质性显著, 2000到2015间其空间格局也有细微的变化。

生态系统服务权衡协同关系具有地域分异特征, 且在不同尺度上表现不同, 在生态区尺度下, 泰加林和温带草原灌木与其他生态区相比差异显著, 如生境质量与碳固定, 生境质量与水产量之间为权衡关系。在国家尺度下, 生态系统服务间的关系在不同国家表现不同, 如俄罗斯4种生态系统系统服务间的协同关系相较于其他国家整体较弱, 俄罗斯、沙特、埃及和乌克兰产水量与碳固定为权衡关系。生态系统服务间权衡与协同关系对人口有的响应程度不同, 随着人口的增加, 碳固定与产水量和土壤保持和水产量的协同关系逐渐减弱, 响应程度较强, 而土壤保持与碳固定间的关系对人口无明显的响应。