CN103942749B - 一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法 - Google Patents

Abstract

本发明属于机器学习技术领域，可用于高光谱地物分类，具体公开了一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，主要解决了高光谱地物分类中由于标记样本少，混合像元现象，地物识别正确率不高等问题。其步骤为：首先对样本进行预处理，并把样本分为训练样本与测试样本；其次计算出训练样本的图Laplacian矩阵；再次随机的产生网络的输入权值与隐藏层偏置，利用标记样本计算网络输出权值的初始值；然后交替计算每个未标记样本属于每类的概率与网络的输出权值矩阵，直到目标函数相邻两次的差值在容许误差内；最后利用优化的模型预测测试样本标签。本方法具有较强实用性与鲁棒性，实验结果说明该发明很适合高光谱地物分类。

Description

一种基于修正聚类假设和半监督极速学习机的高光谱地物分 类方法

技术领域

本发明属于机器学习技术领域，可用于高光谱地物分类。具体涉及一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法。

背景技术

利用高光谱所提供的丰富地面信息实现地物分类是近几年很热的一个研究课题。随着高光谱成像技术的进步，获取海量数据已经变得相当容易，而且代价较低，但是这些样本点都是未标记的。而要获得标记样本则代价较高，需要专家标记和实地考察等。而且有限的标记样本在空间分布上很可能不均匀，因此仅仅利用这些标记样本很难构造准确的分类器。另外，随着成像技术的高度发展，虽然高光谱数据的光谱分辨率得到大大提高，但是其空间分辨率没有得到相应的提高，因此，在高光谱数据中，一个像素点很可能由不同物质组成，即其中存在着混合像元。这就导致了同物异谱和同谱异物现象的出现。所谓的同物异谱，指的是处于不同位置的同一种物质，其光谱响应曲线差别很大；而同谱异物指的是，不同种类的物质有可能具有相同的光谱响应曲线。因此对于确定某一像素属于哪一类而不属于其它类变得更困难。

传统的解决高光谱地物分类的方法有监督学习，半监督学习与无监督学习。监督学习的分类策略是通过训练大量的标记样本来构造精确的分类器，因而在标记样本很有限的前提下使用这种方法分类时往往很难正确的得到遥感图像中的地物统计分布。无监督学习在学习的过程中是直接依靠数据本身的特性进行模型的构建，因此没有利用到标记样本中所含的先验信息，所以分类精度也不理想。而半监督学习方法在学习过程中融合了标记样本和未标记样本的信息，既利用了标记样本的先验信息，又利用海量未标记样本所含的信息来修正学习过程，即同时利用标记样本和未标记样本来构造分类器，因此是目前高光谱地物分类的主流技术。常见的高光谱地物分类方法有：直推式支持向量机(TransductiveSupport Vector Machine,TSVM)、基于图的半监督分类方法等。TSVM是支持向量机(SVM)在半监督分类中的推广，但是由于部分训练样本的类标记缺失，所以求解TSVM的最优解是一个NP难问题，实际运用中所得到的常常是它的近似最优解。和TSVM方法相比，基于图的半监督分类方法具有很好的数学解释性和良好的学习性能，并且能够避免收敛到局部最优，但是该类方法没有考虑到前面所提到的混合像元问题。

极速学习(Extreme Learning Machine,ELM)算法是最近几年提出的一种新颖的机器学习算法，是以单隐层前馈神经网络为模型的。由于其输入层权值和隐藏层偏值随机设置，输出层权值能够以解析的形式求出来，所以在大规模数据分类中表现出运行速度快，分类正确率高，泛化性能好等特点。该方法属于监督学习模式。

如何利用极速学习机的优势，并把它很好的应用在当下的高光谱地物分类问题中，是一个很值得研究的内容。但是目前这方面的工作很少。

发明内容

本发明综合利用了上述已有技术的优点，并克服一些技术中存在的缺点，即利用了极速学习机学习速度快并把它扩展到半监督学习领域，利用了半监督学习的优势克服了标记样本少的缺点，并且利用修正聚类假设解决混合像元问题，提出一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，同时提高分类的效率与正确率。

为实现上述目的，本发明的技术方案是：提供一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，具体步骤包括：

(1)预处理训练样本；

用滑窗技术对所有样本进行预处理，窗口大小为5*5，随机选择40％的样本做为训练样本X，剩余的做测试样本X_t，其中X＝{x₁,...,x_l,x_l+1...,x_l+u}，x₁,...,x_l为标记样本，x_l+1...,x_l+u为未标记样本，N＝l+u为总的训练样本数，每个样本x＝[x₁,x₂,…,x_d]，x_i表示每个样本的第i个特征值，其中i＝1,2,…d，d代表样本的特征维数。Y_l＝[y₁,y₂,...,y_l]表示与训练样本中的标记样本x₁,...,x_l对应的标签，其中y_i＝[y_i1,y_i2,...,y_im]^T表示与第i个标记样本对应的标签向量，m是类别数，T表示转置符号。

(2)计算训练样本的图Laplacian矩阵；

(3)计算网络的输出权值矩阵初始值；

(4)计算样本属于每类的概率P；

$p_k\left(x_j\right)=\frac{\mathrm{\λ}\underset{t=1,t\≠j}{\overset{N}{\mathrm{\Σ}}}{p}_k\left(x_t\right)w_{\mathrm{tj}}}{{\left|\right|f\left(x_j\right)-s_k\left|\right|}^2+\mathrm{\λ}\underset{s=1,s\≠j}{\overset{N}{\mathrm{\Σ}}}{w}_{\mathrm{js}}},(k=1,...m;j=l+1,...,N)$

其中p_k(x_j),k＝1,...,m；j＝l+1,...,N表示第j个未标记样本x_j属于第k类的概率，m为类别数目，N＝l+u为总的训练样本数，λ是正则参数，w_js代表像素x_j与像素x_s之间的相似度，f(x_j)表示样本x_j的预测值，是第k类的标准标签，即该向量中第k个元素为1，其余的元素为0。

(5)计算新的输出权值矩阵β；

$\mathrm{\β}=(Y_l{H_l}^T+S\hat{P}{J}^T{H}_u^T){(H_l{H}_l^T+H_{u}J\hat{P}{J}^T{H}_u^T)}^{-1}$

H_l是标记样本的隐藏层输出矩阵，I_u是u×u的单位阵，其中每个S_k是一个m×u的矩阵，该矩阵的k行全为1，其余行全为0；是未标记样本属于各个类别的概率值，是对角矩阵，对角元素为P按行排列时位置值的平方，H_u是未标记样本的隐藏层输出矩阵。

(6)判断目标函数变化是否小于ε；

目标函数为

$\mathrm{tr}(({\mathrm{\βH}}_l-Y_l){({\mathrm{\βH}}_l-Y_l)}^T+({\mathrm{\βH}}_{u}J-S)\hat{P}{({\mathrm{\βH}}_{u}J-S)}^T+{\mathrm{\λPLP}}^T)$

如果相邻目标函数的差值大于ε，则从步骤4开始重复执行步骤4到步骤6；如果小于ε，则往下执行步骤7；

(7)预测测试样本标签；

按照下式计算

$Y_t=H_t\widehat{\mathrm{\β}}$

其中 $Y_t={\left[\begin{array}{c}{y_1}^T\\ \·\\ \·\\ \·\\ {y_t}^T\end{array}\right]}_{t\×m}$ 表示对于测试样本的网络输出矩阵，y_i＝[y_i1,y_i2,...,y_im]^T表示每一个样本的输出向量，表示优化后的输出权值矩阵，H_t表示测试样本的隐藏层输出矩阵。

上述步骤(2)的图Laplacian矩阵的计算步骤具体如下：

2.1)首先计算图的权重矩阵W：

$W=\left\{\begin{array}{cc}{w}_{\mathrm{ij}}& \mathrm{ife}(i,j)\∈E\\ 0& \mathrm{otherwise}\end{array}\right.$

其中E代表边集，表示样本间的某种相似性；边e(i,j)的权重为w(e)＝w_ij，表示顶点i与顶点j，即样本i与样本j之间的相似度，并且w_ij＝w_ji，边的权重计算公式如下：

$w_{\mathrm{ij}}=\left\{\begin{array}{cc}1\sqrt{d_{\mathrm{ij}}}& x_i\∈\mathrm{NB}\left(x_j\right){\mathrm{orx}}_j\∈\mathrm{NB}\left(x_i\right)\\ 0& \mathrm{otherwise}\end{array}\right.$

其中NB(x_i)是以x_i为中心的一个r*r的方形窗中的样本集，d_ij代表样本x_i与样本x_j之间的欧氏距离；

2.2)再计算出图Laplacian矩阵L：

L＝I-D^-1/2WD^-1/2

其中I是单位矩阵,D是由所有结点的度数构成的对角阵，图中结点的度数定义为与该结点相连的所有边的权值之和，即

上述步骤(3)的网络输出权值矩阵初始值的计算步骤具体如下：

3.1)随机产生网络的输入权值w_i,与网络隐藏层偏值b_i，表示隐藏层神经元个数；

3.2)利用sigmoid函数作为极速学习网络的隐藏层神经元激活函数，计算出隐藏层输出矩阵：

$H_l(w_1,...,w_{\tilde{N}},b_1,...,b_{\tilde{N}},x_1,...,x_{\tilde{N}})={\left[\begin{array}{ccc}\mathrm{\φ}(w_1\·x_1+b_1)& \·\·\·& \mathrm{\φ}(w_{\tilde{N}}\·x_1+b_{\tilde{N}})\\ \·& & \·\\ \·& \·\·\·& \·\\ \·& & \·\\ \mathrm{\φ}(w_1\·x_l+b_1)& \·\·\·& \mathrm{\φ}(w_{\tilde{N}}\·x_l+b_{\tilde{N}})\end{array}\right]}_{l\×\tilde{N}}$

其中,w_i·x_j表示w_i与x_j的内积；

3.3)则输出权值矩阵为

Y_l是标记样本标签，表示H_l的Moore-Penrose广义逆运算。

本发明的优点是：首先用5*5的滑窗技术对样本进行预处理，然后随机选择40％的样本做为训练样本，剩余的做测试样本，训练样本中只标记很少的样本；其次，对训练样本构图，算出其Laplacian图矩阵；再次，随机产生网络的输入权值与隐藏层神经元偏置，利用sigmoid函数作为网络隐藏层的激活函数，并利用广义逆(Moore-Penrose)计算输出权值矩阵；然后，根据修正聚类假设，以上一步计算得到的输出权值矩阵做为初始值，交替计算每一样本属于每一类的概率以及输出权值矩阵，直到达到给定阈值为止；最后得到优化好的半监督极速学习机模型，对测试样本进行标签预测。克服了学习速度慢，解决了标记样本少与混合像元问题，同时提高分类的效率与正确率。

以下将结合附图对本发明做进一步详细说明。

附图说明

图1是本发明的总流程图；

图2是本发明输出权值矩阵初始化子流程图；

图3是本发明仿真采用的Indiana高光谱地物分类标记图；

图4是在每类标记的样本数从2个增加到10个时，不同方法分类正确率的变化图；

图5是在每类标记10个样本时，不同方法分类结果的标记图。

具体实施方式

参照附图1，本发明的具体步骤如下：

步骤1.预处理训练样本

用滑窗技术对所有样本进行预处理，窗口大小为5*5，随机选择40％的样本做为训练样本X，剩余的做测试样本X_t；X＝{x₁,...,x_l,x_l+1...,x_l+u}，前l个为标记样本，剩余的u个为未标记样本，N＝l+u为总的训练样本数；Y_l＝[y₁,y₂,...,y_l]，y_i＝[y_i1,y_i2,...,y_im]^T为第i个标记样本对应的标签。

步骤2.计算训练样本的图Laplacian矩阵

2.1)首先计算图的权重矩阵W：

$W=\left\{\begin{array}{cc}{w}_{\mathrm{ij}}& \mathrm{ife}(i,j)\∈E\\ 0& \mathrm{otherwise}\end{array}\right.$

其中E代表边集，表示样本间的某种相似性。边e＝(i,j)的权重为w(e)＝w_ij，表示顶点i,j，即样本i与样本j之间的相似度，并且w_ij＝w_ji。边的权重计算公式如下：

$w_{\mathrm{ij}}=\left\{\begin{array}{cc}1\sqrt{d_{\mathrm{ij}}}& x_i\∈\mathrm{NB}\left(x_j\right){\mathrm{orx}}_j\∈\mathrm{NB}\left(x_i\right)\\ 0& \mathrm{otherwise}\end{array}\right.$

其中NB(x_i)是以x_i为中心的一个r*r的方形窗中的样本集，本发明中r＝3,d_ij代表样本x_i与样本x_j之间的欧氏距离。

2.2)再计算出图Laplacian矩阵L：

L＝I-D^-1/2WD^-1/2

其中I是单位矩阵,D是由所有结点的度数构成的对角阵，图中结点的度数定义为与该结点相连的所有边的权值之和，即

步骤3.计算网络的输出权值矩阵初始值

参照附图2，该步骤具体子步骤如下：

3.1)随机产生网络的输入权值w_i,与网络隐藏层偏值b_i，表示隐藏层神经元个数；

3.2)利用sigmoid函数作为极速学习网络的隐藏层神经元激活函数，计算出隐藏层输出矩阵：

$H_l(w_1,...,w_{\tilde{N}},b_1,...,b_{\tilde{N}},x_1,...,x_{\tilde{N}})={\left[\begin{array}{ccc}\mathrm{\φ}(w_1\·x_1+b_1)& \·\·\·& \mathrm{\φ}(w_{\tilde{N}}\·x_1+b_{\tilde{N}})\\ \·& & \·\\ \·& \·\·\·& \·\\ \·& & \·\\ \mathrm{\φ}(w_1\·x_l+b_1)& \·\·\·& \mathrm{\φ}(w_{\tilde{N}}\·x_l+b_{\tilde{N}})\end{array}\right]}_{l\×\tilde{N}}$

其中,w_i·x_j表示w_i与x_j的内积；

3.3)则输出权值矩阵为

Y_l是标记样本标签，表示H_l的Moore-Penrose广义逆运算。

步骤4.计算样本属于每类的概率P

$p_k\left(x_j\right)=\frac{\mathrm{\λ}\underset{t=1,t\≠j}{\overset{N}{\mathrm{\Σ}}}{p}_k\left(x_t\right)w_{\mathrm{tj}}}{{\left|\right|f\left(x_j\right)-s_k\left|\right|}^2+\mathrm{\λ}\underset{s=1,s\≠j}{\overset{N}{\mathrm{\Σ}}}{w}_{\mathrm{js}}},(k=1,...m;j=l+1,...,N)$

步骤5.计算新的输出权值矩阵β

$\mathrm{\β}=(Y_l{H_l}^T+S\hat{P}{J}^T{H}_u^T){(H_l{H}_l^T+H_{u}J\hat{P}{J}^T{H}_u^T)}^{-1}$

H_l是标记样本的隐藏层输出矩阵，I_u是u×u的单位阵，其中每个S_k是一个m×u的矩阵，该矩阵的k行全为1，其余行全为0。是对角矩阵，对角元素为P按行排列时位置值的平方。

步骤6.判断目标函数变化值是否小于ε

目标函数为

$\mathrm{tr}(({\mathrm{\βH}}_l-Y_l){({\mathrm{\βH}}_l-Y_l)}^T+({\mathrm{\βH}}_{u}J-S)\hat{P}{({\mathrm{\βH}}_{u}J-S)}^T+{\mathrm{\λPLP}}^T)$

如果相邻目标函数的差值大于ε，则从步骤4开始重复执行步骤4到步骤6；如果小于ε，则往下执行步骤7。

步骤7.预测测试样本标签

按照下式计算

$Y_t=H_t\widehat{\mathrm{\β}}$

其中 $Y_t={\left[\begin{array}{c}{y_1}^T\\ \·\\ \·\\ \·\\ {y_t}^T\end{array}\right]}_{t\×m}$ 表示对于测试样本的网络输出矩阵，y_i＝[y_i1,y_i2,...,y_im]^T表示每一个样本的输出向量，表示优化后的输出权值矩阵，H_t表示测试样本的隐藏层输出矩阵。

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

1)实验条件

本实验采用高光谱地物分类Indiana作为实验数据，采用软件MATLAB R2012a作为仿真工具，计算机配置为Intel Core i3/2.27G/2G。

2)实验内容

本实验是对图3所示的Indiana数据进行分类，Indiana数据共21025个样本，包含背景和16类地物，地物共10366个。实验中先对样本进行滑窗预处理，窗口大小选择5*5。对该数据进行16类分类，每类样本标记2～10个，观察其正确率的变化。随机选择40％的样本作为训练样本，剩下的为测试样本。实验中窗口半径r和终止参数ε分别设置为3和10^-3。分别与SELM算法(Semi-supervised Extreme Learning Machine)、SSDR算法(Semi-SupervisedDimensionality Reduction)和SSDRsp算法(Semi-Supervised DimensionalityReduction based on Sparse representation)对比。在SSDR和SSDRsp中，均采用NN分类器。SELM以及本文方法的隐藏层节点数均设为400。在相同条件下，重复进行30次实验，并记录其平均结果以及每种标记下的最优结果。

图4是在每类标记的样本数从2个增加到10个时，不同方法分类正确率的变化图；图5是在每类标记10个样本时，不同方法分类结果的标记图；表1是不同方法在标记10个样本时分类后各种指标数据。

3)实验结果分析

图4给出了随着每类标记样本数的增加，4种方法的测试样本整体分类精度(TestOA)变化曲线图。数据取自独立运行30次后的平均结果。显而易见，随着已标记样本数的增加，四种方法的分类正确率均呈增长趋势。但是我们的方法明显优于其他三种方法，SELM算法次之，SSDRsp算法略优于SSDR算法。在每类仅标记2个样本时我们的方法就能达到很不错的分类结果。

图5给出了在每类标记10个样本时，不同方法的分类结果标记图。也给出标准标记图做参考。从图中可以很容易看出本文方法的优势。

表1从不同衡量指标的数据结果上对四种方法的分类能力进行比较，实验中每类标记10个样本。所使用的衡量指标包括：训练总体精度(train Overall Accuracy,trainOA)，测试平均精度(test Average Accuracy,testAA),测试总体精度(test OverallAccuracy,testOA),测试总体精度的方差(Standard of test Overall Accuracy,Std oftestOA)与Kappa系数(Kappa)。从表中可以看出，本方法不仅在分类正确率上优于其他方法，而且在泛化能力上，以及算法的稳定性上都有很好的表现。

表1：不同方法在每类标记10个样本时分类后各种指标数据

以上例举仅仅是对本发明的举例说明，并不构成对本发明的保护范围的限制，凡是与本发明相同或相似的设计均属于本发明的保护范围之内。

Claims (3)

1.一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，其特征在于：

(1)预处理训练样本

用滑窗技术对所有样本进行预处理，窗口大小为5*5，随机选择40％的样本做为训练样本X，剩余的做测试样本X_t，其中X＝{x₁,...,x_l,x_l+1...,x_l+u}，x₁,...,x_l为标记样本，x_l+1...,x_l+u为未标记样本，N＝l+u为总的训练样本数，每个样本x＝[x₁,x₂,…,x_d]，x_i表示每个样本的第i个特征值，其中i＝1,2,…d，d代表样本的特征维数；Y_l＝[y₁,y₂,...,y_l]表示与训练样本中的标记样本x₁,...,x_l对应的标签，其中y_i＝[y_i1,y_i2,...,y_im]^T表示与第i个标记样本对应的标签向量，m是类别数，T表示转置符号；

(2)计算训练样本的图Laplacian矩阵；

(3)计算网络的输出权值矩阵初始值；

(4)计算样本属于每类的概率P；

$p_k\left(x_j\right)=\frac{\mathrm{\λ}\underset{t=1,t\≠j}{\overset{N}{\Σ}}{p}_k\left(x_t\right)w_{tj}}{\left|\right|f\left(x_j\right)-s_k|{|}^2+\mathrm{\λ}\underset{s=1,s\≠j}{\overset{N}{\Σ}}{w}_{js}},(k=1,\mathrm{...},m;j=l+1,\mathrm{...},N)$

其中p_k(x_j),k＝1,...,m；j＝l+1,...,N表示第j个未标记样本x_j属于第k类的概率，m为类别数目，N＝l+u为总的训练样本数，λ是正则参数，w_js代表像素x_j与像素x_s之间的相似度，f(x_j)表示样本x_j的预测值，是第k类的标准标签，即该向量中第k个元素为1，其余的元素为0；

(5)计算新的输出权值矩阵β；

$\mathrm{\β}=(Y_l{H_l}^T+S\hat{P}{J}^T{H}_u^T){(H_l{H}_l^T+H_{u}J\hat{P}{J}^T{H}_u^T)}^{-1}$

H_l是标记样本的隐藏层输出矩阵，I_u是u×u的单位阵，其中每个S_k是一个m×u的矩阵，该矩阵的k行全为1，其余行全为0；是未标记样本属于各个类别的概率值，是对角矩阵，对角元素为P按行排列时位置值的平方，H_u是未标记样本的隐藏层输出矩阵；

(6)判断目标函数变化是否小于ε；

目标函数为

$tr\left(\right({\mathrm{\βH}}_l-Y_l){\left({\mathrm{\βH}}_l-Y_l\right)}^T+({\mathrm{\βH}}_{u}J-S)\hat{P}{\left({\mathrm{\βH}}_{u}J-S\right)}^T+{\mathrm{\λPLP}}^T)$

如果相邻目标函数的差值大于ε，则从步骤4开始重复执行步骤4到步骤6；如果小于ε，则往下执行步骤7；

(7)预测测试样本标签；

按照下式计算

$Y_t=H_t\widehat{\mathrm{\β}}$

其中表示对于测试样本的网络输出矩阵，y_i＝[y_i1,y_i2,...,y_im]^T表示每一个样本的输出向量，表示优化后的输出权值矩阵，H_t表示测试样本的隐藏层输出矩阵。

2.如权利要求1所述的一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，其特征在于：所述步骤(2)的图Laplacian矩阵的计算步骤具体如下：

2.1)首先计算图的权重矩阵W：

$W=\left\{\begin{array}{cc}{w}_{ij}& ife(i,j)\∈E\\ 0& otherwise\end{array}\right.$

其中E代表边集，表示样本间的某种相似性；边e(i,j)的权重为w(e)＝w_ij，表示顶点i与顶点j，即样本i与样本j之间的相似度，并且w_ij＝w_ji，边的权重计算公式如下：

$w_{ij}=\left\{\begin{array}{cc}1\sqrt{d_{ij}}& x_i\∈NB\left(x_j\right)or{x}_j\∈NB\left(x_i\right)\\ 0& otherwise\end{array}\right.$

其中NB(x_i)是以x_i为中心的一个r*r的方形窗中的样本集，d_ij代表样本x_i与样本x_j之间的欧氏距离；

2.2)再计算出图Laplacian矩阵L：

L＝I-D^-1/2WD^-1/2

其中I是单位矩阵,D是由所有结点的度数构成的对角阵，图中结点的度数定义为与该结点相连的所有边的权值之和，即

3.如权利要求1所述的一种基于修正聚类假设和半监督极速学习机的高光谱地物分类方法，其特征在于：所述步骤(3)的网络输出权值矩阵初始值的计算步骤具体如下：

3.1)随机产生网络的输入权值w_i,与网络隐藏层偏值b_i， 表示隐藏层神经元个数；

3.2)利用sigmoid函数作为极速学习网络的隐藏层神经元激活函数，计算出隐藏层输出矩阵：

$H_l(w_1,...,w_{\tilde{N}},b_1,...,b_{\tilde{N}},x_1,...,x_{\tilde{N}})={\left[\begin{array}{ccc}\mathrm{\φ}(w_1\·x_1+b_1)& \mathrm{...}& \mathrm{\φ}(w_{\tilde{N}}\·x_1+b_{\tilde{N}})\\ .& & .\\ .& \mathrm{...}& .\\ .& & .\\ \mathrm{\φ}(w_1\·x_l+b_1)& \mathrm{...}& \mathrm{\φ}(w_{\tilde{N}}\·x_l+b_{\tilde{N}})\end{array}\right]}_{l\×\tilde{N}}$

其中,w_i·x_j表示w_i与x_j的内积；

3.3)则输出权值矩阵为

Y_l是标记样本标签，表示H_l的Moore-Penrose广义逆运算。