长三角地区生态系统服务价值动态演化及驱动力分析

方 林,蔡 俊①,刘艳晓,袁宏伟,李灿锋 (.安徽农业大学公共管理系,安徽 合肥 30036；.东北师范大学人文地理系,吉林 长春 3004；3.中国地质调查局昆明自然资源综合调查中心,云南 昆明 65000)

2018年习近平总书记提出，支持长江三角洲区域一体化发展并上升为国家战略，2019年12月国务院印发《长江三角洲区域一体化发展规划纲要》，将皖、苏、沪、浙全域纳入规划范围。作为中国经济最活跃、开放水平最高、人口最稠密的区域之一，长江三角洲区域城市拓展、人口涌入所产生的次生环境问题备受关注。践行“两山”理念，正确认识区域生态价值现状，实现区域生态环境共保联治成为亟待思考和完成的课题。

生态系统服务(ecosystem service，ES)指人们利用生态系统自有的结构、过程和功能直接或者间接获得维系生命所必需的物质及服务[1-2]。自1997年COSTANZA等[3]确定生态系统服务价值(ecosystem service value，ESV)分类及评估方法，并对全球ESV进行估算后，这类研究便成为地理、生态和土地等相关学科研究热点，其研究成果是衡量区域生态环境质量变化的重要内容。早期学者主要围绕ES价值定义[4]、分类体系[5]及评估方法[6-7]展开研究，联合国开展的千年生态系统评估项目将生态系统功能分为供给、调节、支持和文化服务，明晰了ES功能分类；我国谢高地等[7]在此分类基础上结合实地调查及专家打分得到中国陆地ESV系数表，为后续不同区域及尺度ESV评估奠定基础，并在中国范围内得到广泛使用和验证[8-10]。随着研究深入及技术进步，相关学者采用CA-Markov[11]，InVEST[12]、ARIES[13]和SoIVES[14]等模型评估和预测多种生态系统服务功能价值动态变化，部分学者通过RS和GIS手段，采用Moran′sI指数、热点探测[15]和LISA集聚图[16]等空间分析方法对不同尺度ESV时空演化特征展开分析。由于生态系统服务种类的复杂性、空间分布的异质性和人类使用的选择性，生态系统服务之间会呈现此消彼长的权衡关系或相互增益的协同关系[17]。系统认识生态系统服务权衡协同关系，对于区域生态治理具有重要现实意义。当前国内外学者围绕生态系统服务之间相互关系[18]、形成机制[19]和尺度效应[20]等不同方向展开了深入研究，在研究方法上多采用相关性分析、均方根偏差和约束线分析，其中，基于线性拟合建立的生态系统服务权衡协同模型能够定量反映不同生态功能之间权衡协同关系的演化趋势及程度，被广泛应用于相关研究[21]中。近年来，有些学者开始关注造成ESV时空演变及空间分异的影响因素及其驱动机制，但多采用地理加权回归[22]、关联度模型[23]、主成分分析和回归分析[24]等传统模型和方法，难以发掘不同因子对ESV空间分异的影响[20]。

综上，已有关于ESV的研究经历了从概念界定、方法探索到动态评估、时空特征分析再到影响因素、驱动机制的纵向深度拓展，取得了较为丰硕的成果[25]，但ESV时空动态研究亟待深化，不同生态功能相互关系尚未厘清，造成ESV空间分异的驱动背因有待发掘。自1978年以来，长三角地区经济持续增长[26]，而区域生态环境却出现一定程度退化[27]，为经济上行与生态下行的矛盾交互带。以长三角地区为研究对象，构建区域ESV估算模型及驱动因素指标体系，围绕ESV时空动态及各生态系统服务功能权衡协同关系进行量化分析，并基于地理探测器分析形成长三角ESV空间分异格局的驱动因素，以便为正确认识长三角地区生态价值现状、准确把握生态价值演变时空规律和实现跨区域生态环境共保联治提供理论参考及科学依据。

1 研究区概况

长三角地区包括苏、浙、皖、沪“三省一市”41个城市，是中国人口集聚和经济活跃程度最高的地区之一。国家统计局资料(http:∥www.stats.gov.cn)显示，截至2018年，域内人口约为2.2亿，占我国总人口的16%；区域GDP达22.12万亿元，占全国GDP的24.17%。区域地处长江中下游平原，地形相对平坦，水热条件优越，年平均气温约为17 ℃，年降水量约为1 368 mm，自然资源十分富足。而受人类活动及经济发展影响，近年来区域生态环境质量有所下降。长三角地区作为生态环境保护与经济社会发展冲突的热点区域，以此为研究对象，对于区域乃至全国生态环境治理具有重要借鉴意义。

2 材料与方法

2.1 数据来源及指标体系选取

以长三角地区为研究对象，基于中国科学院资源环境科学数据中心(https:∥www.resdc.cn)1990、2000、2010和2018年4期30 m×30 m分辨率的土地利用栅格数据，参照中国科学院土地利用变化遥感监测(China land-use and land-cover change，CNLUCC)分类体系，采用ArcGIS 10.2软件Reclass tools工具将4期土地利用数据分为林地、草地、耕地、水域、建设用地和未利用地6类(图1)。

图1 1990—2018年长江三角洲地区土地利用现状

根据已有研究成果[28-30]及研究区实际情况，构建包含自然地理和社会经济2个准则层16个指标因子的驱动力指标体系(表1)，以此探索长三角地区ESV空间分异的驱动因素，指标数据与ESV数据匹配均为2018年。对所得自然地理因素数据进行插值，重采样至30 m×30 m分辨率；社会经济因素通过数据清洗及空间匹配，形成基于长三角41个地级市数据集，进一步通过K-means分析按照自然断裂点法将原始数据做离散化处理，作为地理探测器模型输入数据。

表1 长三角ESV空间分异影响因素指标体系

2.2 研究方法

2.2.1ESV核算

ESV核算在谢高地等[7]研究成果基础上进行。为减少区域差异性对ESV评估结果的影响，对全国ESV系数表进行修正，构建长三角地区ESV估算模型。ESV当量因子系数是基于不同生态系统能够产生的生态系统服务的相对贡献率判定，其值约为研究区域平均粮食单产市场价值的1/7[7]。《中国农产品价格调查年鉴》数据显示，作为稻米主产区，1990—2018年研究区平均稻米产量为5 701.44 kg·hm-2，按同时期稻米平均市场价格为1.89元·kg-1计算，研究区修正后的ESV单位因子价值量为1 539.39元·hm-2。在此基础上计算可得长三角地区单位面积ESV系数(表2)，结合式(1)可得长三角地区ESV。由于建设用地ESV较低，故该研究未考虑。引入权衡协同度模型[31]〔式(2)～(3)〕，探究1990—2018年不同生态系统服务功能间的相互关系。为确保ESV系数修正结果对于ESV估算的可靠性，对长三角地区各土地利用类型ESV及其价值系数的相关性进行敏感性检验[32]。

(1)

式(1)中，VES为生态系统服务价值，元；Ai为第i种土地利用类型面积，hm2；CVi为第i种土地利用类型单位面积ESV系数，元·hm-2。

(2)

(3)

式(2)～(3)中，Sia和Sib分别为研究期初和期末第i种生态系统服务功能价值，元；Ii为第i种生态系统服务功能价值变化率，%；Dij为第i种与第j种生态系统服务功能的权衡协同程度。Dij<0，表示第i种与第j种生态系统服务功能呈负相关，表现为权衡关系；反之，Dij>0，表示第i种与第j种生态系统服务功能表现为同增或同减的协同关系[33]。

表2 长三角地区单位面积ESV系数表

2.2.2ESV空间分析

为探究长三角地区ESV在空间上的分布特征和演化规律，采用ArcGIS 10.2软件Create Fishnet功能创建3 km×3 km、5 km×5 km、7 km×7 km、10 km×10 km、15 km×15 km和20 km×20 km等不同尺度格网，用以表征ESV空间分布情况。经实验，研究区采用10 km×10 km格网尺度既能保证研究需求，且运算量相对较低[27]。采用ArcGIS 10.2软件分析长三角地区ESV空间相关性[34]，空间自相关分析用全局Moran′sI指数衡量。Moran′sI指数取值范围为[-1，1]，当Moran′sI指数大于0时，表示空间呈正相关，指数越高，ESV空间分布相关性越强，空间集聚越明显；当Moran′sI指数小于0时，表示空间呈负相关，指数越低，相关性越弱，空间差异性越大[2]。在此基础上，采用ArcGIS 10.2软件Getis-Ord Gi*热点探测工具分析1990—2018年长三角地区ESV冷热点动态演化规律，并绘制区域ESV热点和冷点空间分布图。

2.2.3ESV空间分异驱动因素

地理探测器作为探究事物空间分异特征及其驱动因素的一种新型统计学方法，相对于传统计量方法，因子探测可以通过q值来度量单个因素对区域ESV空间分异的驱动力度，交互因子探测则可进一步识别不同因素对区域ESV空间分异的交互影响[35]。目前地理探测器已被广泛应用于土地利用[36]、区域规划[37]和生态环境[25]等诸多领域，该模型在解释空间分异和探讨驱动机制方面具有独特优势，适用于从不同角度及深度探索造成长三角地区ESV空间分异的驱动因素，模型表达式为

(4)

式(4)中，q为影响因素对ESV空间分异特征的解释力度，其取值范围为0～1，值越大，自变量X对因变量Y的解释力度就越大，反之则越小；L为ESV空间分异影响因素X的分类或分区数，该研究X为5层；Nh和N分别为层h和研究区样本数，该研究样本数为41；σh2和σ2分别为层h和研究区ESV的方差。

3 结果与分析

3.1 研究区ESV变化分析

如表3～4所示，就1990—2018年长三角地区总ESV而言，1990—2000年小幅度增加，2000年之后逐步减少，呈先增后减、总体减少的变化趋势。其中，耕地ESV持续降低；林地ESV在1990—2000年增幅为3.42%，之后逐年减少，但减少趋势放缓，且总量有所增加；草地ESV在1990—2010年大幅减少，2010年之后小幅增加；水域ESV持续增加；未利用地ESV变幅在2000—2010年间最大。尽管长三角地区ESV整体呈下降趋势，但降幅已有明显变缓趋势，草地和水域ESV出现一定增长，这说明得益于近年来环境保护力度加强，研究区生态系统服务价值整体呈增长趋势。

1990—2018年长三角地区各项生态系统服务功能价值均有不同程度下降(表5)。

表3 长三角地区生态系统服务总价值(ESV)

表4 长三角生态系统服务价值(ESV)变化

表5 长三角单项生态系统服务功能价值(ESV)变化情况

保持土壤和废物处理服务作为区域主要生态系统服务功能，其价值下降最为明显，提供美学景观和原材料生产功能服务价值下降相对较少。1990—2010年研究区单项生态系统服务功能价值减少幅度较大，2010年之后，降幅减缓，这与长三角地区总ESV变化趋势一致。1990—2000、2000—2010和2010—2018年长三角地区各生态系统服务功能权衡协同比分别为14∶22、18∶18和0∶36(图2)。1990—2010年各生态系统服务功能权衡度下降，协同程度上升，2010年之后协同程度进一步上升，1990—2018年长三角生态系统服务功能间逐渐由冲突、竞争的权衡关系转变为同进同退的协同关系，协同关系占据完全主导地位，区域生态系统服务功能相互关系发生改变[38]。长三角区域单项ESV由权衡向协同的转变，并不意味着区域生态环境进一步恶化，相反，结合研究区整体ESV下降变缓趋势，可推测未来长三角地区ESV可能增长，生态状况进一步好转。

单位面积ESV系数可靠性对ESV评估结果具有重要作用[39]，为此对价值系数展开敏感性分析。如表6所示，长三角地区不同年份不同土地类型价值系数敏感性均小于1，即ESV结果对系数不敏感，上述研究结果可靠。

图2 长三角地区生态系统服务功能权衡协同关系

表6 长三角地区ESV敏感性检验结果

3.2 ESV变化空间自相关分析

基于10 km×10 km格网尺度，采用自然间断点法将长三角地区单位面积ESV划分为低值区(0～1 500元)、较低值区(>1 500～2 300元)、中值区(>2 300～3 200元)、较高值区(>3 200～4 700元)和高值区(>4 700～10 000元)(图3)。低值区主要分布在研究区北部以耕地为主的平原地带，较高值区主要分布在皖西、皖南和浙南等山地丘陵地带，高值区集中在江河湖泊地带。总体而言，1990—2018年长三角地区ESV呈现南高北低的分布格局。

全局空间自相关分析结果(表7)表明，各年份Moran′sI值均在0.6以上，说明长三角地区ESV空间分布呈显著空间正相关关系，且保持较高空间集聚度。

图3 1990—2018年长三角地区单位面积ESV空间分布

表7 长三角地区ESV空间自相关分析结果

热点探测结果(图4、表8)显示，冷点集聚区分布在皖北、苏北和上海周边地区，热点集聚区主要集中在大别山地区、江河湖泊附近和南部山区丘陵地带。1990—2018年冷、热点集聚区变化不显著，就冷点数量而言，冷点区数量呈先增后减趋势，2010年为冷点数量增长拐点，之后开始减少；就热点数量而言，1990—2000年呈增长趋势，2000年之后热点数量持续减少，但减幅变小。

图4 1990—2018年长三角地区ESV冷热点分布

表8 1990—2018年长三角地区ESV冷、热点数量变化

3.3 ESV空间分异驱动力分析

在95%置信区间下，长三角地区ESV空间分异不同因子解释力度从大到小排序为海拔(0.768)、日照时长(0.694)、年降水量(0.647)、土壤有机质含量(0.611)、年平均气温(0.443)、人口密度(0.407)和归一化植被指数(0.389)。海拔、日照和降水量为长三角地区ESV空间分异主要影响因子，土壤有机质含量、年平均气温和人口密度为次要影响因子；自然地理因素为区域ESV空间分异特征主要驱动因素，而社会经济因素影响程度不高，仅人口密度因子有一定影响。

从交互探测结果(图5)中剔除生态探测不显著部分，在95%置信区间下，各因子交互解释力度多高于0.500。与单因子对长三角地区ESV空间分异的影响相比，因子交互作用的影响更大。一产占比与海拔交互作用的影响最为明显，达到0.883。

自然地理和社会经济要素的内部交互作用弱于外部交互，即社会经济与自然地理要素间表现出极强的交互性，社会经济因素与日照时长、年降水量和海拔的交互解释力度大于0.700。因子交互类型多为双因子增强，表现为交互解释力度大于任意单因子的解释力度[35]，区域经济生产总值和道路密度与年平均气温、年日照时长、年降水量和归一化植被覆盖指数间呈非线性增强，交互解释力度大于单因子解释力度之和[35]，这说明区域经济生产总值和道路密度与其他因子产生交互作用时的解释力度比自身单独作用更大。

RH为相对湿度，AAT为年平均气温，ASD为年日照时长，AP为年降水量， NDVI为归一化植被指数，DEM为海拔，SOM为土壤有机质含量， DW为水系密度，GDP为经济生产总值，PPI为一产占比， DP为人口密度，UR为城市化率，PFE为财政支出， RD为道路密度，PSI为二产占比，PTI为三产占比。

4 讨论与结论

4.1 讨论

(1)1990—2018年长三角地区ESV呈现先增后减、整体下降的趋势，这与史慧慧等[40]研究结果基本一致。研究初期长三角腹地经济尚未完全起步，区域生态环境受扰动程度较小。2000年之后，在上海及周边城市的辐射带动下，区域整体经济发展迅速，ESV系数高的土地利用类型被ESV系数低的土地利用类型所侵占，建设用地面积显著增长，耕地流失，林、草地退化，进而导致区域生态质量下降。而近年来得益于地方政府切实管护，区域生态有好转迹象。1990—2010年长三角地区生态系统服务功能权衡度降低，协同度升高。2010年之后协同关系占据主导地位，当前长三角各生态系统功能间表现出“同进同退”“相互增益”的协同关系，这从侧面说明在当前区域ESV整体下降的背景下，加强区域生态环境保护的必要性与迫切性。

(2)长三角地区ESV空间分布呈现南高北低的分布格局，这是由于研究区北部以耕地为主，南部则多林、草地，一般而言，耕地ESV系数比林、草地小，易出现低值地带。长三角地区ESV呈现显著正向自相关，说明区域ESV空间分布并非随机，具有明显空间集聚效应，高值趋近，低值相邻。一般而言，ESV空间分布愈集聚，生态空间破碎化程度愈低，较高的正向集聚性有助于生态空间划分，也有助于长三角生态环境共保联治的实现。长三角地区ESV冷点数量先增后减，说明区域人类活动、经济发展所造成的高值损失正逐步降低，但受限于环境恢复的滞后性，当前热点数量仍在减少，但减幅变小，冷、热点数量年际变化在空间上符合区域ESV动态变化，也印证了区域生态环境总体向好的趋势。

(3)在ESV空间分异驱动力探测方面，单因子探测结果表明海拔、日照时长和年降水量是长三角地区ESV空间分异的主要驱动力，相对而言，社会经济因素对区域ESV空间分异的驱动解释力度较小，这与朱治州等[27]关于长三角核心城市群的研究结果一致，而与王波等[30]的研究结果有所差异，这主要是由于研究区域不同造成。银川市地处中西部地区，生态环境脆弱，相较于长三角地区更易受到人为社会经济因素影响。社会经济与自然地理因素的交互解释力度明显高于单因子，说明长三角ESV空间分异是自然地理因素与社会经济因素共同作用的结果[41]。区域自然禀赋造就ESV空间分异，而社会经济与自然地理要素结合会加速区域环境及生态系统结构与功能的改变，进而改变ESV空间分异，跨区域生态环境共保联治的实现需要统筹考量自然与社会2个方面因素。

该研究以长三角地区为研究对象，探讨了1990—2018年区域ESV动态演化、空间分异及其驱动因素，研究结果对长三角生态环境共保联治、区域一体化发展具有一定参考价值。但由于生态系统的复杂性，关于时空演变规律及其驱动机制的挖掘是项循序渐进的工作，地理探测器可以在空间尺度上较好地揭示区域ESV空间分异现象的驱动因素，但在时间尺度上对驱动因素的挖掘仍有待进一步深入，未来需对现有指标体系进行拓展，围绕时间及空间从多尺度探讨区域ESV演变的驱动机制。

4.2 结论

通过对1990—2018年长三角地区ESV动态演化、空间分异及其驱动因素展开分析，形成以下结论：

(1)1990—2018年长三角地区ESV呈现先增后减、总体减少的趋势，区域各生态系统服务功能权衡度明显降低，2010年之后协同关系成为区域生态系统服务功能间的主导关系。

(2)长三角地区ESV呈现南高北地的空间分布格局，呈现明显正向空间自相关；2000年之后区域热点数量持续减少，但减幅变小，冷点数量呈先增后减的演化态势。

(3)自然地理要素对长三角地区ESV空间分异的影响程度大于社会经济因素，海拔、日照时长和年降水量为主要驱动因子；社会经济因素与自然地理因素交互作用增强了各因子对ESV空间分异的解释力度。