劉浩然,范偉偉,徐永勝,林文樹
(東北林業(yè)大學(xué) 工程技術(shù)學(xué)院,黑龍江 哈爾濱 150040)
森林作為陸地生態(tài)系統(tǒng)中最大的碳儲庫在維持全球碳平衡過程中起著重要作用[1-2]。因此,如何快速準(zhǔn)確的估算森林生物量及碳儲量成為亟待解決的問題。目前林業(yè)研究所計算的生物量通常指森林地上生物量,即包括樹干、樹枝、樹葉等所有的地上活物質(zhì)[3]。然而,傳統(tǒng)測定森林地上生物量的方法耗時費(fèi)力、效率低下,還會對森林造成一定的損害。
近年來,激光雷達(dá)作為一種新興的主動遙感技術(shù)在林業(yè)中的應(yīng)用越來越廣泛。國內(nèi)外學(xué)者開展了大量關(guān)于基于激光雷達(dá)技術(shù)提取單木結(jié)構(gòu)參數(shù)的研究,為估算不同尺度下森林生物量提供了數(shù)據(jù)基礎(chǔ)。Gleason 等[4]基于單木樹冠采用支持向量機(jī)、隨機(jī)森林以及多元線性回歸等多種方法分別估計了森林單木生物量,結(jié)果表明支持向量機(jī)法估計的生物量模型精度最高,同時不同方法產(chǎn)生精度差異的原因可能是點云密度、模型差異及單木結(jié)構(gòu)參數(shù)提取精度的影響。Kankare 等[5-6]基于LiDAR 數(shù)據(jù)和單木組分生物量構(gòu)建了新的生物量模型,當(dāng)引入算法提取的單木樹高、冠幅等結(jié)構(gòu)參數(shù)作為輸入變量時,該模型取得了較高的估計精度,表明LiDAR 數(shù)據(jù)對于模型的構(gòu)建以及精度的提升具有重要作用。另外,作者又通過實驗證明了相比于地基激光雷達(dá),機(jī)載激光雷達(dá)在單木生物量預(yù)測方面精度略低,但應(yīng)用范圍更廣、適用性更高,同時利用樹冠結(jié)構(gòu)參數(shù)和幾何特征可對預(yù)測精度進(jìn)行提升。Usoltsev 等[7]引入虛擬變量構(gòu)建新的異速生長方程,以樺樹為例估計其單木和林分的生物量,經(jīng)試驗表明:該模型具有較高的實用性和移植性。同時,國內(nèi)利用LiDAR點云數(shù)據(jù)反演單木生物量的研究也越來越多。馮仲科等[8]將樹冠體積、表面積等樹冠因子引入生物量模型,獲得了相較于傳統(tǒng)CAR(Constant allometric ratio)模型精度更高的生物量模型。劉東起等[9]利用點云數(shù)據(jù)提取樹高、冠幅等結(jié)構(gòu)參數(shù),基于異速生長方程估測了小興安嶺涼水地區(qū)的森林生物量;劉峰等[10]利用條件隨機(jī)場和最優(yōu)化方法實現(xiàn)LiDAR 點云的單木分割,基于單木點云空間結(jié)構(gòu)、回波特征以及地形特征等作為特征變量,采用回歸模型估測喬木層單木地上生物量。李桂林等[11]以黃河三角洲孤島林場刺槐林為研究對象,利用無人機(jī)激光雷達(dá)數(shù)據(jù)從單木尺度提取刺槐的樹高和冠幅,結(jié)合背包移動雷達(dá)數(shù)據(jù)提取的單株刺槐胸徑,構(gòu)建了刺槐生物量估測模型,并對模型估算過程中存在的不確定性進(jìn)行了討論。
綜上所述,森林生物量與森林的垂直結(jié)構(gòu)有很緊密的相關(guān)關(guān)系,且易于被激光雷達(dá)獲取[12]。目前國內(nèi)外學(xué)者大多基于機(jī)載激光雷達(dá)或地基激光雷達(dá)獲取的單木結(jié)構(gòu)參數(shù)作為變量,通過一系列反演推算最終得到森林生物量。而無人機(jī)激光雷達(dá)作為一種新型獲取點云數(shù)據(jù)設(shè)備,具有體積小、操作性強(qiáng)、作業(yè)周期短、能快速獲取高精度點云數(shù)據(jù)等優(yōu)點[13]。本研究基于無人機(jī)激光雷達(dá)點云數(shù)據(jù),通過聚類算法獲取單木結(jié)構(gòu)參數(shù)作為變量,并基于多元回歸統(tǒng)計關(guān)系建立生物量與單木結(jié)構(gòu)參數(shù)之間的模型,探討模型的可行性以及多組變量與生物量之間的聯(lián)系,通過精度對比篩選出最佳因子估測森林單木生物量。
研究區(qū)位于哈爾濱市城市林業(yè)示范基地,地理坐 標(biāo) 為127°35′~127°39′E,45°42′~45°44′N,總占地面積為43.95 hm2。研究區(qū)海拔134~163 m,屬溫帶季風(fēng)氣候,夏季溫暖多雨,冬季寒冷多雪,年平均氣溫為3.5℃。實驗林場現(xiàn)有主要樹種為落葉松Larix gmelinii、核桃揪Juglans mandshurica、水曲柳Fraxinus mandshurica、白樺Betula platyphylla、樟子松Pinus sylvestris、蒙古櫟Quercus mongolica等[14]。
本研究在獲取研究區(qū)內(nèi)水曲柳和樟子松樣地點云數(shù)據(jù)時使用的無人機(jī)型號為大疆M600(深圳大疆創(chuàng)新科技有限公司),搭載得激光雷達(dá)系統(tǒng)為SZT-250。數(shù)據(jù)采集時間為2018年9月27日,天氣狀況良好,晴朗無風(fēng);平均飛行地速28.8 km/h,相對飛行高度50.2 m,飛行總時長6 min 56 s,掃描速度為10~100 線/秒,平均點密度為40 點/m2,共獲取西北-東南方向7 個條帶的數(shù)據(jù)。數(shù)據(jù)采用標(biāo)準(zhǔn)LiDAR 存儲格式(.las),根據(jù)需要可轉(zhuǎn)換為(.pcd)、(.txt)等格式。
水曲柳和樟子松樣地的野外調(diào)查于2018年10月3—5日期間開展。結(jié)合林分密度,分別在水曲柳和樟子松樣地內(nèi)設(shè)計并隨機(jī)選取了3 塊邊長15 m×15 m 的樣方進(jìn)行數(shù)據(jù)實測。采用皮尺測量單木東西、南北2個方向的冠幅并求算冠幅平均值,利用測高儀獲取單木樹高。
由于本研究所選林分的郁閉度較高,無人機(jī)激光雷達(dá)很難將激光點打到樹干,無法準(zhǔn)確提取單木胸徑。因此選擇樹高、冠幅、樹冠投影面積和樹冠體積作為反演生物量的單木結(jié)構(gòu)參數(shù)。
2.1.1 預(yù)處理
基于LiDAR360 軟件(北京數(shù)字綠土科技有限公司),對無人機(jī)點云數(shù)據(jù)進(jìn)行了3 步預(yù)處理:1)點云數(shù)據(jù)中常見的噪聲主要為孤立點。本研究采用基于距離的孤立點算法進(jìn)行去噪處理[15-16],以一定球形區(qū)域內(nèi)的點數(shù)少于既定值的點視為噪聲點。根據(jù)本研究樣地點云數(shù)據(jù)的分布情況,設(shè)置鄰域點個數(shù)為3,鄰域搜索半徑為5 m。2)采用漸進(jìn)加密三角網(wǎng)濾波算法[17]將點云劃分為地面點和非地面點。當(dāng)角度為8°,距離為1.4 m 時地面點分離效果最佳。3)將分離得到的地面點以0.5 m×0.5 m 為基本單元建立規(guī)則格網(wǎng),采用反距離權(quán)重法[18]進(jìn)行局部內(nèi)插處理,生成數(shù)字高程模型(Digital elevation model,DEM)和數(shù)字表面模型(Digital surface model,DSM)。由DSM 減去DEM 將得到冠層高度模型(canopy height model,CHM)[19]。
2.1.2 單木冠幅與樹高提取
k 均值聚類算法(k-means clustering algorithm)是一種通過迭代求解的聚類分析算法,首先隨機(jī)選擇K 個元素作為初始聚類中心,然后計算數(shù)據(jù)集所有元素與每個聚類中心之間的距離,并把每個元素分配給距離它最近的聚類中心。全部元素都被分配后,每個聚類的聚類中心會根據(jù)聚類中現(xiàn)有元素重新計算,直到目標(biāo)函數(shù)收斂為止[20]。因此,通過迭代把點云中的每一個點劃分到不同的聚類中,可以實現(xiàn)對單木樹冠的識別與分割。
在軟件Visual Studio 2017 中,利用k-means算法對冠層點云進(jìn)行分割,依據(jù)單木樹冠分割結(jié)果提取東西、南北冠幅,并求算平均冠幅值?;诰植孔畲笾捣ㄕ业絾文緲涔诘臉O大值點,該點的z值即為樹高。
2.1.3 單木樹冠投影面積及體積計算
二維凸包算法是基于離散的點云集構(gòu)建樹冠投影面積主要算法:樹冠投影面為一定范圍內(nèi)所有點的集合,二維凸包就是將點集中最邊緣的點連接起來構(gòu)建凸多邊形,根據(jù)凸多邊形的面積估算樹冠的投影面積[21]。而樹冠體積的求算方法則基于積分的思想,首先將樹冠分為若干層,對每一層實施二維凸包算法求算該層的橫截面積,利用分層高度和橫截面積計算各層體積,最后把所有的切片體積累加求和得到樹冠體積,該方法稱為漸進(jìn)凸包算法[22]。
非線性模型在森林單木生長模型中使用最為普遍。其中,相對生長模型(異速生長關(guān)系)作為最具代表性的非線性模型,大致分為CAR(Constant allometric ratio)和VAR(Variable allometric ratio)兩種類型。胥輝采用再抽樣的方法證明在相對生長模型中CAR 模型參數(shù)估計值更加準(zhǔn)確[23]。
本研究以CAR 模型為基礎(chǔ),引入冠幅、樹冠體積、樹冠投影面積構(gòu)建生物量模型,并與國家林業(yè)局公布水曲柳以及樟子松的異速生長方程估算得到的生物量(表1)進(jìn)行對比[24]。由于各樹冠因子均滿足相對生長關(guān)系,所以將CAR 模型表示為各自變量或組合變量冪函數(shù)相乘的形式[25]。
表1 國家林業(yè)局公布的水曲柳和樟子松的異速生長方程Table 1 Allometric equationof Fraxinus mandshurica and Pinus sylvestris published by State Forestry Administration
據(jù)此可以得到包含樹高、冠幅、樹冠體積、樹冠投影面積的生物量模型如下:
式(1)~(4)中:W為生物量;H為樹高;C為冠幅平均值;S為樹冠投影面積;V為樹冠體積;a1、a2、a3、a4為模型參數(shù)。
本研究使用決定系數(shù)(R2)、均方根誤差(RMSE)、總相對誤差(SRE)、平均相對誤差(MRE)、平均相對誤差絕對值(MREA)5 個統(tǒng)計指標(biāo)評估模型的擬合精度以及預(yù)測能力[26]。各指標(biāo)公式如下所示:
式(4)~(8)中:n為樣本單元數(shù);為模型的預(yù)測值;yi為實測值;為樣本平均值。
利用k-means 算法對樣地點云數(shù)據(jù)進(jìn)行聚類,對聚類后的單木樹冠進(jìn)行冠幅提取,依據(jù)研究區(qū)方位,把單木樹冠點云中x方向的最大值與最小值之差作為東西冠幅值,y方向的最大值與最小值之差作為南北冠幅值,基于東西和南北方向的冠幅數(shù)據(jù)計算平均冠幅值。分別用水曲柳樣地和樟子松樣地中單木平均冠幅實測數(shù)據(jù)與k-means 算法提取的平均冠幅值進(jìn)行回歸分析,并繪制散點圖。由圖1和圖2可知,基于算法提取的水曲柳樣地單木平均冠幅值在3.47~8.10 m 之間,均值為5.54 m;樟子松樣地單木平均冠幅值在2.27~5.78 m 之間,均值為4.11 m。從表2結(jié)果可以看出,k-means 算法提取平均冠幅值x與實測平均冠幅值y的存在較好的線性關(guān)系,水曲柳樣地中平均冠幅的R2為0.82,RMSE 均值為0.61 m;樟子松樣地中平均冠幅的R2分別為0.80,RMSE 均值為0.43 m。結(jié)果表明:基于k-means 算法提取的單木冠幅精度具有較高的可靠性。
表2 提取與實測平均冠幅的回歸關(guān)系Table 2 Regression relationship between extractedand measured mean crown width
圖1 水曲柳樣地中提取冠幅平均值與實測冠幅平均值的回歸關(guān)系Fig.1 Regression relationship between extracted and measured average crown width in the Fraxinus mandshurica plot
圖2 樟子松樣地中提取冠幅平均值與實測冠幅平均值的回歸關(guān)系Fig.2 Regression relationship between extracted and measured average crown width in the Pinus sylvestris plot
基于分割后的點云數(shù)據(jù),利用局部最大值法探測樹冠頂點,該點的Z值即為單木的樹高。分別用水曲柳樣地和樟子松樣地中樹高實測數(shù)據(jù)與算法提取樹高值進(jìn)行回歸分析,并繪制提取值與實測值的散點圖。由圖3~4所示,基于算法得提取的水曲柳單木樹高在17.97~23.27 m 之間,均值為20.78 m;樟子松單木樹高在15.54~19.12 m 之間,均值為17.25 m。從表3所示的提取樹高和實測樹高回歸關(guān)系結(jié)果可以得出:水曲柳樣地中樹高的R2為0.86,RMSE 為0.49 m;樟子松樣地中樹高的R2為0.84,RMSE 為0.40 m。結(jié)果表明算法提取的單木樹高與實測值之間具有較好的相關(guān)性。
圖3 水曲柳樣地中提取樹高與實測樹高的回歸關(guān)系Fig.3 Regression relationship between extracted and measured tree height in the Fraxinus mandshurica plot
圖4 樟子松樣地中提取樹高與實測樹高的回歸關(guān)系Fig.4 Regression relationship between extracted and measured tree height in the Pinus sylvestris plot
表3 提取與實測樹高回歸關(guān)系Table 3 Regression relationship between extractedand measured tree height
對分割后的點云數(shù)據(jù)采用二維凸包算法分別計算兩塊樣地的46 株水曲柳和43 株樟子松的樹冠投影面積,并基于漸進(jìn)凸包算法分別計算水曲柳和樟子松兩塊樣地的樹冠體積,得到的水曲柳和樟子松樹冠投影面積和體積數(shù)據(jù)統(tǒng)計分別如表4~5所示:
表4 兩塊樣地樹冠投影面積統(tǒng)計結(jié)果Table 4 Statistical results of tree crown projection area in two plots
表5 兩塊樣地樹冠體積統(tǒng)計結(jié)果Table 5 Statistical results of tree crown volume in two plots
基于國家林業(yè)局公布的水曲柳和樟子松異速生長方程計算樣地生物量,結(jié)合無人機(jī)激光雷達(dá)提取的單木結(jié)構(gòu)參數(shù),根據(jù)生物量模型1~4,采用非線性最小二乘回歸方法分別對生物量進(jìn)行擬合,建立生物量估算模型,得出模型參數(shù)估計值以及評價指標(biāo),結(jié)果如表6所示。
由表6可知,當(dāng)生物量模型中引入樹冠投影面積、樹冠體積后,相關(guān)系數(shù)R2增大,同時均方根誤差RMSE 減小。結(jié)果表明樹冠投影面積、樹冠體積和樹高、冠幅一樣,與生物量之間具有一定的相關(guān)性。模型4 是兩塊樣地生物量擬合最優(yōu)的,其中,水曲柳樣地模型擬合優(yōu)度R2=0.834,RMSE=18.912 1 kg/株;樟子松樣地模型擬合優(yōu)度R2=0.792,RMSE=8.120 kg/株,與國家林業(yè)局頒布的兩種樹木計算生物量的異速生長方程具有較高的相關(guān)性。綜上所述,通過對擬合優(yōu)度總結(jié)分析,得到兩塊樣地生物量最優(yōu)模型:
表6 生物量模型參數(shù)估計值與擬合優(yōu)度Table 6 Parameters estimation and goodness of fit of biomass models
水曲柳生物量最優(yōu)模型為W=0.14H1.546C0.385S0.247V0.426;樟子松生物量最優(yōu)模型為W=0.193H1.912C0.329S-0.18V0.104。
由于在模型中引入了樹冠投影面積和樹冠體積變量,使得生物量模型預(yù)估精度比僅有樹高、冠幅作為變量的生物量模型有明顯提高。由表7中生物量模型評價指標(biāo)可以看出,兩塊樣地生物量模型的SRE 和MRE 分別在1.5%和0.7%以內(nèi),而MARE 在6.95%~12.59%之間;當(dāng)選用樹高、冠幅、樹冠投影面積和體積作為變量時,兩塊樣地生物量模型的各項評價指標(biāo)均最小。
表7 生物量模型評價指標(biāo)Table 7 Biomass models evaluation index
總之,模型4 在利用最小二乘法擬合生物量時的擬合程度最好、預(yù)測精度最高,表明無人機(jī)激光雷達(dá)數(shù)據(jù)可以通過提取樹高、冠幅、樹冠投影面積及樹冠體積等單木結(jié)構(gòu)參數(shù)對研究區(qū)樣地單木生物量進(jìn)行反演,能取得較高的反演精度。
本研究利用無人機(jī)激光雷達(dá)對水曲柳和樟子松兩塊樣地進(jìn)行掃描,利用編制算法及數(shù)據(jù)處理軟件提取單木結(jié)構(gòu)參數(shù)?;谔崛〉膯文窘Y(jié)構(gòu)參數(shù),建立包含樹高、冠幅、樹冠投影面積及樹冠體積的新生物量模型,并將新模型擬合結(jié)果與國家林業(yè)局頒布水曲柳、樟子松生物量的異速生長方程估算的生物量進(jìn)行對比。本研究具體結(jié)論如下:
1)將算法提取的單木冠幅、樹高與實測值進(jìn)行對比,其中水曲柳樣地平均冠幅在3.47~8.10 m之間,樹高均值為20.78 m;樟子松樣地平均冠幅在2.27~5.78 m 之間,樹高均值為17.25 m,這與研究區(qū)內(nèi)樣地樹木實際情況相符,且單木因子均與實測數(shù)據(jù)相關(guān)性良好,故本研究基于無人機(jī)激光雷達(dá)測得的單木參數(shù)可作為變量應(yīng)用在生物量模型中。
2)當(dāng)選用樹高、冠幅、樹冠投影面積和樹冠體積作為變量時,水曲柳和樟子松生物量模型擬合效果最優(yōu),R2分別為0.834、0.792,RMSE 分別為18.912、8.120 kg/株,且各項評價指標(biāo)均較優(yōu),符合上述擬合結(jié)果。
3)當(dāng)生物量模型中加入樹冠因子時,模型擬合效果和評價指標(biāo)均較優(yōu),尤其當(dāng)加入樹冠體積時更為明顯:對比樹高、冠幅作為自變量的模型,樹冠體積的引入使得水曲柳、樟子松生物量模型的R2分別提高了0.076、0.060,RMSE 分別下降了6.759、1.386 kg/株;同理,對比樹高、冠幅和樹冠投影面積作為自變量的模型,加入樹冠體積則使得兩種生物量模型R2分別提高了0.082、0.060,RMSE 分別下降了5.329、1.085 kg/株。因此,通過分析結(jié)果發(fā)現(xiàn)將樹冠因子引入生物量模型是提高森林單木生物量估計精度的一種可行方法。
1)本研究基于無人機(jī)激光雷達(dá)數(shù)據(jù)提取的樣地單木樹高、冠幅、樹冠投影面積及樹冠體積4種結(jié)構(gòu)參數(shù)擬合單木生物量模型,通過與國家林業(yè)局頒布的水曲柳和樟子松異速生長方程估算的生物量進(jìn)行對比,證明了在沒有樹木胸徑數(shù)據(jù)作為參數(shù)進(jìn)行生物量模型構(gòu)建,依然能夠取得較好的擬合優(yōu)度和較高的預(yù)測精度,結(jié)合無人機(jī)激光雷達(dá)技術(shù)的優(yōu)勢,表明無人機(jī)激光雷達(dá)數(shù)據(jù)具有獨(dú)立、高效、高精度反演生物量的潛力。龐勇、李增元等人利用機(jī)載激光雷達(dá)獲取了研究區(qū)內(nèi)針葉林、闊葉林等多種森林類型的點云數(shù)據(jù),結(jié)果表明:各樣地中除根生物量外各組分生物量模型的R2均高于0.8。這與本研究得到的結(jié)果基本一致,表明了基于無人機(jī)激光雷達(dá)估測生物量的可行性。
2)無人機(jī)激光雷達(dá)數(shù)據(jù)在聚類過程中存在“欠分割”或“過分割”的現(xiàn)象,是因為樣地林間距過小導(dǎo)致樹木存在樹冠重疊現(xiàn)象,因此難以正確地分割出重疊樹冠的邊界,造成冠幅數(shù)據(jù)提取誤差增大;由于樣地林分的高郁閉度,激光雷達(dá)打在地面上的點減少,會造成樹高提取值的誤差增大,選擇更為適合的算法提取單木結(jié)構(gòu)參數(shù)能更好的反演生物量。
3)本研究構(gòu)建的4 個生物量模型,通過對比擬合程度和預(yù)測精度可以得出:引入樹冠投影面積和樹冠體積2 個樹冠因子能明顯提高生物量模型的質(zhì)量,其中樹冠體積對生物量模型的貢獻(xiàn)度更為顯著,這是因為樹冠部分的生物量受樹冠形狀、大小、飽滿程度以及樹木長勢影響很大,與樹冠投影面積相比樹冠體積能更好的反應(yīng)單木樹冠的生長態(tài)勢,因此和生物量是緊密相關(guān)的。
4)選取水曲柳和樟子松2 塊樣地作為研究對象,其中水曲柳屬于闊葉木,樟子松屬于針葉木。通過數(shù)據(jù)分析對比發(fā)現(xiàn):水曲柳生物量擬合模型的更高,但樟子松生物量模型的RMSE 明顯低于水曲柳;在評價指標(biāo)方面,樟子松生物量模型的各項指標(biāo)均優(yōu)于水曲柳生物量模型。分析原因主要是闊葉木樹冠冠幅、體積相對較大,在單木生物量上占比更多,因此引入樹冠因子時,相比于樟子松生物量模型,水曲柳生物量模型R2提升更大。而闊葉木冠型多變,針葉木樹冠結(jié)構(gòu)單一、冠型整齊,造成了水曲柳樣地單木結(jié)構(gòu)參數(shù)獲取精度的降低,因此產(chǎn)生較大的誤差值,符合上述統(tǒng)計結(jié)論。
綜上所述,本研究基于無人機(jī)激光雷達(dá)數(shù)據(jù)引入樹冠因子作為生物量建模的變量,既保證了模型精度又比傳統(tǒng)的通過人工測量單木結(jié)構(gòu)參數(shù)來估測單木生物量成本小。需要注意的是,在進(jìn)行單木樹冠點云分割時可能會產(chǎn)生漏檢、過檢等錯誤分割現(xiàn)象,因此在今后的研究中需要提高分割算法的精度和自動化程度。另外,本研究選擇的樣地是人工林,樣本數(shù)據(jù)相對較少,在后續(xù)研究中應(yīng)需加大樣本數(shù)據(jù)的容量,同時選取結(jié)構(gòu)更加復(fù)雜、樹種更為豐富的天然林林分進(jìn)行生物量估測研究;最后,聯(lián)合多源遙感技術(shù)手段協(xié)同作業(yè)將成為獲取森林高精度單木結(jié)構(gòu)參數(shù)以及優(yōu)化生物量模型的關(guān)鍵。