李魁明，王晓燕，姚罗兰，石云

1.防灾科技学院应急管理学院

2.首都师范大学资源环境与旅游学院

3.北京农业职业学院国际教育学院

4.防灾科技学院生态环境学院

在社会经济发展加速转型时期，人类活动导致土地利用格局的改变，产生直接或间接的生态环境效应。1997年Costanza提出了生态系统服务价值（ESV）估算原理、方法后，生态系统服务概念和价值研究得到广泛重视[1]。近年来，区域土地利用变化对ESV的影响研究，已成为生态学、地理学和环境科学等学科的研究热点：从研究空间尺度上看，有的以行政区[2-4]、城市群[5-6]、经济区[7-8]等社会经济功能区为评价对象；有的以生态功能区[9-10]为评价对象，有的以流域[11]、高原[12]、山地[13]、盆地[14]等自然地理单元为评价对象。从研究方法上看，ESV的评估一般将物质量评价法和价值量评价法结合进行研究[15-16]，其中以基于土地利用类型的生态系统服务价值系数法[3-8]应用最广泛。关于ESV影响因素方面，或是基于土地利用方式的转换分析ESV的变化特征[17]，或是建立指标与ESV的回归结果解释影响因素[18]，2类研究方法均忽视了因子作用强度空间差异。地理探测器能够定量地探测单因子驱动力和多因子交互驱动力，定量识别因子强度空间差异。在城市群尺度上应用地理探测器方法开展生态服务价值的空间分异研究并不多见。

京津冀城市群是我国经济发展格局中最具活力和潜力的前沿阵地[19]，但近年来快速城镇化和工业化进程导致土地利用结构的改变也给该区域带来了严重的生态环境问题。加强对京津冀地区土地利用变化、ESV时空演变及驱动因子的研究，对于京津冀协同发展具有重要意义。汪东川等[20]利用县域尺度的粮食单产对单位面积农田食物生产服务价值进行修正，探讨了ESV的空间异质性，选取3 km×3 km网格尺度揭示京津冀地区2005—2015年ESV的空间动态演变特征。娄佩卿等[21]基于谷歌地球引擎(GEE)云平台，利用当量估算方法定量估算京津冀地区的ESV并结合15 km×15 km网格尺度探明其时空动态变化。唐秀美等[22]在分析生态系统服务供需逻辑关系的基础上选择指标，评估京津冀地区生态系统服务需求状况，从县域尺度分析ESV的时空变化。盛晓雯等[23]采用土地利用动态度和ESV模型等方法，对京津冀地区土地利用变化影响下的ESV在时间和空间上的差异进行了分析。以上研究多从整体上评估京津冀地区ESV，定量分析城市群内部地市的ESV时空演变的较少，并且鲜见ESV影响因子强度空间差异化的系统分析。

鉴于此，笔者针对京津冀地区土地利用变化对生态系统造成的影响，根据谢高地等[24]提出的ESV评估体系，分析2000—2018年京津冀地区土地利用变化引起的ESV时空规律，运用地理探测器分析京津冀地区影响因子强度空间差异，以期为研究区土地利用规划、环境修复与保护提供参考。

1 研究区概况与研究方法

1.1 研究区概况

京津冀地区位于华北平原北部，地势西北高东南低，地貌以平原为主，包括北京、天津两大直辖市以及河北省的11个地级市，共13个地级市以上行政单元，总面积为21.72万km2，占全国总面积的2.26%。在社会经济迅速发展的同时，太行山燕山土壤侵蚀和坝上高原荒漠化加剧，平原洼淀萎缩及消失等也对京津冀原本脆弱的生态环境造成了巨大压力。

1.2 数据来源与处理

京津冀地区土地利用数据源自中国科学院资源环境科学数据中心（https://www.resdc.cn）的2000年、2010年和2018年3期Landsat TM影像；根据中国科学院资源环境信息数据库分类法，对3期土地利用数据进行重分类得到耕地、林地、草地、建设用地、水域和未利用地六大类土地利用类型。粮食单产数据来源于3省（市）相应年份统计年鉴，粮食价格根据《中国农产品价格调查年鉴（2019年）》得出，驱动因子指标中，数字高程模型（DEM）数据来源于美国对地观测全球高程数据；降水量、土壤侵蚀强度来源于中国科学院资源环境科学数据中心；GDP、农民人均纯收入、社会消费品零售额、人均耕地面积、人口密度等来源于相关省（市）统计年鉴。

1.3 研究方法

1.3.1 土地利用动态分析法

为定量分析京津冀地区土地利用变化，引入单一土地利用类型动态度和综合土地利用动态度，计算公式如下：

1.3.2 生态服务价值测算

结合谢高地等[24-25]的研究成果，同时借鉴已有研究[26-28]，通过查阅统计年鉴，计算得出京津冀地区2000年、2010年和2018年平均粮食单产量（5 253 kg/hm2），通过查阅《中国农产品价格调查年鉴（2019年）》计算出京津冀 2018年粮食平均单价（2.65元/kg），确定京津冀地区单位面积生态系统服务价值系数（表1）。由此计算京津冀地区ESV，其计算公式见式（3）。

表1 京津冀地区单位面积生态系统服务价值系数Table 1 Ecosystem service value coefficient per unit area in Beijing-Tianjin-Hebei Region 元/hm2

考虑到京津冀地区13个地级以上城市大小不一，引入地均生态系统服务价值〔式（4）〕和地均生态系统服务价值变化率〔式（5）〕，来分析京津冀地区生态系统服务价值的时空演变。使用敏感性指数[27]分析生态系统服务价值对于生态系统价值系数的依赖程度〔式（6）〕。

1.3.3 地理探测器模型

地理探测器是用来探测地理要素空间分异性，揭示地理要素背后驱动因子的一种方法[29-30]。采用地理探测器模型对京津冀地区影响ESV空间分异的自然因子、社会因子和经济因子进行探测，公式如下：

2 结果与分析

2.1 土地利用动态变化

从土地利用结构上看（图1），2000—2018年，京津冀地区土地利用结构以耕地、林地、草地为主，2000年，这3种土地类型的占比分别为50.81%、20.72%和16.32%，合计占比为87.84%；2018年的占比分别为46.14%、21.13%和15.75%，合计占比为83.02%。2000—2018年，各类型土地面积变化最大的为建设用地，其次是耕地，建设用地面积增加了9 992 km2，耕地面积减少了9 742 km2。而林地、草地、水域、未利用土地4种土地利用类型面积变化幅度不大。

图1 2000—2018年京津冀地区不同土地利用类型面积变化Fig.1 Area change of different land use types in Beijing-Tianjin-Hebei Region from 2000 to 2018

2000—2018年京津冀地区土地利用变化动态度如表2所示。从单一土地利用动态度来看，2000—2010年，京津冀地区建设用地变化最剧烈，其次是水域和未利用土地，其动态度分别为1.19%、-0.38%和-0.37%；2010—2018年，建设用地变化幅度依然最大，其次是水域和未利用土地，其中建设用地变化幅度明显上升，水域呈增长趋势，耕地较上一阶段减少幅度提升。总的来看，京津冀地区近18年来，土地利用幅度变化剧烈的是建设用地、未利用土地，变化幅度大于1%；耕地因面积占比最大，变化不明显；草地和林地变化较小；水域面积波动大，不稳定。

表2 2000—2018年京津冀地区土地利用变化动态度Table 2 Land use change dynamics degree in Beijing-Tianjin-Hebei Region from 2000 to 2018%

2000—2018年，综合土地利用动态度为0.30%，其中 2000—2010年仅为 0.10%，但 2000—2018年为0.57%，这表明京津冀地区在2010年后土地利用变化速率加快，也表明区域社会经济活动对土地影响更加剧烈。土地利用动态度与人类活动强度密切相关，例如耕地被侵占，建设用地面积扩张，未利用地逐渐被开发、实施退耕还林还草等生态修复工程，导致京津冀各地市各时期土地利用动态度出现不同幅度的变化。

2.2 生态系统服务价值变化

2.2.1 京津冀地区生态系统服务价值时间序列分析

根据式（3）得到京津冀地区2000年、2010年、2018年的 ESV分别为 4 346.78亿、4 295.36亿、4 413.11亿元，呈先下降后上升趋势，总共增加了66.33亿元，其中2000—2010年减少了51.42亿元，2010—2018年增加了117.75亿元（表3）。

表3 京津冀地区2000—2018年生态系统功能单项服务价值Table 3 Ecosystem individual service value of Beijing-Tianjin-Hebei Region from 2000 to 2018

从生态系统功能构成一级内容来看，调节服务与支持服务始终是京津冀地区生态系统服务的核心功能；从占比分析，各类服务占比变化幅度较小，说明2000—2018年京津冀地区各类生态系统服务较为稳定，且其中调节服务占比一直保持在71%左右，其次为支持服务，占比一直保持在18%左右，文化服务与供给服务占比均较小。从生态系统功能构成二级内容来看，水文调节的贡献率最大，2000—2018年均在35%左右，其次是气候调节、气体调节、土壤保持，贡献率均在10%左右。

2.2.2 地均生态系统服务价值的时空变化

根据式（4）得到京津冀地区2000年、2010年、2018年的AESV分别为201.67万、199.27万、204.10万元/km2。2000年相较2010年降低了1.19%，2010相较2018年增加了2.42%，呈现先下降后上升的趋势，上升幅度大于下降幅度，说明京津冀地区生态系统服务趋于好转。

将AESV大于200万元/hm2的城市划分为高值区，100万～200万元/hm2为中值区，小于100万元/hm2为低值区，作出京津冀地区13个地级以上城市的AESV空间演变趋势，结果如图2所示。由图2可知，京津冀地区AESV空间分布差异明显，高值区主要分布在冀北和京津地区。具体来看，2000年和2010年的高值区主要集中在北京市、天津市、承德市和秦皇岛市；低值区主要分布在中南部的平原地区，包括衡水市、沧州市、廊坊市。2018年相对于2000年和2010年，高值区增加了唐山市，开始连片分布，低值区少了沧州市。承德市和北京市土地利用类型以林地为主，近年来造林工程的实施和山区草地转换为林地，对ESV贡献较大；天津市东邻渤海地区，地势较低，水域面积广大，在土地利用类型中，水域生态系统服务价值系数最大；秦皇岛市地处燕山山脉东段丘陵地区与山前平原地带，林地、草地占比较大。衡水市、沧州市、廊坊市地处华北平原腹地，耕地面积占比特别大，故AESV较低；唐山市和沧州市AESV的提升得益于2市沿海滩涂滩地面积的扩大。

图2 2000—2018年京津冀地区AESV空间分布Fig.2 Spatial distribution of AESV in Beijing-Tianjin-Hebei Region from 2000 to 2018

根据式（5）得出京津冀地区各地级以上城市AESV变化率的空间分布，结果如图3所示。由图3可知，近18年京津冀地区大部分地市AESV呈减少趋势，表明生态系统仍然面临较大压力。具体来看，2000—2010年，84.62%的地市AESV呈减少趋势，但减幅不大，仅天津市减幅超过5%，15.38%的地市呈增加趋势，但是增幅不大，均小于5%。2010—2018年，38.46%的地市AESV呈增加趋势，61.54%的地市趋于减少，其中天津市、廊坊市和衡水市减幅超过5%，该阶段呈增加趋势的地市较上一阶段增加3个，分别是张家口市、沧州市和唐山市。天津市近18年水域面积大幅转换为建设用地，而水域的生态系统服务价值系数最高，因此AESV减幅较大；2010年以来张家口市的“三北防护林”“退耕还林还草”以及“京津风沙源治理工程”等工程建设成效显著，AESV开始呈现增长趋势；廊坊市紧邻京津地区，受地理位置和环京楼市影响，2010年以后大量耕地转换为建设用地，AESV减幅较大。衡水市土地利用类型单一，耕地面积占比特别大，2010年后大量耕地转换为建设用地，AESV减幅较大。

图3 2000—2018年京津冀地区AESV变化率空间分布Fig.3 Spatial distribution of AESV change rate in Beijing-Tianjin-Hebei Region from 2000 to 2018

2.3 生态系统服务价值对价值系数的敏感性分析

根据式（6），通过上下调整50%的生态系统服务价值系数来计算京津冀地区敏感性指数，结果如表4所示。

表4 京津冀地区生态系统服务价值敏感性指数Table 4 Sensitivity index of ecosystem services value in Beijing-Tianjin-Hebei Region

由表4可见，不同土地利用类型的敏感性指数差别较大，但不同年份之间差别较小。2000年各土地利用类型的敏感性指数高到低依次是林地、水域、耕地、草地、未利用土地；2010年依次为林地、耕地、草地、水域、未利用土地；2018年依次为林地、水域、草地、耕地、未利用土地。各土地利用类型的敏感性指数均小于1，表明京津冀的ESV是缺乏弹性的，相对于价值系数来说是相对稳定的，因此本研究所选取的价值系数适用于该研究区。其中林地的敏感性指数最高，说明林地对ESV影响最大，这是由于京津冀地区林地面积占比大，且林地生态系统服务价值系数也大。林地的敏感性指数为0.411 2～0.415 8，表明当林地的价值系数增加1%时，ESV将增加0.411 2%～0.415 8%。水域面积虽然占比不大，但其当量因子最大，故而敏感性指数较高。耕地的敏感性指数从0.197 6减小到0.177 3，这是因为耕地面积在逐年减少。

2.4 生态系统服务价值驱动因子分析

为了探测京津冀地区ESV的驱动因子，参照前人成果[2,11,28]，基于数据可获取性，选取年降水量（X1）、DEM（X2）、土壤侵蚀强度（X3）3个自然因子，GDP总量（X4）、地均 GDP（X5）、农民人均纯收入（X6）、社会消费品零售额（X7）4个经济因子，人均耕地面积（X8）、人口密度（X9）2个社会因子对京津冀地区的ESV进行探测。

基于京津冀地区13市研究尺度，以AESV作为变量，以上述影响因子为因变量，基于地理探测器（GeoDetector）工具的“因子探测”和“交互作用探测”功能进行驱动因子及驱动因子之间交互作用的定量分析，探讨京津冀地区ESV空间分异的主导因子，以及各驱动因子间交互作用的特征，最终得到地理探测器对京津冀地区AESV空间分异的因子探测结果（表5）和交互探测结果（表6）。

表5 京津冀地区生态系统服务价值驱动因子Table 5 Driving factors of ecosystem services value in Beijing-Tianjin-Hebei Region

表6 地理探测器的交叉探测结果Table 6 Cross detection results of geodetectors

分析显示，京津冀地区AESV在空间上的分布差异受到区域自然因子、经济因子和社会因子共同作用，其中高程、地均GDP、人口密度等因子的贡献率较大，而土壤侵蚀力、社会消费品零售额、人均耕地面积等对AESV空间分异影响的贡献率相对较小。从驱动因子对AESV的影响来看，DEM的解释力最为显著，为27.2%；其次是地均GDP和人口密度，其影响力均在20%左右，是影响区域AESV空间分异的重要因素；另外，年降水量、GDP总量和农民人均纯收入等因子影响力均在14%以上，是影响区域AESV空间分异较为重要的因素；土壤侵蚀力、社会消费品零售额、人均耕地面积的解释在10%以下，对ESV的空间分异影响较小。

交互作用探测器用来检验2种影响因子是自身独立起作用还是相互作用[30]。从表6可以看出，任意2个因子的交互作用均大于单个因子的影响，因此京津冀地区AESV空间分异是由不同驱动因子共同作用的结果。从交互探测的结果来看，DEM（X2）与地均 GDP（X5）、DEM（X2）与年降水量（X1）的交互作用程度达到40%左右，不同驱动力之间的交互作用对京津冀地区AESV的空间分异影响明显增强，不同因子之间的复杂耦合作用所形成的协同增强效应共同影响了京津冀地区AESV在空间上分异效果。

因此京津冀地区在生态系统优化与生态风险管控实践中，应考虑不同驱动因子的作用特点和各驱动因子交互协同增强的效果，采取差异化多元调控策略，选择与区域自然条件、社会经济发展水平相适应的土地利用开发模式，避免不合理或强烈人为土地利用干扰与自然、社会经济因子协同作用增强对区域生态系统的压力。

3 结论

（1）2000—2018年，耕地是京津冀地区最主要的土地利用类型，土地利用类型面积变化最大的为耕地和建设用地。耕地减少了9 742 km2，建设用地增加了9 992 km2。从单一土地利用动态度来看，土地利用幅度变化剧烈的是建设用地，变化幅度高达3.14%，变化幅度最小的是林地。综合土地利用变化度呈增长趋势。

（2）从数值分析，2000—2018年ESV呈先下降后上升趋势，18年间增加了66.33亿元。从生态功能构成来看，调节服务与支持服务始终是京津冀地区生态系统服务的核心功能，其中水文调节的贡献率最大。从时间上看，京津冀地区2000年、2010年、2018年的AESV分别为201.67万、199.27万、204.10万元/km2，京津冀地区生态系统服务趋于好转。从空间上看，京津冀地区AESV空间分布差异明显，高值区主要分布在冀北和京津地区，低值区主要分布在中南部的平原地区。京津冀各地市AESV的变化较为明显，大多数地市呈减少趋势，2010—2018年较前一阶段情况稍有好转。敏感性结果表明京津冀的ESV是缺乏弹性的，其中林地的敏感指数最高，对ESV影响最大。

（3）生态系统本身具有多样性与复杂性，加上ESV计算模型和方法多样，其价值评估结果存在较大差异，因此其计算方法及参考指标数据仍有完善的空间。对比京津冀地区已有的相关研究，本文所采用的当量因子测算的ESV与已有研究结果大致相同，研究结论合理。

（4）地理探测器分析表明，京津冀地区AESV在空间上的分布差异受到区域自然因素、经济因素和社会因素共同作用，其中DEM、地均GDP、人口密度等因子的贡献率较大，不同驱动力之间的交互作用对京津冀地区AESV的空间分异影响明显增强。