夏 旺，趙 展

（1.武漢大學(xué) 測繪學(xué)院，湖北 武漢430079）

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

夏 旺1，趙 展1

（1.武漢大學(xué) 測繪學(xué)院，湖北 武漢430079）

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

變化圖斑；分水嶺分割；區(qū)域合并；變化檢測

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

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

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

1 基于MAD的差異影像獲取

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

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

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

2 差異影像的分水嶺分割

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

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

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

1）將Seg所有像素設(shè)為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 分割區(qū)域合并

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

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

設(shè)區(qū)域i與區(qū)域j的光譜差異RegionDifij為：

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

3.1 根據(jù)面積合并

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

3.2 根據(jù)光譜差異合并

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

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

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

3.3 變化圖斑提取

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

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

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

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

4 實驗與結(jié)果分析

4.1 實驗數(shù)據(jù)

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

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

4.2 實驗分析

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

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

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

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

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

表1 變化檢測精度評定

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

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

5 結(jié) 語

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

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

[2] 吳芳,劉榮,田維春,等.遙感變化檢測技術(shù)及其應(yīng)用綜述[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].武漢大學(xué)學(xué)報(信息科學(xué)版),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] 陳波,張友靜,陳亮.標記分水嶺算法及區(qū)域合并的遙感圖像分割[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。

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

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