基于多源數據的土地利用變化檢測

趙展， 夏旺， 閆利

(武漢大學測繪學院，武漢 430079)

0 引言

我國每年都進行年度土地利用變更調查，以保證我國土地利用調查數據成果具有很好的現勢性。目前的土地利用變化信息獲取還是以人工解譯遙感影像為主，任務繁重，效率較低。利用計算機自動分類技術快速獲得土地利用變化信息已成為當前的迫切需求。

根據變化檢測前是否獲得土地利用/覆蓋的類別信息，可以將變化檢測方法分為分類后比較法和直接比較法。直接比較法是直接對同一區域不同時相遙感影像的光譜和其他特征信息進行比較，利用差異信息直接確定變化的位置與范圍，常用方法有影像代數法[1]、變化矢量分析法[2-3]、多元變化檢測(multivariate alteration detection，MAD)變換法[4]以及光譜梯度差分法[5]等。在這類方法中，除了光譜特征，距離測度[6]、自相關性[7-8]和建筑物指數[9-11]等也被廣泛使用，在直接比較結果的基礎上設置閾值或進行二類分類從而獲得變化區域[2-3，12-13]。分類后比較法是先獲得前后時相的土地利用/覆蓋類別，通過對比前后類別發現變化，常用的方法有基于像元的隨機森林分類[14]、貝葉斯分類[15]、馬爾科夫隨機場集成多特征分類[16-18]，以及面向對象分類[19-20]等。分類后比較法的關鍵問題在于分類本身是一個比較繁瑣的過程，且易于將分類誤差引入變化檢測結果并放大。目前也有學者開始研究利用土地利用數據庫輔助變化檢測，如將后一時相分類結果與土地利用矢量對比發現變化[19,21]、從土地利用數據中學習規則知識[22]等。

盡管這些方法在國內外取得了較好的變化檢測效果，但往往過程較復雜，需要人工設置一定的參數或規則，并存在較多的漏提。土地利用變更調查的原則是盡量不遺漏任何疑似變化區域，同時由于作業范圍大，要求方法具有一定的自動化程度和較快的運行速度。針對這些需求，本文提出一種基于多源數據的土地利用變化自動檢測方法。該方法利用多源數據實現信息互補，將分類后比較法和直接比較法相結合，以提高變化檢測的自動化程度、效率和可靠性。

1 研究區概況與數據源

本文研究區選擇了我國南方和北方2個生態經濟區的典型區域。一個位于黑龍江省杜爾伯特蒙古族自治縣，屬于東北平原地區。該縣近些年來以“畜牧強縣”政策為主發展經濟，大力發展畜牧業，并且鼓勵牧民進城鎮生活，出現了大量新增建設用地。研究區域為該縣政府駐地的康泰鎮及其周圍區域，面積約為150 km2。另一個位于湖南省長沙縣，屬于南方丘陵地區。該縣處于長株潭“兩型社會”(資源節約型和環境友好型社會)綜合配套改革試驗區的核心地帶，近年大力發展工業、建立物流基地，同樣出現大量土地利用變化。區內包括黃興鎮和干杉鎮范圍，面積約為160 km2。本文使用的數據為國家年度變更調查標準數據，包括前、后2時相的遙感影像，以及前一時相的土地利用數據。同時采用后一時相國家土地利用變更調查中人工解譯提取的變化圖斑結果作為對比數據，該結果是作業人員在ArcGIS軟件中通過目視對比前后時相影像提取的。其中，杜爾伯特蒙古族自治縣研究區數據為2012和2013年度的數據，遙感影像是SPOT5影像(空間分辨率為2.5 m)；長沙縣研究區數據為2013和2014年度的數據，遙感影像為國產高分衛星影像(空間分辨率為2 m)。所使用的遙感影像均按照國家年度變更調查標準技術規范進行了幾何糾正和輻射校正處理，與土地利用數據精確配準。研究區域遙感影像分別如圖1和圖2所示。

(a) 2013年影像 (b) 2014年影像疊加人工解譯變化圖斑

2 基于多源數據的土地利用變化檢測

本文提出利用多源數據進行自動的土地利用變化檢測。首先,利用前一時相的土地利用數據，自動獲取后一時相影像上的訓練樣本，進行影像分類，將后一時相分類結果與前一時相土地利用數據對比，自動地發現變化區域；然后，對前、后2時相影像做MAD變換，利用MAD變換影像自動設置閾值去除偽變化；最后，通過目視判別確定提取的疑似變化圖斑是否為真實變化，對于真實變化區域采用人工繪制邊界方法提取變化圖斑。整體技術流程如圖3所示。

圖3 變化檢測與提取整體技術流程

2.1 矢量數據輔助影像自動分類

影像分類是一個比較繁瑣的過程，特征選擇、樣本提取和分類器參數設置等方面均影響分類結果，對操作人員經驗要求較高。本文基于第二次土地利用調查中農村土地利用類別體系建立分類類別體系，研究從已有土地利用矢量數據中自動獲得樣本的影像全自動分類方法，并在此基礎上實現自動的變化檢測。

2.1.1 自動分類類別體系設計

年度土地利用變更調查重點監測的是非建設用地到建設用地的變化，包括已建成區域和在建區域，其規定的變化類型與土地利用調查土地分類體系中的類別并不完全對應一致。在內業變化圖斑提取過程中，只需要發現變化區域，最重要的需求是不遺漏真正變化區域和提取效率盡可能高，最終的變化后地物類別是通過外業核查確定的。針對這一特點，本文將土地利用調查土地分類體系中的地物類別簡化為5大地類進行影像分類，將“耕地”、“園地”、“林地”和“草地”簡化為植被類；“交通運輸用地”中除“機場用地“和“港口碼頭用地”外，均簡化為道路類；“水域及水利設施用地”中除“水工建筑用地”外，均簡化為水體類；“其他土地”中除“設施農用地”外，均簡化為裸地類；“城鎮村及工礦用地”、“機場用地”、“港口碼頭用地”、“水工建筑用地”和“設施農用地”合并為建筑物類。

2.1.2 樣本的自動選擇與精化

樣本的選擇對分類結果有較大的影響。傳統方法中分類樣本是通過人工選取的，選取的結果嚴重依賴于作業員的專業知識和經驗，使分類結果受作業員的主觀性影響較大。在本文方法中，利用已有土地利用數據自動獲得訓練樣本并精化，從而實現全自動的影像分類，在提高效率的同時保證了分類結果的客觀性。

在土地利用變更調查中已有前一時相土地利用矢量數據，盡管后一時相部分區域土地利用會發生變化，但在未發生變化的區域，矢量數據的類別和影像是正確對應的。實際情況中，1 a內變化區域相對只是少數，土地利用數據中圖斑的類別信息可認為絕大部分都是準確的。因此，將前一時相土地利用數據中的所有圖斑作為初始樣本，采用迭代過程剔除少數錯誤樣本，精化初始樣本。其步驟為：首先，計算每個類別的類別特征中心，計算每個樣本與各個類別特征中心的特征距離，選擇所有特征距離最小類別與土地利用數據類別一致的樣本，組成初始的精化樣本集；然后，利用初始的精化樣本集訓練分類器對每個樣本進行分類，對比每一個樣本分類類別和在土地利用數據中的類別，將不一致的樣本剔除，得到新的精化樣本集；最后，利用新的樣本集重新訓練分類器，如此循環迭代，直至所有樣本的分類類別與土地利用數據中的類別一致，得到最終的精化樣本集，并采用訓練得到的最終分類器對影像進行分類，得到最終分類結果。方法的基本流程如圖4所示。

圖4 樣本精化流程

2.1.3 分類器設計與分類特征計算

支持向量機(support vector machine，SVM)遵循結構風險最小化的原則，具有強大的非線性和高維處理能力及泛化能力，非常適合于遙感影像的分類。大量實驗表明，在SVM的各種核函數中，徑向基函數(radial basis function,RBF)具有較好的穩定性。因此，本文選擇基于RBF的SVM作為影像分類器。通過大量實驗，設置RBF的γ參數為2.0，懲罰系數C為1.0。

本文采用的分類特征是影像的光譜特征和光譜特征的歸一化四則運算特征。歸一化四則運算特征計算方式為

f=(pixelValuei-pixelValuej)/(pixelValuei+pixelValuej)

(1)

式中：i,j表示波段，滿足i,j=1,2，...,n(n為影像波段數),且i

完成影像分類后，將分類結果與土地利用矢量對比即可獲得變化檢測結果。土地利用變更調查主要關注非建設用地變化為建設用地。因此，將類別不一致且分類類別為建筑物、道路或裸地的像素標記為變化像素，將相鄰的變化像素聚集成疑似變化圖斑。盡管在分類過程中不可避免會存在分類誤差，但本文方法的主要錯誤分類是將非建設用地誤分為建設用地，因此，變化檢測的結果一般不會出現漏提現象，滿足變更調查不遺漏變化的要求。另一方面，變化檢測結果會存在較多的偽變化，需要進一步去除，本文利用前后時相影像的MAD變換影像來去除偽變化。

2.2 基于MAD變換的偽變化去除

2.2.1 MAD變換原理

MAD變換的基本思想是認為前后時相影像在未變化的區域存在強相關，在變換區域存在弱相關。通過典型相關分析可以分析前后時相影像的整體相關性突出變化區域。

設隨機向量X=[X1,X2,...，Xp]T和Y=[Y1,Y2,...，Yp]T分別是擁有p個波段的不同時相影像。將X和Y分別做零均值歸一化，然后做線性變換，即

(2)

式中：a=[a1,a2,...,ap]∈Rp,b=[b1,b2,...,bp]∈Rp為線性變換系數矩陣。則MAD定義為

(3)

根據典型相關分析原理，2個新的變量U和V之間的相關性越強，兩者之間包含的變差信息越少；反之，相關性越弱，所包含的變差信息越多。所以在歸一化約束條件Var{aT·X}=Var{bT·Y}=1下，以Var{aT·X-bT·Y}最大化為準則，找到變換系數矩陣a和b即可完成MAD變換；然后根據以方差倒數為權值的各波段均方根可以得到最終的MAD變換影像。在MAD影像上變化區域像素一般具有較高的灰度值，而非變化區域具有較低的灰度值。

2.2.2 利用MAD變換去除偽變化

盡管現有研究表明，MAD變換相比于其他正交變換(如主成分分析變換)具有更好的突顯變化的能力[4],但當影像覆蓋范圍大、地物復雜時，MAD影像上變化區域和非變化區域灰度值區間仍然存在一定重疊，設置閾值提取變化區域時，要滿足變更調查“不遺漏變化”的要求，會造成較多的誤提偽變化。因此，本文并不直接利用MAD影像提取變化區域，而是在自動分類變化檢測基礎上，利用MAD影像剔除偽變化區域。

大范圍檢測任務下，影像上變化區域所占比例較小。因此本文通過灰度直方圖自動確定閾值，即將直方圖中一定百分比對應的灰度值作為閾值，小于閾值的像素被認為是一定沒有發生變化的像素，如果一個疑似變化區域所有像素的灰度值均小于閾值，則認為是偽變化區域，將其剔除。另外，高空間分辨率影像上地物局部細節的變化可能會造成細小的偽變化，以及前后時相存在幾何配準誤差也會造成細小的偽變化，這些偽變化可以通過設置面積閾值濾除。

2.3 變化圖斑繪制與提取

實際作業中，一般在偽變化去除后仍然會有一定的偽變化存在。這些偽變化除了來源于分類錯誤之外，更多的是客觀地物與土地利用定義不一致造成的。例如農田由于搭建農業大棚等呈現類似建筑物的特征被誤分為建筑物，鹽堿地等也會呈現類似建筑物的光譜特征；而田間的道路由于一般寬度較窄在土地利用分類中仍被歸為農田類(旱地、水澆地和水田等)的一部分。這些偽變化很難通過自動方法濾除。因此本文方法采用目視觀察剔除這類偽變化，然后通過目視對比確定變化區域的變化類型，最后人工繪制準確的變化圖斑邊界。

3 實驗與分析

采用本文方法對杜爾伯特蒙古族自治縣研究區域數據進行處理，結果如圖5所示。

(a) 類別簡化的2012年土地利用數據 (b) 2013年影像自動分類結果

(c) 初步變化檢測結果 (d) MAD變換影像

(e) 去除偽變化后的變化檢測結果 (f) 最終提取變化圖斑

首先將2012年度土地利用數據進行類別簡化(圖5(a))，然后通過樣本精化自動獲得訓練樣本，利用自動訓練的分類器對2013年影像進行分類，結果如圖5(b)。從圖中可以看出，自動分類結果存在一定的誤差，主要的誤分區域是部分亮水體和道路被誤分為建筑物。其中道路和建筑物都屬于建設用地，在土地利用變化檢測過程中，道路的誤分不會造成影響，而水體的誤分會造成誤提的偽變化，需要通過后續處理去除。另一方面值得注意的是，2013年影像成像時間在雨季，地面積水比較嚴重，在2012年影像上部分草地和灘涂裸地被水淹沒，2013影像的分類結果上表現出較多的水體，這并非真實的變化，而本文變化檢測方法中只關注新增建設用地，這一現象不會造成建設用地變化檢測的漏提或誤提。將分類結果與2012年土地利用結果對比，得到的初步的變化檢測結果如圖5(c)所示。從圖5(c)中可以看出，盡管自動分類存在一定誤差，但其錯誤主要造成偽變化，而不會漏提真實變化區域。初步提取的疑似變化圖斑數量較多(23 735個)，需要進一步去除偽變化。首先利用變化區域的面積大小去除過小的偽變化。我國的土地利用規范要求平原地區最小圖斑面積為600 m2，影像的空間分辨率為2.5 m，因此將面積閾值設為90個像元。然后計算2012年和2013年影像的MAD變換影像，如圖5(d)所示。從圖中可以看出發生變化的區域在MAD影像上亮度值都相對較高，可以通過設置閾值去除一定的偽變化。如上文所述采用累積直方圖確定閾值，本研究中累積直方圖比例設為98%。進行偽變化去除后結果如圖5(e)所示，去除偽變化后剩余疑似變化圖斑減少到512個。剩下的疑似圖斑中仍有比較多的偽變化，這些偽變化除水體外，主要是裸露的農田和農田間的道路，如上文所述，大部分偽變化是客觀地物與土地利用定義不一致造成的，只能通過目視觀察分辨。這些偽變化往往比較集中，可以通過目視觀察快速剔除。經過目視判別后，最終提取并繪制變化圖斑37個，結果如圖5(f)所示。

采用同樣的處理流程對長沙縣數據進行處理，取得類似的處理結果，該研究區域地物細碎且變化較多，最終提取變化圖斑156個，如圖6所示。

(a) 自動變化檢測結果 (b) 最終提取變化圖斑

將本文方法提取的結果與國家年度土地利用變更調查結果對比。對于杜爾伯特蒙古族自治縣，國家土地利用變更調查共提取了24個變化圖斑，本文方法提取了37個變化圖斑；對于長沙縣，國家土地利用變更調查共提取了139個變化圖斑，本文方法提取了146個變化圖斑。相比之下，本文方法能夠提取更多的變化圖斑，更不易遺漏變化區域。圖7示出本文方法與傳統方法提取結果的局部對比。

(a) 2012年影像1 (b) 2013年影像1疊加本文提取變化圖斑 (c) 2012年影像2疊加2012提取變化圖斑 (d) 2013年影像2疊加本文年土地利用信息

如圖7(a)和(b)所示影像1中，一個新增加的綠色建筑物在年度變更調查中被遺漏，而本文方法能夠將其提取出來。除此之外，年度土地利用變更調查方法只對比前后時相影像，對于歷史原因造成的土地利用數據與實際情況不一致則無法發現。如圖7(c)和(d)所示影像2中，土地利用矢量中類別是旱地，但在2012年影像上已經表現出建設用地特征，二者不一致，這是歷史原因造成的土地利用矢量與實際情況不符。由于該區域在2012和2013年影像上均表現為建設用地，僅目視對比影像認為沒有變化發生。而本文方法同時利用土地利用矢量和前后時相影像，因此能夠發現這一歷史原因造成的沒有及時更新的變化。

在實際應用效率方面，本文方法在自動變化檢測基礎上提取變化圖斑，相比于現有目視判讀整幅影像的方法能夠顯著減少工作量，提高作業效率。由于無法得知實際變更調查中研究區域數據的實際作業時間，因此本文設計了對比實驗衡量本文方法效率，由相同熟練程度的操作人員按照年度土地利用變更技術規程要求目視解譯提取研究區域變換圖斑，將其效率與本文方法對比。對于杜爾伯特蒙古族自治縣研究區，本文方法所需時間約為23 min，而目視解譯作業時間約為50 min；對于長沙縣研究區，本文方法處理所需時間約為70 min，而目視解譯作業時間約為2.5 h。與目視判別相比，本文方法所需時間縮短了一半。

4 結語

針對土地利用變更調查提出了一種基于多源數據的土地利用變化檢測方法，能夠實現自動化、高效率的土地利用變化檢測與提取。研究取得如下結論：

1)利用已有前一時相土地利用數據自動獲得后一時相影像的訓練樣本，并通過迭代分類過程自動剔除錯誤樣本，能夠實現全自動的影像分類，在提高效率的同時保證了分類結果的客觀性；在此基礎上與前一時相土地利用數據對比可以實現全自動的變化檢測。

2)對前、后時相影像做MAD變換，利用MAD變換影像可以通過設置簡單閾值有效地剔除大部分偽變化區域，提高工作效率。

3)研究證明相比于現有作業方法，本文方法的自動化程度更高，作業速度提高一倍以上，且更不易于遺漏真實變化區域。

4)本文方法的不足之處在于變化圖斑的最終提取方面還需要一定的人工干預，需要人工判斷變化類型以及繪制變化圖斑準確邊界，如何進一步提高這方面的自動化程度，是下一步研究的方向。