CN104463848A - 采用局部因果窗的高光谱图像实时探测方法 - Google Patents

Abstract

本发明涉及一种采用局部因果窗的高光谱图像实时探测方法，其特征在于:步骤1：读入高光谱图像数据，并进行初始化；步骤2：建立局部因果滑动窗模型，设定局部因果滑动窗宽度，所说局部因果滑动窗内仅含有待检测像元之前的已知像元信息；步骤3：利用递归方程得出当前窗口局部相关矩阵逆矩阵；步骤4：利用RX异常算子构造局部实时探测算子，得出当前检测值；步骤5：扫描是否超出图像范围，如果没有超出，则返回步骤3；否则，停止探测，输出检测结果。

Description

采用局部因果窗的高光谱图像实时探测方法

技术领域

本发明涉及一种采用局部因果窗的高光谱图像实时探测方法。

背景技术

目标检测技术是高光谱图像处理的一个重要研究领域。根据有无先验知识条件，目标检测技术可分为已知背景或目标信息的目标检测算法和无需任何先验知识的异常目标检测算法。由于在很多实际情况下，研究者总是没有足够的先验知识来表征目标类别的统计信息，因此，无需先验知识的高光谱图像异常检测得到了越来越广泛的关注。

RX算子是Reed和Xiaoli于1990年提出的一种广义似然比算子，被人们称之为异常检测算法的基准算法而得到了国内外学者的广泛学习和改进。由已有文献可知，异常目标检测可分为两大类，全局异常探测和局部异常探测。全局异常探测是指利用全局统计信息作为依据构造检测算子，比较典型的有经典的全局RX算法，简称为RXD。Chang等人将经典的基于协方差矩阵的RXD进行了改进，提出了基于相关矩阵的R-RXD算子。然而，当异常目标比较微弱或者仅仅在局部范围异常而淹没于全局背景时，高光谱全局异常探测模型将失效。局部异常检测是指利用局部滑动窗口计算局部背景统计信息构造检测算子，例如局部RX算子，又称为LRXD，它是利用一个以待检测像元为中心点的同心双层滑动窗来采集统计信息，实现局部异常检测。

在实际应用中，实时探测尤为重要。这是因为许多异常目标，例如运动目标，它的停留时间是非常短的，会突然出现然后很快就消失。因此，要想实现这种异常目标的检测，其检测过程必须是实时的。另一方面，随着现代遥感技术的飞速发展，高光谱遥感在获取更加丰富地物信息的同时，其大量的数据量也给数据存储、卫星下行传输和后续处理带来了巨大的压力。实现机载实时处理可以减少数据传输量，也就相应减少了数据压缩下传而导致的信息损失，减少了地面数据存储和处理工作量，提高了图像分类和目标检测的效率和实效性。但现有技术中，能够有效实现高光谱图像局部异常探测的方法尚属空白。

发明内容

本发明目的在于提供一种采用局部因果窗的高光谱图像实时探测方法，能够实时、准确地实现高光谱图像局部异常探测。

实现本发明目的技术方案：

一种采用局部因果窗的高光谱图像实时探测方法，其特征在于:

步骤1：读入高光谱图像数据，并进行初始化；

步骤2：建立局部因果滑动窗模型，设定局部因果滑动窗宽度，所说局部因果滑动窗内仅含有待检测像元之前的已知像元信息；

步骤3：利用递归方程得出当前窗口局部相关矩阵逆矩阵；

步骤4：利用RX异常算子构造局部实时探测算子，得出当前检测值；

步骤5：扫描是否超出图像范围，如果没有超出，则返回步骤3；否则，停止探测，输出检测结果。

步骤3中，通过两次利用Woodbury恒等式定理，建立局部相关矩阵的递归方程。

局部因果滑动窗为局部因果阵列滑动窗。

局部相关矩阵表示如下：

$\begin{array}{c}{\tilde{R}}_w(n+1)=\frac{1}{w}\underset{\mathrm{window}}{\underset{i\∈\mathrm{current}}{\mathrm{\Σ}}}{r}_i{r}_i^T=\frac{1}{w}[\underset{\mathrm{window}}{\underset{i\∈\mathrm{previous}}{\mathrm{\Σ}}}{r}_i{r}_i^T-r_{n-w}{r}_{n-w}^T+r_n{r}_n^T]\\ ={\tilde{R}}_w\left(n\right)-\frac{1}{w}{r}_{n-w}{r}_{n-w}^T+\frac{1}{w}{r}_n{r}_n^T\end{array}$

式中，w是局部因果阵列窗的窗宽，表示窗宽为w的局部相关矩阵，是上一状态滑动窗内所有像元的相关矩阵，r_n-w是当前窗口溢出的像元信息，r_n是新进入当前窗口的像元信息。

步骤2中，局部因果窗为局部因果矩阵滑动窗。

本发明具有的有益效果：

本发明为采用局部因果窗的高光谱图像实时探测方法，不同与传统局部异常探测算子，本发明提出了局部因果滑动窗，使得传统不具有因果特性的局部窗模型变成具有因果特性的局部阵列窗，滑动窗口能够只用当前待检测像元之前的的已知像元信息构造局部背景特性，而不需未知信息，实现数据处理的因果化。然后通过两次利用Woodbury恒等式定理，得出局部相关矩阵逆矩阵的递归表达，来更新局部相关矩阵的逆矩阵，而无需重新反复计算矩阵求逆，大大提高了算法运行速度，更有力地实现实时探测。本发明不仅能与经典的RX算子结合使用，同时也可以和高光谱图像异常检测的其他算子结合使用，具有很强的可移植性，更易满足高光谱检测的需求。利用合成高光谱数据进行仿真实验，本发明方在取得准确检测效果的同时，有效实现了高光谱图像数据的实时探测。

本发明首次提出了局部因果滑动窗，滑动窗口能够只用当前待检测像元之前的的已知像元信息构造局部背景特性，而不需要用到以后信息，实现数据处理的因果化；本发明局部相关矩阵求逆只需当前溢出像元r_n-w以及进入像元r_n的信息，而没有用其他像元信息，因此，可大大减少算子所需的存储空间；本发明无需重复计算及可大大降低算法的运算时间，提高算子速度，实现算子实时探测。

附图说明

图1为本发明方法流程图；

图2为传统的局部双层窗模型；

图3为本发明局部因果矩阵窗模型；

图4为本发明局部因果阵列窗模型；

图5为真实的高光谱图像数据及五种地物位置信息；

图6为基于图5合成的高光谱图像数据；

图7为采用经典全局K-RXD和R-RXD异常探测算子后的检测结果灰度图及其相应的3D图；

图8为不同窗口大小的传统LRXD算子检测结果灰度图；

图9为不同因果阵列窗宽度的LRTCARXD算子检测结果灰度图；

图10为LRXD算子和LRTCARXD算子运算时间对比图。

具体实施方式

如图1所示，本发明采用局部因果窗的高光谱图像实时探测方法包括如下步骤：

步骤1：读入高光谱图像数据，并进行初始化；

读入高光谱数据这一图像的大小为m×n，每个像素有L个波段特征，x_i,j为样本数据集合R^L中的样本，代表L维特征空间。

初始化主要是初始化局部相关矩阵。为了避免相关矩阵为病态矩阵，一般应有L个以上像元后才开始计算初始化局部相关矩阵，并且局部窗宽度要大于L。

步骤2：建立局部因果滑动窗模型，设定局部因果滑动窗宽度，所说局部因果滑动窗内仅含有待检测像元之前的已知像元信息；

传统高光谱局部异常探测主要是采用经典局部双层滑动窗，如图2所示，该滑动窗是一个以待检测像元为中心的同心双层窗，其中外层窗是用于提取背景的统计信息，其大小应远远大于异常目标大小，而内层窗则是可能的目标区域，其大小取决于异常目标的大小。由于实时检测一个必要条件就是算子应是因果的，因此这种传统局部双层窗模型由于不满足必要条件，而不能实现实时探测。

如图3所示，本发明局部因果滑动窗可采用局部因果矩阵滑动窗，窗口内仅含有待检测像元之前的像元信息，而不包括未知像元信息，实现窗口的因果特性。图3中，(a)是其演化过程，黑色中心点表示当前检测像元，白色部分是背景信息的组成，灰色部分将不用于局部因果矩阵窗；(b)是窗口滑动时，窗内像元信息的变化。由于所提出的因果矩阵滑动窗口在滑动时会有多个像元信息发生变化(溢出当前窗口或者进入当前窗口)，因此每进入到下一个像元检测而滑动窗口时，都需要一个追踪模板记录当前溢出和进入的像元信息，这使得实时探测变得比较复杂。

基于这一问题，如图4所示，本发明的最佳实施例采用局部因果阵列滑动窗，以下方法步骤均基于局部因果阵列滑动窗。图4中，(a)是其演化过程，黑色中心点表示当前检测像元，白色部分是背景信息的组成；(b)是窗口滑动时，窗内像元信息的变化，当窗口滑动时，局部因果阵列窗只有两个像元发生变化，即溢出一个最远像元，进入一个最新像元。

步骤3：利用递归方程得出当前窗口局部相关矩阵逆矩阵；

1)建立局部相关矩阵的状态方程。

根据图4所示的因果阵列窗口模型，局部相关矩阵可以表示如下：

$\begin{array}{c}{\tilde{R}}_w(n+1)=\frac{1}{w}\underset{\mathrm{window}}{\underset{i\∈\mathrm{current}}{\mathrm{\Σ}}}{r}_i{r}_i^T=\frac{1}{w}[\underset{\mathrm{window}}{\underset{i\∈\mathrm{previous}}{\mathrm{\Σ}}}{r}_i{r}_i^T-r_{n-w}{r}_{n-w}^T+r_n{r}_n^T]\\ ={\tilde{R}}_w\left(n\right)-\frac{1}{w}{r}_{n-w}{r}_{n-w}^T+\frac{1}{w}{r}_n{r}_n^T\end{array}---\left(1\right)$

其中上标用于区分全局相关矩阵，w是局部因果阵列窗的窗宽，表示窗宽为w的局部相关矩阵，是上一状态滑动窗内所有像元的相关矩阵，r_n-w是当前窗口溢出的像元信息，r_n是新进入当前窗口的像元信息。

2)Woodbury恒等式引理

在数学中，尤其是线性代数中，Woodbury恒等式指出，某些具有k阶修正因子的矩阵的逆可以通过原始矩阵逆矩阵及其k阶修正因子来表示。该准则又称为矩阵逆矩阵引理，Sherman-Morrison-Woodbury引理或者直接称之为Woodbury引理。该引理可表示如下：

(A+UCV)^-1＝A^-1-A^-1U(C^-1+VA^-1U)^-1VA^-1       (2)

其中，A,U,C和V都是具有特定大小的矩阵，A是n×n矩阵，U是n×k矩阵，C是k×k矩阵，V是k×n矩阵。在此，如果矩阵C的维数远远小于矩阵A，则该等式的效率要远远高于直接计算(A+UCV)的逆阵。尤其是当C是1×1的单位矩阵式，等式(2)可简化成Sherman-Morrison恒等式，表示如下:

(A+UV)^-1＝A^-1-A^-1U(1+VA^-1U)^-1VA^-1       (3)

这个等式在已知A^-1且需要计算(A+UCV)^-1的时候是很有用的。在后续异常算子的计算中，利用该恒等式可有效的避免相关矩阵和协方差矩阵求逆的计算，实现矩阵求逆的递归更新。式(3)也可表示成如下形式：

(A-UV)^-1＝A^-1+A^-1U(1-VA^-1U)^-1VA^-1      (4)

3)通过两次利用Woodbury恒等式引理，推出局部相关矩阵的递归更新公式

由等式(1)，将前两部分看作是一个整体，记为等式(5),从而局部相关矩阵(1)可简化为(6):

${\widetilde{\stackrel{~}{R}}}_w\left(n\right)={\tilde{R}}_w\left(n\right)-\frac{1}{w}{r}_{n-w}{r}_{n-w}^T---\left(5\right)$

${\tilde{R}}_w(n+1)={\widetilde{\stackrel{~}{R}}}_w\left(n\right)+\frac{1}{w}{r}_n{r}_n^T---\left(6\right)$

如果想求局部相关矩阵的逆矩阵我们将第一次利用Woodbury恒等式引理,将等式(3)用于(6)可得:

$\begin{array}{c}{{\tilde{R}}_w}^{-1}(n+1)={[{\widetilde{\stackrel{~}{R}}}_w\left(n\right)+\frac{1}{w}{r}_n{r}_n^T]}^{-1}\\ ={{\widetilde{\stackrel{~}{R}}}_w}^{-1}\left(n\right)-\frac{(1/w){{\widetilde{\stackrel{~}{R}}}_w}^{-1}\left(n\right)r_n{r_n}^T{{\widetilde{\stackrel{~}{R}}}_w}^{-1}\left(n\right)}{1+(1/w){r_n}^T{{\widetilde{\stackrel{~}{R}}}_w}^{-1}\left(n\right)r_n}\end{array}---\left(7\right)$

由等式(7)可知,只要得到的逆矩阵我们便可根据等式(7)递推更新得到局部相关矩阵的逆矩阵在此通过第二次利用Woodbury恒等式引理,将等式(4)用于(5)可得:

$\begin{array}{c}{{\widetilde{\stackrel{~}{R}}}_w}^{-1}\left(n\right)={[{\tilde{R}}_w\left(n\right)-\frac{1}{w}{r}_{n-w}{r}_{n-w}^T]}^{-1}\\ ={{\tilde{R}}_w}^{-1}\left(n\right)+\frac{(1/w){{\tilde{R}}_w}^{-1}\left(n\right)r_{n-w}{r}_{n-w}^T{{\tilde{R}}_w}^{-1}\left(n\right)}{1-(1/w)r_{n-w}^T{{\tilde{R}}_w}^{-1}\left(n\right)r_{n-w}}\end{array}---\left(8\right)$

由此，通过两次利用Woodbury恒等式引理，可以得到局部相关矩阵的递归表达。

步骤4：利用RX异常算子构造局部实时探测算子，得出当前检测值；

利用RX异常检测算子得到基于局部滑动因果窗的新型高光谱局部实时探测算法，简记为LRTCARXD(Local Real Time Causal Array RX Detector)：

${\mathrm{\δ}}^{\mathrm{LRTCARXD}}\left(r_n\right)={r_n}^T\tilde{R}{\left(n\right)}^{-1}{r}_n$

步骤5：扫描是否超出图像范围，如果没有超出，则返回步骤3；否则，停止探测，输出检测结果。

下面结合仿真实验进一步说明本发明的有益效果。

图5是真实的高光谱图像数据及五种地物位置信息。实验数据为内华达州的某一矿区，拍摄于1997年。该AVIRIS高光谱成像仪有224个波段，图像数据大小为350×350。在去除1-3，105-115以及150-170等水吸收波段已经低信噪比波段，189-波段数据用于后续的试验中。在图中，我们给出了五种地物的地理分布：A,B,C,K,M。利用这五种地物光谱特性，设计了合成数据如图6所示。图6是合成的高光谱图像数据，合成数据大小为200×200，有5×5＝25个异常目标，其成分如图6所示，该图可用于检测不同大小的异常目标探测以及亚像元探测。

图7是采用经典全局K-RXD和R-RXD异常探测算子后的检测结果灰度图及其相应的3D图，其中(a)是K-RXD检测结果灰度图，(b)是R-RXD检测结果灰度图，(c)和(d)分别是(a)和(b)对应的检测结果3D图；图8是不同窗口大小的传统LRXD算子检测结果灰度图，其中图(a)双层窗口尺寸为3/15(内窗尺寸/外窗尺寸)，图(b)为5/17；图9是不同因果阵列窗宽度的LRTCARXD算子(本发明)检测结果灰度图，其中图(a)窗宽度为15*15，图(b)窗宽度为19*19；

图10是LRXD算子和LRTCARXD(本发明)算子运算时间对比图。从仿真实验对比看出，本发明能够实时、准确地实现高光谱图像局部异常探测。

Claims (5)

1.一种采用局部因果窗的高光谱图像实时探测方法，其特征在于:

步骤1：读入高光谱图像数据，并进行初始化；

步骤2：建立局部因果滑动窗模型，设定局部因果滑动窗宽度，所说局部因果滑动窗内仅含有待检测像元之前的已知像元信息；

步骤3：利用递归方程得出当前窗口局部相关矩阵逆矩阵；

步骤4：利用RX异常算子构造局部实时探测算子，得出当前检测值；

步骤5：扫描是否超出图像范围，如果没有超出，则返回步骤3；否则，停止探测，输出检测结果。

2.根据权利要求1所述的采用局部因果窗的高光谱图像实时探测方法，其特征在于：步骤3中，通过两次利用Woodbury恒等式定理，建立局部相关矩阵的递归方程。

3.根据权利要求2所述的采用局部因果窗的高光谱图像实时探测方法，其特征在于：局部因果滑动窗为局部因果阵列滑动窗。

4.根据权利要求3所述的采用局部因果窗的高光谱图像实时探测方法，其特征在于：局部相关矩阵表示如下：

$\begin{array}{c}{\tilde{R}}_w(n+1)=\frac{1}{w}\underset{\mathrm{window}}{\underset{i\∈\mathrm{current}}{\mathrm{\Σ}}}{r}_i{r}_i^T=\frac{1}{w}[\underset{\mathrm{window}}{\underset{i\∈\mathrm{previous}}{\mathrm{\Σ}}}{r}_i{r}_i^T-r_{n-w}{r}_{n-w}^T+r_n{r}_n^T]\\ ={\tilde{R}}_w\left(n\right)-\frac{1}{w}{r}_{n-w}{r}_{n-w}^T+\frac{1}{w}{r}_n{r}_n^T\end{array}$

式中，w是局部因果阵列窗的窗宽，表示窗宽为w的局部相关矩阵，是上一状态滑动窗内所有像元的相关矩阵，r_n-w是当前窗口溢出的像元信息，r_n是新进入当前窗口的像元信息。

5.根据权利要求2所述的所述的采用局部因果窗的高光谱图像实时探测方法，其特征在于：步骤2中，局部因果窗为局部因果矩阵滑动窗。