劉關琴,施凱澤,萬玉潔,胥 輝
(1. 西南林業(yè)大學 西南地區(qū)生物多樣性保育國家林業(yè)局重點實驗室,云南 昆明 650224;2. 云南省林業(yè)調查規(guī)劃院,云南 昆明 650051)
森林是陸地生態(tài)系統(tǒng)中最大的有機碳儲庫,對全球碳循環(huán)至關重要[1],評價森林生態(tài)系統(tǒng)碳匯能力的其中一個重要指標是森林生物量的分布和變化[2],森林碳匯估算的關鍵就在于如何準確估算生物量。傳統(tǒng)的野外測量不能提供大規(guī)模連續(xù)觀測數(shù)據(jù),且費時費力,遙感模型估測法可以實時、動態(tài)、無破壞地對森林地上生物量進行估算,尤其適用于大區(qū)域的生物量估算[3]。植被在紅外波段和可見光的光譜特征具有明顯區(qū)別,可利用植被的光譜特征來監(jiān)測植被生長狀況和估測森林生物量,其中,Landsat影像已成為眾多學者使用最普遍的遙感數(shù)據(jù)源[4],但仍具有不確定性,不確定性的來源主要有生物量模型預估表現(xiàn)、野外樣地數(shù)據(jù)獲取、遙感數(shù)據(jù)源及其空間分辨率、森林的空間異質性、森林生境信息因子等,其中光飽和點的不確定性尤為突出[5]。
近年來,諸多學者對如何解決森林生物量數(shù)據(jù)飽和帶來的問題進行了研究,歐光龍等[6]通過構建不同模型對森林生物量遙感進行估測以探索減少飽和現(xiàn)象造成的估測精度問題。周律等[7]基于遙感變量構建空間回歸模型,避免了模型自相關性,提高了思茅松林地上生物量估測精度。許振宇[8]等結合主被動遙感數(shù)據(jù),采用4種模型估測分析桂東縣森林生物量,為亞熱帶森林生物量建模分析提供了新思路。
金沙江流域位于長江上游,在云南分布廣泛,生態(tài)系統(tǒng)十分脆弱[9]。按森林植被類型劃分可將此流域分為3個區(qū),東部中山山原峽谷區(qū)、中部中山峽谷和滇中高原區(qū)、西部高山峽谷區(qū)[10]。中部中山峽谷和滇中高原區(qū)的主要植被類型為云南松、高山松、針葉林和旱冬瓜、櫟類等形成的混交林,此區(qū)域缺少相關優(yōu)勢樹種光學遙感估測的飽和點數(shù)據(jù)。目前,其他地區(qū)已有對這些優(yōu)勢樹種光學遙感估測飽和點的探索結果,盧騰飛等[11]利用三次模型擬合曲靖市云南松生物量飽和值為167 t/hm2。趙盼盼[12]估算浙江省中西部針葉林的生物量飽和值得到馬尾松與杉木遙感生物量估測飽和值分別為159 t/hm2和143 t/hm2。吳勇等[13]通過皮爾遜相關系數(shù)比較篩選出Landsat 8 OLI Band 4波段變量與AGB進行曲線擬合,求得高山松林與云冷杉林遙感生物量估測飽和值分別為149.09 t/hm2和162.3 t/hm2。
基于此,將以云南省金沙江流域中段的5類主要樹種為研究對象,基于Landsat 8 OLI遙感影像,結合地面二調小班數(shù)據(jù),構建生物量估測模型,探索光學遙感數(shù)據(jù)源的飽和點閾值確定方法,并反演研究區(qū)5類主要樹種的地上生物量。實現(xiàn)對研究區(qū)生物量估算精度的提升,減小因存在數(shù)據(jù)飽和問題所造成的不確定性。為進一步提高大尺度流域生物量光學遙感估測精度提供理論依據(jù),為今后金沙江流域森林碳匯計算提供理論基礎,也為其他流域的森林生物量光學遙感估測及飽和點分析研究提供參考。
金沙江流域發(fā)源于青藏高原,流域全長1 560 km,介于25°~36°N和90°~103°E范圍,流域內地形極為復雜,眾多高山峽谷相間并列,地勢北高南低,落差約5 100 m。各區(qū)段自然環(huán)境差異較大,有明顯的垂直氣候特征[14]。研究區(qū)域位于大理、麗江、楚雄3個州市所在的金沙江流域,本文中統(tǒng)一稱為金沙江流域中段(圖1)。
2.1.1小班數(shù)據(jù)
1)小班單位面積地上生物量計算
選取研究區(qū)2016年森林資源二類調查小班數(shù)據(jù),包括金沙江流域流經(jīng)的楚雄、大理、麗江3個州市。根據(jù)實際情況,本研究選取云南松林、常綠櫟類林、高山松林、其他闊葉林、其他針葉林進行分析建模(表1)。參考胥輝等[15]的蓄積量—生物量轉換模型計算各小班單位面積森林生物量[16]。計算公式為:
B=V×SVD×BEF
(1)
式中:B為小班單位面積地上生物量(t/hm2);V為小班單位面積蓄積量(m3/hm2);SVD為基本木材密度(t/m3);BEF為生物量轉換因子(無量綱)。
2)樣本小班的確定
優(yōu)勢樹種或樹種組在研究區(qū)分布如圖1B所示,將小班中的異常數(shù)據(jù)剔除,篩選出平均胸徑大于等于5 cm的小班;利用 ArcGIS“子集要素”工具,等比例保留小班數(shù)據(jù);再利用三倍標準差法刪除異常離群值,最終選出4 144個小班樣地,其中分別取70%小班進行建模,30%小班作為檢驗樣本,具體信息如表2所示。
(a)金沙江流域中段Landsat 8 OLI合成遙感影像 (b)金沙江流域中段各樹種或樹種組分布
表1 生物量轉換因子法計算模型Tab.1 Calculation model of biomass conversion factor method
2.1.2遙感影像數(shù)據(jù)
在地理空間數(shù)據(jù)云網(wǎng)站下載與二類調查數(shù)據(jù)同時期的流域區(qū)域的Landsat 8 OLI衛(wèi)星影像數(shù)據(jù)(表3),利用ENVI 5.3進行輻射定標、FLAASH大氣校正、地形校正預處理,進而通過影像融合、裁剪及鑲嵌獲得研究區(qū)預處理后的遙感影像數(shù)據(jù)。
2.1.3遙感因子提取與篩選
以云南松林、常綠櫟類林、高山松林、其他闊葉林和其他針葉林的小班面狀數(shù)據(jù)為單位,將小班面狀矢量數(shù)據(jù)和134個遙感特征因子數(shù)據(jù)層相疊加,統(tǒng)計小班內各遙感因子的平均值[18]。本研究提取的4類遙感因子分別如下:
表2 小班樣地的選擇Tab.2 Selection of subcompartment plot
表3 研究區(qū)Landsat 8 OLI影像基本信息Tab.3 Basic information of Landsat 8 OLI images in study area
1)原始單波段[17](7):B1、B2、B3、B4、B5、B6、B7;
2)植被指數(shù)[18](9):ARVI、B、DVI、G、MSAVI、NDVI、PVI、RVI、SAVI;
3)K-T和K-L變換因子[18](6):KT1、KT2、KT3;PCA1、PCA2、PCA3;
4)紋理特征[18](3*3、7*7,112):Mean、Variance、Homogeneity、Contrast、Dissimilarity、Entropy、Second Moment、Correlation。
研究采用皮爾遜相關系數(shù)分析,SPSS 軟件對生物量與各遙感變量進行相關性分析,選擇相關顯著性(P≤0.05)的遙感因子作為建模備選變量。
根據(jù)地上生物量與各波段光譜反射率相關性,選擇生物量與相關性最為顯著的波段反射率做散點圖,發(fā)現(xiàn)波段反射率會隨生物量的增大逐漸減小,當生物量達到一定水平時波段反射率不再發(fā)生變化,該臨界值即為光學遙感估測生物量時的飽和點。為進一步量化云南松林、常綠櫟類林、高山松林、其他闊葉林和其他針葉林的生物量飽和值,通過擬合生物量與波段反射率的函數(shù)關系,計算函數(shù)所對應的拐點值確定生物量遙感估測的光飽和點值。
2.4.1線性回歸模型構建
運用IBM SPSS Statistics 22對遙感因子進行相關性分析,篩選出顯著性較好的因子,自變量為遙感因子,因變量為單位生物量,以IBM SPSS Statistics 22為平臺,建立線性逐步回歸模型,其公式為:
Y=β0+β1x1+β2x2+…+βnxn+ε
(2)
式中:Y為生物量;β0為常數(shù)項;β1,β2…βn為模型系數(shù);x1,x2…xn為相關遙感因子;ε為模型滿足隨機正態(tài)分布下的隨機殘差;n為自變量個數(shù)。
運用SPSS進行多元線性回歸,輸入觀測值即單位生物量、自變量即相關性較顯著的遙感因子值,輸出結果包括非標準化系數(shù)、標準誤差、標準系數(shù)、t、Sig、VIF。在變量間的多重共線性問題極其常見,且對模型造成嚴重影響,用VIF函數(shù)對自變量進行共線性檢驗。
2.4.2KNN模型構建
KNN算法能用于森林參數(shù)的估計和各種空間數(shù)據(jù)的因變量估測,通過計算K個樣本點特征空間中最鄰近參考像元的屬性并計算其反距離權重平均值得到未知像元點的屬性。一般情況下,訓練樣本點的影像隨距離的增加逐漸減小。本文主要采用歐氏距離,計算公式為:
(3)
式中:ρ為點(x2,y2)與點(x1,y1)之間的歐氏距離;|X|為點(x2,y2)到原點的歐氏距離。
2.4.3模型獨立性檢驗
選取決定系數(shù)(R2)、平均絕對誤差(MAE)、平均相對誤差(MRE)、均方根誤差(RMSE)和相對均方根誤差(rRMSE)5個指標作為模型的獨立性檢驗指標。公式為:
決定系數(shù)(R2):
(4)
平均絕對誤差(MAE):
(5)
平均相對誤差(MRE):
(6)
均方根誤差(RMSE):
(7)
相對均方根誤差(rRMSE):
(8)
2.4.4“刀切法”殘差分析
殘差是觀測值與回歸值的差,其通過模型的擬合效果對模型預測結果產(chǎn)生影響。通過“刀切法”殘差分析對研究區(qū)5類優(yōu)勢樹種或樹種組和總體樣本分別進行劃分區(qū)間,計算分析區(qū)間模型殘差情況,確定模型在遙感估測中是否存在誤差。采用3個指標作為“刀切法”殘差結果評價指標,分別是平均殘差(ME)和平均相對殘差(MRE),計算公式為:
(9)
式中:n為樣本總數(shù);yi為實測值;y^i為估測值。
圖2表示5類優(yōu)勢樹種或樹種組的遙感因子散點圖和飽和點擬合曲線。由曲線趨勢可以看出,5類優(yōu)勢樹種或樹種組的單波段遙感因子值隨其單位生物量的增加而逐漸趨于一個固定值,這個值所對應的單位生物量即為各樹種的飽和點。利用飽和點曲線的二項式方程計算得出云南松林、常綠櫟類林、高山松林、其他闊葉林、其他針葉林的飽和點分別為 98.870 t/hm2、101.220 t/hm2、 127.938 t/hm2、88.432 t/hm2、160.201 t/hm2;R2分別為0.471、0.325、0.232、0.336、0.312。
3.2.1生物量與遙感因子相關性分析
篩選出5類優(yōu)勢樹種或樹種組,并將這5類優(yōu)勢樹種或樹種組的遙感因子的均值統(tǒng)計量與生物量進行Pearson’s相關性分析,將相關性高的遙感因子進行多重共線性診斷,選取方差膨脹因子VIF<10的遙感特征變量,最終得到的結果見表4。
3.2.2模型擬合
1)線性回歸模型構建
5類樹種或樹種組參與逐步線性回歸模型構建的遙感因子之間的線性關系見表5。
從表5可以看出,各樹種或樹種組與參與模型構建的遙感因子之間的線性關系,眾多遙感因子中紋理特征被納入建模次數(shù)較多,存在較好的相關性,其中COR7_1、SM7_4分別被不同的3個樹種或樹種組納入3次,而K-T和K-L變換因子僅被其他闊葉林和其他針葉林納入一次建模,表明這類遙感因子與這5類樹種或樹種組的生物量相關性不強。不同樹種或樹種組的R2也各不相同,其他針葉林的R2最高,為0.356;其次是其他闊葉林的R2,其值為0.344,且RMSE最低,為5.190,表明此模型對于這兩類樹種的地上生物量具有最好的模型擬合效果;常綠闊葉林、高山松林的R2均在0.2以上,且RSME值均相對較低,表明模型對這幾類森林的生物量具有較好的模型預估精度;云南松林的R2只有0.112,模型擬合精度最低,表明多元逐步線性回歸模型對于計算該地區(qū)的云南松林生物量的適用性不強。
圖2 5類優(yōu)勢樹種或樹種組飽和點擬合曲線Fig.2 Fitting curve of saturation point for five dominant tree species or species groups
表4 5類優(yōu)勢樹種或樹種組的相關性分析及共線性診斷Tab.4 Correlation analysis and collinearity diagnosis results of five dominant tree species or species groups
表5 5類樹種或樹種組逐步線性回歸模型構建Tab.5 Construction of stepwise linear regression model for five dominant tree species or species groups
2)KNN模型構建
KNN模型以線性回歸篩選的遙感因子為自變量,以單位生物量為因變量,對各樹種進行建模,就決定系數(shù)(R2)來講,有4個模型達到0.3以上,其他針葉林最高為0.416;從均方根誤差(RMSE)來講,其他闊葉林模型最低為5.110,總體來看KNN模型的擬合效果更好一些(表6)。
表6 KNN模型精度評價結果Tab.6 KNN model accuracy evaluation results
3.2.3模型檢驗
總體來看,估測模型中總體平均相對誤差、平均絕對誤差和相對均方根誤差基本上是KNN算法表現(xiàn)的最低,而決定系數(shù)KNN算法表現(xiàn)的最高(表7)。因此,選擇KNN作為云南省金沙江中段遙感估測的最終模型。
表7 線性逐步回歸模型和KNN模型獨立性檢驗Tab.7 Independence test results of stepwise linear regression model and KNN model
3.2.4“刀切法”殘差分析
通過分段殘差分析,整體上KNN模型的預估精度較高,誤差相對逐步回歸模型較小。從生物量分段看,兩個模型均存在生物量低值高估和高值低估情況,其他闊葉林只有<50 t/hm2生物量段,且均為高值低估;在低生物量段(<50 t/hm2)時,除其他闊葉林和逐步回歸模型估測的其他針葉林的MRE為正值外,其余均為負值。在(>50 t/hm2)生物量段時,均為正值,且ME和MRE隨生物量段的增加呈增大趨勢,直到(>150 t/hm2)生物量段時,常綠櫟類林和高山松林的ME均大于90 t/hm2,且MRE值均大于1,出現(xiàn)顯著的高值低估情況,說明在高生物量段兩模型估測能力均顯不足(表8)。
遙感影像的各特征變量均是對森林生物量特征及其變化的響應,但普遍存在當生物量增大到一定程度時遙感信息發(fā)生飽和的現(xiàn)象,這一問題已成為影響生物量估測精度的主要原因,進一步研究解決遙感數(shù)據(jù)估測森林生物量普遍存在的數(shù)據(jù)飽和問題十分迫切。由于林木生物量中58%由樹干部分構成,而植被在生長最為繁茂的季節(jié),其樹干部分的信息被表層樹冠掩蓋,導致傳感器記錄的信息不完整,因此出現(xiàn)了以Landsat 8_OLE數(shù)據(jù)為數(shù)據(jù)源的生物量估測在高密度地區(qū)以及森林結構復雜地區(qū)容易出現(xiàn)飽和點問題[19]。
表8 “刀切法”殘差分析Tab.8 Residual analysis using Jackknife method
本文以Landsat8_OLI遙感影像為主要數(shù)據(jù)源,結合森林資源規(guī)劃設計調查統(tǒng)計數(shù)據(jù),探索參數(shù)模型(逐步線性回歸)和非參數(shù)模型(KNN)兩種算法在云南省金沙江流域中段森林生物量遙感估測中的不同效果,通過對比兩種建模方法在生物量估測中的效果,選擇適合的模型,實現(xiàn)對云南省金沙江流域中段森林中5類優(yōu)勢樹種或樹種組的生物量遙感估測,并就生物量估測中普遍存在的飽和點問題進行定量分析比較。研究區(qū)樣本生物量飽和值為:其他針葉林(160.201 t/hm2)>高山松林(127.938 t/hm2)>常綠櫟類林(101.220 t/hm2)>云南松林(98.870 t/hm2)>其他闊葉林(88.432 t/hm2),該值低于趙盼盼[12]研究所得的森林生物量光飽和值(馬尾松林159 t/hm2;杉木林143 t/hm2;闊葉樹林123 t/hm2),這可能是因為生物量光飽和點受多種因素影響,云南省金沙江流域中段山地地形復雜,立體氣候明顯,垂直高差大,且林分結構繁復,森林異質性水平更高。
合適的建模方法可以在一定程度上減小數(shù)據(jù)飽和的影響,逐步性線回歸模型和KNN模型已成為遙感估測森林生物量常用方法,從模型獨立性檢驗來看,KNN模型下各優(yōu)勢樹種或樹種組的R2都比逐步性線回歸模型下的高,表明KNN模型的估測精度高于逐步線性回歸模型的估測精度,這與戚玉嬌[20]、郭穎[21]、韓宗濤[22]等的研究結果一致。從不同的生物量分段結果來看,兩個模型均存在低值高估和高值低估情況,但總體上,不同優(yōu)勢樹種或樹種組在所有生物量分段上的ME和MRE均表現(xiàn)為KNN模型優(yōu)于一般逐步線性回歸模型,可見KNN模型有著更為精確的模型估測能力,且在一定程度上能夠減小數(shù)據(jù)飽和引起的估測誤差,但當生物量段>150 t/hm2時,常綠櫟類林和高山松林的ME均大于90 t/hm2,且MRE值均大于1,出現(xiàn)顯著的高值低估情況,說明在高生物量段兩模型估測能力均顯不足。
本研究利用半變異函數(shù)求解森林生物量光飽和值為流域森林生物量光飽和值求解提供了一定的參考,且對比分析了逐步線性回歸和KNN模型,結果表明,KNN模型在一定程度上提高了生物量的估測精度,但與實際值對比仍然具有明顯的高值低估情況。由于本研究所選取的遙感變量因子較為單一,模型的選取也只考慮了KNN模型和逐步線性回歸模型,然而實際的生物量估測還受其他環(huán)境變量因素的影響。建模過程中影響因素眾多,各種方法各有優(yōu)勢,每種方法的適用均存在一定的局限性,只能根據(jù)數(shù)據(jù)客觀條件優(yōu)選某一種或多種方法,今后的研究可以結合多元遙感數(shù)據(jù),考慮更多的遙感變量和環(huán)境因子以及森林信息來估測森林生物量進而提高估測精度。
通過對5類樹種或樹種組進行函數(shù)求解,參數(shù)方程組合,求出各樹種或樹種組的光學遙感估測生物量光飽和值:其他針葉林160.201 t/hm2,高山松林127.938 t/hm2,常綠櫟類林101.220 t/hm2,云南松林98.870 t/hm2,其他闊葉林88.432 t/hm2。
逐步線性回歸模型和KNN模型對森林生物量均具有較好的估測能力,其中,KNN模型的估測精度及擬合效果明顯優(yōu)于逐步線性回歸模型。生物量殘差分段分析表明,兩模型均存在低值高估和高值低估的情況,且在低生物量段(<50 t/hm2)時,除其他闊葉林和逐步線性回歸模型估測的其他針葉林的MRE為正值外,其余均為負值。在(>50 t/hm2)生物量段時,均為正值。在中、低生物量段(<100 t/hm2)時,KNN模型預估精度高于逐步線性回歸模型,在高生物量段(>150 t/hm2)時,兩種模型的預估精度均不高。整體來看,KNN模型的預估精度遠高于逐步線性回歸模型。