CN109598701A - 一种基于信息扩展的多光谱遥感影像非监督变化检测方法 - Google Patents

Abstract

本发明涉及一种基于信息扩展的多光谱遥感影像非监督变化检测方法，包括以下步骤：1)基于非线性波段生成的波段光谱信息扩展：分别对前一时相和后一时相的多光谱影像的任意两个原始波段进行波段光谱信息扩展，获取新扩展的非线性人工光谱波段；2)基于多尺度形态学重构的波段空间信息扩展：分别对原始波段和非线性人工光谱波段通过空间重构运算获取空间扩展波段特征集后，进行差分运算得到差分波段集；3)非监督多类变化检测。与现有技术相比，本发明具有实现有限原始波段扩展和空间信息重构、实现多类变化的自动化、非监督检测等优点。

Description

一种基于信息扩展的多光谱遥感影像非监督变化检测方法

技术领域

本发明涉及相多光谱遥感影像处理领域，尤其是涉及一种基于信息扩展的多光谱遥感影像非监督变化检测方法。

背景技术

地表的动态变化从时间维度上揭示了环境的演变过程。地球观测卫星的持续对地观测，可探测人类活动和自然因素对环境及地表覆盖的影响，这对全球变化研究有着极其重要的意义。变化检测技术是对地观测研究中最为重要的技术之一，能够自动或半自动地识别不同卫星观测时间段中，同一地理区域内，两幅甚至更多幅影像之间的实际发生土地覆盖/土地利用变化。在卫星上搭载的众多传感器中，多光谱扫描仪可以在较大的地理范围内获取较高空间分辨率的遥感影像。在过去几十年中，考虑到卫星遥感数据的可获得性，多光谱影像是地表要素识别、提取和变化检测研究与应用中的主要数据源。

常规在使用多时相、多光谱影像进行变化检测应用时，由于数据源本身的限制，如影像波段的质量较差(包括为校准的水汽波段、影像坏道、波段噪声等)、波段间的空间位置不一致性、多时相数据间的不一致性(如云、阴影等干扰)、以及波段光谱范围的限制(如只覆盖可见光波段，无近红外、中场波红外)等，进而在使用极少数波段进行多类变化检测时，原始光谱波段不能充分的反映真实的地物变化信息，从而增加了变化检测的难度。虽然基于核函数或深度神经网络的方法在一定程度上可以提供解决此类问题的途径，但都需要长时间的调参过程和庞大的运算量，且这些方法所需的大量训练样本在实际变化检测应用中都很难获得。因此，需要一种简单、高效且实用的基于波段光谱-空间信息扩展的多光谱遥感影像非监督变化检测方法。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于信息扩展的多光谱遥感影像非监督变化检测方法。

本发明的目的可以通过以下技术方案来实现：

一种基于信息扩展的多光谱遥感影像非监督变化检测方法，包括以下步骤：

1)基于非线性波段生成的波段光谱信息扩展：分别对前一时相和后一时相的多光谱影像的任意两个原始波段进行波段光谱信息扩展，获取新扩展的非线性人工光谱波段；

2)基于多尺度形态学重构的波段空间信息扩展：分别对原始波段和非线性人工光谱波段通过空间重构运算获取空间扩展波段特征集后，进行差分运算得到差分波段集；

3)非监督多类变化检测。

所述的步骤1)中，非线性人工光谱波段的计算式为：

其中，X_1,α、X_1,β分别为B维的前一时相多光谱影像X₁的任意两个原始波段，M₁、R₁分别为前一时相多光谱影像X₁扩展后的非线性人工光谱波段，X_2,α、X_2,β分别为B维的后一时相多光谱影像X₂的任意两个原始波段，M₂、R₂分别为后一时相多光谱影像X₂扩展后的非线性人工光谱波段，α、β为波段号。

所述的步骤1)中，在进行除法运算时，通过添加一常数t来避免出现分母为极小非零值的情况。

所述的步骤2)具体包括以下步骤：

21)给定圆形结构元素尺度i的范围，并对影像X₁进行形态学开重建和闭重建得到形态学剖面结构特征，即：

其中，MPⁱ(X₁)为影像X₁在尺度i下进行开闭运算得到的特征集合，为影像X₁在尺度i下进行开运算得到的特征，为影像X₁在尺度i下进行闭运算得到的特征，为影像X₁第b个波段在尺度i下的形态学开重建结果，为影像X₁第b个波段在尺度i下的形态学闭重建结果；

22)分别对影像X₁、影像X₂以及非线性人工光谱波段M₁、R₁、M₂和R₂采用与步骤21)相同的空间重构运算生成空间扩展波段特征集F_X₁、F_X₂、F_M₁、F_M₂、F_R₁和F_R₂；

23)对空间扩展波段特征集进行差分运算得到最终经过波段光谱-空间信息扩展后的差分波段集E_D。

述的步骤22)中，影像X₁和影像X₂对应的空间扩展波段特征集F_X₁和F_X₂的表达式为：

F_X₁＝[MP²(X₁),...,MP⁶(X₁)]

F_X₂＝[MP²(X₂),...,MP⁶(X₂)]。

所述的步骤23)中，差分波段集E_D的计算式为：

E_D＝[X_D,F_D]

X_D＝X₂-X₁

F_D＝[F_X₂,F_M₂,F_R₂]-[F_X₁,F_M₁,F_R₁]

其中，X_D为原始光谱变化表达，F_D为扩展后的增强变化表达。

所述的步骤3)中，采用压缩变化矢量分析、序列光谱变化矢量分析或迭代加权多变量变化检测进行多类变化检测。

当采用压缩变化矢量分析和序列光谱变化矢量分析时，通过压缩变化矢量表达描述光谱变化，定义变量变化强度ρ和变化方向θ，根据变量变化强度ρ和变化方向θ构建二维极坐标变化表达域，ρ值用以表示变化的可能性，即ρ值越大，越有可能发生变化，θ值表示可能的不同类别的变化，所述的变化强度ρ和变化方向θ的表达式为：

其中，E_D,n为E_D中第n个元素，N为E_D中的元素总数，即扩展波段总维度，r_n为单位常数向量中的第n个元素。

与现有技术相比，本发明具有以下优点：

1)实现有限原始波段的光谱信息扩展：利用邻近波段迭代乘、除非线性函数对原有波段光谱信息进行计算，生成新的伪波段以实现光谱信息的增强与扩展；

2)实现有限原始波段的空间信息重构：利用多尺度形态学重构运算，对原始光谱波段和扩展的光谱波段进行多尺度空间信息重构，以综合多尺度地物目标的形态学表达，实现空间信息的增强与扩展。

3)实现多类变化的自动化、非监督检测，无需训练样本数据和先验知识：采用三种主流非监督多类变化检测方法，包括压缩变化矢量分析、序列光谱变化矢量分析或迭代加权多变量变化检测方法，对波段扩展叠加后的特征影像进行自动化变化检测。

附图说明

图1为本发明的整体方法流程图。

图2为撒丁岛数据集，其中，图(2a)为撒丁岛1995年获取的Landsat TM彩红外合成影像，图(2b)为撒丁岛1996年获取的Landsat TM彩红外合成影像，图(2c)为变化参考图。

图3为盐城市数据集，其中，图(3a)为盐城市2006年获取的模拟波段Hyperion彩红外合成影像，图(3b)为盐城市2007年获取的模拟波段Hyperion彩红外合成影像，图(3c)为变化参考影像。

图4为撒丁岛数据集中不同方法获得的变化检测图，其中，图(4a)为S²CVA基于X_D的变化检测图，图(4b)为IR-MAD基于X_D的变化检测图，图(4c)为M²C²VA基于X_D的变化检测图，图(4d)为S²CVA基于E_D的变化检测图，图(4e)为IR-MAD基于E_D的变化检测图，图(4f)为变化参考图。

图5为盐城市数据集中不同方法获得的变化检测图，其中，图(5a)为S²CVA基于X_D的变化检测图，图(5b)为S²CVA基于E_D的变化检测图，图(5c)为IR-MAD基于X_D的变化检测图，图(5d)为IR-MAD基于E_D的变化检测图，图(5e)为M²C²VA基于X_D的变化检测图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。

如图1所示，本发明提供一种基于信息扩展的多光谱遥感影像非监督变化检测方法，主要由以下三个步骤组成：

(1)基于非线性波段生成的波段光谱信息扩展

原始波段描述影像的一阶统计量，二阶统计信息可用波段间的相关性来表示，从而作为原始波段光谱信息的有效补充。假设X_1,α和X_1,β是B维的前一时相多光谱影像X₁的任意两个波段，两个新扩展的波段集M₁和R₁可通过在两个波段间逐像素进行乘法和除法运算得到，即：

其中，M₁(或R₁)的维度为B²/2-B/2。对后一时相影像X₂进行相同的操作，可得到新扩展的波段集M₂和R₂。对X₁和X₂中的波段进行归一化运算，使得数据的范围保持在[0,1]。在进行除法运算的时候，添加一个较小值常数t，即以避免出现分母为极小非零值的情况。

(2)基于多尺度形态学重构的波段空间信息扩展

第二步旨在对原始光谱波段和步骤(1)中生成的非线性人工光谱波段的基础上进行空间信息的重构，以实现波段空间信息的扩展。通过建立形态学剖面(MP)，考虑多尺度几何结构信息的光谱-空间一致性表达。假设i为给定的圆形结构元素的尺度，对影像X₁可通过形态学开重建(O_R)和闭重建(C_R)得到形态学剖面结构特征，即：

其中MPⁱ(X₁)为在尺度i下进行开闭运算得到的特征集合。为了融合从多尺度形态学剖面特征中提取的信息，结构元素的大小i被设定为2至6连续区间。

F_X₁＝[MP²(X₁),...,MP⁶(X₁)] (6)

利用相同的空间重构运算分别生成空间扩展波段特征集：F_X₂,F_M₁,F_M₂,F_R₁,F_R₂。并分别对两个时相间的原始影像、合成后的光谱-空间波段扩展影像进行差分运算。

X_D＝X₂-X₁ (7)

F_D＝F₂-F₁＝[F_X_D,F_M_D,F_R_D]＝[F_X₂,F_M₂,F_R₂]-[F_X₁,F_M₁,F_R₁] (8)

E_D＝[X_D,F_D] (9)

E_D为最终经过波段光谱-空间信息扩展后的差分波段集，其中包含了原始光谱变化表达X_D以及新增的经过两次波段扩展后的增强变化表达F_D。其中F_D集成了从光谱相异性角度估计的新增光谱变化(即M)，强化了消除光照影响后的光谱变化(即R)，以及变化目标在空间信息重构后的变化表征能力(即F_X_D,F_M_D,F_R_D)。

(3)非监督多类变化检测对前两步中生成的增强后的波段集E_D进行多类变化检测。采用最新的三种主流非监督方法，压缩变化矢量分析(C²VA)、序列光谱变化矢量分析(S²CVA)以及迭代加权多变量变化检测(IR-MAD)进行多类变化检测。

其中，C²VA和S²CVA通过压缩变化矢量表达来描述光谱变化。通过在增强后的波段集上，定义变化强度和变化方向两变量：

其中，E_D,n表示为E_D中第n个元素，N为E_D中的元素总数(即扩展波段总维度)。在C²VA方法中，r_n表示某个单位常数向量中的第n个元素。通过压缩变量ρ和θ可构建一个二维极坐标变化表达域，即：{ρ∈[0,ρ_max]andθ∈[0,π]}，其中ρ值代表了变化的可能性(ρ值越大，越有可能发生变化)，θ值表示可能的不同类别的变化。ρ_max为ρ可能取值的最大值。变化部分中同质聚类块的数量代表了存在的变化类别数。不同于C²VA方法，S²CVA通过序列分析来实现自适应且动态的变化表达。根据奇异值分解将第一特征向量定义为参考向量r，在不同检测的层次只考虑部分的变化矢量，所以通过迭代投影可以将这局部的矢量投影至一个能代表变异信息量最大化的参考方向上去，并且较好地保留不同变化之间的可区分信息。在多光谱数据集中，层次结构的层数相对较少，取决于变化类别的数目和可分离度。对于二值变化信息的检测，沿着变化强度分量方向定义阈值T_ρ，该阈值可通过贝叶斯理论框架下的最大期望算法获得。对于多类变化信息，利用模糊c均值聚类方法对沿着变化方向分量进行聚类，通过50次随机初始化迭代结果进行平均以减少初始化对于聚类结果影响的不确定性。IR-MAD方法通过计算单时相扩展波段特征集的正则变量，并计算它们之间的差异特征，而差异特征彼此之间相互正交。其中，为了优化不变化背景与变化目标之间的分离度，运用了一个迭代定权方法。所有正则差异变量的平方和，构成了一幅描述变化强度的卡方影像。最后，基于该卡方影像用模糊c均值算法聚成两类进行二值变化检测，基于正则差异变量在二值检测的基础上用模糊c均值聚K类进行多类变化检测。

实施例：

1、实验数据

撒丁岛数据集：

第一组真实遥感数据集为1995年9月(X₁)和1996年6月(X₂)分别在意大利撒丁岛获取的两时相Landsat5TM影像。影像大小为412×300像素，经过辐射校正及几何配准，试验用六个光谱波段(波段1-5，7)。该数据包含的变化有露天采石场(C₁)和焚烧场(C₂)面积的扩大以及由于Mulargia湖(C₃)水位上涨引起的水域面积扩张三种变化。前后两时相的假彩色影像如图(2a)和(2b)所示。

盐城市数据集：

第二组真实遥感数据集为2006年5月(X₁)和2007年4月(X₂)分别在中国江苏省盐城市获取的沿海湿地区域农田的高光谱EO-1Hyperion遥感影像。选取原始影像中的一个子区域，大小为220×430像素。从原始242个波段中经过预处理(未校正和噪声波段的去除、坏条带修复、大气校正、几何配准等)选出128波段(即13-53，85-96，103-118，135-164，188-199，和202-218)用于后续实验。为模拟多光谱变化检测情景，在计算R²相关系数矩阵后，将具有相似光谱信息的波段按波长范围分成五个组，每组通过均值合成仿真波段。该数据集包含植被、裸地、水体和土壤间的五种相互变化。前后两时相的假彩色影像、变化参考影像如图(3a)，(3b)和(3c)所示。

2、实验结果

撒丁岛数据集

采用C²VA，S²CVA，IR-MAD三种方法分别对原始影像波段和本发明提出方法进行波段扩展后的影像进行非监督多类变化检测，并与最新的多尺度形态学压缩变化矢量分析(M²C²VA)方法进行对比实验。为评估统计结果，以S²CVA方法检测结果为例，计算原始影像及波段扩展后影像变化检测结果(即Kappa系数)间的z检验值。若计算的|z|值超过1.96，表明在95％的置信水平下，两个检测结果间具有显著差异。表1为差异性统计结果，其中有13个|z|值大小超过1.96(表中加粗表示)，表明不同变化检测结果间具有显著性差异，同样说明融合不同扩展波段的变化检测结果可以进一步优化变化检测结果。

表1利用Kappa系数进行z检验的差异性统计结果(S2CVA方法)

表2为使用不同变化检测方法获得实验结果的总体精度(OA)、Kappa系数、错检、漏检及总误差的像元数。从表中可以看出，相比在原始影像上进行变化检测，文章提出的基于波段扩展的变化检测方法在进行多类变化检测时整体精度有了显著的提升。C²VA方法中整体精度从95.81％提高至98.34％，错误像元个数从5338降至2400个；S²CVA方法中整体精度从95.89％提升至98.35％，错误像元个数从5123降至2395个；IR-MAD方法中整体精度从96.00％提升至97.98％，错误像元个数从4984降至2909个。S²CVA方法在三个方法中获得最高的精度。通过和对比实验结果进行比较，M²C²VA方法尽管也考虑了多尺度光谱-空间信息，文章提出的基于波段扩展的方法依旧取得了更高的整体精度、Kappa系数和较少的错误数，这表明了本文提出方法的有效性。

表2不同变化检测方法的精度及误差统计

盐城市数据集

在盐城市数据集上进行相同的实验，从表3和图5可以看出，使用基于波段扩展的影像进行多类变化检测时，其结果要极大地优于在原始影像上进行变化检测的结果，其在三种变化检测方法中均取得更高的整体精度和更低的总体误差数。其中，相比较于用原始波段，C²VA和S²CVA的实验结果整体精度提高了约3％，IR-MAD的实验结果整体精度提高了约5％，而他们的错误像元个数则显著降低，这表明了使用主流变化检测方法对本文提出的基于波段扩展的影像进行变化检测的有效性。此外，通过对比实验得到的结果可发现，尽管M²C²VA对原始影像进行变化检测时结果优于文中使用的三种方法，但本发明提出的基于波段扩展的影像变化检测结果，比M²C²VA获得了更高的整体精度、Kappa系数和更少的错检测误差。

表3不同变化检测方法的实验精度及错误数

本专利提出了一种全新的基于波段光谱-空间信息扩展的非监督多类变化检测方法，主要包含：1)使用邻近波段迭代乘、除非线性函数对波段光谱信息进行扩展；2)利用多尺度形态学重构对波段空间进行扩展；3)采用三种主流非监督多类变化检测方法对波段扩展后的影像进行变化检测。在两个真实的双时相遥感数据集上的实验结果表明了提出方法的有效性。

Claims (8)

1.一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，包括以下步骤：

1)基于非线性波段生成的波段光谱信息扩展：分别对前一时相和后一时相的多光谱影像的任意两个原始波段进行波段光谱信息扩展，获取新扩展的非线性人工光谱波段；

2)基于多尺度形态学重构的波段空间信息扩展：分别对原始波段和非线性人工光谱波段通过空间重构运算获取空间扩展波段特征集后，进行差分运算得到差分波段集；

3)非监督多类变化检测。

2.根据权利要求1所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤1)中，非线性人工光谱波段的计算式为：

其中，X_1,α、X_1,β分别为B维的前一时相多光谱影像X₁的任意两个原始波段，M₁、R₁分别为前一时相多光谱影像X₁扩展后的非线性人工光谱波段，X_2,α、X_2,β分别为B维的后一时相多光谱影像X₂的任意两个原始波段，M₂、R₂分别为后一时相多光谱影像X₂扩展后的非线性人工光谱波段，α、β为波段号。

3.根据权利要求2所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤1)中，在进行除法运算时，通过添加一常数t来避免出现分母为极小非零值的情况。

4.根据权利要求2所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤2)具体包括以下步骤：

21)给定圆形结构元素尺度i的范围，并对影像X₁进行形态学开重建和闭重建得到形态学剖面结构特征，即：

其中，MPⁱ(X₁)为影像X₁在尺度i下进行开闭运算得到的特征集合，为影像X₁在尺度i下进行开运算得到的特征，为影像X₁在尺度i下进行闭运算得到的特征，为影像X₁第b个波段在尺度i下的形态学开重建结果，为影像X₁第b个波段在尺度i下的形态学闭重建结果；

22)分别对影像X₁、影像X₂以及非线性人工光谱波段M₁、R₁、M₂和R₂采用与步骤21)相同的空间重构运算生成空间扩展波段特征集F_X₁、F_X₂、F_M₁、F_M₂、F_R₁和F_R₂；

23)对空间扩展波段特征集进行差分运算得到最终经过波段光谱-空间信息扩展后的差分波段集E_D。

5.根据权利要求4所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤22)中，影像X₁和影像X₂对应的空间扩展波段特征集F_X₁和F_X₂的表达式为：

F_X₁＝[MP²(X₁),...,MP⁶(X₁)]

F_X₂＝[MP²(X₂),...,MP⁶(X₂)]。

6.根据权利要求4所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤23)中，差分波段集E_D的计算式为：

E_D＝[X_D,F_D]

X_D＝X₂-X₁

F_D＝[F_X₂,F_M₂,F_R₂]-[F_X₁,F_M₁,F_R₁]

其中，X_D为原始光谱变化表达，F_D为扩展后的增强变化表达。

7.根据权利要求6所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，所述的步骤3)中，采用压缩变化矢量分析、序列光谱变化矢量分析或迭代加权多变量变化检测进行多类变化检测。

8.根据权利要求7所述的一种基于信息扩展的多光谱遥感影像非监督变化检测方法，其特征在于，当采用压缩变化矢量分析和序列光谱变化矢量分析时，通过压缩变化矢量表达描述光谱变化，定义变量变化强度ρ和变化方向θ，根据变量变化强度ρ和变化方向θ构建二维极坐标变化表达域，ρ值用以表示变化的可能性，即ρ值越大，越有可能发生变化，θ值表示可能的不同类别的变化，所述的变化强度ρ和变化方向θ的表达式为：

其中，E_D,n为E_D中第n个元素，N为E_D中的元素总数，即扩展波段总维度，r_n为单位常数向量中的第n个元素。