基于多源遙感數據提高山地森林識別精度
——以祁連山國家公園肅南縣段為例

宋潔，劉學錄*

（1.甘肅農業大學資源環境學院，甘肅 蘭州730070；2.甘肅農業大學土地利用研究所，甘肅 蘭州730070）

森林作為陸地生物圈的主體，因其在全球碳循環中占據重要地位，并能夠為人類社會提供如生物多樣性保護、氣候調節和水資源涵養等諸多生態系統服務而受到廣泛的重視［1］。森林植被分布監測，不僅為區域乃至全球的環境變化研究提供了基礎信息［2］，也為林業規劃和生態保護政策的制定提供了基本的決策依據。在此背景下，探索低成本、高時效、操作方便并具有一定精度保證的森林識別方法對森林的監測與保護具有重要的意義。

近年來，遙感技術已成為大尺度森林識別的首選［3-4］。被動光學遙感數據包括高光譜［5-6］以及多光譜遙感影像［7-8］在目前森林信息獲取研究中的應用較為普遍。然而，多光譜和高光譜數據都只能提供物體表面的光譜信息，可能無法分離出光譜相似但結構不同的物體，例如灌木與喬木［9］。激光雷達是一種基于可見光、近紅外或短波紅外激光束測距的主動遙感技術，它可以通過物體返回脈沖所記錄的三維位置提供關于物體的垂直結構信息［10］。也正是基于此特點，激光雷達數據已逐步被單獨、或與光譜數據結合使用于土地覆被分類的研究中：Sasaki 等［11］結合機載激光雷達數據與高分辨率影像，采用面向對象的分類方法獲得97.5%的土地覆被分類精度以及48.4%的樹種分類精度；Zhou 等［12］通過星載激光雷達ICESat/GLAS 數據的波形特征將研究區域分類為樹木、建筑物和平坦地面并獲得87.2%的總體精度；Liu 等［13］結合星載激光雷達ICESat/GLAS 數據與Landsat 數據，將研究區分類為農田、森林、灌木、水體和不透水表面，總體精度達到91%。在分類過程中，以上研究均發現與單獨使用激光雷達數據進行分類相比，將主動式激光雷達數據與被動式光學遙感數據相結合，綜合依據地物的垂直結構特征以及光譜特征進行分類，能夠對具有相似垂直結構但光譜特征不同的地物進行更詳細的區分，從而獲得更高的分類精度。

雖然機載激光雷達能夠獲得更為準確的地物垂直結構信息，但其昂貴的價格以及較大的數據處理量使其較難在大尺度尤其是地形復雜區域應用［14］。與之相比，星載激光雷達如搭載在冰、云和陸地高程衛星（ICESat）上的地學激光測高系統（geoscience laser altimeter system，GLAS）因其全球性的采樣策略、較為準確的地理定位以及能夠免費獲取的數據已成為大尺度森林信息獲取的有力工具［14-15］。作為第一顆大腳印激光雷達衛星，ICESat/GLAS 于2003年1月由NASA 發射，到2009年10月終止運行。7年中，GLAS 上的激光器通過40 Hz·s-1的頻率發射信號，以直徑約為65 m 的橢圓光斑（腳印）照亮地球表面并記錄返回的激光能量來獲取測高信息，收集森林地區的激光雷達回波數據超過2.5 億束［16］。作為近年來最常用的激光雷達數據［12］，GLAS 數據被廣泛應用于諸如植被高度測算、植被冠層結構研究、地上生物量估計、林木蓄積量估計、建筑物高度測算以及土地覆被分類等諸多研究領域［14］。然而，其在土地覆被分類方向的研究主要集中在地形較為平坦的地區，而針對山地的相關研究較少。這主要是因為GLAS 數據對地面坡度較為敏感，地形起伏會導致其波形失真從而引起對地物垂直結構的誤判［17］。另外，在結合光學影像進行分類時，山地地區由于山體掩映造成的地形陰影會導致嚴重的“同物異譜”和“異物同譜”的現象，也給山地森林范圍及類型的識別帶來了較大困難［18-19］。而山地森林由于其所受人為干擾較少，往往在區域甚至更大尺度扮演著重要的生態屏障角色，尤其在目前人為土地利用對熱帶森林的影響日益增大的背景下［20］，對山地森林尤其是北方山地森林范圍及類型進行監測顯得尤為重要且迫切。

本研究以位處西北干旱地區的祁連山國家公園肅南縣段為例，結合GLAS 星載激光雷達數據、Landsat OLI影像、Google Earth 高分辨率影像、DEM 數據以及樣地調查數據，綜合利用各數據提供的垂直結構特征、光譜特征、季相特征和地形特征，探索基于多源遙感數據的山地森林識別精度提升方法，在將森林與其他土地利用類型進行區分的基礎上，進一步識別不同的森林類型（針葉林、闊葉林和針闊混交林），以期為今后山地森林識別提供借鑒。本研究旨在回答以下3 個方面的問題：1）GLAS 數據提供的垂直結構信息能否提升山地森林范圍的識別精度；2）不同的地形特征能否提升山地森林類型的識別精度；3）分類時采用不同的波段運算組合其結果是否會有不同。

1 材料與方法

1.1 研究區概況

研究區位于祁連山國家公園肅南縣段，位于97°23′-102°21′E，37°37′-39°37′N，總面積約1.26×104km2（圖1）。祁連山地處青藏高原東北部，甘肅與青海省交界處，是我國主要山脈之一。該地區分布著大量的森林、灌木、草原和冰川，是整個河西走廊的“水庫”，在水土保持和生物多樣性保護中發揮著重要作用。

研究區海拔1770～5740 m，坡度0°～58°。屬典型的溫帶大陸性山地氣候，晝夜溫差大，年平均降水量450～650 mm，年平均蒸發量2134 mm 左右，年平均氣溫-3～1 ℃，近10年來逐漸升高。植被分布受海拔高度影響從低到高可分為荒漠草原、森林草原、亞高山灌叢草甸、高寒荒漠和冰雪帶。現有森林植被主要由青海云杉（Picea crassifolia）、祁連圓柏（Juniperus przewalskii）、山楊（Populus davidiana）和白樺（Betula platyphylla）組成。其他散生樹種有油松（Pinus tabuliformis）、蒙古櫟（Quercus mongolica）和側柏（Platycladus orientalis）。青海云杉主要分布在海拔2500～3300 m 的陰坡和半陰坡上，祁連圓柏主要分布在海拔2700～3300 m 的陽坡、半陰坡或半陽坡上，常與高山云杉混交。山楊和白樺是次生林的主要營建樹種，山楊主要分布在溝壑、溪流邊和海拔2700 m 以下的洪泛區。

1.2 遙感數據

1.2.1 Landsat OLI 數據 搭載在Landsat 8 上的陸地成像儀（operational land imager，OLI），通過一系列如帶寬增強、信噪比提高和附加光譜帶等技術改進，目前已成為Landsat 系列中最強的衛星［21］，且近年來已被廣泛用于區分不同的土地覆被類型甚至不同的森林類型。本研究通過美國地質調查局國家地球資源觀測與科學中心網站（http：//GLOVIS.usgs.gov/）免費獲取覆蓋研究區域的4 景2018年L1T 級Landsat8 OLI 影像，影像基本參數見表1。

表1 Landsat OLI 影像基本參數Table 1 The basic parameters of Landsat OLI images

1.2.2 ICESat/GLAS 數據 ICESat/GLAS 發布了15 個一級和二級數據產品（標簽為GLA1 到GLA15），這些產品可在美國國家冰雪數據中心網站（https：//NSIDC.org/data/icesat）免費獲取。本研究使用了二級GLA14產品，它提供了有關激光腳印的地理位置、表面高程、波形參數和反射比（如信號的起始和結束位置）等信息。同時，雖然最初設計的GLAS 腳印直徑約為65 m，但在實際任務中它們的大小因激光器不同而不同，激光器1、2 和3 的實際腳印直徑分別約為110、90 和55 m［10］。因為3 號激光器腳印的橢圓度最小且接近圓形［22］，本研究僅使用來自3 號激光器的數據，以便量化復雜地形對GLAS 波形的影響。最終在研究區獲取的GLAS 數據包括L3F（2006年5-6月）、L3G（2006年10-11月）、L3H（2007年3-4月）、L3I（2007年10-11月）、L3J（2008年2-3月）和L3K（2008年10月），共計215 個腳印點（圖1）。

圖1 研究區位置示意圖Fig.1 The location of study area

1.3 外業調查數據

本研究的外業調查于2018年10月進行，主要調查內容包括GLAS 腳印點內植被高度測量和分類精度驗證樣本實地采樣兩部分。由于本研究使用的激光雷達GLAS 數據采樣時間為2006-2008年，且研究區內起伏的地形會對GLAS 波形產生較大影響從而引起對腳印點內植被高度的誤判，因此，首先對部分GLAS 腳印點內實際植被高度進行測量以對GLAS 數據提取的植被高度進行校正和驗證。鑒于不同地形坡度對GLAS 波形的影響可能不同［17］，采樣前首先將研究區坡度分為5 級：0°～5°、5°～15°、15°～25°、25°～35°、35°～58°，并通過綜合考慮GLAS腳印點的可及性等相關因素，分層隨機抽樣位于植被區域的65 個GLAS 腳印點進行野外測量。其中12 個腳印點位于0°～5°，18 個腳印點位于5°～15°，13 個腳印點位于15°～25°，17 個腳印點位于25°～35°，5 個腳印點位于35°～58°。在GLAS 腳印點使用精度約為3 m 的手持全球定位系統（global positioning system，GPS）接收器進行空間定位，以所選GLAS 足跡中心為圓心建立直徑為55 m 的圓形采樣區域，采用樣方調查法測量腳印點內植被的高度：將每個樣地分為4 個圓形子樣地，一個位于主樣地的中心，另3 個子樣地分別位于中心樣地360°、240°和120°的方位角，每個子樣地的直徑為15 m，每個子樣地圓心之間的距離為20 m。直接采用手持激光高度計（高度分辨率為0.01 m）測量每個子樣地內植被高度，并計算子樣地所測植被高度平均值。同時記錄腳印點內植被類型、森林類型以及地形特征。腳印點植被平均高度定義為各子樣地植被平均高度的平均值。同時在調查過程中拍攝了腳印點周邊具有代表性的景觀，并記錄其坐標以及土地覆被類型、植被類型等信息，以補充土地覆蓋分類的精度驗證數據。

對用于分類精度驗證的實地調查樣本，綜合研究區各類別土地分布狀況、面積比例等因素，采用分層抽樣的方法獲取典型地類和主要森林類型樣本共計72 個。采用GPS 記錄樣本中心坐標，詳細記錄其土地覆被類型、植被類型等信息并拍攝照片，同時對森林范圍的樣地記錄其森林類型、樹種組成、海拔、坡度、坡向等主要調查因子。確保外業調查所采集樣本涵蓋研究區域典型地類和所有主要森林類型，并盡量在空間上均勻分布。

1.4 輔助數據

Google Earth 上的Quick Bird 冬季影像被用作分類輔助數據，近年來，Google Earth 高清影像資源一直被用作基于遙感的土地覆蓋分類的可靠參考數據源［23］。

本研究采用的DEM 數據是用研究區1∶5 萬地形圖數據生成，利用ArcGIS 10.2 軟件生成研究區的海拔、坡度及坡向分布圖，為GLAS 數據的地形校正以及森林類型識別提供輔助信息。此外，研究區2018年度更新的森林資源二類清查數據（包括林分的詳細信息以及土地覆被類型信息，以下簡稱森林資源清查數據）也被用來輔助建立本次分類的訓練及驗證樣本。

1.5 Landsat OLI 數據預處理

對光學遙感影像進行預處理主要是為了保證影像光譜信息真實且空間位置精確［18］，主要包括遙感影像的輻射定標、大氣校正、正射校正、地形輻射校正、拼接、配準和裁剪。對4 景Landsat OLI 影像做輻射定標和大氣校正，以消除或減少大氣分子和氣溶膠的散射和吸收對地物反射率的影響，使植被波譜曲線趨于正常。因為L1T文件已經使用地面控制點和數字高程模型進行了幾何精校正，所以坐標精度基本能滿足中小比例尺的要求。但由于山區地形復雜，為防止影像發生幾何畸變，通過DEM 數據對影像進行正射校正，以改正影像傾斜和投影差，消除因山區地形起伏和傳感器系統誤差等而引起的像點位移。

在復雜的山區地形中，受地形遮蔽的影響，對應遙感圖像每一像元的地面所接受的有效光照有很大差別，增加了地物光譜信息提取的難度。因此，為提取地物的真實光譜值，必須先進行地形輻射校正以減弱遙感影像地形效應。本研究選擇C 校正模型［24］基于研究區DEM 數據以及影像太陽高度角與太陽方位角信息進行地形輻射校正。隨機選擇研究區域10%的像元作為樣本，用來擬合未校正前地表反射率和局部太陽入射角余弦值之間的線性回歸方程，得到斜率和截距，進一步得到經驗因子c。采用目視分析［25］和去相關分析［26］評價地形校正的結果，發現此次地形校正在坡度平緩地區校正效果較為明顯，校正后陰影區減少，陰坡和陽坡的亮度差異減小。但在高山及丘陵集中區域，在進行地形輻射校正之前，地表反射率與局部太陽入射角余弦值之間的相關系數為0.7149，校正后相關系數為0.5694，說明出現了較為嚴重的過校正現象［26］。因此就整個研究區域而言，地形校正效果不盡理想，仍需就“同物異譜”以及“同譜異物”的地物光譜值的提取進行進一步的研究。

對經過校正后的4 景Landsat OLI 影像通過直方圖匹配算法進行鑲嵌，以Google Earth 高清影像作為參考，選取位于道路、河流交叉處以及建筑物點的79 個控制點，采用二次多項式模型進行相對配準，總體均方根誤差（total root mean square error）為0.481。最后采用研究區域的矢量范圍文件對影像進行裁剪。

1.6 ICESat/GLAS 數據預處理

1.6.1 無效數據過濾 由于GLAS 數據質量易受到云污染、信號飽和以及大氣前向散射的影響，因此，依照以下兩個標準首先對無效數據進行過濾：1）使用GLA14 產品中記錄的檢測標志（i_satNdx=0，i_FRir_qaFlag=15）識別無飽和以及無云的GLAS 數據并保留［27］；2）將GLA14 數據中記錄的高程（d_elev）與DEM 數據進行比較，由于研究區域內均為低矮建筑物，因此將高差閾值設為大于該地區最大冠層高度（50 m），當GLAS 數據與DEM 數據的高差超過70 m 時，將其視為潛在的低云數據，并將其從分析中刪除［28］。

1.6.2 地理定位精度驗證 GLAS 數據在研究區的地理定位精度對其能否作為本次研究分類依據至關重要。因此，參照前人研究［29］將GLAS 腳印點的地面高程與DEM 高程值進行比較以評估其定位精度。

首先利用GDAL 開源庫，依據GLA14 產品中記錄的每個腳印點的地理位置與腳印長短軸及方位角，將GLA14 的HDF5 格式的數據轉換為ESRI ShapeFile 矢量格式，并轉換為WGS84 坐標系。對每個腳印內的DEM數據所有像素的高程進行加權平均，與GLA14 數據記錄的高程值進行比較，根據距腳印中心的距離，使用以下公式對DEM 中每個像素的值進行加權［30］。

式中：W是GLAS 腳印點內每個DEM 像素的權重；a和b分別是GLAS 腳印的半長軸和半短軸，記錄在GLA14產品中；（x，y）和（x0，y0）分別是DEM 像素和GLAS 腳印點的中心坐標；x′和y′是DEM 像素沿腳印點長軸和短軸的坐標；α是足跡長軸的方位角，同樣記錄在GLA14 產品中。在加權平均之前，將腳印點內所有DEM 像素權重之和歸一化為1。

共計隨機選取30 個GLAS 腳印點對其地理位置精度進行驗證，相關系數R2為0.964，均方根誤差RMSE（root mean square error）為1.49 m，因此認為本研究區域GLAS 數據的地理定位精度誤差較小［31］，能夠滿足分類的需要。

1.6.3 GLAS 腳印點提取地物高度的校正 在GLA14 產品中，記錄了GLAS 原始波形通過高斯濾波器擬合的1 個波形信號開始范圍增量（d_SigBegOff）以及最多6 個高斯波峰質心范圍增量（d_gpCntRngOff），本研究將波形信號開始范圍增量作為地物頂端波形信號，最低兩個高斯波峰中較強峰質心范圍增量近似表示地表波形信號。對其求差所得的波形幅度（Wext）可以表示從GLAS 中提取的地物高度，這種方法也被廣泛應用于從GLAS 數據中提取森林冠層高度［32］。

在平坦地區，Wext可以近似表示實際的地物高度，但在山區，地形起伏會對GLAS 波形形狀產生影響，導致從中提取地物高度出現誤差。因此，本研究參考前人研究［32-33］，采用一種基于幾何模型的地形校正方法，引入平移系數a和b，得到經過地形校正后從GLAS 數據中提取的地物高度為：

式中：d表示GLAS 腳印點直徑，本研究中取值為55 m；θ表示腳印點所在的坡度等級的平均坡度值（°）；c表示光速（m·ns-1）；FWHM（full width at half maximum）表示激光脈沖寬度，為6 ns；a，b為平移系數，采用通用全局優化算法（levenberg marquard，LM）通過實地調查GLAS 腳印點的實際測量植被高度擬合得到；HGLAS表示地形校正后從GLAS 數據中提取的地物高度。

由于所使用的GLAS 數據與待分類影像有近10年的時間差距，因此在利用GLAS 數據分類之前首先通過咨詢當地林業及土地資源管理部門將10年間新增的植樹造林地區、林木砍伐地區、火燒跡地以及道路及房屋建設或拆遷區域的GLAS 腳印點進行剔除，以免造成土地覆被類型的誤分。另外，由于GLAS 腳印點提取的地物高度主要用于分類波形結構相似但垂直高度不同的植被（森林、灌木、草地和農田），因此，GLAS 腳印點實地調查主要在植被區域進行。而10年間喬木的生長量會對平移系數a和b的擬合精度產生影響，從而影響從GLAS 數據中提取的除喬木外其他植被高度的精度。通過分析研究區森林資源清查數據可知，研究區喬木林齡范圍在9～100年，其中林齡在15年以下的林分僅占總林分的0.29%，且林齡較小的樹種多為楊類，區內中齡林及近熟林為區內森林主要齡組，占總林分的97%以上。因此，認為10年間喬木的高度增長量對灌木與喬木的區分影響較小。并在采用GLAS 實地調查數據擬合平移系數之前，根據樹齡及樹種生長曲線將10年間的樹高增加值從實際測量的喬木高度中減去［34］，以減小因時間差距導致的GLAS 數據提取植被垂直高度的誤差對分類精度的影響。

1.7 森林范圍識別精度比較

為明確GLAS 激光雷達數據提供的地物垂直結構特征對森林范圍識別精度的影響，在分類過程中，本研究僅依據光譜特征以及綜合依據光譜及垂直結構特征分別建立分類訓練樣本。根據研究區地表覆蓋狀況和應用需求，參考《國家森林資源連續清查技術規定》地類劃分標準，將項目區土地分類為森林、灌木、草地、農田、裸地、建設用地和水域。

以研究區GLAS 腳印點為分類訓練樣本綜合依據光譜及垂直結構特征識別研究區不同地類。由于本研究使用的GLA14 產品經過統一的濾波和高斯分解擬合處理，其波形形狀主要取決于腳印點內物體的垂直結構，因此，首先根據不同波形形狀將研究區分為3 類土地覆被類型：平坦表面、建筑物和植被［12］。一般將GLA14 產品僅具有一個峰值（i_numPk=1）的擬合波形定義為平坦表面［28］。將具有兩個峰值（i_numPk=2）的GLAS 擬合波形定義為建筑物［12-13］。值得注意的是，有研究證明［12］，兩個峰值的GLAS 波形一般為平頂的建筑物，而坡頂的建筑物由于其屋頂表面未在同一水平面會導致GLAS 波形更加復雜。由于研究區內建筑物多為平頂民居或低層的工業建筑，因此本研究將i_numPk=2 的GLAS 數據作為建筑物樣本。植被由于其通常有多個目標在同一GLAS 腳印點中或其本身垂直結構（葉片、枝干、冠層等）的復雜性，其波形往往表現出比以上地物更為復雜的特征，因此將具有3 個以上且數量不定的波峰（3≤i_numPk≤6）的擬合波形定義為植被［35］。

根據腳印點光譜特征對以上3 大類土地覆被類型進行進一步細分，其中植被由于其光譜差異較小，但不同植被類型其垂直高度不同，因此依據高度對其進行分類［36］。參考項目區實際情況及前人研究［35］，當HGLAS≥6.0 m，將其分類為森林，當HGLAS為1.5～6.0 m，將其分類為灌木，當HGLAS≤1.5 m，將其歸類為草地或農田。

從研究區森林資源清查數據中獲取不同類型的土地覆被樣本數據，建立僅依據光譜特征的訓練樣本對研究區影像再次進行分類并將結果與上述分類結果進行對比，以總體精度（計算混淆矩陣）和Kappa 系數作為分類精度評價指標。為避免分類精度是由于增加了分類樣本提升的，保持兩次分類的訓練樣本像素數量一致。另外，為驗證影像分類時Landsat OLI 影像參與運算的波段對分類結果產生的影響，本研究在兩次分類時分別按照OLI 影像分類時識別植被常用的不同波段組合參與運算，并對結果進行比較。

1.8 森林類型識別精度比較

由于不同的森林類型很難依靠GLAS 數據中記錄的波形信息對其進行區分［12，36］，而在依據光譜信息區分森林類型時，由于研究區地形復雜，雖然進行了正射校正以及地形輻射校正等預處理，仍難以完全消除地形影響導致的光譜異常、信息畸變和類別混淆。而山地不同坡向往往分布著不同類型的森林，更使得“同物異譜、同譜異物”現象對不同森林類型的識別造成困難。因此，本研究引入季相及地形特征［37-38］對不同類型的森林進行分類，并比較相關特征的加入對分類精度的提升程度。

首先將獲取于冬季的Google Earth 高清影像與Landsat OLI 影像進行對比，借助季相特征，根據森林資源清查數據中不同森林類型的面積比例，采集森林類型單一集中連片，面積約近100 m2的林分，隨機建立針葉林、闊葉林以及針闊混交林的分類訓練樣本。然后將通過研究區DEM 數據生成的海拔及坡向分布圖與影像疊加并設置一定的透明度，對同一森林類型根據不同海拔及坡向的影像光譜特征進行取樣，建立分別位于不同海拔及坡向的針葉林、闊葉林以及針闊混交林的分類訓練樣本。最后，分別依據隨機建立的森林類型訓練樣本，位于不同海拔的森林類型訓練樣本、不同坡向的森林類型訓練樣本以及不同海拔及坡向的森林類型訓練樣本進行分類以識別不同的森林類型，以總體精度和Kappa 系數作為評價指標對各分類結果的精度進行評價。

2 結果與分析

2.1 GLAS 腳印點提取植被高度校正結果

按照不同的坡度級別，依據調查GLAS 腳印得到的實測植被高度利用LM 通用全局優化算法對方程（5）進行參數擬合，得到各坡度等級的平移系數a和b的最佳估計值分別為：0°～5°，a=1.15，b=2.91；5°～15°，a=0.58，b=1.13；15°～25°，a=0.62，b=0.60；25°～35°，a=0.38，b=0.52；35°～58°，a=0.29，b=0.51。基于Google Earth高分辨率影像隨機選取51 個位于森林的GLAS 腳印點（其中6 個腳印點位于0°～5°，14 個腳印點位于5°～15°，11個腳印點位于15°～25°，14 個腳印點位于25°～35°，6 個腳印點位于35°～58°），依據擬合的平移系數計算GLAS 數據提取的森林冠層高度HGLAS。以研究區2018年更新的森林資源清查數據作為獨立驗證數據，將10年間高度增加值從其記錄的森林冠層高度中減去并分別與地形校正前直接從GLAS 數據中提取的波形寬度Wext以及地形校正后GLAS 數據提取的冠層高度HGLAS進行比較，得到各坡度等級地形校正前后的相關系數（R2）和均方根誤差（root mean square error，RMSE）（表2）。

表2 不同坡度等級地形校正前后GLAS 提取植被高度的精度比較Table 2 Accuracy comparison of GLAS derived canopy heights before and after topographic correction in different slope gradient classes

地形校正前，隨著坡度的增加，GLAS 數據提取樹木高度的RMSE 從2.83 m 增加到12.59 m。而地形校正后，不同坡度等級GLAS 數據提取樹木高度的RMSE 穩定在1.84～3.91 m，在很大程度上消除了地形坡度對GLAS 波形的影響，從而使得從GLAS 數據中所提取的地物高度可用于下一步的分類研究。

2.2 森林范圍識別結果

按照光譜特征對GLAS 腳印點中被分類為平坦地面的道路、裸地和水域進行進一步區分，3 種地物在各波段光譜特征趨勢均具有明顯差異（圖2），因此結合光譜特征對其進行分類，并將道路與建筑物樣本合并為建設用地分類依據。

圖2 道路、裸地以及水域的平均光譜曲線Fig.2 Average spectral curves of road，bare land and water

可知，OLI 影像上草地和農田在紅波段（Band4）和近紅外波段（Band5）表現出明顯的植被光譜響應“峰谷”特征，且峰值有一定差異（圖3），另外，本研究還參照冬季的Google Earth 高清影像，依據項目區耕作特點，對草地和農田的分類進行輔助。

圖3 草地和農田的平均光譜曲線Fig.3 Average spectral curves of grassland and farmland

由于GLAS 腳印點內可能存在不止一類土地覆被類型（如混合有建筑物、農田及樹木的腳印），因此只有腳印內某地類覆蓋其足跡面積50%以上才被視為該地類。按照OLI 影像常用于植被識別的波段組合參與運算，采用目前已被證實分類精度較好的基于支持向量機的監督分類方法分別依據光譜特征以及光譜和垂直結構綜合特征對研究區進行分類。基于實地調查數據以及森林資源清查數據，通過分層隨機抽樣生成在空間上均勻分布的389 個驗證樣本（共3627 個像元），分別對分類結果進行精度驗證（表3）。

表3 不同波段組合及分類方法基于不同特征的分類結果對比Table 3 Comparison of classification results based on different features of different band combination and classification methods

可以看出，不同的光學波段組合運算對分類精度并無影響。同時，相比僅依據光譜特征進行分類，加入垂直結構特征后，分類精度明顯提升，總體分類精度提高了10.67%（表3）。為了進一步明確GLAS 垂直結構特征對各地類分類精度的影響，本研究以NIR、Red、Green 波段組合為例比較了基于光譜特征和基于光譜與垂直結構綜合特征分類的不同地類的分類精度（表4）。

由表4可知：1）加入垂直結構信息后，森林、灌木及草地的分類精度明顯提升，漏分誤差（omission error）分別減少了3.38%，24.82%以及31.05%，錯分誤差（commission error）分別減少了21.47%，14.48%以及5.30%。說明垂直結構信息對光譜特征相似但垂直結構不同的植被的區分具有明顯的協助作用。2）建設用地以及水域的制圖精度和用戶精度均有所提升。這主要是因為研究區建設用地面積通常較小，僅依靠光譜信息分類時，容易與周邊地類光譜信息混淆，而GLAS 腳印點垂直結構信息的加入增加了建設用地的位置精度，從而提高了分類精度。垂直結構信息的加入同時消除了由于光線反射或地形陰影造成的水域的錯分和漏分現象。3）加入垂直結構信息后，裸地的用戶及制圖精度均略有降低，這可能是因為項目區的裸地多為礫石地面，同時山地地形的復雜性也增加了GLAS 數據波形的復雜度，因此加入垂直結構信息后反而使分類精度略有降低。

表4 依據不同分類特征的各土地類型分類精度比較Table 4 Comparison of classification accuracy of different land types based on different classification characteristics

為了進一步驗證方法的應用精度，在研究區選擇面積為3671.14 km2試驗區（占總面積的29%），將基于光譜與垂直結構綜合特征分類的面積統計結果與研究區2018年度地理國情及影像解譯面積統計結果進行對比。

從表5可以看出，基于光譜與垂直結構綜合特征分類結果中各地類所占總面積比重與調查及解譯結果各地類比重基本一致，7 種土地覆被類型面積平均相對精度為92.09%（表5），說明將Landsat OLI 影像的光譜特征以及經地形校正后的GLAS 數據的垂直結構特征相結合用于山區土地分類能夠提高山區森林范圍識別的精度。

表5 各地類分類面積統計與對比Table 5 Area statistics and comparative analysis of different classified land types

2.3 森林類型識別結果

依次根據不同森林類型的隨機樣本，加入海拔信息特征樣本、坡向信息特征樣本和海拔及坡向信息特征樣本分別對所提取的森林按照針葉林、闊葉林及針闊混交林進行森林類型分類。采用基于支持向量機的監督分類方法以NIR、Red、Green 波段組合參與運算，保持每次分類的訓練樣本數量一致。

根據不同森林類型的面積比例，基于外業調查數據及森林資源清查數據，選取在不同海拔級別（1770～2770 m、2770～3770 m 和3770～4770 m）以及坡向（陰坡、陽坡、半陰坡和半陽坡）上均勻分布的305 個樣本點（共2594個像元）對各分類結果進行驗證。由于研究區內3770～4770 m 海拔上僅分布少量的針葉林，而4770～5740 m 海拔沒有森林分布，因此僅對1770～2770 m 和2770～3770 m 海拔上的森林類型分類精度進行了驗證（表6）。

可以看出，地形信息的加入能夠顯著提升森林類型分類精度，而相比海拔信息，坡向信息的加入對提升分類精度效果更為顯著（表6）。這主要是因為地形信息的加入能夠盡量全面的考慮到不同地形特征上各森林類型的不同光譜特征，且山體掩映造成的“同物異譜”以及“同譜異物”現象主要由光照區和陰影區造成。另外，雖然森林類型的空間分布格局在區域尺度上主要受到海拔以及緯度的影響，但森林類型分布格局在局部尺度上主要受到諸如微氣候、光以及可獲取水分等因素的影響，而在受水分條件限制的半干旱和干旱地區森林類型分布受這些因素的影響尤為顯著。因此，就本研究區域而言，加入坡向信息比海拔信息更能夠提高森林類型的分類精度。

表6 地形信息對森林類型分類精度的影響比較Table 6 Comparative analysis of the effect of terrain information on the accuracy of forest types classification

坡向信息加入后對各坡向森林類型分類精度的提升作用由高到低依次為半陰坡（21.76%）、陰坡（17.78%）、半陽坡（12.67%）以及陽坡（9.52%）。當同時加入坡向和海拔信息時，分類精度比單獨使用海拔或坡向信息時均有提升，但并不是兩者提升精度的線性加和，這表明海拔和坡向從不同的空間角度上提高了分類精度，同時二者提供的信息可能存在重疊。本研究進一步對不同的森林類型分類精度進行了對比，當依據海拔及坡向特征分類時，針葉林的相對精度最高，主要因為其存在大面積純林，相對最容易提取，而混交林由于其通常以小片面積分布于山地陰影地區的闊葉林及針葉林過渡地帶，所以受地形因素影響較大，相對分類精度最低。將以上兩層次的最優精度分類結果合并，得到研究區森林信息及土地覆被（圖4）。

圖4 研究區森林類型及土地覆被分布Fig.4 Forest types and land cover distribution of study area

3 討論

本研究綜合利用多源遙感數據提供的垂直結構信息、光譜信息、季相和地形信息探索山地森林識別精度提升方法，研究的不確定性主要來自以下幾個方面：

由于本研究獲取的星載激光雷達GLAS 數據的采樣時間與其他數據（樣地調查、遙感影像和森林資源清查數據）的采樣時間有約10年的差距，因此，在利用樣地調查數據擬合地形校正模型的平移系數時，首先將10年間的樹高增加值從實測的喬木樹高中減去，以避免造成GLAS 數據提取的除喬木外其他植被（灌木、草地、農田）高度的誤差。雖然區內中齡林及近熟林占總林分的97%以上，但這一方法仍會對部分幼齡林的識別產生較大影響，給森林范圍的識別帶來不確定性。另外，在利用GLAS 數據分類之前，為避免10年間新增的植樹造林地區、林木砍伐地區、火燒跡地以及道路及房屋建設或拆遷區域的GLAS 腳印點影響分類，對這一部分的腳印點進行剔除，這使得能夠用于分類訓練樣本的GLAS 腳印點減少，也是影響整體分類精度的一個重要因素。

為了便于量化地形坡度對GLAS 波形產生的影響，本研究選擇僅來自3 號激光器的GLAS 數據以統一腳印大小。為獲取盡可能多的數據以用于分類，導致所獲GLAS 數據的季節性差異較大（數據采樣時間為2-6月和10-11月）。由于葉片是激光能量作用的主要介質，冬季落葉以及冠層積雪可能引起植被高度估計的誤差［39］，從而降低分類精度。另外，為避免云覆蓋對土地覆被分類的影響，所下載Landsat OLI 影像的季相差異也會降低不同森林類型的識別精度。

作為第一顆大腳印星載激光雷達數據，GLAS 數據在本研究中的腳印點被認為近似圓形且直徑為55 m，其內存在不止一類土地覆被類型的腳印點，往往將其識別為腳印點中的主要土地覆被類型。由此造成的混合像元問題也成為影響分類精度的重要原因之一。

2018年9月，NASA 成功發射了延續早先的ICESat 任務工作的ICESat-2，其上搭載的先進地形激光測高儀系統（advanced topographic laser altimeter system，ATLAS）相比GLAS 具有不同的傳感器規格以及不同的采樣模式，激光雷達采樣強度和空間分辨率均大大提高，對植被結構的監測更加精確［40］。其提供的ATL08 產品，已于2019年春季開始陸續向公眾提供2018年10月至今的全球陸地與植被高度數據［40］。雖然本次研究在制定外業調查方案時未能獲取到該產品，但在未來的研究中，將利用ATL08 更加密集的采樣點與更精確的植被測高數據，與遙感圖像光譜信息以及紋理信息相結合，進一步提升山地森林類型乃至樹種的識別精度。

4 結論

本研究以祁連山國家公園肅南縣段為例，綜合利用GLAS 星載激光雷達數據、Landsat OLI 影像、Google Earth 高分辨率影像、DEM 數據以及樣地調查數據對山地森林信息進行分層次提取，比較了加入垂直結構特征、季相特征和地形特征后對山地森林范圍及類型識別精度的提升作用，結果表明：1）相比僅依據光譜特征進行分類，依據光譜及垂直結構綜合特征分類時具有相似光譜特征但不同垂直結構的地物分類精度明顯提升，森林范圍的識別精度提高。說明GLAS 數據提供的垂直結構信息的加入能夠增強不同類別之間的可分性，使土地類型分類及森林范圍識別更為有效。2）地形信息的加入能夠顯著提升森林類型識別精度，而相比海拔信息，坡向信息對提升分類精度效果更為顯著。3）多源多時相遙感影像提供的季相特征能夠對不同森林類型的識別提供幫助，而不同波段組合雖然對地物增強的效果不同，但其對分類精度幾乎沒有影響。

另外，雖然GLAS 數據提供的垂直結構信息與光譜信息結合能夠提高山區復雜地形條件下中等空間分辨率遙感影像森林范圍識別的精度，但由于GLAS 數據對地形較為敏感，在坡度大于10°的情況下，波形失真會導致地物高度的估計精度降低，因此在山區應用GLAS 數據時，必須首先對其進行校正以消除地形對波形產生的影響。同時，由于GLAS 數據能夠提供2003-2009年間的地球表面測高數據，因此其對歷史影像的分類精度也具有較好的提升作用，能夠為分析森林變化情況以及制定自然資源管理政策提供依據。