劉美爽 邢艷秋 李立存 楊 超 王 蕊
(東北林業(yè)大學,哈爾濱,150040)
森林生態(tài)系統(tǒng)作為陸地生態(tài)系統(tǒng)的主體,在維持氣候穩(wěn)定、調(diào)節(jié)全球碳平衡等方面有著不可替代的作用。森林生物量,是評估森林碳貯量的重要參數(shù),其精確測量可以為森林經(jīng)營與資源管理者制定政策時提供參考依據(jù)。森林生物量的估算,主要通過森林垂直結(jié)構(gòu)參數(shù)與水平結(jié)構(gòu)參數(shù)來實現(xiàn)。目前,星載激光雷達(ICESat-GLAS)數(shù)據(jù)已經(jīng)成功用于反演森林垂直結(jié)構(gòu)信息;而在反演森林水平結(jié)構(gòu)信息方面,特別是森林類型,因方法不成熟等原因,結(jié)果不很理想[1-2]。
支持向量機(SVM),是一種基于結(jié)構(gòu)風險最小化的新型機器學習方法,具有推廣性能優(yōu)越、全局收斂、不依賴經(jīng)驗信息等優(yōu)點,在模式識別領(lǐng)域取得了一定的研究進展[3-5]。支持向量機,能有效的解決小樣本條件下高維數(shù)據(jù)模型的構(gòu)建問題,在解決有限樣本、高維、非線性模式識別問題有更好的性能[6-7]?;菸娜A[8]提出了一種基于支持向量機的遙感圖像分類方法,結(jié)果表明,支持向量機的遙感圖像分類精度明顯高于傳統(tǒng)最大似然法的分類精度;劉志剛等[9]提出了一種基于加權(quán)無標識支撐向量機,并在此基礎(chǔ)上提出了一種新的不完全監(jiān)督分類方法,該方法在數(shù)據(jù)模擬和遙感影像分類中均達到了理想的分類精度;員永生[10]基于支持向量機研究了面向?qū)ο蟮男滦屯恋馗脖粓D像分類,并在此基礎(chǔ)上提出了5類新型的遙感圖像分類方法,研究結(jié)果顯示,新型方法分類結(jié)果好于傳統(tǒng)的K最鄰近節(jié)點算法(KNN)面向?qū)ο蟮谋O(jiān)督分類識別方法。對于線性不可分的問題,引進了線性不可分支持向量分類機(C-SVC),該方法的思想是通過引進一個非線性映射,將低維線性不可分問題轉(zhuǎn)化為高維線性可分問題。
為提高星載激光雷達識別森林類型的精度,本文提出一種基于線性不可分支持向量分類機的星載激光雷達數(shù)據(jù)識別森林類型方法,并采用K-折交叉驗證方法對核函數(shù)選擇進行評價。
研究區(qū)為長白山汪清林區(qū)。該區(qū)處于寒溫帶森林生態(tài)系統(tǒng),位于長白山系的中低山區(qū)(43°05'~43°40'N,129°56'~131°04'E),總面積 30.4 萬 hm2。地面高程為360~1 477 m,坡度變化范圍為0~45°。該區(qū)域植物種類繁多,結(jié)構(gòu)復雜。針葉樹主要有紅松(Pinus koraiensis)、云杉(Picea)和臭松(Symplocarpus salisb);闊葉樹多為椴樹(Tilia)、蒙古櫟(Quercus monglica)、楓樺(Betula davuric)、色木(Acermono)和白樺(Betula platyphylla)等。
ICESat-GLAS是極地軌道星載激光雷達設(shè)備,飛行高度600 km,用于全球地面的連續(xù)觀測。GLAS脈沖激光器每秒發(fā)射40個激光脈沖,脈沖寬度為4 ns,每個激光脈沖在地面上覆蓋的光斑直徑大約為70 m,相鄰光斑間距為170 m。ICESat-GLAS共提供15個數(shù)據(jù)產(chǎn)品,即 GLA01、GLA02、…、GLA15。其中,全球測高數(shù)據(jù)產(chǎn)品GLA01記錄的完整波形數(shù)據(jù),反映了對應(yīng)地面激光光斑內(nèi)的地物信息;與波形數(shù)據(jù)相應(yīng)的地面光斑的地理位置和高程數(shù)據(jù),由高程數(shù)據(jù)產(chǎn)品GLA14記錄。本研究使用的GLA01和GLA14數(shù)據(jù),由美國國家冰雪數(shù)據(jù)中心(http://nsidc.org/data/ice-at)提供,共獲取了該研究區(qū)的2003—2009年的波形數(shù)據(jù)3 167組。
在2006、2007、2010年9月份,3次在長白山汪清林區(qū)根據(jù)預(yù)設(shè)計的采樣計劃,采用分層隨機采樣法,針對2種森林類型(針葉林和闊葉林),隨機選取287個GLAS激光光斑進行相應(yīng)的地面數(shù)據(jù)調(diào)查(其中:針葉林地96個,闊葉林地159個,其余樣地為裸地、農(nóng)田、草甸和濕地等)。以針葉林和闊葉林為研究對象,對255個樣地進行分析;剔除8個無效樣地后,在剩余的247個樣地中,隨機選取164個樣地數(shù)據(jù)作為訓練集(其中針葉林地60個,闊葉林地104個),其余83個樣地數(shù)據(jù)作為驗證集(其中針葉林地35個,闊葉林地48個),利用GLAS波形進行森林類型識別研究并進行驗證。在樣地調(diào)查過程中,利用GPS對已選定的激光光斑采樣點進行定位。依據(jù)森林調(diào)查的統(tǒng)計原理,為有效地對樣地內(nèi)林木進行調(diào)查,首先測量樣地的坡度;再以光斑中心點為圓心,建立水平投影面積為500 m2的圓形樣地;記錄樣地內(nèi)植被分布情況、植被類型和植被覆蓋度。結(jié)合我國森林資源調(diào)查主要技術(shù)規(guī)定,將針葉林分布比例≥60%的樣地定義為針葉林,而闊葉林分布比例≥60%的樣地定義為闊葉林。
C-SVC具體算法如下:
(1)給定樣本訓練集。T={(x1,y1),(x2,y2),… ,(xn,ym)},其中n維輸入向量xi∈Rn,yi∈{1,-1},i=1、…、n。
(2)選擇適當?shù)暮撕瘮?shù)K(xi,xj)和適當?shù)膽土P參數(shù)C,構(gòu)造并求解最優(yōu)化問題。目標函數(shù):
約束條件:
其中,核函數(shù)K(xi,xj)是滿足Mercer條件的內(nèi)積函數(shù),采用不同的核函數(shù),就可以構(gòu)造不同的分類機。常用的核函數(shù)包括:q次多項式核函數(shù)、徑向基核函數(shù)、Sigmoid核函數(shù)、Fourier級核函數(shù)和B-樣條核函數(shù)等。前3種核函數(shù),針對不確定性對象[11]、無背景數(shù)據(jù)[12]等有著相對優(yōu)良的性能??紤]到研究區(qū)森林水平結(jié)構(gòu)的非均質(zhì)性,分別采用q次多項式、徑向基以及Sigmoid三種不同的核函數(shù)構(gòu)造分類機,具體函數(shù)表達式如公式(3)~公式(5)所示。
q次多項式核函數(shù):
徑向基核函數(shù):
Sigmoid核函數(shù):
參數(shù)C決定了誤判樣本的懲罰程度以及學習機器的復雜性;取值越小,表明對經(jīng)驗誤差的懲罰越小,而且學習機器越簡單;但是,增大了經(jīng)驗風險值。隨著C的增加經(jīng)驗風險值會逐漸減小,當達到某個適當值時,訓練誤差達到最小,精度達到最大。參數(shù)的選擇問題,本質(zhì)是一個優(yōu)化問題,主要包括交叉驗證、貝葉斯法、梯度下降等,其中交叉驗證是目前比較普遍的一種應(yīng)用方法。采用K-折交叉驗證方法進行K-交叉驗證,是在實驗過程中,將數(shù)據(jù)集分成K份,輪流將其中K-1份做訓練1份做測試,K次結(jié)果的最佳值作為對算法精度的估計[13]。算法中K取10,采用十折交叉驗證方法來計算準確率,進而確定合適的懲罰參數(shù)C及核函數(shù)中的系數(shù)r的值。
為進一步驗證所選核函數(shù)C-SVC的分類效果,利用驗證集的83組數(shù)據(jù)對其分類結(jié)果進行了分析評價。
ICESat-GLAS波形數(shù)據(jù)采用專門定義的二進制格式,包括元數(shù)據(jù)信息和其它數(shù)據(jù)信息。首先在IDL平臺上將原始二進制數(shù)據(jù)轉(zhuǎn)化成ASCII格式數(shù)據(jù);為了有效比較波形數(shù)據(jù),再對波形數(shù)據(jù)進行標準化處理,并利用高斯算法對原始ICESat-GLAS回波進行分解和擬合處理,在Matlab平臺上輸出波形數(shù)據(jù),詳細過程參考文獻[14-16]。
從森林自身特征分析,由于存在葉片大小,形狀和層片結(jié)構(gòu)的差異,針葉林和闊葉林分別具有各自典型的垂直結(jié)構(gòu),會分別對應(yīng)不同的典型激光雷達(LIDAR)回波波形特征;從數(shù)學角度分析,LIDAR回波波形可以看作是n個高斯曲線的組合[17],分別對應(yīng)著森林不同層片及其分布的回波能量特征。據(jù)此可以推斷,通過分析LIDAR回波高斯曲線特征的差異,可以實現(xiàn)森林類型的識別。
針葉林的垂直分層很明顯地多于闊葉林,而針葉林層厚則明顯小于闊葉林。因此,當LIDAR接收針葉林層間返回的能量時,與其相關(guān)的回波寬度要小于闊葉樹的相應(yīng)值。這樣,針葉林對應(yīng)的波形分解高斯曲線斜率就要大于闊葉林[18]??紤]到以上因素,基于擬合后的GLAS波形,提取并分析與高斯曲線斜率相關(guān)的波形特征參數(shù),包括曲線斜率、曲線斜率的均值和曲線斜率的標準均方根誤差,進行森林類型區(qū)分的研究。
除對應(yīng)地面回波的最后一個高斯曲線外,對其余所有高斯曲線提取有效回波點處對應(yīng)的時間值、波峰處標準化的能量值和波峰處對應(yīng)的回波時間值,其中,有效回波的臨界值由背景噪聲峰值與4倍的標準方差之和確定(以圖1所示的GLAS波形為例)。
圖1 參數(shù)提取
當原始曲線分解為n個高斯曲線時,曲線的斜率Ki如公式(6)所示。
式中:Qi為第i個高斯曲線標準化能量值;t2i-1為第i個高斯曲線有效回波點對應(yīng)的時間值;t2i為第i個高斯曲線波峰處對應(yīng)的時間值。
曲線斜率均值ˉK和曲線斜率標準均方根誤差ΔK分別如公式(7)和公式(8)所示。
為克服各指標之間量綱的影響,首先將提取的波形參數(shù)進行標準化。采用Z-score標準化(標準差標準化)方法,將需要標準化的參數(shù),包括高斯分解后第一個高斯曲線的斜率K1、冠層高斯曲線斜率的均值ˉK和冠層高斯曲線斜率均值的標準差ΔK等3個變量如式(9)所示進行標準化。
式中:x*為標準化后的波形特征參數(shù);x值為需標準化的波形特征參數(shù);μ為針葉林和闊葉林各類分類樣本對應(yīng)的該波形特征參數(shù)均值;σ為各類分類樣本數(shù)據(jù)對應(yīng)的該波形特征參數(shù)標準差。在此后提到的K1、ˉK和ΔK均為標準化后的參數(shù)。
基于2.3所述的算法及其運算步驟,以K1、ˉK、ΔK三個變量的組合為輸入量,利用3種核函數(shù)(多項式核函數(shù)(Poly)、徑向基核函數(shù)(Rbf)和Sigmoid核函數(shù))及不同的懲罰因子(0.01、0.10、1.00、10.00、100.00、1 000.00)來構(gòu)造不同的分類機 C-SVC,對164組訓練集樣本進行森林類型分類訓練,并用十折交叉驗證法對分類結(jié)果進行精度評價。
基于 ICESat-GLAS波形數(shù)據(jù)處理,ICESat-GLAS波形數(shù)據(jù)在Matlab平臺上輸出波形數(shù)據(jù)如圖2所示。
GLAS回波經(jīng)高斯分解后得到的第一個高斯曲線反映了森林冠層的信息,考慮到闊葉林和針葉林的冠層結(jié)構(gòu)存在著一定的差異,因此,利用第一個高斯曲線斜率值K1對森林的冠層進行分析,進而發(fā)現(xiàn)針葉林和闊葉林的不同特征。
圖3、圖4、圖5分別為針葉林地和闊葉林地所對應(yīng)的GLAS波形第一個高斯曲線的斜率K1、森林冠層高斯曲線斜率均值和森林冠層高斯曲線斜率均方差ΔK分布情況。綜合分析圖3、圖4、圖5及表1發(fā)現(xiàn),針葉林的K1、和ΔK三個特征參數(shù)均值都大于闊葉林;表明這三個參數(shù)能夠在一定意義上體現(xiàn)兩種森林類型冠層的不同特征。因此,利用K1、、ΔK三個變量的組合作為輸入特征參數(shù),基于C-SVC方法進行森林類型識別。
圖2 GLAS波形數(shù)據(jù)分解和擬合結(jié)果
圖3 第一個高斯曲線斜率K1
圖4 植被冠層高斯曲線斜率均值
表1 不同森林類型冠層K1、、ΔK值比較
表1 不同森林類型冠層K1、、ΔK值比較
注:表中“-”為負號。
波形參數(shù) 森林類型 最大值 最小值 平均值 標準差K1 闊葉 2.21 -0.90 -0.31 0.64針葉 4.04 -0.82 0.49 1.19ˉK闊葉2.33-1.23-0.270.71針葉 3.99 -1.17 0.37 1.19 ΔK 闊葉 4.91 -1.26 -0.08 0.99針葉 2.58 -1.27 0.22 1.00
圖5 高斯曲線斜率均方差ΔK
以K1、、ΔK三個變量的組合為輸入量,利用上述3種核函數(shù)及不同的懲罰因子來構(gòu)造不同的分類機C-SVC,對164組訓練集樣本進行森林類型識別訓練,并用十折交叉驗證方法對識別結(jié)果進行精度評價(見表2)。
表2 不同懲罰參數(shù)(C)及核函數(shù)對應(yīng)的交叉驗證分類精度對比
分析表2的結(jié)果,發(fā)現(xiàn)徑向基核函數(shù)構(gòu)建的CSVC分類效果最好,精度可達80.95%;多項式核函數(shù)次之;Sigmoid核函數(shù)最差(見表3)。
表3 不同懲罰因子(C)及核函數(shù)系數(shù)(r)所對應(yīng)的交叉驗證分類精度對比
利用徑向基核函數(shù)構(gòu)建的C-SVC對研究區(qū)森林類型進行識別,并進一步確定分類機最優(yōu)懲罰因子(C)和核函數(shù)系數(shù)(r)。由表3可見,當C=100.00、r=1.000時,分類效果最好,交叉驗證分類精度可達85.24%。
利用驗證集83組數(shù)據(jù),對徑向基核函數(shù)構(gòu)建的C-SVC分類結(jié)果(見表4)進行分析評價。闊葉林樣地總數(shù)為48個,其中46個正確分類,分類精度為95.83%;針葉林總數(shù)為35個,正確分類21個,分類精度為 60.00%;總體分類精度為 85.24%,kappa系數(shù)為0.426 8。說明應(yīng)用徑向基核函數(shù)構(gòu)建的CSVC對闊葉林的分類效果,明顯比對針葉林的分類效果要好。這一結(jié)果表明,確定的三個森林類型分類變量K1、ˉK、ΔK對闊葉林反應(yīng)敏感,對針葉林反應(yīng)不甚理想。造成這一分類結(jié)果的原因可能有二:其一,相對闊葉林地而言,針葉林地采樣點偏少,未能全面客觀地表現(xiàn)出分類規(guī)律,尚有待增加此類采樣點,以進一步深入研究;其二,研究區(qū)內(nèi)闊葉林地分布較廣,且很多為純闊葉林,但針葉林地內(nèi)卻常常混有一定比例的闊葉樹,影響了GLAS回波的波形特征,為針葉林的正確識別造成了一定的困難。
表4 分類驗證數(shù)據(jù)的分類結(jié)果
本文提取了GLAS波形參數(shù)(即第一個高斯曲線的斜率K1、森林冠層高斯曲線斜率均值ˉK和森林冠層高斯曲線斜率均方差ΔK)作為輸入量,應(yīng)用CSVC分類方法和十折交叉驗證方法進行森林類型識別并驗證其精度。結(jié)果顯示:應(yīng)用C-SVC分類方法能較好地對闊葉林和針葉林兩種林分進行識別,總體精度為85.24%。研究結(jié)果表明:目前構(gòu)建的分類機C-SVC所選取的三個森林類型分類變量(K1、ˉK、ΔK),對針葉林識別表現(xiàn)不理想,而且缺乏對針闊混交林的識別研究分析。因此,建議在今后的研究中補充對針闊混交林的采樣,分析提取更適合用于森林類型識別的波形特征參數(shù),結(jié)合其他分類方法對ICESat-GLAS識別森林類型的能力進行更加深入的探討。
[1] Duong V H,Lindenbergh R,Pfeifer N,et al.Single and two epoch analysis of ICESat full waveform data over forested areas[J].International Journal of Remote Sensing,2008,29(5):1453-1473.
[2] Duong H,Pfeifer N,Lindenbergh R.Analysis of repeated ICESat full waveform data:methodology and leaf-on/leaf-offcomparison[C]//Proceedings:Workshop on 3D Remote Sensing in Forestry.Available online at:http://www.rali.boku.ac.at/3drsforestry.html,2006:239-248.
[3] 巫兆聰,歐陽群東,胡忠文.應(yīng)用分水嶺變換與支持向量機的極化SAR圖像分類[J].武漢大學學報:信息科學版,2012,37(1):7-10.
[4] 張淑芬,邢艷秋,吳紅波,等.基于GIS和RS技術(shù)的木材運輸線路優(yōu)化研究:以吉林省汪清林區(qū)為例[J].森林工程,2011,27(2):48-51.
[5] 臧淑英,張策,張麗娟,等.遺傳算法優(yōu)化的支持向量機濕地遙感分類:以洪河國家級自然保護區(qū)為例[J].地理科學,2012,32(4):434-441.
[6] 蕭嶸,王繼成,張福炎.支持向量機理論綜述[J].計算機科學,2000,27(3):1-3.
[7] 李立存,張淑芬,邢艷秋.全站儀和測高儀在樹高測定上的比較分析[J].森林工程,2011,27(4):38-41.
[8] 惠文華.基于支持向量機的遙感圖像分類方法[J].地球科學與環(huán)境學報,2006,28(2):93-95.
[9] 劉志剛,史文中,李德仁,等.一種基于支撐向量機的遙感影像不完全監(jiān)督分類新方法[J].遙感學報,2005,9(4):363-373.
[10] 員永生.基于支持向量機分類的面向?qū)ο笸恋馗脖粓D像分類方法研究[D].楊凌:西北農(nóng)林科技大學,2010.
[11] 陳佳,顏學峰,鐘偉民,等.基于多項式核RVM的非線性模型預(yù)測控制[J].控制工程,2008,15(2):158-160.
[12] 林茂六,陳春雨.基于傅立葉核與徑向基核的支持向量機性能之比較[J].重慶郵電學院學報:自然科學版,2005,17(6):647-650.
[13] 朱向榮,李娜,史新元,等.最小二乘支持向量機算法與紫外光譜法用于鑒別清開靈注射液四混中間體[J].分析化學,2008,36(6):770-774.
[14] 邢艷秋,王立海.基于ICESat—GLAS完整波形的坡地森林冠層高度反演研究:以吉林長白山林區(qū)為例[J].武漢大學學報:信息科學版,2009,34(6):696-700.
[15] 邱賽,邢艷秋,李立存,等.基于小波變換的ICESAT-GlAS波形處理[J].森林工程,2012,28(5):33-35.
[16] 李俊明,邢艷秋,楊超.基于森林類型光譜特征的最佳波段選擇研究:以 HJ/1A高光譜影像為例[J].森林工程,2013,29(4):42-46.
[17] Wagner W,Ullrich A,Ducic V,et al.Gaussian decomposition and calibration of a novel small-footprint full-waveform digitising airborne laser scanner[J].ISPRS Journal of Photogrammetry and Remote Sensing,2006,60(2):100-112.
[18] Zhang J,De gier A,Xing Y,et al.Full Waveform-based analysis for forest type information derivation from large footprint spaceborne lidar data[J].Photogrammetric Engineering and Remote Sensing,2011,77(3):281-290.