刘应帅,余 瑞,郑彬彬,刘嘉慧,宋 奇,陈荣昊,严 哲
(1. 海南大学 生态与环境学院,海口 570228; 2. 海南集思勘测规划设计有限公司,海口 570203)
定量描述陆地生态系统碳动态变化及影响因素对区域生态系统碳循环的研究很重要[1−2]。而森林净生态系统生产力(Net Ecosystem Productivity,NEP)的大小能够表征生态系统的固碳能力,它最初由WOODWELL 等[3]在分析陆地生物圈源、汇问题时提出,其生态学的含义表示为生态系统总的生产力与生态系统呼吸之差。通过NEP 来指示生态系统的固碳状态,NEP 为正值,说明生态系统呈现碳汇状态,反之则为碳源。陆地生态系统一直扮演着强大的碳汇角色[4],而陆地碳汇其中的很大一部分位于热带森林区域。热带森林覆盖面积约占全球表面积的10%,其森林碳汇约占陆地总碳库的25%,并且占全球植被碳储存的50%[5−6]。由于热区异常气候频发,热带森林生态系统碳汇的估算对于气候的敏感性的分析研究存在较大的不确定性[7]。IPCC 的报告指出,气候变化会对全球大部分区域的水文循环造成影响,进而导致区域干旱频发,这在亚马逊地区已经得到证实[8−10]。因此,对关于亚马逊热带森林的干旱敏感性问题产生了一些争论。一些学者认为,亚马逊热带森林旱季呈碳汇状态,雨季则为碳源[11];反驳者认为,热带森林长期的维持碳汇功能,但在遭遇干旱事件时会逆转为碳源[8,12]。从这些争论中可以看出,气候因素会对热带森林植被碳汇大小产生影响。因此,对热带地区森林植被碳汇进行气候驱动因素分析极为重要。温度和降雨是影响生态系统净
生产力的最主要的2 个气候因素[13−15],它能通过影响植物的光合作用和呼吸作用[16−17],进而影响到NEP。当然还有其他的因素在调控,如土地利用变化、土壤水分、CO2施肥效应以及氮沉降效应等[18−21]。例如CO2升高对生产力的潜在影响可能受到全球变暖和降水模式改变的影响,显示出了水的可用性和温度在全球植被光合作用和呼吸作用中承担着强大的驱动因素的角色[22]。一方面,热带林木长期处于热稳定的环境下,未来气温的持续升高,很有可能导致热带树木的热不育状态[15];另一方面,LIU 等[23]的研究表明了降水阈值对于生态系统生产和呼吸的调节作用。因此,了解热带森林生态系统净生产力对于气温和降雨的响应至关重要。当然,也有研究表明,地形因素会影响到气候因素,同一块区域较高海拔的地方温度和降雨会与低海拔地方呈现不一样的模式,植被类型也会有很大的区别[24−25]。
本研究结合多种遥感数据产品,首先利用时间序列分析和一元线性回归的方法,通过季节性分解得到了NEP 的年际趋势、年内趋势以及季节性趋势。其次基于降水的年内变化划分了干季(1−3 月),湿季(7−9 月),并从时空层面分析了不同季节的NEP 变化,以掌握NEP 的时空动态特征,并采用增强回归树的方法,探索了气候因子(气温、降水)、地形因子(海拔、坡度)对于NEP 的贡献程度,在此基础上,分海拔讨论不同季节NEP 趋势和降水趋势、温度趋势的相关性,进一步确定海南岛森林NEP 对气候因素的响应,旨在为理解全球气候变化背景下区域碳循环研究提供资料。
1.1 研究区概况海南岛位于中国的南部,热带北缘,地理位置介于108°03′~111°03′E 和18°10′~20°10′N 之间。海南岛的气候类型属热带海岛季风气候,长夏无冬,1 年分干湿2 季,全年平均温度22.5~25.6 ℃,年均降雨量900~2 500 mm[26−27]。海南岛拥有较高的森林覆盖度,其中,热带山地雨林是海南岛热带森林植被中面积最大、分布集中的垂直自然地带性的植被类型[28],主要分布在吊罗山、五指山、霸王岭、尖峰岭以及鹦哥岭等600 m以上海拔的山地。而低地森林主要分布在600 m以下中低海拔区域,包括红树林、橡胶林等[29]。
1.2 土地利用数据集及海南岛森林提取土地利用数据集源自中国研制的30 m 空间分辨率全球地表覆盖数据GlobalLand30[30],共包含2000、2010、2020 年3 期的土地利用数据(http://www.globallandcover.com/)。基于海南岛的边界图裁出海南岛的土地利用数据集,然后统计近20 a 海南岛的林地区域,绘制出稳定的林地区域图像,即在以上3 组土地利用数据中均为林地的区域会被提取出来,形成稳定的林地分布图(图1),再将其分辨率重采样至0.01°以匹配NEP 数据的空间分辨率。
图 1 研究区概况图(a)2020 年海南岛土地利用图;(b)基于3 期的土地利用提取的林地区域图;(c)海南岛高程图。
1.3 NEP 数据集本研究使用的NEP 是来自全球环境研究中心(Center for Global Environmental Research,CGER),该产品提供了1999—2018 年空间分辨率为0.01°,时间分辨率为10 d 的全球净生态系统交换(Net Ecosystem Exchange,NEE)数据[31](https://db.cger.nies.go.jp/DL/10.17595/20200227.001.html.en)。初始NEE 数据格式被转换为栅格,取相反数即为本研究所使用的NEP 数据。 全年共分为36 幅栅格图像,再分别将每10 d 的NEP均值数据合成月尺度数据和年尺度数据。
1.4 气候因子数据集气温和降水数据来源于国家科技基础条件平台—国家地球系统科学数据共享服务平台−黄土高原科学数据中心(http://loess.geodata.cn),该数据集是根据CRU 发布的全球0.5°气候数据集以及WorldClim 发布的全球高分辨率气候数据集,通过Delta 空间降尺度方案在中国地区降尺度生成的,并使用496 个独立气象观测点数据进行了数据验证[32]。
1.5 Digital Elevation Model(DEM)数据集
DEM 数据集来自中科院资源环境科学与数据中心(https://www.resdc.cn/),该数据集为基于最新的SRTM V4.1 数据经整理拼接生成的90 m 的分省数据。数据采用WGS84 椭球投影。随后将其空间分辨率采样至0.01°,并基于1.2 中提取的林地区域进一步提取林地所处的海拔高度。笔者使用ArcGIS 中的空间分析工具根据海拔计算出坡度和坡向。
1.6 NEP 和降水的时间序列分解针对存在季节性变化的月度NEP,采用季节性分解,将其分解为趋势因子(Trend Component,TC)、季节性因子(Seasonal Component,SC)和 误 差 因 子(Error Component,EC)。TC 能准确把握数据的长期变化;SC 能捕捉到数据1 年内的周期性变化;而EC能反映那些不能被趋势或季节效应解释的变化[33]。通过相加模型可以表示为:
式中: NEPt指 每月的NEP,T Ct、S Ct、 E Ct分别指趋
势、季节性以及随机误差。
线性滤波器是估计时间序列趋势常用的方法,最常见的线性滤波器之一是滑动平均。针对NEP 月度数据,选择了12 点移动平均法,线性滤波器如下:
式中: TCt指去除季节性效应的NEP 趋势值;1999—2018 共20 年,每年12 个月,合计240 个月,t 从第7 月开始取值。
在获取趋势因子之后,通过下式可得到季节性因子:
式中:S Ct为第t月时算得的季节性因子。
式中:Slope 是NEP 逐像元线性回归方程的斜率;代表年数,n 为时间跨度。当Slope>0 时,NEP呈增加趋势;当Slope=0 时,NEP 基本稳定,无明显变化;当Slope<0 时,NEP 呈减少趋势。
这个季节性因子的估计也包含了每个t 时间的随机误差因子 ECt,通过计算每个月季节性效应估计的平均值并在所有年份重复此序列来获取整体的季节性效应。
式中: i 为年份,j 为月份, j 取1~12, N EPij为第i 年第 j 月的平均NEP,N EPj为20 年的j 月平均NEP。
1.7 基于降雨效应的干湿季划分及趋势分析
基于像元尺度的趋势分析法能模拟研究区中每个栅格单元的变化趋势,从而反映NEP 变化的方向和速率。计算公式为:
1.8 基于BRT 的相对贡献率分析使用BRT 来评估1999—2018 年气候因子与地形因子对于海南岛森林NEP 变化的相对影响程度。BRT 分析具有容纳任何数据分布的能力,因此,在分析过程中无需进行数据的转换。在进行BRT 分析之前,对森林NEP 与气候因子以及地形因子进行皮尔逊相关分析和显著性分析。使用R 语言中的GBM包进行BRT 分析,将1999—2018 年平均NEP 的逐像元数据作为响应变量,同一时期的气候因子和地形因子作为解释变量。BRT 的本质就是1 个加法模型,每次建立模型是在之前建立模型损失函数的梯度下降方向。
1.9 海拔划分及趋势相关性分析基于海南岛森林的分布情况[29],以600 m 海拔为界限,将海南岛森林分为低海拔森林和高海拔森林分别研究。笔者计算了近20 年NEP 趋势与气候因子趋势的相关性,如下所示:
式中: xi为1999—2018 年温度或降雨的逐像元趋势, yi为1999—2018 年NEP 的逐像元趋势,为像元数。随后对NEP 趋势、温度趋势以及降雨趋势构建了多元线性回归模型,并以回归系数来估计各因子的贡献度,进而确定不同条件下对NEP 影响的主导因素。
2.1 NEP 和降雨的时间序列分解对1999—2018 年月度NEP 数据进行时间序列分析,发现海南岛森林NEP 随时间的变化趋势及季节性波动。通过使用滑动平均滤波器,得到了近20 年去除季节性影响的较为平稳的NEP 变化趋势(图2)。总体来看,NEP 在20 年间变化起伏波动,在2010 年前后出现最大趋势,而在2016 年出现最低趋势。笔者进一步分析了NEP 年内的季节性效应,发现年内NEP 呈现先增后减的季节影响,1−5 月,NEP 逐渐增大,在5 月达到最大值,随后的6−12 月,NEP 逐渐减小。同理,可以得到降雨的年内季节性效应。1−9 月,降雨随月份逐渐增加,在9 月达到最大值,随后的10−12 月,降雨逐渐减少。基于此可以将海南岛1−3 月划分为干季,7−9 月划分为湿季。
2.2 海南岛森林NEP 的年度及季节变化通过计算每年海南岛森林NEP 的均值,得到 1999—2018 年逐年NEP 的年际变化趋势图(图3-a),整体来看NEP 随时间呈现不显著的下降趋势,有机碳的变化率为−0.57 g·m−2·a(P>0.05),分段来看,前10 年NEP 呈现不显著的增长趋势,有机碳的变化率为3.3 g·m−2·a(P>0.05);后10 年NEP 呈现显著下降趋势,有机碳的变化率为−8.4 g·m−2·a(P<0.05)。近20 年NEP 的均值(有机碳)为483.23 g·m−2·a,其中,有9 年的NEP 高于均值。在2010 年,NEP 达到最大值,为534.68 g·m−2·a。在2016 年,NEP 有最小值,为439.47 g·m−2·a。
同样,在1999—2018 年逐年干季和湿季的NEP 趋势变化图中(图3-b)可以看到湿季NEP 有不显著的下降趋势,而干季NEP 有不显著的上升趋势。湿季NEP 的均值为136.76 g·m−2(有机碳),其中有11 年NEP 超过均值,最大值出现在2005年,为156.21 g·m−2,最小值出现在2015 年,为117.10 g·m−2。干季NEP 的均值为102.09 g·m−2,仅有9 年NEP 超过均值,最大值出现在2009 年,为129.54 g·m−2,最小值出现在2005 年,为81.24 g·m−2。干、湿2 季对比,除2009 年干季NEP 超过湿季NEP,其余年份湿季NEP 的值均高于干季。
图 2 时间序列趋势图(a)1999—2018 年NEP 月趋势图;(b)经过时间序列分解去除季节性影响后的趋势图;(c)NEP 的季节性变化图;(d)随机误差;(e)季节性影响下NEP 的年内变化;(f)季节性影响下降雨的年内变化。
图 3 不同年份和季节下NEP 的变化(a)NEP 的年际变化;(b)干季和湿季NEP 的变化。
2.3 海南岛森林NEP 的空间分布及趋势变化
通过绘制全年、干季以及湿季的NEP 空间分布图(图4-a、b、c)来观察NEP 在空间方位上的分布。年际的NEP 最高值为825.30 g·m−2·a(有机碳),最低值为−54.64 g·m−2·a,其中NEP 的高值集中在海南岛中部和东南部,而在海南岛北部以及海岸线附近多为NEP 低值聚集。此外,大部分地区的NEP 均大于0。干季的NEP 最高值为236.03 g·m−2,最低值为−69.89 g·m−2,与年际NEP 相似,干季NEP 的高值多分布在海南岛中部五指山、霸王岭以及尖峰岭一带,NEP 低值则聚集在海南岛东北部。湿季NEP 的分布较为破碎化,高值和低值相间,其中最高值为211.84 g·m−2,最低值为−19.74 g·m−2,湿季NEP 最高值较干季减少了24.19 g·m−2,但最低值较干季增加了50.15 g·m−2。
图 4 1999—2018 年NEP 的空间分布及趋势变化(a)、(b)、(c)分别为年际、干季、湿季NEP 的空间分布;(d)、(e)、(f)分别为年际、干季、湿季NEP 的空间趋势变化。
应用一元线性回归对海南岛NEP 近20 年的年际趋势变化进行分析,进而绘制了NEP 的趋势变化空间分布图(图4-d, 图4-e, 图4-f)。年际的趋势变化如图4 所示,NEP 呈现显著下降趋势的区域占比为18.49%(P<0.05),NEP 呈现显著增长趋势的区域占比为17.99%(P<0.05),NEP 无显著变化的区域占比为63.52%。从分布来看,NEP 显著增长区域主要集中在海南岛东北部,而显著下降区域主要集中在海南岛中部。干季的趋势变化:NEP 呈现显著下降趋势的区域占比为6.39%(P<0.05),NEP 呈现显著增长趋势的区域占比为23.05%(P<0.05),NEP 无显著变化的区域占比为70.56%。在空间分布上,干季的NEP 趋势变化与年际的趋势变化相似,但是呈显著下降趋势的区域明显减少。湿季的趋势变化:NEP 呈现显著下降趋势的区域占比为23.66%(P<0.05),NEP呈现显著增长趋势的区域占比为12.65(P<0.05),NEP 无显著变化的区域占比为63.69%。与年际趋势相比,湿季的NEP 显著增长趋势占比降低了5.43%,显著下降趋势占比提高了5.17%;与干季趋势相比,湿季的NEP 显著增长趋势占比降低了10.40%,显著下降趋势占比提高了17.27%。
2.4 基于BRT 的各因素贡献度分析通过对不同季节的NEP 与气候因子、地形因子进行皮尔逊相关分析以及显著性检验,发现年际NEP 与气温、降雨、海拔以及坡度均有极显著的相关性(P<0.01),相关系数分别为−0.53、−0.03、0.57、0.50,而与坡向无显著相关性(图5-a);干季NEP 与气温、降雨、海拔以及坡度均有极显著的相关性(P<0.01),相关系数分别为−0.54、−0.07、0.73、0.66,而与坡向无显著相关性(图5-b);湿季NEP 与气温、降雨、海拔、坡度以及坡向均有极显著的相关性(P<0.01),相关系数分别为0.07、−0.04、−0.09、−0.10、−0.03(图5-c)。基于以上分析,选择气温、降雨、海拔和坡度作为自变量因子用于BRT 分析,以便计算各因素在不同季节对NEP 的贡献程度。在进行参数调优后,在R 语言中使用GBM 包进行BRT 分析,得到不同季节不同因子对NEP 的相对贡献率。如图6-a 所示,对年际NEP 影响程度最大的是海拔,其相对贡献率为45.46%,其次分别为降雨、气温和坡度,相对贡献率依次为25.78%、16.17%、12.59%;如图6-b 所示,对干季NEP 影响程度最大的是海拔,其相对贡献率为40.58%,其次为坡度、降雨和温度,相对贡献率依次为28.13%、16.12%、15.17%;如图6-c 所示,对湿季NEP 影响程度最大的是降雨,其相对贡献率为30.75%,其次为温度,海拔和坡度,相对贡献率依次为26.77%、21.88%、20.60%。除湿季外,海拔在年际NEP 和干季NEP 中扮演着重要的角色。
图 5 相关性热力图pre、tem、dem、aspect、slope 分别表示降雨、温度、海拔、坡向以及坡度;(a)、(b)、(c)分别表示年际、干季和湿季3 个时期。
图 6 基于BRT 的贡献度分析pre、tem、dem、slope 分别表示降雨、温度、海拔以及坡度;(a)、(b)、(c)分别表示年际、干季和湿季3 个时期。
2.5 不同海拔区域NEP 的驱动因子分析为进一步说明海南岛森林NEP 的主要影响因素,分别对高海拔区域和低海拔区域进行了NEP 趋势和气候因子(气温、降雨)趋势的相关性分析。结果显示,在海拔600 m 以上的森林生态系统中,全年的NEP 趋势和湿季的NEP 趋势与降雨趋势均有显著的负相关性,相关系数分别为−0.11 和−0.19(图7-d, 图7-f),而干季的NEP 趋势则与降雨趋势没有显著的相关性(图7-e)。不论是年际的NEP趋势还是干季、湿季的NEP 趋势均与温度的趋势没有显著的相关性(图7-a, 图7-b, 图7-c)。而在
海拔600 m 以下的森林生态系统中,则呈现出十分不同的模式。年际的NEP 趋势对于温度趋势有微弱的显著负相关性,相关系数为−0.02(图8-a),对于降雨趋势有极显著的正相关性,相关系数为0.46(图8-d);干季的NEP 趋势与温度趋势没有显著的相关性(图8-b),与降雨趋势有极显著的负相关性,相关系数为−0.06(图8-e);湿季的NEP 趋势与温度趋势、降雨趋势均有极显著的正相关关系,相关系数分别为0.11 和0.30(图8-c,f)。如图9 所示,多元线性回归结果表明,在干季600 m 以下的森林生态系统中,降水趋势对NEP 趋势有显著的负贡献,贡献率为−53%。在湿季600 m 以下的森林生态系统中,温度趋势对NEP 趋势有显著的正贡献,贡献率为90%。
图 7 海拔600 m 以上森林NEP 趋势与温度、降雨趋势的相关性散点图(a)、(b)、(c)分别表示年际、干季、湿季与温度趋势的相关性;(d)、(e)、(f)分别表示年际、干季、湿季与降雨趋势的相关性。
图 8 海拔600 m 以下森林NEP 趋势与温度、降雨趋势的相关性散点图(a)、(b)、(c)分别表示年际、干季、湿季与温度趋势的相关性;(d)、(e)、(f)分别表示年际、干季、湿季与降雨趋势的相关性。
图 9 基于多元线性回归的不同时期不同海拔下气候因子的贡献度分析图
红色五角星表示方程具有显著性,黑色星号表明回归系数或残差通过显著性P=0.05 水平检验。
3.1 海南岛森林稳定的碳汇功能海南岛森林NEP 表现出了季节性变化,但是与亚马逊地区不同的是,亚马逊地区在季节性变换时会存在碳动态的转变,如碳源和碳汇的相互转换[8,11−12,34−35],而海南岛热带森林不论是从年际角度出发,还是从干季、湿季角度出发,均呈现稳定碳汇状态。所不同的是,海南岛森林NEP 在湿季时普遍要比干季高。地理位置差异所带来的气候类型差别可能是海南岛热带森林和亚马逊热带雨林呈现不同碳动态格局的原因。亚马逊地区属热带雨林气候,长年高温多雨,对干旱的敏感性更高[36];而海南岛则属热带海岛季风气候,受季风影响分明显的干季和湿季。相对固定的季节模式,使得海南岛森林固碳动态虽有季节性变化,但总体趋势保持稳定。更为重要的一点是,不只是气候变化的影响,农业扩张所带来的森林砍伐,火灾和干旱之间的相互作用也是造成亚马逊森林碳损失的原因[36−37]。相比而言,海南岛的森林保存的较为完整。此外,本研究采用基于多期土地利用遥感数据提取稳定林地的方法,一定程度上规避了森林变化或者损失而带来的净生态系统生产力的损失。
3.2 不同海拔不同季节的NEP 与气候因素的相关性分析从趋势的相关性上来看,相比气温,降雨与NEP 有着更为显著的相关性。但是在不同海拔下,降雨与NEP 的相关性表现出了相反的模式。在海拔600 m 以下,年际、湿季降雨趋势分别与年际、湿季NEP 趋势呈现出显著的正相关,但是在海拔600 m 以上,又呈现出截然相反的显著负相关性结果。群落结构和生态系统过程往往随海拔梯度变化[24],一般而言对海拔的响应都是由于温度变化所驱动,但也不绝对,降雨等因素也会随海拔梯度变化[24,38]。海拔梯度下,海南岛降雨趋势的变化不是十分明显,所以这种结果可能是因为群落结构发生变化。海南岛海拔600 m 以上保存较为原始的热带山地雨林,在年际和湿季雨量较为充沛的阶段,随着降雨趋势的增加,附生植物旺盛生长会导致NEP 下降。而在低海拔区域,由于人为干扰,林型更为复杂,人工林的固碳能力相较高海拔原始林更为强大。
3.3 降水与温度的主导性分析仅从气候因素考虑,在旱季时,降雨对低海拔森林NEP 变化有显著的负向作用,而在湿季时,温度对低海拔森林NEP 变化有积极作用。气温和降雨是影响生态系统生产力的主要的因素[39],但是随着全球变化加速,陆地生态系统的碳交换对于气候变化的响应和反馈仍旧存在不确定性[40−41]。尤其在热带森林生态系统,降雨量在年内的波动导致季节性干旱的产生,干旱事件对热带森林生态系统树木的生长以及生态系统的功能有着严重的影响[42]。低海拔区域更多为人工林,以橡胶林为例,干旱会导致其落叶甚至死亡,光合作用受到抑制[43]。另一方面,旱季的降雨会增加土壤水分[44],进一步影响了土壤呼吸[45]。雨季时,温度呈现出显著的正向作用。降水的充沛使得水分并不是NEP 的主要限制因子,而此时的温度却有着积极的调节能力,与NEP 变化呈现出一致性。研究表明,温度的升高是利于生态系统呼吸[16]、生态系统总产量[46]以及净碳吸收[47]。可能在雨季时,由于不再受到水分条件的限制,海南岛森林植被得到较大恢复,处于生长阶段的旺盛期,此时温度对于整个生态系统的光合作用影响远远大于呼吸作用。
近20 年海南岛森林NEP 的年际变化波动较大,存在明显的季节性。湿季的NEP 高于干季的NEP。空间上,NEP 呈显著增长趋势主要分布在海南岛东北部,显著下降趋势主要分布在海南岛中部。不同季节NEP 与气候因子的关系有所差异,在干季低海拔区域,降雨对NEP 有显著的负向贡献,贡献度为−53%(P<0.05);在湿季低海拔区域,温度对NEP 有显著的正向贡献,贡献度为90%(P<0.05)。