徐 超，詹天明

南京審計大學信息工程學院，南京211815

1 引言

隨著成像光譜技術的快速發(fā)展，高光譜遙感在近二十年內取得了重大突破。高光譜遙感是以數(shù)百個窄頻率且連續(xù)的光譜通道記錄反射光成像，這些窄頻率覆蓋可見光、近紅外和短波紅外波段，光譜分辨率可達納米數(shù)量級[1-3]。高光譜圖像是三維的立方體數(shù)據(jù)，其每個像素點是一個由向量表示的光譜，每個光譜波段代表每個光譜帶上反射值的輻射度[4-6]。三維的立方體則包含在這些連續(xù)波段上的一系列二維空間圖像[7-8]。因此，高光譜圖像較好地體現(xiàn)了場景的空譜聯(lián)合信息?；谶@一特點，高光譜圖像可以更好地對地物中的不同種類進行區(qū)分以及對特定目標進行識別。因此，高光譜圖像在地質學、生態(tài)學、地貌仿真、農林科技、巖石分類、大氣科學、環(huán)境監(jiān)測等研究領域得到了廣泛的應用，為政府相關部門以及企業(yè)提供了有力的決策支持[9-10]。

圖像中，有些具有特定形狀的像素點與其周圍環(huán)境具有較大差異的目標被稱為異常目標[11]。由于高光譜圖像中每個像素點包含豐富的光譜信息，與彩色圖像和多光譜圖像相比，高光譜圖像更適用于異常檢測，這在許多軍事和民用應用領域中引起了人們的極大興趣。由于異常目標往往非常小且以一種低概率形式出現(xiàn)，因此異常目標很難通過任何的監(jiān)督手段或者肉眼觀察來檢測。從理論角度分析，異常檢測可以看作是將圖像中的像素點分為異常和背景兩個類別的問題。但是在沒有先驗光譜信息的前提下準確檢測到與周邊背景有一定差異的目標是非常具有挑戰(zhàn)性的工作[12]。

近幾十年來，許多專家學者對異常檢測問題進行深入分析，提出了許多異常檢測方法。在高光譜異常檢測研究領域，應用最為廣泛的是由Reed和Yu提出的Reed-Xiaoli(RX)方法[13]。該方法假設背景服從多維高斯分布，RX算子利用多維高斯分布的概率密度函數(shù)測量待測像素屬于背景的概率值。所得的廣義似然比測試結果其實是輸入測試像素的光譜矢量與其鄰域內背景均值的馬氏距離。但在實際應用場景中，由于噪聲和異常像素的存在，用多維高斯分布描述復雜的背景存在較大的誤差。因此許多改進算法紛紛被提出，例如基于正則化的RX算法[14]、基于分割的RX算法[15]、基于線性濾波器的RX算法[16]、核RX算法[17]等。此外，一些沒有基于RX算法的異常檢測方法也被廣泛提出，其中基于表示學習的異常檢測方法在近幾年受到了廣泛關注。這類方法假設高光譜特征可由字典線性表示。通過對表示系數(shù)進行約束，可以得到不同的檢測算子。其中，比較知名的表示學習檢測算法有稀疏表示的異常檢測算法[18]、協(xié)同表示的異常檢測算法[19]，以及低秩表示的異常檢測算法[20]。

基于低秩表示的異常檢測算法將高光譜圖像分為背景和異常兩個部分，并對背景部分進行低秩約束。異常區(qū)域則可通過原始圖像和恢復背景區(qū)域的殘差求得。但是，高光譜圖像中存在大量的混合像元以及噪聲的干擾，該方法會將部分屬于背景的像素點錯誤劃分為異常目標。為此，本文提出一種基于低秩和全變差（total variation，TV）正則化約束的高光譜異常檢測方法。首先，針對混合像元的問題，利用線性和非線性解混方法提取高光譜圖像的豐度圖，將豐度圖與原始圖像進行融合得到待檢測的圖像。其次，在融合數(shù)據(jù)中構建背景區(qū)域的字典，并構建圖像的低秩表示模型。然后，對背景的表示系數(shù)進行低秩約束，對表示異常的光譜維進行TV正則化約束，建立異常檢測的低秩和TV正則化模型。最后，通過最優(yōu)化目標函數(shù)，求解分類結果得到異常檢測結果。在三組典型的高光譜圖像中進行對比實驗，結果表明本文方法與現(xiàn)有方法相比可獲得更好的高光譜異常檢測性能。

2 基于低秩表示的異常檢測

設X∈?m×n×B表示高光譜圖像，對空間域進行拉升并轉置，其可以表示為I∈?B×N，其中N=m×n。在高光譜圖像中，隸屬于異常的像素點與隸屬于背景的像素點是不同的，同時隸屬于背景的像素點之間應存在較強的相關性。因此，文獻[11]將矩陣I分解為背景和異常兩部分：

I=DS+E（1）

其中，DS表示背景部分，D是由隸屬于背景的像素點構成的字典，S是表示系數(shù)，而E則表示異常部分。

在已知I和D的前提下，求解S和E即可實現(xiàn)異常檢測的目的。文獻[11]認為背景像素點與其鄰域的像素點相關性較強，應隸屬于同一子空間，因此表示系數(shù)S應該服從低秩約束。而異常像素點與其他像素點有較強的差異性，因此可對E加上稀疏約束。高光譜異常檢測轉化為求解如下優(yōu)化問題：

其中，rank(·)表示秩函數(shù)，λ是調節(jié)參數(shù)，||·||2,1表示矩陣l21范數(shù)。

式（2）中的優(yōu)化問題是NP難問題，文獻[11]指出該優(yōu)化問題求解等價于求解如下的凸優(yōu)化問題：

其中，||·||*表示矩陣核范數(shù)。當目標函數(shù)達到最優(yōu)化，異?？捎蒃的最優(yōu)值的二范數(shù)求得：

其中，E*表示E的最優(yōu)值。若T(Ii)大于設定的閾值時，當前像素點Ii則隸屬于異常。

3 本文方法

Fig.1 Flow chart of proposed method圖1 本文方法流程圖

為克服混合像元以及噪聲對異常檢測的影響，本文提出基于低秩與TV正則化的異常檢測方法，圖1顯示了本文方法的整體流程。首先，利用光譜線性和非線性解混方法對高光譜圖像進行解混，得到的豐度圖像和原始圖像進行融合形成待檢測混合數(shù)據(jù)。其次，對挖掘背景圖像特征，構建背景字典。然后，建立異常檢測模型，并進行優(yōu)化求解。最后，求解異常部分的二范數(shù)得到異常檢測結果。

3.1 解混

高光譜圖像的空間分辨率不高導致圖像中存在很多由不同地物組成的混合像元，影響后續(xù)的分類和識別精度。為提升高光譜圖像處理的性能，可通過光譜解混獲得亞像素級的信息。目前，高光譜解混通過對每個像素的光譜進行分解，將混合像元分解為多個純像元以及所占比例。高光譜解混主要分為線性解混和非線性解混。線性解混假設混合像元是由多個純像元即端元的線性組合而成。而非線性解混則假設混合像元是由多個端元非線性組合而成。在已知端元的前提下，線性解混和非線性解混的目標都是求解端元對應的豐度，從而達到解混目的。目前，線性解混和非線性解混在高光譜圖像應用領域各有優(yōu)勢，為綜合利用它們這些優(yōu)勢，本文對高光譜圖像進行線性和非線性解混，利用解混得到的豐度信息對原始高光譜圖像進行信息補償以提升后續(xù)異常檢測性能。

本文利用文獻[21]中的高光譜解混方法對圖像進行線性和非線性解混。設Bl=[b1,b2,…,bC]∈?B×C表示端元矩陣，表示雙線性端元矩陣，Al,An∈?m×n×w分別是線性和非線性解混的豐度矩陣，文獻[21]中的解混目標方程為：

其中，?和χ是非負常數(shù)，O是稀疏噪聲，Z(i,j),k=

在線性和非線性解混結束后，將得到的豐度圖Al,An∈?m×n×w疊加到原始高光譜圖像X中，得到待檢測的混合數(shù)據(jù)Xl=[X;Al;An]∈?m×n×(B+2w)。該混合數(shù)據(jù)與原始的高光譜圖像相比，在光譜維上增加了亞像素級的特征，可以更好地用于不同地物的鑒別和分類中。

3.2 構建字典

字典的構建在異常檢測任務中非常重要。在基于稀疏表示的目標檢測方法中，由背景字典和目標字典組成的全局字典是已知的。但在異常檢測中，字典是未知的，需要分析原始數(shù)據(jù)的背景和異常差異，構建完備的背景區(qū)域字典。由于異常部分在圖像中所占的比例非常小，因此本文首先用文獻[22]中所提的ERS（entropy rate superpixel）超像素分割方法對混合數(shù)據(jù)Xl進行分割，得到多個同質區(qū)域。ERS方法對高維數(shù)據(jù)分割速度快，且能設置生成的超像素個數(shù)，超像素邊緣貼合度較好，在高光譜處理領域得到了廣泛的應用。其次，將像素個數(shù)較小的超像素從數(shù)據(jù)中剔除，并將剩余的超像素按照其在圖像中的順序進行排列。然后，利用K均值算法對排列好的數(shù)據(jù)進行分類，得到K個類別。借鑒文獻[11]中的預檢測方法，首先將所有類別中像素個數(shù)小于一定閾值的直接剔除，然后對每個類別計算其均值μi和協(xié)方差Σi，再對每個類別中的像素計算預檢測值：

3.3 建模與優(yōu)化

對Xl進行展開得到其二維矩陣的表示形式Il，與文獻[11]類似，Il可由背景部分和異常部分分別表示：

其中，Sl是混合數(shù)據(jù)背景部分字典的表示系數(shù)，El是混合數(shù)據(jù)的異常部分。

為增強異常檢測結果，在模型中對異常部分的光譜維引入TV正則項，結合對Sl的低秩約束，本文的變化檢測模型如下：

其中，模型中的第二項是對表示系數(shù)進行稀疏性約束，第四項是光譜維TV正則項，dh表示垂直方向的梯度算子。該模型結合了背景表示系數(shù)的低秩全局結構信息和稀疏局部結構信息，同時對異常部分進行了TV正則化約束，可以更好地對異常進行檢測。

式（8）的最優(yōu)化可用交替迭代算法[23-24]進行求解。因為式（8）對每個變量都是凸函數(shù)，因此其迭代解是收斂的。首先對目標方程進行分裂，引入中間變量J、P1和P2，滿足J=Sl，P1=dhEl。原始優(yōu)化問題轉化為求解如下最優(yōu)化問題：

3.4 算法步驟

本文提出的基于低秩與TV正則化高光譜異常檢測算法總體流程如下所示：

算法1基于低秩與TV正則化高光譜異常檢測算法

4 實驗結果與分析

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

為驗證本文方法對高光譜異常檢測的效果，使用三組真實高光譜圖像中的包含小目標的圖像塊進行仿真實驗。第一組高光譜圖像是由機載可見光/紅外成像光譜儀在美國加利福尼亞州圣地亞哥上空拍攝的。該圖像空間分辨率為3.5 m/pixel，尺寸為400×400，成像波長范圍從370到2 610 nm，共包含224個光譜波段，去掉38個水汽、低信噪比和信息缺失波段后，剩余186個波段用于實驗[11]。圖2（a）顯示了該數(shù)據(jù)的偽彩色圖像。第二組高光譜圖像是選自HYDICE高光譜數(shù)據(jù)集。該圖像空間分辨率為1 m/pixel，尺寸為307×307，光譜分辨率為10 nm，共包含227個光譜波段，去掉45個水汽和低信噪比波段后，剩余182個波段用于實驗[20]。圖2（b）顯示了該數(shù)據(jù)的偽彩色圖像。

Fig.2 Real hyperspectral image圖2 真實高光譜圖像

在第一組真實高光譜圖像中，選擇其左上角包含飛機的圖像塊驗證本文方法。圖3（a）顯示了左上角圖像塊，其空間大小為100×100。圖3（b）顯示了其對應的真實異常圖，包含3個共57個像素點的異常目標。在第二組真實高光譜圖像中，選擇右上角包含汽車和屋頂?shù)膱D像塊驗證本文方法。圖4（a）顯示了右上角圖像塊，其空間大小為80×100。圖4（b）顯示了其對應的真實異常圖，其中共包含9種異常目標。圖5（a）顯示了圣地亞哥機場中的停機坪仿真?zhèn)尾噬珗D像。空間分辨率為100×100。在該仿真圖像中植入了16個異常目標，圖5（b）顯示了對應的真實異常圖像。

Fig.3 Image block in upper left corner of the first hyperspectral image圖3 第一個高光譜圖像左上角圖像塊

Fig.4 Image block in upper right corner of the second hyperspectral image圖4 第二個高光譜圖像右上角圖像塊

Fig.5 Simulated image block of airport in the first hyperspectral image圖5 第一個高光譜圖像停機坪仿真圖像塊

4.2 參數(shù)分析

本文方法在字典構建時主要參數(shù)包括預檢測時的類別個數(shù)K，構建每個類別子字典的p值參考文獻[11]分別取K=15,p=20。在預檢測前進行分類并將類別中像素個數(shù)小于一定閾值的直接剔除，主要是提高字典構建以及后續(xù)異常檢測的效率。該閾值根據(jù)圖像大小設置，針對本文中的實驗數(shù)據(jù)，閾值可設為20至50。本文實驗中分別設置20、30、40和50。發(fā)現(xiàn)不同值對檢測結果影響不大，主要原因是在分類結果中每個類別的像素個數(shù)是固定的，去除像素個數(shù)較小的那些類別即可。因此，可通過分類結果將每個類別中的像素個數(shù)進行排序，再按照排序結果動態(tài)設置閾值。

本文實驗發(fā)現(xiàn)除K值和p值外，超像素的個數(shù)對字典構建效果也有較大影響。圖6～圖8分別顯示了圖3（a）、圖4（a）和圖5（a）選取100、200、300、400，300、400、500、600以及300、400、500、600個超像素的分割結果。從分割結果中看出，圖3（a）選取300個超像素時異常目標的分割結果較好，圖4（a）選取500個超像素時異常目標的分割結果較好，而圖5（a）選取500個超像素時異常目標的分割結果較好。

Fig.6 Different scales superpixel image in Fig.3(a)圖6 圖3(a)不同尺度超像素圖

Fig.7 Different scales superpixel image in Fig.4(a)圖7 圖4(a)不同尺度超像素圖

Fig.8 Different scales superpixel image in Fig.5(a)圖8 圖5(a)不同尺度超像素圖

本文異常檢測模型中包含參數(shù)α、β和λ。參考文獻[11]α和β值分別設為α=0.1，β=0.1。固定α和β值后，用AUC（area under curve）曲線顯示不同λ值對檢測結果的影響。圖9顯示了兩組數(shù)據(jù)集中選取不同λ值的AUC值。從圖中可以看出，當值λ選取過大時，檢測精度有所下降，這是因為λ值較大，模型中光譜鄰域信息的約束增強。而當λ值較小時，TV項的約束對檢測結果起不到應有影響。因此，λ值的合理范圍設定在[0.000 1,0.001 0] 。本文實驗中，設λ=0.000 5。

Fig.9 AUC curve obtained by different λ values圖9 不同λ 值得到AUC曲線圖

4.3 對比分析

為驗證本文方法在高光譜圖像中的異常檢測效果，在兩個數(shù)據(jù)集中與global-RX[13]、SRD（sparse representation-based detector）[25]、RPCA-RX（robust principal component analysis-RX）[13,26]以及LRASR（low-rank and sparse representation）[11]算法所取得的檢測結果進行對比。同時為分析線性解混和非線性解混的性能，用線性解混和非線性解混分別與本文3.3節(jié)所提異常檢測模型結合進行異常檢測（簡稱LU-LRTV（linear unmixing-low rank and total variation）和NLU-LRTV（nonlinear unmixing-low rank and total variation）方法）。圖10、圖11和圖12分別顯示了不同方法在圖3（a）、圖4（a）和圖5（a）上的異常檢測結果。從圖中可以看出本文所提方法由于利用了解混后的豐度圖以及在檢測模型中增加了低秩和TV正則化約束，其結果不僅可以較為準確地檢測出異常區(qū)域，此外從結果中可以看出若只用線性解混或者非線性解混，其檢測結果都沒有本文所用線性-非線性融合的檢測效果好。說明線性解混和非線性解混在高光譜異常檢測任務中各有優(yōu)勢，將其融合使用則可以得到最高的檢測精度。

Fig.10 Anomaly detection results of different methods in Fig.3(a)圖10 不同方法對圖3(a)異常檢測結果圖

Fig.11 Anomaly detection results of different methods in Fig.4(a)圖11 不同方法對圖4(a)異常檢測結果圖

Fig.12 Anomaly detection results of different methods in Fig.5(a)圖12 不同方法對圖5(a)異常檢測結果圖

為定量分析不同方法的檢測精度，表1顯示了不同算法在兩組數(shù)據(jù)集上檢測結果的AUC值，AUC值越高說明檢測精度越好。從表中結果可知本文方法在實驗所用三組高光譜圖像異常檢測中都取得了比其他方法更優(yōu)的結果。同時，對異常目標較大的檢測，本文方法的檢測精度最高，AUC值可以達到0.992 0。表2列出了本文方法在三組數(shù)據(jù)集中的檢測速度，在100×100×186以及80×100×172大小的高光譜圖像變化檢測任務中其速度分別在30 s和15 s左右，可以看出本文方法的檢測速度較快，可滿足實際應用需求。

Table 1 AUC values of different methods on three datasets表1 不同方法在三組數(shù)據(jù)集上的AUC值

Table 2 Detection time of proposed method on three datasets表2 本文方法在三組數(shù)據(jù)集上的檢測時間s

5 結論

針對高光譜圖像異常檢測問題，本文提出了一種基于低秩正則化約束的高光譜異常檢測方法。該方法主要分為四個步驟：（1）端元提取后對高光譜圖像進行線性和非線性解混，將得到的兩組豐度圖像與原高光譜圖像進行融合；（2）根據(jù)背景區(qū)域在融合數(shù)據(jù)中的特征并利用超像素分割結果構建圖像背景的字典，同時構建圖像的低秩表示模型；（3）由背景和異常目標各自特點，建立異常檢測正則化模型；（4）對模型進行優(yōu)化求解，得到異常檢測結果。在兩組真實高光譜圖像數(shù)據(jù)集中進行實驗，并與經典的高光譜圖像異常檢測方法的檢測結果進行對比分析，結果表明本文方法可獲得較優(yōu)的高光譜異常檢測的性能。