生态系统服务(ecosystem service, ES)是指通过生态系统的结构、过程和功能直接或间接得到的生命支持产品和服务^[1]，而居民福祉(Human well-being, HWB)则是表征人类生活状态及水平的多维概念^[2]。生态系统服务桥接自然系统与居民福祉，其通过生态过程和功能为人类提供各种产品和服务，而人类则通过对生态系统服务的消费来满足和提高自身福祉^[3]。进入21世纪以来，生态系统的惠益为人类生活水平与幸福感的提升提供了物质基础^[4]，然而随着社会经济快速发展，生态系统退化与人类需求提升之间的矛盾不断加剧，生态退化、环境污染、资源紧缺等问题日益凸显^[5]。探讨生态系统服务与居民福祉关联的科学认知已成为生态学、地理学等学科重要的研究领域。随着“人—地”系统耦合研究的不断深化，学者们的研究从生态系统服务定量测度、生态系统服务权衡协同开始扩展到生态系统服务与居民福祉的关联、耦合^[6]。从研究尺度上看，学者们的研究揭示了从县域、市域尺度到省域尺度上生态系统服务与居民福祉的关联^[7-9]；从生态系统视角看，草地生态系统、河流生态系统、城市生态系统与居民福祉关联关系差异显著^[9-11]；从研究方法上看，生态系统服务价值当量表和InVEST模型等被用来测度生态系统服务价值，而人类发展指数^[12]、生活满意度^[13]、国家福祉指数^[14]及全球幸福指数等^[15]则被用来测度居民福祉水平。总体来看，学者们的研究揭示了不同区域生态系统服务与居民福祉的关系，为区域生态系统管理服务于社会经济的可持续发展提供了参考。然而生态系统服务的水平取决于生态系统功能，而生态系统服务的动态变化则受制于非线性的生态系统结构和功能的变化。因此，不同类型生态系统及其与居民福祉之间的动态关联特征具有明显的地域特征^[16-17]。中国北方农牧交错带是保护中国中东部地区免受风沙侵扰的重要生态安全屏障，同时也是区域重要的涵养水带。农牧交错带具有经济形态独特、生态环境脆弱、土地利用受人类活动影响显著等特征^[18]。由于气候变化、人类不合理开发以及大型生态修复工程(防沙治沙、退耕还林和禁牧封育等)的实施，农牧交错带生态系统经历了显著变化，该区“人—地”互作模式也随社会经济发展显著改变。然而，目前对中长时间序列的中国北方农牧交错带的生态系统服务与居民福祉的研究尚未见诸报道。鉴于此，本文以中国北方农牧交错带的桑干河流域(山西段)为研究对象，运用单位面积价值当量因子法、熵权法、聚类分析等方法，探究桑干河流域(山西段)2010—2020年内生态系统服务与居民福祉的空间关联特征。研究结果可为推动桑干河流域(山西段)生态系统服务从科学认知走向管理决策，加强生态系统服务管理和制定合理的生态补偿方案等提供决策辅助。

1 研究区概况

桑干河是京津冀重要的水源涵养区和防风固沙生态屏障带。桑干河发源于山西省北部，属海河水系，为永定河的上游。桑干河自忻州市宁武县管涔山庙儿沟自西南向东北流向大同市阳高县南徐，并在该地出山西省境进入河北省，于官厅水库上游与洋河汇合后注入永定河。桑干河主要支流有源子河、黄水河、浑河、御河等。桑干河流域(山西段)介于东经112°35′—114°25′，北纬39°23′—40°15′之间，面积为16 767 km²，干流长261 km，年降水量为350~520 mm，天然河川径流量多年平均值为5.39×10¹² m³。区内主要山脉有恒山、洪涛山、管涔山等。本研究按照流域界线，将宁武县、朔城区、山阴县、应县、怀仁县、平鲁区、南郊区、新荣区、大同县、浑源县、左云县、广灵县、天镇县、阳高县作为研究区。桑干河(山西段)地处干旱—半干旱过渡带，区内生态环境敏感且脆弱。随着社会经济发展及气候变化，区内生态系统胁迫特征凸显，主要表现为：河流水生态退化、沙化土地面积和土壤侵蚀面积大、历史遗留矿山点多面广、盐碱地范围广且治理难度大。近年来，随着河流生态环境综合治理的实施，区内生态系统服务功能明显向好。然而中长时间尺度上桑干河流域(山西段)生态系统服务与居民福祉关联的空间分异规律尚不明晰。

2 数据来源与研究方法 2.1 数据来源

本文的数据主要包括桑干河流域(山西段)土地利用变化数据和社会经济统计数据。桑干河土地利用变化数据(2010，2015和2020年)来自中国科学院资源环境科学与数据中心。社会经济统计数据包括人口、主要农作物播种面积、粮食产量、居民消费价格指数、城镇和农村居民人均可支配收入和蔬菜、水果、水产品占有量等，该类数据来自《2011年山西省统计年鉴》《2016年山西省统计年鉴》《2021年山西省统计年鉴》《全国农产品成本收益资料汇编》及各县(区)政府统计公报。

2.2 研究方法 2.2.1 桑干河流域(山西段)生态系统服务价值评估方法

目前生态系统服务价值评估主要有两种方法，一种是利用InVEST模型对生态系统服务价值进行计算^[19]，另一种是谢高地在Costanza^[20]的基础上确定的生态系统功能当量因子法。InVEST模型在单项生态系统服务价值的测算及可视化表达上有优势，而生态系统功能当量因子法则在大尺度生态系统服务价值评估上占优。本文利用谢高地等最新的生态系统功能当量表计算桑干河流域(山西段)生态系统服务价值^[21]。根据研究区的土地利用状况，本文对谢高地等提出的生态系统服务价值当量表进行重新修正(表 1)^[22-23]。具体修正方法为：由于研究区耕地以旱地为主，因此耕地的生态系统服务价值当量取旱地的生态系统服务价值当量值；林地的生态系统服务价值当量取针叶、针阔混交、阔叶、灌木的平均值；草地取草原、灌草丛、草甸的平均值；水域取湿地和水系的平均值；未利用地取荒漠和裸地的平均值；建设用地生态系统服务的价值当量设为0。

桑干河流域(山西段)生态系统服务价值测算方法为：

(1) 根据研究区各县(区)的农林业生产特征，选取玉米、大豆和谷子作为主要粮食作物，并根据研究区各县(区)2010，2015和2020年主要粮食作物的平均粮食单产、粮食播种面积、粮食总播种面积和粮食单价，计算研究区2010，2015和2020年的粮食单产市场价值(Q)。计算公式为：

$ Q=\frac{m_i p_i q_i}{M} $

(1)

式中：m_i为第i种粮食作物的单价；p_i为第i种粮食作物的平均粮食单产；q_i为第i种粮食作物的播种面积；M为粮食作物的总播种面积。

(2) 令单位面积(hm²)的生态系统服务价值当量因子的经济价值等于当年平均粮食单产市场价值的1/7，根据谢高地等提出的中国陆地生态系统单位面积生态服务价值当量和研究区粮食单产及其经济价值，可计算出研究区生态服务单位当量价值(E_a)、单位面积价值(V_ci)以及生态系统服务价值(ESV)。计算公式为：

$ E_a=Q / 7 $

(2)

$ V_{c i}=E_a a_{c i} $

(3)

$ \mathrm{ESV}=A_i V_{c i} $

(4)

式中：Q为研究区平均粮食单产市场价值; a_ci为第c类生态系统第i中生态服务类型的生态服务价值当量; V_ci为第c类生态系统第i中生态服务类型的单位面积生态服务价值; A_i为第i类土地利用类型的面积。

(3) 由于经济发展和通货膨胀等因素，根据研究区2010，2015和2020年的平均消费价格指数计算出当年的购买力指数(以2010年作为基期年进行折算)。利用购买力指数对生态系统服务价值进行核算^[24]，得到桑干河流域(山西段)各县(区)2010，2015和2020年核算后的生态系统服务价值单价(表 2)。

2.2.2 桑干河流域(山西段)居民福祉的评价方法

基于《千年生态系统评估报告》中的福祉内涵及相关研究^[25-26]，考虑数据的可得性，本文从收入、物质、安全3个维度出发构建了桑干河流域(山西段)居民福祉评价指标体系(利用购买力指数对居民人均可支配收入和城镇居民人均可支配收入进行核算)。在评估中，首先，通过利用极差标准化法对各指标进行标准化处理，并运用熵权法确定各指标权重(表 3)^[27]；然后，通过加权求和，测算桑干河流域(山西段)各地区的居民福祉指数。

2.2.3 桑干河流域(山西段)生态系统服务与居民福祉空间关系分析方法

本文用K-means聚类分析法来划分生态系统服务与居民福祉的空间关联模式。聚类分析法是基于聚类样本自身属性，按照属性间的相似性或差异性特征，定量识别对象之间的亲疏关系，并以此定量的亲疏关系程度对样本进行聚类。算法流程为：

(1) 随机选取K个样本点作为初始聚类中心；

(2) 计算样本特征向量x_i和每个聚类中心c_j的距离，计算公式为：

$ \begin{aligned} & d\left(x_i, c_j\right)=\sum\limits_{k=1}^m\left(x_{k i}-c_{k j}\right)^2 \\ & x_i=\left[x_{1 i}, x_{2 i} \cdots x_{m i}\right]^T, c_j=\left[c_{1 j}, c_{2 j} \cdots c_{m j}\right]^T \end{aligned} $

(5)

i∈{1, 2…n}, j∈{1, 2…K}, n为数据集的样本数，K为聚类的类别数；

(3) 将样本x_i划分到离他最近的聚类中心所对应的类中；

(4) 使用样本特征向量的均值，更新各个聚类中心，均值的计算公式为：

$ \bar{x}_G=\frac{1}{|G|} \sum\limits_{i=1}^{|G|} x_i, x_i \in G $

(6)

式中：x_i为样本特征向量；|G|为类别G的样本个数；

(5) 重复步骤(2)—(4)直到聚类中心不再变化。

3 结果 3.1 桑干河流域(山西段)生态系统服务价值与居民福祉的空间动态特征

(1) 桑干河流域(山西段)生态系统服务价值空间动态特征。本文利用ArcGIS中的自然间断点分级法将桑干河流域(山西段)14个县(区)2010，2015和2020年人均生态系统服务价值划分为低、中、高3个等级(图 1)。从空间分布看，研究区人均生态系统服务价值的低值区分布在流域的中部、西南部，而其他区域尤其是平鲁区、山阴县、天镇县、阳高县、广灵县等县(区)人均生态系统服务价值较高。流域人均生态系统服务高值区林草地比重大；桑干河干流及其支流(恢河等)流经其中，形成了大片的洪积冲积扇，该区土壤肥沃，湖沼湿地及地下水资源丰富；桑干河自然保护区等高服务生态系统分布其中。2010—2020年，桑干河流域(山西段)人均生态系统服务价值发生变化的有9个县(区)，占研究区全部县(区)的64%(图 2)。2010—2020年研究区人均生态系统服务价值总体上呈上升趋势，这与地方政府对生态环境保护和修复的高度重视有关。生态修复前，桑干河流域(山西段)沙化土地面积大(沙化面积占流域总面积的22.2%)、河流水生态退化、水资源总量匮乏、地表径流量持续减少、流域土壤侵蚀面积占流域总面积的55%。历史遗留矿山点多面广，对生态系统破坏严重、盐碱地范围广，治理难度大。近年来，大同、朔州、忻州三市联动，共管共治，全力保障桑干河水生态安全。三市加强专项规划与河流保护相关规划的衔接融合，全面启动实施《桑干河流域生态修复与保护规划》，统筹流域山水林田湖草沙系统治理。区内水环境、水土流失、土壤沙化、土地盐碱化等生态胁迫现象得到显著改善。同时该区大力整顿矿山开采、医药化工等污染环境、破坏生态、过度消耗资源的项目。上述举措显著提高了桑干河流域(山西段)生态系统服务价值。

(2) 桑干河流域(山西段)居民福祉空间动态特征。利用自然间断点分级法将桑干河流域(山西段)14个县(区)2010，2015和2020年居民福祉划分成低、中、高3个等级(图 3)。2010—2020年桑干河流域(山西段)居民福祉的高值区分布在沿主干道的西南—东北方向上。山阴县、应县、阳高县、大同县较高，朔城区、怀仁县等次之，其余区域较低。从居民福祉变化看，2010—2020年桑干河流域(山西段)发生变化的区域有9个县(区)，占全部研究区域的64%(图 4)。10 a间研究区居民福祉整体呈现上升趋势。朔州市具有典型农牧交错带特征的社会经济结构，其是中国重要的奶牛和肉羊养殖基地，同时也是全国重要的动力煤基地、华北重要的电力基地及全国较大的日用瓷生产基地，区内朔城区、山阴县、怀仁县等县(区)经济发展水平较高，同时，该区大力推广的现代设施蔬菜、羔羊养殖加工园区建设和设施农业建设，逐渐形成“黑色煤炭、白色乳品、绿色生态”三色循环经济发展模式，充分挖掘生态系统服务价值，惠益了区内居民生活水平与幸福感提升。桑干河流域(山西段)煤炭等矿产资源丰富，区内多为资源主导型经济结构，产业结构发展不均衡，近年来区内社会经济转型发展带来的“阵痛”降低了区内居民福祉。然而，随着区域现代农业的快速发展，以及大同火山、大同黄花等为主线的特色文化旅游带动了桑干河流域(山西段)经济高质量发展，特别是现代农业、商贸、康养、体育、生态经济等相关产业的相互渗透，促进了区域居民福祉提升。

3.2 桑干河流域(山西段)生态系统服务与居民福祉的空间关联模式及其变化

(1) 桑干河流域(山西段)生态系统服务与居民福祉空间关联模式的划分。2010，2015和2020年桑干河流域(山西段)生态系统服务与居民福祉空间关联聚类结果所对应参数详见表 4。根据聚类分析结果，本文将研究时段内桑干河流域(山西段)生态系统服务与居民福祉存在的空间关联模式划分为8种：高服务高福祉(HH)、高服务低福祉(HL)、高服务中福祉(HM)、低服务高福祉(LH)、低服务中福祉(LM)、低服务低福祉(LL)、中服务高福祉(MH)、中服务中福祉(MM)(图 5)。2010年HH主要分布于山阴县和应县，HL主要分布于宁武县，HM主要分布于天镇县和广灵县，LH主要分布于怀仁县和南郊区，MH主要分布于朔城区、阳高县和大同县，MM主要分布于平鲁区、左云县、新荣区和浑源县；2015年HH主要分布于平鲁区、山阴县、应县、天镇县、阳高县和广灵县，LL主要分布于朔城区、南郊区、怀仁县和大同县，MM主要分布于宁武县、左云县、浑源县和新荣区；2020年HH主要分布于天镇县、阳高县和广灵县，HL主要分布于宁武县，HM主要分布于新荣区和浑源县，LH主要分布于怀仁县，LM主要分布于朔城区和南郊区，MH主要分布于山阴县、应县和大同县，MM主要分布于平鲁区和左云县。

(2) 桑干河流域(山西段)生态系统服务与居民福祉空间关联模式的变化。2010—2015年研究区内各县(区)空间关联模式变化较大(图 6)，主要变化出现在朔州南部和大同市，其中朔城区和大同县的变化方向为MH→LL，南郊区和怀仁县的变化方向为LH→LL，天镇县和广灵县的变化方向为HM→HH，阳高县的变化方向为MH→HH，平鲁区的变化方向为MM→HH，宁武县的变化方向为HL→MM。2015—2020年，研究区空间关联模式变化主要发生在朔州市、大同市大部分地区，其中山阴县和应县的变化方向为HH→MH，新荣区和浑源县的变化方向为MM→HM，朔城区和南郊区的变化方向为LL→LM，平鲁区的变化方向为HH→MM，怀仁县的变化方向为LL→LH，宁武县的变化方向为MM→HL。桑干河流域(山西段)10 a来生态系统服务价值的变化总体上呈增高的趋势，居民福祉的变化大部分地区呈增长趋势，少数地区则呈降低趋势。研究区生态系统服务价值升高与流域水环境治理、水土流失治理、土壤风沙盐碱治理等系统生态修复工程的实施以及医药化工等污染环境、破坏生态、过度消耗资源的项目的整顿有关。

研究区居民福祉的变化则与区域经济结构调整关系密切，区内产业结构发展不均衡，煤炭等资源产业在区域经济结构中占主导地位，资源型产业的发展对区内耕地的自然资源造成胁迫，导致居民福祉降低，然而随着近年来区内资源型经济转型发展，特别是特殊农业、康养、生态经济等相关产业的发展，区内居民福祉提升。

(3) 桑干河流域(山西段)生态系统服务与居民福祉不同关联模式的发展建议。根据桑干河流域(山西段)生态系统服务与居民福祉的空间关联模式，结合研究区生态环境和社会经济的特征，对不同的空间关联模式提出生态系统管理建议。高服务高福祉、高服务中福祉、中服务高福祉及中服务中福祉区主要有山阴县、应县、天镇县、广灵县、浑源县、左云县等，区内分布有桑干河国家湿地公园、大同火山群国家地质公园、应县木塔等风景名胜区，该区人均生态系统服务价值较高。区内沿桑干河及其支流分布的盆地区则为重要的粮食产区，近年来桑干河(山西段)水生态治理、风沙治理、水土流失治理等举措直接惠益居民福祉的改善和提高。该区宜不断加强旅游与现代农业、商贸、康养、生态经济等相关产业的相互渗透，积极发展“绿色生态、黑色煤炭、白色乳品”产业模式。而宁武县是高服务低福祉，宁武县作为桑干河的源头，生态系统服务价值高，因其交通不便，经济欠发达，且离经济发达地区较远，导致其居民福祉处于低水平状态，该地区拥有丰富的旅游资源，可加强对旅游资源的宣传，加强基础设施建设，发展生态旅游带动经济发展。桑干河流域(山西段)低服务高福祉、低服务低福祉、低服务中福祉主要分布于南郊区、怀仁县、朔城区、大同县等县(区)。该区以资源型经济结构为主，受经济结构及自然条件胁迫，该区风沙侵蚀、水土流失、水生态及水环境问题较突出。该区宜在积极推进资源型经济转型发展的同时加强国土空间生态修复，推动污染型企业整改或转型，推进现代农业发展力度，加大如怀仁羔羊肉等地理标志产品宣传推广，产业结构不断优化、生态环境质量不断提升。

4 讨论与结论 4.1 讨论

本文将货币购买力变化纳入生态系统服务价值的测算中，并尝试用购买力指数对人均生态系统服务价值基础单价及居民福祉进行修正，研究结果可在一定程度消除因不同年份比较基准的不同而产生的偏差。本文运用单位面积价值当量因子法测算生态系统服务价值，与InVEST模型相比，该方法在大尺度生态系统服务价值评估上占优，然而在测算精度上尚存在一定劣势。本文结合研究区特点及数据的可得性，构建了居民福祉评价指标体系。然而鉴于居民福祉的表达兼具主观、客观双重特性，本文评价尚需要进一步融合居民福祉主观表达指标。因此，未来可结合专家咨询、入户调查、部门访谈等方式对居民福祉进行主客观综合测度。桑干河流域(山西段)是京津冀地区重要的水源涵养区和防风固沙生态屏障带，探讨桑干河流域(山西段)生态系统服务与居民福祉的空间关系模式，能够找到改善生态环境，提高居民福祉的关键点，为京津冀地区类似地区生态和经济协调发展提供决策依据。

4.2 结论

(1) 2010—2020年桑干河流域(山西段)人均生态系统服务价值呈现不断增长的趋势，居民福祉整体也呈增长趋势。人均生态系统服务价值呈现中部和西南部低的空间分布格局，而居民福祉水平的高值区分布在西南—东北方向上，特别是山阴县、应县、阳高县等地区的居民福祉10 a来保持较高水平。

(2) 研究区生态系统服务和居民福祉的空间关联模式有：高服务高福祉、高服务低福祉、高服务中福祉、低服务高福祉、低服务中福祉、低服务低福祉、中服务高福祉、中服务中福祉。2010—2020年，研究区内各县(区)空间关联模式变化主要出现在朔州和大同市部分地区。

(3) 2010—2020年间桑干河流域生态系统服务价值与居民福祉均呈增长趋势，两者的空间关联模式呈现明显的空间分异规律。桑干河流域(山西段)宜通过积极推进资源型经济转型发展及加强国土空间生态保护与修复等途径提高区域生态系统服务与居民福祉。