基于融合差異圖的變化檢測(cè)方法及其在洪災(zāi)中的應(yīng)用

黃平平 段盈宏 譚維賢 徐 偉

(內(nèi)蒙古工業(yè)大學(xué)信息工程學(xué)院 呼和浩特 010051)

(內(nèi)蒙古自治區(qū)雷達(dá)技術(shù)與應(yīng)用重點(diǎn)實(shí)驗(yàn)室 呼和浩特 010051)

1 引言

由于惡劣的氣候和天氣影響自然災(zāi)害頻發(fā)，其中，洪災(zāi)嚴(yán)重破壞了生態(tài)環(huán)境，損害了私人和公共財(cái)產(chǎn)，對(duì)人類社會(huì)的生存發(fā)展有重大影響[1]。因此，對(duì)洪害前后進(jìn)行變化檢測(cè)，確定洪災(zāi)的淹沒(méi)范圍，對(duì)災(zāi)情的評(píng)估和災(zāi)后的管理規(guī)劃都有一定的指導(dǎo)作用。合成孔徑雷達(dá)(Synthetic Aperture Radar,SAR)可以全天時(shí)、全天候地獲取遙感數(shù)據(jù)，不受光照和天氣的影響，彌補(bǔ)了光學(xué)和紅外遙感的不足。隨著遙感技術(shù)的發(fā)展，不同時(shí)相、不同頻帶以及不同極化方式的多種SAR圖像資源為SAR圖像變化檢測(cè)提供了數(shù)據(jù)基礎(chǔ)，并廣泛應(yīng)用在災(zāi)區(qū)定位和災(zāi)害評(píng)估等方面[2]。

SAR圖像的變化檢測(cè)方法一般可分為兩種[3]：監(jiān)督和非監(jiān)督。監(jiān)督類變化檢測(cè)方法需要大量的先驗(yàn)知識(shí)，而實(shí)際中往往缺乏真實(shí)的參考信息，適用范圍較窄。非監(jiān)督類變化檢測(cè)方法不需要先驗(yàn)信息，可以直接對(duì)獲得的兩時(shí)相圖像進(jìn)行變化檢測(cè)，在實(shí)際應(yīng)用中更為廣泛，且非監(jiān)督類變化檢測(cè)技術(shù)多在差異圖生成和變化信息提取兩方面進(jìn)行創(chuàng)新。傳統(tǒng)方法多使用不同的變化檢測(cè)算子構(gòu)建差異圖，如經(jīng)典的差值法和比值法[4,5]。相對(duì)熵，又被稱為KL散度(Kullback-Leibler divergence)或信息散度，是兩個(gè)概率分布間差異的非對(duì)稱性度量，以此構(gòu)造的差異圖可衡量?jī)煞跋耖g的相似程度。考慮到散斑噪聲影響以及檢測(cè)結(jié)果的精度，Gong等人[6]利用鄰域像元灰度值和中心像元灰度值的相關(guān)性，提出一種基于鄰域比值算子的差異圖生成方法。單一類型的差異圖像有自身的局限性，如均值比差異圖對(duì)目標(biāo)區(qū)域的群分性較好，但對(duì)紋理信息感知較弱，且存有大量噪聲，對(duì)數(shù)比差異圖在一定程度上削弱了背景部分的斑點(diǎn)噪聲，但殘留的噪聲仍對(duì)后續(xù)處理有很大的影響[7–9]?；诖耍恍┭芯咳藛T選擇使用圖像融合的方法，利用不同差異圖像的互補(bǔ)信息以及不同傳感器的并行信息，在城市化和河口區(qū)域建筑的變化檢測(cè)中展現(xiàn)了較好的結(jié)果[10–12]。在不同雷達(dá)圖像數(shù)據(jù)中，通過(guò)融合多時(shí)相雷達(dá)數(shù)據(jù)以及相干和非相干性分析，對(duì)氣旋、滑坡等災(zāi)害區(qū)域的評(píng)估也有一定的效果[13–15]。小波變換相對(duì)于輪廓波等多尺度變換而言，具有更小的計(jì)算復(fù)雜度，同時(shí)，小波變換可以將圖像的時(shí)空信息轉(zhuǎn)換為頻域信息，進(jìn)而可以方便地提取圖像的詳細(xì)信息，在圖像融合中有廣泛的應(yīng)用[16–18]。在變化信息提取方面，Sumaiya等人[19]使用加權(quán)閾值分割對(duì)數(shù)均值比差異圖，獲得變化檢測(cè)結(jié)果，計(jì)算過(guò)程快速簡(jiǎn)便。Zhang等人[20]利用廣義高斯分布對(duì)差異圖像的每一類進(jìn)行建模，采用圖割算法提取空間先驗(yàn)信息并通過(guò)模糊C均值算法很好地初始化模型的參數(shù)，在最大后驗(yàn)的貝葉斯準(zhǔn)則下提取變化信息，其精度與算法速度都有一定的提升。Li等人[21]提出了一種稱為SDPCANet的監(jiān)督學(xué)習(xí)網(wǎng)絡(luò)模型，該模型將顯著性檢測(cè)與主成分分析網(wǎng)絡(luò)相結(jié)合，并在多時(shí)相SAR圖像中驗(yàn)證了方法的有效性。Li等人[22]提出了一種基于卷積神經(jīng)網(wǎng)絡(luò)(Convolutional Neural Networks,CNN)的SAR圖像變化檢測(cè)新方法，該方法通過(guò)CNN直接從原始的兩個(gè)SAR圖像中生成分類結(jié)果，無(wú)需進(jìn)行任何預(yù)處理操作，消除了生成差異圖像的過(guò)程，從而減少了差異圖對(duì)源圖像的影響，并且在異質(zhì)圖像的變化檢測(cè)中也有較好的效果。

洪災(zāi)發(fā)生時(shí)往往伴隨著多云多雨，而光學(xué)傳感器無(wú)法穿透云層，無(wú)法獲取有效的地面信息，因此，SAR圖像的變化檢測(cè)更有效地應(yīng)用在洪災(zāi)區(qū)域中。由于洪災(zāi)發(fā)生時(shí)降雨影響了地物散射特性[23]，在整幅圖像中許多沒(méi)有受洪災(zāi)影響的區(qū)域也呈現(xiàn)為偽暗區(qū)。在傳統(tǒng)的差異圖中，偽暗區(qū)與背景信息差異較大，而與變化信息更為接近，一般的變化檢測(cè)算法更容易將偽暗區(qū)作為變化信息處理，這也導(dǎo)致了檢測(cè)結(jié)果中較高的虛警率和較低的檢測(cè)精度[24]。由于偽暗區(qū)的存在，對(duì)差異圖的直接聚類易導(dǎo)致較多的錯(cuò)誤信息，故可利用皮爾遜相關(guān)系數(shù)對(duì)模糊局部信息C均值(Fuzzy Local Information C-Means,FLICM)聚類方法的初始聚類結(jié)果進(jìn)行二次分類，減小對(duì)檢測(cè)結(jié)果的影響。一般而言，洪水淹沒(méi)區(qū)域的相關(guān)系數(shù)更低，而標(biāo)準(zhǔn)未變化區(qū)域的相關(guān)系數(shù)更高，所以通過(guò)比較待定區(qū)域的相關(guān)系數(shù)，待定區(qū)域的類別，更容易劃分為變化區(qū)域和未變化區(qū)域中兩者相關(guān)系數(shù)均值更小的類別[25–28]。

本文結(jié)構(gòu)安排如下：本文第2部分提出了一種利用融合差異圖的SAR圖像變化檢測(cè)方法，然后在此差異圖下得到基于相關(guān)系數(shù)FLICM算法二次分類后的檢測(cè)結(jié)果，最后經(jīng)迭代條件模型和馬爾科夫隨機(jī)場(chǎng)(Iterative Conditional Model-Markov Random Field,ICM-MRF)算法優(yōu)化得到變化檢測(cè)最終結(jié)果；第3部分首先利用Bern地區(qū)的ERS-2實(shí)測(cè)SAR數(shù)據(jù)以及加拿大Ottawa地區(qū)的Radarsat實(shí)測(cè)SAR數(shù)據(jù)進(jìn)行結(jié)果分析并驗(yàn)證本文方法的有效性，之后利用該方法對(duì)中國(guó)鄱陽(yáng)湖區(qū)域2020年6月和7月洪災(zāi)前后的Sentinel-1-A實(shí)測(cè)SAR數(shù)據(jù)進(jìn)行分析，進(jìn)而對(duì)受災(zāi)面積和其后的災(zāi)害趨勢(shì)進(jìn)行評(píng)估；第4部分對(duì)全文做出總結(jié)。

2 SAR圖像變化檢測(cè)算法

本文算法流程主要由融合差異圖構(gòu)造和ICMMRF分割[29,30]兩部分組成，算法流程如圖1所示。

融合差異圖構(gòu)造：首先得到兩幅SAR圖像的均值比差異圖和利用鄰域異質(zhì)性信息的相對(duì)熵差異圖，然后經(jīng)小波變換，通過(guò)加權(quán)平均的方式將兩幅圖像的低頻系數(shù)重構(gòu)，與各自的高頻系數(shù)，經(jīng)小波逆變換得到兩幅融合低頻系數(shù)后的差異圖，通過(guò)局部能量法融合為最終的差異圖。

ICM-MRF分割：在ICM-MRF算法中，使用FLICM算法將變化檢測(cè)結(jié)果分為3類，分別為變化區(qū)域、未變化區(qū)域和待定區(qū)域，對(duì)待定區(qū)域通過(guò)洪災(zāi)前后圖像的皮爾遜相關(guān)系數(shù)(后簡(jiǎn)稱為相關(guān)系數(shù))進(jìn)一步分為變化區(qū)域和未變化區(qū)域，并以此作為初始分割，得到最終的變化檢測(cè)結(jié)果。

2.1 融合差異圖構(gòu)造

假設(shè)X1,X2為已經(jīng)過(guò)配準(zhǔn)預(yù)處理的兩時(shí)像SAR洪災(zāi)前后變化區(qū)域圖像，則離散圖像數(shù)據(jù)X2對(duì)X1的相對(duì)熵可表示為

其中，N表示像元總個(gè)數(shù)，xn表示第n個(gè)像元的值，Ys表示X2對(duì)X1的相對(duì)熵。

圖1 基于融合差異圖的變化檢測(cè)算法流程圖Fig.1 Flow chart of change detection algorithm in this paper

在相對(duì)熵的計(jì)算中引入像元空間鄰域的異質(zhì)性信息，突出圖像邊緣，增強(qiáng)連續(xù)性，同時(shí)抑制背景信息。鄰域異質(zhì)性信息可表示為

?表示圖像像元間的異質(zhì)度，其定義為像元鄰域的方差δ(x)和其鄰域的均值μ(x)的比值。則引入鄰域異質(zhì)性信息后的相對(duì)熵計(jì)算公式可定義為

其中，Yd為X2對(duì)X1的鄰域異質(zhì)性相對(duì)熵，?n1,?n2分別為圖像X1,X2同一位置像元對(duì)應(yīng)的異質(zhì)性系數(shù)歸一化后的結(jié)果，μ1,μ2為其分別對(duì)應(yīng)的去除中心像元點(diǎn)的鄰域均值

式(4)為相對(duì)熵計(jì)算的對(duì)稱性處理，Yk為本文改進(jìn)的相對(duì)熵差異圖的最終結(jié)果。使用鄰域異質(zhì)信息來(lái)調(diào)節(jié)圖像相對(duì)熵的計(jì)算，熵值越小表示兩幅圖像同一位置對(duì)應(yīng)像元的相似程度越高，反之則相似程度越低。

引入鄰域異質(zhì)性信息的相對(duì)熵計(jì)算方法，相對(duì)于傳統(tǒng)的單一像元值的熵值計(jì)算方法，有較好的背景抑制和邊緣連續(xù)性效果，但目標(biāo)區(qū)域的完整性和平滑性一般，故將目標(biāo)區(qū)域完整性和平滑性較好的均值比算子差異圖與本文引入鄰域異質(zhì)性信息的差異圖，通過(guò)小波變換進(jìn)行融合。均值比算子可定義為

其中，μ1(n)和μ2(n)分別表示在圖像X1和X2第n個(gè)像元其鄰域(一般取3×3大小)內(nèi)的平均值，Ym表示均值比差異圖。

本文圖像通過(guò)小波變換進(jìn)行融合，小波變換一般先將圖像分解為更低分辨率水平上的低頻輪廓信息，以及在水平垂直和對(duì)角線方向的高頻細(xì)節(jié)信息，再經(jīng)過(guò)小波系數(shù)的替換、選擇權(quán)值或者疊加運(yùn)算進(jìn)行融合。本文選擇的融合規(guī)則如下：(1)取兩圖像小波分解后的加權(quán)低頻系數(shù)均值作為新的低頻系數(shù)；(2)將新的低頻系數(shù)與各自的高頻系數(shù)經(jīng)小波逆變換得到融合低頻系數(shù)后的差異圖；(3)將兩幅新的差異圖通過(guò)局部能量法得到目標(biāo)差異圖。經(jīng)過(guò)本文規(guī)則得到的差異圖中，具有均值比差異圖中較好的區(qū)域完整性和平滑性，也具有改進(jìn)相對(duì)熵差異圖中背景抑制性較強(qiáng)、邊緣連續(xù)性較好的優(yōu)點(diǎn)。

2.2 ICM-MRF分割

FLICM是一種無(wú)監(jiān)督的模糊聚類算法，通過(guò)圖像的統(tǒng)計(jì)特性可以較好地完成分類，進(jìn)而提取圖像的變化信息。使用FLICM算法將融合差異圖分為兩類的目標(biāo)函數(shù)為

若以μn1和μn2分別表示為第n個(gè)像元對(duì)于兩類的隸屬度，通過(guò)比較兩類的隸屬度大小，則可以確定初始分割結(jié)果。

通過(guò)FLICM算法將差異圖分為3類，分別代表變化區(qū)域、待定區(qū)域以及未變化區(qū)域。計(jì)算洪災(zāi)前后兩SAR圖像的相關(guān)系數(shù)，如式(10)所示

式中，Ic(i,j)表示為(i,j)像元點(diǎn)對(duì)應(yīng)的相關(guān)系數(shù)值；X1(i,j)和X2(i,j)均表示圖像的(i,j)像元點(diǎn)對(duì)應(yīng)的像元值；表示X(i,j)像元點(diǎn)周圍鄰域像元值的均值。計(jì)算各類區(qū)域相關(guān)系數(shù)均值，比較待定區(qū)域相關(guān)系數(shù)與變化區(qū)域和未變化區(qū)域相關(guān)系數(shù)的均值，進(jìn)行二次分類。

融合差異圖通過(guò)以上算法二次分類后，其結(jié)果可作為ICM-MRF算法的初始分割。假設(shè)融合差異圖像為Y={yη=yij}，初始分割圖像X={xη=xij}為定義在Y上的馬爾科夫隨機(jī)場(chǎng)，是對(duì)圖像所做的標(biāo)記；xη是點(diǎn)η的類別標(biāo)記，xη ∈I(I=1,2,···,k),k為類別數(shù)，本文中k=2。設(shè)點(diǎn)η的鄰域值記為sη,s={sη,η ∈Y}是η的鄰域系統(tǒng)。求解最優(yōu)圖像標(biāo)記問(wèn)題的關(guān)鍵在于求解先驗(yàn)概率P(X)和概率密度函數(shù)P(Y |X)?；谪惾~斯決策理論以及馬爾科夫隨機(jī)場(chǎng)和吉布斯分布的等價(jià)性，選取2階鄰域系統(tǒng)和合適的模型，用β表示馬爾科夫參數(shù)，則局部條件先驗(yàn)概率可表示為

Vh為使用Potts數(shù)學(xué)模型的基團(tuán)勢(shì)能函數(shù)，H為所有基團(tuán)的集合。由于要得到后驗(yàn)概率P(X|Y)最大，即等效于求解U(X|Y)最小，采用迭代條件模型(Iterative Conditional Model,ICM)求解，它的核心思想是用逐像元點(diǎn)最大化條件概率，迭代地實(shí)現(xiàn)每個(gè)像元點(diǎn)的更新，如式(15)所示

在馬爾科夫隨機(jī)場(chǎng)(Markov Random Field,MRF)算法中，作為初始分割的標(biāo)簽圖對(duì)MRF的分割結(jié)果有較大的影響，因此將利用圖像相關(guān)性對(duì)FLICM算法二次分類后的結(jié)果作為MRF算法的初始標(biāo)簽圖，以期得到更好的初始分割結(jié)果，則基于融合差異圖的ICM-MRF算法如表1所示。

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

本文選擇的驗(yàn)證數(shù)據(jù)為瑞士Bern地區(qū)在1999年4月和5月的ERS-2數(shù)據(jù)，兩時(shí)相遙感圖像為圖2(a)、圖2(b)，大小為301×301，變化檢測(cè)參考圖為圖2(c)；加拿大Ottawa地區(qū)在1997年5月和8月的Radarsat遙感數(shù)據(jù)，兩時(shí)相遙感圖像為圖3(a)、圖3(b)，大小為290×350。兩數(shù)據(jù)集產(chǎn)生圖像變化的主要原因均為洪水淹沒(méi)陸地所致[2]。

表1 基于融合差異圖的ICM-MRF圖像分割計(jì)算過(guò)程Tab.1 ICM-MRF image segmentation calculation process based on fusion difference map

使用上述數(shù)據(jù)進(jìn)行仿真試驗(yàn)，通過(guò)對(duì)Bern地區(qū)數(shù)據(jù)以及Ottawa地區(qū)數(shù)據(jù)均值比、相對(duì)熵方法得到的差異圖，以及本文融合差異圖進(jìn)行比較，得到了上述所有方法在兩數(shù)據(jù)集下不同差異圖的期望最大化(Expectation-Maximization,EM)算法、K-means算法、FLICM算法以及本文ICM-MRF算法的檢測(cè)結(jié)果，并分別在兩數(shù)據(jù)集下進(jìn)行主觀和客觀的評(píng)價(jià)。其后通過(guò)本文方法對(duì)2020年中國(guó)鄱陽(yáng)湖地區(qū)數(shù)據(jù)進(jìn)行變化檢測(cè)，進(jìn)而估計(jì)由于鄱陽(yáng)湖湖水覆蓋而受災(zāi)的區(qū)域面積和變化趨勢(shì)。

本文采用虛警率(False Positive,FP)、總錯(cuò)誤率(Overall Error,OE)、準(zhǔn)確率(Probability Correct Classification,PCC)和Kappa系數(shù)作為量化評(píng)估的指標(biāo)。其中，PCC表示總體的分類正確性，Kappa系數(shù)綜合考慮了像元點(diǎn)的檢測(cè)狀況，是衡量檢測(cè)精度的關(guān)鍵指標(biāo)，其值越大表示檢測(cè)精度越高。

3.1 差異圖2分類和3分類影響分析

Bern地區(qū)數(shù)據(jù)集的不同算法差異圖如圖4(a)—圖4(c)所示，并以A表示均值比差異圖，B表示相對(duì)熵差異圖，C表示本文差異圖，在此3種差異圖下，以EM算法、K-means算法和FLICM算法為例，分析了2分類和3分類后基于相關(guān)系數(shù)二次分類對(duì)結(jié)果的影響，如圖5、圖6和圖7所示。對(duì)比各算法的差異圖可得，均值比算法在背景中存在大量噪聲，相對(duì)而言，相對(duì)熵法雖大量抑制了噪聲，但模糊了目標(biāo)區(qū)域的結(jié)構(gòu)。本文構(gòu)造的差異圖綜合了相對(duì)熵和均值比差異圖的優(yōu)勢(shì)，邊緣連續(xù)性和區(qū)域完整性較相對(duì)較好，與背景的對(duì)比度較兩差異圖更為明顯，且噪聲點(diǎn)較均值比差異圖更少。從以下變化檢測(cè)結(jié)果圖中不難看出，當(dāng)僅將差異圖劃分為2類時(shí)，噪聲點(diǎn)和一些區(qū)域也將直接定義為變化區(qū)域或者為變化區(qū)域，直接影響了結(jié)果的區(qū)域劃分和檢測(cè)精度。而將差異圖先劃分為變化區(qū)域，未變化區(qū)域以待定區(qū)域3類區(qū)域，再將其根據(jù)不同劃分區(qū)域的相關(guān)系數(shù)二次分類后，將有效地減小噪聲的影響，并加強(qiáng)了區(qū)域的完整性，進(jìn)而提升檢測(cè)結(jié)果的精度。相對(duì)于均值比和異質(zhì)相對(duì)熵差異圖，在不同算法下，本文差異圖也顯示了更好的分類效果。

圖2 Bern地區(qū)數(shù)據(jù)Fig.2 Bern region data

圖4 Bern數(shù)據(jù)差異圖Fig.4 Bern data difference map

圖5 基于均值比差異圖2分類和3分類后二次分類結(jié)果Fig.5 Based on the difference map in mean ratio 2 and 3 classification results after secondary classification

圖6 基于相對(duì)熵差異圖2分類和3分類后二次分類結(jié)果Fig.6 Based on the difference map in relative entropy 2 and 3 classification results after secondary classification

3.2 Bern和Ottawa地區(qū)實(shí)驗(yàn)結(jié)果與分析

Bern地區(qū)數(shù)據(jù)集對(duì)應(yīng)的每種差異圖的不同算法變化檢測(cè)結(jié)果的相關(guān)指標(biāo)分析如表2所示。對(duì)于不同差異圖Bern數(shù)據(jù)集的本文算法變化檢測(cè)結(jié)果如圖8(a)—圖8(c)所示。檢測(cè)結(jié)果中白色點(diǎn)為變化區(qū)域的像元點(diǎn)。從變化檢測(cè)結(jié)果圖中也可看出，在本文ICM-MRF算法下，融合差異圖的檢測(cè)結(jié)果顯示了更好的變化區(qū)域邊界，與實(shí)際的變化情況更為接近。表2中也顯示了針對(duì)不同差異圖下同一種算法的變化檢測(cè)結(jié)果，本文差異圖相對(duì)于其他差異圖而言，具較低的虛警率與總錯(cuò)誤率，更高的正確率和Kappa系數(shù)。針對(duì)本文差異圖而言，由于EM算法檢測(cè)結(jié)果中存在較多噪聲，故Kappa系數(shù)較低，但仍可看清變化區(qū)域的完整結(jié)構(gòu)，并且本文ICM-MRF算法具有更高的Kappa系數(shù)和較低的錯(cuò)誤率，相對(duì)于檢測(cè)效果較好的FLICM算法，Kappa系數(shù)提高至0.837。

圖7 基于融合差異圖2分類和3分類后二次分類結(jié)果Fig.7 Based on the fusion difference and the results of secondary classification after classification 2 and 3

表2 Bern數(shù)據(jù)集不同差異圖不同算法變化檢測(cè)指標(biāo)分析Tab.2 Analysis of change detection indicators of different algorithms in Bern dataset

圖8 Bern數(shù)據(jù)本文ICM-MRF算法變化檢測(cè)結(jié)果Fig.8 Change detection results of this paper algorithm in Bern data

圖9 Ottawa數(shù)據(jù)差異圖及其變化檢測(cè)結(jié)果Fig.9 Ottawa data difference map change detection result

Ottawa地區(qū)數(shù)據(jù)集的不同算法差異圖如圖9(a)—圖9(c)所示，并以A表示均值比差異圖，B表示相對(duì)熵差異圖，C表示本文差異圖，其對(duì)應(yīng)的每種差異圖的不同算法變化檢測(cè)結(jié)果的相關(guān)指標(biāo)分析如表3所示。對(duì)于不同差異圖Ottawa數(shù)據(jù)集的本文ICM-MRF算法變化檢測(cè)結(jié)果如圖9(d)—圖9(f)所示，EM算法變化檢測(cè)結(jié)果如圖9(g)—圖9(i)所示，K-means算法變化檢測(cè)結(jié)果為圖9(j)—圖9(l)所示，F(xiàn)LICM算法變化檢測(cè)結(jié)果為圖9(m)—圖9(o)所示。檢測(cè)結(jié)果中白色點(diǎn)為變化區(qū)域的像元點(diǎn)。通過(guò)對(duì)比各算法的差異圖不難看出，均值比算法、相對(duì)熵算法與Bern區(qū)域結(jié)果分析一致，且本文的差異圖構(gòu)造方法在Ottawa數(shù)據(jù)集內(nèi)也體現(xiàn)出邊緣連續(xù)性，區(qū)域完整性以及背景對(duì)比度均相對(duì)較好的效果。不同差異圖的變化檢測(cè)結(jié)果中，構(gòu)造的融合差異圖在本文ICM-MRF算法下變化檢測(cè)結(jié)果相對(duì)于均值比差異圖的檢測(cè)結(jié)果具有更好的細(xì)節(jié)表現(xiàn)，而相對(duì)于相對(duì)熵差異圖的檢測(cè)結(jié)果具有更好的區(qū)域完整性，且Kappa系數(shù)也有一定提升。表3中也顯示了與Bern數(shù)據(jù)集基本相似的指標(biāo)分析結(jié)果，由此可見本文算法在兩數(shù)據(jù)集均具有較好的適用性。整體而言，本文算法可以提高變化檢測(cè)結(jié)果的精度，減小錯(cuò)誤率。

3.3 鄱陽(yáng)湖區(qū)域洪災(zāi)前后變化檢測(cè)結(jié)果分析

由于連續(xù)多日暴雨的影響，中國(guó)鄱陽(yáng)湖區(qū)域發(fā)生洪災(zāi)，淹沒(méi)了周邊大片陸地區(qū)域，人工無(wú)法有效地實(shí)時(shí)實(shí)地測(cè)量受災(zāi)區(qū)域而制訂善后處理計(jì)劃，利用SAR衛(wèi)星遙感圖像，可以快速準(zhǔn)確地提取受災(zāi)區(qū)域，評(píng)估災(zāi)區(qū)的空間分布和動(dòng)態(tài)變化。本文選取中國(guó)鄱陽(yáng)湖地區(qū)2020年6月和7月Sentinel-1-A數(shù)據(jù)進(jìn)行實(shí)驗(yàn)來(lái)評(píng)估此次洪災(zāi)。由于完整的鄱陽(yáng)湖區(qū)域分布在Sentinel-1-A數(shù)據(jù)的兩景圖像內(nèi)，故本文使用對(duì)兩景數(shù)據(jù)分別檢測(cè)的方法進(jìn)行實(shí)驗(yàn)。具體實(shí)驗(yàn)數(shù)據(jù)為2020年6月26日數(shù)據(jù)，拼合后完整鄱陽(yáng)湖區(qū)域圖像如圖10(a)所示，7月8日數(shù)據(jù)如圖10(b)所示，7月20日數(shù)據(jù)如圖10(c)所示，3日包含鄱陽(yáng)湖整體區(qū)域的Sentinel-1-A數(shù)據(jù)。圖10(a)、圖10(b)和圖10(c)大小均為999×651，像素大小為243×243 m2，產(chǎn)生圖像變化的主要原因是鄱陽(yáng)湖水淹沒(méi)附近區(qū)域陸地所致。

通過(guò)EM算法、FLICM算法、K-means算法和本文ICM-MRF算法，在融合差異圖下對(duì)鄱陽(yáng)湖Sentinel-1-A衛(wèi)星遙感圖像進(jìn)行變化檢測(cè)，EM算法變化檢測(cè)結(jié)果如圖11(a)、圖11(b)所示，F(xiàn)LICM算法變化檢測(cè)結(jié)果如圖11(c)、圖11(d)所示，K-means算法變化檢測(cè)結(jié)果如圖11(e)、圖11(f)所示，本文ICM-MRF算法變化檢測(cè)結(jié)果如圖11(g)、圖11(h)所示。通過(guò)對(duì)比6月26日與7月8日的變化檢測(cè)結(jié)果，明顯看到在本文差異圖下兩種算法均有較好的檢測(cè)結(jié)果，得到的變化區(qū)域也較為相似，但不同的是，在圖11(c)、圖11(d)的FLICM算法中，黃色方框區(qū)域?yàn)樘摼瘏^(qū)域，此區(qū)域在真實(shí)雷達(dá)圖中并未被洪水淹沒(méi)，圖11(a)、圖11(b)EM算法和圖11(e)、圖11(f)K-means算法中，黃色方框區(qū)域?yàn)槁┚瘏^(qū)域，此區(qū)域?yàn)楸缓樗蜎](méi)的區(qū)域，而本文算法均未將兩區(qū)域錯(cuò)誤劃分，整體來(lái)看，本文算法下，變化區(qū)域結(jié)構(gòu)更為完整，檢測(cè)結(jié)果更為準(zhǔn)確。

基于本文算法，圖11(g)—圖11(j)變化檢測(cè)結(jié)果為3幅鄱陽(yáng)湖洪災(zāi)前后變化檢測(cè)結(jié)果。其中，圖11(g)為6月26日與7月8日的變化檢測(cè)結(jié)果，圖11(h)為其對(duì)應(yīng)在實(shí)際中的變化區(qū)域。圖11(i)為7月8日與7月20日的變化檢測(cè)結(jié)果，圖11(j)為其對(duì)應(yīng)在實(shí)際中的變化區(qū)域。白色、紅色和藍(lán)色像元點(diǎn)代表變化區(qū)域。由此可見，本文算法較完整且清晰地檢測(cè)出鄱陽(yáng)湖區(qū)域的洪災(zāi)前后變化范圍。選取的包含完整鄱陽(yáng)湖區(qū)域的SAR圖像實(shí)際面積約為38402.46 km2，由于受暴雨影響，鄱陽(yáng)湖水位上升，7月8日—7月20日，鄱陽(yáng)湖西南部區(qū)域受災(zāi)嚴(yán)重，湖水淹沒(méi)了湖周邊大片陸地，包括附近水域淹沒(méi)的水中小島和陸地等。由表4鄱陽(yáng)湖區(qū)域變化檢測(cè)結(jié)果量化分析可知，此期間變化像素?cái)?shù)為12589，占總像素點(diǎn)的1.94%，總體受災(zāi)面積大約為 743.37 km2。7月8日—7月20日，在表4中仍可看到5475個(gè)變化像素點(diǎn)，占總像素點(diǎn)的0.84%，這表明總體受災(zāi)面積仍在提升，但變化像素點(diǎn)減少了大約56.5%，洪水淹沒(méi)區(qū)域較之前明顯減少，期間總受災(zāi)面積大約為323.29 km2，洪災(zāi)擴(kuò)散勢(shì)頭逐漸減小。值得注意的是，圖11(j)中藍(lán)色區(qū)域雖為變化區(qū)域，但不是洪水淹沒(méi)區(qū)域，而是由水體轉(zhuǎn)為陸地的區(qū)域，這也說(shuō)明部分區(qū)域的洪水被引流或者逐漸退去，總體形趨于穩(wěn)定。

表3 Ottawa數(shù)據(jù)集不同差異圖不同算法變化檢測(cè)指標(biāo)分析Tab.3 Analysis of change detection indicators of different algorithms in Ottawa dataset

圖10 鄱陽(yáng)湖地區(qū)數(shù)據(jù)Fig.10 Poyang Lake area data

圖11 鄱陽(yáng)湖區(qū)域EM,FLICM,K-means和本文算法洪災(zāi)前后變化檢測(cè)結(jié)果Fig.11 Change detection results of EM,FLICM,K-means and this paper algorithm before and after the flood disaster in Poyang Lake area

表4 鄱陽(yáng)湖區(qū)域變化檢測(cè)結(jié)果量化分析Tab.4 Quantitative analysis of the detection results of Poyang Lake area change

4 結(jié)論

針對(duì)洪災(zāi)區(qū)域的變化檢測(cè)，本文提出一種基于融合差異圖的變化檢測(cè)方法，利用鄰域異質(zhì)性信息保持相對(duì)熵差異圖中變化區(qū)域的連續(xù)性，抑制背景信息，并通過(guò)小波變換，構(gòu)造與均值比差異圖融合的最終差異圖。通過(guò)FLICM聚類算法與圖像的相關(guān)系數(shù)結(jié)合而得到更接近于洪災(zāi)變化區(qū)域的二次分類結(jié)果，然后將其經(jīng)過(guò)ICM-MRF算法優(yōu)化得到最終的變化檢測(cè)結(jié)果。在Bern地區(qū)ERS-2遙感數(shù)據(jù)以及Ottawa地區(qū)Radarsat遙感數(shù)據(jù)實(shí)驗(yàn)結(jié)果均表明，本文方法對(duì)洪災(zāi)區(qū)域的變化檢測(cè)結(jié)果有較好的檢測(cè)精度，變化區(qū)域結(jié)構(gòu)更為完整準(zhǔn)確。通過(guò)本文算法得到2020年6,7月份中國(guó)鄱陽(yáng)湖地區(qū)的3景Sentinel-1-A遙感圖像變化檢測(cè)結(jié)果，并對(duì)其受災(zāi)區(qū)域面積以及變化趨勢(shì)進(jìn)行估計(jì)，具有一定的實(shí)用價(jià)值。今后將繼續(xù)研究洪災(zāi)區(qū)域產(chǎn)生虛警、漏警的原因以及在變化檢測(cè)中自動(dòng)區(qū)別洪水淹沒(méi)區(qū)域和水體轉(zhuǎn)為陸地區(qū)域的方法。