CN113255698B - 用于高光谱影像空间特征提取的超像素级自适应ssa方法 - Google Patents

Abstract

本发明公开了用于高光谱影像空间特征提取的超像素级自适应SSA方法，包括以下步骤：S1、对原始高光谱影像进行PCA降维，得到第一主成分图；S2、在第一主成分图上进行熵率分割，得到超像素图；S3、确定超像素图中每个波段上所有超像素区域的尺寸大小，对较大的超像素区域使用带有自适应窗口的2DSSA来进行特征提取；对较小的超像素区域使用SSA来进行特征提取；S4、将2DSSA和SSA处理后的区域或向量，重新投影到原来的超像素区域；最终得到所有超像素区域的特征图像；S7、将特征图像进行堆叠操作，得到与原始高光谱影像相同尺寸大小的特征影像；S8、使用支持向量机分类器对特征影像进行精度评价并获得影像分类结果。

Description

用于高光谱影像空间特征提取的超像素级自适应SSA方法

技术领域

本发明涉及图像处理领域，尤其涉及用于高光谱影像空间特征提取的超像素级自适应SSA方法。

背景技术

随着传感器技术的快速发展，高光谱影像在20世纪80年代之后，开始在环境监测、国土资源调查评估、城市规划等领域发挥重要的作用。在此过程中，高光谱的分类精度直接决定了其应用的有效程度。然而，高光谱影像分类还存在诸多挑战。高光谱数据的光谱维数，即波段数，通常比其标记类别中的样本数大得多。这导致了模型的复杂度与训练样本的数量严重不匹配，即Hughes现象。此外，光谱和空间域的数据冗余以及环境因素和传感器引起的噪声需要考虑到。因此，有效的特征提取和降维是高光谱遥感土地覆盖分类的关键。

在过去的几十年里，研究人员提出一系列的特征提取方法，如主成分分析(PCA)、线性判别分析(LDA)、最大噪声分数(MNF)等，已经被设计用于特征变换。经过特征变换后，高光谱数据的维度得到显著降低，并且大部分的信息得到保留，去除了噪声等影响成分。但是这些方法作用于光谱域时，只考虑光谱信息，缺乏空间分析，分类中存在严重的点状噪声现象。

空间特征提取方法可以有效利用像素的空间关系，有效的提高分类精度。基于形态轮廓(MPs)的空间特征提取方法是目前应用最广泛的空间特征提取方法之一，该方法利用一系列结构元素(SEs)来提取图像的空间几何细节。通过从高光谱影像的前几个主成分图像中提取MPs，得到了扩展的MPs(EMPs)，能够实现空间光谱特征提取。考虑到MPs的关联结构元素(SEs)的影响，研究人员进一步提出了形态属性轮廓(APs)，该方法将一组形态属性滤波器应用于减少的HSI成分，获得其特征的多级特征，用于建模不同类型的结构信息。此外，近年来，一些基于图学习的光谱空间特征提取方法也被提出，这种方法用于揭示高光谱影像复杂的高阶关系，降低数据维数，用简单的分类器获得较好的分类结果。

近年来，一种时间序列分析技术——奇异谱分析(SSA)被引入到HSI中进行有效的特征提取。对于给定的一个基于光谱像素向量，SSA可以将其分解成几个子向量，每个子向量的大小与原始向量相同。通过选择子分量重建谱轮廓，可以增强有用信息，有效抑制噪声或不具有代表性的信号，提高分类精度。与主成分分析相比，SSA保留了更多的光谱信息，在后续的数据分类任务中具有更好的可分性。针对二维数据阵列，特别是图像的情况，扩展的二维SSA(2DSSA)被提出，并被认为其具有与SSA相似的能力。Zabalza等进一步将2DSSA作用于高光谱影像的每个波段图像进行空间特征提取。与传统空间技术对固定结构元素或小邻域(窗口)敏感不同，2DSSA可以通过嵌入嵌套的Hankel轨迹矩阵来利用HSI的全局空间相关性，优于许多传统方法。

然而，SSA和2DSSA在高光谱空间特征提取方面仍有一定的局限性。传统的SSA只提取光谱趋势特征，而忽略了丰富的空间上下文信息，导致分类图中出现误分类和噪声，将SSA应用于空间域似乎更能提高分类性能。至于2DSSA，它有几个主要的缺点。首先，考虑到遥感图像地物形状不规则，大小不一致，采用固定矩形窗口提取空间特征不适合2DSSA，这可能导致分类结果较差，特别是小区域和地物边界。此外，固定大小的嵌入窗口或提取尺度不能充分利用丰富的空间特征。如何选择合适的尺度进行特征提取是一个值得探索的问题。最后，传统的2D-SSA在内存和处理方面存在较大的计算代价，尤其是对于大型嵌入窗口和大型图像。为了更加高效的利用SSA技术解决空间特征提取的问题，需要对该技术进行改进与提高。

发明内容

本发明目的是针对高光谱影像空间特征提取中的地物尺度复杂、噪声问题，以及传统SSA方法的不足，本发明的目的是提供一种用于高光谱影像空间特征提取的超像素级自适应SSA方法，该方法通过合理评价超像素尺寸大小，在超像素区域上进行SSA或2DSSA特征提取，可以解决传统SSA分类精度较低、大窗口下2DSSA计算复杂度高、无法提取联合谱空间特征等缺点，深入利用不同的均匀区域的空间特性。

为了实现上述目的，本发明的技术方案是：

用于高光谱影像空间特征提取的超像素级自适应SSA方法，包括以下步骤：

S1、对原始高光谱影像进行PCA降维，得到第一主成分图；

S2、在第一主成分图上进行熵率分割，得到超像素图；

S3、确定超像素图中每个波段上所有超像素区域的尺寸大小，并将其划分为两个阈值：T1，T2，满足T2>T1；对于大于T1的超像素区域，将超像素区域填补为规则的区域，使用带有自适应窗口的2DSSA来进行特征提取；对于小于T1的超像素区域，将超像素区域中的所有像素拉伸为一维向量，并使用SSA来进行特征提取；

S4、将2DSSA和SSA处理后的区域或向量，重新投影到原来的超像素区域；最终得到所有超像素区域的特征图像；

S7、将特征图像进行堆叠操作，得到与原始高光谱影像相同尺寸大小的特征影像；

S8、在原始高光谱影像与特征影像中，随机挑选一定数量的样本作为训练集，剩余样本作为测试集，使用支持向量机分类器对特征影像进行精度评价并获得影像分类结果。

进一步的，所述步骤S2中进行熵率分割包括以下步骤：

S21、将第一主成分图映射为图G，G＝(V,E)，其中V是图G的顶点，代表像素向量，E是图G的边缘集合，表示相邻像素之间的相似性；

S22、选择子边缘集A使得图G恰好包含S个离散的子图，每个子图代表一个超像素；

S23、使用构造的熵率项函数H(·)和平衡项函数B(·)得到目标函数，通过最优化目标函数来获得紧凑且均匀的S个超像素，目标函数如下：

其中，λ≥0是平衡项函数的权重调节参数。

进一步的，所述步骤S2中，对于分割后的超像素区域，其尺寸大小由参数f来衡量，记为：

f＝min(col,row)；

其中，row和col为超像素区域的高度和宽度。

进一步的，所述步骤S3中需要对超像素区域特征提取时的处理窗口大小进行判别，其判别公式为：

其中，L_1D是SSA的处理窗口；L_x和L_y是2DSSA的处理窗口，L_SpaSSA表示SpaSSA的处理窗口大小。

进一步的，所述步骤S3中使用带有自适应窗口的2DSSA进行特征提取包括以下步骤：

S31、定义处理窗口的尺寸为L_x×L_y，将处理窗口从超像素区域的左上角滑动到右下角，不同位置处处理窗口内所有像素按列拉伸为向量

将其作为轨迹矩阵P_2D的列向量，可得轨迹矩阵P_2D：

S32、计算协方差矩阵P_2DP_2D ^T的特征值(λ₁≥λ₂≥...≥λ_LxLy)和相应的特征向量(U₁,U₂,...,U_LxLy)，轨迹矩阵P_2D可以写成如下形式：

其中U_i和V_i分别为轨迹矩阵P_2D的左奇异向量和右奇异向量；

S33、通过两步对角平均处理操作，可将矩阵P_2D ¹再次转换为row×col大小的规则区域。

进一步的，所述步骤S3中使用SSA进行特征提取包括以下步骤：

S34、将超像素区域内的所有像素拉伸为一个一维向量，设其长度为N，使用一维处理窗口L_1D将其嵌入为一个轨迹矩阵P_1D：

其中p_i代表该超像素内的所有像元值，该轨迹矩阵的元素在反对角线上是相同的，也称为Hankel矩阵；

S35、对轨迹矩阵P_1D进行奇异值分解，得到多个子矩阵相加的形式：

其中κ为奇异值，u_i和v_i分别为轨迹矩阵P_1D的左奇异值向量、右奇异值向量；

S36、选取最大奇异值对应的子矩阵来代替原始矩阵P_1D，通过对角平均过程，将选取的子矩阵重新转化为长度为N的向量。

与现有技术相比，本发明具有的优点和积极效果是：

1、SpaSSA可以自适应地提取高光谱影像的空间特征。对于每个波段的大部分超像素区域采用2D-SSA来获取的局部空间特征，而对于超像素过小的情况则采用SSA来获取空间特征，通过2D-SSA和SSA自适应地应用于局部均匀区域，以增强不同目标之间的差异，保持边缘特征。

2、本发明可以根据超像素的大小自适应地确定相关联的2D-SSA的嵌入窗口大小(或提取尺度)，降低了2DSSA的计算复杂度，提高了计算效率。

3、本发明中的SpaSSA以超像素的方式操作，并不是在完整图像上执行，SVD的轨迹矩阵的大小和计算复杂度显著降低，因此大大降低了计算成本，提高了工作效率。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。

图1为本发明SpaSSA方法用于高光谱空间特征提取和分类的总体结构图；

图2为超像素2DSSA的空间特征提取处理示意图；

图3为超像素SSA的空间特征提取处理示意图；

图4为SpaSSA与传统SSA和2DSSA的分类对比图；其中，a为传统SSA方法，b为传统2DSSA方法，c为本发明提出的SpaSSA方法。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

本发明提供了一种用于高光谱影像空间特征提取的超像素级自适应SSA方法，可用于由AVIRIS(AirborneVisibleInfraredImagingSpectrometer)传感器获取的IndianPine地区高光谱影像的处理。AVIRIS传感器的空间分辨率和波谱分辨率分别为10米和10纳米，波谱范围400-2400nm，共包括224个波段，其中4个波段的值均为0，所以被去除。此外，IndianPine影像上有20个波段，包括波段[104-108]、[150-163]和220的值容易受吸水波段的影响，这些波段在本案例中也被去除。因此，本案例采用的IndianPine影像共有200个波段，图像尺寸为145×145像元大小，共包括16种地物类型(主要是不同农作物)，其余像素被认为是背景。

如图1所示，本发明提供的用于高光谱影像空间特征提取的SpaSSA方法包括以下详细步骤：

一、对于原始高光谱影像数据，使用主成分分析(PCA)对影像进行降维，得到的第一主成分图作为分割的基础图像，以降低后续的计算量。

二、确定合适的超像素数量，基于第一主成分图，采用熵率分割(ERS)方法获取超像素图。具体来讲，ERS将基准图像映射为一个图G＝(V,E)，其中V是图的顶点，代表像素向量，E是图的边缘集合，表示相邻像素之间的相似性。选择子边缘集A使得图G恰好包含S个离散的子图，每个子图代表一个超像素。然后，使用构造的熵率项函数H(·)和平衡项函数B(·)得到目标函数，通过最优化目标函数来获得紧凑且均匀的S个超像素，目标函数如下：

式中，λ≥0是平衡项函数的权重调节参数。通过基于稀疏表示的距离度量为每个超像素指定一个公共标签。

三、对于高光谱每个波段，由分割得到S个超像素。每个超像素的尺寸大小由参数f来衡量，记为：

f＝min(col,row)

其中row和col为超像素的高度和宽度。该参数可用于判断确定SpaSSA方法中处理不同大小超像素的方式。

四、定义两个阈值：T₁，T₂，满足T₂>T₁。较小的阈值用以确定选择SSA还是2DSSA作用于超像素区域，即大于T₁的超像素区域使用2DSSA来处理，反之使用SSA处理；较大的阈值T₂用来决定2DSSA的嵌入处理窗口的大小，以保证2DSSA的处理效率。具体的判别方式如下：

其中L_1D是SSA的处理窗口，考虑到影响较小，将其固定为10；L_x和L_y是2DSSA的处理窗口，L_SpaSSA表示SpaSSA的处理窗口大小。经过实验与验证，我们将阈值T₁和T₂分别固定为3和11，以保证2DSSA处理尽可能多的超像素和尽可能高的处理效率。

五、对于较大尺寸的超像素，使用自适应窗口的2DSSA来进行处理，如图2所示。对于任一较大的超像素，考虑到2DSSA通常作用于规则的区域，因此首先为超像素构建一个规则的作用区，即根据超像素的行数和列数，结合邻域像素来补充超像素区域，得到的超像素尺寸为row×col。根据上述原则，定义处理窗口L_x×L_y，该窗口从超像素区域的左上角滑动到右下角，不同位置处窗口内所有像素按列拉伸为向量

将其作为轨迹矩阵P_2D的列向量，可得：

得到的轨迹矩阵P_2D有一个名为HbH的结构，即Hankel-block-Hankel矩阵。计算协方差矩阵P_2DP_2D ^T的特征值(λ₁≥λ₂≥...≥λ_LxLy)和相应的特征向量(U₁,U₂,...,U_LxLy)，轨迹矩阵可以写成如下形式：

其中U_i和V_i分别为轨迹矩阵的左奇异向量和右奇异向量。选择P_2D ¹成分作为P_2D的近似值，主要因为它包含了最主要的空间信息。通过两步对角平均处理操作，可将矩阵P_2D ¹再次转换为row×col大小的规则区域。将原始超像素区域的像素值提取并重新投影到超像素位置，可以得到特征提取之后的超像素。2DSSA的处理窗口时自适应的，在其不大于阈值T₂的情况下，它的大小取与该区域大小相适应的最大窗口，而在大于T₂的情况下，则使用T₂作为最大处理窗口，以防处理窗口过大而造成处理效率降低。

六、对于尺寸小于阈值T₁的超像素，使用SSA方法进行特征提取，如图3所示。超像素区域的所有像素，首先被拉伸为一个一维向量，假设其长度为N，使用一维处理窗口L_1D将其嵌入为一个轨迹矩阵P_1D：

其中p_i代表该超像素内的所有像元值，该轨迹矩阵的元素在反对角线上是相同的，也称为Hankel矩阵。然后对矩阵P_1D进行奇异值分解，得到多个子矩阵相加的形式：

其中κ为奇异值，u_i和v_i分别为左右奇异值向量。选取最大奇异值对应的子矩阵P_1D ¹来代替原始矩阵P_1D。最后，通过对角平均过程，将选取的子矩阵重新转化为长度为N的向量。该向量被重新投影到原始的超像素区域内，得到SSA特征提取后的超像素。

七、对高光谱影像的每个波段内，所有不同尺寸的超像素区域分别进行步骤五或步骤六的处理，可以得到每个经过特征提取的特征图像，将所有波段处理后的特征图进行堆叠，得到与原始高光谱影像相同尺寸大小的特征影像。

八、对得到的特征影像，随机选取一定数量的像素样本作为训练，剩余样本作为测试，使用SVM分类器来对特征进行精度评估，并获得影像分类结果。

为了验证空间特征提取的分类效果，分别从影像中随机选取10％的样本作为训练，剩余的样本用于测试，采样方式为分层抽样，训练集和测试集中的样本互不重叠。精度评价指标包括三个：总体精度(OA)，平均精度(AA)和kappa系数。此外，为避免系统误差和减少随机偏差，实验统计了在10次独立训练和测试下的平均分类精度。

表1给出了本发明中的方法和原始的SSA、2DSSA方法用于IndianPine影像特征提取，并经过SVM分类的总体分类精度。

表1特征提取后分类的精度对比

评价指标	SSA+SVM	2DSSA+SVM	SpaSSA+SVM
				OA	84.02	95.76	97.97
AA	80.37	94.63	97.79
				kappa	81.78	95.16	97.69

如表1所示，相比于传统的SSA和2DSSA算法，本发明中的方法能够取得更优的分类精度。分类结果对应的分类图如图4所示，其可以看出，SpaSSA的分类结果图很好的消除了SSA中存在的分类“椒盐噪声”现象，并且降低了2DSSA中的条带误分类，在地物边缘处的误分类较少。

本发明相较于传统SSA具有以下特点：

1、SpaSSA可以自适应地提取高光谱影像的空间特征。对于每个波段的大部分超像素区域采用2D-SSA来获取的局部空间特征，而对于超像素过小的情况则采用SSA来获取空间特征，通过2D-SSA和SSA自适应地应用于局部均匀区域，以增强不同目标之间的差异，保持边缘特征。

2、本发明可以根据超像素的大小自适应地确定相关联的2D-SSA的嵌入窗口大小(或提取尺度)，降低了2DSSA的计算复杂度，提高了计算效率。

3、本发明中的SpaSSA以超像素的方式操作，并不是在完整图像上执行，SVD的轨迹矩阵的大小和计算复杂度显著降低，因此大大降低了计算成本，提高了工作效率。

Claims (4)

1.用于高光谱影像空间特征提取的超像素级自适应SSA方法，其特征在于：包括以下步骤：

S1、对原始高光谱影像进行PCA降维，得到第一主成分图；

S2、在第一主成分图上进行熵率分割，得到超像素图；

S3、确定超像素图中每个波段上所有超像素区域的尺寸大小，并将其划分为两个阈值：T1，T2，满足T2>T1；对于大于T1的超像素区域，将超像素区域填补为规则的区域，使用带有自适应窗口的2DSSA来进行特征提取；对于小于T1的超像素区域，将超像素区域中的所有像素拉伸为一维向量，并使用SSA来进行特征提取；

S4、将2DSSA和SSA处理后的区域或向量，重新投影到原来的超像素区域；最终得到所有超像素区域的特征图像；

S7、将特征图像进行堆叠操作，得到与原始高光谱影像相同尺寸大小的特征影像；

S8、在原始高光谱影像与特征影像中，随机挑选一定数量的样本作为训练集，剩余样本作为测试集，使用支持向量机分类器对特征影像进行精度评价并获得影像分类结果；

所述步骤S3中需要对超像素区域特征提取时的处理窗口大小进行判别，其判别公式为：

其中，L_1D是SSA的处理窗口；L_x和L_y是2DSSA的处理窗口，L_SpaSSA表示SpaSSA的处理窗口大小；

所述步骤S3中使用带有自适应窗口的2DSSA进行特征提取包括以下步骤：

S31、定义处理窗口的尺寸为L_x×L_y，将处理窗口从超像素区域的左上角滑动到右下角，不同位置处处理窗口内所有像素按列拉伸为向量

将其作为轨迹矩阵P_2D的列向量，可得轨迹矩阵P_2D：

S32、计算协方差矩阵P_2DP_2D ^T的特征值(λ₁≥λ₂≥…≥λ_LxLy)和相应的特征向量(U₁,U₂,…,U_LxLy)，轨迹矩阵P_2D可以写成如下形式：

P_2D＝P_2D ¹+P_2D ²+…+P_2D ^LxLy

其中U_i和V_i分别为轨迹矩阵P_2D的左奇异向量和右奇异向量；

S33、通过两步对角平均处理操作，可将矩阵P_2D ¹再次转换为row×col大小的规则区域。

2.如权利要求1所述的用于高光谱影像空间特征提取的超像素级自适应SSA方法，其特征在于：所述步骤S2中进行熵率分割包括以下步骤：

S21、将第一主成分图映射为图G，G＝(V,E)，其中V是图G的顶点，代表像素向量，E是图G的边缘集合，表示相邻像素之间的相似性；

S22、选择子边缘集A使得图G恰好包含S个离散的子图，每个子图代表一个超像素；

S23、使用构造的熵率项函数H(·)和平衡项函数B(·)得到目标函数，通过最优化目标函数来获得紧凑且均匀的S个超像素，目标函数如下：

其中，λ≥0是平衡项函数的权重调节参数。

3.如权利要求2所述的用于高光谱影像空间特征提取的超像素级自适应SSA方法，其特征在于：所述步骤S2中，对于分割后的超像素区域，其尺寸大小由参数f来衡量，记为：

f＝min(col,row)；

其中，row和col为超像素区域的高度和宽度。

4.如权利要求1所述的用于高光谱影像空间特征提取的超像素级自适应SSA方法，其特征在于：所述步骤S3中使用SSA进行特征提取包括以下步骤：

S34、将超像素区域内的所有像素拉伸为一个一维向量，设其长度为N，使用一维处理窗口L_1D将其嵌入为一个轨迹矩阵P_1D：

其中p_i代表该超像素内的所有像元值，该轨迹矩阵的元素在反对角线上是相同的，也称为Hankel矩阵；

S35、对轨迹矩阵P_1D进行奇异值分解，得到多个子矩阵相加的形式：

其中κ为奇异值，u_i和v_i分别为轨迹矩阵P_1D的左奇异值向量、右奇异值向量；

S36、选取最大奇异值对应的子矩阵来代替原始矩阵P_1D，通过对角平均过程，将选取的子矩阵重新转化为长度为N的向量。

Images