基于遙感數(shù)據(jù)與森林詳查估算林齡的空間分布

李凡，李明澤，史澤林，姜宏宇，安金鵬

0 引言

林齡結(jié)構(gòu)是分析干擾對(duì)森林碳循環(huán)影響的重要替代因子，同時(shí)也是評(píng)估區(qū)域尺度森林生態(tài)系統(tǒng)碳匯潛力的重要參數(shù)［1-3］。林齡不僅涉及森林生態(tài)系統(tǒng)不同生長(zhǎng)階段的碳源/碳匯功能，同時(shí)也反映了森林受到的自然干擾和森林管理等影響［4-5］。研究林齡的空間分布可以為林業(yè)、企業(yè)、事業(yè)單位制定經(jīng)營(yíng)計(jì)劃提供依據(jù)；為生態(tài)方面碳循環(huán)的研究提供幫助；為全面開(kāi)展以碳匯為目的的森林經(jīng)營(yíng)打好基礎(chǔ)；同時(shí)還可以協(xié)助檢查、評(píng)價(jià)森林經(jīng)營(yíng)效果和計(jì)劃執(zhí)行情況。

在已有的研究中，森林詳查和遙感是現(xiàn)今獲取林齡數(shù)據(jù)的主要手段。森林詳查通過(guò)有限的樣地?cái)?shù)據(jù)獲取較小尺度的林齡數(shù)據(jù)，很難獲取大尺度區(qū)域上森林林齡空間分布。遙感具有較高的空間異質(zhì)性監(jiān)測(cè)能力，是監(jiān)測(cè)區(qū)域尺度上森林變化信息的有效技術(shù)。戴銘等利用森林詳查數(shù)據(jù)結(jié)合降尺度技術(shù)估算林齡的空間分布［1］；Pan Y等利用森林清查數(shù)據(jù)、歷史火災(zāi)數(shù)據(jù)，結(jié)合衛(wèi)星遙感數(shù)據(jù)估算并繪制了北美森林林齡空間分布圖［6］；Zhang C等利用森林遙感高度和類型信息，繪制了分辨率為1 km的全國(guó)森林林齡分布圖［7］。目前對(duì)森林林齡的研究幾乎都是結(jié)合森林詳查數(shù)據(jù)與遙感數(shù)據(jù)，但是鮮有通過(guò)直接構(gòu)建林齡與其他林分因子的統(tǒng)計(jì)模型并結(jié)合多源遙感數(shù)據(jù)來(lái)估測(cè)森林林齡分布的研究。

本研究以森林清查數(shù)據(jù)、資源三號(hào)影像數(shù)據(jù)以及激光雷達(dá)數(shù)據(jù)為數(shù)據(jù)源，通過(guò)構(gòu)建林齡與樹(shù)高的模型關(guān)系，利用遙感數(shù)據(jù)獲取區(qū)域尺度林分平均高和類型信息，通過(guò)模型擬合估算出區(qū)域尺度的林齡空間分布，并探討其空間分布特征。

1 研究地區(qū)與研究方法

1.1 研究區(qū)概況

研究地區(qū)為黑龍江省東北林業(yè)大學(xué)帽兒山實(shí)驗(yàn)林場(chǎng)中林施業(yè)區(qū)。該林場(chǎng)地理坐標(biāo)為東經(jīng)127°30′～ 127°34′,北緯 45°20′～ 45°25′。總面積為26 496 hm2,南北長(zhǎng)30 km,東西寬20 km,全場(chǎng)劃分為10個(gè)施業(yè)區(qū),即尖砬溝、老爺嶺、新墾、躍進(jìn)、中林、北林、東林、太平、密峰和老山施業(yè)區(qū)。植被屬長(zhǎng)白植物區(qū)系，原始地帶性植被為闊葉紅松林，現(xiàn)為東北東部山區(qū)較典型的天然次生林區(qū)［8］，次生林類型多樣且具有代表性，群落類型有硬闊葉林、軟闊葉林、針葉林和針闊混交林，平均林齡為50 a。研究區(qū)主要闊葉樹(shù)種有水曲柳、胡桃楸、椴樹(shù)、榆樹(shù)、色木、柞樹(shù)、白樺、楊樹(shù)和楓樺等，主要針葉樹(shù)種為落葉松等。本次采用的是2014年國(guó)家森林資源二類調(diào)查的帽兒山林場(chǎng)中林施業(yè)區(qū)小班調(diào)查數(shù)據(jù)。

1.2 數(shù)據(jù)處理與研究方法

1.2.1 數(shù)據(jù)處理

選取資源三號(hào)衛(wèi)星影像面積為2 500 km2單景全色波段影像以及面積為2 704 km2單景多光譜影像，獲取時(shí)間為2014年9月。

根據(jù)影像數(shù)據(jù)特點(diǎn)對(duì)其做預(yù)處理，依次為：依靠RPC文件和DEM數(shù)據(jù)進(jìn)行圖像正射校正；使用FLAASH大氣校正方法對(duì)多光譜影像進(jìn)行大氣校正；使用image to image方式進(jìn)行圖像配準(zhǔn)，以全色波段影像為基準(zhǔn)對(duì)多光譜圖像進(jìn)行配準(zhǔn)，共生成匹配點(diǎn)76個(gè)，RMSE為0.237；采用適用于資源三號(hào)數(shù)據(jù)的Gram-Schmidt Pan Sharpening(GS)融合方法進(jìn)行圖像融合。

1.2.2 平均樹(shù)高提取方法

激光雷達(dá)可以提供直接、準(zhǔn)確的森林垂直結(jié)構(gòu)測(cè)量數(shù)據(jù)。利用機(jī)載激光雷達(dá)點(diǎn)云數(shù)據(jù)可以獲取研究區(qū)域范圍內(nèi)的樹(shù)高值［9］。研究結(jié)果表明，雷達(dá)數(shù)據(jù)獲取樹(shù)高可以達(dá)到實(shí)地調(diào)查的精度，其誤差絕對(duì)平均值均在1.4 m之內(nèi)［10］。

提取樹(shù)高方法往往與DEM和DSM的生成緊密相關(guān),兩者差分而生成的高度即地物絕對(duì)高度,是樹(shù)高提取的直接依據(jù)［11］。冠層高度模型（Canopy Height Model，簡(jiǎn)稱CHM）表達(dá)了林木高度和冠層的空間分布狀態(tài)［12］，數(shù)字表面模型（Digital Surface Model，簡(jiǎn)稱DSM）是對(duì)地物表面的模擬，包括建筑和植被等的表面。由地面點(diǎn)內(nèi)插生成DEM，非地面點(diǎn)生成DSM。在森林區(qū)域由DSM和DEM相減就能得到CHM，并可用于樹(shù)高的提取［13-16］。

1.2.3 森林類型提取方法

ZY-3比較適合于森林植被信息的提取，特別是針葉林和闊葉林［17］。由于遙感圖像大都是多波段圖像，在提取信息特征值時(shí)，就會(huì)遇到波段的選取問(wèn)題。如每個(gè)波段都計(jì)算其特征值并應(yīng)用于分類，會(huì)導(dǎo)致信息冗余，反而使分類精度下降，而且時(shí)間開(kāi)銷比較大?？梢圆捎脤⒁恍╆P(guān)鍵波段融合成一張含有多波段信息的合成圖像并從圖像中尋找特征值的方法［18］。本研究合成圖像由資源三號(hào)影像數(shù)據(jù)四個(gè)波段、ISODATA非監(jiān)督分類結(jié)果、NDVI和DEM組成。其中ISODATA非監(jiān)督分類分類類別數(shù)量設(shè)置為10～20，迭代次數(shù)為10，并將分類結(jié)果通過(guò)對(duì)比遙感影像合并類別。NDVI（歸一化植被指數(shù)）可以很好的分離出林地和非林地。計(jì)算公式為：

NDVI值的范圍為-1～1，一般綠色植被區(qū)范圍為0.2～0.8。

2 結(jié)果和分析

2.1 樹(shù)高與林齡相關(guān)性分析

2014年國(guó)家森林資源二類調(diào)查的帽兒山林場(chǎng)中林施業(yè)區(qū)219個(gè)小班數(shù)據(jù)，提取每個(gè)小班平均林齡和平均樹(shù)高數(shù)據(jù)做林齡-樹(shù)高的相關(guān)性分析。經(jīng)分析，林齡和樹(shù)高的相關(guān)系數(shù)為0.852，可以看出，林齡-樹(shù)高的相關(guān)性較強(qiáng)（圖1）。

圖1 林齡-樹(shù)高相關(guān)性圖Fig.1 Correlation between stand age and tree height

2.2 樹(shù)高與林齡模型擬合

通過(guò)統(tǒng)計(jì)分析結(jié)合資料去除了其中的45個(gè)異常點(diǎn)，最終用171個(gè)小班數(shù)據(jù)擬合模型。將數(shù)據(jù)按照森林類型分為針葉混交林、闊葉混交林和針闊混交林三類，分別擬合模型。其中，由于針闊混交林的樣本量較小，采用交叉驗(yàn)證的方法來(lái)對(duì)數(shù)據(jù)進(jìn)行擬合，擇優(yōu)選擇，保證精度。交叉驗(yàn)證是沒(méi)有任何前提假定直接估計(jì)泛化誤差的模型選擇方法，由于沒(méi)有任何假定，可以應(yīng)用于各種模型選擇中［19］。具體做法是將數(shù)據(jù)隨機(jī)等分為四份，選其中三份作為訓(xùn)練樣本進(jìn)行模型擬合，剩下一份作為檢驗(yàn)樣本進(jìn)行精度檢驗(yàn)，之后選訓(xùn)練樣本中的一份數(shù)據(jù)作為檢驗(yàn)樣本，其他三份數(shù)據(jù)作為訓(xùn)練樣本再次進(jìn)行模型擬合。這個(gè)過(guò)程一直進(jìn)行，直到所有的樣本都有且僅有一次作為檢驗(yàn)樣本。最后，對(duì)比四個(gè)模型的精度，選擇R2最大的模型。

針葉混交林，闊葉混交林，針闊混交林?jǐn)M合模型的R2分別為0.83、0.81、0.90，均接近于1，模型精度分別為85.1%、91.2%、80.3%。說(shuō)明林齡與樹(shù)高模型擬合精度較高，可以通過(guò)樹(shù)高較高精度的估算林齡。擬合模型見(jiàn)表1。

表1 林齡-樹(shù)高模型擬合表Tab.1 Model of stand age and tree height

2.3 提取小班平均樹(shù)高

進(jìn)行點(diǎn)云數(shù)據(jù)的濾波處理，提取地面點(diǎn)，分離非地面點(diǎn)。利用地面點(diǎn)數(shù)據(jù)生成DEM，非地面點(diǎn)數(shù)據(jù)生成DSM。生成CHM:使用ARCGIS柵格計(jì)算器，用DSM減去DEM得到研究區(qū)域的CHM（圖2）。利用CHM提取研究區(qū)域中的單木樹(shù)高，然后計(jì)算小班的平均高。

單木樹(shù)高的提取使用的是數(shù)字綠土公司的LiDAR360軟件，主要步驟為：使用CHM Segmentation工具，設(shè)置樹(shù)高閾值為5～40。處理 得 到 包 含 TreeLocationX，TreeLocationY、TreeHeight、CrownDiameter和 CrownArea 5個(gè)字段的.csv格式的表格文件。將表格文件加載到ArcGIS中，添加X(jué)Y數(shù)據(jù)，X字段為TreeLocationX，Y字段為TreeLocationY，Z字段為TreeHeight，導(dǎo)出為shp文件。然后根據(jù)小班來(lái)剪裁生成的shp文件，再根據(jù)屬性表計(jì)算出每個(gè)小班的平均樹(shù)高。利用統(tǒng)計(jì)后的53個(gè)小班的實(shí)測(cè)平均樹(shù)高對(duì)估算的平均樹(shù)高進(jìn)行精度檢驗(yàn)，結(jié)果：RMSE為3.69，精度為80.46%。

圖2 中林施業(yè)區(qū)冠層高度模型Fig.2 Canopy height mode of research area

2.4 基于資源三號(hào)數(shù)據(jù)提取森林類型

2.4.1 監(jiān)督分類

根據(jù)衛(wèi)星影像的可判性，結(jié)合森林二類清查資料將研究區(qū)劃分為6個(gè)類別：闊葉混交林、針葉混交林、道路、水體、建筑物、耕地或裸地，其中研究區(qū)域針闊混交林極少，并不單獨(dú)劃分一類。生成一幅多元數(shù)據(jù)集的合成圖像。

分別應(yīng)用最小距離法、最大似然法和神經(jīng)網(wǎng)絡(luò)法進(jìn)行監(jiān)督分類。通過(guò)在圖像上疊加分類結(jié)果，目視判斷出最小距離法不能很好的區(qū)分道路、建筑物、水體與耕地，但對(duì)林地與非林地的區(qū)分效果好。最大似然法能很好的區(qū)分建筑物、水體、耕地和林地，但是對(duì)林地與云、陰影以及林地與道路的區(qū)分效果不是很理想。神經(jīng)網(wǎng)絡(luò)法運(yùn)算速度慢，對(duì)林地與非林地的分離較為理想，但對(duì)林分類型的區(qū)分效果差。綜合考慮土地類型的提取與林分類型的區(qū)分效果，最終選擇采用最大似然法的分類結(jié)果，如圖3所示。

圖3 土地利用分類結(jié)果Fig.3 Land use classification results

2.4.2 分類精度評(píng)價(jià)

對(duì)最大似然法分類結(jié)果精度進(jìn)行評(píng)價(jià)，使用的是根據(jù)地面實(shí)際ROIS來(lái)求混淆矩陣的方法。地面實(shí)際ROIS是根據(jù)2016年5月份的高分辨率航拍相片劃定的。該相片空間分辨率為0.2 m，地物信息清晰可見(jiàn)。

評(píng)價(jià)結(jié)果為：林地與非林地能得到很好的分離，針葉混交林與闊葉混交林也可以較高精度的分離，針葉混交林生產(chǎn)者精度不高，主要原因是針葉混交林分類結(jié)果存在“胡椒鹽”現(xiàn)象。通過(guò)對(duì)像元的手動(dòng)歸類消除“胡椒鹽”。在本次分類中，很好的提取出林地，且對(duì)林地類型的分類精度也較為理想。評(píng)價(jià)檢驗(yàn)見(jiàn)表2。

表2 研究區(qū)土地利用分類結(jié)果Tab.2 Land use classification results in the study area

2.4.3 估算林齡

將機(jī)載雷達(dá)數(shù)據(jù)提取的平均樹(shù)高按照林班進(jìn)行劃分，計(jì)算出每個(gè)林班的平均樹(shù)高。共獲得18個(gè)林班的平均樹(shù)高數(shù)據(jù)，根據(jù)森林分類信息將其分別帶入統(tǒng)計(jì)模型，得到林齡的估算值。并將其劃分成齡組。估計(jì)結(jié)果見(jiàn)表4。

利用實(shí)測(cè)林齡數(shù)據(jù)對(duì)估測(cè)的林齡進(jìn)行精度檢驗(yàn)，檢驗(yàn)結(jié)果：平均相對(duì)誤差的絕對(duì)值為0.125，估測(cè)精度為85.4%。在研究區(qū)的矢量圖上添加估測(cè)的齡組信息，顯示研究區(qū)林齡的空間分布，如圖4所示。

圖4 中林施業(yè)區(qū)齡組空間分布Fig.4 Spatial distribution of age group of research area

表3 林齡估計(jì)概況Tab.3 Introduction of estimates stand age

3 討論

中林施業(yè)區(qū)林齡空間分布差異較小。森林平均林齡為51 a左右，從整體上來(lái)看林齡普遍偏低?？拷缆芬约案嘏缘臉?shù)種多為落葉松人工林，其林齡結(jié)構(gòu)偏年輕，這是人為活動(dòng)干擾較大造成的。

本研究通過(guò)有限的森林詳查數(shù)據(jù)，構(gòu)建了林齡與樹(shù)高的統(tǒng)計(jì)模型。通過(guò)機(jī)載雷達(dá)數(shù)據(jù)和資源三號(hào)衛(wèi)星影像在區(qū)域尺度上提取出森林的平均樹(shù)高數(shù)據(jù)以及森林的類型信息。通過(guò)遙感獲取的森林?jǐn)?shù)據(jù)應(yīng)用統(tǒng)計(jì)模型獲得區(qū)域尺度上的林齡空間分布。事實(shí)證明這種研究方法是可行的。當(dāng)林分平均高到達(dá)20 m及以上時(shí)，林齡與樹(shù)高關(guān)系出現(xiàn)飽和。林齡增加樹(shù)高不再發(fā)生變化。通過(guò)計(jì)算，估測(cè)林齡的一元線性模型中：針葉混交林，闊葉混交林，針闊混交林對(duì)應(yīng)模型適用林齡范圍分別為：50、70、100 a。

本研究存在的不足之處：未采用多元線性以及非線性模型，今后工作中應(yīng)該嘗試構(gòu)建林齡與多個(gè)森林調(diào)查因子的統(tǒng)計(jì)關(guān)系，使模型的精確度更高；在提取樹(shù)高的過(guò)程中，今后應(yīng)該進(jìn)一步提高CHM的提取精度；在遙感分類方面，可以加入更多輔助信息作為評(píng)判參數(shù)，或者采用決策樹(shù)，面向?qū)ο蠓诸惖确椒?，提高森林類型劃分的精度?/p>

［1］戴銘,周濤,楊玲玲,等.基于森林詳查與遙感數(shù)據(jù)降尺度技術(shù)估算中國(guó)林齡的空間分布［J］. 地理研究,2011, 30(1):172-184.

［2］Deng F, Chen J M, Pan Y, et al. The use of forest stand age information in an atmospheric CO2 in-version applied to North America［J］. Biogeosciences, 2013, 10(8), 5335-5348.

［3］Christiansen J R, Gundersen P. Stand age and tree species affect N2O and CH4 exchange from afforested soils［J］.Biogeosciences, 2011, 8(9), 2535-2546.

［4］周蕾,王紹強(qiáng),周濤,等. 1901～2010年中國(guó)森林碳收支動(dòng)態(tài):林齡的重要性［J］.科學(xué)通報(bào),2016,61(18):2064-2074.

［5］徐冰,郭兆迪,樸世龍,等. 2000～2050年中國(guó)森林生物量碳庫(kù):基于生物量密度與林齡關(guān)系的預(yù)測(cè)［J］. 中國(guó)科學(xué):生命科學(xué),2010,40(7):587-594.

［6］Pan Y, Chen J M, Birdsey R, et al. Age structure and disturbance legacy of North American forests［J］. Biogeosciences, 2011, 8(3),715-732.

［7］Zhang C, Ju W, Chen J M, et al. Mapping forest stand age in China using remotely sensed forest height and observation data［J］.Journal of Geophysical Research Biogeosciences, 2014, 119(6):1163-1179.

［8］豐緒霞,劉兆剛,張海玉,等. 基于RS和GIS帽兒山林場(chǎng)森林立地分類及質(zhì)量評(píng)價(jià)［J］. 東北林業(yè)大學(xué)學(xué)報(bào),2010,38(8):27-30+56.

［9］吳迪,范文義.激光雷達(dá)協(xié)同多角度光學(xué)遙感數(shù)據(jù)反演樹(shù)高［J］.北京林業(yè)大學(xué)學(xué)報(bào),2014,36(4):8-15.

［10］劉東起,范文義,李明澤.利用小光斑激光雷達(dá)估測(cè)林分參數(shù)和生物量［J］. 東北林業(yè)大學(xué)學(xué)報(bào),2012,40(1):39-43.

［11］龐勇,趙峰,李增元,等. 機(jī)載激光雷達(dá)平均樹(shù)高提取研究［J］. 遙感學(xué)報(bào),2008(1):152-158.

［12］邢艷秋,尤號(hào)田,霍達(dá),等. 小光斑激光雷達(dá)數(shù)據(jù)估測(cè)森林樹(shù)高研究進(jìn)展［J］. 世界林業(yè)研究,2014,27(2):29-34.

［13］范鳳云,陳爾學(xué),劉清旺,等. 基于低密度機(jī)載LiDAR和CCD數(shù)據(jù)的林分平均高提?。跩］.林業(yè)科學(xué)研究,2010,23(2):151-156.

［14］Zhen Z, Quackenbush L J, Zhang L. Impact of tree-oriented growth order in marker controlled region growing for individual tree crown delineation using airborne laser scanner(ALS) data［J］.Remote Sensing, 2014, 6(1): 555-579.

［15］Zhen Z, Quackenbush L J, Stehman S V, et al. Agent-based region growing for individual tree crown delineation from airborne laser scanning ALS data［J］. International Journal of Remote Sensing, 2015, 36(7): 1965-1993.

［16］Zhen Z, Quackenbush L, Zhang L. Trends in automatic individual tree crown detection and delineation evolution of LiDAR data［J］.Remote Sensing, 2016, 8(4): 333.

［17］張雨,林輝. 資源3號(hào)衛(wèi)星影像的植被信息提?。跩］. 中南林業(yè)科技大學(xué)學(xué)報(bào),2014,34(11):130-134.

［18］陳敏. 基于紋理特征的遙感圖像分類技術(shù)研究［D］.福州：福州大學(xué),2005.

［19］范永東. 模型選擇中的交叉驗(yàn)證方法綜述［D］.太原：山西大學(xué),2013.