基于光谱角累加的高光谱图像异常检测算法

2013-05-06杜小平夏鲁瑞

中国光学 2013年3期

关键词：波段光谱向量

杜小平，刘明，夏鲁瑞，陈杭

（中国人民解放军装备学院，北京101416）

基于光谱角累加的高光谱图像异常检测算法

杜小平*，刘明，夏鲁瑞，陈杭

（中国人民解放军装备学院，北京101416）

针对传统方法中用作背景的像素中存在干扰像素的问题，提出一种基于光谱角累加的高光谱图像异常检测算法。通过计算测试区域待测像素光谱向量与其他光谱向量之间的夹角，并将其累加，得到图像中每个像素的异常程度;然后使用波段选择预处理方法进一步提高检测性能。HyMap高光谱数据验证表明，在虚警概率设为0.008时，检测概率达到0.73，即在提高异常检测可靠度的同时，降低了虚警概率。

高光谱图像；异常检测；光谱角；波段选择

1 引言

近年来，高光谱遥感引起人们的广泛关注，利用高光谱图像进行目标检测与识别逐渐成为该领域的研究热点。高光谱传感器能够在上百个波段上对地物进行成像，精细记录地物的光谱信息，因此具有光谱鉴别能力。由于反射率反演算法尚不完善，而对于很多应用领域获得先验信息也比较困难，所以不需要任何先验信息的异常检测算法具有重要的研究价值［1-3］。

用传统方法完成高光谱图像异常检测是根据图像中的异常点与背景的统计特征的差异来进行的，如Reed和Xiaoli Yu提出的RX算法［4］，此类方法需要假设背景统计模型，一般采用高斯模型，此外还需要求协方差估计矩阵的逆矩阵，计算量较大且精度受到影响［5］。Harsanyi提出的LPD异常检测算法［6］利用子空间理论，将抑制背景后的数据各波段进行叠加，从而凸显没有被抑制的异常信号，此方法对剩余背景信息非常敏感，受背景的复杂度影响较大。Banerjee提出的基于支持向量域的异常检测算法［7］，用待测像素周围的像素描述背景的几何边界，然后考察待测像素是否为异常，此方法也受背景区域中奇异像素“污染”的影响，窗口的选择需要一定的先验知识［8］。为了解决高光谱图像异常检测中多类背景和背景中奇异像素的影响，本文提出了一种基于光谱角累加的异常检测算法，该方法利用像素之间的相似信息构造异常度衡量因子，不需要任何假设和先验信息即可完成图像异常检测，文中实验验证了算法的有效性。

2 异常检测算法

2.1 光谱角匹配

光谱角用于描述光谱之间的相似程度，在高维空间中，两个向量的夹角越小，两条光谱越相似。将l个波段的光谱响应作为l维空间的矢量，光谱角定义为两个光谱向量的夹角，在高维空间中两矢量T=［t1，t2，…，tl］与R=［r1，r2，…rl］的夹角用反余弦表示为：

上式可进一步写为：

式中，θ值越小，T和R的相似性越大。

从上式可以看出，光谱角与矢量的模无关，即不受向量本身长度的影响。因此，用光谱角度量两个光谱向量的相似性，可以消除或减弱因太阳入射角、地形、坡向和观测角等因素引起的光谱变异的同物异谱现象，此外，光谱角对于乘性干扰具有良好的抗干扰性，不受照度变化的影响［9］。由于所有波段参与计算，虽能够充分利用数据所包含的信息，但也受到来自每个波段噪声的影响，因此有必要剔除受噪声影响严重的波段，提高匹配的精度。

2.2 根据噪声水平选择波段

成像光谱仪在成像过程中会引入许多噪声，用噪声分解的办法对噪声强度进行分析具有一定的难度。然而信号在各波段间具有强烈的相关性，所以可以用去相关法将噪声分离出来［10-11］，通过线性回归去除图像中具有高相关性的信号。回归公式为：

式中：xi，j，k为第i列，第j行，第k波段的像素灰度值，为xi，j，k的线性拟合值，a，b，c，d为线性回归系数。对于xp，k：

去相关后的残差图像为：

式中，ri，j，k为残差值，有：

第k波段图像噪声方差：

式中：M-4为自由度，w，h分别为每一小块图像的宽度和高度。

去相关法既适用于场景均匀的图像，也适用于不均匀的情况。用去相关法估计出每个波段的噪声水平后，选择噪声水平较小的波段参与计算，能够在一定程度上提高计算精度。

2.3 构造异常度因子

在高光谱图像异常检测中，异常定义为与周围背景像素存在光谱差异的像素。高光谱图像中异常通常表现为出露概率低，所占像素少，而背景为大量存在的信号。异常检测的目的是在检测区域中，根据光谱信息找出少量存在的光谱作为异常。

由以上分析，可以用待测像素与检测区域内像素之间的相似性构造异常度因子，考察与待测像素相似的像素数量，并以此作为待测像素的异常程度。如果待测像素为检测区域内大量存在的像素，将待测像素判为背景，反之，则判为异常。

光谱角匹配方法只能衡量两个光谱向量之间的相似度，对于多个光谱之间的关系无法描述，一般目标检测与识别中用光谱角匹配的方法搜索相似像素。现将光谱角匹配的方法用于异常检测中，定义异常度因子为检测区域内待测像素与所有非待测像素夹角之和，即光谱角累加：

式中，θ为待测像素光谱向量与检测区域内其它像素的光谱向量的夹角，n为检测区域内像素数目，l为波段数。T为待测光谱向量，R为检测区域内其他光谱向量。

在式（8）异常度因子计算过程中有以下3点说明：

（1）Ri（i=1，2，…，n-1）中包含检测区域内所有像素。因为待测像素光谱与自身的夹角为零，所以计算时不必除去待测像素;

（2）Ri（i=1，2，…，n-1）来源于局部采样窗口，窗口选择应使窗口内包含足够多的样本像素又不至使计算量过大;

（3）待测像素T周围不必设置安全窗口。

如果式（8）中待测像素T取到检测区域中大量存在的像素（背景），由于类内光谱夹角很小，此时异常度因子γ的各分量中将包含较多的小量，累加结果仍然较小。如果待测像素T取到少量存在的像素，由于类间光谱夹角比较大，此时异常度因子γ的各分量中包含较多的“大量”，此时累加结果较大。γ值会随着待测像素所属类别的像素数占检测区域内总像素数的比例增大而减小，当待测区域内不存在与待测像素相似的像素时，γ值最大，随着检测区域内与待测像素相似的像素的增多，γ逐渐减小，即待测像素的异常程度随之降低。由此可见，以γ作为像素异常度的衡量标准是合理的。

图1 高光谱图像背景与异常在三维空间的分布Fig.1 Distribution of background and abnormality for hyperspectral imagery in 3-dimensional space

噪声和物质本身的纯度影响，使得背景信号在高维空间中表现为一组聚集的点集，而非一个确定的点，异常信号同样如此，如图1所示，图中背景像素与异常像素分布相对聚集，但两组点集中心相距较远。对于异常检测而言，噪声的影响要小于背景的影响，类间光谱向量的夹角α远大于类内光谱向量的夹角β，这样在所构造的异常度因子中会突出类间差异而减弱类内差异的影响。而且，上述类间夹角和类内夹角之差在一个数量级内，这可以减小背景样本中离群像素的影响。例如，图2中将图像的中心点（中心点为背景）作为待测像素，在背景样本中存在两个明显的奇异像素，计算γ时，分量中含有两个类间夹角分量α1，α2。假设窗口内有1 000个样本像素，那么α1，α2至多只占γ值的0.2%，而不会影响待测像素属性的判断，即仍会将中心像素判为背景（非异常）。

图2 包含异常像素的单波段图像Fig.2 Single band imagery with anomaly pixel

2.4 分割异常目标

将计算所得的异常度因子数据做适当的调整，在图3中画出异常度因子直方图，横坐标为异常度因子值，纵坐标为异常度因子值对应的像素数目。

图3 异常度因子直方图Fig.3 Histogram of anomaly degree factor

异常度因子沿横轴增大，较小的异常度因子对应的像素为背景，占有较多像素;较大的异常度因子对应异常，占有较少的像素。多数像素聚集在异常度因子较小的区域，随着异常度因子的增大，像素数目突然变得很小，认为此时的异常度因子已经不在背景区域。一般情况下，异常像素的异常度因子与背景像素的异常度因子差异较大，二者之间的异常度因子值对应背景和异常的混合像素，在高维空间中处于两种地物点集之间的区域，此区域的像素相对稀疏，因此，选择背景右侧第一个极小值处的异常度因子值作为异常检测的阈值，进行异常目标分割。

2.5 异常检测算法具体实现步骤

基于以上的讨论，算法采用局部矩形窗口进行像素采样，窗口尺寸的选择主要依据潜在异常目标的大小来确定［12-14］，具体实现步骤为：

（1）估计高光谱图像每个波段的噪声水平，按从大到小的顺序重新排列各波段图像，剔除受噪声影响较大的波段，余下的波段参与计算。

（2）计算检测区域内的每个像素的异常度因子，从而生成一幅“像素异常度图像”。在这幅图像上，每个像素的灰度值表示像素的异常程度，因此像素灰度值越大，表示像素的异常度越高。

（3）设定阈值分割异常像素，得到异常的二维空间分布图。

3 实验分析

图4 第30波段图像和目标分布图Fig.4 Imagery and targetmaps of the 30th band image

为了验证算法的有效性，利用HyMap高光谱数据进行仿真实验，该图像覆盖了从短波到红外的光谱范围，剔除信噪比低的波段后，余下的126波段参与验证。从HyMap数据中截取部分图像，图像的大小为100 pixel×100 pixel，共126波段，图像的第30波段如图4（a）所示，人为插入待检测异常目标，异常目标光谱曲线从图像的其他区域的典型目标上提取，然后以一定的比例与背景混合置换图像中的像素，目标分布如图4（b）所示，第一行为1像素目标，目标间隔为10像素，第二行为1像素目标，目标间隔为5像素，第三行为2像素目标，第四行为4像素目标，第3、4行的目标间隔为4像素，图像中共包含18个目标。分别执行RX算法和本文算法，实验中将3种算法的检测窗口大小都设置为30 pixels×30 pixels。图5为2种算法检测结果。

图5 RX算法和本文算法检测结果Fig.5 Test results by RX algorithm and proposed algorithm in this paper

从实验结果对比可以看出，本文算法所获得的检测结果要好于RX算法。这是因为RX算法中需假设高光谱数据服从高斯分布，而实际数据往往不能很好地符合该假设。此外，如果在RX算法中用于估计背景协方差的像素中存在着异常像素，使估计的精度降低，会导致检测效果不佳。而本文算法不需要任何假设，用异常度因子作为像素异常程度的量度，能够适应多类背景的情况，同时根据以上讨论所得结论，本文算法在一定程度上克服了背景数据中异常像素的影响。为了更加直观地表示检测结果，图6给出了相同虚警率（0.008）下的分割结果。从图6（a）与图6（b）的对比中可以看出（b）中探测到的目标明显较多。经计算，RX方法在该虚警率下的检测率为0.33，本文方法的检测率为0.73。通过比较可以发现，本文方法抑制了RX方法中存在的很多分散的虚警，但虚警的面积增大。

算法特性曲线（ROC曲线）用于描述不同检测阈值下检测概率Pd与虚警概率Pf之间的变化关系，提供算法检测性能的定量分析。将检测概率Pd定义为检测到的真实目标像素数目与地面真实目标像素数目的比值，虚警概率Pf定义为检测到的虚警像素数目与整幅图像像素数目总和的比值。利用不同的阈值得到不同条件下的Pd和Pf绘制算法的ROC曲线。

图7给出了2种算法的ROC曲线，本文算法ROC曲线全部位于RX算法ROC曲线之上，通过ROC曲线进一步说明了本文算法有较好的性能。由于不需要高斯分布的假设以及使用波段选择和抑制背景的处理方法，本文的算法具有较强的适应性。

图7 RX算法和本文算法ROC曲线Fig.7 ROC curves using RX algorithm and algorithm proposed in this paper

4 结论

本文提出了一种基于光谱角累加的高光谱异常检测算法，用待测像素与检测区域内样本像素的光谱角累加作为异常度因子，用异常度因子衡量待测像素的异常程度，然后设定合理的阈值进行异常检测。构造异常度因子过程中，针对此种方法的局限，用基于噪声水平的波段选择方法加以改进。实验证明，本文算法为一种有效的高光谱图像异常检测算法，在虚警概率设为0.008时，检测概率达到0.73。与高光谱图像异常检测领域的标准检测算法RX算法相比，本文的算法具有较高的可靠性以及较强的适应性。

［1］赵慧洁，秦宝龙，贾国瑞.高光谱遥感系统调制传递函数的在轨测试［J］.光学精密工程，2011，19（6）：1235-1243. ZHAO H J，QIN B L，JIA G R.On-orbit detection formodulation transfer function of hyperspectral remote sensing system［J］.Opt.Precision Eng.，2011，19（6）：1235-1243.（in Chinese）

［2］董超，赵慧洁，王维，等.局部正交子空间投影高光谱图像异常检测算法［J］.光学精密工程，2009，17（8）：2004-2010. DONg CH，ZHAO H J，WANGW，etal..Hyperspectral image anomaly detection algorithm based on local orthogonal subspace projection［J］.Opt.Precision Eng.，2009，17（8）：2004-2010.（in Chinese）.

［3］李欢，周峰.星载超光谱成像技术发展与展望［J］.光学与光电技术，2012，10（5）：38-44. LIH，ZHOU F.Developments of spaceborne hyperspectral imaging technique［J］.Opt.Optoelectronic Technology，2012，10（5）：38-44.（in Chinese）.

［4］REED IS，YU X.Adaptivemultiple-band CFAR detection of an optical pattern with unknown spectral distribution［J］. IEEE T.Acoust.，Speech，1990，38（10）：1760-1770.

［5］ADLER-GOLDEN SM.Improved hyperspectral anomaly detection in heavy-tailed backgrounds［C］//Hyperspectral Image and Signal Processing：Evolution in Remote Sensing，2009WHISPERS′09.Firstworkshopon，2009：1-4.

［6］HARSANYI JC，CHANG C I.Hyperspectral image classification and dimensionality reduction：an orthogonal subspace projection approach［J］.IEEE T.Geosci.Remote，1994，32（4）：779-785.

［7］BANERJEE A，BURLINA P，DIEHLC.A supportvectormethod for anomaly detection in hyperspectral imagery［J］.IEEE T.Geosci.Remote，2006，44（8）：2282-2291.

［8］KHAZAIS，HOMAYOUNIS，SAFARIA，et al..Anomaly detection in hyperspectral images based on an adaptive support vectormethod［J］.IEEE Geosci.Remote Sensing Lett，，2011，8（4）：646-650.

［9］杜培军，谭琨，夏俊士.高光谱遥感影像分类与支持向量机应用研究［M］.北京.科学出版社，2012：29-31. DU P J，TAN K，XIA JSH.Research on Application of Support Vector Machine and classification of Hyperspectral Remote Sensing Imagery［M］.Beijing：Science Press，2012：29-31.

［10］ROGER R E，ARNOLD JF.Reliably estimating the noise in AVIRIS hyperspectral images［J］.International J.Remote Sensing，1996，17（10）：1951-1962.

［11］AIAZZIB，ALPARONE L，BARONTIF，etal..Hyperspectral noisemodelling and estimation hyperspectral data from airborne imaging spectrometers［J］.Anuals of Geophysics，2006，49（1）：3141

［12］穆治亚，魏仲慧，何昕，等.基于NIOSⅡ的红外弱小目标检测算法实现［J］.液晶与显示，2012，（6）：814-819. MU ZH Y，WEIZH H，HE X，etal..Design of dim-small target detection on infrared image based on NIOSⅡ［J］.Chinese J.Liquid Crystals abd Displays，2012，（6）：814-819.（in Chinese）

［13］吴君钦，刘昊，罗勇静.态背景下的运动目标检测算法［J］.液晶与显示，2012，（5）：682-686. WU JQ，LIU H，LUO Y.Algorithm ofmoving object detection in static background［J］.Chinese J.Liquid Crystals abd Displays，2012，（5）：682-686.（in Chinese）

［14］刘翔，周桢.基于分块背景建模的运动目标检测技术［J］.液晶与显示，2011，26（6）：831-835. LIU X，ZHOU ZH.Moving objectives detecting technology based on block processing of background model establishment［J］.Chinese J.Liquid Crystals abd Displays，2011，26（6）：831-835.（in Chinese）

Anomaly detection algorithm for hyperspectral imagery based on summation of spectral angles

DU Xiao-ping*，LIU Ming，XIA Lu-rui，CHEN Hang
（Academy of Equipment，Chinese People′s Liberation Army，Beijing 101416，China）
*Corresponding author，E-mail：dxp8600@163.com

As interference pixels influence the background features in traditionalmethods，a new anomaly detection algorithm for hyperspectral imagery is proposed based on the summation of spectral angles.The anomaly degree of each pixel in the imagery is acquired by calculating the angles between the spectral vector of test pixel and other spectral vectors in the test region，and accumulating the angles.Then，the pretreatmentmethod of band selection is used to further improve detection performance.HyMap hyperspectral data experiments show that the detection probability reaches0.73 when the probability of false alarm is set to be0.008.The reliability of anomaly detection is improved，while the probability of false alarm is reduced.

hyperspectral imagery;anomaly detection;spectral angle;band selection

TP751.1

10.3788/CO.20130603.0325