黃屹杰，張加龍，胡耀鵬，程 滔

（1. 西南林業(yè)大學(xué) 林學(xué)院，云南 昆明 650233；2. 國(guó)家基礎(chǔ)地理信息中心 調(diào)查監(jiān)測(cè)部，北京 100830）

森林生物量是森林生態(tài)系統(tǒng)的最基本數(shù)量特征，是研究林業(yè)和生態(tài)問(wèn)題的基礎(chǔ)[1]。地上生物量(aboveground biomass， AGB)是森林生態(tài)系統(tǒng)生產(chǎn)力的重要指標(biāo)和質(zhì)量的綜合體現(xiàn)，對(duì)研究全球氣候變化具有重要意義。目前，森林地上生物量的調(diào)查方法有地面調(diào)查和遙感監(jiān)測(cè)2種。傳統(tǒng)的地面調(diào)查法耗時(shí)長(zhǎng)、成本高，而遙感監(jiān)測(cè)方法快速、無(wú)損，具有在更大尺度上應(yīng)用的優(yōu)勢(shì)，已成為地上生物量估測(cè)的主流[2-3]。在森林調(diào)查中，常使用統(tǒng)計(jì)模型來(lái)預(yù)測(cè)樣地內(nèi)個(gè)體樹(shù)木的生物量，然后對(duì)同一樹(shù)種的單株地上生物量進(jìn)行匯總，來(lái)作為遙感估算生物量的訓(xùn)練數(shù)據(jù)和精度評(píng)估數(shù)據(jù)[4]。用于地上生物量估測(cè)的遙感數(shù)據(jù)源包括光學(xué)遙感、機(jī)載雷達(dá)、激光測(cè)量數(shù)據(jù)等[5]。與成本較高的機(jī)載雷達(dá)數(shù)據(jù)相比，光學(xué)影像穿透森林冠層獲取樹(shù)干信息的能力較差，但隨著新算法的不斷改進(jìn)，使用光學(xué)遙感影像進(jìn)行地上生物量估算的精度也在不斷提高[6]。

目前，遙感估算地上生物量在因子選取、建模方法、數(shù)據(jù)飽和等方面存在較多不確定性問(wèn)題[7-8]。區(qū)域尺度森林地上生物量遙感估算包含測(cè)量和抽樣、單株生物量模型、遙感因子和地形因子等不確定性來(lái)源[9-10]。其中，與測(cè)量和抽樣誤差相比，模型的不確定性對(duì)地上生物量估算的影響較大，如SHETTLES等[11]利用激光雷達(dá)和地面數(shù)據(jù)估測(cè)地上生物量發(fā)現(xiàn)：模型不確定性占總不確定性的55%。同時(shí)，用于建模的生物量數(shù)據(jù)是使用該區(qū)域樹(shù)種的異速生長(zhǎng)方程計(jì)算得出的。因此，單株生物量模型的誤差會(huì)傳播到樣地級(jí)的生物量估測(cè)，從而影響遙感估算生物量的不確定性[12-13]?，F(xiàn)階段不確定性的量化方法主要分為3種[14]：過(guò)程模型分析法、隨機(jī)誤差傳遞法和Monte Carlo模擬法。

本研究基于云南省香格里拉市Landsat 8影像和外業(yè)調(diào)查數(shù)據(jù)，建立了生物量遙感估測(cè)模型，分析了遙感估算樣地高山松Pinus densata地上生物量的總不確定性，以期為提高森林地上生物量估算精度提供參考。

1 研究地區(qū)與數(shù)據(jù)

1.1 研究區(qū)概況

研究區(qū)位于云南省西北部的迪慶藏族自治州香格里拉市 (26°52′～28°52′N(xiāo)，99°20′～100°19′E)，該區(qū)總面積為11 613 km2，平均海拔為3 459 m。香格里拉市的森林覆蓋率較高，達(dá)75%，主要植被類(lèi)型為寒溫性針葉林，優(yōu)勢(shì)樹(shù)種有云杉Picea asperata、冷杉Abies fabri、高山松、云南松Pinus yunnanensis、高山櫟Quercus semicarpifolia等。

1.2 數(shù)據(jù)來(lái)源

于2015年11月和2016年3月，在每個(gè)鄉(xiāng)鎮(zhèn)高山松純林區(qū)域隨機(jī)布設(shè)60塊樣地，大小為30 m×30 m，且每塊樣地間隔3 km以上。樣地高山松林分的齡級(jí)主要為近熟林和成熟林。外業(yè)調(diào)查記錄胸徑大于5 cm樣木的樹(shù)高、胸徑，同時(shí)對(duì)每塊樣地的林分密度進(jìn)行了計(jì)算(表1)。

表1 外業(yè)調(diào)查樣地統(tǒng)計(jì)Table 1 Basic statistics of the sample plots of field surveys

Landsat 8影像數(shù)據(jù)來(lái)源于美國(guó)地質(zhì)調(diào)查局(http://glovis.usgs.gov/)，3景數(shù)據(jù)為2015年11—12月的成像，云量均小于2%。

2 研究方法

2.1 樣地地上生物量計(jì)算

采用每木檢尺，通過(guò)胸徑、樹(shù)高和單株生物量模型[15]計(jì)算每株地上生物量，進(jìn)而得出每塊樣地的地上生物量。單株生物量模型()如下：

2.2 Landsat數(shù)據(jù)處理

對(duì)影像數(shù)據(jù)進(jìn)行輻射定標(biāo)，采用FLAASH方法進(jìn)行大氣校正。以研究區(qū)SPOT-5影像為參考數(shù)據(jù)，選取100個(gè)地面控制點(diǎn)，采用二項(xiàng)式方法對(duì)影像進(jìn)行幾何校正，使用雙線性?xún)?nèi)插法將影像重新采樣為30 m×30 m，誤差控制在1個(gè)像元內(nèi)[8]；接著，采用坡度匹配模型[16]進(jìn)行地形校正。具體操作方法參考文獻(xiàn)[17]。

對(duì)預(yù)處理后的影像數(shù)據(jù)提取4類(lèi)遙感因子：①原始波段因子，分別為C、B1、B2、B3、B4、B5、B7；②植被指數(shù)因子，分別為B43、B42、B53、B54、B57、B73、B74、B3Albedo、 B473、NDVI、ND32、ND53、ND54、ND57、ND452、DVI；③信息增強(qiáng)因子，分別為VIS123、Albedo、MID57；④紋理信息因子[18]，分別為均值(ME)、方差(VA)、同質(zhì)性(HO)、反差(CO)、相異(DI)、熵(EN)、角二階矩(SM)、相關(guān)性(CR)、偏度(SK)。用R5和R9分別代表5×5和9×9窗口，如R9B5CR代表9窗口第5波段的相關(guān)性紋理。

岳彩榮[19]提取了5×5、7×7、9×9、15×15的4種窗口紋理信息因子對(duì)香格里拉森林生物量進(jìn)行了遙感估測(cè)；張加龍等[20]基于遙感影像和連續(xù)清查固定樣地對(duì)高山松地上生物量進(jìn)行估算時(shí)，使用了5×5、9×9的2種窗口紋理信息因子。以上研究表明：高山松地上生物量與5×5和9×9窗口的紋理信息因子相關(guān)性更強(qiáng)。因此，本研究的紋理信息因子選用5×5和9×9窗口。4類(lèi)遙感因子的提取和計(jì)算參考文獻(xiàn) [8， 17， 20]。

2.3 因子篩選方法

基于SPSS軟件對(duì)樣地地上生物量與備選遙感因子進(jìn)行相關(guān)性分析，結(jié)果按照Pearson系數(shù)進(jìn)行排序。本研究因子入選的顯著性水平為P≤0.05，因子剔除的水平為P≥0.10，得到相關(guān)性排前14個(gè)的遙感因子。為避免解釋變量之間的多重共線性問(wèn)題，將各個(gè)遙感因子之間相關(guān)性較高的因子剔除，最終篩選9個(gè)因子參與模型建立，分別為B74、9窗口第5波段相關(guān)性(R9B5CR)、9窗口第6波段方差(R9B6VA)、9窗口第3波段熵(R9B3EN)、5窗口第5波段角二階矩(R5B5SM)、9窗口第6波段相關(guān)性(R9B6CR)、5窗口第3波段熵(R5B3EN)、5窗口第4波段熵(R5B4EN)、5窗口第4波段偏度(R5B4SK)。

2.4 建模方法

遙感估算高山松地上生物量的3種方法分別為多元線性回歸(MLR)、梯度提升回歸樹(shù)(GBRT)以及隨機(jī)森林(RF)。從60個(gè)樣地調(diào)查數(shù)據(jù)中隨機(jī)選取48個(gè)(80%)進(jìn)行建模，剩余12個(gè)(20%)用來(lái)獨(dú)立性檢驗(yàn)。建模和檢驗(yàn)樣本的地上生物量統(tǒng)計(jì)結(jié)果見(jiàn)表2。

表2 建模和檢驗(yàn)樣本的地上生物量實(shí)測(cè)值Table 2 Measured values of aboveground biomass of modeling and testing samples

2.4.1 多元線性回歸 多元線性回歸模型可以同時(shí)對(duì)多個(gè)解釋變量與1個(gè)因變量進(jìn)行擬合，用回歸方程來(lái)表示擬合關(guān)系。本研究多元線性回歸建模在SPSS軟件中實(shí)現(xiàn)。

2.4.2 梯度提升回歸樹(shù) 該方法通過(guò)構(gòu)建多個(gè)弱分類(lèi)器，經(jīng)過(guò)多次迭代最終組合成1個(gè)強(qiáng)分類(lèi)器[21]。本研究基于Python語(yǔ)言的Sklearn工具包提供的梯度提升回歸樹(shù)算法，對(duì)高山松地上生物量與遙感因子進(jìn)行建模分析。在建模分析前需要對(duì)弱分類(lèi)器的最大迭代次數(shù)(n_estimators)、學(xué)習(xí)速率(learning_rate)、最大深度(max_depth)等參數(shù)進(jìn)行確定。參數(shù)確定往往需要遵循一定的經(jīng)驗(yàn)，比如最大迭代次數(shù)通常在預(yù)測(cè)值收斂的情況下越小越好，學(xué)習(xí)速率通常小于0.10，最大深度通常不應(yīng)大于15。依據(jù)這些經(jīng)驗(yàn)設(shè)置上述參數(shù)的范圍，比較參數(shù)值組合下模型學(xué)習(xí)的效果，以得到模型輸出預(yù)測(cè)結(jié)果[8]。本研究根據(jù)擬合優(yōu)度的最佳參數(shù)組合，選取最大迭代次數(shù)為60，子采樣為0.5，學(xué)習(xí)速率為0.05，決策樹(shù)最大深度為7，葉子節(jié)點(diǎn)最少樣本數(shù)為3。

2.4.3 隨機(jī)森林 隨機(jī)森林適用于多數(shù)分類(lèi)與回歸的問(wèn)題，由一系列決策樹(shù)組成?；赑ython語(yǔ)言的Sklearn工具包進(jìn)行高山松地上生物量建模分析。程序使用bootstrap重采樣方法從樣本集中提取多個(gè)重采樣樣本進(jìn)行建模，之后對(duì)多個(gè)模型值進(jìn)行預(yù)測(cè)并組合，再通過(guò)投票得出最終值，一般會(huì)在隨機(jī)產(chǎn)生的分類(lèi)樹(shù)中選出重疊次數(shù)最多的決策樹(shù)作為最終模型[6]。本研究初步設(shè)置決策樹(shù)數(shù)量為50，逐步增加決策樹(shù)數(shù)量模擬建模過(guò)程，回歸誤差趨于穩(wěn)定，最終確定決策樹(shù)數(shù)量為300。

2.5 精度評(píng)價(jià)方法

采用的評(píng)價(jià)指標(biāo)包括決定系數(shù)(R2)、均方根誤差(RMSE)、預(yù)測(cè)精度(F)。各指標(biāo)計(jì)算如下：

2.6 不確定性分析與度量

2.6.1 模型殘差變異不確定性計(jì)算 樣地尺度上單株生物量模型的殘差不確定性()計(jì)算如下：

采用六步法[23]計(jì)算地上生物量遙感估算模型殘差變異產(chǎn)生的誤差：第1步升序排列樣地地上生物量實(shí)測(cè)值(y)；第2步利用生物量模型的預(yù)測(cè)值()計(jì)算殘差()，殘差為樣地地上生物量實(shí)測(cè)值與預(yù)測(cè)值的差值；第3步將60塊樣地進(jìn)行分組，每6塊地為1組，共10組；第4步計(jì)算每組樣地地上生物量預(yù)測(cè)值的平均值()及殘差標(biāo)準(zhǔn)差。每組樣地地上生物量預(yù)測(cè)值的平均值、殘差、殘差標(biāo)準(zhǔn)差計(jì)算如下：

第6步將各塊樣地地上生物量預(yù)測(cè)值代入擬合后得出的公式，對(duì)所有樣地的殘差標(biāo)準(zhǔn)差求和，除以地上生物量實(shí)測(cè)值的和，從而得到模型殘差的不確定性。

2.6.2 模型參數(shù)不確定性計(jì)算 對(duì)于單株生物量模型和多元線性回歸模型的參數(shù)不確定性，可通過(guò)泰勒級(jí)數(shù)一階展開(kāi)式進(jìn)行量化計(jì)算。生物量模型經(jīng)泰勒級(jí)數(shù)一階展開(kāi)如下：

3 結(jié)果與分析

3.1 樣地尺度上地上生物量的不確定性

單株生物量模型參數(shù)的協(xié)方差矩陣為：

將協(xié)方差矩陣帶入式(11)，得出單株生物量模型的參數(shù)不確定性為12.42%。

由式(5)計(jì)算模型的殘差不確定性為10.76%。對(duì)殘差不確定性與參數(shù)不確定性合成，可得出單株生物量模型的不確定性為16.43%。當(dāng)模型誤差疊加到樣地尺度時(shí)，通過(guò)式(14)得出樣地尺度不確定性為7.07%。

3.2 多元線性回歸 (MLR)模型的不確定性

本研究建立了地上生物量與遙感因子的多元線性回歸模型，根據(jù)擬合優(yōu)度選出最佳模型為：

其中：y表示地上生物量預(yù)測(cè)值，xR5B4EN為5窗口第4波段熵，xR9B5CR為9窗口第5波段相關(guān)性，xR9B3EN為9窗口第3波段熵。模型R2為0.326，RMSE為26.12，預(yù)測(cè)精度為64.91%(圖1A)。

對(duì)殘差標(biāo)準(zhǔn)差(y)和樣地地上生物量預(yù)測(cè)分組平均值(x)進(jìn)行相關(guān)性擬合發(fā)現(xiàn)：擬合效果最佳的模型為計(jì)算后得出多元線性回歸模型的殘差不確定性為34.86%(圖 1B)。

圖1 基于多元線性回歸的地上生物量預(yù)測(cè)效果和殘差不確定性Figure 1 Prediction effect and the residual uncertainty of aboveground biomass based on multiple linear regression

多元線性回歸模型參數(shù)的協(xié)方差矩陣為：

將協(xié)方差矩陣代入式(11)可得出參數(shù)變異不確定性為21.30%。對(duì)殘差不確定性和參數(shù)不確定性合成，通過(guò)多元線性回歸模型，估算出高山松地上生物量的不確定性為40.84%。

3.3 梯度提升回歸樹(shù) (GBRT)模型的不確定性

圖2A所示：梯度提升回歸樹(shù)模型的R2為0.815，RMSE為14.24，預(yù)測(cè)精度為74.72%。圖2B所示：殘差標(biāo)準(zhǔn)差(y)與樣地地上生物量預(yù)測(cè)分組平均值(x)擬合的最佳模型為 y =0.260 1x-2.127 3，R2=0.670，得出殘差不確定性為22.01%。非參數(shù)建模方法不考慮參數(shù)不確定性，基于梯度提升回歸樹(shù)模型估算得出高山松地上生物量的不確定性為22.01%。

圖2 基于梯度提升回歸樹(shù)的地上生物量預(yù)測(cè)效果和殘差不確定性Figure 2 Prediction effect and the residual uncertainty of aboveground biomass based on gradient boost regression tree

3.4 隨機(jī)森林 (RF)模型的不確定性

圖3 A所示：隨機(jī)森林模型的R2為0.889，RMSE為11.02，預(yù)測(cè)精度為76.47%。圖3B所示：殘差標(biāo)準(zhǔn)差(y)與樣地生物量預(yù)測(cè)分組平均值(x)擬合后的最佳模型為線性關(guān)系，模型為 y =0.214 3x-1.475 7，R2=0.863。與多元線性回歸和梯度提升回歸樹(shù)模型的殘差不確定性計(jì)算方法相同，在不考慮參數(shù)不確定性的情況下，得出隨機(jī)森林估算高山松地上生物量的不確定性為18.09%。

圖3 基于隨機(jī)森林的地上生物量預(yù)測(cè)效果和殘差不確定性Figure 3 Prediction effect and the residual uncertainty of aboveground biomass based on random forest

3.5 遙感估測(cè)生物量總不確定性

由表3可見(jiàn)：在樣地尺度上，單株生物量模型的不確定性為16.43%，總不確定性為7.07%；在3種遙感估算模型中，多元線性回歸模型的不確定性為40.84%，梯度提升回歸樹(shù)模型的不確定性為22.01%，隨機(jī)森林模型的不確定性為18.09%?；诙嘣€性回歸、梯度提升回歸樹(shù)、隨機(jī)森林模型，估算出高山松地上生物量的總不確定性分別為41.45%、23.12%和19.42%。

表33 種地上生物量估測(cè)模型的不確定性Table 3 Uncertainty results of the three aboveground biomass estimation models

4 結(jié)論與討論

本研究結(jié)合誤差傳遞的方法，量化了3種生物量遙感估測(cè)模型與樣地尺度單株生物量模型產(chǎn)生的誤差，并將各個(gè)來(lái)源的不確定性合成，得出遙感估算高山松地上生物量的總不確定性。遙感估算地上生物量的不確定性主要來(lái)源于遙感估測(cè)模型的不確定性，其中多元線性回歸為40.84%，梯度提升回歸樹(shù)為22.01%，隨機(jī)森林為18.09%。樣地尺度不確定性相對(duì)較小，為7.07%。從預(yù)測(cè)效果來(lái)看，非參數(shù)的隨機(jī)森林和梯度提升回歸樹(shù)模型的不確定性要優(yōu)于參數(shù)的多元線性回歸模型，隨機(jī)森林模型的不確定性低于梯度提升回歸樹(shù)模型。單株生物量模型參數(shù)變異為12.42%，殘差不確定性為10.76%，模型的不確定性主要來(lái)源于參數(shù)變異，殘差不確定性相對(duì)較低。

3種遙感估算模型的不確定性研究表明：多元線性回歸通過(guò)回歸方程直觀反映遙感因子和地上生物量的關(guān)系，但不能有效描述它們之間復(fù)雜的非線性關(guān)系，要提高參數(shù)模型的預(yù)測(cè)精度則需要大量樣地觀測(cè)數(shù)據(jù)。張加龍等[8]基于外業(yè)調(diào)查和Landsat 8 OLI影像，建立了多元線性回歸、偏最小二乘、隨機(jī)森林和梯度提升回歸樹(shù)4種模型，估測(cè)了香格里拉市高山松地上生物量，結(jié)果表明：隨機(jī)森林和梯度提升回歸樹(shù)2種非參數(shù)模型的預(yù)估精度優(yōu)于其他2種參數(shù)模型；陳蜀蓉等[25]運(yùn)用4種模型對(duì)浙江縉云縣公益林生物量進(jìn)行建模估算的結(jié)果同樣表明：隨機(jī)森林和Erf-BP神經(jīng)網(wǎng)絡(luò)2種非參數(shù)模型優(yōu)于參數(shù)模型。本研究建立的梯度提升回歸樹(shù)和隨機(jī)森林2種非參數(shù)模型有較好的擬合優(yōu)度和預(yù)測(cè)精度，不確定性?xún)?yōu)于參數(shù)的多元線性回歸模型。遙感估測(cè)模型的不確定性對(duì)森林地上生物量的估算精度占主要的影響，因此，提升模型精度仍是今后生物量估算研究的重要方向。

已有不少學(xué)者進(jìn)行了單株生物量模型不確定性的研究。秦立厚等[10]分別用30、42、48株杉木Cunninghamia lanceolata進(jìn)行單株生物量建模時(shí)發(fā)現(xiàn)：一元生物量模型的殘差變異引起的不確定性分別為15.2%、12.3%、11.7%，二元生物量模型的不確定性為13.3%、9.4%、8.7%，表明二元單株生物量模型的不確定性要低于一元模型，且建模樣本的增加可以顯著降低模型殘差變異不確定性。與其相比，本研究在對(duì)樣地尺度生物量不確定性進(jìn)行度量時(shí)，引用了基于冪函數(shù)形式的二元單株生物量估算模型[15]。由于缺少實(shí)際伐倒的樣木生物量數(shù)據(jù)，運(yùn)用模型預(yù)測(cè)精度和均方根誤差間接計(jì)算模型不確定性，結(jié)果可能存在偏差。CHAVE等[26]計(jì)算的熱帶雨林地區(qū)生物量模型殘差變異誤差為31.3%，而本研究模型殘差變異誤差為10.76%。這可能是熱帶樹(shù)種之間的特征差異大，導(dǎo)致產(chǎn)生的殘差變異較大，而本研究建模對(duì)象僅為亞熱帶高山松單一樹(shù)種，未對(duì)抽樣、測(cè)量誤差以及遙感影像坐標(biāo)校正誤差進(jìn)行分析。這在后續(xù)研究中還需進(jìn)一步完善。