CN102208034B - 基于半监督维数约减的高光谱图像分类方法 - Google Patents

Abstract

本发明公开了一种基于半监督维数约减的高光谱图像分类方法，主要解决高光谱图像数据维数过高而导致的计算量大以及现有方法的分类正确率低的问题。其步骤包括：将高光谱图像的每一个像素点用特征向量表示，选出有标记训练集、测试集以及总训练集；分别构造有标记训练集的局部类间和局部类内不相似性矩阵，得到总的局部不相似性矩阵；构造特征值方程并求解，得到投影矩阵；分别将有标记训练集和测试集投影到低维空间，得到新的有标记训练集和测试集；将新的有标记训练集和测试集输入支撑矢量机进行分类，得到测试集的类别信息。本发明采用了半监督的思想，能够获得较高的分类正确率，可用于地图制图，植被调查，军事情报获取领域。

Description

基于半监督维数约减的高光谱图像分类方法

技术领域

本发明属于图像处理技术领域，涉及高维数据的降维，用于对高光谱遥感图像的分类。

背景技术

高光谱遥感是遥感领域在20世纪80年代起最重要的发展之一，也是当今及今后几十年内遥感的前沿技术。高光谱遥感技术利用成像光谱仪以纳米级的光谱分辨率，以几十或几百个波段同时对地表物成像，能够获得地物的连续光谱信息，实现地物空间信息、辐射信息、光谱信息的同步获取，具有“图谱合一”的特性，从而使人类对地观测和信息获取能力向前迈进了一大步。

常用的高光谱图像数据包括由美国宇航局NASA喷气推进实验室的空载可见光/红外成像光谱仪AVIRIS获得的Indian Pine数据集和Kennedy Space Center(简称KSC)数据集，以及NASA的EO-I HYPERION光谱仪获得的Botswana数据集等。

遥感图像分类就是将遥感图像中的每个像元划归到类别中去的过程。高光谱图像分类是建立在遥感图像分类的基础上，结合高光谱遥感图像特点，对高光谱图像进行像元的识别和分类。由于高光谱遥感技术获取的图像包含了丰富的空间、辐射和光谱三重信息，这些信息特别适合分类，但仍存在巨大的挑战和困难：(1)数据量大，至少几十个波段，导致计算复杂度很高，也给数据的存储、传递和显示带来了挑战；(2)维数过高，存在冗余数据，会降低分类精度；(3)波段多，且波段间相关性高，导致所需训练样本数目增多，如果训练样本不足，将导致从训练样本得到的参数不可靠。因此如何有效的降低高光谱图像的维数，减少数据量是高光谱图像分类中的一个重要问题。

目前，已经有许多有效的维数约简方法，如主成分分析PCA，核主成分分析KPCA，线性判别分析LDA，核判别分析KDA，局部线性嵌入LLE等。在这些已有的降维方法中，PCA和LDA是最著名，也是应用最广的线性降维方法。PCA的主要目标是寻找在最小均方意义下最能代表原始数据的投影方向，该投影方向是通过线性变换得到的一组最优的单位正交向量基，这些向量的线性组合可以重构原始样本，并且重构后的样本和原样本之间的误差最小。PCA在许多模式识别应用中取得了较好的效果，但是由于它是一种无监督的降维方法，不适用于反映样本之间的差异。相对于PCA方法，LDA是一种有监督的线性降维方法，主要目标就是寻找在最小均方意义下，最能够分开各类数据的投影方向，实现上是通过寻找一组线性变换以达到类内散度最小且类间散度最大的目的。目前，LDA是最基本也是应用最广的降维方法。在此基础上，已经有学者利用核方法将PCA和LDA推广到非线性形式，得到核主成分分析KPCA和核判别分析KDA。

上述PCA和KPCA属于无监督降维方法，因为没有有效的监督信息即类别信息的指导，往往难以得到较好的后续分类结果；LDA和KDA属于有监督降维方法，其虽然可以充分利用大量的监督信息而取得较为理想的结果，但是获取大量监督信息需要付出很大的代价。尤其对于高光谱图像数据，其类别信息的获取是耗时耗力且非常困难的。

发明内容

本发明的目的在于同时利用少量的有标记数据以及大量的无标记数据，提出一种基于半监督维数约减的高光谱图像分类方法，以减小计算复杂度，并在有标记样本较少的情况下提高分类正确率。

实现本发明目的的技术方案是，在利用少量有标记样本提供判别信息的基础上，通过引入一个正则项，加入由大量无标记样本得到的几何结构信息，获得更优的低维空间，将原始高维数据转换到低维空间中再进行分类，从而获得更好的分类性能。其实现步骤包括如下：

(1)输入高光谱图像，该高光谱图像包含c类像素点，其中有N个有标记像素点，m个无标记像素点，图像的每一个像素点为一个样本，每个样本用一个特征向量表示，样本的特征维数为D；

(2)用N个有标记样本构成有标记训练集其对应的类别标记集为

用m个无标记样本构成测试集其中，x_i表示有标记训练集的第i个有标记训练样本，y_j表示测试集的第j个测试样本，每个样本用一个列向量表示，l_i是第i个有标记训练样本所属的类别标号，N是有标记训练样本总个数，m是测试样本总个数，D表示样本的特征维数，R表示实数域；

(3)从m个无标记样本中选出η个，与N个有标记样本组合，得到总训练集

其中，s_k表示总训练集的第k个训练样本，总训练集的样本个数为N+η，D表示样本的特征维数，R表示实数域；

(4)构造有标记训练集X的局部类间不相似性矩阵和局部类内不相似性矩阵

$\tilde{C}=\underset{p=1}{\overset{c}{\mathrm{\Σ}}}\underset{h=1}{\overset{n_p}{\mathrm{\Σ}}}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{N}_{\mathrm{hi}}^b(x_h-x_i){(x_h-x_i)}^T$

$\tilde{A}=\underset{p=1}{\overset{c}{\mathrm{\Σ}}}\underset{h=1}{\overset{n_p}{\mathrm{\Σ}}}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{N}_{\mathrm{hi}}^w(x_h-x_i){(x_h-x_i)}^T$

其中，c表示样本的类别数，n_p表示有标记训练集中第p类的有标记训练样本个数，N表示有标记训练样本总个数，x_h表示第p类的第h个有标记训练样本，x_i表示有标记训练集中的第i个有标记训练样本，T表示矩阵的转置，

为调节类间近邻数的系数，

为调节类内近邻数的系数，

式中，q_b表示类间近邻数，即从第p类以外的有标记训练样本中选出的与x_h距离最近的近邻样本个数，N_b(x_h)定义为有标记训练样本x_h的类间q_b近邻域，下标b均表示类间，q_w表示类内近邻数，即从第p类有标记训练样本中选出的与x_h距离最近的近邻样本个数，N_w(x_h)定义为有标记训练样本x_h的类内q_w近邻域，下标w均表示类内；

(5)对所述

和

进行求和，得到有标记训练集总的局部不相似性矩阵： $\tilde{T}=\tilde{C}+\tilde{A};$

(6)构造总训练集S的相似度矩阵M，计算拉普拉斯矩阵L＝D′-M，其中，D′是对角矩阵，对角元素为

M_kt为矩阵M的第k行t列元素；

(7)根据有标记训练集X的局部类间不相似性矩阵

和总的局部不相似性矩阵

拉普拉斯矩阵L，以及总训练集S，构造特征值方程

e＝1，Λd，其中，β是权重系数，v_e表示第e个特征向量，λ表示特征值，T表示矩阵的转置，d表示降维后的维数；

(8)设置降维后的维数d的值，求解特征值方程

e＝1，Λd，得到d个最大特征值及其对应的特征向量v₁，v₂，Λv_d，构成投影矩阵V＝[v₁，v₂，Λv_d]∈R^D×d，D表示样本的原始特征维数，R表示实数域；

(9)分别将有标记训练集

和测试集

投影到由投影矩阵V所张成的低维空间，得到投影后的新的有标记训练集和新的测试集

其中，x′_i＝<V^T·x_i>为第i个有标记训练样本的新的特征向量，x_i为原有标记训练集的第i个训练样本，y′_j＝<V^T·y_j>为第j个测试样本的新的特征向量，y_j为原测试集的第j个测试样本，N表示有标记训练样本总个数，m表示测试样本总个数，d表示新的样本集的特征维数，R表示实数域；

(10)将新的有标记训练集

和类别标记集

以及新的测试集

输入到支撑矢量机SVM中进行分类，得到测试集的分类结果

其中，l′_j表示第j个测试样本所属的类别标号。

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

1、本发明由于采用维数约减算法对高光谱图像数据进行降维后，再进行分类，大大减少了计算量，提高了分类的速度。

2、本发明由于同时利用了有标记样本和无标记样本，能够获得更理想的低维投影空间，从而获得高维数据在低维空间中更精确的表示，提高了后续分类的准确率；

对比实验表明，本发明有效降低了计算的复杂度，提高了高光谱遥感图像的分类准确率。

附图说明

图1是本发明的流程示意图；

图2是本发明仿真采用的Indian Pine图像；

图3是本发明与现有方法在有标记训练样本个数不同时，采用不同的特征维数得到的分类结果：

图4是本发明与现有方法对Indian Pine图像的分类结果图。

具体实施方式

参照图1，本发明的具体实施步骤包括：

步骤1，输入高光谱图像，从高光谱图像的每一类像素点中选出等量的像素点，作为有标记像素点，有标记像素点的总个数为N，该高光谱图像的剩余m个像素点作为无标记像素点，用每一个像素点的各个波段的灰度值作为该像素点的特征向量，每个像素点就是一个样本，样本的原始特征维数为D。

步骤2，选取有标记训练集X、测试集Y以及总训练集S。

2a)用N个有标记样本构成有标记训练集

其对应的类别标记集为

其中，x_i表示有标记训练集的第i个有标记训练样本，每个有标记训练样本用一个列向量表示，l_i是第i个有标记训练样本所属的类别标号，N是有标记训练样本总个数，D表示样本的特征维数，R表示实数域；

2b)用m个无标记样本构成测试集

其中，y_j表示测试集的第j个测试样本，每个测试样本用一个列向量表示，m表示测试样本个数，D表示样本的特征维数，R表示实数域；

2c)从m个无标记样本中选出η个，作为无标记训练样本，与N个有标记样本组合，得到总训练集

其中，s_k表示总训练集的第k个训练样本，总训练集的样本个数为N+η，D表示样本的特征维数，R表示实数域。

步骤3，构造有标记训练集X的局部类间不相似性矩阵

和局部类内不相似性矩阵

3a)计算有标记训练集中两两样本之间的欧氏距离：d_ig＝||x_i-x_g||，其中，x_i和x_g分别为有标记训练集中第i个和第g个有标记训练样本，d_ig表示第i个和第g个有标记训练样本之间的欧式距离；

3b)设置类间近邻数q_b和类内近邻数q_w的值，根据有标记训练样本之间的欧氏距离大小，对有标记训练集中每一个属于第p类的有标记训练样本x_h，从第p类以外的有标记训练样本中选出与x_h距离最近的q_b个类间近邻样本，从第p类有标记训练样本中选出与x_h距离最近的q_w个类内近邻样本，将N_b(x_h)定义为有标记训练样本x_h的类间q_b近邻域，将N_w(x_h)定义为有标记训练样本x_h的类内q_w近邻域，下标b和w分别用来区分类间与类内，p的取值范围是1到c，c表示样本的类别数，即原始像素点的类别数，x_h表示第p类的第h个有标记训练样本，n_p表示有标记训练集中第p类的有标记训练样本个数；

3c)根据步骤3b)中确定的参数，分别计算调节类间近邻数的系数

和调节类内近邻数的系数

3d)根据步骤3b)和步骤3c)中确定的参数，分别计算每一类有标记训练样本的局部类间不相似性矩阵

及局部类内不相似性矩阵

${\tilde{C}}_p=\underset{h=1}{\overset{n_p}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{N}_{\mathrm{hi}}^b(x_h-x_i){(x_h-x_i)}^T,$

${\tilde{A}}_p=\underset{h=1}{\overset{n_p}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{N}_{\mathrm{hi}}^w(x_h-x_i){(x_h-x_i)}^T;$

其中，N是有标记训练样本总个数，T表示矩阵的转置；

3e)对每一类有标记训练样本的局部类间不相似性矩阵

进行求和，得到有标记训练集的局部类间不相似性矩阵： $\tilde{C}=\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\tilde{C}}_p=\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}\underset{p=1}{\overset{n_p}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{N}_{\mathrm{hi}}^b(x_h-x_i){(x_h-x_i)}^T;$

3f)对每一类有标记训练样本的局部类内不相似性矩阵

进行求和，得到有标记训练集的局部类内不相似性矩阵： $\tilde{A}=\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}{\tilde{A}}_p=\underset{p=1}{\overset{c}{\mathrm{\&Sigma;}}}\underset{h=1}{\overset{n_p}{\mathrm{\&Sigma;}}}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{N}_{\mathrm{hi}}^w(x_h-x_i){(x_h-x_i)}^T.$

步骤4，对所述

和

进行求和，得到有标记训练集总的局部不相似性矩阵： $\tilde{T}=\tilde{C}+\tilde{A}.$

步骤5，构造总训练集S的相似度矩阵M，并计算拉普拉斯矩阵L＝D′-M。

5a)计算总训练集S中任意两个样本之间的相似度：

其中，s_k和s_t分别是总训练集S中的第k个和第t个训练样本，σ为相似度矩阵的核参数，exp{}为自然底数e的指数函数；

5b)由所述M_kt得到相似度矩阵M，即M_kt作为矩阵M的第k行第t列元素，矩阵M的大小为N+η行N+η列，N+η为总训练集的样本个数；

5c)由所述M_kt计算对角矩阵D′，对角元素计算公式为：

5d)由相似度矩阵M和对角矩阵D′，得到拉普拉斯矩阵L＝D′-M。

步骤6，构造特征值方程。

6a)设置权重系数β，根据总训练集S和拉普拉斯矩阵L，构造正则项矩阵βSLS^T，该正则项矩阵大小为D×D，D为原始样本特征维数；

6b)根据有标记训练集的局部类间不相似性矩阵

和总的局部不相似性矩阵

以及正则项矩阵βSLS^T，得到特征值方程：

e＝1，Λd，其中，v_e表示第e个特征向量，λ表示特征值，T表示矩阵的转置，d表示降维后的维数。

步骤7，求解特征值方程，得到投影矩阵。

7a)设置降维后的维数d的值，对特征值方程

e＝1，Λd进行变形，得到矩阵

再对该矩阵进行特征分解，得到特征值

以及对应的特征向量其中，λ_e表示矩阵S的第e个特征值，v_e为特征值λ_e对应的特征向量，r表示特征值的个数；

7b)取特征值中d个最大特征值，用这d个特征值对应的特征向量v₁，v₂，Λv_d，构成投影矩阵V＝[v₁，v₂，Λv_d]∈R^D×d，D表示样本的原始特征维数，R表示实数域。

步骤8，分别将有标记训练集X和测试集Y投影到由投影矩阵V所张成的低维空间，得到投影后的新的有标记训练集和新的测试集。

8a)将原始的有标记训练集

投影到由投影矩阵V张成的空间中，得到新的有标记训练集

其中，x′_i＝<V^T·x_i>为第i个有标记训练样本的新的特征向量，x_i为原有标记训练集的第i个训练样本，N表示有标记训练样本总个数，d表示新的有标记训练集的样本特征维数，R表示实数域；

8b)将原始的测试集

投影到由投影矩阵V张成的空间中，得到新的测试集

其中，y′_j＝<V^T·y_j>为第j个测试样本的新的特征向量，y_j为原测试集的第j个测试样本，m表示测试样本个数，d表示新的测试集的样本特征维数，R表示实数域。

步骤9，将新的有标记训练集

和类别标记集

以及新的测试集输入到支撑矢量机SVM中进行分类，得到测试集的识别结果其中，l′_j表示第j个测试样本所属的类别标号。

本发明的效果可以通过以下仿真实验进一步说明：

1.仿真条件：

仿真实验采用美国宇航局NASA喷气推进实验室的空载可见光/红外成像光谱仪AVIRIS于1992年6月在印第安纳西北部获取的Indian Pine图像，如图2所示，图像大小为145×145，共220个波段，去除噪声以及大气和水域吸收的波段还有200个波段，共16类地物信息，由于部分类别的数据个数非常少，在仿真实验中，只考虑表1所示的9类数据，整幅图像被分为9类。

仿真实验在CPU为Intel Core(TM)2Duo、主频2.33GHz，内存为2G的WINDOWSXP系统上用MATLAB 7.0.1软件进行。

表1IndianPine图像中的9类数据

类别	类别名称	个数
			1	Corn-no till	1434
2	Corn-min	834
			3	Grass/Pasture	497
4	Grass/Trees	747
			5	Hay-windrowed	489
6	Soybeans-no till	968
			7	Soybeans-min	2468
8	Soybean-clean	614
			9	Woods	1294

2.仿真内容及分析：

使用本发明与现有三种方法对高光谱图像Indian Pine进行分类，现有三种方法分别是：核判别分析结合支撑矢量机KDA+SVM，主成分分析结合支撑矢量机PCA+SVM，核主成分分析结合支撑矢量机KPCA+SVM。本发明基于半监督维数约减结合支撑矢量机的高光谱图像分类方法的缩写为SSLGC+SVM。

各分类方法KDA+SVM、PCA+SVM、KPCA+SVM和本发明SSLGC+SVM中的分类器SVM核参数γ分别设置为10^-7、10^-8、10^-18和10^-8，惩罚因子C一律设置为10000。本发明的类间近邻数q_b设置为10，类内近邻数q_w设置为7，相似度矩阵的核参数σ设置为200，权重参数β设置为1，无标记训练样本个数η固定为2000。

A.从表1所示9类数据中选取部分像素点作为有标记像素点，这9类数据中的剩余像素点作为无标记像素点，用本发明与现有三种方法对9类数据进行20次分类实验，取分类结果的平均值，作为最终分类正确率，如图3所示，其中图3a是四种方法在每类有标记训练样本个数为8时的分类正确率与特征维数的曲线关系图，图3b是四种方法在每类有标记训练样本个数为12时的分类正确率与特征维数的曲线关系图，横坐标为特征维数，维数范围是1至60，纵坐标为分类正确率。从图3a可以看出当特征维数大于10时，本发明的分类正确率达到65％以上，明显高于现有方法；从图3b中可以看出当特征维数大于10时，本发明的分类正确率达到70％以上，明显高于现有方法。从图3a和图3b中还可以看出，在维数大于15后，本发明的结果趋于稳定，因此只需要采用15维特征，即可以得到较高的识别率，从而大大减少了计算量。

B.从表1所示的9类数据的每一类像素点中选出10个作为有标记像素点，整幅Indian Pine图像的剩余像素点作为无标记像素点，用四种方法对整幅Indian Pine图像的所有像素点进行分类，各方法中降维后的特征维数设置为15，结果如图4所示，其中图4a是KDA+SVM的分类结果图，图4b是PCA+SVM的分类结果图，图4c是KPCA+SVM的分类结果图，图4d是本发明的分类结果图，从图中可以看出，本发明的结果图最为清晰，边缘一致性比现有方法更好。

综上，本发明在基于半监督维数约减的基础上结合支撑矢量机对高光谱图像进行分类，充分利用有标记与无标记信息，能够减少计算量，并得到较高的识别率，与现有的方法相比具有一定的优势。

Claims (3)

1.一种基于半监督维数约减的高光谱图像分类方法，包括如下步骤：

(1)输入高光谱图像，该高光谱图像包含c类像素点，其中有N个有标记像素点，m个无标记像素点，图像的每一个像素点为一个样本，每个样本用一个特征向量表示，样本的特征维数为D；

(2)用N个有标记样本构成有标记训练集 

其对应的类别标记集为 

用m个无标记样本构成测试集 

其中，x_i表示有标记训练集的第i个有标记训练样本，y_j表示测试集的j个测试样本，每个样本用一个列向量表示，l_i是第i个有标记训练样本所属的类别标号，N是有标记训练样本总个数，m是测试样本总个数，D表示样本的特征维数，R表示实数域；

(3)从m个无标记样本中选出η个，与N个有标记样本组合，得到总训练集 其中，s_k表示总训练集的第k个训练样本，总训练集的样本个数为N+η，D表示样本的特征维数，R表示实数域；

(4)构造有标记训练集X的局部类间不相似性矩阵 

和局部类内不相似性矩阵 

其中，c表示样本的类别数，n_p表示有标记训练集中第p类的有标记训练样本个数，N表示有标记训练样本总个数，x_h表示第p类的第h个有标记训练样本，x_i表示有标记训练集中的第i个有标记训练样本，T表示矩阵的转置， 

为调节类间近邻数的系数， 

为调节类内近邻数的系数， 式中，q_b表示类间近邻数，即从第p类以外的有标记训练样本中选出的与x_h距离最近的近邻样本个数，N_b(x_h)定义为有标记训练样本x_h的类间q_b近邻域，下标b均表示类间，q_w表示类内近邻数，即从第p类有标记训练样本中选出的与x_h距离最近的近邻样本个数，N_w(x_h)定义为有标记训练样本x_h的类内q_w近邻域，下标w均表示类内；

(5)对所述 

和 

进行求和，得到有标记训练集总的局部不相似性矩阵：

(6)构造总训练集S的相似度矩阵M，计算拉普拉斯矩阵L=D′-M，其中，D′是对角矩阵，对角元素为 

M_kt为矩阵M的第k行t列元素；

(7)根据有标记训练集X的局部类间不相似性矩阵 

和总的局部不相似性矩阵 

拉普拉斯矩阵L，以及总训练集S，构造特征值方程 

其中，β是权重系数，v_e表示第e个特征向量，λ表示特征值，T表示矩阵的转置，d表示降维后的维数；

(8)设置降维后的维数d的值，求解特征值方程 得到d个最大特征值及其对应的特征向量v₁，v₂，…v_d，构成投影矩阵V=[v₁，v₂，…v_d]∈R^D×d，D表示样本的原始特征维数，R表示实数域；

(9)分别将有标记训练集 

和测试集 

投影到由投影矩阵V所张成的低维空间，得到投影后的新的有标记训练集 

和新的测试集 

其中，x′_i=<V^T·x_i>为第i个有标记训练样本的新的特征向量，x_i为原有标记训练集的第i个训练样本，y′_j=<V^T·y_j>为第j个测试样本的新的特征向量，y_j为原测试集的第j个测试样本，N表示有标记训练样本总个数，m表示测试样本总个数，d表示新的样本集的特征维数，R表示实数域； 

(10)将新的有标记训练集 

和类别标记集 

以及新的测试集 

输入到支撑矢量机SVM中进行分类，得到测试集的分类结果 

其中，l′_j表示第j个测试样本所属的类别标号。

2.根据权利要求1所述的高光谱图像分类方法，其中步骤(6)所述的构造总训练集S的相似度矩阵M，按如下矩阵元素公式计算： 

其中，M_kt表示矩阵M的第k行第t列元素，s_k和s_t分别是总训练集S中的第k个和第t个训练样本，σ为核参数，exp{}为自然底数e的指数函数，矩阵M的大小为N+η行N+η列，N+η为总训练集S的样本个数。

3.根据权利要求1所述的高光谱图像分类方法，其中步骤(7)所述的构造特征值方程

按如下步骤进行：

7a)设置权重系数β，根据总训练集S和拉普拉斯矩阵L，构造正则项矩阵βSLS^T,该正则项矩阵大小为D×D，D为原始样本特征维数；

7b)根据有标记训练集的局部类间不相似性矩阵 

和总的局部不相似性矩阵 

以及正则项矩阵βSLS^T，得到特征值方程： 

其中，v_e表示第e个特征向量，λ表示特征值，T表示矩阵的转置，d表示降维后的维数。 

Images