1 引言

高光谱遥感技术已经普遍用于地表材料识别和地球观测,与传统的全色和多光谱遥感图像相比,高光谱图像(HSI)可以为地表物质提供连续的光谱曲线^[1-3]。高光谱异常通常指观测的目标或图像与同质背景具有明显不同的光谱特征。一般来说,异常目标在整幅图像中出现的概率比较低,而且和周围像素的光谱差异较大。近年来,高光谱异常检测蓬勃发展,已成功应用于矿产勘察、边界监测、搜救等领域^[4]。

近几十年,大量的异常检测算法被提出。Reed等^[5]提出的RX算法是高光谱异常目标检测领域的一个基石工作,它假设背景模型服从多元正态分布,并通过计算待检测像素和周围的背景像素之间的马氏距离来判断该像素是否异常。根据估计的统计信息是全局还是局部,RX算法可以分为:全局RX(Global RX,GRX)和局部RX(Local RX,LRX)。RX系列的算法存在以下两个缺点:1) 大多数实际获取的数据并不满足多元正态分布统计假设;2) 背景一般会受到异常目标的污染。

为了解决不能准确估计背景统计信息的问题,基于表达类的算法被提出,可以分为两大类:基于稀疏表达^[6]和基于协同表达^[7]。基于稀疏表达和基于协同表达的算法都没有对数据的统计分布进行假设。基于稀疏表达的算法首先是学习背景字典,然后通过背景字典表达待检测像素,如果重构误差很大则证明待检测像素是异常的,否则为背景像素。基于协同表达的算法和基于稀疏表达的算法相比,它只是用待检测像素周围的像素来表达待检测像素。如果背景中存在噪声,或者被异常目标所污染,这类算法就存在很高的虚警率。

为了解决上述问题,文献[ 8-10]中采用低秩稀疏矩阵分解算法(LRaSMD)将矩阵分解为低秩矩阵、稀疏矩阵和噪声矩阵。其中低秩矩阵对应背景矩阵,稀疏矩阵可以用来检测异常目标,而且高光谱数据的低秩属性可以通过低秩稀疏矩阵分解算法更好地挖掘。通过低秩稀疏矩阵分解算法,可以得到比较干净的背景,因为它已经将噪声和异常矩阵都从原始数据中分离出去了,这就很好地解决了基于稀疏表达算法遇到的问题。利用LRaSMD进行高光谱异常检测,由于检测决策依赖于稀疏分量,稀疏分量可能包含稀疏参数的较大值处的一些非异常像素,所以它们可能会使检测性能变差,并具有更多的误报。

因此,本文提出了一种基于低秩稀疏矩阵分解和稀疏字典表达(LRaSMD-SR)的高光谱异常目标检测算法。LRaSMD-SR算法,通过LRaSMD算法获取比较干净的背景,通过稀疏表达的方式从背景中构建背景字典模型,最后通过计算重构误差来检测异常目标。

2 低秩稀疏矩阵分解

LRaSMD技术可将一个矩阵分解为低秩矩阵、稀疏矩阵和噪声矩阵。由于HSI中的异常目标是低概率出现的,故可认为这些异常目标是稀疏的。因此,Sun等^[8]采用LRaSMD分解出来的稀疏矩阵进行异常目标检测,并取得了非常不错的效果。对HSI数据而言,相邻像素具有相似的光谱特性,其光谱特征的相关度很高,即可认为HSI的背景部分具有低秩特性^[11]。因此,通过LRaSMD分解还可得到低秩矩阵,可认为其对应背景矩阵。通过获取的背景矩阵,可对背景模型进行精确地估计,从而有效解决了噪声和异常对背景模型的污染问题。

高光谱图像X∈R^M×P可通过LRaSMD分解为低秩的背景矩阵B∈R^M×P、稀疏的异常矩阵S∈R^M×P和具有独立分布的高斯噪声矩阵N∈R^M×P,其中M代表整幅HSI中的总像素数,P代表HSI中的波段数。因此LRaSMD在HSI领域的分解公式为

X=B+S+N。(1)

虽然近几年有很多的LRaSMD优化算法被提出来,但是绝大多数的方法依然采用GoDec算法来求解低秩背景分量和稀疏分量。GoDec是一种快速近似算法,它通过约束背景矩阵B的秩和稀疏矩阵S的稀疏度来控制重建模型的复杂度。GoDec算法通过最小化公式(2)中的分解误差来求解背景矩阵B和稀疏矩阵S:

minB,SX-B-SF2,s.t.rank(B)≤r,card(S)≤kM。(2)

式中: ·F代表范数;rank表示矩阵秩,card表示矩阵分量个数;r表示背景矩阵B秩的最大值;k表示稀疏矩阵S的稀疏度,反映了图像中的稀疏分量,通常将其定义为S的l₀范数。背景可近似地表示为几个基向量的线性组合,基向量的数量等于背景分量矩阵的秩。一般来说,随着迭代次数的增加,分解误差会单调减少。因此,(2)式中的优化问题可以转换为交替求解以下两个参数的子问题:

Bt=argminrank(B)≤rX-B-St-1F2,(3)St=argmincard(S)≤kMX-Bt-1-SF2,(4)

式中t为矩阵维数。利用基于双边随机投影算法(BRP)低秩近似理论^[12]来求解(3)式,假设

Y1=XA1Y2=XTA2,(5)

式中:A₁∈R^B×r和A₂∈R^M×r都是随机矩阵。A₁可以通过MATLAB中randn函数获得,randn函数可以产生服从标准正态分布的随机矩阵;A₂可以通过A₂=Y₁=XA₁求得。基于BRP算法可以得到:

B=Y1(AT2Y1)-1YT2。(6)

(4)式中S_t取决于X-B_t-1的硬阈值,即

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

式中:P_Ω(·)是矩阵在输入集Ω上的投影,Ω是 X-Bt-1的前kM个最大项的非零子集。

Godec算法步骤如下:

1) 输入数据矩阵X∈R^M×B,容错系数ε,最大迭代次数M_iter;

2) 初始化B₀=X,S₀:=0,t:=0,A₁=randn(B,r);

3) 重复进行M_iter次直到 X-Bt-StF2/ XF2<ε

(a) t:=t+1;

(b) Y₁=(X-S_t-1)A₁,A₂=Y₁,Y₂=(X-S_t-1)^TA₂;

(c) 如果r( A2TY₁)<r 则r:=r( A2TY₁),重复(b)步骤;

(d) B_t=Y₁( A2TY₁)^-1Y2T;

(e) S_t=P_Ω(X-B_t-1)。

3 基于低秩稀疏矩阵分解和稀疏字典表达的高光谱异常目标检测

在高光谱异常目标检测任务中,由于异常情况的多样性,学习异常目标的字典存在很大挑战,并且难以实现,因此采用背景的字典进行目标重构。稀疏字典表达模型认为一个待检测样本可以近似用背景训练样本和目标训练样本来线性表达^[13]。将由背景训练样本和目标训练样本组成的矩阵作为字典D。对于待检测样本可根据目标和背景类别的训练样本来对其进行如下建模:

X≈Dα,(8)

式中:X=[x₁x₂… x_P]^T为具有P个光谱带的测试像素;D为由背景和目标的训练样本组成的字典;α为稀疏分量,其对应于D中各样本的权重。稀疏分量α可表示为

α=argminX-Dα2+λα1,(9)

其中λ为重构误差项和正则项之间的平衡参数。在高光谱异常检测的问题中,由于目标字典很难获得,因此,(9)式中的字典,只采用背景部分,而此部分可以通过LRaSMD算法获得。如果一个待检测像素不能用背景字典表达,这个像素就是异常目标;如果可以用背景字典表达,那就是背景。因此,可以通过求解重构误差的算法来进行异常目标检测。重构误差计算公式可表示为

r(X)=X-Dα2。(10)

通过LRaSMD去掉大部分噪声和异常,可以获得相对比较干净的背景。基于已经获得的干净背景,进行背景字典学习,具体的背景字典学习算法如下。

输入:数据矩阵X,背景矩阵B,字典基的个数N,字典个数T;

输出:学习到的背景字典D_T;

初始化:λ=0.005,Q=200,A₀:=0,B₀:=0,随机从背景矩阵B中选取Q个样本作为初始背景字典A₀,B₀;

for t=1:T

求解稀疏系数:α_t=arg min xt-Dt-1α2+λα1;

更新字典:A:=A_t-1+α_tαtTB:=B_t-1+XαtT;

for j=1∶N

更新第j列:

μ_j:= 1Ajj(b_j-Dα_j)+d_j

d_j:= 1max(μj2,1)μ_j

end

Dt=Δargmin1t∑i=1t12xi-Dαi2+λαi1

end

4 实验及分析

从定性和定量两个角度对算法的性能进行评价,实验使用一组模拟和两组真实的高光谱图像对算法进行验证。

4.1 数据描述

第1组是模拟数据,模拟数据集有105个光谱带,大小为200 pixel×200 pixel。背景包括5种不同的材料,这些材料选自美国地质调查局植被光谱数据,分别是风滚草、干长草、草坪草、鼠尾草和黑色的刷叶。模拟场景的假彩色图像的上半部分由干燥的长草和草坪草组成,而图像的下部则是黑毛叶和鼠尾草的混合而成。在合成过程中,每种材料的混合比率是随机的。模拟场景上半部分的异常目标是鼠尾草,下半部分的异常目标是干长草,并且风滚草作为异常目标散布在整幅图。这3种异常目标是纯像元而非混合像元。异常目标的大小分别是1 pixel×1 pixel、2 pixel×2 pixel和3 pixel×3 pixel。模拟数据由12个稀疏分布在背景中的异常目标组成,模拟数据如图1所示。

图 1. 模拟数据。(a)伪彩色图像;(b)真实图像

Fig. 1. Simulated data. (a) Pseudo color image; (b) ground truth

下载图片 查看所有图片

第2组是真实数据:AVIRIS数据库。AVIRIS飞机数据集由美国加州圣地亚哥的机载可见红外成像光谱仪(AVIRIS)收集。AVIRIS高光谱图像是一个大小为100 pixel×100 pixel的区域,光谱波段数为224。光谱波长范围从370~2510 nm。考虑到吸水区域、低信噪比和坏波段,以下波段(1~6,33~35,97,107~113,153~166和221~224)被去除,只有189个波段信息被保留,图像中的3架飞机被认为是异常。AVIRIS飞机数据及其相应的真实图像如图2所示。

第3组数据也是真实数据,使用的数据是美国德克萨斯州某城郊住宿区的HYDICE高光谱图像。原始高光谱图像的尺寸为307 pixel×307 pixel×210 pixel,光谱波长范围为400~2500 nm。考虑到吸水区域、低信噪比和坏波段,以下波段(1~4,76,87,101~111,136~153和198~210)被删除,剩余162个波段。在本文中只截取右上角150 pixel×150 pixel的高光谱图像用于实验,其中汽车和屋顶等目标被视为异常现象。HYDICE城市数据及其相应的真实图像如图2所示。

图 2. 真实数据。(a) AVIRIS飞机数据伪彩色图像;(b) AVIRIS飞机数据真实图像;(c) HYDICE城市数据伪彩色图像;(d) HYDICE城市数据真实图像

Fig. 2. Real data. (a) Pseudo color image on AVIRIS airplane data; (b) ground truth on AVIRIS airplane data; (c) pseudo color image on HYDICE urban data; (d) ground truth on HYDICE urban data

下载图片 查看所有图片

为了验证LRaSMD-SR算法的有效性,需要对算法的性能进行定性和定量的验证。接收机工作特性(ROC)曲线可以定性地表征算法的性能,而AUC(Area Under Curve)为接收机工作特性曲线下所覆盖的面积,使用AUC可以定量对检测性能作评价,AUC的数值越大,表明分类器的工作性能越好。本研究采用ROC曲线和AUC值从定性和定量的角度对LRaSMD-SR算法作有效性评价。为了证明LRaSMD-SR算法的性能,选取4个算法进行对比实验,分别是:全局RX(GRX)、局部RX(LRX)^[5]、低秩稀疏矩阵分解算法(LRaSMD)和稀疏表达算法(SRD)^[14]。

4.2 模拟数据实验

模拟数据检测结果的可视化参见图3。从图3可以看出,LRaSMD-SR算法能够检测出所有的异常目标,而且异常目标的强度值都很大。GRX也可以检测出所有的异常目标,但是很明显它所标记的异常目标的强度值比LRaSMD-SR低,这也证明了LRaSMD-SR算法的优越性。SRD算法的性能是比较差的,因为它只检测出一部分的异常目标。LRX的性能比SRD好,但是又比LRaSMD和LRaSMD-SR算法差。

图 3. 模拟数据集上检测结果的可视化图。(a) GRX;(b) LRX;(c) LRaSMD;(d) SRD;(e) LRaSMD-SR

Fig. 3. Visualization of detection results on simulated datasets. (a) GRX; (b) LRX; (c) LRaSMD; (d) SRD; (e) LRaSMD-SR

下载图片 查看所有图片

可视化视图虽然可以直观地看出算法性能的优劣,但是它不能精确地评价算法,因此ROC曲线被用来作进一步的对比分析,如图4所示。从图4(a)可以明显看出,GRX曲线被LRaSMD-SR曲线覆盖,LRX和LRaSMD曲线与LRaSMD-SR曲线部分重合,LRaSMD-SR算法性能比其他所有算法都好,因为其ROC曲线在其他算法的ROC曲线上面。不同算法在模拟数据上的AUC值可参考表1的第2列。GRX和LRaSMD-SR算法的AUC值都是1,这个结果和观察所得的可视化图像和ROC曲线的结果一致。而且,所有算法在模拟数据上的AUC值都比较高,这是因为模拟数据中异常目标和背景的差异性较大。

图 4. ROC曲线。(a)模拟数据;(b) AVIRIS飞机真实数据;(c) HYDICE城市数据

Fig. 4. ROC curves. (a) Simulated data; (b) AVIRIS airplane true data; (c) HYDICE urban data

下载图片 查看所有图片

4.3 AVIRIS飞机真实数据实验

所有算法在真实的AVIRIS飞机数据上的可视化如图5所示。由图5可以看出,GRX漏检测的异常目标很多。SRD虽然比GRX检测出的异常目标要多,但是还是存在一定的漏检情况。LRX、LRaSMD和LRaSMD-SR相比于其他2个算法都取得了不错的效果,基本都能将异常目标检测出来。但是LRaSMD-SR算法检测出来的异常目标的形状和真实图像最接近。LRaSMD检测的异常目标强度值比较低,而LRX稍微存在一点漏检情况。所有算法在真实的AVIRIS飞机数据上的ROC曲线如图4(b)所示。可以明显看出,LRaSMD-SR算法ROC曲线在其他算法ROC曲线的上面,这可以证明LRaSMD-SR算法性能比其他算法好。但是ROC曲线有交叉部分,所以难以判断其他对比算法的优劣性。因此,为了更加准确地对算法进行评价,计算了所有算法AUC值见表1的第3列。LRaSMD-SR的AUC值是0.9908,这是所有算法中最高的。LRaSMD-SR算法结合利用的LRaSMD算法产生的背景矩阵作为背景,有效地去除了异常和噪声对背景字典构建的影响。可以从表1中看出,LRaSMD-SR算法的AUC值比SRD提高了不少,SRD算法的AUC值只有0.9389。

图 5. AVIRIS数据集上检测结果的可视化图。(a) GRX;(b)LRX;(c) LRaSMD;(d) SRD;(e) LRaSMD-SR

Fig. 5. Visualization of detection results on AVIRIS dataset. (a) GRX; (b) LRX; (c) LRaSMD; (d) SRD; (e) LRaSMD-SR

下载图片 查看所有图片

4.4 HYDICE城市真实数据实验

所有算法在真实的HYDICE城市数据上的可视化如图6所示。由图6可以看出GRX和LRaSMD存在漏检情况。SRD和LRX虽然比GRX检测出的异常目标要多,但是检测的异常目标强度值比较低。LRaSMD算法检测出来的异常目标的形状和真实图像最接近。所有算法在真实的HYDICE数据上的ROC曲线如图4(c)所示,可以明显看出LRaSMD-SR算法ROC曲线在其他算法ROC曲线的上面,证明LRaSMD-SR检测性能比其他算法都好。不同算法在HYDICE城市数据上的AUC值可参考表1的第4列,LRaSMD-SR算法AUC值为0.9705,是所有算法中最高的。

图 6. HYDICE数据集上检测结果的可视化图。(a) GRX;(b) LRX;(c) LRaSMD;(d) SRD;(e) LRaSMD-SR

Fig. 6. Visualization of detection results on HYDICE dataset. (a) GRX; (b) LRX; (c) LRaSMD; (d) SRD; (e) LRaSMD-SR

下载图片 查看所有图片

4.5 显著性分析

为了分析LRaSMD-SR算法检测效果的显著性,对HYDICE真实数据进行了5组实验,每种算法的AUC结果如表2所示。SRD算法由于在构建背景的过程中是随机选择背景像素,有可能将异常像素作为背景像素加入字典中,因此SRD算法5次实验结果的方差比较大,而LRaSMD-SR算法在构建字典时采用的是LRaSMD方法获得背景集合,因此5次实验结果都比SRD方法好,而且方差小。

F分布的临界值为F_0.05(4,20)=2.87,在表3中通过对5组数据的方差分析计算出F=450>2.87,从而得出LRaSMD-SR算法效果是显著的。总体而言,通过分析ROC曲线和AUC值,LRaSMD-SR和其他算法相比具有出色的检测性能。

5 结论

提出一种基于低秩稀疏矩阵分解和稀疏字典表达(LRaSMD-SR)的高光谱异常目标检测算法,通过LRaSMD获得背景样本矩阵,通过稀疏表达的方式将背景分量作为字典,采用重构误差的算法检测异常点。最后对高光谱图像进行实验验证,证明了LRaSMD-SR算法可以取得非常好的异常检测结果。