基于變權(quán)馬爾可夫隨機(jī)場的遙感圖像變化檢測

高 鵬 翔

(中國工程物理研究院計算機(jī)應(yīng)用研究所 四川 綿陽 621900)

0 引 言

隨著人類科技水平的迅猛發(fā)展，人類在地球上的活動范圍愈發(fā)廣泛。因此及時、精確地開展對地觀測來實(shí)時監(jiān)測地表的覆蓋變化情況，對人類的發(fā)展具有重要意義[1]。當(dāng)代航空、航天和傳感器等技術(shù)的進(jìn)步，促使了遙感圖像獲取手段的日益成熟[2]。當(dāng)下采用遙感數(shù)據(jù)進(jìn)行變化檢測已是監(jiān)測地表覆蓋變化最精確、高效的技術(shù)手段，亦是遙感研究的熱點(diǎn)。遙感圖像變化檢測是指通過分析一對或者一系列在相同地區(qū)不同時間點(diǎn)獲取的遙感圖像，來確定地表是否發(fā)生改變的技術(shù)[3]。

變化檢測已在城市規(guī)劃與建設(shè)農(nóng)業(yè)、牧業(yè)、林業(yè)、工業(yè)、軍事等領(lǐng)域內(nèi)得到了普遍地運(yùn)用，如土地利用/土地覆蓋情況監(jiān)測、地理數(shù)據(jù)庫更新、災(zāi)害評估、森林植被覆蓋監(jiān)測、土地沙漠化監(jiān)測、水系變化監(jiān)測、城市建設(shè)規(guī)劃、軍事目標(biāo)監(jiān)測等[4-8]。

馬爾可夫隨機(jī)場模型(MRF)由于其完備的數(shù)學(xué)理論，以及能夠充分利用圖像中的空間上下文信息[9]，在遙感圖像變化檢測上取得了可喜的成就。Bruzzone等[10]首先使用MRF進(jìn)行遙感圖像變化檢測，對差值圖像使用EM算法，分別估計不變類和變化類概率密度，之后引入MRF，通過條件迭代式(ICM)算法[11]求解。Benedek等[12]采用條件混合馬爾可夫模型來進(jìn)行變化檢測，這是一種多層MRF模型，結(jié)合了混合馬爾可夫模型和信號的條件獨(dú)立隨機(jī)場。Chen等[13]提出了一種上下文敏感的弱監(jiān)督變化檢測技術(shù)，該方法通過MRF模型分析概率高斯過程分類器的后驗(yàn)概率實(shí)現(xiàn)。Chen等[14]針對傳統(tǒng)的MRF方法在不連續(xù)區(qū)域(如邊界、山脊、山谷等)經(jīng)常得到錯誤的變化檢測結(jié)果，提出了一種無監(jiān)督變化檢測算法。Gong等[9]提出了一種基于模糊聚類和改進(jìn)MRF能量函數(shù)的SAR變化檢測算法，該方法在計算能量函數(shù)時同時考慮了像素點(diǎn)的隸屬度和同種類別的鄰域像素數(shù)目改變。Subudhi等[15]提出了一種基于模糊聚類和MRF的變化檢測算法，使用模糊吉布斯馬爾可夫隨機(jī)場(GMRF)來對差值圖像的空間灰度特性建模，采用基于可變鄰域搜索的全局收斂準(zhǔn)則來迭代估計模糊GMRF模型的參數(shù)，避免了傳統(tǒng)算法容易得到局部最優(yōu)值問題。范奎奎等[16]將雙樹復(fù)小波變換和MRF相結(jié)合，解決了變化檢測中多尺度圖像降噪時丟失高頻信息與單像素孤立問題。

然而，傳統(tǒng)的基于馬爾可夫隨機(jī)場的變化檢測算法，對于中心像素和鄰域像素使用固定的權(quán)重，會造成空間鄰域信息的過度使用[14-15]，使得傳統(tǒng)的標(biāo)記場不能精確地確定鄰域像素間的空間關(guān)系，導(dǎo)致變化檢測結(jié)果過度平滑，影響變化檢測的精度。

針對上述問題，本文提出一種基于變權(quán)馬爾可夫隨機(jī)場的遙感圖像變化檢測算法。首先通過兩時圖像生成差值圖像；隨后使用模糊C均值聚類算法對差值圖像上的像素點(diǎn)進(jìn)行分類；接著使用一個改進(jìn)的MRF模型來對差值圖像建模；最后使用ICM算法求解能量函數(shù)最小化得到變化檢測結(jié)果。

1 傳統(tǒng)的基于MRF的變化檢測

定義大小為M×N的兩時差值圖像為XD，其滿足XD={x1,x2，…，xn}，n=M×N。變化類別與未變化類別為ω={ωch,ωun}，Y為類別標(biāo)記場，對于二值變化檢測問題，將圖像上所有像素點(diǎn)分為變化與未變化兩個類別。根據(jù)MRF理論有：

Y=arg max{P(Y)P(XD|Y)}

(1)

式中：P(Y)為類別標(biāo)記的先驗(yàn)概率；P(XD|Y)是數(shù)據(jù)集中像素點(diǎn)的條件概率密度。

根據(jù)Hammersley-Clifford定理，可知吉布斯隨機(jī)場(GRF)與馬爾可夫隨機(jī)場(MRF)等價，則式(1)值最大時等同于MRF能量函數(shù)最小。對于每個像素點(diǎn)xi，其能量函數(shù)UMRF(xi)定義為兩部分組成：

UMRF(xi)=Uspectral(xi)+Uspatial(xi)

(2)

式中：Uspectral(xi)是光譜能量函數(shù)，描述了像素點(diǎn)xi的光譜特性，是該像素點(diǎn)的灰度統(tǒng)計描述；Uspatial(xi)是空間能量函數(shù)，描述了像素點(diǎn)xi與其所在鄰域內(nèi)其余像素點(diǎn)之間的類別依賴。

通常假定兩時差值圖像的灰度統(tǒng)計特性滿足高斯分布，于是采用高斯分布來描述光譜能量函數(shù)Uspectral(xi)：

(3)

對于空間能量函數(shù)Uspatial(xi)，其根據(jù)中心像素點(diǎn)xi的標(biāo)記與xi的鄰域像素的標(biāo)記情況來確定。且Uspatial(xi)滿足：

(4)

(5)

式中：β為懲罰系數(shù)，通常β>0，其值由用戶設(shè)定，通過取不同的β值可以控制鄰域像素對中心像素的影響，調(diào)整空間上下文信息對變化檢測過程的影響；Ni是像素點(diǎn)xi的鄰域，通常采用二階鄰域，即對于每一個中心像素點(diǎn)周圍有8個像素點(diǎn)；l(xi)和l(xj)是像素點(diǎn)類別標(biāo)記，采用Potts模型來定義中心像素點(diǎn)與其鄰域像素點(diǎn)類別標(biāo)記關(guān)系I(l(xi),l(xj))；對于二值指示函數(shù)I(l(xi),l(xj))，當(dāng)中心像素點(diǎn)標(biāo)記與其鄰域像素的標(biāo)記相同時其值取為-1，否則取為0。

變化檢測可視為典型的組合優(yōu)化問題，其結(jié)果可通過迭代搜索算法所獲取的全局或局部最優(yōu)解來得到。使用ICM算法便可以獲得式(2)的最優(yōu)值。

2 基于變權(quán)MRF的變化檢測

2.1 變權(quán)懲罰系數(shù)

在傳統(tǒng)的基于MRF模型的變化檢測算法中，Potts模型使用一個固定的懲罰系數(shù)β來計算空間能量函數(shù)。對于差值圖像中所有的像素點(diǎn)來說，β的取值是相同的，這沒有考慮差值圖像中像素灰度值的分布特性，易造成差值圖像中灰度值極大或者極小的像素點(diǎn)對于圖像鄰域信息的使用過度，從而影響變化檢測精度。

對于差值圖像首先使用FCM算法，將圖像上的像素點(diǎn)初步分為兩個類別，這會得到兩個聚類中心。基于此結(jié)果，可以得到兩個閾值H1、H2，因此可以將差值圖像分為三部分。其中H1左側(cè)為未變化類別，H1與H2之間為不確定類別，H2右側(cè)為變化類別。H1與H2的計算方式如下：

(6)

式中：C1為未變化類別聚類中心；C2為變化類別聚類中心；Mmid為模糊C均值聚類后具有相同變化與未變化的隸屬度的像素點(diǎn)所具有的灰度值，即該像素點(diǎn)屬于變化類別與未變化類別的可能性相同，兩種隸屬度值均為0.5；α1和α2是兩個常數(shù)，可以用來調(diào)節(jié)三個部分的面積。該劃分如圖1所示。

圖1 變化、不確定、未變化類別劃分示意圖

在此基礎(chǔ)上，根據(jù)差值圖像中像素點(diǎn)的灰度值所在的區(qū)間，采用三個不同的表達(dá)式來計算所在區(qū)間的懲罰系數(shù)。于是修正式(4)中的固定懲罰系數(shù)β。其計算公式為：

(7)

式中：βm(X(i,j))是像素點(diǎn)位置為(i,j)的懲罰系數(shù)；Xmin和Xmax分別是差值圖像中灰度最小和最大的值。

當(dāng)圖像中灰度值小于H1時，隨著灰度值從H1減少到Xmin時，懲罰系數(shù)線性減小，這減小了具有相對較小的灰度值的像素對于空間鄰域信息的影響，這是因?yàn)檫@些像素點(diǎn)有極大的概率保持不變。同理，當(dāng)灰度值大于H2時，隨著灰度值由H2增加到Xmax，懲罰系數(shù)線性減小，這減小了具有相對較大的灰度值的像素對于空間鄰域信息的影響，這是因?yàn)檫@些像素點(diǎn)有極大的概率保持變化。對于灰度值處于H1與H2中間時，其屬于不確定的類別，對于這類像素點(diǎn)采用一個固定的懲罰系數(shù)。以此重新定義了使用的空間信息的權(quán)重。

2.2 變權(quán)指示函數(shù)

在傳統(tǒng)的基于MRF模型的變化檢測算法中假設(shè)某個像素點(diǎn)常常與其所在鄰域中的像素點(diǎn)具有相同的變化類別標(biāo)記。然而由于外部條件所限，獲得的遙感數(shù)據(jù)中包含大量混合像素。此外由于圖像數(shù)據(jù)和算法本身的固有不確定性，差異圖像的初始標(biāo)記類別并不是確定的，這尤其體現(xiàn)在圖像細(xì)節(jié)豐富的區(qū)域中。所以純粹地將指示函數(shù)的值取為-1和0是不合適的，該操作會顯得太絕對和不精確。本文對Potts模型進(jìn)行修正，來重新定義鄰域像素間的空間關(guān)系。

Foody[17]提出了一種基于熵的方法來度量鄰域空間的像素標(biāo)記的不確定性。采用該方法在模糊C均值聚類后，對于任意像素點(diǎn)可以得到其屬于變化與未變化兩個類別的隸屬度。故可得像素點(diǎn)xj的熵，定義如下：

Ej=-uj,ch×log2uj,ch-uj,un×log2uj,un

(8)

式中：uj,ch和uj,un分別是該像素點(diǎn)屬于變化與未變化類別的隸屬度值。于是可以將式(5)修正為：

(9)

(10)

這樣在計算空間能量函數(shù)時，可以充分考慮到鄰域像素的類別標(biāo)記的不確定性，使變化檢測結(jié)果的過度平滑問題得到了很好的控制。故式(4)可以重新改寫為：

(11)

2.3 算法描述

該變化檢測算法的主要步驟概括如下：

1) 對配準(zhǔn)、校正后的兩時圖像做差獲得差值圖像。

3) 使用變權(quán)MRF模型來修正空間能量函數(shù)Uspatial(xi)的計算，使用改進(jìn)的變權(quán)懲罰系數(shù)和變權(quán)指示函數(shù)替換原始計算方式，重新定義空間上下文信息。

4) 使用ICM算法求解能量函數(shù)最小值，得到最終的變化檢測結(jié)果。

該算法的示意圖如圖2所示。

圖2 基于變權(quán)MRF變化檢測算法示意圖

3 實(shí) 驗(yàn)

3.1 實(shí)驗(yàn)數(shù)據(jù)及環(huán)境

實(shí)驗(yàn)使用的數(shù)據(jù)集為AICDDataset數(shù)據(jù)集和Szada數(shù)據(jù)集。AICDDataset為光學(xué)遙感模擬數(shù)據(jù)集，該數(shù)據(jù)集包含1 000組大小為800×600的兩時圖像及對應(yīng)的變化檢測參考圖像，空間分辨率為0.5 m×0.5 m。Szada為真實(shí)遙感影像數(shù)據(jù)集，包含7組大小為952×640的兩時圖像以及對應(yīng)的變化檢測參考圖像，空間分辨率為1.5 m×1.5 m。

分別進(jìn)行如下實(shí)驗(yàn)：(1) 從2個數(shù)據(jù)集中隨機(jī)選取2組數(shù)據(jù)，分別進(jìn)行2組模擬實(shí)驗(yàn)、2組真實(shí)實(shí)驗(yàn)進(jìn)行算法效果展示，所用的數(shù)據(jù)如圖3所示。(2) 從AICDDataset數(shù)據(jù)集中隨機(jī)抽取100組數(shù)據(jù)，比較本文算法與對比實(shí)驗(yàn)算法的效果。變化參考圖中白色區(qū)域代表兩時圖像發(fā)生變化的位置，黑色區(qū)域代表未發(fā)生變化的位置。

圖3 部分實(shí)驗(yàn)數(shù)據(jù)

實(shí)驗(yàn)所使用的數(shù)據(jù)均為經(jīng)過精確配準(zhǔn)、幾何校正、輻射校正的數(shù)據(jù)，消除了由于外部不良因素導(dǎo)致的偽變化，通過MATLAB R2015b編制程序求解。

3.2 實(shí)驗(yàn)結(jié)果及分析

將本文算法與文獻(xiàn)[10]和文獻(xiàn)[15]的變化檢測算法進(jìn)行對比，以上4組實(shí)驗(yàn)結(jié)果分別如圖4-圖7所示。本文采用Kappa系數(shù)來評價變化檢測的結(jié)果，其取值范圍在[-1,1]之間，Kappa系數(shù)的大小反映了變化檢測結(jié)果圖與變化參考圖之間一致性的程度，其值越接近于1，表示變化檢測結(jié)果圖與變化參考圖的匹配程度越高。

圖4 實(shí)驗(yàn)(1)結(jié)果圖

圖5 實(shí)驗(yàn)(2)結(jié)果圖

圖6 實(shí)驗(yàn)(3)結(jié)果圖

圖7 實(shí)驗(yàn)(4)結(jié)果圖

上述實(shí)驗(yàn)的Kappa系數(shù)統(tǒng)計結(jié)果如表1所示。

表1 Kappa系數(shù)統(tǒng)計表

從AICDDataset數(shù)據(jù)集中隨機(jī)抽取100組數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，分別統(tǒng)計本文算法與對比實(shí)驗(yàn)算法的Kappa系數(shù)，結(jié)果如圖8所示。

圖8 Kappa系數(shù)對比圖

上述兩次實(shí)驗(yàn)結(jié)果表明，在絕大部分情況下，本文提出的變化檢測算法與變化參考圖具有最高的相似度，取得了比對比算法更好的變化檢測結(jié)果。本文算法實(shí)現(xiàn)了對變化信息的有效提取，其含有較少的漏檢和誤檢像素，對噪聲等偽變化信息實(shí)現(xiàn)了良好的甄別，噪聲點(diǎn)更少，實(shí)驗(yàn)結(jié)果也更平滑。

4 結(jié) 語

針對傳統(tǒng)基于MRF模型變化檢測算法在能量函數(shù)的計算中使用固定的懲罰系數(shù)和固定指示函數(shù)的問題，提出了一種基于變權(quán)MRF模型的無監(jiān)督的遙感變化檢測算法。該算法通過使用變權(quán)懲罰系數(shù)和基于熵的局部不確定鄰域標(biāo)記的變權(quán)指示函數(shù)，重新定義了空間信息權(quán)重和MRF模型中的鄰域像素的空間關(guān)系，進(jìn)而實(shí)現(xiàn)對空間能量函數(shù)的改造，優(yōu)化了傳統(tǒng)MRF模型。實(shí)驗(yàn)結(jié)果表明本文算法提供了比傳統(tǒng)MRF模型更精確的變化檢測，避免了變化檢測的過度平滑問題。因此本文算法是一種可行、有效的變化檢測算法。