廣東西樵山國家森林公園森林碳儲量空間分布研究

張凌宇,趙慶,吳曉君,許東先,羅皓,謝進金

(廣東省森林培育與保護利用重點實驗室/廣東省林業科學研究院,廣州 510520)

0 引言

全球氣候變化已經成為21世紀人類社會面臨的最大挑戰,由于社會經濟發展、人口增加以及能源消耗等,二氧化碳排放量逐年增加是當代氣候變化的根本原因,并對碳循環造成重大影響[1]。基于此,2020年9月,中國政府向國際社會宣布,將提高國家自主貢獻力度,采取更有力的政策和措施,努力爭取2030年前二氧化碳排放達到峰值,2060年前實現碳中和。這是中國首次向全球明確實現碳中和的時間點,也是迄今為止世界各國中做出的最大減少全球變暖預期的氣候承諾[2]。

森林在全球碳循環過程占據著重要位置,特別是在增加碳匯、維持碳平衡以及氣候調節等方面發揮著重要作用[3]。作為發展碳匯林業、減排增匯的基礎,精準獲取區域尺度森林碳儲量已成為現如今全球氣候變化與碳循環研究的重要內容,森林碳儲量空間分布格局的研究也可以從布局和決策上指導森林經營方案的制定[4-5]。

目前,有關森林碳匯估算方法主要分為3種,分別為樣地清查法、利用變分原理建立多物理量模型的渦度相關法以及模型模擬法[6-9],然而無論是以換算因子法為基礎的樣地調查法,還是以地統計學和遙感估算技術為基礎的碳儲量估測模型進行模擬,都存在一定的局限性,實際上都沒有考慮森林碳儲量數據收集過程中的空間異質性和空間自相關性[10]。由于林業數據收集過程中,處于不同位置上的林分存在差異,其森林碳儲量數據普遍存在空間異質性,這類空間非平穩性數據,容易影響基于傳統統計方法所得到結果的準確性[11],而空間自相關性則隨著空間距離的變化而變化,在大尺度研究森林碳儲量的空間分布時,往往忽略這種相關性的影響,使研究結果的精度大幅度降低,一般采用Moran'sI表示研究區整體和局域的空間自相關性[12]。

近年來大量研究表明,地理加權回歸模型(Geographically weighted regression model, GWR)在解決空間異質性和空間自相關性問題上具有明顯的優勢。GWR模型最早在國外應用在環境科學和社會科學等研究領域[13-14]。引入國內后,GWR模型在林業領域做出諸多應用。例如王海賓等[15]以浙江省內一景Landsat 8影像覆蓋范圍內的喬木林為研究對象,采用GWR模型和協同克里格插值法分別構建了喬木林地上碳密度估算模型,結果顯示GWR模型構算的模型精度要明顯高于協同克里格插值法的模型精度,同時GWR模型較好地保留了估算變量的空間異質性;李明澤等[16]通過將遙感因子和林分因子相結合的方法建立了基于GWR模型的大興安嶺森林立地質量估算模型,結果顯示加入空間樣地信息的GWR模型有效地降低了空間自相關性。除了對比分析傳統GWR模型和其他模型在估算精度上的差異以外,也有學者采用GWR模型和其他模型相結合的方法對研究目標進行模擬和預測,如Zhen等[17]以黑龍江涼水自然保護區各小班風倒木數據為研究對象,分別利用地理加權邏輯回歸(Geographic weighted logistic regression,GWLR)、地理加權泊松模型(GWPR)和GWR模型對風倒木發生概率、數量和蓄積進行了研究,結果顯示局域模型在模型擬合優良性上要優于相應的全局模型,并形成了良好的模型殘差分布;戚玉嬌[18]采用最小二乘模型(Ordinary Least Squares, OLS)、空間誤差模型(Spatial error model, SEM)、GWR模型3種空間模型和普通最小二乘模型對森林的碳儲量進行了估測,研究結果顯示,GWR模型的擬合效果要明顯優于其他3種模型,相比于OLS模型和SEM模型,GWR模型在最大程度上消除了模型殘差的空間自相關性。

林地“一張圖”是以“互聯網+”、云計算、大數據為戰略背景,在一定程度上實現了森林資源管理的精細化與科學化,本研究以2020年廣東西樵山國家森林公園森林資源管理“一張圖”數據為基礎,運用OLS模型、SEM模型和GWR模型分別對該地區森林碳儲量的空間分布進行模擬,分析全局模型和不同參數配置下的空間模型在解決空間自相關性上的優劣,為解決空間非平穩性問題提供了一種思路和方法,此外,通過分析研究區不同區域的空間異質性,也可以為該地區建立基于區域尺度的城郊森林公園景觀健康評價指標體系以及森林景觀異質性研究提供理論支持[19-21]。

1 研究地區與研究方法

1.1 研究地區概況

研究區廣東西樵山國家森林公園位于佛山市南海區的西南部,屬南亞熱帶季風氣候區,土壤主要為赤紅壤,植物種類豐富,其中主要的喬灌草植物有380多種,代表性植物群落為南亞熱帶季風常綠闊葉林,主要喬木由殼斗科(Fagaceae )、樟科(Lauraceae)、山茶科(Theaceae )、梧桐科(Sterculiaceae)和大戟科( Euphorbiaceae )等種類的樹種組成;灌木主要有青江藤(Celastrushindsii)、五指毛桃(Ficushirta)、鬼燈籠(Clerodendrumfortunatum)和車輪梅(Rhaphiolepisindica)等;草本植物主要有禾本科和蕨類植物火炭母(Polygonumchinense)、露籽草(Ottochloanodosa)、芒萁(Dicranopterispedata)和三叉蕨(Tectariasubtriphylla)等。

1.2 森林碳儲量數據來源

森林碳儲量數據來源于2020年廣東西樵山國家森林公園森林資源管理“一張圖”數據,本研究中的森林碳儲量指森林植被碳儲量,包括喬木碳儲量、林下植被碳儲量、灌木碳儲量和草本碳儲量。單位面積森林的碳儲量,即森林碳密度,是單位面積森林生物量與碳轉換系數的乘積,國際常用的碳轉換系數為0.5,單位面積生物量的估算采用廣東省森林資源監測中心提供的生物量擴展因子方程,利用連續清查數據中的測樹因子,分樹種推算出樣地實測生物量,而對于一些沒有對應生物量方程的樹種,選擇相近樹種的生物量方程代替計算[22-23]。

1.3 模型變量的選擇

本研究采用逐步選擇法(引入和剔除的顯著性標準ɑ =0.05)對模型變量進行選取,通過引入和剔除交替進行自變量檢查,直到無統計學意義的新變量可以引入也無自變量可以剔除時為止。結合平均胸徑、平均樹高、林分蓄積量、林分密度、郁閉度、坡度、坡向、海拔、下木蓋度、土壤厚度和枯枝落葉厚度11個初始變量進行選擇,在對各參數進行相關性檢驗的基礎上,得到林分平均胸徑、林分平均高、海拔和坡度4個影響森林碳儲量分布的變量,本研究中所有的獨立變量都進行了標準化處理。

1.4 最小二乘模型(OLS)

最小二乘模型是利用多組觀察值求得p個自變量X與因變量Y之間殘差平方和最小的一種擬合方法,其表達式如下。

(1)

式中:β為估算的未知固定效應的模型系數向量;ε為模型殘差,服從N(0,σ2)。

1.5 空間誤差模型(SEM)

空間誤差模型表示在空間上鄰域的Y對該點沒有直接影響,空間相關性來自其他因素影響,其表達式如下。

Y=Xβ+λWε+ξ。

(2)

式中:Wε為空間誤差項;ξ是真正的模型誤差項,服從N(0,σ2I);λ為空間自相關參數。

1.6 地理加權回歸模型(GWR)

地理加權回歸模型是將距離權重作為一種解決空間異質性的方法納入到模型中,從而在很大程度上減少甚至消除由于空間異質性的作用導致回歸模型參數估計不準確的影響,其表達式如下。

(3)

式中:β0(ui,vi)為回歸模型的截距;(ui,vi)為第i取樣點的地理坐標;βk(ui,vi)為第i取樣點上第k個參數的估計值,屬于一個地理函數坐標,ε為模型殘差,服從N(0,σ2)。

1.7 模型殘差的空間自相關性

本研究分別采用全局和局域的MoranI來評價模型殘差的空間分布,全局MoranI計算公式如下。

(4)

(5)

在本研究中,局域MoranI用來檢測每一個獨立樣地點模型殘差的局部聚集狀況,當局域MoranI<0時,說明在該樣地點周圍聚集著不同的模型殘差,此時模型殘差的空間分布最好,當該點的局域MoranI>0時,表示在該樣地點周圍聚集著相似的模型殘差[24]。

1.8 模型擬合優度比較

1.9 數據處理

本研究分別采用SAS9.3、GeoDa、GWR4.0軟件建立OLS模型、SEM模型和GWR模型,采用Excel下宏文件ROOKCASE對全局和局域MoranI值進行計算,最后基于ArcGIS 軟件平臺,采用反距離權重法(Inverse Distance Weight,IDW)對GWR模型下的空間分布圖進行繪制。

2 結果與分析

2.1 模型擬合

由表1可知,對于OLS模型而言,所有變量的方差膨脹因子(Variance Inflation Factor,VIF)均小于5,表明模型各自變量間不存在多重共線性,并且在ɑ =0.05水平上所有變量的參數估計值均表現出統計顯著性,可以對因變量進行很好的解釋。在OLS模型各變量中,林分平均胸徑、林分平均高和海拔的參數估計值為正,說明這3個變量與碳儲量呈正相關,坡度的參數估計值為負,說明其與碳儲量呈負相關,所有變量對森林碳儲量影響由大到小順序為:林分平均高、林分平均胸徑、坡度、海拔。

表1 最小二乘模型各變量的參數估計值、標準誤差和顯著性檢驗Tab.1 Parameter estimates, standard error and significance test of variables in OLS model

由表2可知,SEM模型各變量的參數估計值同樣通過了顯著性檢驗,且與森林碳儲量在正負相關性上與OLS模型表現相同,但是在對森林碳儲量的影響大小上與OLS模型存在差異,由大到小順序為:林分平均高、林分平均胸徑、海拔、坡度。

表2 空間誤差模型各變量參數估計值、標準誤差和顯著性檢驗Tab.2 Parameter estimates, standard error and significance test of variables in SEM model

相對于OLS模型和SEM模型得到的一組參數,GWR模型給出了5組模型參數估計值,分別為最小值、25%分位值(Q1)、中值、75%分位值(Q3)和最大值,見表3。

表3 地理加權回歸模型各變量參數估計值的描述性統計Tab.3 Descriptive statistics of parameter estimates of variables in GWR model

由表3可知,GWR模型產生了較大范圍的參數估計值,且GWR模型各參數估計值的范圍均包括了OLS模型和SEM模型的參數估計值,GWR模型在中值位置下的各參數估計值與OLS模型和SEM模型的參數估計值在數值大小上表現出相近性,且與森林碳儲量在正負相關性上與二者表現一致。與OLS模型和SEM模型不同的是,GWR模型隨著位置的變化,各參數估計值對森林碳儲量的影響大小也時刻發生著變化,如GWR模型在最小值位置下對森林碳儲量的影響由大到小順序為:海拔、坡度、林分平均胸徑、林分平均高,而在最大值位置下對森林碳儲量的影響由大到小順序為:林分平均高、林分平均胸徑、海拔、坡度。

為了更清楚地表現GWR模型各參數估計值的空間分布情況,圖1給出了GWR模型各變量對應參數的空間分布圖,林分平均胸徑參數估計值除了在研究區南部極小范圍內為負值外,在其他地區參數估計值均為正值,說明在這些區域內,森林碳儲量隨著林分平均胸徑的增大而增大,此外,從研究區中部向東西兩側逐漸延伸,林分平均胸徑對森林碳儲量的影響逐漸增大(圖1(a));在整個研究區范圍內,林分平均高參數估計值均為正值,且自西南向東北方向逐漸增大(圖1(b));在研究區中部較大范圍內海拔參數估計值為正值,說明在該區域范圍內,隨著海拔的逐漸增高,森林碳儲量逐漸增大(圖1(c));坡度方面,從總體上看,坡度參數估計值為正的區域主要集中在研究區的東北部和北部小范圍地區,參數估計值為負值的區域約占研究區總面積的85%以上,說明在本研究區內坡度對于森林碳儲量增加起到抑制作用的區域要明顯大于起到促進作用的區域(圖1(d))。

圖1 地理加權回歸模型各變量對應參數的空間分布Fig.1 Distributions of variable parameters in GWR model

2.2 模型擬合效果評價

表4 3種模型評價的統計量Tab.4 Statistics for the evaluation of the three models

2.3 模型殘差的空間自相關性

圖2給出了3種模型在不同空間距離上模型殘差全局Moran'sI的變化情況,由圖3可知,從總體上看SEM模型和OLS模型在空間距離達到3 900 m以后逐漸趨于穩定,GWR模型的全局Moran'sI基本上都在0上下浮動,變化幅度較小,且各空間距離下的全局Moran'sI值均小于其他2種模型,模型的穩定性最好,SEM模型變化幅度略小于OLS模型,OLS模型的穩定性最差。

圖2 3種模型殘差的全局Moran's IFig.2 Global Moran's I correlogram of three models

圖3給出了3種模型殘差局域Moran'sI的空間分布圖,由圖3可知,在OLS模型和SEM模型中,研究區中北部、東北部及南部的部分區域產生了較大范圍且取值為正的局域Moran'sI值,說明OLS模型和SEM模型的殘差在這些區域傾向于呈現出相似聚集狀態。而在GWR模型中,局域Moran'sI取值為正的區域明顯減少,模型殘差相似聚集狀態大幅度下降,此時,GWR模型殘差的空間自相關性明顯低于其他2種模型,產生了最好的局域化空間分布效果。

2.4 森林碳儲量的空間分布

圖4給出了3種模型對研究區森林碳儲量模擬的空間分布圖,由圖4可知,3種模型在總體范圍上對于森林碳儲量空間分布的模擬趨勢基本相同,即研究區的中部森林碳儲量較多,在該地區3種模型的碳儲量含量均在35 t/hm2以上,其他區域碳儲量相對較少,但是大部分含量也在30 t/hm2以上,其中OLS模型和SEM模型在模擬森林碳儲量的空間分布上呈現出明顯的相似性,GWR模型模擬的森林碳儲量在35 t/hm2以上的面積要高于OLS模型和SEM模型對應的模擬面積,此外,3種模型的擬合偏差存在差異,其中OLS模型的擬合偏差為2.14 t/hm2,SEM模型的擬合偏差為2.11 t/hm2,GWR模型的擬合偏差為1.36 t/hm2。

圖4 3種模型森林碳儲量的空間分布Fig.4 Spatial distribution of forest carbon storage of three models

3 討論

本研究采用了OLS模型、SEM模型和GWR模型對西樵山森林碳儲量的空間分布情況進行了模擬,對影響森林碳儲量分布的地形因子和林分因子進行了分析,并對3種模型的擬合效果和空間自相關性進行了分析,結果表明,3種模型下各變量參數在對研究區森林碳儲量的影響上基本一致,其中林分平均胸徑對森林碳儲量的影響最大,兩者存在正相關關系,這也符合森林生長規律,隨著林分平均胸徑的增加,林分整體固碳能力增強,森林碳儲量逐漸增大,這與陳科屹等[5]、劉正顯[25]的研究結果一致。

GWR模型在各參數對森林碳儲量靈敏性的描述上也好于其他模型,如處于不同位置下各參數對森林碳儲量影響的大小順序存在差異,這是由于GWR模型在運算過程中,隨著模型所處位置的變換,在每個位置上都會形成一個局域空間模型,該模型在運行過程中充分考慮了地形、林分因子變化對森林碳儲量的影響,通過獨立運算降低或消除空間自相關性的影響,從而形成了不同位置下各參數對森林碳儲量影響大小和方向的不同,而OLS模型和SEM模型在運行過程中只形成了一個模型[12, 17-18],因此在各參數對森林碳儲量影響的靈敏性描述上存在一定的局限性。在模型擬合方面,GWR模型擬合精度要優于OLS模型和SEM模型,原因在于GWR模型充分考慮了變量的空間異質性和空間相關性對擬合結果產生的影響,考慮了森林碳儲量和相關變量在研究范圍內的特定關系[3]。在模型殘差空間自相關性方面,分析了3種模型殘差的全局空間自相關性及空間分布,結果發現,GWR模型的全局自相關性和全局Moran'sI的起伏程度都明顯小于OLS模型和SEM模型,其局域Moran'sI形成了不同觀測值聚集的較好的空間分布效果,表明GWR模型的穩定性較好,這與其他學者的研究結果一致[17, 25-27]。

在本研究中,無論從全局Moran'sI在各步長的起伏跨度上,還是局域Moran'sI的可視化空間分布上,OLS模型和SEM模型都表現出極大的相似性,這是由于研究區面積較小,并不能很好地反映SEM模型相對于OLS模型的優越性;在模型擬合方面,SEM模型的擬合精度略高于OLS模型,這也說明SEM模型雖然在一定程度上可以提高模型的精度,事實上SEM模型在運算過程中考慮了變量對模型誤差產生的空間依賴性的影響,但模型誤差項有時并不能完全反映出實際的空間自相關和空間異質性情況[10],而GWR模型在研究區內則表現出較好的優越性,有效地消除空間自相關性和解決空間異質性問題,從而不受研究區面積大小的限制,同時,GWR模型更容易消除尺度效應的影響,近年來研究顯示,尺度效應在森林更新[28]、景觀健康評價[21]、物種豐富度[29]、森林景觀類型空間關聯性[30]以及優勢樹種混交林的空間格局[31]等方面都起著重要的作用。

森林植被、土壤屬性、立地條件、氣候條件、凋落物和根系輸入等都影響著森林碳儲量的分布,本研究在借助逐步回歸方法的基礎上,只選取了對森林碳儲量分布影響最大的林分因子和地形因子進行分析,由于研究區面積限制,并未考慮氣象、土壤和環境等因子對森林碳儲量的影響。此外,研究區以國家森林公園為主體,兼具游憩和自然教育等功能,人類活動可能對森林碳儲量分布存在一定的影響,如研究區中軸線多以景區為主,中部地區碳儲量高于其他區域,其樹干形態包括群落結構等都經過大規模人工干預;其他區域雖然人工林占比也較大,但受人工干預相對較小,因此人類活動是否影響碳儲量的分布還有待進行后續研究。

4 結論

2)林分因子對西樵山森林碳儲量的影響最大,從總體上看,林分平均胸徑和林分平均高與該地區森林碳儲量呈正相關關系。

3)SEM模型在研究區面積較小的情況下對空間自相關性的反映并不靈敏,而GWR模型則不受研究區大小的限制。

4)相對于OLS模型和SEM模型,GWR模型的全局Moran'sI基本上都在0上下浮動,起伏程度較小,形成更好的空間分布效果,有效地降低了空間自相關性。

5)研究區西樵山森林碳儲量整體上表現為中部區域高于周邊區域的特點,其中中部區域的碳儲量分布范圍在35～40 t/hm2,周邊區域碳儲量分布范圍在30～35 t/hm2。