薛傳平， 高志海， 孫斌， 李長龍， 王燕， 張媛媛

(中國林業(yè)科學(xué)研究院資源信息研究所，北京 100091)

0 引言

渾善達(dá)克沙地是內(nèi)蒙古4大沙地之一，是中國京津冀經(jīng)濟(jì)區(qū)重要的生態(tài)屏障。沙地榆樹疏林是一種極具區(qū)域特色的植被類型，由于它主要分布于溫帶地區(qū)，并且與分布于熱帶、亞熱帶地區(qū)的稀樹草原(savanna)非常相似，一些研究就稱之為溫帶savanna[1-3]。榆樹疏林是渾善達(dá)克沙地生態(tài)系統(tǒng)的重要組成部分，它的退化會直接導(dǎo)致沙地生態(tài)系統(tǒng)的退化和荒漠化的擴(kuò)張，因此保護(hù)好榆樹疏林對于防治沙地荒漠化擴(kuò)張很有必要[4]。準(zhǔn)確的榆樹疏林覆蓋數(shù)據(jù)是其他科學(xué)合理研究的前提條件，其準(zhǔn)確的監(jiān)測結(jié)果是管理決策者重要參考依據(jù)。由于常規(guī)的野外樹木調(diào)查入庫花費(fèi)人力、財力很大，尋找一種替代的收集榆樹疏林覆蓋數(shù)據(jù)的方法很有必要。近年來，遙感技術(shù)發(fā)展十分迅速，影像空間分辨率有了很大提高，高空間分辨率遙感影像的可選擇性越來越高，遙感已成為一種分析大范圍林地生態(tài)系統(tǒng)的重要工具。用高空間分辨率遙感影像提取樹冠輪廓，可有效提高調(diào)查效率，減少外業(yè)工作強(qiáng)度，對摸清榆樹疏林的規(guī)模和分布狀況，從而制定合理的保護(hù)管理對策有實(shí)際意義。但如何有效提取樹冠仍是一個難題，也是近幾年來國內(nèi)外遙感領(lǐng)域研究的熱點(diǎn)。

目前，國內(nèi)對榆樹疏林的研究多局限于生態(tài)方面[5-6]，運(yùn)用遙感技術(shù)對榆樹疏林樹冠提取的研究較少。面向地理對象影像分析技術(shù)(geographic object based image analysis，GEOBIA)是將影像分割為有意義的地理對象，利用對象的光譜、幾何和紋理等信息，發(fā)展與應(yīng)用自動化或半自動化遙感影像解譯分析的理論、方法和工具，提高影像解譯準(zhǔn)確度及效率，節(jié)省人力、物力和財力，被認(rèn)為是處理遙感中各種復(fù)雜特征識別問題的一種有效方法[7]。GEOBIA已在非洲、澳大利亞和美國等稀樹草原樹冠提取上得到很好應(yīng)用。Laliberte等使用QuickBird數(shù)據(jù)，采用多尺度分層分割和CART決策樹方法提取干旱區(qū)灌木[8]；Gibbes等使用IKONOS數(shù)據(jù)，基于面向?qū)ο蠓诸惖姆椒ㄌ崛∠∈铇淠綶9]；Boggs使用SPOT-5和QuickBird數(shù)據(jù)，利用基于像素的歸一化植被指數(shù)(normalized difference vegetation index，NDVI)閾值方法和面向?qū)ο蠓诸惙椒ǚ謩e提取研究區(qū)稀疏樹冠，結(jié)果表明面向?qū)ο蠓椒ㄒ乳撝堤崛〗Y(jié)果精度普遍高[10]；施敏燕等基于QuickBird數(shù)據(jù)，采用面向?qū)ο蠖鄬哟畏指罘椒ㄌ崛☆~濟(jì)納胡楊林樹冠信息，結(jié)果表明該方法可準(zhǔn)確示出植株的分布情況[11]。

國外關(guān)于Savanna的研究可為榆樹疏林的研究提供有益的參考價值，本文探討使用國產(chǎn)高空間分辨率數(shù)據(jù)利用遙感技術(shù)提取榆樹覆蓋的方法。該方法將很大程度降低外業(yè)調(diào)查的工作量，拓寬國產(chǎn)高空間分辨率影像在疏林草原生態(tài)系統(tǒng)上的應(yīng)用前景。

1 研究區(qū)概況及數(shù)據(jù)源

1.1 研究區(qū)概況

研究區(qū)位于內(nèi)蒙古自治區(qū)正藍(lán)旗，該區(qū)為天然沙地榆集中分布區(qū)，以“藍(lán)旗榆”著稱，其北部為渾善達(dá)克山地腹地，南部為低山丘陵，該區(qū)地理坐標(biāo)為：E115°00′～116°42′，N41°56′～43°11′。正藍(lán)旗氣候?qū)儆诟珊?、半干旱大陸性季風(fēng)氣候，年平均氣溫為1.7 ℃，年平均降水量為350 mm左右，年平均蒸發(fā)量為2 700 mm左右，呈現(xiàn)“收入少，支出多”的特點(diǎn)。研究區(qū)屬于歐亞草原區(qū)，植物物種豐富，是一個含有散生沙地榆樹林的溫帶型干旱地區(qū)，榆樹疏林群落結(jié)構(gòu)和組成均單一，很少有其他樹種。根據(jù)錫林郭勒盟森林資源統(tǒng)計，1976年正藍(lán)旗內(nèi)榆樹疏林面積達(dá)2 482 hm2[12]。沙地榆樹多呈叢狀分布，每叢3～5株，最多達(dá)20余株，根據(jù)樹齡不同可分為5 a以下幼苗、5～30 a幼樹、30～70 a成熟樹及70 a以上老齡樹，最大高度約9～10 m，最大胸徑約90 cm，最大冠幅約9～12 m[13]。沙地榆樹在每年5月初開始萌動，于10月底落葉，另外考慮到8月底—9月初是該地區(qū)打草季節(jié)，為減少草地對榆樹信息提取的影響，本研究首選影像成像時間是9月中旬—10月底或5月底—7月初。研究區(qū)地理位置如圖1所示。

圖1 研究區(qū)地理位置示意圖

1.2 數(shù)據(jù)源

1)遙感數(shù)據(jù)及預(yù)處理。本研究選取了一景GF-2數(shù)據(jù)，獲取時間為2016年9月16日。數(shù)據(jù)的預(yù)處理包括輻射校正、幾何糾正和影像融合等。輻射校正是在ENVI5.3軟件中FLAASH模型下進(jìn)行，以消除因輻射誤差而引起的影像畸變；幾何糾正是在ERDAS9.2軟件下進(jìn)行，參考影像是空間分辨率為2 m高精度正射校正后ZY-3數(shù)據(jù)，幾何糾正的絕對誤差小于1個像元。對幾何糾正后影像進(jìn)行投影轉(zhuǎn)換為WGS84 UTM Zone 50N。影像融合算法為NNDiffuse Pan Sharpening，該算法能很好地保留影像的色彩、紋理和光譜信息。GF-2衛(wèi)星傳感器具體參數(shù)如表1所示。

表1 GF-2衛(wèi)星傳感器參數(shù)

2)地面調(diào)查數(shù)據(jù)。自2011年開始連續(xù)6 a對渾善達(dá)克沙地及其周邊開展野外調(diào)查，獲得了大量野外樣地實(shí)測資料。其中，2014—2016年間重點(diǎn)對正藍(lán)旗境內(nèi)榆樹疏林進(jìn)行3次地面調(diào)查，設(shè)置大小為50 m×50 m樣地，調(diào)查內(nèi)容包括樣地內(nèi)每木經(jīng)緯度、胸徑、冠幅(包括東西、南北2個方向)和樹高，同時拍攝樣地照片作為補(bǔ)充。本研究所用的GF-2影像區(qū)域內(nèi)含8個實(shí)測榆樹樣地，140棵單木。為彌補(bǔ)GPS定位誤差，結(jié)合外業(yè)調(diào)查數(shù)據(jù)與目視解譯判讀結(jié)果，在融合后影像上目視手動勾繪150個圖斑作為參考數(shù)據(jù)，用于榆樹區(qū)粗提取過程閾值及條件的設(shè)置，其中有24個榆樹群和126棵榆樹單木。

2 研究方法

為提高識別精度和數(shù)據(jù)處理速度，本研究以融合后的GF-2影像為數(shù)據(jù)源，先建立NDVI閾值初步判斷，再采用GEOBIA方法深度分析研究區(qū)的典型沙地榆樹，這對北方干旱、半干旱地區(qū)典型樹種識別具有重要參考意義；同時對評價該地區(qū)的生態(tài)環(huán)境也具有重要的現(xiàn)實(shí)意義。

2.1 榆樹區(qū)粗提取

由于研究區(qū)榆樹大多分布在沙質(zhì)土壤上，并且分布稀疏，可用NDVI初步掩模非榆樹區(qū)。NDVI能夠很好地區(qū)分植被區(qū)與非植被區(qū)，但閾值太小會保留大面積林下植被，閾值過大會排除一些小榆樹。樹冠反射率受光合物質(zhì)葉子和非光合物質(zhì)枝干的共同影響，與其他植被類型相比，尤其是草地，它們有較明顯的光譜變化。因此，基于影像對象的空間和光譜信息提取榆樹，可通過對NDVI設(shè)置閾值生成對象，使用對象的NIR標(biāo)準(zhǔn)差作為區(qū)分樹冠的一個指標(biāo)。為防止大面積背景參與分類，結(jié)合統(tǒng)計數(shù)據(jù)，設(shè)置面積閾值大于20 m2且小于3 000 m2的限制條件。單棵榆樹或小面積榆樹群的形狀接近橢圓形，利用其圓度和橢圓度幾何特征可以很好地將其識別。NDVI計算公式為

(1)

式中：RNIR為近紅外波段的反射率；RRED為紅光波段的反射率。

人工設(shè)置NDVI閾值及條件A，粗提取榆樹分布區(qū)。NDVI閾值從最小值0.1開始，以0.05為逐步增量，直到NDVI達(dá)最大值0.4為止。只有滿足條件A，即面積介于20～3 000 m2之間、NIR標(biāo)準(zhǔn)差大于450，且橢圓度大于0.4的對象為潛在榆樹分布區(qū)P。否則，小于NDVI最小值的對象為其他地類，這里統(tǒng)一歸并為背景B，此過程是迭代進(jìn)行的。圖2中B為背景，P為潛在榆樹區(qū)，條件A為面積、NIR標(biāo)準(zhǔn)差和橢圓度3個限制條件的交集。具體榆樹區(qū)粗提取過程如圖2所示。

圖2 榆樹粗提取過程

根據(jù)統(tǒng)計目視勾繪的參考樣本特征值設(shè)置NDVI閾值及條件A。本文參考數(shù)據(jù)特征值統(tǒng)計結(jié)果見表2。

表2 參考數(shù)據(jù)特征值統(tǒng)計

①：b為NDVI均值；n為對象內(nèi)像素數(shù)量；bi為NDVI值；σ為NIR標(biāo)準(zhǔn)差；c為NIR均值；ci為NIR值；a為面積；p為像素大?。籸為圓度；s為最小包圍對象橢圓半徑；l為對象內(nèi)最大封閉橢圓半徑。

2.2 基于GEOBIA方法的榆樹區(qū)準(zhǔn)確提取

本研究使用GEOBIA方法對粗提取結(jié)果進(jìn)一步準(zhǔn)確提取。此方法從提出到現(xiàn)在10余a內(nèi)，已廣泛應(yīng)用于城市監(jiān)測、災(zāi)害監(jiān)測等很多領(lǐng)域，但在干旱區(qū)植被提取方面應(yīng)用很少[14-15]，其優(yōu)勢是以對象為處理單元克服像素級分類的椒鹽現(xiàn)象，且能夠充分利用影像對象的光譜、幾何和紋理等特征，但特征選擇和規(guī)則集構(gòu)建的不確定性是制約其發(fā)展的關(guān)鍵因素。因此，本研究基于分離閾值(separability and thresholds，SEaTH)算法優(yōu)選特征及自動構(gòu)建規(guī)則集。該算法由Nussbaum等人提出，是一種基于J-M(Jeffries-Matudita)距離評價2個類別在某特征上的關(guān)聯(lián)程度和根據(jù)高斯概率分布模型計算特征閾值的方法[16]。

1)分割方法。本文采用的分割方法為Multiresolution segmentation 多尺度算法。它是一種自下而上的區(qū)域生長算法，通過合并相鄰像元或者對象，保證對象內(nèi)部像元之間同質(zhì)性最大、異質(zhì)性最小。多尺度分割主要設(shè)置參數(shù)包括分割尺度、波段權(quán)重及勻質(zhì)性因子。分割尺度的設(shè)置至關(guān)重要，如何選擇最佳分割尺度也是很多學(xué)者研究的熱點(diǎn)，目前已出現(xiàn)了一些確定最佳分割尺度的方法或工具，如最大面積法、目標(biāo)函數(shù)法和面積比均值法等。但很多研究仍采用多次試驗(yàn)比較分割結(jié)果的方法來確定最佳分割尺度，因?yàn)殡S著分割尺度增大分割對象也增大，但是兩者之間并沒有直接關(guān)系，因此，仍需要多次嘗試來確定分割尺度。

2)特征篩選。目視查看潛在榆樹區(qū)發(fā)現(xiàn)，除榆樹外還包括灌木、草地、草本濕地、人工楊樹林和沙地5種類別。根據(jù)樣本的代表性和分布均勻性，本文選取470個分割對象作為樣本，其中灌木47個、草地52個、草本濕地78個、榆樹152個、沙地108個和人工楊樹林33個。分別統(tǒng)計其特征值，估計其概率分布，計算2種類別C1和C2之間的J-M距離來評價其可分離度[17]。J-M距離取值范圍為[0，2]，其值為0時表示2個類別在某特征上完全不能區(qū)分，其值趨向于2時表示2個類別在某特征上的分離度較大，因此選擇J-M距離較大的特征組成最優(yōu)特征集。具體計算公式為

J=2(1-e-B)

(2)

(3)

式中：B為巴氏距離；m1和m2表示類別的某特征值；σ1和σ2表示類別的某特征標(biāo)準(zhǔn)差。

影像對象常用特征包括光譜特征、幾何特征和紋理特征。本文考慮的光譜特征包括各波段及NDVI的均值、標(biāo)準(zhǔn)差、比率、最大方差及亮度等；幾何特征包括面積、長度、寬度、長寬比、圓度、橢圓適合性和矩形適合性等；紋理特征包括各個波段圖像灰度共生矩陣(gray-level co-occurrence matrix，GLCM)的勻質(zhì)性、反差、熵、均值、標(biāo)準(zhǔn)差、相關(guān)性、相性和角二矩陣等。

(4)

(5)

式中n1和n2是類別C1和C2的樣本數(shù)目。由式(4)和(5)計算出2個閾值T1和T2，選取位于m1和m2之間的閾值作為最佳閾值T。但分割時會存在一些類別之間因J-M距離值偏低而不易區(qū)分，若此時直接使用T作為最佳閾值分割結(jié)果會不準(zhǔn)確。針對這個問題，Marpu等[18]提出閾值調(diào)整規(guī)則，T′為調(diào)整后閾值，J為J-M距離值，調(diào)整規(guī)則如下：

如果J<0.5，那么忽略該特征；

如果0.5≤J≤1.25，那么T′=m2；

如果1.25

如果J≥1.75，那么T′=T。

4)分類識別方法。根據(jù)SEaTH算法篩選的特征及閾值結(jié)果對榆樹類別進(jìn)行類描述，使用eCogni-tion中規(guī)則分類的classification算法進(jìn)行榆樹分類識別，將其他地類歸為背景; 最終將影像分為榆樹和背景2個類別。

3 結(jié)果分析

3.1 榆樹區(qū)粗提取

局部榆樹粗提取結(jié)果如圖3所示，其中，藍(lán)色為背景，紅色為粗提取榆樹分布區(qū)，可以發(fā)現(xiàn)大多數(shù)非榆樹區(qū)被歸為背景類中，榆樹區(qū)被初步識別出來。目視查看榆樹粗提取效果，發(fā)現(xiàn)冠幅小的榆樹單木能夠被很好地識別，并且提取的冠幅接近其真實(shí)冠幅，而聚集分布的多株榆樹除了榆樹之外，其林下植被、沙地和灌叢等也會被粗提取出來。此外，粗提取結(jié)果還包括一些草本濕地、人工楊樹林等非榆樹區(qū)，仍需要進(jìn)一步準(zhǔn)確提取。

(a) GF-2融合影像 (b) 粗提取結(jié)果

3.2 多尺度分割

利用NDVI閾值粗提取榆樹分布區(qū)只是一個初步的提取過程，提取結(jié)果仍包括一些灌叢、草地、草本濕地、人工楊樹林及沙地等非榆樹區(qū)，需進(jìn)一步對粗提取結(jié)果進(jìn)行準(zhǔn)確提取，首先是影像分割。本文選取多尺度分割方法，分割尺度及勻質(zhì)性因子是通過多次試驗(yàn)比較分割結(jié)果來確定。由于分割過程中主要利用影像的光譜因子，所以光譜因子相對更為重要。波段權(quán)重的設(shè)置參考Ren等[19]提出的最佳指數(shù)因子(optimum index factor，OIF)，該方法是一種簡單有效篩選波段的方法，通過計算各波段的標(biāo)準(zhǔn)差和波段間的相關(guān)系數(shù)的比值選擇最佳的波段組合，其具體計算公式為

(6)

式中:Si為第i波段的標(biāo)準(zhǔn)差；Rj為任意2個波段的相關(guān)系數(shù)。OIF值越大表示其波段組合信息量越大。

通過計算GF-2影像波段間最佳指數(shù)因子篩選出最佳波段組合，進(jìn)而設(shè)置多尺度分割的波段權(quán)重。本文GF-2影像各波段間相關(guān)系數(shù)及OIF值分別見表3和表4。

表3 GF-2各波段的標(biāo)準(zhǔn)差及波段間相關(guān)系數(shù)

①：B4為近紅外波段，B3為紅光波段，B2為綠光波段，B1為藍(lán)光波段。

表4 GF-2各波段組合OIF值

需要注意分割尺度不易設(shè)置過小，防止單棵樹冠被過度分割。經(jīng)過多次嘗試，最終確定分割尺度為40，各波段權(quán)重分別為0(B1)、1(B2)、2(B3)、1(B4)、1(NDVI)，光譜因子為0.7，形狀因子為0.3，光滑度和緊密度為0.5。局部分割效果如圖4所示。

(a) 原始影像 (b) 分割效果

3.3 基于SEaTH算法的榆樹疏林提取規(guī)則集構(gòu)建

SEaTH算法可根據(jù)統(tǒng)計的樣本特征值，自動篩選出區(qū)分兩兩類別之間的最優(yōu)特征及其對應(yīng)的最佳特征閾值，然后根據(jù)篩選結(jié)果構(gòu)建提取某地類的規(guī)則。本文主要目的是提取榆樹的空間分布，不考慮其他地類，因此重點(diǎn)是構(gòu)建提取榆樹的規(guī)則。表5為榆樹與其他地類區(qū)分的規(guī)則，其中，榆樹與草地、人工楊樹林、沙地之間的J-M距離值均大于1.75，尤其是與沙地間J-M距離值大于1.9，表示榆樹與這些類別之間容易區(qū)分，特征閾值直接取閾值T。而榆樹與灌木、草本濕地之間J-M距離值約1.6，表示榆樹與灌木、草本濕地之間易混淆，此時，根據(jù)Marpu等[18]提出的閾值調(diào)整規(guī)則，特征閾值取(T+m2)/2。

表5 榆樹的提取規(guī)則

①：RatioG為綠光波段比值;StandarddeviationNIR為近紅外波段標(biāo)準(zhǔn)差;GLCMEntropy(90)為灰度共生矩陣熵特征;RatioR為紅光波段比值;HSITransformationIntensity(R=r，G=g，B=b)為HSI變換強(qiáng)度。

3.4 榆樹分類識別結(jié)果及精度評價

根據(jù)表5采用SEaTH算法構(gòu)建的榆樹準(zhǔn)確提取規(guī)則，進(jìn)一步掩模非榆樹區(qū)，將能夠同時滿足其中5條規(guī)則的對象被識別為榆樹，其他地類則為背景。

圖5為研究區(qū)局部地區(qū)榆樹覆蓋信息識別提取結(jié)果，選取2種不同生境(沙丘、塔拉)榆樹種群并對比其提取結(jié)果，其中藍(lán)色為背景，紅色為最終提取的榆樹。

(a) 沙丘生境榆樹影像 (b) 識別結(jié)果

(c) 塔拉生境榆樹影像 (d) 識別結(jié)果

本研究使用混淆矩陣對分類識別結(jié)果進(jìn)行精度評價。在榆樹識別結(jié)果上對榆樹與背景2個類別分別生成隨機(jī)點(diǎn)各300個，然后在融合后的GF-2影像上目視判讀各個點(diǎn)的類別作為參考，通過對比榆樹分類識別結(jié)果，逐一判斷這600個隨機(jī)點(diǎn)的分類正誤，構(gòu)建混淆矩陣進(jìn)行評價，結(jié)果見表6。

表6 榆樹識別精度評價

本方法提取的總體精度達(dá)88.17%，Kappa系數(shù)為0.76，其中，背景的用戶精度達(dá)99.33%，榆樹的制圖精度達(dá)99.14%。研究區(qū)榆樹種群復(fù)雜，不同的發(fā)育階段和生境對榆樹疏林識別有不同的影響。其中，冠幅較小的幼樹容易被忽略。由于與周圍背景有明顯差異，沙丘生境中榆樹(圖5(a))很容易被提取，并且準(zhǔn)確率較高，而塔拉(草原)生境中榆樹(圖5(b))由于受背景的影響，其提取精度明顯不如沙丘生境中的榆樹。此外，高大灌木與榆樹易混淆，選取的300個榆樹類別樣本中有69個為誤分，其中46個樣本實(shí)際地類為灌木，其高度約0.5～2.5 m，主要包括榆樹幼苗、紅柳、錦雞兒和柴樺等。

4 結(jié)論與討論

目前對沙地榆樹疏林的研究多是采用樣地調(diào)查、實(shí)地走訪和查閱歷史文獻(xiàn)等傳統(tǒng)方法，并且研究內(nèi)容多局限于生態(tài)方面，運(yùn)用遙感技術(shù)去探討榆樹疏林空間分布的研究很少。本文嘗試通過遙感手段使用GF-2數(shù)據(jù)采用基于GEOBIA技術(shù)提取研究區(qū)榆樹，為沙地榆樹信息提取提供有效決策支撐，獲得主要結(jié)論如下：

1)NDVI閾值可快速掩模大多數(shù)非榆樹區(qū)，準(zhǔn)確提取榆樹的潛在分布區(qū)，減少后續(xù)榆樹準(zhǔn)確提取的工作量，適用于榆樹疏林區(qū)的粗提取。

2)本文使用GF-2基于GEOBIA方法進(jìn)一步識別榆樹，具有一定的可行性。經(jīng)驗(yàn)證，其總體精度為88.17%，榆樹提取的制圖精度達(dá)99.14%，基本可以滿足榆樹疏林區(qū)識別的要求。

3)采用單景GF-2影像，結(jié)合研究區(qū)榆樹疏林的分布特性，提出先使用NDVI粗提取榆樹分布區(qū)，與直接采用GEOBIA方法提取榆樹分布區(qū)相比，可提高識別精度并加快數(shù)據(jù)處理速度。

4)考慮到特征選擇和規(guī)則集構(gòu)建的不確定性對GEOBIA方法的制約，提出基于SEaTH算法優(yōu)選特征及自動計算特征閾值，減少人為的主觀干擾，取得了較好的效果，使渾善達(dá)克沙地榆樹疏林信息提取進(jìn)入了新的階段，對于該地區(qū)榆樹疏林管理有一定的幫助。但是渾善達(dá)克沙地面積廣闊，情況復(fù)雜多樣，本研究提出的方法需要進(jìn)一步驗(yàn)證。

5)本文方法也存在一些不足，仍有待進(jìn)一步改善。例如塔拉生境中榆樹由于受背景影響，使用本文方法提取效果比較一般。下一步研究可考慮結(jié)合多個時相的遙感數(shù)據(jù)來提高塔拉生境中榆樹的識別準(zhǔn)確率，改進(jìn)本研究方法的不足。