高光譜協(xié)同稀疏與非局部低秩張量變化檢測(cè)

詹天明，宋 博，孫 樂(lè)，萬(wàn)鳴華，楊國(guó)為

1.南京審計(jì)大學(xué) 信息工程學(xué)院，南京211815

2.南京審計(jì)大學(xué) 審計(jì)信息工程與技術(shù)協(xié)同創(chuàng)新中心，南京211815

3.南京信息工程大學(xué) 計(jì)算機(jī)學(xué)院，南京210094

變化檢測(cè)是地球科學(xué)以及遙感領(lǐng)域研究最多的問(wèn)題之一。它的任務(wù)是利用同一場(chǎng)景的兩個(gè)或多個(gè)時(shí)間點(diǎn)的圖像，以檢測(cè)不同時(shí)間點(diǎn)內(nèi)發(fā)生的變化。其應(yīng)用范圍包括自然災(zāi)害監(jiān)測(cè)、城市和森林的長(zhǎng)期跟蹤、農(nóng)作物的生長(zhǎng)狀態(tài)觀測(cè)、土地資源開(kāi)發(fā)利用審計(jì)等。變化檢測(cè)從任務(wù)角度主要分為：判斷是否發(fā)生變化、變化位置定位、變化發(fā)生的種類(lèi)等。而變化檢測(cè)過(guò)程主要包括預(yù)處理、變化檢測(cè)算法和效果評(píng)估這三個(gè)步驟。其中，變化檢測(cè)算法是其核心步驟。目前，常被用于變化檢測(cè)的主要有合成孔徑雷達(dá)圖像、全色圖像、多光譜圖像以及高光譜圖像。具有高光譜分辨率的高光譜圖像可以提供比其他遙感圖像更豐富的光譜信息，因此高光譜圖像常被用于識(shí)別細(xì)微變化，反映出不同物體的詳細(xì)成分等問(wèn)題。近年來(lái)，基于高光譜圖像的變化檢測(cè)被很多學(xué)者廣泛研究，但其較低的空間分辨率、光譜維冗余信息以及較少的數(shù)據(jù)集都給高光譜變化檢測(cè)方法的發(fā)展帶來(lái)了較大影響。

近年來(lái)，基于二分類(lèi)的變化檢測(cè)方法較為常見(jiàn)，它的目的是將多時(shí)相遙感影像劃分為變化類(lèi)和非變化類(lèi)。Malila提出了變化向量分析（change vector analysis，CVA）方法，該方法計(jì)算兩個(gè)時(shí)相的光譜變化向量，再根據(jù)變化幅度和變化方向判定每個(gè)像素點(diǎn)是否發(fā)生變化。在此基礎(chǔ)上，Baisantry 等人提出了基于增強(qiáng)主成分分析（principal component analysis，PCA）和反三角函數(shù)閾值的變化矢量分析方法（PCA-CVA）。該方法核心思想是用兩時(shí)相圖像的主成分差分和反三角函數(shù)合理選擇閾值，進(jìn)而判定當(dāng)前像素點(diǎn)是否發(fā)生變化。Nielsen 等人利用典型性相關(guān)分析的降維屬性提出多元變化檢測(cè)（multivariate alteration detection，MAD）方法，該方法通過(guò)求解多元變化值以判斷各個(gè)像素是否在時(shí)間維上發(fā)生變化。此后，Nielsen 在多元變化檢測(cè)方法基礎(chǔ)上提出了迭代加權(quán)多元變化檢測(cè)方法（iteratively reweighted MAD，IR-MAD），該方法在迭代步驟中對(duì)每個(gè)像素賦予的權(quán)重進(jìn)行更新，迭代結(jié)束后以權(quán)重作為判斷像素是否變化的標(biāo)準(zhǔn)。

基于監(jiān)督的遙感影像的變化檢測(cè)方法也被廣泛研究。其中一類(lèi)方法的主要思想是通過(guò)監(jiān)督學(xué)習(xí)對(duì)各時(shí)相遙感圖像進(jìn)行分類(lèi)，然后根據(jù)每個(gè)像素點(diǎn)在不同時(shí)間點(diǎn)上的類(lèi)別信息進(jìn)行變化檢測(cè)。Bovolo 等人提出的變化檢測(cè)框架中先對(duì)兩時(shí)相的遙感影像分別進(jìn)行逐像素分類(lèi)，如果一個(gè)像素點(diǎn)在兩時(shí)相的影像中分類(lèi)不同，就可以認(rèn)為這個(gè)像素點(diǎn)產(chǎn)生了變化。另一類(lèi)方法是對(duì)差異圖像進(jìn)行監(jiān)督分類(lèi)，將差異圖像分為變化和非變化兩類(lèi)，從而完成變化檢測(cè)任務(wù)。Nemmour 等人提出了多支持向量機(jī)（multiple support vector machine，Multiple SVM）的變化檢測(cè)方法，該方法先對(duì)變化圖像進(jìn)行二分類(lèi)，再用模糊積分和吸引子動(dòng)力學(xué)指標(biāo)對(duì)分類(lèi)結(jié)果進(jìn)行優(yōu)化。隨著深度學(xué)習(xí)方法的發(fā)展，基于深度學(xué)習(xí)的高光譜變化檢測(cè)方法也被廣泛研究。Wang 等人提出端到端2D 卷積神經(jīng)網(wǎng)絡(luò)框架的高光譜變化檢測(cè)方法。該方法利用兩時(shí)相的光譜信息以及解混結(jié)果為每個(gè)像素點(diǎn)構(gòu)建一個(gè)混合矩陣，再用卷積神經(jīng)網(wǎng)絡(luò)進(jìn)行分類(lèi)。

雖然遙感變化檢測(cè)方法被廣泛提出，但受限于高光譜圖像的低空間分辨率和高光譜分辨率，直接應(yīng)用于高光譜變化檢測(cè)時(shí)其精度都不是太高。而現(xiàn)有基于深度學(xué)習(xí)方法的高光譜變化檢測(cè)方法需要大量的訓(xùn)練樣本才能保證較高的檢測(cè)精度。此外，這些方法需要將三維高光譜圖像轉(zhuǎn)換為二維矩陣，一定程度上破壞了高光譜圖像空間關(guān)聯(lián)性。為解決該問(wèn)題，將三維高光譜圖像表示為一個(gè)張量的形式進(jìn)行處理，在高光譜圖像降維、解混、分類(lèi)和融合中引起了廣泛關(guān)注。張量模型利用多線(xiàn)性代數(shù)的強(qiáng)大框架，可以聯(lián)合利用光譜和空間信息，有助于提升后續(xù)高光譜圖像處理的性能。

目前，張量建模在高光譜圖像融合以及分類(lèi)領(lǐng)域已得到較為廣泛的研究，但在高光譜變化檢測(cè)問(wèn)題中利用張量進(jìn)行建模的方法較少。為充分利用張量模型中空譜聯(lián)合特征以提高高光譜變化檢測(cè)精度，本文提出一種基于協(xié)同稀疏與非局部低秩張量的高光譜變化檢測(cè)方法。該方法首先求得兩相位的高光譜差分圖像，再根據(jù)差分圖像中圖像塊非局部特點(diǎn)，提取不同的非局部相似張量簇。然后基于協(xié)同稀疏正則化和低秩正則化建立變化檢測(cè)模型，利用極小化方法求解得到表示系數(shù)。最后根據(jù)表示系數(shù)求出張量在變化和非變化兩類(lèi)中的投影殘差，由投影殘差判斷該張量是否發(fā)生變化。該方法主要有如下兩個(gè)優(yōu)點(diǎn)：（1）將非局部圖像塊組合在一起形成一系列低階張量，可以深挖高光譜空譜聯(lián)合信息，同時(shí)降低維度；（2）協(xié)同稀疏表示分類(lèi)可以在較少的訓(xùn)練樣本下保證一定的分類(lèi)精度，且在高光譜圖像稀疏分解方面比傳統(tǒng)的低階正則化分類(lèi)器具有更好的性能，從而克服高光譜變化圖像數(shù)據(jù)集不足的問(wèn)題。在仿真和真實(shí)數(shù)據(jù)集中進(jìn)行對(duì)比實(shí)驗(yàn)，結(jié)果表明本文方法可以取得較好的高光譜變化檢測(cè)精度。

1 高光譜圖像的張量表示

設(shè)χ ∈R為-階張量，各個(gè)I表示模式，張量中的每個(gè)元素表示為x。階張量的F 范數(shù)和范數(shù)定義為：

階張量沿模展開(kāi)表示為∈R，將還原成張量則可表示為χ=fold()，其中fold是張量沿模展開(kāi)的逆操作。設(shè)∈R表示大小為J×I的矩陣，和χ 的階積運(yùn)算定義為：

設(shè)J∈R表示高光譜圖像，∈R是字典，S ∈R是系數(shù)張量，N∈R是噪聲。它們之間的關(guān)系可以表示為：

其中，、表示高光譜圖像的空間維的高和寬，表示高光譜圖像的波段數(shù)，表示字典中的原子個(gè)數(shù)。

2 本文方法

本文方法首先將兩時(shí)相圖像J和J做差分取絕對(duì)值得到J=|J-J|∈R，然后對(duì)差異絕對(duì)值圖建立字典={,}∈R，其中、分別表示變化和不變區(qū)域的子字典。建立基于協(xié)同稀疏與非局部低秩張量正則化模型，求出J 在字典中的表示系數(shù)，最后根據(jù)其在兩個(gè)子空間的投影殘差判斷每個(gè)像素是否發(fā)生變化?？傮w流程圖如圖1 所示。

圖1 本文方法總體流程圖Fig.1 General flow chart of proposed method

2.1 協(xié)同稀疏模型

由式（4）高光譜圖像在字典下的建模方式，兩時(shí)相圖像差異絕對(duì)值圖像可建模為：

其中，S ∈R是未知的表示系數(shù)的張量形式，N∈R是殘差。根據(jù)表示學(xué)習(xí)理論，通過(guò)求解如下優(yōu)化問(wèn)題即可求得表示系數(shù)S ：

其中，(S)是S 的先驗(yàn)函數(shù)。在表示學(xué)習(xí)理論中，將系數(shù)的張量形式乘以字典，必然導(dǎo)致系數(shù)的張量形式在空間維上僅包含少量非零值，但在最近鄰理論中，每個(gè)像素點(diǎn)與字典中的少數(shù)原子具有較強(qiáng)的相關(guān)性，因此每個(gè)像素點(diǎn)在字典下的表示系數(shù)應(yīng)該具有稀疏性。為此，基于上述先驗(yàn)，建立如下協(xié)同稀疏模型：

根據(jù)文獻(xiàn)[15]，協(xié)同稀疏正則項(xiàng)||S||定義如下：

2.2 協(xié)同稀疏與非局部低秩張量分解模型

將上述局部張量塊應(yīng)用于整個(gè)圖像域中，表示系數(shù)張量形式S 的核范數(shù)可表示為：

將式（9）與式（6）相結(jié)合，得到如下的協(xié)同稀疏與非局部低秩張量模型：

其中，是非負(fù)常數(shù)，第一項(xiàng)是求得J 在字典中的一個(gè)投影系數(shù)，而第二項(xiàng)則體現(xiàn)了該系數(shù)的協(xié)同稀疏先驗(yàn)，第三項(xiàng)S 的核范數(shù)即為圖像的非局部低秩約束，體現(xiàn)了該系數(shù)的空間非局部相關(guān)性先驗(yàn)。對(duì)該模型極小化即可得到J在字典中最優(yōu)表示系數(shù)。

2.3 優(yōu)化算法

模型（10）無(wú)法直接求解，為此采用交替方向乘子法將目標(biāo)函數(shù)轉(zhuǎn)化為如下多個(gè)子問(wèn)題分別求解：

其中，U、U、是引入變量。式（11）的拉格朗日方程可表示為：

其中，是拉格朗日懲罰因子，V、V、V、V是拉格朗日乘子，是懲罰函數(shù)，當(dāng)≥0 時(shí)，()=0，當(dāng)＜0 時(shí)，()=+∞。

由交替方向乘子法，上述變量可通過(guò)交替迭代方式求解。求解變量S 的目標(biāo)函數(shù)為：

其迭代解為：

其中，是單位矩陣，是的轉(zhuǎn)置。

求解變量U、U、U、U的目標(biāo)函數(shù)為：

其中，函數(shù)表達(dá)式為：

函數(shù)表達(dá)式為：

(,)=diag(max{()-,0})(·)表示奇異值分解。

拉格朗日乘子V、V、V、V更新步驟為：

本文方法的總體算法步驟如下：

基于協(xié)同稀疏與非局部低秩張量的高光譜變化檢測(cè)算法

輸入：第一時(shí)相高光譜圖像J，第二時(shí)相高光譜圖像J。

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

3.1 數(shù)據(jù)集

實(shí)驗(yàn)采用了3 組數(shù)據(jù)，第一組是Farmland 數(shù)據(jù)集，第二組和第三組是UASFC（Urban area in San Francisco City）數(shù)據(jù)集。

Farmland數(shù)據(jù)集是由EO-1（Earth observing-1）Hyperion設(shè)備采集的2006年5月3 日 和2007 年4 月23 日江蘇省鹽城市一個(gè)農(nóng)田的高光譜圖像。其光譜范圍為0.4～2.5 μm，光譜分辨率10 nm 左右，空間分辨率大約為30 m。該圖像共包含155 個(gè)波段，空間大小為450×140 像素。兩相位偽彩色圖及真實(shí)變化圖如圖2 所示。

圖2 Farmland 數(shù)據(jù)集Fig.2 Farmland dataset

UASFC 數(shù)據(jù)集是由無(wú)人機(jī)機(jī)載傳感器采集的2009年9月18日和2015年5月11日舊金山市市區(qū)的高光譜圖像。拍攝頻率為1.26 GHz，光譜空間分辨率為1.66 m，方位分辨率為1 m，入射角為25°～65°。共兩組數(shù)據(jù)：第一組數(shù)據(jù)（UASFC1）空間大小為200×200 像素，含138 波段；第二組數(shù)據(jù)空間（UASFC2）大小為100×100 像素，含138 波段。兩相位偽彩色圖及真實(shí)變化圖如圖3 和圖4 所示。

圖3 Urban area in San Francisco City 第一組數(shù)據(jù)集Fig.3 The first dataset of Urban area in San Francisco City

圖4 Urban area in San Francisco City 第二組數(shù)據(jù)集Fig.4 The second dataset of Urban area in San Francisco City

3.2 模型評(píng)價(jià)指標(biāo)

本文使用總體準(zhǔn)確度（overall accuracy，OA）和Kappa 系數(shù)兩個(gè)指標(biāo)評(píng)估模型的變化檢測(cè)效果。為了計(jì)算OA 和Kappa 值，需要統(tǒng)計(jì)模型對(duì)每個(gè)像素的預(yù)測(cè)情況：模型正確分類(lèi)了變化的像素?cái)?shù)量記為真陽(yáng)性（true positive，TP）；模型正確分類(lèi)了未變化的樣本數(shù)量記為真陰性（true negative，TN）；模型錯(cuò)誤地將未變化的像素分類(lèi)為變化像素的數(shù)量記為假陽(yáng)性（false positive，F(xiàn)P）；模型錯(cuò)誤地將變化的像素分類(lèi)為未變化像素的數(shù)量記為假陰性（false negative，F(xiàn)N）。記為像素總數(shù)，則有：

用于衡量模型正確的分類(lèi)像素比例，計(jì)算公式為：

Kappa 系數(shù)衡量模型預(yù)測(cè)結(jié)果和地面真值圖的一致性，計(jì)算公式為：

其中，p的計(jì)算公式為：

3.3 參數(shù)設(shè)置及分析

本文的迭代終止條件為=10，若迭代次數(shù)超過(guò)30 次尚未達(dá)到該條件則終止迭代。本文模型中的兩個(gè)參數(shù)和分別控制聯(lián)合稀疏約束項(xiàng)和非局部低秩約束項(xiàng)，其值較低時(shí)約束較小，空譜聯(lián)合特征不明顯，檢測(cè)結(jié)果中噪聲過(guò)多會(huì)導(dǎo)致精度下降；而當(dāng)其值較高時(shí)約束較大，檢測(cè)結(jié)果過(guò)多地受到空譜聯(lián)合約束的影響，數(shù)據(jù)本身對(duì)變化檢測(cè)的決策能力下降，變化檢測(cè)精度也不高。為此，根據(jù)各項(xiàng)的數(shù)值量級(jí)關(guān)系取∈[0.000 5,0.003 0],∈[0.001 0,0.005 0] 進(jìn)行實(shí)驗(yàn)。圖5 顯示了不同的和值下，本文方法的檢測(cè)精度。從圖中結(jié)果可以看出，本文方法在=0.001 0 和=0.002 5 時(shí)取得最優(yōu)的檢測(cè)精度，且在∈[0.001 0,0.002 5]，∈[0.001 0,0.005 0]的范圍內(nèi)都保持一個(gè)較高的檢測(cè)精度。

圖5 選取不同參數(shù)時(shí)三組數(shù)據(jù)集的平均OA 值Fig.5 Average OA values of three datasets with different parameters

對(duì)于監(jiān)督學(xué)習(xí)的方法，訓(xùn)練樣本個(gè)數(shù)對(duì)變化檢測(cè)精度的影響也是比較敏感的。為此，本文在每個(gè)數(shù)據(jù)集中的訓(xùn)練樣本總數(shù)取1 000 到8 000 進(jìn)行實(shí)驗(yàn)，其中正負(fù)樣本個(gè)數(shù)各占50%。圖6 顯示了不同訓(xùn)練樣本個(gè)數(shù)下本文方法的變化檢測(cè)精度。從圖中可以看出，本文方法在Farmland 數(shù)據(jù)集中選取2 000 個(gè)訓(xùn)練樣本即可得到較高的檢測(cè)精度，在Urban area in San Francisco City 兩組數(shù)據(jù)集中選取4 000 個(gè)訓(xùn)練樣本得到較高的檢測(cè)精度。

圖6 三組數(shù)據(jù)集取不同訓(xùn)練樣本數(shù)量時(shí)的OA 值Fig.6 OA values of three datasets with different number of training samples

本文方法需要選擇圖像塊構(gòu)造非局部張量，為選擇最優(yōu)尺度的圖像塊，分別選取3×3，5×5 和7×7 三種尺度在Farmland、UASFC1 和UASFC2 數(shù)據(jù)集上進(jìn)行變化檢測(cè)。對(duì)應(yīng)的OA 值如表1 所示。從表中可以看出，圖像塊尺度為3×3 大小時(shí)，本文方法在3 個(gè)數(shù)據(jù)集上都能取得較好的OA 值。這是因?yàn)樵诟吖庾V變化檢測(cè)問(wèn)題中，變化區(qū)域包含有很多較小且不規(guī)則的圖像塊。若尺度較大，會(huì)對(duì)這些圖像塊進(jìn)行平滑，從而影響變化檢測(cè)的精度。在本文實(shí)驗(yàn)中，選取3×3 大小的圖像塊構(gòu)造非局部張量進(jìn)行變化檢測(cè)。

表1 不同尺度圖像塊在三組數(shù)據(jù)集上的OA 值Table 1 OA values of different scales image blocks on three datasets

3.4 對(duì)比實(shí)驗(yàn)

為定量分析本文方法的變化檢測(cè)效果，選取基于SVM 的變化檢測(cè)方法、基于卷積神經(jīng)網(wǎng)絡(luò)（convolutional neural network，CNN）的變化檢測(cè)方法、基于解混和CNN 的變化檢測(cè)方法（general end-to-end 2-D CNN，GETNET）、基于低秩約束的變化檢測(cè)方法（non-local low-rank sparse unmixing，NLLRSU）、基于張量回歸的變化檢測(cè)方法（generalized tensor regression，GTR）以及本文方法去除非局部張量約束的變化檢測(cè)方法（collaborative sparsity tensor，CST）與本文方法CSNLRT（collaborative sparsity non-local low-rank tensor）進(jìn)行比較。圖7～圖9 分別顯示了這些方法在Farmland 數(shù)據(jù)集和Urban area in San Francisco City（UASFC）兩組數(shù)據(jù)集上的變化檢測(cè)結(jié)果。從圖中可以看出，基于SVM 的變化檢測(cè)方法和基于CNN 的變化檢測(cè)方法在空間上未進(jìn)行任何約束，檢測(cè)結(jié)果中噪聲較多。而GETNET 方法、基于低秩空間約束方法NLLRSU 以及本文方法因使用了空譜聯(lián)合的約束，檢測(cè)結(jié)果具有較強(qiáng)的抗噪性。表2 給出了這些方法在各個(gè)數(shù)據(jù)集中的OA 值和Kappa 值，從結(jié)果中可以看出，本文方法在3 組數(shù)據(jù)集中都能取得最優(yōu)的OA 值和Kappa 值。從OA 值和Kappa 值中可以看出，本文所提的張量表示方法可以有效地從兩時(shí)相高光譜圖像中準(zhǔn)確檢測(cè)出變化的部分。此外，本文所提CSNLRT 方法比CST 方法在OA 以及Kappa指標(biāo)上都取得較好的結(jié)果，這說(shuō)明本文所提變化檢測(cè)模型中的非局部低秩的約束項(xiàng)在抑制噪聲提高檢測(cè)精度上起到了非常好的效果。由于本文方法利用了非局部張量低秩建模以及迭代求解，本文方法的運(yùn)行時(shí)間較長(zhǎng)。

圖7 不同方法在Farmland 數(shù)據(jù)集上的檢測(cè)結(jié)果Fig.7 Detection results of different methods on Farmland dataset

圖8 不同方法在UASFC1 數(shù)據(jù)集上的檢測(cè)結(jié)果Fig.8 Detection results of different methods on UASFC1 dataset

圖9 不同方法在UASFC2 數(shù)據(jù)集上的檢測(cè)結(jié)果Fig.9 Detection results of different methods on UASFC2 dataset

表2 不同方法在3 組數(shù)據(jù)集上的實(shí)驗(yàn)結(jié)果Table 2 Experimental results of different methods on 3 datasets

4 總結(jié)

本文針對(duì)高光譜變化檢測(cè)問(wèn)題，提出一種基于協(xié)同稀疏與非局部低秩張量的高光譜圖像變化檢測(cè)方法。該方法首先求得前后時(shí)間點(diǎn)的高光譜差分圖像，再根據(jù)差分圖像中圖像塊的非局部分布特點(diǎn)，提取不同的非局部相似張量簇。然后基于協(xié)同稀疏正則化和低秩正則化建立變化檢測(cè)模型，利用極小化方法對(duì)模型求解得到表示系數(shù)。最后根據(jù)表示系數(shù)求得張量在不同類(lèi)別中的投影殘差，進(jìn)而判斷該張量是否發(fā)生變化。在Farmland 和UASFC 數(shù)據(jù)集中進(jìn)行實(shí)驗(yàn)，本文方法的OA 和Kappa 值都優(yōu)于其他方法所得結(jié)果，表明本文方法在兩時(shí)相高光譜數(shù)據(jù)變化檢測(cè)任務(wù)中可以取得較好的檢測(cè)精度。