一种波段聚类和多尺度结构特征融合的高光谱图像分类模型

王彩玲,张 静,王洪伟,宋晓楠,纪 童

1. 西安石油大学计算机学院,陕西 西安 710065 2. 西北工业大学光电与智能研究院,陕西 西安 710072 3. 甘肃农业大学草业学院,甘肃 兰州 730070

引 言

高光谱遥感技术利用光谱探测与成像技术对感兴趣物体成像,形成了高光谱图像立方体。 它包含丰富的光谱和空间信息,是一个三维数据[1],因此被广泛应用于生态保护[2]、 灾害监测[3]、 医疗卫生[4]等领域。 在高光谱(hyperspectral image,HSI)数据处理任务中,分类仍然是一项具有挑战性的任务,旨在为每个像元分配唯一标签,对许多实际应用非常关键。 一开始,我们通过人眼进行分类,取得了不太理想的效果。 到后来,利用机器学习方法虽然提高了分类精度,但仍存在许多问题: (1)样本稀少引发的“Hughes”现象; (2)高冗余信息以及噪声问题; (3)“同物异谱”和“同谱异物”现象。 这些问题导致分类精度较低,不利于后续的图像分析。

为了克服这些问题,多种特征提取方法用于保持光谱分辨能力。 一些基于特征变换的方法如主成分分析[5],独立主成分分析,它们在降低光谱维度的同时保留了数据的原有信息,但却改变了波段本身的物理意义。 另外一种基于波段选择的方法由于保留了高光谱数据的固有属性而备受关注。 它在不损失物理意义和信息的前提下,选择了低相关性,信息量大的子集。 波段相似性度量方式包括: 信息散度(information divergence,ID)、 欧氏距离(Euclidean distance,ED)[6]等。

研究表明,仅提取数据内部的光谱信息对高光谱分类并不充分。 近年来,利用HSI光谱和图像信息的光谱-空间分类方法已被广泛研究。 例如,为了得到图像内部的纹理特征,1973年,Haralick等引入了灰度共生矩阵(gray level co-ccurennce matrix,GLCM)的概念,通过计算像元之间的联合概率密度得到灰度关系。 2019年,Tan等提出了三维GLCM光谱-空间分类方法,该模型利用GLCM特征图像提取异质性组织内部的纹理特征,取得了不错的效果[7]。 Ojala等提出了局部二值模式(local binary patterns,LBP),由于对灰度变化敏感度低且计算速度快,也被用于图像的纹理分析。 再后来,形态学(morphology profile,MP)、 扩展形态学[8]、 多属性轮廓特征[9]等方法被提出,使得分类精度进一步提高。 2015年,Dalla等将独立主成分分析(independent component analysis,ICA)和扩展形态学相结合对高光谱遥感影像进行分类,实验分析证明了此方法在图像空间信息建模的有效性[10]。 Hou等提出了一种基于三维形态剖面(3D-MP)的光谱空间分类方法,该模型充分利用了三维光谱空间特征,在样本有限的情况下提高了分类精度[11]。 近些年,众多学者开展基于滤波的分类技术研究。 典型的有小波变换、 Gabor滤波[12]等。 Pablo等提出一种基于小波变换和数学形态学的光谱空间技术,该方法不仅提取到了数据的深层空谱特征,还证明了对图像噪声问题的鲁棒性[13]。 Jia等提出了基于光谱-空间Gabor表面特征(GSF)融合的空谱分类模型,该方法保留了数据内部丰富的光谱-图像信息,进一步提升了分类性能[14]。 慕彩红等提出了一种基于空间坐标的光谱-空间特征融合的HSI分类算法,将HSI在空间维度上划分为若干小图像,根据空间坐标特征和光谱特征分别得到各类别样本的概率,SVM进一步对概率特征进行分类,得到最终的分类结果[15]。 目前,基于空间信息和光谱信息相结合的方法已经成为高光谱图像分类的主流方法,研究证明,空谱特征联合是提高分类性能的有效方法,能降低椒盐噪声的影响,取得较好的区域一致性结果。 然而,现有的光谱空间特征提取方法仍面临以下几个挑战:

(1)这些方法引入空间信息的方式更复杂,在一定程度上增加了计算成本。

(2)地物含有复杂的纹理、 结构特征,现有的模型大多只考虑了单一尺度下的光谱特征与空间特征的融合,而没有考虑多尺度下光谱特征与空间特征的丰富相关性。

(3)现有的特征提取方法大多没有考虑到由于弱光或其他因素可能产生的图像噪声。 在这种不完善的情况下,现有的光谱空间特征提取方法的性能往往会严重下降。

为了解决以上这些问题,从数据中提取足够的空谱联合特征,提出了一种波段聚类和多尺度结构特征融合的高光谱图像分类模型,其新颖性和贡献主要包括以下两个方面:

(1)考虑到相邻波段的相关性高于非相邻波段的相关性,使用基于光谱散度信息的聚类方法将高光谱数据分为有限的、 波段数量不等的子立方体,在每个子立方体中求取平均值叠加后作为光谱特征输出,从而实现高质量带图像的选择。

(2)图像的整体结构信息才是人类视觉感知系统的主要感知数据,而纹理等“细节特征”往往是不需要的。 在不同的参数设置下去除图像的纹理细节部分,将主结构从复杂的纹理图像中提取出来,得到多尺度结构特征信息,可以有效的去除噪声,增加不同类别之间的差异。

实验结果表明,该模型不仅充分提取了高光谱数据内部的空谱特征,还对小样本问题非常有效。

1 研究方法

1.1 归一化

为了消除指标间量纲以及奇异样本数据导致的不良影响,需要对数据进行归一化处理,如式(1)所示

(1)

Indian Pines数据内部光谱波段差别较大,因此我们对其进行归一化。 设X为空间中某一像素点,表示为X={x1,x2,…,x200},由于它有200个光谱波段,对其归一化处理的步骤如式(2)所示

(2)

式(2)中,xmin为X中的最小光谱值,xmax为X中的最大光谱值。Xnew(i)为X归一化后新的像素值。

University of Pavia数据内部空间物质相似而像素值差别很大,因此对空间维进行归一化,该数据空间维大小为610×340,将其拉伸为一维向量,可以表示为Y={y1,y2,…,y207 400},对其归一化处理的步骤如式(3)所示

(3)

式(3)中,ymin为Y中的最小光谱值,ymax为Y中的最大光谱值。Ynew(i)为Y归一化后新的像素值。

1.2 光谱特征提取框架-自适应子空间划分策略

在高光谱数据中,波段之间的相似性与其距离相关。 设原始高光谱图像立方体为L,I为光谱波段总数。 为了快速对高光谱数据进行分割,首先,我们采取粗划分策略,将数据L均分为n个波段数为Z的子空间,如式(4)和式(5)。

(4)

L→{l1,l2,…,ln}

(5)

其中,li(i≤n)为粗划分后的子空间。 为了对每个粗子空间li进行更加准确的表示,我们采取细划分策略来获取新的子空间。 由于相距远的波段相关性极小,因此我们仅考虑相邻子空间li,li+1(i+1≤n)内的波段。 具体过程如下: (1)将每个空间波段拉伸为一维向量; (2)用信息散度构造相邻波段的相似性矩阵; (3)对相邻子空间按照最大化类间距离与类内距离比值的方法进行自适应划分。 细划分策略后目标函数的形式如式(6)所示

argmin(dintra-dinter)

(6)

式(6)中,dinter是相邻两个子空间的类间距离,定义如式(7)所示

dinter=max|dID(m,n)|

(7)

式(7)中,dID(m,n)表示第i(i≤n)个子空间中的波段m和第i+1个子空间中的波段之间的信息散度值。

dintra是子空间类内距离,定义如式(8)—式(10)所示

dintra=Mi+Mi+1

(8)

(9)

(10)

式中,dID(p,q)为同一个子空间内第p和第q个波段之间的信息散度值。Mi和Mi+1表示相邻子空间散度和的平均值,t为第i个子立方体中包含的频带总数。 通过粗细划分策略后,子空间表示如式(11)。

B={b1,b2,b3,…,bn}

(11)

由于自适应划分后的子空间内包含不同数目的波段,我们对每个子空间内的光谱波段求取平均值,具体如式(12)所示

Avgb={avgb1,avgb2,avgb3,…,avgbn}

(12)

式(12)中,Avgb表示每个子空间光谱平均值avgbi叠加后的结果。 图1为自适应子空间划分的光谱特征提取框架。

图1 自适应子空间划分策略-光谱特征提取框架Fig.1 Adaptive subspace partitioning strategy-spectral feature extraction framework

1.3 空间特征提取框架-多尺度相对全变分

为了从高光谱图像数据中获取有用的空间特征,结构特征提取扮演着越来越重要的角色。 相对全变分(relative total variation,RTV)方法作为其中的一种,可以有效地减小类内差异,增大类间差异,并且对于椒盐噪声和纹理特征的去除也至关重要。 然而,由于参数σ和λ直接影响着RTV的效果,因此,单一尺度下的结构特征无法满足高光谱数据结构丰富,物质多样的特性。

针对RTV的缺点,通过设置多参数,对光谱特征数据进行不同尺度下的相对全变分提取,形成多尺度相对全变分(multi-scale RTV,MRTV)空间特征提取框架,如式(13)和式(14)所示。

fi=RTV(I,λi,σi) (i≤c)

(13)

F={f1,f2,f3,…,fc}

(14)

式中,I为原始高光谱数据,假设I的波段数为m,对I进行c次不同参数设置下的RTV计算,可以得到c个波段数为m的特征子集{f1,f2,f3,…,fc}。λi和σi分别表示控制平滑度和空间尺度的参数。fi为第i个特征子集,波段数为m。F表示c个特征子集按光谱维堆叠后的立方体数据,波段数为m×c。

虽然MRTV有效地提取了高光谱数据内部的结构特征,但由于堆叠操作导致了高维性,增加了计算成本。 由于核主成分分析(kernel PCA,KPCA)对于提高像素的光谱可分离性和降低光谱维度非常有效,因此在数据堆叠之后进行KPCA处理,具体过程如式(15),N为P中保留的主成分数量,P为最终特征输出。 特征提取过程如图2所示。

图2 空间特征提取框架-多尺度相对全变分Fig.2 Spatial feature extraction framework-multi-scale relative total variation

P=KPCA(F,N)

(15)

1.4 空谱特征提取框架-波段聚类和多尺度结构特征融合

提出了一种波段聚类和多尺度结构特征融合的高光谱图像分类模型(ASPS-MRTV),为了充分提取高光谱数据内部的光谱特征和空间特征,分别构建了自适应子空间划分策略(adaptive subspace partitioning strategy,ASPS)光谱特征提取框架和多尺度相对全变分空间特征提取框架。 过程如下:

(1) 对高光谱数据进行归一化处理;

(2) 根据式(4),将原始高光谱数据L均分为n个波段数为Z的子空间;

(3) 将每个空间波段拉伸为一维向量后用ID构造相邻波段的相似性矩阵;

(4) 根据式(7)—式(10)对细划分策略后的目标函数式(6)进行求解;

(5) 通过粗细划分策略后得到自适应子空间bi(i

(6) 对每个子空间bi求取光谱波段平均值得到avgbi;

(7) 将每个自适应子空间内的平均值avgbi进行叠加,得到Avgb;

(8) 利用MRTV空间特征提取框架对Avgb进行多尺度特征提取,将得到的多个子空间按光谱维堆叠形成立方体数据F;

(9) 对立方体数据F进行KPCA处理后得到空谱特征P。

流程图如图3所示。

图3 空谱特征提取框架-波段聚类和多尺度结构特征融合Fig.3 Frame-band clustering and multi-scale structural feature fusion for spatial spectrum feature extraction

2 实验部分

2.1 数据集

Indian Pines数据(简称为IP)由红外成像光谱仪在印第安拍摄所获取的。 数据大小为145×145×220,具有20 m的空间分辨率,光谱范围为0.2～0.4 μm。 实验中去除20个水吸收区域的光谱带,保留了200个波段。 图4为IP数据的伪彩图,地面真值图及标签。

图4 Indian Pines数据(a): 伪彩图; (b): 地面真值图; (c): 标签Fig.4 Indian Pines data(a): False color map; (b): Ground truth map; (c): Label

University of Pavia数据(简称为UP)是由ROSIS-3在帕威亚大学拍摄所获取的。 数据大小为610×340×115,实验中去除12个水吸收区域的光谱带,保留了103个波段,空间分辨率为1.3 m,光谱覆盖范围为0.43～0.86 μm。 图5为UP数据的伪彩图,地面真值图及标签。

图5 University of Pavia数据(a): 伪彩图; (b): 地面真值图; (c): 标签Fig.5 University of Pavia data(a): False color map; (b): Ground truth map; (c): Label

2.2 实验设置

实验环境为Windows10和Matlab2019a,处理器为Inter Core I7-8750,GPU为NVIDIA GeForce GTX1050Ti,内存为16 GB。 实验过程中,训练样本为5%的IP数据以及1%的UP数据,测试样本为剩余的样本。 为了验证波段聚类和多尺度结构特征融合框架(ASPS-MRTV)的有效性,我们使用惩罚参数C和核参数σ都为24.5的SVM进行分类,并采用OA、 AA和 Kappa 系数来评价分类效果的好坏。 实验中各参数及其含义如表1所示。

表1 各参数及其含义Table 1 Parameters and their meanings

2.3 参数优化

结构特征质量的好坏是由RTV的参数λ和σ决定的。 通常情况下,将λ固定在[0,0.05]之间,σ固定在[1,3]之间。 为了避免高光谱数据结构特征的单一性,考虑到时间复杂度和空间复杂度的影响,设置了包含大小尺度的5个λ,λ1=0.003,λ2=0.015,λ3=0.01,λ4=0.02,λ5=0.05,实验过程中,采用控制变量的方法对σ调优,设置K=N=30,σ1-σ5为1,σ在IP上的寻优的结果如表2所示,最优精度加粗表示。

表2 σ在IP上的调优过程Table 2 Process of optimizing σ for IP

参数K的优化也至关重要。 一开始,将K设置为30,由于IP数据包含200个光谱波段,所以将K固定在[10,100]之间,实验结果如表3所示。 同样地,也将N固定在[10,100]之间,具体结果见表4,最优结果加粗表示。

表3 参数K对IP的分类效果Table 3 IP classification effect of parameter K

表4 参数N对IP的分类效果Table 4 IP classification effect of parameter N

从表2的实验结果可知,当σ1=1,σ2=1,σ3=1,σ4=2,σ5=3时IP数据的分类效果最好,此时OA=96.76%,Kappa=96.28%,AA=85.66%。 由表3可以得知,当K=30时,OA,Kappa,AA都取得了最高的分类精度,分别达到了96.76%,96.28%,85.66%。 另外我们也发现,当参数K不断增大时,OA呈现出先增大后减小的趋势。 这是由于当K较小时包含的有用信息也少,但随着K的不断增大,有用信息变多,分类精度也就提高了,然而当K过大时信息的冗余又导致分类精度下降,所以才会出现先增大后减小的现象。 由表4可知,当N=20时,OA,Kappa都达到了最高,分别为97.06%,96.64%。 随着N的增大,分类精度也呈现出先增大后减小的趋势。 实验过程中以OA的好坏为第一标准。 对于IP数据,考虑到训练时间的影响,综合各指标,设置λ1=0.003,λ2=0.015,λ3=0.01,λ4=0.02,λ5=0.05,σ1=1,σ2=1,σ3=1,σ4=2,σ5=3,K=30,N=20。

对于UP数据,σ调优过程的实验结果如表5所示。 由于UP数据包含103个光谱波段,所以将K和N分别固定在[10,60],具体结果见表6与表7。 最优结果加粗表示。

表5 σ在UP上的调优过程Table 5 Optimization process for σ on UP

表6 参数K对UP的分类效果Table 6 UP classification effect of parameter K

表7 参数N对UP的分类效果Table 7 UP classification effect of parameter N

由表5可知,当σ1=2,σ2=1,σ3=2,σ4=2,σ5=1时,UP数据的分类精度最高,OA达到98.90%,Kappa达到98.54%,AA达到88.61%。 同样的,对于UP数据,也研究了参数K和N对分类结果的影响。 从表6可以看出,随着K增大,OA,AA和Kappa大体上呈现出先增大后减小的趋势,当K=20时,各项指标最高,此时OA=98.95%,Kappa=98.61%,AA=88.65%。 由表7可知,当N=50时,分类精度最高,OA达到98.98%,Kappa达到98.62%,AA达到88.76%,除此之外,OA随着N的增加也呈现出先增大后减小的趋势。 综上,对于UP数据,设置λ1=0.003,λ2=0.015,λ3=0.01,λ4=0.02,λ5=0.05,σ1=2,σ2=1,σ3=2,σ4=2,σ5=1,K=20,N=50。

2.4 结果对比

为了验证本文提出的波段聚类和多尺度结构特征融合框架(ASPS-MRTV)的效果,分别与SVM、 ASPS(ED)、 ASPS(ID)、 ASPS-LBP、 ASPS-GlCM、 ASPS-BF等方法进行对比。 实验中统一使用C和σ都为24.5的SVM进行分类。

(1) SVM: 直接使用SVM对高光谱数据进行分类;

(2) ASPS(ED): 仅使用光谱特征提取框架-自适应子空间划分策略获取光谱维度信息,波段相似性度量使用欧式距离;

(3) ASPS(ID): 仅使用光谱特征提取框架-自适应子空间划分策略获取光谱维度信息,波段相似性度量使用信息散度;

(4) ASPS-LBP: 使用光谱特征提取框架-自适应子空间划分策略获取光谱维度信息,再利用局部二值模式(LBP)获取空间信息,设置B=10,R=1;

(5) ASPS-GlCM: 使用光谱特征提取框架-自适应子空间划分策略获取光谱维度信息,再利用灰度共生矩阵(GLCM)获取空间信息,选择4个不同的方向(0°,45°,90°,135°);

(6) ASPS-BF: 使用光谱特征提取框架-自适应子空间划分策略获取光谱维度信息,再利用双边滤波器(BF)获取图像的空间信息,设置w=15、δα=8、δr=3;

(7) ASPS-MRTV(AM): 使用我们提出的波段聚类和多尺度结构特征融合的方法获取光谱-空间信息。

首先对IP数据使用上述7种方法进行实验,得到的结果如表8所示。

表8 IP上不同方法的分类精度Table 8 Classification accuracy of different methods on IP

从表8中可以看出,与其他方法相比,我们提出的波段聚类和多尺度结构特征融合的空谱特征提取框架有着最高的分类性能,OA=97.06%,Kappa=96.64%,AA=85.95%。 由分类结果可知,SVM的分类精度高于ASPS(ED),验证了欧式距离(ED)虽然对计算波段间相似性有效,但对高维空间并不适合。 与ASPS(ED)方法相比,ASPS(ID)的OA提高了28.44%,验证了信息散度不仅对波段信息相似度有效,也可以在降低高光谱图像立方体数据维度的同时,保留更完整和有用的信息。 与SVM相比,ASPS(ID)的OA提高了9.8%,验证了光谱特征提取框架-自适应子空间划分策略对IP数据的类别可分性更加有效。 与ASPS-GlCM、 ASPS-LBP相比,ASPS-BF、 ASPS-MRTV的OA有了大幅度的提高,验证了对于高光谱数据而言,结构特征比纹理特征更能提高分类精度,这是由于结构特征提取技术更能去除噪声,纹理等无用的空间信息。 与ASPS-BF相比,ASPS-MRTV的分类精度提高了5.24%,验证了空间特征提取框架-多尺度相对全变分(MRTV)更能充分利用高光谱数据内部的空间信息,有效减少椒盐噪声和纹理特征的影响。 图6为七种方法对IP数据的分类结果图,可以看出与其他六种方法相比,ASPS-MRTV的错误分类标签最少,与地面真实图更相近。

图6 不同特征提取对IP分类图(a): 地面真实图; (b): SVM; (c): ASPS(ED); (d): ASPS(ID); (e): ASPS-LBP; (f): ASPS-GlCM; (g): ASPS-BF; (h): AMFig.6 IP classification maps by different feature extraction methods(a): Ground map; (b): SVM; (c): ASPS(ED); (d): ASPS(ID); (e): ASPS-LBP; (f): ASPS-GlCM; (g): ASPS-BF; (h): AM

对于UP数据,表9给出了不同方法获得的分类精度。 最优结果如表中加粗值所示,可以看到,OA=98.98%,Kappa=98.62%,AA=88.76%。 与SVM相比,其余六种方法都获得了较高的分类精度,验证了通过提取高光谱图像的空间和光谱信息,可以提高图像特征分类的准确性。 与ASPS(ED)相比,ASPS(ID)的OA提高了0.91%,再次验证了在细划分策略中,使用信息散度构造波段间相似矩阵的有效性。 与ASPS(ID)和ASPS(ED)相比,ASPS-LBP、 ASPS-GlCM、 ASPS-BF、 ASPS-MRTV的分类精度更高,这是因为前两种方法仅考虑了高光谱数据内部的光谱特征,忽略了图像内部相邻像素之间的空间相关性,导致分类精度较低,也再次说明了空谱特征提取方式对于识别不同地物类别非常有效。 与ASPS-GlCM和ASPS-LBP相比,ASPS-BF和ASPS-MRTV分类效果更好,再次验证了对于高光谱图像立方体而言,充分利用数据内部的结构特征对分类精度提升更有效。

表9 UP不同方法的分类精度Table 9 Classification accuracies of different methods for UP

图7为七种方法对UP数据的分类效果图。 可以看出,SVM与仅使用光谱特征提取的ASPS(ID),ASPS(ED)存在更多的错误分类标签,相比之下,ASPS-BF和ASPS-MRTV分类误点较少,证明了结构特征对严重噪声问题的鲁棒性。

3 结 论

提出了一种波段聚类和多尺度结构特征融合的高光谱图像分类模型(ASPS-MRTV)。 一方面,为了充分提取图像内部的光谱信息,构造了基于自适应子空间划分策略的特征提取框架,首先通过粗划分对原始数据进行等分,之后通过细划分策略对等子空间进行自适应,最后将每个自适应子空间的平均值进行叠加形成光谱特征。 另一方面,构造了多尺度相对全变分方法,用于提取高光谱图像的结构特征,通过设置多参数,对粗细划分策略后的光谱数据进行不同尺度下的RTV提取,达到了信息互补的效果,最后使用KPCA对堆叠数据进行处理,得到高光谱图像立方体的空谱特征。 与其余六种方法相比较,ASPS-MRTV获得了更高的精度,但也存在一些不足,比如运行时间过长,在后续的研究中,将进一步优化模型的时间复杂度,以提取更具判别力的空谱特征。