CN104504686A - 一种采用局部自适应阈值分割的高光谱图像异常探测方法 - Google Patents

Abstract

本发明属于高光谱图像异常目标检测领域，具体涉及一种采用局部自适应阈值分割的高光谱图像异常探测方法。读入高光谱数据；采用非线性KRX算子对高光谱数据处理得到检测结果灰度图；对检测的灰度图像进行预处理消除背景干扰；将图像分成若干个m×n的子图像对每一个子图像进行阈值迭代法选取阈值；用得到的自适应阈值T_i对子图像二值化，遍历整幅图像进而得到最终检测结果。本发明提出了针对大面积背景干扰抑制的算法，利用形态学滤波开运算的结构元素将背景干扰有效的提取和消除；采用迭代法计算出局部最佳阈值，无需大量试验验证求取阈值可大大减少实际处理中的工作量，提高了异常目标检测的效率和准确率。

Description

一种采用局部自适应阈值分割的高光谱图像异常探测方法

技术领域

本发明属于高光谱图像异常目标检测领域，具体涉及一种采用局部自适应阈值分割的高光谱图像异常探测方法。

背景技术

目标检测是高光谱遥感图像应用的一个重要方向。根据是否需要先验知识，现有的目标检测技术可分为需要先验目标信息或已知背景的目标检测算法和无需任何先验信息的异常目标检测算法。由于在很多实际处理中缺少足够的光谱先验信息，因此，无需任何先验信息的高光谱图像异常检测更符合实际应用的需求。

RX算子是目前使用最广泛的异常探测算法之一，作为经典算法，它源自多光谱图像检测，并由Reed和Xiaoli Yu在1990年最早提出来。该算法是在一些简化的假设条件下，构造似然比检测算子直接对高光谱图像进行处理，这使得检测结果会产生较高的虚警。传统的RX算子是基于线性统计模型推到得到的，它忽略了高光谱数据上百个波段间丰富的非线性信息，导致在实际应用中检测效果不理想。为此，Kwon等人进一步地提出了一种基于核的非线性RX改进算法，该算法利用核函数将原始输入空间的非线性分析映射为高维特征空间的线性分析，更好地利用高光谱图像波段之间的非线性统计特性。然而，当存在较大面积不感兴趣的背景干扰时，KRX算子往往会出现将干扰误判成目标的现象。新型背景抑制核RX(Opening-operation Kernel RX,OKRX)算法，该算法针对KRX检测的灰度值图像大面积亮背景干扰的提取，从根本上解决这种干扰引起的虚警问题。

在实际异常检测应用中，阈值的自适应选取尤为重要。这是由于对检测结果的获得必须是具有实效性的，通过人工分割时，是对全局进行最佳阈值的选取，需要多次实验验证才能确定。随着现代遥感技术的飞速发展，高光谱遥感在处理中的地物信息更加丰富，这使得通过人工分割方法获取最佳阈值变得非常困难。另一方面，对于仅存在于局部范围内而淹没于全局背景中的异常目标或一些较弱的异常目标，全局阈值分割方法将失效，给异常目标检测性能带来了巨大的影响。采用局部自适应阈值分割的方法可以减少实际处理中的工作量，提高了异常目标检测的效率和准确率。为了准确地进行高光谱遥感图像的异常探测，同时解决仅存在于局部范围内而淹没于全局背景中的异常目标检测，全局阈值分割将失效这一问题，本发明提出了一种局部自适应阈值分割的高光谱异常目标检测方法(Opening-operation localadaptive threshold Kernel RX,OAKRX)。

发明内容

本发明的目的在于提供一种无需通过反复大量试验获取最佳阈值，降低算法在实际处理中的工作量，具有更强的实用性的一种采用局部自适应阈值分割的高光谱图像异常探测方法。

本发明的目的是这样实现的：

(1)读入高光谱数据；

(2)采用非线性KRX算子对高光谱数据处理得到检测结果灰度图：

(2.1)选择同心双层窗，对原始高光谱图像数据归一化，然后根据合成图像大小和目标分布，选择目标检测窗口设为3×3像素，背景信息提取窗口大小设为11×11像素；

(2.2)选择现有的核函数中高斯径向基核函数RBF，选用径向基核，其中x_i为核函数中心，即目标窗内的待检测像元，σ²为核函数的宽度参数，x为背景窗内待检测像元；

k(x·x_i)＝exp(||x-x_i||²/σ²)

(2.3)采用非线性KRX算子对原始高光谱数据处理得到检测结果灰度图F；

(3)对检测的灰度图像进行预处理消除背景干扰：

(3.1)采用一个直径r大于最大异常目标的圆形结构元素B对检测结果灰度值图像进行形态学开运算，对检测结果灰度图F进行邻域连接，将相邻的噪声连接起来构成一块区域，消去图像中的异常目标，获得图像背景干扰的灰度图像G：

G＝FοB

(3.2)灰度图像G中的背景干扰体现为高亮度信息，根据图像上的点p通过对灰度图像进行灰度分解得到的一个的二值图像存在灰度为m∈[0,M-1]使得M为灰度图像中的最大灰度级：

${\hat{G}}_p=\left\{\begin{array}{c}1,m\left(p\right)<m\\ 0,m\left(p\right)\&GreaterEqual;m\end{array}\right.$

利用分解得出的二值图像与原灰度图像F矩阵的Hadamard积，消去图像背景：

$\tilde{F}={(F*\hat{G})}_{\mathrm{ij}};$

(4)将图像分成若干个m×n的子图像

(5)对每一个子图像进行阈值迭代法选取阈值：

将图像中的像素点最大值和最小值的均值作为初始阈值T_i0：

$T_{i0}=(\max \left({\tilde{F}}_i\right)+\min \left({\tilde{F}}_i\right))/2;$

将初始阈值T_i0作为阈值把图像的像素点分成两部分，即图像和图像

$\left\{\begin{array}{cc}{\tilde{F}}_{i1}={\tilde{F}}_{i1}^1,& {\tilde{F}}_{i1}^1>T_{i0}\\ {\tilde{F}}_{i2}={\tilde{F}}_{i2}^1,& {\tilde{F}}_{i2}^1\&le;T_{i0}\end{array}\right.;$

再取两部分图像的灰度均值作为新的阈值：

$T_{i1}=(\mathrm{mean}\left({\tilde{F}}_{i1}\right)+\mathrm{mean}\left({\tilde{F}}_{i2}\right))/2,$

计算新的阈值T_i1代替T_i0，重复迭代过程直至T_ik+1＝T_ik为止；

对阈值进行预判决，即存在一个为大于0的实数ε，当ε→0，有或判决子图像属于背景，求取的阈值设置为子图像中最大灰度值，否则取迭代停止时的T_ik作为该子图像的分割阈值，其中T_i为对应子图像的分割阈值：

(6)用得到的自适应阈值T_i对子图像二值化，遍历整幅图像进而得到最终检测结果。

本发明的有益效果在于：：

1)提出了针对大面积背景干扰抑制的算法，利用形态学滤波开运算的结构元素将背景干扰有效的提取和消除；

2)采用局部阈值分割的方法，有效的解决了仅存在于局部范围内而淹没于全局背景中的异常目标的分离问题；

3)采用迭代法计算出局部最佳阈值，无需大量试验验证求取阈值可大大减少实际处理中的工作量，提高了异常目标检测的效率和准确率。

附图说明

图1采用局部自适应阈值分割的高光谱图像异常探测算法流程图；

图2真实的高光谱图像数据，其中(a)是原始高光谱的第20波段灰度显示，(b)是真实地物分布；

图3模拟的高光谱图像数据，其中(a)是模拟数据的第20波段图像，灰白色方块部分表示待检测异常目标，不规则灰白色部分代表大面积的背景干扰，灰黑色部分表示背景信息；(b)是模拟数据的真实目标分布；

图4真实高光谱数据和模拟数据分别采用KRX和结合形态学背景抑制的OKRX检测结果，其中(a)是真实数据的KRX检测灰度图像，(b)利用形态学来运算提取大面积异常干扰，(c)消除干扰的OKRX检测灰度图，(d)是模拟数据数据的KRX检测灰度图像，(e)利用形态学来运算提取大面积异常干扰，(f)消除干扰的OKRX检测灰度图；

图5真实高光谱数据和模拟数据的三种检测算法的检测结果二值化对比图，其中(a)是真实数据KRX检测结果，(b)消除干扰的OKRX检测结果(c)是采用局部自适应阈值分割方法OAKRX检测结果，(d)、(e)、(f)分别是(a)、(b)、(c)对应的模拟数据的检测结果图；

图6真实数据采用KRX和结合形态学开运算背景抑制局部自适应阈值分割OAKRX的检测灰度峰度3D图，其中(a)是图5(a)对应的检测结果3D图，(b)是图5(c)对应的检测结果3D图；

图7模拟数据采用KRX和结合形态学开运算背景抑制局部自适应阈值分割OAKRX的检测灰度峰度3D图，其中(a)是图5(d)对应的检测结果3D图，(b)是图5(f)对应的检测结果3D图。

具体实施方式

下面结合附图对本发明做进一步描述。

本发明针对现有的高光谱异常检测算法无法消除大面积异常背景干扰这一问题，首先采用了形态学开运算对背景干扰提取预处理，利用矩阵的Hadamard积消除干扰，然后引入局部窗将经过处理后的灰度图像分为若干个小图像，对每一个子图像采用阈值迭代法进行自适应阈值选取和判决，使得处理时无需反复大量试验，降低算法在实际处理中的工作量，从而实现了阈值的自适应选取的高光谱异常探测。

本发明不仅适用于核RX算子的使用，同时适用于高光谱图像异常检测的其他算子。并具有很强的实用性。

为实现上述的发明目的，本发明采用下述的技术方案：

1、读入高光谱数据。

2、采用非线性KRX算子对高光谱数据处理得到检测结果灰度图。具体步骤如下：

1)选择同心双层窗，对原始高光谱图像数据归一化，然后根据合成图像大小和目标分布，选择目标检测窗口设为3×3像素，背景信息提取窗口大小设为11×11像素。

2)选择适合的核函数。从已有的文献可知，现有的核函数中高斯径向基核函数(RBF)具有较好的正态分布特性，因此选用径向基核，其表达式如公式(1)，其中x_i为核函数中心，即目标窗内的待检测像元，σ²为核函数的宽度参数，x为背景窗内待检测像元。

k(x·x_i)＝exp(||x-x_i||²/σ²)    (1)

3)采用非线性KRX算子对原始高光谱数据处理得到检测结果灰度图F。

3、对检测的灰度图像进行预处理消除背景干扰。具体步骤如下：

1)研究高光谱异常探测的局部双层窗模型。从已有的文献可知，现有的高光谱局部异常探测主要是通过一个同心双层窗的滑动进行探测，该滑动窗是以待检测像元为中心的，其内层窗和外层窗分别对应可能的目标区域和背景信息区域，其中内层窗的大小取决于异常目标的大小，外层窗的尺寸应远远大于异常目标大小。

2)根据上述的局部双层窗检测模型，设计可用于背景干扰提取的结构元素。其设计目的主要是使得经过形态学开运算后仅留下含大面积干扰的背景信息。采用一个直径r略大于最大异常目标的圆形结构元素B对检测结果灰度值图像进行形态学开运算，即对F进行邻域连接，将相邻的噪声连接起来构成一块面积较大区域，消去图像中的异常目标，从而获得图像背景干扰的灰度图像G。

G＝FοB(2)

3)G中的背景干扰体现为高亮度信息，根据图像上的点p我们定义通过对灰度图像进行灰度分解得到的一个的二值图像存在其灰度为m∈[0,M-1]使得。M为灰度图像中的最大灰度级。

${\hat{G}}_p=\left\{\begin{array}{c}1,m\left(p\right)<m\\ 0,m\left(p\right)\&GreaterEqual;m\end{array}\right.---\left(3\right)$

利用分解得出的二值图像与原灰度图像F矩阵的Hadamard积，消去图像背景，即：

$\tilde{F}={(F*\hat{G})}_{\mathrm{ij}}---\left(4\right)$

4、考虑到KRX检测的灰度图像存在各处的对比度不同、背景变化的情况，很多时候利用全局阈值难以将目标和背景分开，通过局部窗将图像分成若干子图像，即将图像分成若干个m×n的子图像

5、以及建立阈值T的迭代方程，可以实现灰度图像局部自适应阈值的二值分割，具体分析如下：

1)将中的像素点最大值和最小值的均值作为初始阈值T_i0，如公式(5)所示。

$T_{i0}=(\max \left({\tilde{F}}_i\right)+\min \left({\tilde{F}}_i\right))/2---\left(5\right)$

将Ti0作为阈值把的像素点分成两部分，即和

$\left\{\begin{array}{cc}{\tilde{F}}_{i1}={\tilde{F}}_{i1}^1,& {\tilde{F}}_{i1}^1>T_{i0}\\ {\tilde{F}}_{i2}={\tilde{F}}_{i2}^1,& {\tilde{F}}_{i2}^1\&le;T_{i0}\end{array}\right.---\left(6\right)$

再取两部分图像的灰度均值作为新的阈值。

$T_{i1}=(\mathrm{mean}\left({\tilde{F}}_{i1}\right)+\mathrm{mean}\left({\tilde{F}}_{i2}\right))/2---\left(7\right)$

计算新的阈值T_i1代替T_i0，重复迭代过程直至T_ik+1＝T_ik为止。

考虑到其中子图像可能会出现完全属于背景，采用对阈值进行预判决，即存在一个为大于0的实数ε，当ε→0，有或时，判决子图像属于背景，此时求取的阈值应设置为子图像中最大灰度值，否则取迭代停止时的T_ik作为该子图像的分割阈值。如公式(8)所示，其中T_i为对应子图像的分割阈值。

6、用得到的自适应阈值T_i对子图像二值化，遍历整幅图像进而得到最终检测结果。得到基于局部自适应阈值分割的高光谱图像异常探测算法。

本发明采用了形态学开运算对背景干扰提取预处理，不同于传统异常检测中直接进行二值分割的方法，结合背景抑制的核方法可以有效的去除来源于不感兴趣的大面积背景干扰引起的虚警，具有更好的抗噪干扰性能。

提出了局部自适应阈值，不同于传统异常检测的需要多次实验验证确定最佳阈值的选取方法，自适应阈值仅利用待分割像元之前的像元灰度信息差异进行阈值迭代求解，针对每一个局部窗中的像元建立T的迭代方程，通过反复迭代至收敛，使得算法可针对较弱异常有较好的检测效果，大大提高了算法检测性能。

结合核异常检测算子对高光谱图像进行异常探测，即将上述得到的迭代方程用于KRX算子中，得到采用局部自适应阈值分割的高光谱图像异常探测。

本发明所述的采用局部自适应阈值分割的高光谱图像异常探测算法，不同与传统的异常探测算法，该方法提出了局部自适应阈值分割的方法。结合形态学滤波方法提取大面积背景杂波干扰并对其进行抑制，利用Hadamard积滤除背景干扰，然后将处理后的灰度值图像利用局部窗分成若干个子图像，利用阈值迭代方程求得的每一个子图像的阈值进行判决得到二值结果，遍历整幅图像。有效的解决了仅存在于局部范围内而淹没于全局背景中的异常目标的分离问题，无需大量实验选取最佳阈值，降低了数据处理的工作量。利用合成高光谱数据和真实数据进行了异常检测仿真实验，在降低虚警和取得较好的检测效果的同时，实现了高光谱图像异常检测最佳阈值的自适应阈值选取。

下面结合附图对本发明的具体实现步骤做更详细地描述，结合图1,采用局部自适应阈值分割的高光谱图像异常探测方法的算法流程图,其具体实施过程如下：

1.读入高光谱数据，为了验证本发明提出的方法的有效性，分别选用了真实数据和模拟数据进行验证，如图2和图3所示。其中真实的高光谱图像数据取自美国圣地亚哥机场AVIRIS高光谱图像截取的100×100空间大小的图像，该图像可覆盖可见光到近红外波段的连续光谱范围，有224个波段，去除水的吸收带和信噪比较低的波段后，余下的126个波段用于检测，该图像中含有38个异常目标，且每个目标所占的像素数较少。图2即为原始高光谱的第20波段灰度显示和真实地物分布。基于本文算法针对检测地物中存在大面积干扰消除，是取用真实数据AVIRIS高光谱图像中的部分含大面积干扰的背景信息，采用不删除背景像素的叠加像素插入方法进行目标插入，合成40×40空间大小的图像，该图像有126个波段，含有12个异常目标，且每个目标所占的像素数较少，植入加性高斯噪声，信噪比分别为20:1，合成数据图像与目标分布如图3所示。

2.选择高斯径向基(RBF)核函数k(x·x_i)＝exp(||x-x_i||²/σ²)作为核核函数，将原始数据映射到高维特征空间，充分利用非线性统计特性对原始高光谱数据处理得到检测结果灰度图F，采用真实数据和模拟数据进行实验的结果灰度图分别为图4(a)和(d)所示。

3.采用形态学开运算提取大面积背景干扰，利用矩阵的Hadamard积运算消除干扰，得到去噪后的灰度图像图4(c)和(f)分别为真实数据和模拟数据经背景干扰消除后的对应的灰度图。

4.将处理后的灰度图像分成若干子图像

5.针对每一个子图像利用阈值迭代法选取阈值，选择一个初始阈值T_i0，将灰度图像的像元分为两部分，然后计算两部分的平均灰度值：令T_min为小于T_i0的部分均值，T_max为大于T_i0部分的均值。取T_i1＝(T_min+T_max)/2作为新的阈值代替T_i0，如此迭代，直至T_ik收敛，即T_ik+1＝T_ik，停止迭代。

6.判断子图像是否完全属于背景，对迭代出的阈值进行预判决，即当使得或求取的阈值设置为子图像中最大灰度值，否则令T_i＝T_ik。

用得到的阈值T_i对子图像二值化，遍历整幅图像进而得到最终检测结果。图5为KRX、OKRX及OAKRX三种算法检测结果对比图，图5(c)和(f)所示分别为真实数据和模拟数据的检测结果图，为了更直观的看出背景抑制的效果和局部自适应阈值选取的有效性，实验对比KRX和OAKRX两种算法绘制出了检测灰度图像的3D峰度对比图，图6和图7分别为真实数据和模拟数据的3D峰度对比图，从检测结果图和峰度对比图可以看出，本发明的方法有效的实现了异常探测。

Claims (1)

1.一种采用局部自适应阈值分割的高光谱图像异常探测方法：

(1)读入高光谱数据；

(2)采用非线性KRX算子对高光谱数据处理得到检测结果灰度图：

(2.1)选择同心双层窗，对原始高光谱图像数据归一化，然后根据合成图像大小和目标分布，选择目标检测窗口设为3×3像素，背景信息提取窗口大小设为11×11像素；

(2.2)选择现有的核函数中高斯径向基核函数RBF，选用径向基核，其中x_i为核函数中心，即目标窗内的待检测像元，σ²为核函数的宽度参数，x为背景窗内待检测像元；

k(x·x_i)＝exp(||x-x_i||²/σ²)

(2.3)采用非线性KRX算子对原始高光谱数据处理得到检测结果灰度图F；

(3)对检测的灰度图像进行预处理消除背景干扰：

(3.1)采用一个直径r大于最大异常目标的圆形结构元素B对检测结果灰度值图像进行形态学开运算，对检测结果灰度图F进行邻域连接，将相邻的噪声连接起来构成一块区域，消去图像中的异常目标，获得图像背景干扰的灰度图像G：

G＝FоB

(3.2)灰度图像G中的背景干扰体现为高亮度信息，根据图像上的点p通过对灰度图像进行灰度分解得到的一个的二值图像存在灰度为m∈[0,M-1]使得M为灰度图像中的最大灰度级：

${\hat{G}}_p=\left\{\begin{array}{c}1,m\left(p\right)<m\\ 0,m\left(p\right)\&GreaterEqual;m\end{array}\right.$

利用分解得出的二值图像与原灰度图像F矩阵的Hadamard积，消去图像背景：

$\tilde{F}={(F*\hat{G})}_{\mathrm{ij}};$

(4)将图像分成若干个m×n的子图像

(5)对每一个子图像进行阈值迭代法选取阈值：

将图像中的像素点最大值和最小值的均值作为初始阈值T_i0：

$T_{i0}=(\max \left({\tilde{F}}_i\right)+\min \left({\stackrel{\text{~}}{F}}_i\right))\text{/2;}$

将初始阈值T_i0作为阈值把图像的像素点分成两部分，即图像和图像

$\left\{\begin{array}{cc}{\tilde{F}}_{i1}={\tilde{F}}_{i1}^1,& {\tilde{F}}_{i1}^1>T_{i0}\\ {\tilde{F}}_{i2}={\tilde{F}}_{i2}^1,& {\tilde{F}}_{i2}^1\&le;T_{i0}\end{array}\right.;$

再取两部分图像的灰度均值作为新的阈值：

$T_{i1}=(\mathrm{mean}\left({\tilde{F}}_{i1}\right)+\mathrm{mean}\left({\tilde{F}}_{i2}\right))/2,$

计算新的阈值T_i1代替T_i0，重复迭代过程直至T_ik+1＝T_ik为止；

对阈值进行预判决，即存在一个为大于0的实数ε，当ε→0，有或时，判决子图像属于背景，求取的阈值设置为子图像中最大灰度值，否则取迭代停止时的T_ik作为该子图像的分割阈值，其中T_i为对应子图像的分割阈值：

(6)用得到的自适应阈值T_i对子图像二值化，遍历整幅图像进而得到最终检测结果。