CN104268561A - 基于结构先验低秩表示的高光谱图像解混方法 - Google Patents

Abstract

本发明属于图像处理技术领域，涉及一种基于结构先验低秩表示的高光谱图像解混方法，包括原始高光谱数据进行光谱矫正、估计矫正后的高光谱数据的信号子空间、将校正后的数据以及光谱库作为低秩表示模型的输入数据和字典，通过解此优化问题得到整个光谱库的丰度矩阵、剔除光谱库中非真正端元的光谱、将修剪过的光谱库作为端元重新代入低秩表示模型计算所对应的丰度矩阵，反复迭代最终得到包含真正的端元矩阵以及所对应的丰度矩阵。利用光谱库作为端元，避免了直接从原始数据中提取端元，且由于对字典进行了修剪更新，降低了光谱库中物质光谱自相干性对解混效果的不利影响，能有效地识别出端元信号，提高了丰度估计的准确率。

Description

基于结构先验低秩表示的高光谱图像解混方法

技术领域

本发明属于图像处理技术领域，涉及目标识别方法，可应用于矿物识别，具体涉及一种基于结构先验低秩表示的高光谱图像解混方法。

背景技术

高光谱成像技术是20世纪80年代初在遥感领域发展起来的新兴成像技术，典型的硬件设备为成像光谱仪。成像光谱仪通过分光技术，将电磁波信号分解为许多微小、相邻的波段，对应波段上的能量被不同的传感器捕获，因此形成高光谱遥感图像波段数目多、光谱分辨率高的数据特点。高光谱图像的另外一个显著特点为图谱合一，即同时记录了地物的图像和光谱信息，各个波段分别成像，因此各个波段的图像叠合在一起形成数据立方体，每个像元在各个波段图像上的辐射值可以形成一条连续的光谱曲线。目前，高光谱遥感已经得到广泛的应用，研究表明，高光谱遥感在专题图和地形图的测制和更新、环境监测与灾害评估、岩矿的探测和识别、农作物和植被的精细分类等领域具有良好的应用前景。

然而，混合像元的广泛存在已经成为制约高光谱遥感应用的突出问题。遥感图像像元记录的是探测单元的瞬时视场角所对应的地面范围内目标的辐射能量总和，如果探测单元的瞬时视场角所对应的地面范围仅包含了同一类性质的目标，则该像元记录的是同一性质的地面目标的辐射能量总和，这样的像元称为纯像元；如果探测单元的瞬时视场角所对应的地面范围包含了多类不同性质的目标，则该像元记录的是多类不同性质地面目标的辐射能量总和，这样的像元称为混合像元。由于地面的复杂多样性及传感器空间分辨率的限制，混合像元在高光谱图像上更为广泛存在，它不仅影响了基于高光谱图像的地物识别精度，而且已经成为高光谱遥感向定量化方向深入发展的主要障碍。因此，如何有效地解决混合像元问题是高光谱图像处理技术面临的难题之一。

目前解决混合像元问题最为有效的方法为混合像元分解。它是将混合像元的测量光谱分解为一组组成光谱(端元)和相应的丰度的过程，丰度表示了每种端元在混合像元中所占的比例。混合像元分解能够突破传感器空间分辨率的限制，在亚像元精度上获取混合像元的真实属性信息，提高影像分类的精度。

传统的混合像元分解方法主要分为端元提取和丰度估计这两个步骤。研究者们提出了许多基于统计学和几何学的端元提取方法，这些方法基于一个假设即原始高光谱数据中存在纯象元。但是由于实际这样的假设不一定成立，因此从原始数据中识别端元并不是一个理想的方法。于是近几年研究者们提出了基于稀疏回归的高光谱解混方法，该方法将光谱库作为端元，由于与整个光谱库中的信号相比，一幅图像中的端元数目是很小的，则与之对应的丰度矩阵是稀疏的，但是由于光谱库中光谱的自相干性较高，使得不能有效地识别出端元信号，降低了丰度估计的准确率，致使最终的解混效果不是很理想。

发明内容

本发明的目的是克服现有的高光谱图像解混方法由于光谱库中光谱的自相干性较高，使得不能有效地识别出端元信号，降低了丰度估计的准确率，致使最终的解混效果不是很理想的问题。

为此，本发明提供了一种基于结构先验低秩表示的高光谱图像解混方法，包括如下步骤：

(1)结合光谱库A∈R^L×m对原始高光谱数据进行光谱校正，得到校正后的数据Y∈R^L×N，

其中，L表示高光谱数据的波段数，m表示光谱库中的物质光谱总个数，N表示高光谱数据样本总个数，R表示实数域；

(2)使用基于最小错误高光谱信号子空间识别法估计校正后的数据Y的信号子空间，得到信号子空间维数p，即估计出的端元数目；

(3)设当前光谱库中物质光谱个数为r，初始化r＝m；

(4)将数据Y和光谱库A作为低秩表示模型的输入信号和字典，求解得到对应的系数矩阵X即丰度矩阵：

$\underset{X,E}{\min }{\left|\right|X\left|\right|}_{*}+\mathrm{\λ}{\left|\right|E\left|\right|}_{\mathrm{2,1}}$

约束条件：Y＝AX+E,X≥0，1^TX＝1^T

其中，E∈R^L×N是噪声矩阵，X∈R^r×N是求解得到的丰度矩阵，||X||_*是矩阵X的核范数，||E||_s,1是矩阵E的混合l_2,1范数，X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束；

(5)判断当前光谱库A中的物质光谱个数r与端元数目的估计值p之差是否小于预设的阈值n，对于真实数据，n在10-20之间取值；

如果是，则终止，丰度矩阵X即为最终输出结果；

如果否，则根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，其中m′是光谱库经修剪后保留的物质光谱个数且m′＞p,然后将A_S代替A且r＝m′，重复步骤(4)和(5)，直至输出最终结果。

上述步骤(2)中使用的基于最小错误高光谱信号子空间识别法具体实施方法如下：

2a)输入数据Y＝[y₁,y₂,...,y_N]∈R^L×N，其中y_i为数据Y的第i列，i＝1,2,...,N，并计算

2b)估计噪声 $\widehat{\mathrm{\ϵ}}=[{\widehat{\mathrm{\ϵ}}}_1,{\widehat{\mathrm{\ϵ}}}_2,...,{\widehat{\mathrm{\ϵ}}}_N]\∈R^{L\×N},$ 其中为的第i列，i＝1,2,...,N，具体计算步骤如下：

第一步：计算Z＝Y^T， $\hat{K}\left(Z^{T}Z\right),K^{\′}={\hat{K}}^{-1};$

第二步：设置i＝1,2,...,L，分别计算 ${\widehat{\mathrm{\β}}}_i=({\left[K^{\′}\right]}_{{\∂}_i,{\∂}_i}-{\left[K^{\′}\right]}_{{\∂}_i,i}{\left[K^{\′}\right]}_{i,{\∂}_i}/{\left[K^{\′}\right]}_{i,i}){\left[\hat{K}\right]}_{{\∂}_i,i},{\widehat{\mathrm{\ϵ}}}_i=z_i-Z_{{\∂}_i}{\widehat{\mathrm{\β}}}_i,$ 其中表示矩阵K′去掉第i行第i列后的矩阵，表示矩阵K′去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵[K′]_i,i表示矩阵K′的第i行第i列元素，表示矩阵去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵z_i表示矩阵Z的第i列，表示矩阵Z去掉第i列后的矩阵；

第三步：输出 $\widehat{\mathrm{\ϵ}}=[{\widehat{\mathrm{\ϵ}}}_1,{\widehat{\mathrm{\ϵ}}}_2,...,{\widehat{\mathrm{\ϵ}}}_N]\∈R^{L\×N};$

2c)估计噪声相关矩阵

2d)估计信号相关矩阵 ${\hat{K}}_x=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left((y_i-{\widehat{\mathrm{\ϵ}}}_i)(y_i-{\widehat{\mathrm{\ϵ}}}_i^T)\right);$

2e)计算其中W_p＝[w₁,w₂,...,w_L]，w_i为W_p的第i列，i＝1,2,...,L，是的特征向量；

2f)计算 $p=\arg \underset{p}{\min }\{\mathrm{tr}\left(U_p^{\⊥}{\hat{K}}_y\right)+2\mathrm{tr}\left(U_p{\hat{K}}_n\right)\},$ 其中 $U_p^{\⊥}=I-U_p,$ I为单位矩阵，tr()表示求矩阵主对角线上元素之和。

上述的步骤(4)中，||X||_*是矩阵的核范数，具体定义如下：

${\left|\right|X\left|\right|}_{*}=\mathrm{trace}\left(\sqrt{X^{T}X}\right)=\underset{i=1}{\overset{\min \{r,N\}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i$

其中σ_i(0≤i≤min{r,N})是矩阵X的奇异值。

上述的步骤(4)中，||E||_2,1是矩阵的混合l_2,1范数，定义如下：

${\left|\right|E\left|\right|}_{\mathrm{2,1}}={\mathrm{\Σ}}_{i=1}^L\sqrt{{\mathrm{\Σ}}_{j=1}^N{\left(e_{\mathrm{ij}}\right)}^2}$

其中e_ij表示矩阵E的第i行第j列元素。

上述的步骤(4)中，X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束，其中X≥0是丰度矩阵的“非负”约束，即矩阵X中的所有元素必须大于0；1^TX＝1^T为丰度的“和为一”约束，即矩阵X中的每一列之和必须为1，其中矩阵X中的每一列代表每个像元对应不同端元的丰度值，具体如下列公式所示：

$x_i^j\≥0,i=\mathrm{1,2},...,r,j=\mathrm{1,2},...,N$

$\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_i^j=1,j=\mathrm{1,2},...,N$

其中是X中第i行第j列的元素，表示第j个样本的第i个端元的丰度值。

上述的步骤(5)中，根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，包括如下步骤：

5a)设定修剪阈值t＝T×d，其中T为预设的定值，d为迭代次数，迭代次数即修剪次数，初始值为1，随着迭代次数的增加修剪阈值t也在不断增加，有利于加快光谱库的修剪；

5b)计算丰度矩阵X的第i行向量x_i∈R^1×N,1≤i≤r中小于第一步中定义的修剪阈值的元素的个数k；

5c)判断k是否等于N，如果是，则将其对应的A中第i个物质光谱剔除掉，否则，保留该物质光谱。

本发明的有益效果：

1、本发明由于采用基于结构先验的低秩表示模型进行高光谱图像混合像元分解，与现有的基于稀疏表示的解混方法相比，能够利用丰度矩阵的全局结构进行混合像元分解。

2、本发明与简单的基于低秩表示的解混方法相比，将光谱库作为端元矩阵，避免了直接从原始数据中提取端元，且由于对字典进行了修剪更新，降低了光谱库中信号自相干性对解混效果的不利影响。

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

附图说明

图1是发明的流程示意图。

图2是本发明仿真采用的模拟数据高光谱图像及其所含五种端元的真实丰度图，其中(a)为模拟数据高光谱图，(b)-(f)分别为五种端元的真实丰度图。

图3是真实数据高光谱图像。

图4是高光谱图像在不同解混算法下的丰度图。

具体实施方式

现有的高光谱图像解混方法由于光谱库中光谱的自相干性较高，使得不能有效地识别出端元信号，降低了丰度估计的准确率，致使最终的解混效果不是很理想，为了解决上述问题，本发明提出一种如图1所示的高光谱图像解混方法。

这种基于结构先验低秩表示的高光谱图像解混方法，包括如下步骤：

(1)结合光谱库A∈R^L×m对原始高光谱数据进行光谱校正，得到校正后的数据Y∈R^L×N，

其中，L表示高光谱数据的波段数，m表示光谱库中的物质光谱总个数，N表示高光谱数据样本总个数，R表示实数域；

(2)使用基于最小错误高光谱信号子空间识别法估计校正后的数据Y的信号子空间，得到信号子空间维数p，即估计出的端元数目，具体实施方法如下：

2a)输入数据Y＝[y₁,y₂,...,y_N]∈R^L×N，其中y_i为数据Y的第i列，i＝1,2,...,N，并计算

2b)估计噪声 $\widehat{\mathrm{\ϵ}}=[{\widehat{\mathrm{\ϵ}}}_1,{\widehat{\mathrm{\ϵ}}}_2,...,{\widehat{\mathrm{\ϵ}}}_N]\∈R^{L\×N},$ 其中为的第i列，i＝1,2,...,N，具体计算步骤如下：

第一步：计算Z＝Y^T， $\hat{K}\left(Z^{T}Z\right),K^{\′}={\hat{K}}^{-1};$

第二步：设置i＝1,2,...,L，分别计算 ${\widehat{\mathrm{\β}}}_i=({\left[K^{\′}\right]}_{{\∂}_i,{\∂}_i}-{\left[K^{\′}\right]}_{{\∂}_i,i}{\left[K^{\′}\right]}_{i,{\∂}_i}/{\left[K^{\′}\right]}_{i,i}){\left[\hat{K}\right]}_{{\∂}_i,i},{\widehat{\mathrm{\ϵ}}}_i=z_i-Z_{{\∂}_i}{\widehat{\mathrm{\β}}}_i,$ 其中表示矩阵K′去掉第i行第i列后的矩阵，表示矩阵K′去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵[K′]_i,i表示矩阵K′的第i行第i列元素，表示矩阵去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵z_i表示矩阵Z的第i列，表示矩阵Z去掉第i列后的矩阵；

第三步：输出 $\widehat{\mathrm{\ϵ}}=[{\widehat{\mathrm{\ϵ}}}_1,{\widehat{\mathrm{\ϵ}}}_2,...,{\widehat{\mathrm{\ϵ}}}_N]\∈R^{L\×N};$

2c)估计噪声相关矩阵

2d)估计信号相关矩阵 ${\hat{K}}_x=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left((y_i-{\widehat{\mathrm{\ϵ}}}_i)(y_i-{\widehat{\mathrm{\ϵ}}}_i^T)\right);$

2e)计算其中W_p＝[w₁,w₂,...,w_L]，w_i为W_p的第i列，i＝1,2,...,L，是的特征向量；

2f)计算 $p=\arg \underset{p}{\min }\{\mathrm{tr}\left(U_p^{\⊥}{\hat{K}}_y\right)+2\mathrm{tr}\left(U_p{\hat{K}}_n\right)\},$ 其中 $U_p^{\⊥}=I-U_p,$ I为单位矩阵，tr()表示求矩阵主对角线上元素之和。

(3)设当前光谱库中物质光谱个数为r，初始化r＝m；

(4)将数据Y和光谱库A作为低秩表示模型的输入信号和字典，求解得到对应的系数矩阵X即丰度矩阵：

$\underset{X,E}{\min }{\left|\right|X\left|\right|}_{*}+\mathrm{\λ}{\left|\right|E\left|\right|}_{\mathrm{2,1}}$

约束条件：Y＝AX+E,X≥0，1^TX＝1^T

其中，E∈R^L×N是噪声矩阵，X∈R^r×N是求解得到的丰度矩阵，||X||_*是矩阵X的核范数，||E||_2,1是矩阵E的混合l_2,1范数，X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束；

||X||_*是矩阵的核范数，具体定义如下：

${\left|\right|X\left|\right|}_{*}=\mathrm{trace}\left(\sqrt{X^{T}X}\right)=\underset{i=1}{\overset{\min \{r,N\}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i$

其中σ_i(0≤i≤min{r,N})是矩阵X的奇异值。

||E||_2,1是矩阵的混合l_2,1范数，定义如下：

${\left|\right|E\left|\right|}_{\mathrm{2,1}}={\mathrm{\Σ}}_{i=1}^L\sqrt{{\mathrm{\Σ}}_{j=1}^N{\left(e_{\mathrm{ij}}\right)}^2}$

其中e_ij表示矩阵E的第i行第j列元素。

X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束，其中X≥0是丰度矩阵的“非负”约束，即矩阵X中的所有元素必须大于0；1^TX＝1^T为丰度的“和为一”约束，即矩阵X中的每一列之和必须为1，其中矩阵X中的每一列代表每个像元对应不同端元的丰度值，具体如下列公式所示：

$x_i^j\≥0,i=\mathrm{1,2},...,r,j=\mathrm{1,2},...,N$

$\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_i^j=1,j=\mathrm{1,2},...,N$

其中是X中第i行第j列的元素，表示第j个样本的第i个端元的丰度值。

(5)判断当前光谱库A中的物质光谱个数r与端元数目的估计值p之差是否小于预设的阈值n，对于真实高光谱数据，n在10-20之间取值；

如果是，则终止，丰度矩阵X即为最终输出结果；

如果否，则根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，其中m′是光谱库经修剪后保留的物质光谱个数且m′＞p,然后将A_S代替A且r＝m′，重复步骤(4)和(5)，直至输出最终结果。

根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，包括如下步骤：

5a)设定修剪阈值t＝T×d，其中T为预设的定值，d为迭代次数，迭代次数即修剪次数，初始值为1，随着迭代次数的增加修剪阈值t也在不断增加，有利于加快光谱库的修剪；

5b)计算丰度矩阵X的第i行向量x_i∈R^1×N,1≤i≤r中小于第一步中定义的修剪阈值的元素的个数k；

5c)判断k是否等于N，如果是，则将其对应的A中第i个物质光谱剔除掉，否则，保留该物质光谱。

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

1.仿真条件:

仿真实验采用一个模拟数据和一个真实数据。模拟数据由USGS光谱库(共498种光谱信号)随机选择的240种光谱信号作为光谱库A，并从中随机选择5种光谱线性混合并加入不同信噪比的噪声而生成，如图2(a),图2(b-f)显示的是5种端元的真实丰度图。

真实数据是著名的AVIRISCuprite高光谱图像的一部分，其对应的光谱库是USGS光谱库，实验中将该光谱库整个498种光谱信号作为光谱库A，图3显示了该高光谱数据。在CPU为IntelCore(TM)2Duo、主频2.33GHz，内存为2G的WINDOWSXP系统上用MATLAB7.0.1软件进行仿真。

2.仿真内容:

仿真内容主要包括本发明与现有几种方法在模拟数据上的正确率比较实验以及在真实数据上的定性对比实验，现有几种方法包括：非负约束最小二乘法NCLS，稀疏法SUnSAL，稀疏结合全变异空间正则SUnSAL-TV这三种方法；本发明的基于结构先验低秩表示的高光谱图像解混方法缩写为SP-LRR。

3.仿真实验结果分析：

1)模拟数据

表1展示的是模拟数据在NCLS，SUnSAL，SUnSAL-TV和本发明SP-LRR在信噪比分别为20dB、30dB、40dB时的丰度矩阵信号重构误差(SRE)。用dB衡量为：SRE(dB)≡10log₁₀(SRE)，其中x为真实丰度矩阵，为算法估计出的丰度矩阵，E[﹒]为求期望操作,相关参数设置为：T＝0.02,n＝2。

表1:不同方法在不同信噪比时的SRE(dB)

从表1可以得出，本发明比NCLS、SUnSAL和SUnSAL-TV方法能够得到更高的SRE值，即性能更好。由于本模拟数据具有鲜明的空间一致性，所以SUnSAL-TV方法因为加入了全变异正则项考虑到了空间一致性较NCLS和SUnSAL方法性能有了较大的提高。而本发明由于在字典更新的过程中逐渐剔除了非端元信号仅保留端元信号，所以大大下降了错误率，获得了高于SUnSAL-TV方法的SRE值。

2)真实数据

图4展示了AVIRISCuprite高光谱图像在不同解混算法下的丰度图，通过对比可以得到该数据针对不同算法性能的定性对比结果，相关参数设置为：T＝0.01,n＝15。从图4中可以看出，本发明的结果与NCLS、SUnSAL方法估计出的丰度相比，数值较高，形状轮廓更加鲜明，说明了本发明的有效性。

综上，本发明具有利用数据空间结构以及减轻光谱库互相干性不利影响的优点，并能在模拟数据和真实数据上得到较好的丰度分布图，与现有的方法相比具有一定的优势，能有效地识别出端元信号，提高了丰度估计的准确率。

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

Claims (6)

1.基于结构先验低秩表示的高光谱图像解混方法，其特征是：包括如下步骤：

(1)结合光谱库A∈R^L×m对原始高光谱数据进行光谱校正，得到校正后的数据Y∈R^L×N，

其中，L表示高光谱数据的波段数，m表示光谱库中的物质光谱总个数，N表示高光谱数据样本总个数，R表示实数域；

(2)使用基于最小错误高光谱信号子空间识别法估计校正后的数据Y的信号子空间，得到信号子空间维数p，即估计出的端元数目；

(3)设当前光谱库中物质光谱个数为r，初始化r＝m；

(4)将数据Y和光谱库A作为低秩表示模型的输入信号和字典，求解得到对应的系数矩阵X即丰度矩阵：

$\underset{X,E}{\min }{\left|\right|X\left|\right|}_{*}+\mathrm{\λ}{\left|\right|E\left|\right|}_{\mathrm{2,1}}$

约束条件：Y＝AX+E,X≥0，1^TX＝1^T

其中，E∈R^L×N是噪声矩阵，X∈R^r×N是求解得到的丰度矩阵，||X||_*是矩阵X的核范数，||E||_2,1是矩阵E的混合l_2,1范数，X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束；

(5)判断当前光谱库A中的物质光谱个数r与端元数目的估计值p之差是否小于预设的阈值n，对于真实高光谱数据，n在10-20之间取值；

如果是，则终止，丰度矩阵X即为最终输出结果；

如果否，则根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，其中m′是光谱库经修剪后保留的物质光谱个数且m′＞p,然后将A_S代替A且r＝m′，重复步骤(4)和(5)，直至输出最终结果。

2.根据权利要求1所述的基于结构先验低秩表示的高光谱图像解混方法，其特征是：所述步骤(2)中的基于最小错误高光谱信号子空间识别法具体实施方法如下：

2a)输入数据Y＝[y₁,y₂,...,y_N]∈R^L×N，其中y_i为数据Y的第i列，i＝1,2,...,N，并计算

2b)估计噪声其中为的第i列，i＝1,2,...,N，具体计算步骤如下：

第一步：计算Z＝Y^T， $\hat{K}\left(Z^{T}Z\right),K^{\′}={\hat{K}}^{-1};$

第二步：设置i＝1,2,...,L，分别计算 ${\widehat{\mathrm{\β}}}_i=({\left[K^{\′}\right]}_{{\∂}_i,{\∂}_i}-{\left[K^{\′}\right]}_{{\∂}_i,i}{\left[K^{\′}\right]}_{i,{\∂}_i}/{\left[K^{\′}\right]}_{i,i}){\left[\hat{K}\right]}_{{\∂}_i,i},{\widehat{\mathrm{\ϵ}}}_i=z_i-Z_{{\∂}_i}{\widehat{\mathrm{\β}}}_i,$ 其中表示矩阵K′去掉第i行第i列后的矩阵，表示矩阵K′去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵[K′]_i,i表示矩阵K′的第i行第i列元素，表示矩阵去掉第i列后的矩阵，表示矩阵的第i行，表示矩阵z_i表示矩阵Z的第i列，表示矩阵Z去掉第i列后的矩阵；

第三步：输出 $\widehat{\mathrm{\ϵ}}=[{\widehat{\mathrm{\ϵ}}}_1,{\widehat{\mathrm{\ϵ}}}_2,...,{\widehat{\mathrm{\ϵ}}}_N]\∈R^{L\×N};$

2c)估计噪声相关矩阵

2d)估计信号相关矩阵 ${\hat{K}}_x=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left((y_i-{\widehat{\mathrm{\ϵ}}}_i)(y_i-{\widehat{\mathrm{\ϵ}}}_i^T)\right);$

2e)计算其中W_p＝[w₁,w₂,...,w_L]，w_i为W_p的第i列，i＝1,2,...,L，是的特征向量；

2f)计算 $p=\arg \underset{p}{\min }\{\mathrm{tr}\left(U_p^{\⊥}{\hat{K}}_y\right)+2\mathrm{tr}\left(U_p{\hat{K}}_n\right)\},$ 其中 $U_p^{\⊥}=I-U_p,$ I为单位矩阵，tr()表示求矩阵主对角线上元素之和。

3.根据权利要求1所述的基于结构先验低秩表示的高光谱图像解混方法，其特征是：所述的步骤(4)中，||X||_*是矩阵的核范数，具体定义如下：

${\left|\right|X\left|\right|}_{*}=\mathrm{trace}\left(\sqrt{X^{T}X}\right)=\underset{i=1}{\overset{\min \{r,N\}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i$

其中σ_i(0≤i≤min{r,N})是矩阵X的奇异值。

4.根据权利要求1所述的基于结构先验低秩表示的高光谱图像解混方法，其特征是：所述的步骤(4)中，||E||_2,1是矩阵的混合l_2,1范数，定义如下：

${\left|\right|E\left|\right|}_{\mathrm{2,1}}={\mathrm{\Σ}}_{i=1}^L\sqrt{{\mathrm{\Σ}}_{j=1}^N{\left(e_{\mathrm{ij}}\right)}^2}$

其中e_ij表示矩阵E的第i行第j列元素。

5.根据权利要求1所述的基于结构先验低秩表示的高光谱图像解混方法，其特征是：所述的步骤(4)中，X≥0,1^TX＝1^T是丰度的“非负”以及“和为一”约束，其中X≥0是丰度矩阵的“非负”约束，即矩阵X中的所有元素必须大于0；1^TX＝1^T为丰度的“和为一”约束，即矩阵X中的每一列之和必须为1，其中矩阵X中的每一列代表每个像元对应不同端元的丰度值，具体如下列公式所示：

$x_i^j\≥0,i=\mathrm{1,2},...,r,j=\mathrm{1,2},...,N$

$\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{x}_i^j=1,j=\mathrm{1,2},...,N$

其中是X中第i行第j列的元素，表示第j个样本的第i个端元的丰度值。

6.根据权利要求1所述的基于结构先验低秩表示的高光谱图像解混方法，其特征是：所述的步骤(5)中，根据丰度矩阵X的稀疏性，从光谱库中剔除非真正的端元信号得到新的光谱库A_S∈R^L×m′，包括如下步骤：

5a)设定修剪阈值t＝T×d，其中T为预设的定值，d为迭代次数，迭代次数即修剪次数，初始值为1，随着迭代次数的增加修剪阈值t也在不断增加，有利于加快光谱库的修剪；

5b)计算丰度矩阵X的第i行向量x_i∈R^1×N,1≤i≤r中小于第一步中定义的修剪阈值的元素的个数k；

5c)判断k是否等于N，如果是，则将其对应的A中第i个物质光谱剔除掉，否则，保留该物质光谱。