基于近似后驗(yàn)信息的高光譜異常檢測(cè)

王強(qiáng)輝，華文深*，黃富瑜，張 炎，嚴(yán) 陽

1. 陸軍工程大學(xué)石家莊校區(qū)電子與光學(xué)工程系，河北 石家莊 050003 2. 中國(guó)人民解放軍31681部隊(duì)，甘肅 天水 741000

引 言

高光譜圖像(hyperspectral image, HSI)是一種包含空間和光譜信息的三維圖像，其光譜分辨率極高，包含了從可見光到中紅外甚至長(zhǎng)波紅外數(shù)百個(gè)連續(xù)波段的光譜信息，像元光譜曲線近似連續(xù)，利用光譜信息進(jìn)行目標(biāo)檢測(cè)具有廣泛的應(yīng)用。根據(jù)是否提前掌握目標(biāo)光譜作為先驗(yàn)信息，目標(biāo)檢測(cè)分為需要先驗(yàn)信息的光譜匹配檢測(cè)和不需要先驗(yàn)信息的異常檢測(cè)。由于目前缺少完備且實(shí)用的光譜庫，先驗(yàn)信息的獲取極為困難，不需要先驗(yàn)信息的異常檢測(cè)成為研究的熱點(diǎn)[1]。

異常檢測(cè)是指通過挖掘圖像中的光譜和空間信息將圖像中存在的異常像元檢測(cè)出來。經(jīng)典的異常檢測(cè)算法源自多光譜圖像的RX算法[2]，是由Reed和Yu提出的基于廣義似然比檢驗(yàn)(generalized likelihood ratio test, GLRT)的恒虛警率(constant false alarm rate, CFAR)方法。其假設(shè)背景服從多元高斯分布，而異常目標(biāo)偏離這種分布，通過計(jì)算待測(cè)像元與背景的馬氏距離實(shí)現(xiàn)異常檢測(cè)。由于該算法背景的獲取源于整幅圖像，不可避免地將異常成分混入背景中，造成檢測(cè)精度較低、 虛警率較高。為了更為準(zhǔn)確地獲取背景數(shù)據(jù)、 提高檢測(cè)效果，Taitano等提出了局部RX(local RX, LRX)算法[3]。該算法通過建立滑動(dòng)雙窗口模型，使背景的提取更為準(zhǔn)確，提高了檢測(cè)精度[4]。近年來，基于低秩稀疏矩陣分解(low rank and sparse matrix decomposition anomaly detection, LSMAD)的異常檢測(cè)算法[5]有了較大的發(fā)展，該算法將高光譜數(shù)據(jù)矩陣分解為低秩背景矩陣、 稀疏異常矩陣和噪聲矩陣[6-7]，然后利用低秩背景的協(xié)方差矩陣和均值向量構(gòu)造RX算子，提取異常點(diǎn)實(shí)現(xiàn)異常檢測(cè)。

傳統(tǒng)的正交子空間投影(orthogonal subspace projection, OSP)算法[8]是由Harsanyi和Chang提出的一種光譜匹配檢測(cè)算法，通過抑制背景并且最大化信噪比實(shí)現(xiàn)目標(biāo)檢測(cè)。與無異常檢測(cè)算法相比，其檢測(cè)精度更高、 效果更好[9]?？紤]到從光譜庫中獲取目標(biāo)與背景的先驗(yàn)信息較為困難，本文提出了一種基于近似后驗(yàn)信息的異常檢測(cè)(approximate posterior information-based anomaly detection, APIAD)算法，利用矩陣分解算法得到的背景矩陣構(gòu)造近似背景信息與近似目標(biāo)信息，然后將其作為傳統(tǒng)OSP算法的先驗(yàn)信息實(shí)現(xiàn)異常檢測(cè)。充分利用了異常檢測(cè)算法不需要先驗(yàn)信息和光譜匹配檢測(cè)算法精度高的優(yōu)點(diǎn)，提高了對(duì)異常目標(biāo)的檢測(cè)性能。

1 矩陣分解算法

矩陣分解算法將高光譜數(shù)據(jù)矩陣分解為低秩背景矩陣、 稀疏異常矩陣和噪聲矩陣。在高光譜圖像中，異常目標(biāo)出現(xiàn)概率較低，所以異常部分是稀疏的。背景部分具有很大的相關(guān)性，所以是低秩的。矩陣分解模型如式(1)

XM×H=BM×H+SM×H+NM×H

(1)

式(1)中：M為波段總數(shù)，H為像元個(gè)數(shù)。X為高光譜數(shù)據(jù)矩陣，B為低秩背景矩陣，S為稀疏異常矩陣，N為噪聲矩陣。

矩陣分解可以通過去分解(Go Decompose, GoDec)算法實(shí)現(xiàn)：通過對(duì)B的秩和S的稀疏度進(jìn)行約束，限制矩陣分解模型的復(fù)雜度，最終得到低秩背景矩陣和稀疏異常矩陣(包含噪聲矩陣)。GoDec算法通過最小化分解誤差來求解低秩矩陣B和稀疏矩陣S

(2)

式中：‖·‖F(xiàn)代表F范數(shù)，R為低秩矩陣秩的最大值，k為稀疏矩陣的稀疏度，反映了圖像的稀疏分量，可以將其定義為矩陣S的l0范數(shù)，隨著迭代次數(shù)的增加，分解誤差逐漸減小，可以轉(zhuǎn)換為交替求解以下兩個(gè)參數(shù)

(3)

(4)

式中：Bt和St是第t次迭代的最優(yōu)結(jié)果。在運(yùn)算之前，首先對(duì)矩陣進(jìn)行初始化：迭代次數(shù)t設(shè)置為0，Bt設(shè)置為X，St設(shè)置為零矩陣。在迭代過程中，不斷更新Bt和St

St=PΩ(X-Bt),Ω: |(X-Bt)i, j∈Ω|≠0,

(5)

2 基于近似后驗(yàn)信息的異常檢測(cè)算法

2.1 傳統(tǒng)OSP算法

傳統(tǒng)OSP算法是一種基于線性混合模型的光譜匹配檢測(cè)算法[10-11]，需要同時(shí)具備目標(biāo)與背景的先驗(yàn)信息，它能減少背景光譜對(duì)檢測(cè)效果的影響，而且能抑制各種噪聲。對(duì)于任意一個(gè)光譜向量r，可以認(rèn)為其由目標(biāo)、 背景和噪聲三部分組成

(6)

(7)

將投影后的信號(hào)最大化信噪比，得到傳統(tǒng)OSP算法

(8)

2.2 基于近似后驗(yàn)信息的異常檢測(cè)算法

在實(shí)際應(yīng)用中，我們很難獲取目標(biāo)與背景的先驗(yàn)信息，造成檢測(cè)精度更高的OSP算法的應(yīng)用受到很大限制[12]。為了充分利用傳統(tǒng)OSP算法精度高的優(yōu)點(diǎn)，將矩陣分解算法得到的低秩背景矩陣B作為傳統(tǒng)OSP算法的近似背景信息，近似目標(biāo)信息則通過計(jì)算圖像像元與低秩背景矩陣B的馬氏距離得到

(9)

(10)

(11)

式中：r為待測(cè)像元光譜向量，μb和Γb分別為低秩背景矩陣B的均值向量和協(xié)方差矩陣。設(shè)置閾值η，若檢測(cè)器輸出δ(r)>η，則待測(cè)像元r視為目標(biāo)像元，反之視為背景。計(jì)算圖像中所有像元，得到圖像的所有初始異常目標(biāo)像元。取初始異常目標(biāo)像元的光譜均值作為近似目標(biāo)信息

(12)

式中:ri與l表示初始異常目標(biāo)向量與向量個(gè)數(shù)。

將得到的近似背景信息與近似目標(biāo)信息作為傳統(tǒng)OSP算法的先驗(yàn)信息，得到基于近似后驗(yàn)信息的異常檢測(cè)算法

(13)

2.3 算法的步驟與實(shí)現(xiàn)

步驟1 矩陣分解與提取近似背景信息

采用去分解算法將高光譜數(shù)據(jù)矩陣分解為低秩背景矩陣、 稀疏異常矩陣和噪聲矩陣，低秩背景矩陣作為近似背景信息。

步驟2 提取近似目標(biāo)信息

計(jì)算圖像中所有像元與低秩背景矩陣的馬氏距離得到初始異常目標(biāo)像元，取初始異常目標(biāo)像元的光譜均值得到近似目標(biāo)信息。

步驟3 正交子空間投影

將得到的近似背景信息與近似目標(biāo)信息作為傳統(tǒng)OSP算法的先驗(yàn)信息，得到基于近似后驗(yàn)信息的異常檢測(cè)算法，實(shí)現(xiàn)對(duì)整幅圖像異常檢測(cè)。算法流程圖如圖1所示。

圖1 算法流程圖Fig.1 Algorithm flow chart

3 結(jié)果與討論

為驗(yàn)證算法的有效性和可靠性，采用一組仿真數(shù)據(jù)和一組AVIRIS真實(shí)高光譜數(shù)據(jù)進(jìn)行仿真實(shí)驗(yàn)。仿真環(huán)境：CPU處理器為IntelCore i7-7700HQ，主頻為2.80 GHz，安裝內(nèi)存為8.00 GB，軟件為Matlab R2017a。

3.1 仿真數(shù)據(jù)

3.1.1 實(shí)驗(yàn)過程

首先，利用ENVI軟件從數(shù)據(jù)庫中提取三種不同物質(zhì)進(jìn)行實(shí)驗(yàn)，分別為礦物(minerals2151/PortlanditeGDS525)、 干植被(CDE004)、 濕植被(BlueOak)。將礦物作為異常目標(biāo)，兩種植被作為背景。三組物質(zhì)光譜曲線如圖2所示，光譜范圍為0.4～2.5 μm, 光譜分辨率為4 nm，包含526個(gè)波段。

圖2 三組物質(zhì)光譜曲線Fig.2 Spectral curves of three types of material

圖3 20 dB下第10波段灰度圖Fig.3 10th gray-scale map at 20 dB

圖4 真實(shí)目標(biāo)分布圖Fig.4 Real target distribution map

3.1.2 結(jié)果分析

對(duì)不同信噪比下的仿真數(shù)據(jù)分別運(yùn)用RX算法、 LRX算法、 LSMAD算法和本算法進(jìn)行異常檢測(cè)，得到的檢測(cè)結(jié)果灰度圖像如圖5—圖7所示。

通過對(duì)圖5—圖7三組不同信噪比下的檢測(cè)結(jié)果分析可知：RX算法的檢測(cè)結(jié)果較差，幾乎難以檢測(cè)出目標(biāo)，只有在信噪比較高的情況下檢測(cè)出幾個(gè)混合比例較高的目標(biāo)，受噪聲影響較大。這主要是因?yàn)镽X算子是基于馬氏距離的能量檢測(cè)算子，背景的構(gòu)造過程中混入了大量異常目標(biāo)，對(duì)于與背景差異不大、 混合比例較低的目標(biāo)很難進(jìn)行檢測(cè)。LRX算法與RX算法相比背景的構(gòu)造更加合理，混入的異常目標(biāo)較少，但仍然受信噪比影響較大，在信噪比較低的情況下同樣難以檢測(cè)出異常目標(biāo)。LSMAD算法與前兩種算法相比檢測(cè)效果有明顯的提高，在較低信噪比情況下仍然可以檢測(cè)出混合比例較高的異常目標(biāo)，但受信噪比的影響仍然較大。這是因?yàn)長(zhǎng)SMAD算法進(jìn)行異常目標(biāo)檢測(cè)時(shí)雖然利用GoDec算法得到的低秩背景矩陣構(gòu)造背景，但最后的檢測(cè)結(jié)果仍然是通過RX算子得到，目標(biāo)混合比例較低時(shí)仍然難以檢測(cè)。本文提出的算法利用矩陣分解算法得到的背景矩陣和異常目標(biāo)向量均值作為先驗(yàn)信息，然后通過OSP算法實(shí)現(xiàn)檢測(cè)。與其他算法相比，該算法在信噪比較低的情況下仍然可以將混合比例較低的異常目標(biāo)檢測(cè)出來，受信噪比影響較小，這是因?yàn)樵撍惴ǖ玫降谋尘熬仃嚪浅＜儍?，且利用的異常目?biāo)向量均值和最大信噪比構(gòu)造的檢測(cè)器更進(jìn)一步抑制了噪聲的影響，提高了異常目標(biāo)檢測(cè)精度，達(dá)到了較好的檢測(cè)效果。

圖5 20 dB下4種算法檢測(cè)結(jié)果Fig.5 Test results of 4 algorithms under 20 dB

圖6 15 dB下4種算法檢測(cè)結(jié)果Fig.6 Test results of 4 algorithms under 15 dB

圖7 10 dB下4種算法檢測(cè)結(jié)果Fig.7 Test results of 4 algorithms under 10 dB

接收機(jī)工作特性(receiver operating characteristics, ROC)曲線反映了算法的檢測(cè)效果，它表示檢測(cè)率Pd與虛警率Pfa之間的變化關(guān)系，虛警率越低、 檢測(cè)率越高，即ROC曲線向左上方彎曲越大，檢測(cè)效果越好。利用ROC曲線線下面積(AUC值)定量描述檢測(cè)效果，三種信噪比下ROC曲線如圖8所示，AUC值和運(yùn)行時(shí)間如表1所示。

圖8 三種信噪比下的ROC曲線Fig.8 ROC curves for three signal-to-noise ratios

表1 三種信噪比下AUC值和運(yùn)行時(shí)間Table 1 AUC values and running time under three signal to noise ratios

通過分析ROC曲線、 AUC值及運(yùn)行時(shí)間可得：本算法優(yōu)勢(shì)較為明顯，在虛警率較低情況下仍然具有較高的檢測(cè)率，并且受噪聲影響相對(duì)較小，AUC值更大，其運(yùn)行時(shí)間比LSMAD算法僅增加了1.6 s左右，對(duì)時(shí)間消耗并不算太大。

3.2 真實(shí)高光譜數(shù)據(jù)

3.2.1 實(shí)驗(yàn)過程

該數(shù)據(jù)選取的是AVIRIS傳感器拍攝的洛杉磯機(jī)場(chǎng)的部分?jǐn)?shù)據(jù)，波長(zhǎng)范圍是0.37～2.51 μm，去除信噪比較低和水汽吸收較明顯的波段，選擇剩余205個(gè)波段作為研究對(duì)象。截取的圖像區(qū)域大小為50×50像元，存在2個(gè)異常目標(biāo)，其第10波段灰度圖及真實(shí)目標(biāo)分布圖分別如圖9、 圖10所示。

圖9 第10波段灰度圖Fig.9 10th band gray-scale

圖10 真實(shí)目標(biāo)分布圖Fig.10 Real target distribution map

3.2.2 結(jié)果分析

對(duì)該實(shí)驗(yàn)數(shù)據(jù)分別運(yùn)用RX算法、 LRX算法、 LSMAD算法和本文提出的算法進(jìn)行檢測(cè)，設(shè)檢測(cè)到的目標(biāo)像元數(shù)為50，得到的檢測(cè)結(jié)果圖像及其三維概率圖像如圖11所示。

圖11 四種算法檢測(cè)結(jié)果Fig.11 Detection results obtained by four algorithms

分析以上四種算法檢測(cè)結(jié)果：RX算法存在較大的虛警，而LRX算法虛警率較低，這與LRX算法選取局部而不是全局信息作為背景有關(guān)，檢測(cè)效果較好，但其對(duì)異常目標(biāo)的輪廓識(shí)別較差。LSMAD算法在輪廓識(shí)別上有了較大的改善，但也存在虛警率較高的現(xiàn)象。本算法虛警率較低，并且可以近乎準(zhǔn)確地識(shí)別出異常目標(biāo)的輪廓，檢測(cè)效果相對(duì)較好。

圖11從定性的角度表示了檢測(cè)結(jié)果，下面從定量的角度進(jìn)行分析。算法檢測(cè)到的異常目標(biāo)數(shù)、 目標(biāo)像元數(shù)和虛警數(shù)是衡量算法性能的三個(gè)重要指標(biāo)，檢測(cè)到的總的像元數(shù)為50。RX算法檢測(cè)到2個(gè)異常目標(biāo)，目標(biāo)像元數(shù)為33個(gè)，虛警數(shù)為17個(gè)。LRX算法檢測(cè)到2個(gè)異常目標(biāo)，目標(biāo)像元數(shù)為45個(gè)，虛警數(shù)為5個(gè)。LSMAD算法檢測(cè)到3個(gè)異常目標(biāo)，目標(biāo)像元數(shù)為27個(gè)，虛警數(shù)為23個(gè)。本算法檢測(cè)到2個(gè)異常目標(biāo)，目標(biāo)像元數(shù)為48個(gè)，虛警數(shù)為2個(gè)。ROC曲線及四種算法的定量檢測(cè)結(jié)果分別如圖12、 表2所示。

表2 四種算法的定量檢測(cè)結(jié)果Table 2 Quantitative test results of four algorithms

圖12 真實(shí)數(shù)據(jù)ROC曲線Fig.12 Real data ROC curves

通過對(duì)比分析，我們發(fā)現(xiàn)本算法具有極大的優(yōu)勢(shì)，其漏檢率最低，虛警率最低，AUC值最大，明顯優(yōu)于其他三種算法。與LSMAD算法相比優(yōu)勢(shì)更加明顯，僅犧牲了極少的運(yùn)行時(shí)間，克服了LSMAD算法誤將部分背景檢測(cè)為異常目標(biāo)的問題，極大程度上提高了檢測(cè)效果。

4 結(jié) 論

異常檢測(cè)算法不需要任何先驗(yàn)信息，應(yīng)用范圍較廣，但其檢測(cè)精度與需要先驗(yàn)信息的光譜匹配檢測(cè)算法相比仍有一定的差距。針對(duì)這一問題，提出了一種基于近似后驗(yàn)信息的異常檢測(cè)算法。首先利用矩陣分解算法得到低秩背景矩陣，利用低秩背景矩陣得到近似背景信息與近似目標(biāo)信息，然后利用正交子空間投影算法得到最終的檢測(cè)算法，最后將本算法作用于整幅圖像完成異常檢測(cè)。充分發(fā)揮了光譜匹配檢測(cè)算法精度高和異常檢測(cè)算法不需要先驗(yàn)信息的優(yōu)勢(shì)，與RX算法、 LRX算法、 LSMAD算法相比，本算法極大地改善了檢測(cè)性能，提高了檢測(cè)精度，降低了虛警率，并且對(duì)噪聲的抑制效果更好，取得了很好的檢測(cè)效果。這表明了將異常檢測(cè)算法與光譜匹配檢測(cè)算法相結(jié)合在目標(biāo)檢測(cè)領(lǐng)域可以充分發(fā)揮兩組算法的優(yōu)點(diǎn)，具有較大的發(fā)展空間。