基于模糊标签的高光谱数据子空间投影和分类方法

摘要

本发明公开了一种基于模糊标签的高光谱数据子空间投影和分类方法，主要解决高光谱图像中由于混合像元和噪声引起地物错分和数据判别性差的问题。其步骤为：1.将遥感数据库样本集划分为训练样本和标记样本集；2.计算由子空间投影后的标记样本集生成的判别项；3.构造由训练样本的模糊标签确定的拉普拉斯正则项；4.通过最大化判别项和正则项之差获取最优投影矩阵和模糊标签，以实现有效降维的同时，实现高精度的分类。本发明采用判别子空间投影的方法来构造判别项，将数据投影到低维空间，增强数据的判别性能，进而引入模糊标签来构造拉普拉斯正则，解决了混合像元带来的错分问题，在实现降维的同时，实现高精度的分类。

说明书

技术领域

本发明属于图像处理技术领域，更进一步涉及一种数据降维和分类方法，可 用于遥感影像数据的降维与分类。

背景技术

经上世纪的快速发展，高光谱遥感技术在理论技术和应用上发生了翻天覆地 的变化，广泛应用于农业、林业、国防侦察识别伪装等领域。但是，高光谱数据 处理的技术相对落后，制约了高光谱遥感技术的进一步推广。分类作为高光谱数 据处理的一个重要内容，成为高光谱数据研究领域的一大热点。

高光谱图像可以提供丰富的信息，获得确定物质或地物性质的光谱的同时， 揭示地物之间的空间位置关系，实现了“图谱合一”，进而能显著地提高数据分 析的可靠性和细节性。

虽然，高光谱图像包含丰富的光谱和空间信息，但是，同时也给图像分类算 法带来一系列挑战。一方面，由于空间分辨率的限制及其它因素影响，一个像元 常常由多种地物构成，这种像元被称为混合像元，而混合像元导致了高光谱图像 中“同物异谱(即相同种类地物具有不同的光谱信息)”和“同谱异物(即不同 种类地物具有相同的光谱信息)”现象的存在^[10]，不可避免地在图像分类过程中 引起地物的错分。另一方面，由于高光谱图像中数据维数很高，量化精度随之增 高，所以，在图像分类中，如果有监督信息的训练样本很少，分类精度会有明显 的下降，而且高维数据会带来大量复杂的计算。所以，在高光谱数据处理中，有 效对数据进行降维，并改进混合像元的分解方法，可以提取数据的有效信息，同 时获得更准确的分类结果。

现有的经典的分类方法主要有以下三类：

(一)无监督分类方法：如K均值聚类，是通过最小化聚类中各个点到该 类聚类中心的距离平方和的原则，实现各个点的分类。这种分类方法缺点是不能 自动调节聚类的个数。

(二)监督分类方法：如支撑向量机，是基于结构最小化的分类方法。这种 方法比K均值聚类方法具有更好的泛化能力，但支撑向量机需要有监督信息的 样本参与分类，而获得有监督信息的样本需要耗费大量的人力物力，在有监督信 息样本少的情况下，分类效果变差。

(三)半监督分类方法：这种方法融合了未标记样本和标记样本中所含信息 来改善分类器性能，提高分类精度。但是目前的半监督分类方法往往基于“严格 的聚类假设”，也就是，相似的物质拥有相同的标签的假设，这样的假设不能有 效解决混合像元被错分的问题。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于模糊标签的高光 谱数据子空间投影和分类方法，利用少量的监督信息，同时实现对高光谱遥感数 据的有效降维和分类。

实现本发明目的的技术方案是：通过判别子空间投影的方法将数据投影到低 维空间，增强数据的判别性能，进而通过引入模糊标签来构造拉普拉斯矩阵，解 决混合像元带来的错分问题，在实现降维的同时，实现高精度的分类。具体步骤 如下：

(1)将高光谱遥感影像数据库样本集分为训练样本集X和标记样本集X_l；

(2)计算由标记样本集X_l在子空间投影后生成的判别项：

$L_{dis}=\arg \underset{w}{\max }\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{i_2}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{\mathrm{Wx}}_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{i_1}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{{\mathrm{Wx}}_{i_j}^l}^{\prime }|{|}^2\right)\right),$

其中，L_dis表示判别项，N_l是标记样本的数目，，表示第i个标记样本，是的第k个异类的标记样本，是的第j个同类的标记样本，k_i2是与异 类的标记样本的数目，k_i1是与同类的标记样本的数目,W∈R^D×d是将D维空 间的数据投影到d维空间的投影矩阵，D由高光谱遥感影像自身的性质决定，d 为降维后数据的维数，且d＜＜D，Rⁿ是n维实数空间，||·||²表示两个向量间的距 离的平方；

(3)构造由训练样本集X的模糊标签确定的拉普拉斯正则项：

$R_p=\underset{w}{\arg \min }\left(\underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\right||{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{w}_{st})$

其中，R_p表示由模糊标签确定的拉普拉斯正则项，x_s和x_t分别是训练样本集X 第s和第t个样本，N为高光谱遥感影像数据中样本的个数，w_st表示样本x_s和x_t的相似度，由热核函数w_st＝exp(-||p(x_s)-p(x_t)||²/2σ²)确定，其中，p(x_s)∈R^c×1和p(x_t)∈R^c×1分别是x_s和x_t的模糊标签，p(x_s)和p(x_t)是分别由x_s和x_t依次属 于1到c类的概率组成的c×1的向量，c是高光谱遥感影像的类别数，σ是热核 函数的宽度；

(4)求解投影矩阵W和模糊标签p(x_i),i＝1,...,N

根据判别项和模糊标签构造的拉普拉斯正则项，得到目标函数 L＝L_dis-λR_p，其中，λ是正则项参数，用来平衡判别项和正则项之间的权重； 通过交替迭代求解的方法，求解投影矩阵W和模糊标签p(x_i),i＝1,...,N：

4a)固定模糊标签p(x_i),i＝1,...,N，求解投影矩阵W

此时目标函数表达式可写成：

$\left(\begin{array}{c}{L}_1=\underset{w}{\arg \max }\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s_2}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{\mathrm{Wx}}_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s_1}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{{\mathrm{Wx}}_{i_j}^l}^{\prime }|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{r}_{st}\right)\\ =\underset{w}{\arg \max }\left(tr\left(W^T\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s2}}{\Sigma }}\left|\right|x_i^l-x_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s1}}{\Sigma }}\left|\right|x_i^l-{x_{i_j}^l}^{\prime }|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|x_s-x_t|{|}^2{r}_{st}\right)W\right)\right)\\ =\underset{w}{\arg \max }\left(tr\left(W^{T}SW\right)\right)\end{array}\right)$

投影矩阵W可以通过对上式 $S=\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s2}}{\Sigma }}\left|\right|x_i^l-x_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s1}}{\Sigma }}\left|\right|x_i^l-x_{i_j}^l|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|x_s-x_t|{|}^2{r}_{st}$进行特征分解得到；

4b)固定投影矩阵W，求解模糊标签p(x_i),i＝1,...,N

此时，目标函数表达式可写成：

$\left(\begin{array}{c}{L}_2=\underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{w}_{st}\right)\\ =\underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{e}^{-\frac{\left|\right|p\left(x_s\right)-p\left(x_t\right)|{|}^2}{2{\sigma }^2}}\right)\\ \approx \underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2\left(1-\frac{\left|\right|p\left(x_s\right)-p\left(x_t\right)|{|}^2}{2{\sigma }^2}\right)\right)\end{array}\right)$

通过对L₂关于p(x_s)求导，可得到p的表达式为：

$p_k\left(x_s\right)=\frac{\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_j-{\mathrm{Wx}}_t|{|}^2{p}_k\left(x_t\right)}{\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_j-{\mathrm{Wx}}_t|{|}^2}$

其中，p_k(x_j)表示第j个样本属于第k类的概率，k的取值范围为1到c， p_k(x_t)表示第t个样本属于第k类的概率，N为训练样本的数目；

4c)通过L计算目标函数值，并计算ΔL＝L_n+1-L_n

L_n+1是第n+1次迭代得到的结果，L_n是第n次迭代得到的结果，当ΔL小 于设定的阈值或者迭代次数达到设定的最大迭代次数，则停止迭代转向下一步， 否则转向4a)；

(5)通过对p按列取最大值，找到每列的最大值所在的行数，该行数就是 各训练样本所属的类别数。

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

本发明采用判别子空间投影的方法来构造判别项，通过将数据投影到低维空 间，增强了数据的判别性能，并且引入模糊标签来构造拉普拉斯正则，解决了混 合像元带来的错分问题，在实现降维的同时，实现了高精度的分类。

附图说明

图1是本发明的流程图；

图2是本发明仿真使用的实验高光谱数据IndianPines及其真实标记图；

具体实施方式

参照图1，对本发明做进一步的详细描述。

步骤1：将遥感影像数据库样本集划分为训练数据集X和标记样本集X_l。

1a)在待处理的遥感数据集中，全部数据构成训练样本数据集X∈R^D×N，其 中,D表示训练集样本的维数，Rⁿ表示n维实数空间，N表示训练集样本的总 数；在本发明的实施实例IndianPines数据集中，样本维数D为200，训练集样本 的总数N为10366；

1b)每类随机从训练样本集中选取k个样本作为有监督信息的标记样本 集，其中，N_l＝c×k，c为高光谱图像类别数，在本发明的实施实例 IndianPines数据集中，c为16，k取8；

1c)在标记样本集X_l中，对每个标记样本通过欧氏距离寻找其k_i1个同类近 邻和k_i2个异类近邻，在本发明的实施实例IndianPines数据集中，同类近邻数k_i1为 3，异类近邻数k_i2为6。

步骤2：计算由子空间投影后的标记样本集生成的判别项。

通过对每个标记样本进行判别子空间投影后，使得同类的标记样本间的距 离更近，异类的标记样本的距离更远，因此标记样本集生成的判别项为：

$L_{dis}=\arg \underset{w}{\max }\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{i_2}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{\mathrm{Wx}}_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{i_1}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{{\mathrm{Wx}}_{i_j}^l}^{\prime }|{|}^2\right)\right),$

其中，L_dis表示判别项，N_l是标记样本的数目，，表示第i个标记样本，是的第k个异类的标记样本，是的第j个同类的标记样本，k_i2是与异 类的标记样本的数目，k_i1是与同类的标记样本的数目,W∈R^D×d是将D维空 间的数据投影到d维空间的投影矩阵，D由高光谱遥感影像自身的性质决定，d 为降维后数据的维数，且d＜＜D，Rⁿ是n维实数空间，||·||²表示两个向量间的距 离的平方，在本发明的实施实例IndianPines数据集中，k_i1＝3,k_i2＝6，d＝40。

步骤3：构造由训练样本的模糊标签确定的拉普拉斯正则项。

$R_p=\underset{w}{\arg \min }\left(\underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\right||{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{w}_{st})$

其中，R_p表示由模糊标签确定的拉普拉斯正则项，x_s和x_t分别是训练样本集X 第s和第t个样本，N为高光谱遥感影像数据中样本的个数，w_st表示样本x_s和x_t的相似度，由热核函数w_st＝exp(-||p(x_s)-p(x_t)||²/2σ²)确定，其中，p(x_s)∈R^c×1和p(x_t)∈R^c×1分别是x_s和x_t的模糊标签，p(x_s)和p(x_t)是分别由x_s和x_t依次属 于1到c类的概率组成的c×1的向量，c是高光谱遥感影像的类别数，σ是热核 函数的宽度；

步骤4：求解投影矩阵W和模糊标签p(x_i),i＝1,...,N

根据判别项和模糊标签构造的拉普拉斯正则项，得到目标函数 L＝L_dis-λR_p，其中，λ是正则项参数，用来平衡判别项和正则项之间的权重； 通过交替迭代求解的方法，求解投影矩阵W和模糊标签p(x_i),i＝1,...,N：

4a)固定模糊标签p(x_i),i＝1,...,N，求解投影矩阵W

此时目标函数表达式可写成：

$\left(\begin{array}{c}{L}_1=\underset{w}{\arg \max }\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s_2}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{\mathrm{Wx}}_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s_1}}{\Sigma }}\left|\right|{\mathrm{Wx}}_i^l-{{\mathrm{Wx}}_{i_j}^l}^{\prime }|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{r}_{st}\right)\\ =\underset{w}{\arg \max }\left(tr\left(W^T\left(\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s2}}{\Sigma }}\left|\right|x_i^l-x_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s1}}{\Sigma }}\left|\right|x_i^l-{x_{i_j}^l}^{\prime }|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|x_s-x_t|{|}^2{r}_{st}\right)W\right)\right)\\ =\underset{w}{\arg \max }\left(tr\left(W^{T}SW\right)\right)\end{array}\right)$

投影矩阵W可以通过对上式 $S=\underset{i=1}{\overset{N_l}{\Sigma }}\left(\underset{k=1}{\overset{k_{s2}}{\Sigma }}\left|\right|x_i^l-x_{i_k}^l|{|}^2-\underset{j=1}{\overset{k_{s1}}{\Sigma }}\left|\right|x_i^l-x_{i_j}^l|{|}^2\right)-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|x_s-x_t|{|}^2{r}_{st}$进行特征分解得到；

4b)固定投影矩阵W，求解模糊标签p(x_i),i＝1,...,N

此时，目标函数表达式可写成：

$\left(\begin{array}{c}{L}_2=\underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{w}_{st}\right)\\ =\underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2{e}^{-\frac{\left|\right|p\left(x_s\right)-p\left(x_t\right)|{|}^2}{2{\sigma }^2}}\right)\\ \approx \underset{p}{\arg \max }\left(-\lambda \underset{s=1}{\overset{N}{\Sigma }}\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_s-{\mathrm{Wx}}_t|{|}^2\left(1-\frac{\left|\right|p\left(x_s\right)-p\left(x_t\right)|{|}^2}{2{\sigma }^2}\right)\right)\end{array}\right)$

通过对L₂关于p(x_s)求导，可得到p的表达式为：

$p_k\left(x_s\right)=\frac{\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_j-{\mathrm{Wx}}_t|{|}^2{p}_k\left(x_t\right)}{\underset{t=1}{\overset{N}{\Sigma }}\left|\right|{\mathrm{Wx}}_j-{\mathrm{Wx}}_t|{|}^2}$

其中，p_k(x_j)表示第j个样本属于第k类的概率，k的取值范围为1到c， p_k(x_t)表示第t个样本属于第k类的概率，N为训练样本的数目；

4c)通过L计算目标函数值，并计算ΔL＝L_n+1-L_n

L_n+1是第n+1次迭代得到的结果，L_n是第n次迭代得到的结果，当ΔL小 于设定的阈值或者迭代次数达到设定的最大迭代次数，则停止迭代转向下一步， 否则转向4a)；在本发明的实施实例IndianPines数据集中，设定的阈值为10^-4， 最大迭代次数为40。

步骤5：通过对p按列取最大值，找到每列的最大值所在的行数，该行数就 是各训练样本所属的类别数。

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

1，仿真实验条件。

本实验采用IndianPines数据集作为实验数据，采用软件MATLAB7.10.0作为 仿真工具，计算机配置为IntelCorei5/2.27G/2G。

IndianPines高光谱数据92AV3C：该场景为AVIRIS传感器于1992年6月获 得的印第安纳州西北部的IndianPines测试地，该数据大小为145×145，每个像元 有220个波段，去掉含有噪声的20个波段，仅保留剩下的200个波段，该数据 共包含16类地物，图2(a)给出了IndianPines高光谱数据，图2(b)给出了 IndianPines高光谱数据的真实标记图。

2.仿真实验内容。

仿真1，在图2(a)所给的IndianPines高光谱数据上进行每类8个标记样本 情况下的仿真实验，并将本发明方法与现有的以下四种降维方法进行对比：1) 主成分分析PCA；2)局部fisher判别分析LFDA；3)最大化边际原则MMC； 4)基于成对约束的半监督降维SSDR。

实验中，本发明同类近邻数k_i1＝3，异类近邻数k_i2＝6，降维后的维数d＝40， 正则项参数λ＝0.8，表中OA代表总体分类精度。

表1给出了每类标记样本个数为8，对比方法采用最近邻分类器，每种方法进 行30次仿真时的实验对比结果。

表1：本发明与现有方法在每类8个标记样本个数下的对比结果

方法 本发明 PCA FLDA MMC SSDR OA 83.64％ 65.31％ 78.37％ 65.9％ 62.64％

从表1可见，本发明在每类标记样本个数为8时，精度在表中所列的五种方法 中是最高的，因此具有最好的分类效果。