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

夏 旺，趙 展

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

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

夏 旺1，趙 展1

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

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

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

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

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

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

1 基于MAD的差異影像獲取

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

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

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

2 差異影像的分水嶺分割

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

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

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

1）將Seg所有像素設(shè)為0，標(biāo)記值記為label，初始label記為1；

2）順序遍歷圖像，當(dāng)前遍歷像素為p，若p=0，則標(biāo)記為label；

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

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

5）當(dāng)堆棧為空時(shí)，label加1，返回步驟2）。

3 分割區(qū)域合并

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

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

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

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

3.1 根據(jù)面積合并

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

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

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

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

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

3.3 變化圖斑提取

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

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

計(jì)算合并后剩下的各區(qū)域的mean，對(duì)得到的所有集合作K均值聚類，其中K=2。將mean的集合分為兩類：一類中的mean值較小，代表該類mean所對(duì)應(yīng)的區(qū)域?yàn)榉亲兓瘏^(qū)域，將該類中最大的mean記為Mmax；另一類中的mean值較大，代表該類mean所對(duì)應(yīng)的區(qū)域?yàn)樽兓瘏^(qū)域，將該類中最小的mean記為Mmin，則閾值C為：

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

4 實(shí)驗(yàn)與結(jié)果分析

4.1 實(shí)驗(yàn)數(shù)據(jù)

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

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

4.2 實(shí)驗(yàn)分析

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

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

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

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

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

表1 變化檢測(cè)精度評(píng)定

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

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

5 結(jié) 語

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

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

[2] 吳芳,劉榮,田維春,等.遙感變化檢測(cè)技術(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] 佃袁勇,方圣輝,姚崇懷.一種面向地理對(duì)象的遙感影像變化檢測(cè)方法[J].武漢大學(xué)學(xué)報(bào)(信息科學(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] 陳波,張友靜,陳亮.標(biāo)記分水嶺算法及區(qū)域合并的遙感圖像分割[J].國(guó)土資源遙感,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。

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

夏旺，碩士研究生，主要研究方向?yàn)閿z影測(cè)量與遙感。