基于MODIS時間序列數(shù)據(jù)的竹林地上生物量估算

楊紹欽，王 翔，許 澄，商天其

（1.浙江省森林資源監(jiān)測中心，浙江 杭州 310020；2.浙江省林業(yè)勘測規(guī)劃設計有限公司，浙江 杭州 310020）

森林地上生物量(AGB)是表征森林生產(chǎn)力和陸地生態(tài)碳循環(huán)的重要指標，對研究全球氣候變化、人類經(jīng)營活動等對森林生態(tài)系統(tǒng)的影響具有重要意義[1]。精確估算森林地上生物量對森林可持續(xù)發(fā)展、增加森林碳匯能力及減緩全球氣候變暖具有關鍵作用[2]。樣地調查是精確獲取森林地上生物量的主要手段[3-5]，但該方法耗時費力[6]，難以實現(xiàn)大范圍時空動態(tài)監(jiān)測。近年來，眾多研究者利用遙感變量信息(如光譜特征、植被指數(shù)、冠層結構參數(shù)、紋理等)，結合樣地調查構建遙感估算模型，實現(xiàn)了對地上生物量的時空估算。森林地上生物量遙感估算有參數(shù)和非參數(shù)模型[7-12]，其中，以支持向量回歸(SVR)模型為代表的非參數(shù)模型被廣泛應用于森林地上生物量時空估算研究[13-15]。LI等[16]比較了不同模型(如最小二乘法模型、廣義加性模型、隨機森林、SVR等)在估算地上生物量上的效果后發(fā)現(xiàn)，SVR模型估算效果最好；FENG等[14]發(fā)現(xiàn)隨機森林、SVR模型等非參數(shù)模型估算的地上生物量精度要高于參數(shù)模型(線性模型和乘數(shù)模型)，且SVR模型估算精度最高。

竹林具有“爆發(fā)式生長”的特點。在生長和抽枝展葉期，竹林葉面積指數(shù)(LAI)逐漸增加，地上生物量快速積累[16-17]。LAI時間序列上的變化趨勢能夠反映森林的季節(jié)性變化，并常作為森林生態(tài)系統(tǒng)碳水循環(huán)和地上生物量估算等研究的重要參數(shù)之一[18-20]。相較于歸一化植被指數(shù)(NDVI)，增強型植被指數(shù)(EVI)對高生物量區(qū)域的敏感性更高，受土壤背景和大氣效應的影響更小[21-22]。比值指數(shù)(RVI)是綠色植被的靈敏指示參數(shù)，能夠反映葉干生物量和葉綠素含量的動態(tài)變化。目前，對竹林地上生物量的遙感估算研究多是基于單一時相遙感植被指數(shù)、光譜信息、紋理等特征進行，估算精度較低[23-24]，且不能表征竹林快速生長過程中地上生物量快速積累的動態(tài)變化[6]。本研究以浙江省竹林資源為研究對象，以2014年浙江省中分辨率成像光譜儀(MODIS)的LAI、EVI、RVI時間序列產(chǎn)品為數(shù)據(jù)源，結合2014年竹林地上生物量調查數(shù)據(jù)，構建SVR模型，估算浙江省竹林地上生物量，以期為竹林碳匯遙感監(jiān)測提供參考。

1 研究方法

1.1 研究區(qū)概況

浙江省 (27°06′～31°11′N，118°01′～123°10′E)位于中國東南沿海、長江三角洲南翼，是中國竹林面積最大的省份之一，竹資源十分豐富，竹產(chǎn)業(yè)發(fā)達，素有“世界竹子看中國，中國竹子看浙江”的美譽。目前，浙江省竹林面積為94萬 hm2，其中毛竹Phyllostachys edulis林83萬 hm2，竹林面積占全省森林面積約15%。

1.2 研究數(shù)據(jù)

1.2.1 MODIS 數(shù)據(jù)及預處理 MODIS LAI (MOD15A2H)、地表反射率 (MOD09A1)和 MODIS NDVI(MOD13Q1)數(shù)據(jù)產(chǎn)品在美國國家航空航天局(NASA)官網(wǎng)下載。MODIS LAI和地表反射率的時間和空間分辨率分別為8 d和500 m，MODIS NDVI的時空分辨率分別為16 d和250 m。數(shù)據(jù)利用重投影軟件MRT (MODIS Reprojection Tools)投影，采用最近鄰域法將空間分辨率統(tǒng)一重采樣為 500 m×500 m?；贛ODIS NDVI時間序列產(chǎn)品和地表反射率數(shù)據(jù)提取浙江省竹林豐度信息。首先，利用NDVI數(shù)據(jù)和最大似然法提取林地專題信息；其次，掩膜提取林地的地表反射率，并進行最小噪聲變換；然后，利用影像端元法提取竹林、針葉林和闊葉林端元光譜曲線；最后，利用完全約束最小二乘混合像元分解法，得到浙江省竹林豐度圖。竹林豐度提取方法參考LI等[25]和MAO等[26]。

1.2.2 地面調查數(shù)據(jù) 調查和收集182個不同立竹密度、立地質量和經(jīng)營條件的標準竹林樣地(30 m×30 m)。調查因子包括胸徑(DBH)和年齡(0～1年生的竹子是1度竹，2～3年生是2度竹，4～5年生是3度竹)，記錄樣地中心坐標?；跇拥刂窳中貜胶湍挲g，分別計算毛竹和雷竹Ph.praecox單竹的地上生物量[19]。

其中：M毛和M雷分別表示毛竹和雷竹單竹地上生物量；DBH和A分別為胸徑和竹齡。為與MODIS遙感數(shù)據(jù)匹配，本研究將樣地竹林面積線性擴大為500 m×500 m[11]。

1.3 遙感變量設置

設置141個遙感變量，其中LAI、EVI和RVI時間序列數(shù)據(jù)各46個，各序列年均值各1個。由于MODIS產(chǎn)品易受大氣、氣溶膠、雪覆蓋等因素影響，數(shù)據(jù)存在誤差大、不連續(xù)等問題[20]，本研究利用粒子濾波同化算法提高MODIS LAI時間序列數(shù)據(jù)的精度，主要步驟為：①利用3次樣條帽蓋平滑算法對MODIS LAI產(chǎn)品進行平滑處理[27]；②將平滑的MODIS LAI輸入到LAI動態(tài)模型，模擬葉面積指數(shù)；③將MODIS地表反射率的紅光和近紅外反射率、PROSAIL模型模擬冠層反射率、模擬LAI值輸入到粒子濾波同化算法，得到葉面積指數(shù)同化結果。粒子濾波同化算法詳細參考LI等[17]和MAO等[26]。增強型植被指數(shù)和比值指數(shù)利用MODIS地表反射率的第1 (red)、第2 (nir)和第3 (blue)波段分別計算。

其中：IEVI、IRVI分別表示增強型植被指數(shù)(EVI)和比值指數(shù)(RVI)，ρnir為近紅外波段反射率，ρred為紅光波段反射率，ρblue為藍光波段反射率。利用mSG平滑算法對EVI和RVI進行處理，以減少噪聲影響[28]。

1.4 變量篩選及 SVR 模型構建

變量數(shù)量的增加可能會導致“維數(shù)災難”，造成運算復雜，處理速度下降[19]。為充分利用模型輸入變量的特征，需要選擇最優(yōu)變量特征構建SVR模型。隨機森林算法預測能力強，魯棒性較好，適合處理高維數(shù)據(jù)集，計算自變量對因變量的相對重要性[29]。因此，可采用隨機森林算法選擇SVR模型的輸入變量，即通過計算均方差增量百分數(shù)來度量和評估各變量特征對竹林地上生物量的相對重要性。隨機森林參數(shù)的最優(yōu)ntree值為2 000，mtry采用模型默認值，即變量個數(shù)的1/3[25]。

SVR通過選擇合適的核函數(shù)將樣本數(shù)據(jù)從低維空間轉換到高維空間，使得低維空間內的非線性問題可在高維空間中線性處理，并保證其良好的泛化能力。本研究設置radial、linear、polynomial和sigmoid等4種核函數(shù)[30]，比較分析后選擇最優(yōu)核函數(shù)構建竹林地上生物量估算模型。采用隨機不重復抽樣法從182個竹林樣本中抽取50%作為訓練樣本，剩余的50%作為驗證樣本。

1.5 模型精度評價

利用決定系數(shù)(R)、均方根誤差(RMSE)和P(t檢驗)評價SVR模型模擬地上生物量的精度、誤差和顯著性水平，R較高且RMSE較低，說明模型估算結果較好。計算公式如下：

其中：N為樣本個數(shù)，i=1, 2 ,···,N；ym表示模型模擬值，yo表示實測值；為N個實測值的平均值。

2 結果與分析

2.1 變量篩選

基于隨機森林算法得到EVI、RVI、LAI等141個輸入變量對竹林AGB影響的重要性得分，由圖1可知：不同時間的葉面積指數(shù)、增強型植被指數(shù)和比值指數(shù)對竹林地上生物量影響存在較大的差異，其中LAI時間序列變量占比最大，RVI變量占比最小。在前50個影響較大的變量中，有43個對竹林地上生物量影響的重要性得分大于0.10%，LAI年均值(LAI_Ann)重要性得分最大(2.29%)，LAI時間序列變量中，LAI_153(第153天LAI值)最大，為1.44%，LAI_113最小，為0.10%。因此選取變量重要性得分大于0.10%的43個變量作為SVR模型輸入變量，構建竹林地上生物量估算模型。

圖1 輸入變量對竹林地上生物量影響的重要性得分Figure 1 Importance score of input variables on aboveground biomass of bamboo forest

2.2 SVR 模型構建及分析

利用篩選的變量，分別以radial、linear、polynomial和sigmoid等4種核函數(shù)構建竹林地上生物量估算模型。由圖2可知：基于radial函數(shù)的SVR模型的訓練和測試精度(R)均最高，RMSE最低，估算地上生物量和實測地上生物量的精度(R)分別為0.76和0.72，RMSE分別為5.15和8.03 Mg·hm-2；其次是linear函數(shù)，訓練和測試精度(R)分別為0.64和0.60，RMSE分別為5.95和8.38 Mg·hm-2；sigmoid函數(shù)訓練和測試精度(R)最低，均＜0.15，RMSE最大，均＞9 Mg·hm-2。polynomial函數(shù)構建的SVR模型在估算竹林地上生物量的訓練和測試精度均較高，但測試樣本中估算地上生物量時出現(xiàn)負值，因此本研究選擇基于radial核函數(shù)構建的SVR模型估算浙江省竹林地上生物量。

基于radial核函數(shù)構建SVR模型，分別利用單一變量時間序列數(shù)據(jù)和植被指數(shù)組合(EVI+RVI)測試模型精度。由圖3可見：利用單一變量(EVI或RVI)時間序列數(shù)據(jù)構建的SVR模型估算竹林地上生物量，精度較低(R＜0.65)，RMSE較大；利用單一變量LAI時間序列數(shù)據(jù)和植被指數(shù)組合(EVI+RVI)模擬的AGB的精度(R＞0.65)較高，RMSE較小。結合圖2A可知：基于radial核函數(shù)的SVR模型，利用LAI、EVI和RVI時間序列數(shù)據(jù)估算竹林地上生物量，訓練精度和測試精度均最高，誤差均最低。即利用LAI、EVI、RVI時間序列數(shù)據(jù)，基于radial核函數(shù)構建SVR模型，可有效模擬竹林地上生物量的時空分布。

圖2 不同核函數(shù)構建的SVR模型估算竹林地上生物量的精度Figure 2 Estimating accuracy of 4 kernel functions in SVR model for aboveground biomass of bamboo forest

圖3 基于 SVR 模型的 LAI、EVI、RVI及 EVI+RVI時間序列數(shù)據(jù)估算竹林地上生物量Figure 3 Model accuracy of LAI, EVI, RVI and EVI+RVI time series data for estimating aboveground biomass of bamboo forest using SVR model

2.3 竹林地上生物量空間分布估算

基于radial核函數(shù)構建的SVR模型，結合竹林豐度信息，估算得到浙江省竹林地上生物量統(tǒng)計變化值，由表1可知：浙江省全省竹林地上生物量均值為7.85 Mg·hm-2，總地上生物量為3.31×107Mg；浙江省竹林地上生物量在各市具有明顯的差異性，其中：湖州市、杭州市、金華市、紹興市和寧波市的竹林地上生物量均值均大于全省均值；湖州市竹林地上生物量均值最大為13.56 Mg·hm-2，總地上生物量為2.13×106Mg，占全省總地上生物量的6.44%，但其竹林面積僅占3.73%；杭州市地上生物量均值為9.86 Mg·hm-2，總地上生物量為8.48×106Mg，占全省總地上生物量的25.42%，竹林面積占20.24%；舟山市地上生物量均值最低，為5.72 Mg·hm-2，總地上生物量為6.89×104Mg，占全省總地上生物量的0.21%，竹林面積占0.29%。

表1 基于SVR模型模擬的浙江省竹林地上生物量統(tǒng)計值Table 1 Statistical values of aboveground biomass of bamboo forest in Zhejiang Province based on SVR model

3 討論

本研究基于葉面積指數(shù)、增強型植被指數(shù)和比值指數(shù)等多變量，構建并優(yōu)選SVR模型，實現(xiàn)浙江省竹林地上生物量較高精度的估算。SVR模型耦合LAI、EVI、RVI時間序列數(shù)據(jù)估算竹林地上生物量與實測值間仍存在一定的誤差。當實測值較小時(＜20 Mg·hm-2)，估算值高于實測值，即低值高估；而實測值較高時(＞30 Mg·hm-2)，估算的地上生物量低于實測數(shù)據(jù)，即高值低估。可見本研究構建的SVR模型并不能完全消除地上生物量估算的飽和現(xiàn)象。

本研究發(fā)現(xiàn)：基于單一變量的SVR模型模擬竹林地上生物量效果較差，而多變量的參與尤其是LAI的參與有效提高了模型預測精度。本研究中，植被指數(shù)(EVI+RVI)構建的SVR模型模擬竹林地上生物量的精度分別為0.68和0.71，RMSE分別為5.79和8.13 Mg·hm-2，而加入LAI后，R分別提高了11.76%和1.41%，而RMSE分別降低了11.05%和1.23%，模型預測能力明顯增加。一方面，葉面積指數(shù)是反映植被固碳能力的重要參數(shù)，能夠體現(xiàn)竹林觀測結構及時間動態(tài)；另一方面，本研究采用粒子濾波同化算法對MODIS LAI時間序列數(shù)據(jù)進行同化[17,26]，使得葉面積指數(shù)年際動態(tài)變化趨勢與竹林生長特征相吻合[31]，減弱了信號飽和對地上生物量估算的影響。

已有研究[7,11]認為：調查樣地與遙感影像尺度匹配的差異會加大SVR模型的估算誤差。為與MODIS數(shù)據(jù)匹配，本研究將30 m×30 m竹林樣地進行了線性轉換，但由于沒有考慮混合像元的影響，樣地尺度轉換出現(xiàn)誤差，并由此加大了模型估算誤差。相關研究[6,25]表明：在大尺度森林地上生物量估算中利用森林豐度數(shù)據(jù)能夠減少估算的不確定性。因此，本研究采用混合像元分解方法提取了竹林豐度信息，并乘以模型模擬的竹林地上生物量，進而得到浙江省竹林地上生物量時空分布，盡可能地降低了混合像元對生物量空間估算的影響。

4 結論

本研究基于MODIS LAI、EVI、RVI時間序列數(shù)據(jù)，采用支持向量回歸模型結合隨機森林算法，在變量篩選和模型優(yōu)選的基礎上，實現(xiàn)了浙江省竹林地上生物量較高精度的估算。結果表明：①基于隨機森林模型篩選出43個對竹林地上生物量影響最大的變量，其中LAI時間序列數(shù)據(jù)的貢獻最大；②基于radial核函數(shù)耦合LAI、EVI、RVI時間序列數(shù)據(jù)，構建的SVR模型訓練和測試精度均最高(R分別為0.76和0.72)，均方根誤差均最低(分別為5.15和8.03 Mg·hm-2)；③由SVR模型估算得到浙江省竹林地上生物量，全省竹林地上生物量均值為7.85 Mg·hm-2，總地上生物量為3.31×107Mg；浙江省竹林地上生物量在各市具有明顯的差異性，其中，湖州市、杭州市、金華市、紹興市和寧波市的竹林地上生物量均值均大于全省均值，表明其竹林經(jīng)營效果相對較好。研究結果為大范圍竹林碳匯遙感精準監(jiān)測提供了較好的方法。