基于生态系统服务和生态系统健康的生态风险评价——以长株潭城市群为例

欧阳晓,朱 翔,贺清云

1 湖南财政经济学院 湖南省经济地理研究所，长沙 410205 2 湖南师范大学资源与环境科学学院，长沙 410081

城市化所产生的土地利用和土地覆盖变化给生态系统带来了深刻的影响,逐步改变了生态系统的结构和过程,是生态系统服务丧失的主要驱动力,对生态系统造成重大风险[1- 7]。目前,与土地利用变化相关的生态风险问题已引起世界各国的广泛关注,如何降低城市化对生态风险的压力,实现城市可持续发展,已成为世界各地城市规划者和决策者关注的热点问题[8- 10]。

生态风险评估作为生态环境管理的前奏越来越普遍[11]。1992年由美国环保署提出了生态风险评价的概念,搭建了生态风险评价的研究框架,此后学者们围绕该研究框架在方法和理论等方面进行完善[12- 13]。目前,生态风险的评价主要是为生态环境保护和管理提供政策支撑,传统的研究框架主要研究一个或者多个因素对生态系统产生影响,难以综合的评价生态风险[14- 16],所以,学者们尝试将生态系统服务纳入到生态风险评价的框架中[17- 18],生态系统服务是连接生态系统和人类福祉的纽带[19],从人类福祉的视角结合生态过程和生态风险源进行生态风险的评价可以大幅度地提高时效性[20- 21]。但是,随着研究的深入,一些研究提出生态系统不仅应当提供多样化的生态系统服务,而且需要保持高水平的生态系统健康[22- 23]。健康的生态系统被认为是强调生态系统完整性的生态环境管理的目标,并为生态系统评估提供基础[24]。因此,在区域生态风险评价中,评价终点应将生态系统服务与生态系统健康相结合。然而,现有研究很少尝试建立生态系统服务与生态系统健康相结合的生态风险评估框架。

长株潭城市群是我国新型城镇化重点培育的对象,同时也是长江中游城市群的重要组成部分,是湖南省的经济发展、城市化、人口、环境污染的主要集聚区,当前正处于快速的城市化、工业化发展阶段。通过研究长株潭城市群能够有效解决城市群区域普遍存在的问题。同时,城市群的发展加速了农村景观向城市景观的转变,城市群的自然景观和环境状况变化尤为显著。基于此,以长株潭城市群为例,结合生态系统服务与生态系统健康,构建了生态风险评估框架,分析2000年和2018年长株潭城市群生态系统服务、生态系统健康和生态风险时空变化规律并运用相关分析法剖析生态风险程度与人工表面比率的关系,为实现城市群可持续发展和进一步优化城市群国土生态修复提供新的研究框架和理论参考。

1 研究区概况

长株潭城市群是长江中游城市群的重要组成部分,是湖南省的政治、经济、文化中心,正处于快速的城市化、工业化发展阶段。本文以最新的城市群区域规划以及相关城市群城镇密集区的研究成果为依据,结合区域的完整性,选取城市群都市区为本文的实证研究区,都市区包括长沙、株洲、湘潭三市的市区以及长沙县、株洲县、湘潭县[25- 26],如图1。2018年,都市区的行政区划面积为8629 km2,总人口为912.58 万人,GDP为9175.78 亿元,建成区面积为1038.51 km2,是整个城市群中城市化和生态环境变化最显著的区域,改善生态环境是当务之急,城市群先后实施了湘江流域的重金属治理、生态绿心保护区的“退二进三”等一系列的生态修复工程,取得了一定成效。

图1 研究区Fig.1 Study area

2 数据来源与研究方法

2.1 数据来源及处理

2000年和2018年研究区土地利用数据来源于中国科学院资源环境科学数据中心,分为林地、耕地、草地、水域、建设用地、未利用地等6个类型。归一化植被指数(Normalized Difference Vegetation Index, NDVI)来源于MODIS产品(MOD13Q1)。净初级生产数据(Net Primary Production, NPP)来源于蒙大拿大学数值动力学模拟团队(Numerical Terradynamic Simulation Group, NTSG)。相关数据的类型、分辨率以及来源,见表1。景观指数通过FRAGSTATS 4.2软件基于土地利用数据的栅格图计算得到各网格的景观异质性、景观连通性等一系列景观指数。考虑到本研究使用的数据存在多种空间分辨率,因此,将所有的栅格数据统一到与NDVI一致,对数据进行重采样得到统一的250 m栅格数据。同时,对所有数据利用大小一致的网格全面覆盖研究区,将研究数据分配到对应的网格中,即为直接网格化。本研究设定网格的尺度为1 km×1 km,共9113个网格。然后,按照研究数据的类别进行直接网格化,实现研究数据由行政区划尺度降到网格尺度,能够合理的避免统计指标按行政区划均匀分布的问题。

表1 长株潭城市群基础数据

2.2 研究方法

2.2.1生态风险评估框架

生态风险是人为活动造成的生态系统变化,包括生态系统状况恶化、生态系统服务减少等变化。其中,生态系统服务体现了与生态系统结构和过程相互联系的生态功能;而生态系统健康反映了空间实体的持续状态,涉及生态系统的结构、过程和功能等方面[27- 28]。生态系统服务是功能的体现,生态系统健康是保护的主体,生态系统服务与生态系统健康相结合,能满足生态安全评价终点的要求,可解释保护主体生态系统相关功能[12,29]。景观格局作为不同类型的生态系统镶嵌而成的地理空间单元[30],能够反映出生态系统结构和功能的状况,其变化与生态系统健康密切相关。同时,综合评价终点是生态系统服务与生态系统健康的结合,可以作为衡量生态风险的工具,也可看作是生态安全的表征,加之对生态系统健康评价已成为前沿方向[31- 32]。因此,考虑到生态系统的复杂性,从生态系统的完整性出发,基于土地利用类型的变化,将景观格局与生态学理论相结合,采用机器学习方法来评价生态风险程度的相关因素。

2.2.2生态系统健康评估

生态系统空间实体的生态系统健康(Ecosystem Health, EH)体现了在压力下维持健康结构、自我调节和恢复的能力,可分为三类：活力,组织,弹性[23,33]。计算公式如下：

(1)

式中,EH代表生态系统健康;EV、EO、ER分别代表生态系统活力、组织、弹性。

生态系统活力(Ecosystem Vigor, EV)是指生态系统的新陈代谢或初级生产力。本研究采用NPP对生态系统活力进行量化。已有的研究成果已经证明NPP是评价生态系统初级生产力的有效方法[34]。

生态系统组织(Ecosystem Organization, EO)是指由景观格局决定的生态系统的结构稳定性,特别是景观异质性和景观连通性都会影响生态系统的组织。选取景观异质性、景观连通性等景观指数,利用权重系数模型计算生态系统组织[35],其中,景观指数具体计算公式见表2。

表2 景观格局指数计算公式及解释

EO=0.4×LC+0.4×LH+0.2×IC

=(0.2×SHDI+0.2×SHEI)+(0.1×DIVISION+0.15×IJI+0.15×CONTAG)+0.2×FRAC

(2)

式中，EO指生态系统组织,LC指景观连通性、LH指景观异质性、IC指景观形态。根据现有研究成果[32,36]和长株潭城市群景观格局实际情况,景观连通性和景观异质性在景观格局中的地位相等且处于主导地位,故其权重均取0.4;景观形态在景观格局中的地位处于较低的状态,故其权重取0.2。为量化生态系统组织,进一步将LH、LC和IC分解到具体的景观指数,LH分为SHDI指数和SHEI指数,各占一半的权重;LC分为DIVISION指数、IJI指数以及CONTAG指数,IJI和CONTAG通常被视为景观连通性评价的核心指标,权重要高于DIVISION;IC直接用FRAC指数进行表征。

生态系统弹性(Ecosystem Resilience, ER)是指区域生态系统在受到外界干扰后恢复其原有结构和功能的能力[37]。采用植被覆盖度作为弹性分值来计算区域生态系统弹性,并采用Log函数对土地利用数据和植被指数进行标准化处理[38]。计算公式如下：

(3)

式中,ER代表区域生态系统弹性,ER值越大,区域生态系统的弹性能力越高;i为土地利用类型的总数;Si为第i类土地利用类型的面积;Pi为第i类土地利用类型的弹性分值,即植被指数均值。

2.2.3生态系统服务评估

本文采用谢高地等[39]制定的生态系统服务价值(Ecosystem Serivce Value, ESV)当量表,考虑到区域之间的差异性[40],本文采用2018年长株潭城市群晚稻的平均单产量和平均收购价格。经计算,2018年长株潭城市群晚稻的平均单产量为7146 kg/hm2,平均收购价值为2.83元/kg,基于此制定了长株潭城市群地区不同地类单位面积的生态系统服务价值系数表[41- 42](表3),进而对研究区的ESV进行计算。首先计算ESV当量因子价值量,计算公式如下：

(4)

式中,VCk为ESV当量因子价值量(元 hm-2a-1);P为长株潭城市群晚稻的平均收购价格(元/kg);Q为长株潭城市群晚稻的平均单产量(kg/hm2);n为年份数。

ESV计算公式如下：

ESV=∑Ak·VCk

(5)

式中,ESV为生态系统服务价值;VCk为各地类生态系统单位面积服务价值系数(元 hm-2a-1);Ak是第k类土地利用类型的面积(hm2)。

表3 各地类生态系统单位面积服务价值系数/(元 hm-2 a-1)

2.2.4生态风险特征分析

(1)生态风险特征值(Risk Characterization Value, RCV)计算。生态系统服务与生态系统健康的结合可以拟合出风险评估终点的公式,为平衡各评价单元保护主体与相关功能之间的关系,计算各单元生态风险的公式如下：

(6)

式中,RCV代表生态风险特征值;EH是生态系统健康;ESV是生态系统服务。

(2)生态风险分等定级。采用K-Means聚类方法将生态风险程度分为五个等级。基于平方误差之和,采用了迭代次数为1000次的K-Means方法,将2000年和2018年的风险程度分为5个等级[43]。公式如下：

(7)

式中,SSE表示平方误差的和;k表示k聚类;ci表示k聚类的中心;dis表示k与聚类x之间的欧氏距离。

(3)生态风险空间自相关。考虑区域生态风险的空间异质性,采用Anselin局部Moran′s I进行了生态风险的空间自相关性探讨[44],将区域内风险空间模式分为四种类型。

2.2.5生态风险相关性分析

为了识别人工表面比率(Artificial Surface Ratio,AR)的阈值,采用回归模型分析了人工表面比率与风险表征值之间的相互关系。根据风险评估中的人工表面比率、生态系统服务和生态系统健康三个要素,基于属性将三个要素归一化,归一化范围为0—1。初始时,通过对AR和RCV取平均值来表示RCV与AR之间的线性关系,然后模拟了两者之间的拟合模型,通过推导得到了阈值[45]。

(8)

式中,AR表示网格中人工表面比率;SAR表示网格中人工表面的面积;S网格表示网格的面积。

3 结果分析

3.1 生态系统服务的动态变化特征

2000—2018年期间,网格尺度下的生态系统服务平均值从2000年的1013.71万元下降到2018年的706.28万元,其中超过63%的网格生态系统服务出现下降。图2描述了生态系统服务标准化后的空间格局,可以看出2000年至2018年,生态系统服务高价值(0.8—1.0)区域和中高价值(0.6—0.8)区域的面积都显著下降;生态系统服务低价值(0—0.2)区域的面积显著上升。尽管城市群实施了一系列的生态修复工程,但由于耕地和林地等两大高价值的地类大面积地转化为低价值的人工表面,从而导致生态系统服务价值出现大幅度减少。2000—2018年,生态系统服务高价值区域主要分布在长株潭城市群的南部和北部的长沙县,因为这些区域覆盖了大面积的植被。生态系统服务的低价值区域主要分布在长株潭城市群的中部地区,因为长沙市市区、株洲市市区、湘潭市市区的城区扩张消耗了大量的自然资源,林地、水体等生态系统服务高价值的地类被城市用地侵占,最终形成了生态系统服务的低价值区域。

图2 2000年和2018年生态系统服务空间格局Fig.2 The Spatial pattern of ecosystem services in 2000 and 2018

3.2 生态系统健康的动态变化特征

图3展示了生态系统健康3个组成部分的变化趋势。生态系统活力呈明显下降趋势：2000年生态系统活力值范围为0—0.80,主要分布在0.20—0.70之间;2018年生态系统活力值范围为0—0.68,主要分布在0.10—0.80之间。平均而言,生态系统活力值从2000年的0.40降至2018年的0.30。生态系统组织价值呈轻微下降趋势：2000年生态系统组织价值范围为0—0.53,大部分值分布在0.20—0.34之间;2018年生态系统组织价值范围为0—0.52,大部分值分布在0.19—0.35之间。平均而言,生态系统组织价值从2000年的0.27略微降至2018年的0.26。生态系统弹性值呈轻微下降趋势,2000年的生态系统弹性值高于2018年。相比之下,2018年更多的网格位于较低的取值范围(0—0.60)。平均而言,生态系统弹性值从2000年的0.36降至2018年的0.32。总体而言,3个生态系统健康指数的值均呈现不同程度的下降,其中,生态系统活力值的下降幅度最大,其他两个指数的下降幅度较小。

图3显示了2000年和2018年生态系统健康价值的分布及其差异。2000年的生态系统健康价值主要分布在0.15—0.60之间;2018年的生态系统健康价值主要分布0.31—0.60之间。平均而言,生态系统健康价值从2000年的0.38降至2018年的0.34。

图3 生态系统活力、生态系统组织、生态系统弹性和生态系统健康的直方图Fig.3 The histogram of ecosystem vigor, ecosystem organization, ecosystem resilience, and ecosystem health

3.3 生态风险的动态变化特征

图4揭示了生态风险的空间格局演变特征。2000年,I、II、III、VI及V的风险程度等级网格数所占比例分别为30.43%、25.67%、18.06%、17.14%和8.70%。2018年,这一比例分别为27.12%、21.26%、21.23%、20.37%和10.02%。高风险程度等级(VI级和V级)的网格数比例上升了4.55%。与此相反,低风险程度等级(I级)的网格数比例下降了3.31%,下降区域主要集中株洲县、湘潭县以及长沙县。总体而言,随着区域城市化进程的加快,研究区的生态风险也随之增加。

图4 生态风险程度等级的空间格局Fig.4 The spatial pattern of ecological risk rating Ⅰ级: 低风险程度等级;Ⅱ级: 一般风险程度等级;Ⅲ级: 中等风险程度等级;Ⅳ级: 高风险程度等级;Ⅴ级: 较高风险程度等级

图5所示,两个时期主要存在两种集聚类型,高-高风险集聚区和低-低风险集聚区。高-高风险集聚区分布在城市群的中部,长沙市市区、株洲市市区和湘潭市市区;而低-低风险集聚区分布在城市群的北部和南部。高-高风险集聚区面积呈现一定程度的增加,低-低风险集聚区面积呈现下降趋势。

图5 生态风险等级的空间异质性Fig.5 The spatial heterogeneity of ecological risk rating

3.4 生态风险程度的相关分析与变化

2000—2018年期间,整个城市群人工表面比率存在明显的上升趋势,城市群的平均人工表面比率由2000年的14.96%上升到2018年的23.62%。评价网格的人工表面比率都处于增长阶段,其中2000年,人工表面比率大于10%的网格主要分布在长沙市市区、株洲市市区、湘潭市市区;2018年,人工表面比率大于10%的网格扩张到三个市区的周边地区。

生态指数(EV、EO、ER、EH、ESV、RCV)与人工表面比率的Pearson相关系数均为负值,说明人工表面的扩张是导致生态系统服务和健康下降的主要因素。人工表面和生态指标之间的相关系数分别为生态系统活力(2000年为-0.937,2018年为-0.902);生态系统弹性(2000年为-0.793,2018年为-0.926);生态风险值(2000年为-0.736,2018年为-0.782);生态系统服务价值(2000年为-0.605,2018年为-0.637);生态系统组织(2000年为-0.529,2018年为-0.458);生态系统健康(2000年为-0.173,2018年为-0.542)。

2000年和2018年,RCV与人工表面比率的关系与二次回归方程分别为0.595和0.662(图6)。随着人工表面比率的增加,RCV值迅速下降。基于回归模型的推导,为了进行风险管理,2000年和2018年人工表面比率的阈值分别为20%和36%。

图6 生态风险与人工表面比率的回归关系Fig.6 The regression relationship between ecological risk and artificial surface ratio

2000年AR、ESV、EH的Moran′s I值分别为0.534、0.473、0.456,2018年为0.562、0.452、0.461。结果表明,在2000年和2018年,长株潭城市群存在较强的空间自相关和明显3个主变量的聚类空间格局。MPLE(极大似然估计)的空间自逻辑模型系数如表4所示。在模型构建的基础上,生态系统服务对调节生态风险起到了重要作用,生态系统健康对缓冲风险也起到了一定的作用。然而,AR的压力对风险的引导起着主要作用,空间滞后对风险概率的提高起辅助作用。从系数值的变化来看,只有AR指数呈上升趋势,说明随着人工表面的扩展,对生态系统造成的压力在十年间不断增加。

表4 空间自逻辑模型的回归系数

4 结论与讨论

4.1 结论

基于2000年和2018年长株潭城市群土地利用类型、NPP、NDVI等数据,本研究结合生态系统服务与生态系统健康,构建了生态风险评估框架,对长株潭城市群生态系统服务、生态系统健康和生态风险时空变化规律进行分析,得到以下结论：首先,2000—2018年期间,生态系统服务减少率达到63%;生态系统健康的平均值从0.38下降到0.34。其次,2000—2018年期间,低风险程度等级(I级)的网格数比例下降了3.31%,而高风险程度等级(VI级和V级)的网格数比例上升了4.55%。最后,根据回归分析,人工表面比率的阈值应控制在36%以下,以达到生态风险管理目的。总体而言,该分析框架可以综合评价土地利用过程对生态系统的影响。

4.2 讨论

(1)生态系统健康评价指标的选择

生态系统健康是随着人类主观期望的变化而变化的模糊概念[46]。具体而言,生态系统结构稳定、功能完善,并且能够满足人类生产生活的需求就被认为是健康,反之亦然[32]。对于量化生态系统的属性,相应的指标不是固定的,而是以研究区的生态系统为导向[47]。以生态系统活力因子为例,估算生态系统活力可以用非光合植被生物量、生物土壤结皮、总初级生产量等指标[48],在本研究中,选取净初级生产力是因为它具有更好的可操作性和适用性[49]。因此,生态系统服务价值、生态系统组织、生态系统活力等生态指数被作为定量化指数来使用,以简化生态系统组织过程的复杂性,并从中提取信息[50]。

(2)研究方法的合理性和可行性

生态风险评估是评估人类活动对生态系统或其组成部分造成的潜在不良影响或风险的程度[51]。本研究的一个关键贡献是提出了关联生态系统服务和生态系统健康的生态风险的评估框架。该框架通过关联生态系统服务和生态系统健康,从生态学、遥感和数值模拟等学科的角度出发,选择合适的模型和指标,解决与生态相关的风险评估问题。同时,利用线性相关的方法确定了生态风险与城市化之间的定量关系,对人工表面比率阈值和生态风险模式提出建议可作为区域生态风险管理的补充,有潜力应用于其他地区快速的城市化带来的生态风险评估问题,可以作为区域生态风险的评价终点。然而,“活力-组织-弹性”生态系统健康模型中,涉及多个权重的设置,如何设置权值对最终结果有决定性的影响[52]。学者们根据各自研究区的特征,确定区域具有代表性的景观指数(例如聚集、多样性、破碎化),并以此作为权重设定的依据[32,36,51- 52]。因此,本研究的权重设置可能具有独特性和个体性,其他研究区运用该框架时应该进行调整。

(3)研究的局限性

总体而言,将生态系统服务与生态系统健康相结合可以剖析城市化对生态风险的复杂关系。与之前的研究相比[32,36,51- 52],本研究采用网格化方法,比以行政区划为单元的研究尺度更小。但该研究仍然存在一些局限性。1)对空间异质性明显的研究区域进行生态系统服务量化时,各土地利用类型的生态系统服务平均值的使用存在不确定性。2)虽然生态风险与人工表面比率之间线性关系一定程度上可以反映生态问题,但仍需考虑社会经济驱动力、生态系统功能、社会生态系统环境绩效之间的关系,因为这将有助于量化生态系统服务与生态系统健康之间的关系以及对生态系统的综合影响。3)该框架的应用面临着时空尺度上的限制,需要通过使用现场调查数据和分辨率更高的遥感图像来克服这些限制。