基于優(yōu)化地貌特征和紋理信息的黃土高原溝緣線提取方法

羅志東 劉二佳 齊 實 姚占軍

(1.北京林業(yè)大學水土保持學院， 北京 100083； 2.水利部水土保持監(jiān)測中心， 北京 100053； 3.北京林業(yè)大學林學院， 北京 100083； 4.水土保持國家林業(yè)局重點實驗室， 北京 100083)

0 引言

黃土高原在長期流水侵蝕下，形成支離破碎、溝壑縱橫的黃土地貌景觀。溝緣線作為溝間地(正地形)和溝谷地(負地形)的分界線[1]，是切溝、沖溝最為發(fā)育的部位，影響物質(zhì)與能量傳輸過程，造成正負地貌單元在土壤類型、地表組成、地貌特征、侵蝕營力等方面產(chǎn)生巨大差異[2-3]。研究其空間分布及變化特征有助于全面分析黃土高原地貌演變情況和衡量地表侵蝕狀況[4-5]，為制定區(qū)域水土保持規(guī)劃、水土保持措施布局及土地覆被調(diào)查奠定數(shù)據(jù)基礎。

國內(nèi)溝緣線傳統(tǒng)提取方法主要是利用地形圖或遙感影像為底圖，采取人工目視勾繪方式進行，該方式受人為因素影響多、工作量大。隨著地理信息技術的發(fā)展，基于柵格DEM數(shù)據(jù)自動提取溝緣線技術成為主要研究熱點。目前，基于柵格DEM數(shù)據(jù)自動提取地貌技術的思路，主要是從水文學和地貌形態(tài)學角度出發(fā)的，基于地貌形態(tài)學基本特征單元提出了溝緣線提取方法，利用數(shù)字高程模型單元衍生的坡度、坡向、剖面曲率，匯流路徑、匯水區(qū)域和溝壑分布等地貌特征信息，建立地貌實體形態(tài)組合，判斷其閾值提取規(guī)則，實現(xiàn)溝緣線空間分布的識別[3,6-10]。從目前已有的研究結果來看，該類型方法過于依賴區(qū)域特點，參數(shù)本身的有效性和適用性以及窗口的大小均會對提取結果產(chǎn)生重要影響，且容易出現(xiàn)連續(xù)性差的問題，增大后處理工作量，因此基于此基礎的溝緣線自動提取技術還需進一步改進。另外，王軻等[11]采用地形開度和差值圖像閾值分割原理相結合的方法提取溝緣線，該方法在尋找最佳分割閾值、增強最佳分割閾值，針對不同地域溝緣線形態(tài)的適用性上還存在很大的局限性。晏實江等[12]采用邊緣檢測方法，將數(shù)字高程模型單元看作灰度圖像，通過梯度算子檢測其突變點，實現(xiàn)溝緣線的提取，該方法對灰度變化敏感但是未考慮溝緣線的地理學意義，產(chǎn)生大量隨機噪聲，部分邊緣線條斷裂，且耗時長。面向?qū)ο筇崛〔粌H能充分獲取光譜特征，還可以利用幾何和空間信息，該方法已經(jīng)在典型地形地貌類型的提取上得以應用。齊威[13]利用面向?qū)ο蟮姆椒ǎ诹炼染禈藴什罘ê妥畲竺娣e法確定最優(yōu)分割尺度因子，實現(xiàn)黃土陷穴的自動提取，提取精度達到92.3%，該方法為黃土高原地區(qū)的陷穴提取以及其他典型微地貌的提取提供了技術參考。田丹等[14]以中國1∶1 000 000 DEM為數(shù)據(jù)源，利用相關分析和熵值法篩選確定地貌分類的地形因子組合，基于隨機森林分類樹評價各地形因子的重要性，獲取面向?qū)ο蟾鲌D層的閾值，通過灰度共生矩陣紋理信息構成分類樣本知識庫進行全國地貌分類，分類結果與實際地貌分異基本吻合。SHRUTHI等[15]基于面向?qū)ο蟮姆椒ㄌ崛羡值拇笮　⑿螤詈桶l(fā)生過程，分析溝壑的變遷，認為面向?qū)ο蠓治龇ū葌鹘y(tǒng)的像素方法對地表的反映更客觀。

綜合前人研究方法中存在的問題，本文針對黃土高原溝壑區(qū)，以地表形態(tài)空間分異自然規(guī)律為有效切入點，基于地形空間形態(tài)特征分析思想，兼顧不同尺度的地貌特征，提出一種基于優(yōu)化地貌特征和深度紋理信息的面向?qū)ο蠖喑叨确指詈蜎Q策樹分類的溝緣線自動提取方法，以期進一步改善柵格技術方法的適用性與識別精度，避免地貌像元的信息孤島，優(yōu)化提取結果與效率，為推廣至大范圍宏觀地貌空間分異研究奠定基礎。

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

1.1 研究區(qū)概況

以甘肅省慶陽市西峰區(qū)的硯瓦川典型區(qū)域為研究對象。該區(qū)域?qū)儆跊芎佣壷Я鳎请]東黃土高原董志塬腹地。研究區(qū)面積為24.6 km2，地理位置107°41′15″～107°45′00″E，35°40′00″～35°42′30″N，海拔1 090～1 365 m，地形南北呈扇形，塬面完整，屬溫帶大陸性半干旱氣候，年降水量400～600 mm，且年內(nèi)分配不均，7—9月降水量約占年降水總量的60%，年平均氣溫10℃，年無霜期160～180 d，光照充足，四季分明。主要土壤類型有黑壚土、黃綿土和新積土等。植被屬暖溫性森林草原帶型，基本上無天然森林植被，陰坡生長的優(yōu)勢種有大針茅、鐵桿蒿、茵陳蒿、沙草等，覆蓋度60%～100%，陽坡生長的優(yōu)勢種有白羊草、本氏針茅、麻蒿、茵陳蒿、短花針茅、隱子草等，覆蓋度30%～90%[16]。

1.2 數(shù)據(jù)概況

采用2005年國家測繪局生產(chǎn)的1∶10 000比例尺5 m分辨率的DEM數(shù)據(jù)，數(shù)據(jù)精度滿足國家測繪行業(yè)標準(CH/T 1008—2001)。圖幅號分別為I48G008092和J49G010081，面積為26.02 km2和27.28 km2，其中I48G008092用于探究溝緣線提取方法，J49G010081用于驗證提取方法的適用性。遙感數(shù)據(jù)采用GF-1影像，影像獲取時間為2014年7月11日，數(shù)據(jù)質(zhì)量較好，無條帶和偏色。為了消除影像傳感器和大氣對地物光譜值的影響，首先利用FLAASH模型對多光譜影像進行大氣校正。然后，采用有理函數(shù)模型進行正射校正，并利用最鄰近法進行采樣。最后，采用Pansharp算法將8 m的多光譜數(shù)據(jù)和2 m的全色波段數(shù)據(jù)融合，生成分辨率為2 m的多波段圖像。結合遙感影像和DEM格網(wǎng)數(shù)據(jù)人工目視識別溝緣線，用于與自動提取結果進行對比。

2 研究方法

通過全面分析地形特征因子和紋理因子，構建充分反映溝緣線特征的最優(yōu)因子，并形成“多波段影像”，采用面向?qū)ο蠖喑叨确指罘椒ǎ浞掷蒙疃燃y理信息，通過Jeffries-Matusita距離(J-M距離)法篩選適用于本研究提取目標的多特征指數(shù)，利用決策樹分類方式，生成樹型結構的判別函數(shù)，結合專家知識，確定各特征規(guī)則的參數(shù)，建立下層節(jié)點和分支，實現(xiàn)正負地貌單元的劃分，進而提取溝緣線。

2.1 最優(yōu)地形因子構建

地形因子是地形信息的載體和重要表現(xiàn)形式之一，DEM及其衍生數(shù)據(jù)所構建的數(shù)字矩陣的結合能夠完整地反映地形在某個側(cè)面實體的形態(tài)特征，為地形特征的分析提供了更多的數(shù)據(jù)源。選取數(shù)字地形分析中常用的地形因子：高程、高程變異系數(shù)、地表切割深度、坡度、坡度變率、光照模擬值(圖1)，通過標準差和相關性對適用于地形因子提取的指標進行篩選。其中，標準差反映圖像相對于均值的離散程度，標準差越大則灰度級分布越分散，包含的信息量越趨于最大。由于不同地形因子表達的信息往往存在相關性，物理意義相似，信息具有較大的重疊，因此在溝緣線提取的過程中勢必會造成信息的冗余，需要對相關性弱的地形因子進行篩選。采用信息熵對相關度高的指標進行篩選，該值越大，攜帶的信息量越多，越有利于溝緣線的提取。DEM數(shù)據(jù)紋理特征是反映圖像表面同質(zhì)現(xiàn)象的重要信息，是不同物體差異的重要依據(jù)。研究表明，統(tǒng)計型紋理分析方法更適用于自然紋理的提取[17-18]。本文選擇統(tǒng)計型紋理中的灰度共生矩陣模型，對圖像的紋理特征進行表達。它是通過將高程數(shù)據(jù)和衍生因子進行相應的灰度域映射，統(tǒng)計柵格i和j的圖像灰度，在點對距離為d，柵格點對方向角為θ的情況下同時出現(xiàn)的頻度p(i,j,d,θ)，生成灰度共生矩陣C(d,θ)。從特征值穩(wěn)定性的角度分析，選擇d=3為共生矩陣點對距離參數(shù)，不考慮點對方向改變的影響，將0°、45°、90°、135°共4個方向分別計算的特征參數(shù)取均值。借鑒文獻[18]的研究成果，選取8個紋理特征。

圖1 DEM及其衍生因子Fig.1 DEM and its derivatives diagrams

2.2 地形分割

面向?qū)ο笏枷氲倪b感分類方法的核心技術是影像分割和特征分類算法[19-20]。本文利用分形網(wǎng)絡進化分割算法(Fractal net evolution algorithm，F(xiàn)NEA)，結合溝緣線的分布特征，基于“影像”不同的尺度、不同“波段”組合、同質(zhì)性因子，根據(jù)各地形及其衍生因子的特征產(chǎn)生最初的分割對象，構建在影像上識別和劃分的區(qū)域，選用信息熵作為內(nèi)部同質(zhì)性評價指標，分割強度作為鄰域異質(zhì)性評價指標，通過逐次調(diào)整分割尺度，實現(xiàn)整個區(qū)域的最優(yōu)自動分割[21-23]。

2.3 J-M距離可分離性特征分析

結合正負地貌的地表特征，通過J-M距離變量，獲取正負地貌單元在各特征上的J-M距離可分離性矩陣，選擇區(qū)分性最強的特征進行組合，構建描述某類地物的最優(yōu)特征空間。J-M距離法的取值范圍為[0，2]，J-M距離越大，說明兩類在f特征的可分性越強。其中，0表示兩種類別在某個特征上完全混淆，2表示兩種類別在某個特征上能夠完全分開，通常認為，J-M距離大于1.8的類別分離性較強，而小于1的類別可分性較差，可進行合并處理。計算公式為

J=2(1-e-B)

(1)

(2)

式中J——距離

B——c1、c2兩種類別基于f特征的巴氏距離

m1、m2——c1、c2兩種類別基于f特征的均值

σ1、σ2——c1、c2兩種類別基于f特征的標準差

2.4 決策樹構建和溝緣線信息提取

由于正負地貌單元的某些特征值仍存在交叉的情況，使用單一閾值建立特征規(guī)則會產(chǎn)生混淆或遺漏提取的現(xiàn)象。因此，基于選擇的對象特征的組合，利用深度紋理信息，按照決策樹CART4.5分類的方式，逐層確定各特征的函數(shù)和參數(shù)，形成正負地貌提取的規(guī)則集，識別溝緣線。決策樹信息增益指數(shù)的計算公式為

(3)

(4)

式中N——樣本的集合

Ni——根據(jù)某特征分成的子集

T——結點

pj——j類在T結點中的相對頻率

2.5 精度評價

依據(jù)實踐經(jīng)驗和專業(yè)知識，以2 m分辨率的遙感影像為數(shù)據(jù)源，通過人工識別的方法提取溝緣線，并作為精度評價的基準，評價的指標為提取的溝緣線與人工提取溝緣線的偏移距離，即偏差，計算公式為

(5)

Dmax=max(|xi-yi|)

(6)

式中i——溝緣線柵格編號

xi——編號為i樣本的提取結果

yi——編號為i樣本的目視解譯結果

n——溝緣線柵格點總數(shù)

Dmean、Dmax——絕對誤差均值和最大值

3 結果與分析

3.1 地形因子數(shù)據(jù)集構建

由表1可知，高程和光照模擬值的標準差較大，因此二者可作為后期溝緣線信息提取的主要因子。通過分析表2可知，坡度變率與坡度、坡度變率與地表切割深度之間均存在高度的相關性，因此選擇其中一種因子與其他因子組合以進行正負地貌單元的劃分。坡度、坡度變率、地表切割深度的信息熵分別為0.18、0.22、0.15。結合前人關于地貌特征提取的研究成果[3,6-10,13-14]和本研究的分析，確定了參與目標提取的地形因子為高程、光照模擬值、坡度、坡度變率。為了進一步分析空間特征的合理性，從溝緣線的上下游各選取20個樣點，統(tǒng)計篩選出各個指標值(圖2)。數(shù)據(jù)集曲線表明上述特征因子在溝緣線上下游的正負地貌單元上差異大，利用該數(shù)據(jù)集可作為溝緣線提取的有效分類特征。

表1 各因子統(tǒng)計結果Tab.1 Statistics results of each factor

表2 相關系數(shù)Tab.2 Correlation coefficient

圖2 正負地貌DEM及其衍生因子分布圖Fig.2 Positive and negative topography DEM and its derived factor distribution

3.2 地形紋理數(shù)據(jù)集構建

鑒于各紋理指標的信息具有一定的冗余程度，為避免相似物理含義和相關性強的紋理參數(shù)參與后期提取而影響提取結果，采用DEM及光照模擬數(shù)據(jù)，通過對各特征指標的分類和篩選以及相關性的分析，獲取地形的紋理特征。如表3所示，對比性與均值和、熵與二階矩、相異性與相關性的相關系數(shù)均大于0.80。表4表明熵與二階矩的相關系數(shù)也大于0.80。地形紋理的周期性是描述紋理基本屬性之一，地形形態(tài)具有一定的自相似性，在一定的尺度范圍內(nèi)，地形的特征或結構要素會以一定的規(guī)律和周期頻率重復，其中，熵和二階矩主要是對紋理的周期性進行分析的指標，從圖像上分析，對于提取溝緣線意義不大，應將這兩個因子刪除。黃土高原地區(qū)在一定的尺度上具有空間分異規(guī)律，溝緣線的提取勢必需借助于其上下游或者鄰域的空間差異性。地形特征值對比度和相異性反映鄰近柵格的差異程度，當鄰近柵格的差異度越大，越有可能區(qū)分不同的地形形態(tài)。均值和表示灰度的變化情況，與DEM信息重復，在提取過程中應該刪除。綜上所述，均質(zhì)性、方差、對比性、相關性可用于區(qū)別不同地形特征的差異性，實現(xiàn)目標地物溝緣線的提取。

表3 DEM各地形紋理相關系數(shù)Tab.3 DEM topography texture correlation coefficient

表4 光照模擬值相關系數(shù)Tab.4 Light simulation diagram correlation coefficient

3.3 溝緣線提取模型構建

利用高程、光照模擬值、坡度和坡度變率形成的最優(yōu)“多波段影像”，采用決策樹規(guī)則的面向?qū)ο蟮姆椒ǎ瑯嫿?層結構的面向?qū)ο鬀Q策樹模型，實現(xiàn)溝緣線信息提取(圖3)。通過多次試驗，最終確定的參數(shù)為分割尺度(20)、形狀(0.8)和緊致度(0.5)。節(jié)點層0為總節(jié)點；節(jié)點層1，通過高程框定正地貌和負地貌的范圍；節(jié)點層2，利用光照暈渲圖、坡度和地表切割深度排除負地貌部分的誤分；節(jié)點層3，利用均質(zhì)性、方差、對比性、相關性等紋理特征進一步排除負地貌部分的誤分。

圖3 面向?qū)ο蠓诸悰Q策樹Fig.3 Classification decision tree based on object-oriented method

3.4 精度與對比分析

為進一步探索該方法在區(qū)域上的有效性和適用性，將本算法提取的結果應用于黃土丘陵溝壑區(qū)，結果表明本文分類方法仍具有明顯優(yōu)勢。與人工解譯的溝緣線趨勢基本一致，正負地貌突變明顯的區(qū)域，二者重合度較高。在正負地貌分類圖上隨機抽取均勻分布的20個樣點，以人工目視解譯結果為參照，對比分類結果與真實地表的信息，對溝緣線提取結果進行精度評價，偏移結果在4個像元緩沖范圍內(nèi)為90%，絕對誤差均值為2～3個像元，最大誤差為4～6個像元(圖4a)。誤差大的主要原因有兩個：①DEM數(shù)據(jù)和遙感數(shù)據(jù)時相不一致，DEM數(shù)據(jù)來源于2005年，遙感影像數(shù)據(jù)來源于2014年，由于黃土高原的溯源侵蝕，塬面和梁茆區(qū)域部分懸空，繼而坍塌，溝頭前進，不同時期地形具有一定的差異。

圖4 提取結果對比Fig.4 Comparison of extraction results

②在溝緣線處坡度轉(zhuǎn)折大，坡度小的區(qū)域坡度轉(zhuǎn)折小，坡度大的區(qū)域提取效果有一定誤差。

為比較并驗證本文方法提取結果的準確度，采用一種典型的基于水文學和地貌形態(tài)學的研究方法——基于坡面的形態(tài)特征及匯水過程特點提取算法[9]進行對比試驗分析，綜合評判兩種方法提取結果的準確度及相互關系，結果如圖4b所示。從對比結果可以看出，兩種方法提取的溝緣線在空間分布以及形態(tài)特征等方面基本保持一致，宏觀表現(xiàn)上具有很高的吻合程度，均可用于研究區(qū)的溝緣線提取。采用與本文前述相同的精度分析方法進行精度評價，結果表明，偏移結果在4個像元緩沖范圍內(nèi)為70%，絕對誤差均值為4個像元，最大誤差為8個像元，相對本文方法絕對誤差均值低1～2個像元，最大誤差低2～4個像元，總體精度偏低。同時提取的溝緣線出現(xiàn)不少噪聲，椒鹽現(xiàn)象明顯，線段連續(xù)效果差，并且由于黃土地貌的復雜性，在溝道的內(nèi)部具有少數(shù)散布的偽溝緣線點，后期處理工作量大。通過兩種方法提取結果來看，本文方法較基于坡面的形態(tài)特征及匯水過程特點提取算法具有很好的總體綜合效益。

4 結論

(1)提出了一種基于優(yōu)化地貌特征和深度紋理信息的面向?qū)ο蠖喑叨确指詈蜎Q策樹分類的溝緣線自動提取方法，逐層分類提取溝緣線。

(2)提取結果相對于人工識別的溝緣線，在4個像元緩沖范圍內(nèi)為90%，絕對誤差均值為2～3個像元，最大誤差為4～6個像元，與基于坡面的形態(tài)特征及匯水過程特點的提取算法相比，絕對誤差均值降低了1～2個像元，最大誤差降低了2～4個像元。

(3)本提取方法兼顧不同地貌對象的紋理等特征信息，充分挖掘溝緣線上下游領域像元間的相關性和異質(zhì)性強的特點，有較強的抗噪能力，有效地減弱正負地貌單元混分的現(xiàn)象，解決了連續(xù)性差、噪聲多、耗時長的問題，減少了后處理的數(shù)據(jù)冗余性，實現(xiàn)了分類精度、效率和可操作性的平衡，為推廣至大范圍空間分異研究奠定了基礎。