参考已有研究^[22-23],构建城镇扩张强度指数,见公式(1):

)

式中:P为t至t+n时段内研究单元的城镇扩张强度; U_t,U_t+n分别为研究时段初期与末期的城镇用地面积; n为研究时段长度; A为研究单元总面积。

谢高地等^[19,24]在Costanza等的研究基础上,利用专家知识法,构建并修正了中国生态系统单位面积生态服务价值当量表,将生态系统服务分为供给服务、调节服务、支持服务、文化服务4个一级类型和9个二级类型。1个标准当量因子的定义为1 hm²全国农田年平均自然粮食产量的价值,考虑到人为因素干扰难以消除,在实际计算中可用单位面积农田粮食生产的净利润作为标准当量因子^[21]。根据《中国粮食发展报告》计算得到2000—2015年全国单位面积粮食生产净利润为1 837.46元/hm²,全国粮食平均单产为4 856.97 kg/hm²,考虑到粮食生产能力的空间差异性,通过国家统计局数据计算长江经济带15 a间粮食平均单产为5 058.57 kg/hm²,基于此对全国单位面积粮食生产净利润进行修正,得到15 a间长江经济带标准生态系统服务当量因子价值量为1 913.72元/hm²,并构建长江经济带单位面积生态系统服务价值量表(表1)。

表1 长江经济带单位面积生态系统服务价值量元/(hm2·a)

基于各研究单元土地利用数据与生态系统服务价值量表,计算出研究单元的各项生态系统服务价值量,并构建生态系统服务价值量下降率与变化率指标,公式^[25]如下:

式中:ESV_i,ESV_j,ESV分别为研究单元第i类用地的生态系统服务价值、第j项生态系统服务类型的价值以及生态系统总服务价值; A_i为研究单元第i类用地面积; VC_ij为表1中第i类用地第j项服务的价值系数; D与C分别为研究单元价值量下降率与变化率; ESV_t+n和ESV_t分别表示t+n时刻与t时刻研究单元的生态系统服务价值量。

式中:CS_i为第i种土地利用类型的敏感性指数; ESV₁和ESV₂分别为调整前后的生态系统服务总价值量; VC_1i和VC_2i分别为调整前后第i种土地利用类型的价值系数。若CS_i>1,表明ESV对价值系数具有弹性; 若CS_i<1,则说明ESV对价值系数缺乏弹性,计算过程是可靠的。

Getis-Ord G^*指数可分析研究单元属性的空间分布格局,公式^[26-27]如下:

式中:X_j为第j个研究单元的属性值; n=128; W_ij(d)为以距离d为规则的空间权重; E(G_i^*),var(G_i^*)分别为G_i^*的期望与方差,当Z值显著大于0时,表示第i个研究单元周围都是高值,为热点区,反之则为冷点区。

时空地理加权回归(GTWR)模型在地理加权回归(GWR)模型的基础上引入时间因子,既解决了截面数据样本数量有限性的问题,还考虑了时间与空间非平稳性,能有效估计因子参数,模型如下^[28-29]:

式中:(u_i,v_i,t_i)为第i个样本单元的时空坐标; X,Y分别为解释变量与被解释变量; p为解释变量个数; β₀(u_i,v_i,t_i)为截距项; β_k(u_i,v_i,t_i)为第k个解释变量的估计系数; ε_i为模型残差。本文对GTWR模型的实现主要基于ArcGIS 10.2软件,采用Huang等^[30]制作的GTWR插件,带宽采用AICc优化设置,时空距离参数比值为1。

敏感性指数能验证ESV对价值系数的依赖性,另一方面也能反映各地类对生态系统服务总价值的影响^[25]。将2000—2015年长江经济带各用地类型价值指数VC分别调整50%,根据公式(7)计算其CS值,结果见表2。总体而言,2000—2015年长江经济带各土地利用类型敏感性指数均小于1,说明生态系统服务价值量对于各地类价值系数不敏感,生态系统服务价值量的计算是合理的。其中,林地的敏感性指数相对较大,稳定在0.69,说明林地对于长江经济带整体的生态系统服务价值具有重要的作用。

基于生态系统服务价值当量表,计算出2000—2015年长江经济带各类生态系统服务价值量见表3。在2010—2015年,长江经济带的平均生态系统总服务价值量为73 186.25亿元,其中供给服务、调节服务、支持服务与文化服务分别占比为11.21%,52.59%,29.31%,6.89%。在15 a间,供给服务与支持服务价值量始终呈下降的趋势,而调节服务、文化服务和总服务价值量则在2000—2010年增加,在2010—2015年开始显著下降。总体而言,15 a间长江经济带生态系统的总服务价值量下降了0.27%,其中供给服务、调节服务和支持服务的价值量分别下降了0.71%,0.09%,0.53%,而文化服务则增加了0.12%。

表2 敏感性指数计算结果

对2000—2015年各地区生态系统服务价值量的下降率进行热点分析见图2。供给服务价值量下降的热点区集中于长江经济带下游长三角地区,该区域范围内各研究单元的供给服务价值量均呈现显著下降特征,其中上海市、嘉兴市、苏州市、无锡市分别下降了14.32%,13.81%,13.02%,9.58%; 而供给服务价值量下降的冷点区主要集中在长江经济带上游部分地区。调节服务价值量下降的热点区集中在长三角东南部地区,其中上海市和嘉兴市分别显著降低了9.40%,5.61%; 冷点区集中在长江经济带中游湖北省与湖南省的部分地区。支持服务价值量下降的热点区也集聚在长三角地区,其中上海市、嘉兴市、苏州市分别降低了12.61%,12.32%,8.70%; 冷点区分布在长江经济带上游部分区域。文化服务价值量下降的热点区分布在浙江省与江西省部分地区,其中台州市和绍兴市分别显著降低了2.49%,2.09%; 冷点区则分布在中游湖北省与湖南省部分地区。生态系统服务总价值量下降的热点区集中在长三角地区,其中上海市、嘉兴市、镇江市总价值量分别降低了10.16%,8.00%,5.34%; 冷点区则集中在湖北省与湖南省部分地区。

表3 长江经济带2000-2015年的生态系统服务价值量亿元/a

图2 2000-2015年长江经济带生态系统服务价值量变化率的热点分析

考虑到长江经济带各地区的城镇扩张强度与生态系统服务价值量变化率的时空分布具有非平稳性的特征,采用GTWR模型从局部尺度来分析城镇扩张对价值量变化的影响。利用3个时段的面板数据分别以生态系统供给服务、调节服务、支持服务、文化服务与总服务的价值量变化率为被解释变量,城镇扩张强度为解释变量,进行OLS回归和GTWR回归,并对比其拟合结果。其中,文化服务的回归结果不满足F检验和t检验要求,说明城镇扩张对文化服务价值量变化的线性驱动作用并不显著。其他生态系统服务类型的回归参数见表4。

表4 GTWR和OLS的模型诊断结果

基于R²,RSS与AICc指数,对OLS与GTWR的回归结果进行对比。其中,GTWR模型的拟合度大于OLS的拟合度,而GTWR的残差平方和指数和AICc指数则显著低于OLS,说明采用GTWR模型能比OLS模型更有效地解释城镇扩张对生态系统服务价值量变化的影响作用及其时空异质性特征。

将局部样本的回归系数进行可视化,见图3。总体而言,大部分地区的回归系数为负数,说明随着城镇扩张强度的增大,各地区生态系统的供给服务、调节服务、支持服务及总服务价值量的变化率减小,城镇扩张对长江经济带大部分地区的生态系统服务价值具有显著的负面影响。对于不同的生态系统服务类型,回归系数的大小具有不同的空间分布特征,且在时间序列上发生变化。

对于生态系统供给服务价值而言,回归系数的绝对值较大的区域主要分布在长江经济带的上游和下游地区。在2000—2005年和2005—2010年,城镇扩张对长江经济带上游云南省与贵州省的南部地区和下游长三角地区的负面影响强度较大,主要源于城镇扩张减少了耕地与林地等能提供重要供给服务功能的用地,从而导致生态系统的供给服务价值量显著降低,而中游大部分地区受影响程度较小。2010—2015年,受城镇扩张负面影响较大的区域主要集中在长江经济带的上游地区和中游西部的少数地区,其耕地、林地及草地逐渐向建设用地和水域转化,供给服务价值量逐渐下降,而下游大部分地区受影响较小。

对于生态系统调节服务价值而言,回归系数绝对值较大的区域主要集中在长江经济带上游和中游地区。2000—2005年,长江经济带上游云南省和贵州省的南部地区以及中游湖南省的部分地区受城镇扩张的负面影响较大,主要源于城镇建设对林地、草地与水域的占用,而下游地区受城镇扩张的影响程度相对较小。2005—2010年,城镇扩张对长江经济带上游大部分地区的调节服务价值具有较强的负面影响,而云南省西部德宏傣族景颇族自治州等地区,受益于“三江”流域实施的天然林保护、退耕还林和石漠化治理等重要生态保护工程^[31],其调节服务价值增加,回归系数为正。2010—2015年,由于澜沧江、金沙江等梯形水电站工程^[32]的建设,云南西部部分地区水域面积增大,城镇扩张对该地区的调节服务价值具有正面影响,对长江经济带上游其他地区和中游西部部分地区仍具有较大的负面影响,对下游地区的负面影响较小。

对于生态系统支持服务价值而言,长江经济带上游南部地区回归系数的绝对值较大。2000—2005年,由于耕地和林地、水域等生态用地被破坏,长江经济带上游云南省与贵州省的南部地区和下游北部地区的支持服务价值受城镇扩张的影响较大,中游地区受影响较小。2005—2010年,城镇扩张对支持服务影响较大的地区主要集中在贵州省南部。到了2010—2015年,受城镇扩张负面影响较大的地区主要分布在长江上、中游以南,该地区的耕地、林地、草地和湿地不断向水域和建设用地转化。

对于生态系统总服务价值而言,回归系数绝对值较大的区域主要集中在长江经济带的上游地区。2000—2005年,由于大量耕地和林地等生态用地被建设用地侵占,城镇扩张对总服务价值负面影响较强的地区主要分布在云南省与贵州省南部。2005—2010年,受城镇扩张负面作用较强的地区收缩至贵州省南部。到了2010—2015年,城镇扩张对长江经济带上游大部分地区都具有较强的负面影响,对云南省西部部分地区则具有正面影响,对中游和下游的大部分地区的负面影响较小。