CN113139513B - 基于超像素轮廓和改进pso-elm的空谱主动学习高光谱分类方法 - Google Patents

Abstract

本发明公开一种基于超像素轮廓和改进PSO‑ELM的空谱主动学习高光谱分类方法，属于高光谱图像处理与分类技术领域，该方法利用主成分分析法提取高光谱影像的光谱信息，并通过构造扩展超像素轮廓ESPP来整合主成分图像的光谱和空间信息，再通过构建的ESPP得到的主成分图像空谱特征获取超像素数据，使其适用于主动学习；接着提出了一种改进的极限学习机模型；最后将主成分图像空谱特征的超像素数据输入模型中训练得到最佳模型，将高光谱图像输入训练完的最佳模型中得到分类结果。本发明充分提取了高光谱图像的空谱信息，同时结合改进的粒子群‑极限学习机对输入权值和隐藏偏差进行优化，提高了泛化能力，最终提高了图像的分类精度。

Description

基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类 方法

技术领域

本发明涉及高光谱图像处理与分类技术领域，尤其涉及一种基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法。

背景技术

高光谱遥感图像含有几十或者几百个波段，每个波段都记录着丰富的地物反射能量信息。因此，近年来遥感领域研究者对高光谱遥感图像的关注越来越多，高光谱图像分类的研究也成为热门。

—幅遥感图像中总是包含大量的像素，已标记的像素点只能占极少一部分，存在大量的未标记像素点，要获得这些未标记点的标记需要耗费相当多的人力、物力，同时也会浪费宝贵的时间。可是只在少量的已标记样本上进行学习，很难得到一个泛化性能较好的模型，所以标记未标记的点是一个必要的工作。为了在尽量小的代价下，最大限度提高分类器性能，遥感学者将机器学习领域的主动学习方法引入遥感领域。高光谱图像分类是高光谱遥感影像应用的主要实现技术之一。不同于自然图像的分类，高光谱图像分类的对象是对图像中的每一个像元进行的，目的是给所有像元都预测为一个指定的标签。研究和实践表明，高光谱图像分类实际上跟当前机器学习、计算机视觉和模式识别理论与方法具有紧密连接，在传统的高光谱图像分类中，其研究内容主要集中在两个方面：特征提取以及分类器的设计。

在特征提取方面上，经历了从光谱特征提取到空间-光谱特征提取的阶段。光谱特征提取又根据是否改变原始光谱空间，分为特征提取和波段选择。此外，特征挖掘的概念"也被提出，作为一种新型的重构特征的提取方式，但并没有得到广泛应用。特征提取大多是基于统计的方法，关键是如何构建合适的变换矩阵，以将光谱特征映射到另一个更具表征性的空间上。

在非监督特征提取中，使用最为广泛的有主成分分析(PCA)和独立成分分析(ICA)等。PCA是一个典型的特征降维方法，它是一个线性变换，通过保留数据中对方差贡献最大的特征来实现数据降维。由于高光谱图像所具备的空域特性，仅仅釆用基于光谱特征的提取方法未能考虑空间信息，使得其分类结果往往存在椒盐现象，这与实际的地表分布不相符，为了更好的解决这个问题，很多基于空间-光谱特征提取的方法被提出，并迅速成为当今高光谱图像特征提取的主要趋势。

利用超像素辅助来进行特征提取的方法也成为高光谱图像空谱特征提取的重要手段。在高光谱遥感图像中，由于地物的分布往往是成块存在的，这意味着相邻的像素属于同一类的概率很大。因此，先将图像聚成若干个超像素，使得每个超像素尽可能同质，然后充分利用这种空间局部信息来建立特征提取的方式得到了广泛的应用。现有的分水岭算法和均值漂移方法在处理高光谱图像时，边界吻合程度不高，速度比较慢，占用较大的内存，整体的分割效果也不是特别的理想。

发明内容

针对上述现有技术的不足，本发明提供一种基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法。

为解决上述技术问题，本发明所采取的技术方案是：一种基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，包括如下步骤：

步骤1：输入高光谱遥感影像原始数据集；

步骤2：对输入的高光谱遥感影像原始数据集进行归一化处理，得到归一化后的高光谱遥感影像数据集；

步骤3：利用主成分分析法提取高光谱影像的光谱信息，得到前F个主成分图像作为后续图像，具体方法为：

步骤3.1：计算归一化后高光谱遥感影像数据集的协方差矩阵，过程如下：

步骤3.1.1：将归一化后高光谱遥感影像数据集记为：

其中，m＝A×B为一个A行B列的高光谱影像的像素点个数，n表示高光谱影像的波段数，x_mn表示第n个波段的第m个像素点；n个波段也代表具有n个特征；

步骤3.1.2：计算第i个波段所有像素的均值，如下公式所示：

其中，1≤k≤m，1≤i≤n，为第i个波段所有像素的均值，x_ki表示第i个波段的第k个像素点值；

步骤3.1.3：计算高光谱影像数据集的协方差矩阵C，如下公式所示：

其中，T为矩阵转置运算，x_k表示第k个像素点。

步骤3.2：根据协方差矩阵，求得对应的特征值和特征向量，得到特征值大的前F个主成分图像，并将得到主成分图像作为后续使用图像，过程如下：

步骤3.2.1：对协方差矩阵C进行特征值分解，求得协方差矩阵的特征值λ₁≥λ₂≥λ₃≥λ_i...≥λ_n，及其对应的特征向量α＝[α₁，α₂，α₃，α_i，...α_n]，具体公式如下：

设A是n阶矩阵，如果特征值λ_i和α＝[α₁，α₂，α₃，α_i，...α_n]使关系式Aα＝λ_iα成立，即(C-λ_iE)＝0，那么这样的特征值λ_i称为矩阵A的特征值，α＝[α₁，α₂，α₃，α_i，...α_n]称为A的对应于特征值λ_i的特征向量；其中，E代表单位矩阵；

步骤3.2.2：通过特征值λ₁对应的特征向量α₁得到的第一主成分分量，将得到的第一主成分分量与原始高光谱数据集X相乘，得到第一主成分图像，并将其作为后续图像，同理λ₂对应的特征向量α₂得到的第二主成分分量，λ₃对应的特征向量α₃得到的第三主成分分量，以此类推，如下公式所示：I＝{I₁，I₂，I₃，I_i，...，I_n}

其中，I为图像集合，X为输入的原始数据集图像。

步骤4：通过构造扩展超像素轮廓ESPP来整合主成分图像的光谱和空间信息，过程如下：

步骤4.1：图像I采用简单线性迭代聚类SLIC(simple linear iterativeclustering)算法，根据阈值将图像集合I中的每一个图像划分为指定数量的超像素，并生成平均超像素图像MSPI，具体步骤如下：

步骤4.1.1：SLIC同时考虑像素的颜色和空间属性，在五维空间labxy下实现，其中(l，a，b)是CIELAB颜色空间下的像素点颜色向量，(x，y)是图像平面像素点的位置向量；默认情况下，SLIC算法只有一个参数k_c决定分割结果中超像素的数量；

对于一个像素点数为N_c的图像I，首先从N_c个像素点中以为格网区间间隔，选取k_c个聚类中心C_ic＝[l_ic，a_ic，b_ic，x_ic，y_ic]^T，并将这些中心点移动到其3×3领域内梯度值最小的位置；其次，计算距离，确定每个像素点的最近邻聚类中心，具体公式如下：

其中，d_lab为颜色空间距离，d_xy为平面空间距离，D_s为五维空间距离；m_c为权重；S是格网区间间隔，kc是超像素编号，ic是像素点编号；

SLIC设置搜索范围为聚类中心周围2S×2S邻域范围内的每个像素点分别与当前聚类中心五维向量求距离度量，并将此度量与初始值比较，若比初始值小，则更新当前的度量，并将当前的像素点的标记更新为当前所属的聚类中心；整幅图像扫描之后，每个像素点对应一个聚类中心，相同标记的像素点属于同一个类别；接下来进行迭代，对上一次划分的聚类，求出新的聚类中心，重复上述过程，直到迭代次数达到上限或者分类结果基本不发生改变即划分完成，迭代结束；

步骤4.1.2：划分之后的超像素数量，创建一个有N个超像素的平均超像素图像MSPI，其中每个超像素被分配其覆盖的所有超像素的平均灰度值，公式如下：

其中，每个超像素被分配其覆盖的所有像素的平均灰度值，GV(SP_ω)是第ω个超像素(SP_ω)的灰度值，并且x_i，i＝1，2，...，N_ω是具有N_ω个像素数的SP_ω中第i个像素的灰度值；

步骤4.2：考虑多个阈值，将得到多个MSPI结果，每个结果具有不同的邻域信息；因此，使用递增阈值序列将获得的多个MSPI结果与原始图像叠加，来创建一个数据结构，并将其命名为超像素轮廓SPP(Superpixel Profile)，公式如下：

其中，超像素轮廓SPP表示具有图像I的原始光谱信息及其在邻域边界的空间信息的图像测量；

步骤4.3：将多个不同的主成分的超像素轮廓SPP进行连接构建成扩展超像素轮廓ESPP(Extended Superpixel Profile)，公式如下：

ESPP(H)＝{SPP(PC₁)，SPP(PC₂)，......，SPP(PC_l)}

其中，ESPP(H)为超像素中像素点的光谱特征。

步骤5：根据步骤4构建的ESPP得到的主成分图像空谱特征获取超像素数据D_AL，使其适用于主动学习；

所述获取超像素数据D_AL的过程如下：

其中，n_AL为超像素中像素点个数，ESPP(H)为超像素中像素点的光谱特征，D_AL为适用于主动学习的超像素数据。

步骤6：将自适应粒子群算法与极限学习机ELM相结合，提出了一种改进的极限学习机模型，这种改进的极限学习机ELM模型采用自适应粒子群算法来选择输入权值，以提高泛化性能和单隐层前馈神经网络SLFN的调节效果，具体过程如下：

步骤6.1对于N个任意不同的样本，符号式表达(x_i，t_i)，x_i表示第i个数据实例，t_i表示第i个数据实例对应的标记，集合代指所有训练数据。单隐层前馈神经网络(SLFN)具有H个隐藏层和激活函数g(·)可以零误差逼近这N个样本。

x_i＝[x_i1，x_i2，...，x_in]^T∈Rⁿ，t_i＝[t_i1，t_i2，...，t_im]^T∈R^m

Hwo＝T

其中，wh_i＝[wh_i1，wh_i2，...，wh_in]^T是连接第i个隐藏神经元和输入神经元的权值向量，b_i是第i个隐藏神经元的偏置，wo_i＝[wo_i1，wo_i2，...，wo_im]^T是连接第i个隐藏神经元和输出神经元的权值向量，m为输出层的结点个数。

步骤6.2：求最小化平方差作为评价训练误差(目标函数)，使得该目标函数最小的解就是优化解。即通过最小化近似平方差的方法对连接隐藏层和输出层的权重wo进行求解，目标函数如下：

min||Hwo-T||²，wo∈R^H×m

其中，H是隐藏层和输出层的输出矩阵，T是训练数据的目标矩阵。求得最小二乘解：

wo＝H⁺T

其中，H⁺为矩阵H的Moore-Penrose广义逆矩阵。

步骤6.3：采用奇异值分解SVD计算广义逆矩阵，H的奇异值分解公式如下：

其中

任意矩阵S的广义逆，由于U，V均为正交矩阵，即有：

U^TU＝UU^T＝E，V^TV＝VV^T＝E，U^-1＝U^T，V^-1＝V^T

所以，由此求得任意矩阵S的Moore-Penrose广义逆为

G⁺＝VS⁺U^T

其中，S⁺是对角阵，对角元素为

步骤6.4：首先，群是随机生成的。群中的每个粒子都由一组输入权重和隐藏偏差组成，表达式如下：

P_i＝[wh₁₁，wh₁₂，...，wh_1n，wh₂₁，wh₂₂，...，wh_2n，wh_H1，wh_H2，...，wh_Hn，b₁，b₂，...，b_H]

其中，P_i为第i个粒子，粒子中的所有元素取值范围都在[-1，1]。

步骤6.5：其次，根据步骤6.2计算每个粒子对应的输出权值，然后对每个粒子的适应度进行评估。为了避免SLFN过拟合，我们仅采用每个粒子的适应度作为验证集的均方根误差RMSE，公式如下：

其中，n_v是验证集的数量，t_j表示隐藏层第j个隐藏单元对应的网络输出。

步骤6.6：接下来，通过所有粒子的适应度，计算所有粒子的最佳位置(P_b)和相邻粒子的最佳位置(P_g)。神经网络的权值越小，泛化性能越好。因此，为了进一步提高泛化性能，考虑输出权值的范数和验证集上的(RMSE)来确定所有粒子的最佳位置(P_b)和群的最佳位置(P_g)。具体描述如下：

其中，η为容忍度，f(P_i)，f(P_ib)和f(P_g)分别代表第i个粒子的适应度值，第i个粒子的最佳位置的适应度和所有粒子的全局最佳位置的适应度。和/>分别代表当输入权值设为第i个粒子，第i个粒子的最佳位置和所有粒子全局最佳位置时，MP广义逆得到的相应的输出权值。

步骤6.7：每个粒子采用自适应粒子群优化算法产生新的种群。粒子中的所有元素取值范围都在[-1，1]。该先验信息通过以下公式在改进粒子群算法中：

步骤6.8：最后，重复上述优化过程，直到达到目标或完成最大优化时期。从而得到具有最优输入权值和隐藏偏差的扩展学习模型，并将最优扩展学习模型应用于测试数据。

所述PSO的工作原理是在搜索空间中随机初始化一群鸟，其中每只鸟都被称为“粒子”。这些粒子以一定的速度飞行，经过多次迭代后得到全局最优位置。在每次迭代中，每个粒子根据自身的动量以及最佳位置(P_b)和相邻粒子的最佳位置(P_g)的影响来调整速度矢量，从而得到粒子飞行的新位置。假设搜索空间的维数为D，粒子的总数为n，第i个粒子的位置可表示为向量X_i＝(x_i1，x_i2，...，x_iD)。到目前为止，第i个粒子搜索的最佳位置为P_ib(p_i1，p_i2，...p_iD)。所有搜索粒子的最佳位置为向量P_g(p_g1p_g2，...p_gD)。第i个粒子的速度表示为V_i(v_i1，v_i2，...v_iD)。公式如下：

v_id(t+1)＝v_id(t)+c₁·rand()·[p_id(t)-x_id(t)]+c₂·rand()·[p_gd(t)-x_id(t)]

x_id(t+1)＝x_id(t)+v_id(t+1)1≤i≤n，1≤d≤D

其中，υ_id为粒子速度，p_id为粒子历史最佳位置，x_id为粒子目前位置，p_gd为整个粒子群体历史最佳位置。c₁，c₂是学习因子，通常为常数2。rand()是一个介于0和1之间的随机数。原始算法的实施还需要对速度(v_max)进行限制。

所述自适应粒子群优化算法APSO描述如下：

v_id(t+1)w·v_id(t)+c₁·rand()·[p_id(t)-x_id(t)]+c₂·rand()·[p_gd(t)-x_id(t)]

其中，w叫做惯性因子，其值为非负。目前采用较多的是线性递减权值(LinearlyDecreasing Weight，LDW)策略。公式如下：

w(t)＝w_max-t(w_max-w_min)/iter max

其中，w_maxw_min，itermax依次代表初始惯性权值、最终惯性权值和最大搜索代数。

步骤7：从步骤5中得到的主成分图像空谱特征的超像素数据D_AL中随机选取少数样本并分配类别标签对步骤6中的改进ELM模型进行训练，当前训练样本的标签保留在标签池L中；

步骤8：将已选择的样本从所有样本集中移除，其余的样本保存在一个未标记的样本池U中；在每次迭代中，利用基于不确定性、多样性组合条件设计的查询函数对样本的信息量进行评估，根据样本信息量的大小，从U中选择一批信息量大的样本分配它们的类别标签并被附加到标签池L中来更新L；

步骤9：使用更新后的L对改进的ELM模型进行再训练，直到满足迭代停止条件；

步骤10：使用训练后的模型对超像素进行数据分类，得到每个类别，并将每个类别分别赋值给本类别所包含的每个像素点，达到了高光谱遥感图像分类的目的。

采用上述技术方案所产生的有益效果在于：本发明提供的方法采用的超像素轮廓和改进PSO-ELM的空谱主动学习高光谱图像分类方法，利用超像素分割提取图像的空谱信息，更好的将图像特征信息表示出来；主动学习技术能够在训练样本更少的情况下，提升分类性能。同时结合改进的粒子群-极限学习机对输入权值和隐藏偏差进行优化，采用最小范数和最小二乘法对输出权值进行解析确定。改进粒子群算法在选择输入权值和隐藏偏差的过程中，不仅考虑了验证集上的均方根误差，还考虑了输出权值的范数，并将输入权值和隐藏偏差限制在合理的范围内。提高了泛化能力，最终提高了图像的分类精度。

附图说明

图1为本发明实施例中基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法的流程图；

图2为本发明实施例中提供的SLIC算法聚类中心初始化示意图；其中，(a)为初始化聚类中心示意图；(b)为聚类中心调整示意图；

图3为本发明实施例中提供的高光谱图像主成分MSPI；其中，(a)为初始图像；(b)为超像素分割示意图；(c)为平均超像素图像；

图4为本发明实施例中提供的实验选取的高光谱影像数据集的地面参考数据；

图5为本发明实施例中提供的高光谱影像数据集的伪彩色影像；

图6为本发明实施例中提供的采用四种不同方法进行实验的分类结果图；其中，(a)为EQB方法分类结果图；(b)为MS-ABD方法分类结果图；(c)为MCLU-ECBD方法分类结果图；(d)为采用本发明方法分类结果图。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。

如图1所示，本实施例中基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法如下所述：

步骤1：输入高光谱遥感影像原始数据集；

本实施例中，实验数据集是由ROSIS-03(反射光学系统成像光谱仪)光学传感器在意大利Pavia大学市区采集的高光谱图像，如图4所示。该图像的空间分辨率为1.3m，大小为610×340像素，103个光谱波段，包含9类。图5为其伪彩色图像。

步骤2：对输入的高光谱遥感影像原始数据集进行归一化处理，得到归一化后的高光谱遥感影像数据集；

步骤3：利用主成分分析法提取高光谱影像的光谱信息，得到前F个主成分图像作为后续图像，具体方法为：

步骤3.1：计算归一化后高光谱遥感影像数据集的协方差矩阵，过程如下：

步骤3.1.1：将归一化后高光谱遥感影像数据集记为：

其中，m＝A×B为一个A行B列的高光谱影像的像素点个数，n表示高光谱影像的波段数，x_mn表示第n个波段的第m个像素点；n个波段也代表具有n个特征；

步骤3.1.2：计算第i个波段所有像素的均值，如下公式所示：

其中，1≤k≤m，1≤i≤n，为第i个波段所有像素的均值，x_ki表示第i个波段的第k个像素点值；

步骤3.1.3：计算高光谱影像数据集的协方差矩阵C，如下公式所示：

其中，T为矩阵转置运算，x_k表示第k个像素点。

步骤3.2：根据协方差矩阵，求得对应的特征值和特征向量，得到特征值大的前F个主成分图像，并将得到主成分图像作为后续使用图像，过程如下：

步骤3.2.1：对协方差矩阵C进行特征值分解，求得协方差矩阵的特征值λ₁≥λ₂≥λ₃≥λ_i...≥λ_n，及其对应的特征向量α＝[α₁，α₂，α₃，α_i，...α_n]，具体公式如下：

设A是n阶矩阵，如果特征值λ_i和α＝[α₁，α₂，α₃，α_i，...α_n]使关系式Aα＝λ_iα成立，即(C-λ_iE)＝0，那么这样的特征值λ_i称为矩阵A的特征值，α＝[α₁，α₂，α₃，α_i，...α_n]称为A的对应于特征值λ_i的特征向量；其中，E代表单位矩阵；

步骤3.2.2：通过特征值λ₁对应的特征向量α₁得到的第一主成分分量，将得到的第一主成分分量与原始高光谱数据集X相乘，得到第一主成分图像，并将其作为后续图像，同理λ₂对应的特征向量α₂得到的第二主成分分量，λ₃对应的特征向量α₃得到的第三主成分分量，以此类推，如下公式所示：

其中，I为图像集合，X为输入的原始数据集图像。

本实施例中，特征值大的前5个主成分图像作为后续使用图像。

步骤4：通过构造扩展超像素轮廓ESPP来整合主成分图像的光谱和空间信息，过程如下：

步骤4.1：图像I采用简单线性迭代聚类SLIC(simple linear iterativeclustering)算法，根据阈值将图像集合I中的每一个图像划分为指定数量的超像素，并生成平均超像素图像MSPI，其中阈值集合：{50，100，500，1000，1500}，具体步骤如下：

步骤4.1.1：SLIC同时考虑像素的颜色和空间属性，在五维空间labxy下实现，其中(l，a，b)是CIELAB颜色空间下的像素点颜色向量，(x，y)是图像平面像素点的位置向量；默认情况下，SLIC算法只有一个参数k_c决定分割结果中超像素的数量；

对于一个像素点数为N_c的图像I，首先从N_c个像素点中以为格网区间间隔，选取k_c个聚类中心C_ic＝[l_ic，a_ic，b_ic，x_ic，y_ic]^T，并将这些中心点移动到其3×3领域内梯度值最小的位置；其次，计算距离，确定每个像素点的最近邻聚类中心，具体公式如下：

其中，d_lab为颜色空间距离，d_xy为平面空间距离，D_s为五维空间距离；m_c为权重；S是格网区间间隔，kc是超像素编号，ic是像素点编号；

SLIC设置搜索范围为聚类中心周围2S×2S邻域，如图2中的(a)所示，范围内的每个像素点分别与当前聚类中心五维向量求距离度量，并将此度量与初始值比较，若比初始值小，则更新当前的度量，并将当前的像素点的标记更新为当前所属的聚类中心，如图2中的(b)所示；整幅图像扫描之后，每个像素点对应一个聚类中心，相同标记的像素点属于同一个类别；接下来进行迭代，对上一次划分的聚类，求出新的聚类中心，重复上述过程，直到迭代次数达到上限或者分类结果基本不发生改变即划分完成，迭代结束；

步骤4.1.2：划分之后的超像素数量，以500个超像素为例，如图3中的(b)所示，图3中的(a)为原始数据集。创建一个有N个超像素的平均超像素图像MSPI，如图3中的(c)所示，其中每个超像素被分配其覆盖的所有超像素的平均灰度值，公式如下：

其中，每个超像素被分配其覆盖的所有像素的平均灰度值，GV(SP_ω)是第ω个超像素(SPω)的灰度值，并且x_i，i＝1，2，...，N_ω是具有N_ω个像素数的SP_ω中第i个像素的灰度值；

步骤4.2：考虑多个阈值，将得到多个MSPI结果，每个结果具有不同的邻域信息；因此，使用递增阈值序列将获得的多个MSPI结果与原始图像叠加，来创建一个数据结构，并将其命名为超像素轮廓SPP(Superpixel Profile)，公式如下：

其中，超像素轮廓SPP表示具有图像I的原始光谱信息及其在邻域边界的空间信息的图像测量；

步骤4.3：将多个不同的主成分的超像素轮廓SPP进行连接构建成扩展超像素轮廓ESPP(Extended Superpixel Profile)，公式如下：

ESPP(H)＝{SPP(PC₁)，SPP(PC₂)，......，SPP(PC_l)}

其中，ESPP(H)为超像素中像素点的光谱特征。

步骤5：根据步骤4构建的ESPP得到的主成分图像空谱特征获取超像素数据D_AL，使其适用于主动学习；

所述获取超像素数据D_AL的过程如下：

其中，n_AL为超像素中像素点个数，ESPP(H)为超像素中像素点的光谱特征，D_AL为适用于主动学习的超像素数据。

PSO的工作原理是在搜索空间中随机初始化一群鸟，其中每只鸟都被称为“粒子”。这些粒子以一定的速度飞行，经过多次迭代后得到全局最优位置。在每次迭代中，每个粒子根据自身的动量以及最佳位置(P_b)和相邻粒子的最佳位置(P_g)的影响来调整速度矢量，从而得到粒子飞行的新位置。假设搜索空间的维数为D，粒子的总数为n，第i个粒子的位置可表示为向量X_i＝(x_i1，x_i2，...，x_iD)。到目前为止，第i个粒子搜索的最佳位置为P_ib(p_i1，p_i2，...p_iD)。所有搜索粒子的最佳位置为向量P_g(p_g1，P_g2，...p_gD)。第i个粒子的速度表示为V_i(v_i1，v_i2，...v_iD)。公式如下：

v_id(t+1)＝v_id(t)+c₁·rand()·[p_id(t)-x_id(t)]+c₂·rand()·[p_gd(t)-x_id(t)]

x_id(t+1)＝x_id(t)+v_id(t+1)1≤i≤n，1≤d≤D

其中，υ_id为粒子速度，p_id为粒子历史最佳位置，x_id为粒子目前位置，P_gd为整个粒子群体历史最佳位置。c₁，c₂是学习因子，通常为常数2。rand()是一个介于0和1之间的随机数。原始算法的实施还需要对速度(v_max)进行限制。

自适应粒子群优化算法(APSO)描述如下：

v_id(t+1)w·v_id(t)+c₁·rand()·[p_id(t)-x_id(t)]+c₂·rand()·[p_gd(t)-x_id(t)]

其中，w叫做惯性因子，其值为非负。目前采用较多的是线性递减权值(LinearlyDecreasing Weight，LDW)策略。公式如下：

w(t)＝w_max-t(w_max-w_min)/itermax

其中，w_max，w_min，itermax依次代表初始惯性权值、最终惯性权值和最大搜索代数。

步骤6：将自适应粒子群算法与极限学习机ELM相结合，提出了一种改进的极限学习机模型，这种改进的极限学习机ELM模型采用自适应粒子群算法来选择输入权值，以提高泛化性能和单隐层前馈神经网络SLFN的调节效果，具体过程如下：

步骤6.1对于N个任意不同的样本，符号式表达(x_i，t_i)，x_i表示第i个数据实例，t_i表示第i个数据实例对应的标记，集合代指所有训练数据。单隐层前馈神经网络(SLFN)具有H个隐藏层和激活函数g(·)可以零误差逼近这N个样本。

x_i＝[x_i1，x_i2，...，x_in]^T∈Rⁿ，t_i＝[t_i1，t_i2，...，t_in]^T∈R^m

Hwo＝T

其中，wh_i＝[wh_i1，wh_i2，...，wh_in]^T是连接第i个隐藏神经元和输入神经元的权值向量，b_i是第i个隐藏神经元的偏置，wo_i＝[wo_i1，wo_i2，...，wo_im]^T是连接第i个隐藏神经元和输出神经元的权值向量，m为输出层的结点个数。

步骤6.2：求最小化平方差作为评价训练误差(目标函数)，使得该目标函数最小的解就是优化解。即通过最小化近似平方差的方法对连接隐藏层和输出层的权重wo进行求解，目标函数如下：

min||Hwo-T||²，wo∈R^H×m

其中，H是隐藏层和输出层的输出矩阵，T是训练数据的目标矩阵。求得最小二乘解：

wo＝H⁺T

其中，H⁺为矩阵H的Moore-Penrose广义逆矩阵。

步骤6.3：采用奇异值分解SVD计算广义逆矩阵，H的奇异值分解公式如下：

其中

任意矩阵S的广义逆，由于U，V均为正交矩阵，即有：

U^TU＝UU^T＝E，V^TV＝VV^T＝E，U^-1＝U^T，V^-1＝V^T

所以，由此求得任意矩阵s的Moore-Penrose广义逆为

G⁺＝VS⁺U^T

其中，S⁺是对角阵，对角元素为

步骤6.4：首先，群是随机生成的。群中的每个粒子都由一组输入权重和隐藏偏差组成，表达式如下：

P_i＝[wh₁₁，wh₁₂，...，wh_1n，wh₂₁，wh₂₂，...，wh_2n，wh_H1，wh_H2，...，wh_Hn，b₁，b₂，...，b_H]

其中，P_i为第i个粒子，粒子中的所有元素取值范围都在[-1，1]。

步骤6.5：其次，根据步骤6.2计算每个粒子对应的输出权值，然后对每个粒子的适应度进行评估。为了避免SLFN过拟合，我们仅采用每个粒子的适应度作为验证集的均方根误差RMSE，公式如下：

其中，n_v是验证集的数量，t_j表示隐藏层第j个隐藏单元对应的网络输出。

步骤6.6：接下来，通过所有粒子的适应度，计算所有粒子的最佳位置(P_b)和相邻粒子的最佳位置(P_g)。神经网络的权值越小，泛化性能越好。因此，为了进一步提高泛化性能，考虑输出权值的范数和验证集上的(RMSE)来确定所有粒子的最佳位置(P_b)和群的最佳位置(P_g)。具体描述如下：

其中，η为容忍度，f(P_i)，f(P_ib)和f(P_g)分别代表第i个粒子的适应度值，第i个粒子的最佳位置的适应度和所有粒子的全局最佳位置的适应度。和/>分别代表当输入权值设为第i个粒子，第i个粒子的最佳位置和所有粒子全局最佳位置时，MP广义逆得到的相应的输出权值。

步骤6.7：每个粒子采用自适应粒子群优化算法产生新的种群。粒子中的所有元素取值范围都在[-1，1]。该先验信息通过以下公式在改进粒子群算法中：

步骤6.8：最后，重复上述优化过程，直到达到目标或完成最大优化时期。从而得到具有最优输入权值和隐藏偏差的扩展学习模型，并将最优扩展学习模型应用于测试数据。

本实施例中，最大优化期为20，种群规模为200。在本仿真实验中均选取初始惯性权重w_max和最终惯性权重w_min分别为1.2和0.4。c₁＝c₂＝2。随着优化次数的增加，产生的时间也会增加。同时会影响最优解的寻找。导致分类精度的降低。

步骤7：从步骤5中得到的主成分图像空谱特征的超像素数据D_AL中随机选取少数样本并分配类别标签对步骤6中的改进ELM模型进行训练，当前训练样本的标签保留在标签池L中；

步骤8：将已选择的样本从所有样本集中移除，其余的样本保存在一个未标记的样本池U中；在每次迭代中，利用基于不确定性、多样性组合条件设计的查询函数对样本的信息量进行评估，根据样本信息量的大小，从U中选择一批信息量大的样本分配它们的类别标签并被附加到标签池L中来更新L；

步骤9：使用更新后的L对改进的ELM模型进行再训练，直到满足迭代停止条件；

步骤10：使用训练后的模型对超像素进行数据分类，得到每个类别，并将每个类别分别赋值给本类别所包含的每个像素点，达到了高光谱遥感图像分类的目的。

本实施例中，还给出了采用本发明的方法和其他不同方法对数据集进行分类的对比。表1显示了迭代后PaviaU数据集的各类地物分类精度、总体分类精度(OA)、平均分类精度(AA)、和kappa精度的实验结果。图6给出了分类结果图。本发明方法的分类结果OA，AA和Kappa系数分别约为：98.56％，97.76％和98.10％，相比EQB、MS-ABD和MCLU-ECBD方法的分类结果均有所提高，OA分别提高了约：7.01％、6.1％和2.15％，AA分别提高了约：6.54％、5.93％和2.47％，Kappa系数分别提高了约：9.4％、8.14％和2.87％，同时获取了9个最佳类别的地物分类精度。从表1中可以看出，使用空谱信息的技术比仅基于光谱值的方法表现更好。这显示了空间信息在用主动学习方法对高光谱图像分类中的重要性。此外，可以注意到，基于ESPP的方法比其他所有方法都要执行的更好。本发明的方法能够对大多数类别进行分类，这也表明了本申请的方法在对高光谱图像进行分类时的泛化能力与潜力。

表1四种方法分类结果(％)

传统的方法仅使用了光谱值，没有结合空间信息的提取，严重影响了高光谱遥感影像的分类。采用ESPP的方法，同时获取图像的空间和光谱信息。提高了分类精度。同时高光谱图像在分类器的选择以及其存在大量非线性参数都会导致运行反应时间较长，分类精度不准确。运用改进的粒子群结合极限向量机加快了学习速度和分类的泛化能力。

Claims (5)

1.一种基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，其特征在于，包括如下步骤：

步骤1：输入高光谱遥感影像原始数据集；

步骤2：对输入的高光谱遥感影像原始数据集进行归一化处理，得到归一化后的高光谱遥感影像数据集；

步骤3：利用主成分分析法提取高光谱影像的光谱信息，得到前F个主成分图像作为后续图像，具体方法为：

步骤3.1：计算归一化后高光谱遥感影像数据集的协方差矩阵；

步骤3.2：根据协方差矩阵，求得对应的特征值和特征向量，得到特征值大的前F个主成分图像，并将得到主成分图像作为后续使用图像；

步骤4：通过构造扩展超像素轮廓ESPP来整合主成分图像的光谱和空间信息；

步骤5：根据步骤4构建的ESPP得到的主成分图像空谱特征获取超像素数据D_AL，使其适用于主动学习；

步骤6：将自适应粒子群算法与极限学习机ELM相结合，提出了一种改进的极限学习机模型，这种改进的极限学习机ELM模型采用自适应粒子群算法来选择输入权值，以提高泛化性能和单隐层前馈神经网络SLFN的调节效果；

步骤7：从步骤5中得到的主成分图像空谱特征的超像素数据D_AL中随机选取少数样本并分配类别标签对步骤6中的改进ELM模型进行训练，当前训练样本的标签保留在标签池L中；

步骤8：将已选择的样本从所有样本集中移除，其余的样本保存在一个未标记的样本池U中；在每次迭代中，利用基于不确定性、多样性组合条件设计的查询函数对样本的信息量进行评估，根据样本信息量的大小，从U中选择一批信息量大的样本分配它们的类别标签并被附加到标签池L中来更新L；

步骤9：使用更新后的L对改进的ELM模型进行再训练，直到满足迭代停止条件；

步骤10：使用训练后的模型对超像素进行数据分类，得到每个类别，并将每个类别分别赋值给本类别所包含的每个像素点，达到了高光谱遥感图像分类的目的。

2.根据权利要求1所述的基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，其特征在于，所述步骤3.1的具体方法为：

步骤3.1.1：将归一化后高光谱遥感影像数据集记为：

其中，m＝A×B为一个A行B列的高光谱影像的像素点个数，n表示高光谱影像的波段数，x_mn表示第n个波段的第m个像素点；n个波段也代表具有n个特征；

步骤3.1.2：计算第i个波段所有像素的均值，如下公式所示：

其中，1≤k≤m，1≤i≤n，为第i个波段所有像素的均值，x_ki表示第i个波段的第k个像素点值；

步骤3.1.3：计算高光谱影像数据集的协方差矩阵C，如下公式所示：

其中，T为矩阵转置运算，x_k表示第k个像素点。

3.根据权利要求1所述的基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，其特征在于，所述步骤3.2的具体方法为：

步骤3.2.1：对协方差矩阵C进行特征值分解，求得协方差矩阵的特征值λ₁≥λ₂≥λ₃≥λ_i...≥λ_n，及其对应的特征向量α＝[α₁，α₂，α₃，α_i，...α_n]，具体公式如下：

设A是n阶矩阵，如果特征值λ_i和α＝[α₁，α₂，α₃，α_i，...α_n]使关系式Aα＝λ_iα成立，即(C-λ_iE)＝0，那么这样的特征值λ_i称为矩阵A的特征值，α＝[α₁，α₂，α₃，α_i，...α_n]称为A的对应于特征值λ_i的特征向量；其中，E代表单位矩阵；

步骤3.2.2：通过特征值λ₁对应的特征向量α₁得到的第一主成分分量，将得到的第一主成分分量与原始高光谱数据集X相乘，得到第一主成分图像，并将其作为后续图像，同理λ₂对应的特征向量α₂得到的第二主成分分量，λ₃对应的特征向量α₃得到的第三主成分分量，以此类推，如下公式所示：

其中，I为图像集合，X为输入的原始数据集图像。

4.根据权利要求1所述的基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，其特征在于，所述步骤4的具体方法为：

步骤4.1：图像I采用简单线性迭代聚类SLIC算法，根据阈值将图像集合I中的每一个图像划分为指定数量的超像素，并生成平均超像素图像MSPI，具体步骤如下：

步骤4.1.1：SLIC同时考虑像素的颜色和空间属性，在五维空间labxy下实现，其中(l，a，b)是CIELAB颜色空间下的像素点颜色向量，(x，y)是图像平面像素点的位置向量；默认情况下，SLIC算法只有一个参数k_c决定分割结果中超像素的数量；

对于一个像素点数为N_c的图像I，首先从N_c个像素点中以为格网区间间隔，选取k_c个聚类中心C_ic＝[l_ic，a_ic，b_ic，x_ic，y_ic]^T，并将这些中心点移动到其3×3领域内梯度值最小的位置；其次，计算距离，确定每个像素点的最近邻聚类中心，具体公式如下：

其中，d_lab为颜色空间距离，d_xy为平面空间距离，D_s为五维空间距离；m_c为权重；S是格网区间间隔，kc是超像素编号，ic是像素点编号；

SLIC设置搜索范围为聚类中心周围2S×2S邻域范围内的每个像素点分别与当前聚类中心五维向量求距离度量，并将此度量与初始值比较，若比初始值小，则更新当前的度量，并将当前的像素点的标记更新为当前所属的聚类中心；整幅图像扫描之后，每个像素点对应一个聚类中心，相同标记的像素点属于同一个类别；接下来进行迭代，对上一次划分的聚类，求出新的聚类中心，重复上述过程，直到迭代次数达到上限或者分类结果基本不发生改变即划分完成，迭代结束；

步骤4.1.2：划分之后的超像素数量，创建一个有N个超像素的平均超像素图像MSPI，其中每个超像素被分配其覆盖的所有超像素的平均灰度值，公式如下：

其中，每个超像素被分配其覆盖的所有像素的平均灰度值，GV(SP_ω)是第ω个超像素(SP_ω)的灰度值，并且x_i，i＝1，2，...，N_ω是具有N_ω个像素数的SP_ω中第i个像素的灰度值；

步骤4.2：考虑多个阈值，将得到多个MSPI结果，每个结果具有不同的邻域信息；因此，使用递增阈值序列将获得的多个MSPI结果与原始图像叠加，来创建一个数据结构，并将其命名为超像素轮廓SPP，公式如下：

其中，超像素轮廓SPP表示具有图像I的原始光谱信息及其在邻域边界的空间信息的图像测量；

步骤4.3：将多个不同的主成分的超像素轮廓SPP进行连接构建成扩展超像素轮廓ESPP，公式如下：

ESPP(H)＝{SPP(PC₁)，SPP(PC₂)，......，SPP(PC_l)}

其中，ESPP(H)为超像素中像素点的光谱特征。

5.根据权利要求1所述的基于超像素轮廓和改进PSO-ELM的空谱主动学习高光谱分类方法，其特征在于，所述获取超像素数据D_AL的过程如下：

其中，n_AL为超像素中像素点个数，ESPP(H)为超像素中像素点的光谱特征，D_AL为适用于主动学习的超像素数据。