劉應帥,余 瑞,鄭彬彬,劉嘉慧,宋 奇,陳榮昊,嚴 哲
(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)。