基于無人機多光譜的高山松地上生物量估測

楊正道，舒清態，黃金君，周文武，胥 麗，羅紹龍，王書偉

（1.西南林業大學 林學院，云南昆明 650224；2.廣西壯族自治區中國科學院廣西植物研究所，廣西桂林 541006）

森林是陸地重要的生態系統之一，發揮重要作用，不僅具有豐富的生物多樣性和提供眾多生態服務，還是地球上最大的碳庫[1]。森林生物量是森林生態系統的重要組成部分，包括地上生物量和地下生物量；其中，地上生物量由樹干生物量、樹枝生物量和樹葉生物量組成，一般用單位時間或單位面積內累積的干物質或能量表示[2]，在全球森林生態系統效益評估中發揮主導作用[3-5]。對森林地上生物量進行估測對于了解碳循環、評估碳匯潛力及制定森林資源管理和林業可持續發展策略有重要價值和意義[6-7]。

森林地上生物量估測方法主要包括野外實測和遙感估測[8]。基于樣地調查的野外實測生物量估測方法準確性高，但存在一些缺點，如消耗大量時間、人力和資金，破壞性較大，在偏遠和崎嶇地區實地調查困難等，不適用于大范圍森林地上生物量估測和動態變化監測[9-10]。隨著地理信息技術不斷發展，遙感（Remote Sensing，RS）數據在森林地上生物量估測方面得到廣泛應用。遙感連續動態監測和地理信息系統（Geographic Information System，GIS）空間數據分析能力不斷增強，學者們結合地面數據和衛星遙感數據進行森林地上生物量估測，主要利用的衛星遙感數據為Landsat 8光學遙感影像[11]。盡管衛星遙感技術在森林地上生物量估測方面取得了一定進展，但仍存在一些局限。首先，獲取生物量像元的真實值并不容易；森林地上生物量是通過間接推斷得出的，通常需依賴實地調查數據進行校正和驗證，這一過程需耗費大量時間和人力資源。其次，衛星遙感數據的時間分辨率和空間分辨率相對較低，可能無法捕捉森林地上生物量的細微變化和誤差導致的空間異質性[12]。

近年，無人機遙感技術以其低成本、高安全性、高時效性、高分辨率和低空影像獲取等優勢彌補了衛星遙感技術的局限[13]，在林業方面得到廣泛應用。申鑫等[14]利用無人機獲取的高光譜和高空間分辨率數據提取遙感因子，并結合樣地實測數據，建立亞熱帶天然次生林林分尺度生物量估測模型，進行生物量反演；Jones 等[15]利用無人機影像生成高分辨率三維樹冠模型，在林分尺度上對紅樹林地上生物量進行估測；Tamiminia 等[16]利用無人機多光譜影像進行紐約州柳樹（Salixspp.）林分尺度上森林地上生物量的反演；李濱等[17]基于無人機遙感技術獲取目標樣地樟子松（Pinussylvestrisvar.mongolica）的遙感影像，并結合樣地實測數據，在單株尺度上建立胸徑-生物量模型，用于估測東北林業大學城市林業示范基地內目標樣地的地上生物量；楊雪峰等[18]利用無人機低空遙感和高分辨率衛星影像，獲取塔里木河下游胡楊（Populuseuphratica）林單株尺度的森林結構參數，估測該地區的地上生物量。目前的研究主要集中在單株和林分尺度上的森林地上生物量估測，對于最佳觀測窗口和低緯度高海拔地區森林地上生物量估測，尚缺乏較深入的研究。探索最佳觀測窗口，并在不同地理環境和生態系統中驗證無人機遙感技術的適用性，有助于推動生物量估測方法的發展及其精確性的提高。

本研究通過地統計學變異函數確定高山松（Pinusdensat）地上生物量最佳觀測窗口；采用偏最小二乘法（Partial Least Squares Regression，PLS）和隨機森林（Random Forest，RF）模型對遙感影像中篩選的特征變量和地面實測生物量建模，對比分析兩種模型的精度；采用篩選出的最優模型對無人機飛行區高山松地上生物量進行估測，并分析其空間分布情況，可為進一步開展高寒山區森林生態系統物質與能量循環研究及保護和改造脆弱生態環境提供參考。

1 材料與方法

1.1 研究區概況

研究區位于云南省迪慶藏族自治州香格里拉市小中甸鎮下吉沙林場（99°49′E，27°24′～27°25′N），平均海拔3 300～3 500 m，總面積76.39 hm2，年均氣溫6.3 ℃，年均降水量649.4 mm。林場主要地貌為山地，地形呈東北高、西南低的趨勢。高山松為林場從寒溫帶向亞熱帶過渡的先鋒針葉樹種，適應性廣，再生能力強，具有喜光、耐旱和耐瘠薄等特點。具體樣地分布如圖1所示。

圖1 研究區位置與樣地分布Fig.1 Location of study area and distributions of sample sites

1.2 數據采集

1.2.1 無人機影像數據獲取與預處理

采用大疆四旋翼無人機搭載的航測系統和AURedEdge 多光譜相機（含藍光、綠光、紅光、紅邊和近紅外5 個波段，表1）采集試驗區高山松多光譜影像。航拍過程中，相機面向地面90°俯拍，地面控制點采用RTK 測定，坐標系統為GCS—WGS—1984。航飛時間為2021 年11 月1 日10:00～14:30，天氣晴朗，無風云干擾。通過5 次試飛和航測點調整，確定無人機飛行高度為260 m，航速為8.8 m/s。航向和旁向重疊度分別為75%和70%。地面分辨率為0.128 m，航測區面積為63.73 hm2。

表1 多光譜傳感器波段參數Tab.1 Band parameters of multi-spectral sensors

獲得無人機多光譜影像后，進行預處理。對無人機拍攝到的多光譜影像進行甄別，剔除姿態角不正常、成像有問題的影像。采用大疆智圖軟件（DJI Terra）進行影像數據預處理。將多光譜相機拍攝的原照片導入軟件進行二維多光譜重建，獲得不同波段反射率圖像、數字表面模型（Digital Surface Model，DSM）和數字正射影像（Digital Orthophoto Map，DOM）。采用ENVI 5.3 軟件進行校正，并將數據存儲為柵格格式（圖2）。

圖2 無人機航測區影像Fig.2 Aerial survey area image of UVA

1.2.2 地面調查數據獲取與預處理

根據無人機獲取的多光譜影像及其尺度，采集對應時間內試驗區36 個高山松角規控制檢尺樣地的胸徑、樹高和株數等信息，記錄樣地中心經緯度和高程信息。采用三鼎星T66系列RTK進行坐標點測量，誤差控制在千分之二以內；計算每個樣地的平均胸徑、平均樹高、總蓄積量和斷面積。

參考黃金君等[19]文獻，采用分層貝葉斯法擬合單木生物量模型，通過公式（1）計算不同樣地地上單木平均生物量；通過公式（2）計算每公頃林木株數，并結合地上單木平均生物量，求算樣地公頃生物量；將樣地公頃生物量轉換為與影像像元對應的樣地生物量。36塊樣地數據描述性分析見表2。

表2 高山松樣地地上生物量統計分析Tab.2 Statistical analysis on above-ground biomass of P.densata sample sites

式中，M為地上單木平均生物量（t/hm2）；D為平均胸徑（cm）；H為平均樹高（m）。

式中，N為每公頃林木株數（株/公頃）；Fg為斷面積系數；k為樣地計數木總株數；gi為樣地中第i株計數木的斷面積（cm2）；ki為樣地中第i株計數木的株數。

1.3 研究方法

1.3.1 最佳觀測窗口選擇

最佳觀測窗口常被稱為最佳分辨率，即影像中能識別最小地物的最大空間分辨率，能清楚探索和表述這一尺度層的地理過程和現象內容，同時也能明確地分析這一尺度層地理空間信息的組合規律[20]。

半方差也稱變異函數，是建立在二階平穩假設基礎上的空間變異結構分析工具，也是地統計學的主要理論基礎，在地物空間變異分析方面得到廣泛應用。1997 年，Atkinson[21]將變異函數法應用到遙感影像最佳觀測尺度研究中，通過分析變異函數中變量因子變異的空間程度，確定最佳分辨率。傳統統計學僅對像素值進行統計分析，忽視像素值間的隨機性和相關性；變異函數可以揭示區域化變量在全尺度上的空間變異模式[21]，計算公式為：

式中，γ(h)為區域變量Z（w）的變異函數值；N（h）為空間間隔距離等于h的樣本點對數量；Z（w）i和Z（wi+h）分別為區域變量在點wi和點wi+h位置上的像素值；h為兩個像素點的分隔距離。

在變異函數的計算中，結果值為散點圖。在對整個區域進行定量描述時，需對散點進行擬合。通常采用指數模型、球狀模型和高斯模型3 種方式進行擬合[20，22]。

（1）指數模型表達式為：

h=3α時，γ（h）≈C0+C1，此模型的有效變程為3α。

（2）球狀模型表達式為：

h=α時，γ（h）=C0+C，此模型的有效變程為α。

（3）高斯模型表達式為：

式中，C0為塊金值；C1為偏基臺值；C（C=C0+C1）為基臺值；α為變程或相關距。模型的變程α即最佳觀測窗口。實際應用時，以高山松重采樣尺度大小（即重采樣時像元大小）為h，以對應尺度下的像元光譜值為γ（h），求算出的變程α為最佳觀測窗口。在模型優選基礎上，利用塊金值、基臺值、變程和結構比等模型參數，對區域化變量進行空間分布特征分析，一般選用具有較大決定系數和較小殘差的理論模型擬合變異函數[23]。

變異程度為塊金值與基臺值的比值（C0/C），能反映整個空間的變異狀態，主要包括0～25%、25%～75%和75%以上3 個層次，代表空間自相關程度。空間自相關程度高（變異程度低），表明主要由結構性因素引起變異；空間自相關程度低（變異程度高），表明以隨機部分為主。

1.3.2 最佳觀測窗口設計方案

選取拼接后的無人機多光譜影像，篩選與實測生物量均值最接近的樣方作為中心點，在ENVI 5.3軟件中裁剪出大小為100 m×100 m、質量較好的一景作為觀測尺度試驗影像樣本，利用Band Math 工具對試驗影像樣本進行紋理特征指數計算。

本研究中，無人機多光譜的空間分辨率為0.125 m，數據量大，冗余多，所以利用ENVI 5.3軟件中的Resize Data工具對進行紋理特征指數計算后的影像樣本進行重采樣，將空間分辨率提升至1 m。常見的重采樣內插方法有最鄰近法、雙線性內插法和三次卷積內插法。本研究選擇最鄰近法進行重采樣。最鄰近法將最鄰近的像元值賦予新像元，輸出的圖像仍保持原來圖像的像元值，處理速度快。

采用ArcGIS 軟件將重采樣后柵格數據中每個像元的紋理特征值提取至點，計算每個點的X和Y坐標；借助地統計學GS+軟件，分別采用指數模型、球狀模型和高斯模型作為變異函數擬合模型，對高山松地上生物量最優觀測尺度進行分析，通過模型決定系數（R2）和殘差平方和（Residual Sum of Squares，RSS）對模型進行評價，獲得最佳觀測窗口。R2越大，RSS越小，模型效果越好。

1.3.3 遙感因子提取

研究表明，在高分辨率影像中，生物量與紋理特征的相關性較強，且紋理特征不會出現植被指數的光飽和特性[24-27]。采用二階紋理算法中的灰度共生矩陣（Gray Level Co-occurrence Matrix，GLCM）提取紋理特征，提取均值（Mean，Mea）、相異性（Dissimilarity，Dis）、角二階矩（Second Moment，SM）、方差（Variance，Var）、對比度（Contrast，Con）、熵（Entropy，Ent）、同質性（Homogeneity，Hom）和相關性（Correlation，Cor）8種紋理特征變量，具體計算公式見文獻[28-29]，紋理窗口大小為3×3、5×5、7×7、9×9和11×11。

1.3.4 生物量反演模型構建

隨機抽取70%數據集采用PLS和RF方法建模，剩余30%數據集用于驗證。同一種變化類型中，PLS和RF模型輸入的建模因子相同。

PLS 是一種集合主成分分析、典型相關分析和多元線性回歸分析的數據統計方法，可有效解決多元回歸分析中變量的多重相關性和噪聲，在此基礎上具有預測功能并能消除參數結構不確定和模型無法辨識的影響[30]。

RF 是一種采用多棵決策樹構建的增強型分類（回歸器），當新的數據被帶入隨機森林后，全部決策樹生成分類或預測結果，隨機森林以這些結果的眾數或均值為數據輸出，可在不選取特征的情況下對高維數據進行處理。由于隨機選取樣本，在每次學習中，決策樹均會采用不同的訓練集進行訓練，所以隨機森林能一定程度避免過擬合現象[31]。本研究中，隨機森林模型算法的實現基于R語言平臺。

1.3.5 模型精度評價

選取R2、均方根誤差（Root Mean Squared Error，RMSE）和估測精度（P）為模型精度評價指標。R2越大，RMSE越小，模型擬合度越好。計算公式[19，30]為：

式中，yi為樣本實測值；i為模型預測值；為樣本實測值均值；n為樣本數。

2 結果與分析

2.1 最佳觀測窗口選擇

3 種模型中，球狀模型R2最大（0.889），RSS 最小（0.93×10-6），說明球狀模型擬合結果最佳；其最佳觀測窗口為5.2 m（表3）。球狀模型的變異程度為7.7%，表明高山松無人機多光譜影像空間尺度變異主要由結構性因素引起，其內在空間自相關性強。

表3 變異函數模型參數Tab.3 Parameters of variance function models

2.2 遙感因子相關性分析

以高山松地上生物量為因變量，提取的遙感因子為自變量；采用Pearson 相關系數作為評價指標，其絕對值在0～1 之間，越靠近1，相關性越強。為獲取相關性強的因子，將窗口依次設置為3×3、5×5、7×7、9×9和11×11；優選出每個窗口相關系數較強的前5個因子（表4）。高山松地上生物量模型中，紅邊與近紅外波段敏感性較強，相關系數為0.358～0.442，在提取窗口為9×9時，相關性較強。

表4 遙感因子與生物量的相關系數Tab.4 Correlation coefficients among remote sensing factors and biomass

在以上篩選的基礎上，將相關系數≥0.428 設為閾值，且在0.01水平顯著；5個相關性較高的遙感因子依次為RE9x9_Dis、RE9x9_Con、RE9x9_Var、NIR9x9_Dis和RE5x5_Con，相關系數分別為0.442、0.439、0.430、0.429和0.428。

2.3 飛行區高山松地上生物量模型選擇

基于RF 構建的生物量模型優于PLS，R2較大（0.90），RMSE 較小（17.96 t/hm2），P較大（84.98%）（表5，圖3）。選擇RF 模型進行飛行區高山松地上生物量估算和反演。

表5 模型精度對比Tab.5 Comparisons on model accuracies

圖3 預測值與實測值比較（a：RF；b：PLS）Fig.3 Comparisons between predicted values and measured values(a:RF;b:PLS)

2.4 飛行區高山松地上生物量反演

飛行區高山松地上生物量均值為130.48 t/hm2，總生物量為8 343.53 t（圖4）。反演結果顯示，飛行區高山松地上生物量空間分布呈中間高、四周低的趨勢，生物量集中在123.99～147.32 t/hm2，少數處于172.48～211.13 t/hm2，與實際情況相符，反演結果較好。

圖4 飛行區高山松地上生物量空間分布Fig.4 Spatial distributions of P.densata above-ground biomass in flight area

3 討論與結論

本研究采用無人機多光譜影像為遙感數據源，結合36塊樣地實測數據，首先通過地統計學變異函數確定高山松地上生物量最佳觀測窗口；結果顯示，變異函數3 種模型中，球狀模型擬合效果最佳，其最佳觀測窗口為5.2 m。其次，采用ArcGIS 軟件對無人機多光譜影像進行重采樣至5.2 m，提取200個紋理特征參數，計算其與高山松生物量的Pearson相關系數，并選取相關系數≥0.428 且在0.01 水平顯著的5 個紋理特征因子，建立PLS 和RF 地上生物量模型；結果顯示，RF 模型優于PLS 模型。基于RF 模型，無人機飛行區的高山松地上生物量均值為130.48 t/hm2，總生物量為8 343.53 t。本研究結果可為高寒山地森林生物量便捷、準確估測提供技術參考。本研究中，樣地生物量均值為130.48 t/hm2，超出謝福明[2]研究中香格里拉市高山松地上生物量（17.45～107.46 t/hm2），可能是因為本研究區處于林場，人工撫育強度大，保護措施較好，樹木長勢好。

因條件限制，本研究布設樣地為1 hm2，無人機多光譜影像分辨率為厘米級，將樣地公頃生物量推算至影像像元尺度生物量時，結果存在一定誤差。

確定最佳觀測尺度的方法主要有3 種，包括變異函數、局部方差和離散度。本研究僅探討了基于變異函數的最佳觀測尺度，對基于局部方差和離散度的最佳觀測尺度有待進一步探索。

利益沖突：所有作者聲明無利益沖突。

作者貢獻聲明：楊正道負責試驗調查與設計、數據收集與分析、論文撰寫和文獻檢索；舒清態負責論文指導與修改；黃金君負責試驗調查與設計和數據收集與分析；周文武負責試驗調查與設計；胥麗、王書偉負責文獻檢索；羅紹龍負責數據收集與分析。