土地利用變化圖斑提取方法

夏 旺，趙 展

（1.武漢大學 測繪學院，湖北 武漢430079）

土地利用變化圖斑提取方法

夏 旺1，趙 展1

（1.武漢大學 測繪學院，湖北 武漢430079）

為提高變化圖斑提取效率，減少人工操作難度，在多元變化檢測（MAD）自動發現變化信息的基礎上，提取感興趣區域生成差異影像，并分割差異影像得到完整變化圖斑。利用3種典型變化地物類型進行圖斑提取實驗，結果表明，該方法可以獲取完整﹑連通的變化圖斑，在錯檢率和漏檢率方面優于MAD，且Kappa系數分別從0.770﹑0.810﹑0.729提高到了0.916﹑0.894﹑0.934。

變化圖斑；分水嶺分割；區域合并；變化檢測

獲取兩幅不同時相遙感影像的變化圖斑是土地利用遙感監測的關鍵[1]，但目前通常仍是采用人工勾繪的方式來獲取變化圖斑，我國每年土地變化信息量十分巨大，手工提取變化圖斑耗時耗力，且準確勾繪邊界線也對操作人員的技術水平提出了較高的要求，因此減少人工操作難度具有重要意義。

利用變化檢測方法可以自動識別不同時相遙感影像的變化區域。根據處理單元可將檢測算法分為像素級和對象級兩類[2-3]。對象級變化檢測是將一組有特定關聯的像素集合為一個檢測單元來進行變化檢測，可對上下文信息建模，有效考慮局部鄰接像素之間的關系，能較好地去除由光譜差異或配準誤差造成的偽變化區域[4]，但該算法依賴于分割與分類算法，細小的變化區域可能被湮沒在背景中而漏檢[5]，同時由于光照﹑視角﹑大氣條件等影響，前后兩個時相影像的分割與分類結果難以保持一致，將導致大量偽變化區域[6]，因此在土地覆蓋﹑土地利用調查更新等多個應用領域中，仍廣泛使用像素級變化檢測。像素級變化檢測最大的優點是可以全面檢測影像中任意形狀﹑大小的變化信息，達到很好的效果；但由于算法沒有考慮鄰接像素的信息，容易對由光譜差異和配準造成的誤差敏感，產生離散﹑破碎﹑不連通的檢測結果[7]，且存在大量偽變化區域，不能直接利用變化信息提取變化圖斑。

針對變化檢測算法結果離散﹑破碎﹑不連通的缺點，本文首先利用MAD結果輔助人工選取感興趣區域（ROI），再利用ROI切割前后時相影像獲得影像切片，最后對切片進行歸一化再相減獲得差異影像。差異影像每個像素的像素值表征該像素的變化程度。利用分水嶺算法分割差異影像，并針對分水嶺算法過分割問題，設計了一種二段的區域合并方法。

1 基于MAD的差異影像獲取

MAD是基于典型相關分析進行變化區域檢測的方法。為了有效集中差異信息，提高檢測精度，MAD采用計算一對典型變量并相減的方式，既可在最大程度上消除不同時相﹑不同通道相關信息的影響，又可將所有變化信息分配到互相獨立的MAD變量中去[8]。

對MAD變量進行后處理和閾值化[9-10]可得到一 幅二值圖像。由于MAD結果中存在大量偽變化區域，不能直接用于獲取ROI。本文通過人工目視觀測MAD提取的變化結果，確定所提取的變化區域是否為真實變化，并選取包含變化區域的外圍矩形為ROI，利用ROI切割前后時相影像，獲得影像切片Cut1與Cut2。

雖然用于變化檢測的影像一般會先經過輻射畸變校正，但輻射畸變并不能被完全消除，且由于兩個時相影像可能來自不同傳感器，獲取時刻光照強度也不完全相同，兩個影像仍存在輻射差異。為了盡可能消除輻射差異的影響，本文先將Cut1與Cut2所有波段分別歸一化至0～255，歸一化后的切片分別記為NCut1與NCut2，再將NCut1與NCut2對應波段相減并取絕對值獲取差異影像Dif（Dif的第j個波段記為Difj）：

2 差異影像的分水嶺分割

本文利用分水嶺算法對差異影像進行分割，分水嶺算法對邊緣具有良好的響應，可得到封閉連續邊緣。差異影像通道數與原始影像相同，往往是多通道。本文對差異影像的每個波段分別作分水嶺分割，再將不同波段的分割結果合并在一個波段上。

對于Difj，先利用Sobel算子求得梯度影像Gj。為了便于分水嶺算法操作，需將像素值歸一化至0～255。分水嶺分割的關鍵在于種子點的選取，本文將Gj所有的局部最小值作為種子點，即遍歷Gj每個像素，判斷該像素在其3×3的鄰域內是否為唯一極小值，若是，則該點為種子點。獲取種子點后就可對Gj作分水嶺分割[11]，得到分割結果Segj。Segj是一個標記影像，將Segj像素p的標記值記為Segj(p)。

對Dif各波段都求取相應的分割結果后，需將所有分割結果合并為一個新的分割結果影像，記為Seg，它保留了各波段的分割結果，這是由于變化信息可能分布在不同波譜上，所以需在合并時保留所有變化信息。合并的具體步驟為：

1）將Seg所有像素設為0，標記值記為label，初始label記為1；

2）順序遍歷圖像，當前遍歷像素為p，若p=0，則標記為label；

3）遍歷p的4個鄰接像素，令當前遍歷的鄰接像素為s，若s≠0，且對于任意一個波段j都滿足Segj(s)=Segj(p)，則s被標記為label，并將s壓入堆棧；

4）從堆棧中取出一個像素作為p返回步驟3）；

5）當堆棧為空時，label加1，返回步驟2）。

3 分割區域合并

分水嶺算法對噪聲敏感，同時由于種子點的選取非常密集，分割結果會出現過分割現象。針對這兩個缺點，本文設計了一種二段區域合并的方法。

對Seg的每個分割區域計算區域面積（區域像素個數）以及各波段的像素平均值。若影像有N個波段，則第i個分割區域各波段的像素平均值為meani1, meani2,…, meaniN。

設區域i與區域j的光譜差異RegionDifij為：

若RegionDifij較大，則區域i與區域j的光譜差異較大。

3.1 根據面積合并

本文選取的種子點較為密集，分割結果將存在嚴重的過分割現象。分水嶺算法對誤差十分敏感，由于小面積區域的像素個數較少，噪聲對分割區域像素平均值的影響較大，所以需先將小面積區域與鄰接區域合并，若該區域存在多個鄰接區域，則選擇與它光譜差異最小的鄰接區域合并。具體算法為：遍歷Seg的所有分割區域，令當前遍歷的分割區域為Reg，若Reg的像素個數小于某一閾值A，則根據式（2）計算該區域與它所有鄰接區域的光譜差異RegionDif，選取RegionDif最小的鄰接區域與Reg合并。閾值A的大小與影像大小有關，本文實驗所用A的大小為Seg的像素總數除以500。

3.2 根據光譜差異合并

通過面積合并后，誤差都淹沒在較大面積的分割區域中，所以對分割區域像素平均值影響較小。此時可計算每個分割區域與各自鄰接區域的光譜差異性，若光譜差異較大，則代表兩個區域不屬于同一類別；若較小，則將兩個區域合并。可利用RegionDifij大小判斷區域i與區域j是否應被合并為一個區域，需尋找合適的閾值B，若RegionDifij＜B，則將區域i與區域j合并為一個區域。

不同的差異影像的變化區域與非變化區域的反差并不相同，所以每個Seg的最佳合并閾值也不相同。本文根據K均值聚類算法的思想，尋找自適應的閾值B。先利用式（2）計算面積合并后剩下的所有區域與其鄰接區域的RegionDif，得到所有區域與各自鄰接區域RegionDif的集合，再利用K均值聚類算法處理RegionDif的集合。本文設K=2，將RegionDif的集合分為兩類：一類中的RegionDif值較小，代表該類中RegionDif值對應的兩個分割區域光譜差異較小，應合并為一個區域，將該類中最大的RegionDif值記為Rmax；另一類中的RegionDif值較大，代表該類中RegionDif值對應的兩個分割區域光譜差異大，不應該合并，將該類中最小的RegionDif記為Rmin，則閾值B為：

若RegionDifij＜B，則代表相鄰的兩個區域的光譜值接近，可合并為一個區域。

3.3 變化圖斑提取

根據光譜差異合并后，影像會被劃分為少量幾個區域，每個區域內像素值接近，區域與區域之間像素值差異較大。若區域的像素平均值較大，則該區域為變化區域；反之，則為非變化區域。令區域i所有像素所有波段的均值為：

可以找到合適的閾值C，使得若meani>C，則將區域i判定為變化區域。但是不同的差異影像會有不同的最佳閾值，所以此處仍使用K均值聚類算法來尋找自適應的閾值C。

計算合并后剩下的各區域的mean，對得到的所有集合作K均值聚類，其中K=2。將mean的集合分為兩類：一類中的mean值較小，代表該類mean所對應的區域為非變化區域，將該類中最大的mean記為Mmax；另一類中的mean值較大，代表該類mean所對應的區域為變化區域，將該類中最小的mean記為Mmin，則閾值C為：

若meani>C，則區域i為變化區域；反之，則為非變化區域。

4 實驗與結果分析

4.1 實驗數據

本文實驗采用的遙感影像數據為2012年﹑2013年杜爾伯特蒙古族自治縣地區的SPOT衛星影像，包含紅﹑綠﹑藍3個波段，分辨率為2.5 m。兩期影像的地物類型豐富，主要由水體﹑建筑物﹑道路﹑植被等多類地物構成，且影像的分辨率較高。對實驗數據進行MAD變換，得到MAD結果。目視判別MAD結果，選擇真實變化區域的外圍矩形作為ROI，并利用本文算法獲取最后變化圖斑。

本文選擇了土地利用遙感監測中常見的3種典型變化地物類型，處理結果如圖1～3所示。圖1a為人工勾勒的變化圖斑，作為真實值對算法進行精度評定；圖1b為ROI上生成的差異影像；圖1c為ROI對應的MAD結果；圖1d為利用本文算法得到的結果。對MAD算法和本文算法進行精度評定，結果如表1所示。

4.2 實驗分析

圖1c中MAD檢測的道路不完整，上部分出現斷裂，同時存在離散的噪聲點；圖2c中MAD檢測建筑中存在空洞，同時由于配準誤差，上部分未變化的曲線道路被檢測為變化區域；圖3c中MAD檢測結果更加離散﹑破碎，其中細長部分的建筑斷開導致區域不連通，同時存在空洞，且左上部分的建筑屋頂由于輻射差異導致部分像素被MAD判定為變化區域。綜上可知，MAD可有效發現變化，但檢測的結果離散﹑破碎﹑不連通，形狀不夠完整，同時存在偽變化區域，不能直接用于圖斑提取。

圖1～3的差異影像充分反映了前后時相影像的變化信息，影像中較亮的部分代表變化區域。觀察MAD檢測結果和差異影像可以發現：①由配準誤差造成的偽變化區域在差異影像上的像素值仍較亮，但分布離散，面積較小；由輻射差異造成的偽變化區域像素值雖比背景區域的像素值高，但比真實變化區域的像素值要低。②MAD確定的變化區域內存在空洞部分，在差異影像上偏暗。

圖1 植被-道路類型變化圖斑提取結果

圖2 植被-建筑類型變化圖斑提取結果

圖3 植被-推填土類型變化圖斑提取結果

表1 變化檢測精度評定

從圖1～3的分割結果可以看出，本文算法可剔除MAD上由于輻射差異﹑配準誤差等因素造成的誤差。因為在區域合并時，首先將小面積區域和周圍區域合并，而由于配準誤差形成的偽變化區域分布離散，面積小，所以會被周圍較暗的背景區域合并。輻射差異造成的偽變化區域像素值低于真實變化區域，在利用區域間差異性進行區域合并時，會在K均值聚類中被判定為背景。從分割結果中還可以看出，變化區域內部的空洞被填充，是一個完整的變化圖斑，破碎﹑不連通的部分被合并為一個整體。這是由于在第一次區域合并時，空洞對應的較暗區域與周圍較亮區域合并，整體呈現變化特征，在最后的K均值聚類中被判定為變化區域。

從表1可以看出本文算法漏提與錯提區域明顯小于MAD算法，分割結果與人工提取的結果相似，3種檢測結果的Kappa系數為分別從0.770﹑0.810﹑0.729提高到了0.916﹑0.894﹑0.934。

5 結 語

常用的變化檢測方法檢測得到的變化區域離散﹑破碎﹑不連通，且存在大量偽變化區域，不能直接用于圖斑提取。本文設計的影像分割與合并算法有效地分割了差異影像，得到了完整的變化圖斑，解決了MAD結果離散﹑破碎﹑不連通的問題，在錯檢率和漏檢率方面優于MAD算法。本文算法避免了人工勾勒變化圖斑邊界，可有效提高變化圖斑提取效率，減少人工操作難度。下一步需要研究ROI自動化提取辦法，實現全自動化變化圖斑提取。

[1] 王素敏, 翟輝琴. 遙感技術在我國土地利用/覆蓋變化中的應用[J].地理空間信息,2004,2(2):31-32

[2] 吳芳,劉榮,田維春,等.遙感變化檢測技術及其應用綜述[J].地理空間信息,2007,5(4):57-60

[3] Tewkesbury A P, Comber A J, Tate N J, et al. A Critical Synthesis of Remotely Sensed Optical Image Change Detection Techniques[J].Remote Sensing of Environment,2015(160):1-14

[4] Volpi M, Tuia D, Bovolo F, et al. Supervised Change Detection in VHR Images Using Contextual Information and Support Vector Machines[J].International Journal of Applied Earth Observation and Geoinformation,2013,20(2):77-85

[5] Listner C, Niemeyer I. Object-based Change Detection[J]. Photogrammetrie Fernerkundung Geoinformation,2011(4):233-245

[6] CHEN X, CHEN J, SHI Y，et al. An Automated Approach for Updating Land Cover Maps Based on Integrated Change Detection and Classification Methods[J].ISPRS Journal of Photogrammetry and Remote Sensing,2012,71(326):86-95

[7] 佃袁勇,方圣輝,姚崇懷.一種面向地理對象的遙感影像變化檢測方法[J].武漢大學學報(信息科學版),2014,39(8): 906-912

[8] Bruzzone L, Prieto D F. Automatic Analysis of the Difference Image for Unsupervised Change Detection[J].IEEE Transactions on Geoscience and Remote Sensing,2000,38(3):1 171-1 182

[9] Nielsen A A. The Regularized Iteratively Reweighted MAD Method for Change Detection in Multi- and Hyperspectral Data[J]. IEEE Transactions on Image Processing a Publication of the IEEE Signal Processing Society,2007,16(2):463-478

[10] GONG P. Change Detection Using Principal Component Analysis and Fuzzy Set Theory[J].Canadian Journal of Remote Sensing,1993,19(1):22-29

[11] 陳波,張友靜,陳亮.標記分水嶺算法及區域合并的遙感圖像分割[J].國土資源遙感,2007,19(2):35-38

P237

B

1672-4623（2017）09-0093-04

10.3969/j.issn.1672-4623.2017.09.028

2016-06-28。

項目來源：國土資源部公益性行業科研專項資助項目（201511009-01）。

夏旺，碩士研究生，主要研究方向為攝影測量與遙感。