基于無人機影像天山云杉林主伐跡地提取研究

王雅佩，王振錫，劉夢婷，李 擎，楊勇強

(新疆農業大學林學與園藝學院，新疆教育廳干旱區林業生態與產業技術重點實驗室，烏魯木齊 830052)

0 引 言

【研究意義】天山云杉(PiceaSchrenkianavartianshanica)是新疆山地森林的主要優勢針葉樹種，是北方森林重要的組成部分，天山云杉林在水源涵養、水土保持以及維持新疆綠洲生態平衡等方面發揮著不可替代的重要生態功能[1-2]。天然林資源及其更新對森林的演變及發展具有重要的作用，是長期以來森林生態系統主要研究領域之一，這個過程中其受環境條件、自然干擾和人為干擾類型等因素及其相互作用的影響[3-4]。在林業生產中，森林主伐及其更新對森林資源的有效利用和恢復有著不可分割的作用，并且與森林的結構、功能和動態之間存在著密切的關系[5]。該方法耗費大量的人力物力財力、周期長、時效性低，特別是對主伐跡地而言，采伐地點分散并且不明確，致使野外調查時間延長、工作量增加，不足以滿足當代林業發展的需要[6-7]。遙感技術是森林資源高效快速監測調查的重要技術手段。無人機遙感作為新興遙感技術，能夠在短時間內獲取所需數據，通過自主設置航線更加靈活獲取遙感數據，具有環境適應力強、分辨率高、低成本等特點，使無人機遙感技術成為獲取高分辨率遙感數據的重要技術手段之一[8]。開展森林資源調查和統計對及時獲取林木的生長狀況，和森林未保育有實際意義。【前人研究進展】利用無人機影像對地物信息提取的應用研究已越來越廣泛，由于其高分辨率的特性，使其分類方法主要是以面向對象方法為主。Chenari A[9]基于無人機正射影像，使用面向對象分類方法提取研究區的野生開心果和野生杏仁，并估算平均樹冠面積，結果表明該方法對野生開心果估計精度為0.92，野生杏仁估計精度為0.89。吳金勝等[10]基于無人機影像，采用面向對象對估算研究區水稻面積，結果表明總體分類精度達到93%以上，該方法能夠替代目視解譯法，有效提高了遙感面積估算效率。彭佳憶等[11]基于無人機可見光影像采用面向對象方法提取荒漠地表類型信息，結果表明面向對象方法精度優于最大似然法，面向對象多尺度分割方法提取荒漠地表類型信息具有高精度。【本研究切入點】利用中高分辨率遙感影像在采伐監測中的研究有很多，但基于無人機可見光影像對山區主伐跡地提取及其分布研究方面還未見報道。研究基于無人機可見光影像，采用面向對象多尺度分割方法，利用最鄰近分類方法與SEaTH算法相結合的分類方法提取不同等級的主伐跡地的面積。【擬解決的關鍵問題】基于無人機可將光影像，利用eCognition 9.0軟件，采用面向對象多尺度分割技術，確定各層次最優分割尺度，通過最鄰近分類以及SEaTH算法相結合的分類方法篩選分類特征，并建立類別層次結構及分類規則集，以實現提取不同等級主伐跡地面積，為新疆天保工程實施后山區森林資源更新恢復及評價提供科學依據。

1 材料與方法

1.1 材 料

1.1.1 研究區概況

研究區位于天山山脈中段中山帶的新疆農業大學實習林場內，區域地處天山北麓烏魯木齊西南方向110 km處，位于頭屯河上游，地理坐標范圍為43°16′～43°26′N，86°46′～86°57′E。林場年平均降水量約600 mm，其海拔高度在1 700～3 200 m間。經營面積為1.044 4×104hm2，其中林業用地面積為5.116×103hm2。地形多樣，南部地勢高，北部地勢低，地形切割較為劇烈，坡度主要在10°～40°間。土壤發育較好，林下土壤為普通灰褐色森林土。林區以天山云杉純林為主。

1.2.1 影像數據

采用的無人機搭載的相機型號是SONY DSC-RX1RM2(35 mm)相機，搭載鏡頭焦距為35 mm。為了避免夏季背景植被(草本、灌木等)的影響，無人機數據拍攝于2017年11月15日上午，天氣晴朗，能見度高，風速較小。無人機起飛前，預先設定好飛行的方向、航線和高度載入飛行控制系統，此次飛行路線的航向重疊為80%，旁向重疊為60%，飛行高度為923 m。獲取的無人機影像地面分辨率為0.118 m。

1.2.2 調查數據

數據包括新疆農業大學實習林場二類調查數據，樣地調查數據以及林場伐區設計資料等其他輔助資料。2017年8月，在新疆農業大學實習林場內進行實地調查。根據天然林保護工程實施前的伐區設計材料確定主伐跡地，在主伐跡地范圍內選取典型樣地，設置30 m×30 m的典型樣方，對樣地內天山云杉進行每木檢尺，實測并記錄樣方內全部天山云杉的樹高、胸徑、年齡、冠幅、株數等信息以及樣地的郁閉度。利用亞米級的集思寶GPS記錄樣方的中心點及4個角點的位置信息，以便于將實測樣地與無人機影像進行一一對應。

1.2 方 法

1.2.1 無人機影像處理

采用PIX4D軟件對無人機數據進行處理，主要過程包括：(1)首先篩掉航拍中成像效果差、模糊等不合格的照片；(2)建立測區，導入有效的照片、POS數據和控制點文件，設置參數后提交，軟件會自動提取照片及相機的相關信息；(3)數據的自動拼接，調整拼接線、投影切換和混合影像，最終生成正射影像(DOM)。

1.2.2 多尺度分割

在eCognition Developer 9.0軟件，采用面向對象的方法提取不同等級的主伐跡地。面向對象信息提取中影像最優分割尺度的確定是關鍵技術之一。采用ESP算法[12]對最優分割尺度進行計算。ESP算法完成后，得到橫坐標為尺度參數，縱坐標分別為局部方差(local variance，LV)和局部方差變化率(rate of change，ROC)的折線圖。當ROC呈波峰狀態即出現局部最大值時，則表示此尺度為某種地物對象的最佳分割尺度。

在多尺度分割技術中，影像分割結果的質量與均質性因子的設置密切相關。均質性因子包括顏色因子(Color)和形狀因子(Shape)，其中形狀因子又包含平滑度(Smoothness)以及緊致度因子(Compactness)。其中顏色與形狀因子值是相對的，權值和為1；平滑度與緊致度因子權值和為1。顏色因子提供了對象的光譜信息，形狀因子反映了分割對象的大小，平滑度因子是優化分割對象邊界的平滑程度，緊致度因子使用形狀標準來考慮整體緊致性，并優化影像對象結果的緊湊度。

1.2.3 分類

最鄰近分類算法是一種監督分類方法，是通過選擇訓練樣本，構建特征空間，計算相關特征之間的距離，根據距離確定對象所屬類別[13]，進行分類。如果距離值越小，則表示該對象越可能隸屬于該樣本。為了避免由于采用過多參與分類的特征所帶來的計算量大、特征冗余和分類低等問題，eCognition軟件提供了特征空間優化的模塊(Feature Space Optimazation)，選取預先選擇的特征，通過該模塊計算得到特征優選的結果，可以確定在幾維時類別與類別之間的可分性最大，特征優選的結果可以直接應用于分類體系中執行分類。計算公式如下：

(1)

表1 初選特征因子
Table 1 Primary selection factor

特征Factor特征因子Texturefactor特征Factor特征因子Texturefactor紋理特征Texturefactor均值(Mean)方差(Variance)熵(Entropy)對比度(Contrast)均值性(Homogeneity)相異性(Dissimilarity)二階距(SecondMoment)相關性(Correlation)光譜特征Spectralfactor形狀因子Shapefactor平均值(Mean)亮度(Brightness)標準差(StdDev)波段比(Ratio)模擬植被指數(VI')長寬比(Length/Width)形狀指數(ShapeIndex)

其中，VDVI表示可見光波段差異植被指數，VI' 表示模擬歸一化植被指數提取植被的一個參數，公式如下：

(2)

VI'=(2×G'-R'-B')-(1.4×R'-G') .

(3)

SEaTH(Seperability and Thresholds)算法[14]是基于各類別的樣本特征值，利用J-M(Jeffries-Matusita)距離實現兩兩類別間最優特征的自動選擇，及其最優特征的閾值確定。SEaTH算法主要包括兩部分：

(1)對區分兩兩類別特征進行優選

首先，選擇各個類別的影像對象代表性樣本，選擇預先考慮到的特征，并導出樣本的特征值。SEaTH算法采用J-M距離計算兩兩類別之間在某個特征的分離度，以評估兩個類別的可分離性。計算J-M距離公式如下：

J=2(1-e-B) .

(4)

(5)

(2)最優特征閾值的確定

該算法根據高斯概率分布公式計算兩兩類別在某一個特征的最佳閾值。其公式如下：

p(x)=p(x|C1)p(C1)+p(x|C2)p(C2)

(6)

當p(x|C1)=p(x|C2)時，這兩種類別間的混分最少。C1、C2表示兩個類別，當特征值為X1時，C1和C2這兩個類別的分離性最好，此時X1即為最佳特征閾值T，計算方式如下：

(7)

(8)

式中，T為該特征的閾值，T介于m1和m2之間，m1和m2分別為兩個類別在某個特征上的均值，σ1σ2和σ2σ2分別為兩個類別在某個特征上的標準差，n1和n2為類別C1和C2的樣本數量。圖1

圖1 C1和C2最佳分離特征閾值
Fig.1 The schematic diagram of separation threshold between C1 and C2

1.2.4 分類精度評價

研究選取生產者精度、用戶精度、Kappa系數和總體精度四個指標對不同等級主伐跡地提取的分類結果進行精度評價。

2 結果與分析

2.1 主伐跡地等級劃分

根據實地主伐跡地調查數據包括林分平均胸徑、平均樹高、年齡及林分密度均和目視解譯法勾繪獲得的樣地郁閉度及平均冠幅，利用R軟件，采用中間距離法進行聚類分析，將13塊典型樣地劃分為三個等級。主伐跡地天山云杉林的林分密度隨著主伐跡地等級的增加逐漸提高。主伐跡地林木相較于周圍林木樹冠較小，且數量相對較多，且林木排列較為規律伴隨著周圍有天然更新的林木。不同等級的主伐跡地在圖像上表現為等級Ⅰ主伐跡地相對林木樹冠小，林木株數少以及裸地面積大；等級Ⅱ主伐跡地較等級Ⅰ的樹冠大小適中，林木株數較多；等級Ⅲ主伐跡地樹冠林木樹冠較大，林木株數數量多以及裸地面積較小，且林木大小不一。表2,圖2,圖3

表2 主伐跡地分級及指標
Table 2 Classification and indicators of main cutting areas

等級Grade平均胸徑AverageDBH(cm)平均樹高Averageheight(m)平均年齡Averageage(a)平均冠幅Averagecrownwidth(m)郁閉度Crowndensity林分密度(株/hm2)Standdensity(plant/hm2)Ⅰ5.33～5.883.66～4.3615～221.68～1.860.20～0.54833～2100Ⅱ1.81～6.591.78～5.7517～261.01～2.480.28～0.872267～3367Ⅲ6.67～9.566.93～9.7326～291.13～1.990.41～0.653833～4632

圖2 聚類結果
Fig.2 Clustering result graph

注:(a)(b)為Ⅰ等級，(c)(d)為Ⅱ等級，(e)(f)為Ⅲ等級

Annotation：(a) (b) is grade I, (c) (d) is grade II, (e) (f) is grade III

圖3 不同等級主伐跡地部分樣地示意
Fig.3 Images of parts of main cutting sites of different grades

2.2 最優分割尺度的確定

在eCognition軟件支持下，采用ESP算法對研究區無人機影像數據進行最佳分割尺度計算。一次分割主要目的是為了區分植被和非植被，本次ESP處理計算結果中最佳分割尺度是43、82，113、188、245、274、316、343、404、467和504，使用這些最佳分割尺度分別對影像進行分割，通過觀察其視覺效果，當分割尺度為43時，不論是道路、地面、建筑物還是林地等地物都被分割的十分破碎，未能達到此次分割的目的；當分割尺度為113時，發現不同地物均有合并現象，但分割效果仍然十分破碎；當分割尺度在113至245之間時，發現地面以及林地地物逐漸被合并，而道路和建筑物地物合并效果并不明顯；當分割尺度在245之后開始，可以發現道路地物逐漸被合并，林地地物中也有部分區域進行進一步的合并；當分割尺度為343時，道路、地面、建筑物和林地地物中基本完成了對象的合并，此時分割的效果較為接近第一次分割的目的；當分割尺度為504時，道路地物的線狀信息表達較為清晰，地面以及林地地物中的對象大面積的合并，幾乎不存在破碎現象，達到此次分割目的，因此將分割尺度為504作為第一次分割的最優分割尺度，通過不同形狀因子以及緊湊度因子的測試，最終確定形狀因子為0.4，緊致度為0.5。

在eCognition軟件支持下，采用ESP算法進行二次分割，其目的為在植被類別中提取主伐跡地，經過對比最后確定提取主伐跡地的最優分割尺度為292。對形狀因子以及緊致度因子參數進行調節，最終確定形狀因子以及緊致度因子分別為0.5、0.3。圖4、圖5,表3

圖 4 一次分割ESP最佳分割尺度
Fig. 4 The optimal segmentation scale line chart in the first Segmentation 表3 多層次分割參數Table 3 Scale parameters of the multi-level segmentation

層次等級Hierarchicallevel目標地類Targetlocationclass背景地類Backgroundgroundclass分割尺度Divisionscale形狀因子Shapefactor緊致度因子CompactnessfactorLevel1植被Vegetation非植被Non-vegetation5040.40.5Level2主伐跡地Finalfellingarea非主伐跡地Notfinalfellingarea2920.50.3

圖5 一次分割不同分割尺度對比
Fig.5 The comparison chart of different segmentation scales in the first Segmentation

2.3 面向對象主伐跡地提取

2.3.1 level 1分類

在level 1中，主要區分植被與非植被，提取目標地類為植被。由于研究所使用的影像為無人機可見光影像，通過查閱文獻最終確定VI' 和VDVI兩個植被提取指標，經過多次試驗對比觀察，最終以VI' 大于等于-0.14將植被與非植被區分開。圖6

圖6 level 1層分類結果
Fig.6 The classification result of level 1

2.3.2 level 2分類

在level 2中，對level 1植被類別進行繼承，在對植被類別二次分割前對一次分割結果進行合并，主要是在植被類別中提取主伐跡地，由于主伐跡地是2000年以前進行主伐過后的，相對于其他未主伐過的林地，林木相對于較整齊且樹冠相對于較小，且單元面積內的林木株數相對于較多。通過SEaTH算法，計算主伐跡地和非主伐跡地類別間各種特征的J-M距離及其閾值的確定，選擇J-M距離最大的前5個特征進行分類。表4、圖7

圖7 level 2層分類結果
Fig.7 The classification result of level 2表4 目標地類與背景地類特征J-M距離值及其閾值Table 4 The J-M and threshold value of features between target land covers and background land covers

目標地類Targetlocationclass背景地類Backgroundgroundclass特征FactorJ-M距離J-Mdistance閾值Threshold主伐跡地Finalfellingarea非主伐跡地NotfinalfellingareaHIS顏色特征亮度(Intensity,I)0.97≤0.497R波段均值(Mean_R)0.97≥88.5R波段波段比(Ratio_R)0.93≤0.167灰度共生矩陣均值(GLCM_Mean)0.83≤122.4亮度(Brightness)0.83≤112.79

2.3.3 Ⅰ、Ⅱ、Ⅲ等級的主伐跡地提取

在level 2層次提取的主伐跡地類別中，由于SEaTH算法適用于兩兩類別間的區分，而此次是為了提取Ⅰ、Ⅱ、Ⅲ等級的主伐跡地，根據主伐跡地等級，采用之前篩選的遙感特征因子與SEATH算法計算每兩個等級間的J-M距離最大的特征因子相結合的方式，最終確定提取各等級主伐跡地的特征因子，經過多次試驗對比觀察確定閾值，建立其分類規則，提取Ⅰ、Ⅱ、Ⅲ等級的主伐跡地。表5、圖8

表5 不同等級的主伐跡地分類規則
Table 5 Classification rules for main cutting areas of different grades

目標地類Targetlocationclass特征Factor規則表達式RegularexpressionⅠ色飽和度(Saturation,S)≤0.7476亮度(Brightness)≥133.313G波段平均值(MeanG)≥102.1649ⅡR波段平均值(MeanR)≤157.2225灰度共生矩陣均值(GLCMMean)>141.6125R波段波段比(RatioR)≤0.2273ⅢHIS顏色特征亮度(Intensity,I)≥0.4595亮度(Brightness)≤117.1661灰度共生矩陣方差(GLCMVariance)>21.437

圖8 不同等級的主伐跡地分類結果
Fig.8 Classification results of main cutting areas of different grades

2.4 Ⅰ、Ⅱ、Ⅲ等級的主伐跡地面積精度評價及其分布

所使用無人機影像分辨率較高，達到了0.118 m，因可根據目視解譯結果對分類結果進行驗證，并提取不同等級的主伐跡地面積與目視解譯結果作對比。利用eCognition軟件中提供的基于樣本的誤差矩陣的精度評價方法，通過選取目視解譯與分類結果圖像中對應位置的提取類別進行比較計算，采用混淆矩陣的統計方法對分類結果進行分析以及對比。

從生產者精度來看，主伐跡地Ⅰ、Ⅱ和Ⅲ精度分別為80.77%、83.33%和64.10%；從用戶精度來看，主伐跡地Ⅰ、Ⅱ和Ⅲ精度分別為84.00%、80.81%和71.43%。從總體看，利用最鄰近分類與SEaTH算法相結合的分類方法總體分類精度達到81.82%，Kappa系數為0.74，實驗分類質量取得了較好的效果。表6

表6 分類精度評價結果
Table 6 Classification accuracy evaluation results

地類CategoryⅠⅡⅢ非主伐跡地NotfinalfellingareaⅠ42602Ⅱ48096Ⅲ06254非主伐跡地notfinalfellingarea64587列總和Columnsum52963999生產者精度ProducerAccuracy,PA80.77%83.33%64.10%87.88%用戶精度UserAccuracy,UA84.00%80.81%71.43%85.29%總體精度OverallAccuracy,OA81.82%Kappa系數Kappacoefficient0.74

對于不同等級的主伐跡地面積的提取結果與目視解譯結果對比，主伐跡地Ⅰ的面積相對誤差為12.91%，面積吻合度為87.09%；主伐跡地Ⅱ的面積相對誤差為20.14%，面積吻合度為79.86%；主伐跡地Ⅲ的面積相對誤差為33.67%，面積吻合度為66.33%；提取主伐跡地總面積為72.574 3 hm2，目視解譯面積為92.174 9 hm2，面積相對誤差為21.26%，面積吻合度為78.74%。該方法能夠提取研究區主伐跡地。表7

表7 不同等級主伐跡地面積精度
Table 7 Precision scale of main cutting area of different grades

等級Grade提取面積Extractionarea(hm2)目視解譯面積Visualinterpretationarea(hm2)面積相對誤差ArearelativeError(%)面積吻合度Areacoincidencedegree(%)Ⅰ11.086112.728912.9187.09Ⅱ51.906865.000620.1479.86Ⅲ9.581414.445333.6766.33合計Total72.574392.174921.2678.74

3 討 論

在面向對象分類多尺度分割技術中，最優尺度直接決定影像對象的大小及信息提取的精度[15]。通常情況下最優尺度的選擇需要反復多次對分割尺度進行人工調整，然后通過目視判斷法確定，但運用EPS算法可有效減少嘗試次數，能夠準確、快速的確定最優尺度。

在分類特征優化方面，研究采用最鄰近和SEaTH算法相結合的分類方法提取研究區內不同等級的主伐跡地，其總體分類精度達到了81.82%，取得了相對較好的分類效果，可有效避免了分類特征較多而導致的特征冗余、計算量大和分類精度不高等問題，并為最優特征及其閾值的確定提供了較為便捷、有效的途徑。這與劉鑫[16]2017年基于無人機影像提取四川低丘區耕地信息的研究結論一致。

在分類精度方面，不同等級(Ⅰ、Ⅱ、Ⅲ)的主伐跡地均存在不同程度的漏分和錯分情況，主要原因可能是不同等級(Ⅰ、Ⅱ、Ⅲ)主伐跡地之間的圖像閾值具有一定的重疊度，同時本研究使用的無人機圖像僅是可見光范圍的地物成像，植被信息提取所需的光譜信息量不足。在今后的研究中可以嘗試采用無人機多光譜圖像，持續改進和更新分類算法開展研究，進一步提高分類精度。

從不同等級主伐跡地面積提取精度來看，主伐跡地Ⅰ提取面積吻合度最高，這是由于主伐跡地Ⅰ內的天山云杉林分密度低于主伐跡地Ⅱ和Ⅲ，其林分平均冠幅、平均年齡等因子也都相對較小，因此在影像上更好區分。

4 結 論

4.1 通過調查，研究區主伐跡地分為三個等級，HIS顏色特征亮度(Intensity，I)、R波段均值(Mean_R)、R波段波段比(Ratio_R)、灰度共生矩陣均值(GLCM_Mean)和亮度(Brightness)是主伐跡地提取的重要特征。色飽和度(Saturation，S)、亮度(Brightness)和G波段平均值(Mean G)是主伐跡地Ⅰ提取的重要特征；R波段平均值(Mean R)、灰度共生矩陣均值(GLCM Mean)和R波段波段比(Ratio R)是主伐跡地Ⅱ提取的重要特征；HIS顏色特征亮度(Intensity，I)、亮度(Brightness)和灰度共生矩陣方差(GLCM Variance)是主伐跡地Ⅲ提取的重要特征。

4.2 基于面向對象多尺度分割技術，采用最鄰近以及SEaTH算法相結合的分類方法，借助eCognition軟件獲取最優特征空間，可有效避免“椒鹽”現象，總體分類精度為81.82%，該方法在不同等級主伐跡地分類中具有較高的分類精度，相比僅用單一方法能夠更便捷、更有效的建立類別層次結構及分類規則集，是一種行之有效的方法。