CN108734672B - Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression - Google Patents

Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression Download PDF

Info

Publication number
CN108734672B
CN108734672B CN201810016715.4A CN201810016715A CN108734672B CN 108734672 B CN108734672 B CN 108734672B CN 201810016715 A CN201810016715 A CN 201810016715A CN 108734672 B CN108734672 B CN 108734672B
Authority
CN
China
Prior art keywords
spectral library
spectral
library
hyperspectral
hyperspectral data
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810016715.4A
Other languages
Chinese (zh)
Other versions
CN108734672A (en
Inventor
肖亮
李生富
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Science and Technology
Original Assignee
Nanjing University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Science and Technology filed Critical Nanjing University of Science and Technology
Priority to CN201810016715.4A priority Critical patent/CN108734672B/en
Publication of CN108734672A publication Critical patent/CN108734672A/en
Application granted granted Critical
Publication of CN108734672B publication Critical patent/CN108734672B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/70—Denoising; Smoothing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/50—Image enhancement or restoration using two or more images, e.g. averaging or subtraction
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing
    • G06T2207/10036—Multispectral image; Hyperspectral image
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A40/00—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production
    • Y02A40/10—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production in agriculture

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明公开了一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，包含如下步骤：输入光谱库和待解混高光谱数据；光谱库初始化与光谱数据矩阵构造；构建光谱拟合误差保真项；构建丰度矩阵协同稀疏性约束项；构建光谱库裁剪的稀疏正则化项；建立光谱库裁剪和协同稀疏回归模型；迭代求解光谱库裁剪和协同稀疏回归模型；输出裁剪后的光谱库和端元丰度图。本发明为缓解光谱库端元与实际场景中端元的误匹配问题，建立协同稀疏回归模型，进行目标函数最优化解混，提高了解混的精度，降低了端元误匹配率，增强了对光谱幅度变化和噪声的鲁棒性；可广泛应用于环境监测、矿产勘探和精准农业等领域的高光谱数据解混应用。

Figure 201810016715

The invention discloses a hyperspectral data unmixing method based on spectral library clipping and cooperative sparse regression, comprising the following steps: inputting a spectral library and hyperspectral data to be unmixed; initializing the spectral library and constructing a spectral data matrix; constructing spectral fitting Error fidelity term; construct abundance matrix collaborative sparsity constraint term; construct sparse regularization term for spectral library clipping; build spectral library clipping and collaborative sparse regression model; iteratively solve spectral library clipping and collaborative sparse regression model; output clipped Spectral library and endmember abundance map. In order to alleviate the problem of mismatch between endmembers in the spectral library and endmembers in actual scenarios, the invention establishes a collaborative sparse regression model, performs objective function optimization and unmixing, improves the accuracy of unmixing, reduces the endmember mismatch rate, and enhances the accuracy of the unmixing. Robustness to spectral amplitude changes and noise; can be widely used in hyperspectral data unmixing applications in environmental monitoring, mineral exploration, and precision agriculture.

Figure 201810016715

Description

基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法Hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression

技术领域technical field

本发明涉及遥感高光谱数据处理技术，具体涉及一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法。The invention relates to a remote sensing hyperspectral data processing technology, in particular to a hyperspectral data unmixing method based on spectral library clipping and cooperative sparse regression.

背景技术Background technique

高光谱数据由于其光谱相关性及丰富的空间信息而被广泛应用于军事监测、精准农业和矿物勘探等领域。其中，高光谱数据解混是定量遥感分析的关键技术。高光谱数据解混的基本原理是将单个像元光谱分解成若干个纯净像元光谱的组合。其理论依据是由于成像光谱仪空间分辨率的限制，获得的高光谱图像中存在大量的混合像元，每个混合像元中包含多种纯净物质。Hyperspectral data are widely used in military monitoring, precision agriculture, and mineral exploration due to their spectral correlation and rich spatial information. Among them, hyperspectral data unmixing is a key technology for quantitative remote sensing analysis. The basic principle of hyperspectral data unmixing is to decompose a single pixel spectrum into a combination of several pure pixel spectra. The theoretical basis is that due to the limitation of the spatial resolution of the imaging spectrometer, there are a large number of mixed pixels in the obtained hyperspectral image, and each mixed pixel contains a variety of pure substances.

目前已经提出了许多针对高光谱数据的解混算法，包括纯净像元指数、顶点成分分析、独立成分分析和相关成分分析等。但是上述几类方法都是从高光谱数据中提取端元。经过数十年的研究，人们已经收集许多矿物的标准光谱，形成了庞大的光谱库，完全可以利用光谱库省去端元提取的步骤，直接进行光谱解混。Many unmixing algorithms for hyperspectral data have been proposed, including pure pixel index, vertex component analysis, independent component analysis, and correlated component analysis. However, the above-mentioned methods all extract endmembers from hyperspectral data. After decades of research, people have collected standard spectra of many minerals and formed a huge spectral library. The spectral library can be used to omit the step of endmember extraction and directly perform spectral unmixing.

2014李云松等人提出了基于邻域光谱加权的高光谱图像稀疏解混方法[西安电子科技大学,中国,基于邻域光谱加权的高光谱图像稀疏解混方法[p].发明申请 ,CN103810715A,2014-03-12]，取得了良好的解混效果。但是由于温度、光照等的影响，同一种矿物在实验室收集到的光谱曲线和在真实环境下收集到的光谱曲线有所偏差。上述方法不能有效处理光谱变化现象，当数据存在大量噪声时，算法的解混效果下降。2014 Li Yunsong et al. proposed a hyperspectral image sparse unmixing method based on neighborhood spectral weighting [Xidian University, China, Hyperspectral image sparse unmixing method based on neighborhood spectral weighting [p]. Invention application, CN103810715A, 2014 -03-12], and achieved a good unmixing effect. However, due to the influence of temperature, light, etc., the spectral curve of the same mineral collected in the laboratory deviates from the spectral curve collected in the real environment. The above methods cannot effectively deal with the phenomenon of spectral changes, and when there is a large amount of noise in the data, the unmixing effect of the algorithm decreases.

发明内容SUMMARY OF THE INVENTION

本发明的目的在于提供一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法。The purpose of the present invention is to provide a hyperspectral data unmixing method based on spectral library clipping and cooperative sparse regression.

实现本发明目的的技术解决方案为：一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，包含如下步骤：The technical solution for realizing the purpose of the present invention is: a hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression, comprising the following steps:

步骤1，输入光谱库和待解混高光谱数据；Step 1, input the spectral library and the hyperspectral data to be unmixed;

步骤2，光谱库初始化与光谱数据矩阵构造；Step 2, spectral library initialization and spectral data matrix construction;

步骤3，构建光谱拟合误差保真项；Step 3, construct a spectral fitting error fidelity term;

步骤4，构建丰度矩阵协同稀疏性约束项；Step 4, constructing the synergistic sparsity constraint term of the abundance matrix;

步骤5，构建光谱库裁剪的稀疏正则化项；Step 5, construct the sparse regularization term for spectral library clipping;

步骤6，建立光谱库裁剪和协同稀疏回归模型；Step 6: Establish spectral library clipping and collaborative sparse regression model;

步骤7，迭代求解光谱库裁剪和协同稀疏回归模型；Step 7, iteratively solve the spectral library clipping and collaborative sparse regression model;

步骤8，输出裁剪后的光谱库和端元丰度图。Step 8, output the cropped spectral library and endmember abundance map.

本发明与现有技术相比，其显著优点在于：(1)本发明充分考虑了稀疏解混过程中存在的光谱误匹配现象，利用光谱变化具有稀疏性的先验知识，并结合协同稀疏回归框架，对光谱误匹配现象进行建模，与传统的稀疏解混方法相比，改善了解混的效果，增强了对光谱幅度变化和噪声的鲁棒性；(2)本发明可广泛应用于环境监测、矿产勘探和精准农业等领域的高光谱数据解混应用。Compared with the prior art, the present invention has significant advantages as follows: (1) The present invention fully considers the spectral mismatch phenomenon existing in the sparse unmixing process, utilizes the prior knowledge that spectral changes have sparseness, and combines with collaborative sparse regression Compared with the traditional sparse unmixing method, the unmixing effect is improved, and the robustness to spectral amplitude variation and noise is enhanced; (2) the present invention can be widely used in the environment Hyperspectral data unmixing applications in monitoring, mineral exploration, and precision agriculture.

下面结合附图对本发明作进一步详细描述。The present invention will be described in further detail below with reference to the accompanying drawings.

附图说明Description of drawings

图1是本发明的基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法流程图。FIG. 1 is a flow chart of the hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression of the present invention.

图2(a)是Cuprite数据集上使用Tetracorder软件得到的明矾石矿物丰度图。Figure 2(a) is a map of the mineral abundance of alumite using Tetracorder software on the Cuprite dataset.

图2(b)是Cuprite数据集上使用SUnSAL方法得到的明矾石矿物丰度图。Figure 2(b) is a map of the mineral abundance of alumite obtained using the SUnSAL method on the Cuprite dataset.

图2(c)是Cuprite数据集上使用CLSUnSAL方法得到的明矾石矿物丰度图。Figure 2(c) is the mineral abundance map of alumite obtained using the CLSUnSAL method on the Cuprite dataset.

图2(d)是Cuprite数据集上使用DANSER方法得到的明矾石矿物丰度图。Figure 2(d) is a map of alumite mineral abundance obtained using the DANSER method on the Cuprite dataset.

图2(e)是Cuprite数据集上使用SDCSR方法得到的明矾石矿物丰度图。Figure 2(e) is the alumite mineral abundance map obtained using SDCSR method on the Cuprite dataset.

图2(f)是Cuprite数据集上使用Tetracorder软件得到的水铵长石矿物丰度图。Fig. 2(f) is a map of the mineral abundance of hydroammonium feldspar obtained using Tetracorder software on the Cuprite dataset.

图2(g)是Cuprite数据集上使用SUnSAL方法得到的水铵长石矿物丰度图。Fig. 2(g) is a map of the mineral abundance of hydroammonium feldspar obtained using the SUnSAL method on the Cuprite dataset.

图2(h)是Cuprite数据集上使用CLSUnSAL方法得到的水铵长石矿物丰度图。Fig. 2(h) is a map of the mineral abundance of hydroammonium feldspar obtained using the CLSUnSAL method on the Cuprite dataset.

图2(i)是Cuprite数据集上使用DANSER方法得到的水铵长石矿物丰度图。Fig. 2(i) is a map of the mineral abundance of hydrophosphite obtained using the DANSER method on the Cuprite dataset.

图2(j)是Cuprite数据集上使用SDCSR方法得到的水铵长石矿物丰度图。Fig. 2(j) is a map of the mineral abundance of fenugreek obtained by SDCSR method on the Cuprite dataset.

图2(k)是Cuprite数据集上使用Tetracorder软件得到的玉髓矿物丰度图。Figure 2(k) is the mineral abundance map of chalcedony obtained using Tetracorder software on the Cuprite dataset.

图2(l)是Cuprite数据集上使用SUnSAL方法得到的玉髓矿物丰度图。Figure 2(l) is the mineral abundance map of chalcedony obtained using the SUnSAL method on the Cuprite dataset.

图2(m)是Cuprite数据集上使用CLSUnSAL方法得到的玉髓矿物丰度图。Figure 2(m) is the mineral abundance map of chalcedony obtained using the CLSUnSAL method on the Cuprite dataset.

图2(n)是Cuprite数据集上使用DANSER方法得到的玉髓矿物丰度图。Figure 2(n) is the mineral abundance map of chalcedony obtained using the DANSER method on the Cuprite dataset.

图2(o)是Cuprite数据集上使用SDCSR方法得到的玉髓矿物丰度图。Figure 2(o) is the mineral abundance map of chalcedony obtained using SDCSR method on the Cuprite dataset.

图3(a)是初始光谱库和裁剪后光谱库中的明矾石光谱曲线对比图。Figure 3(a) is a comparison diagram of the spectral curves of alumite in the initial spectral library and the tailored spectral library.

图3(b)是初始光谱库和裁剪后光谱库中的水铵长石光谱曲线对比图。Figure 3(b) is a comparison diagram of the spectral curves of hydroammonium feldspar in the initial spectral library and the tailored spectral library.

图3(c)是初始光谱库和裁剪后光谱库中的玉髓光谱曲线对比图。Figure 3(c) is a comparison diagram of the chalcedony spectral curves in the initial spectral library and the tailored spectral library.

具体实施方式Detailed ways

结合图1，一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，包含如下步骤：Combined with Figure 1, a hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression, including the following steps:

步骤1，输入光谱库和待解混高光谱数据；Step 1, input the spectral library and the hyperspectral data to be unmixed;

步骤2，光谱库初始化与光谱数据矩阵构造；Step 2, spectral library initialization and spectral data matrix construction;

步骤3，构建光谱拟合误差保真项；Step 3, construct a spectral fitting error fidelity term;

步骤4，构建丰度矩阵协同稀疏性约束项；Step 4, constructing the synergistic sparsity constraint term of the abundance matrix;

步骤5，构建光谱库裁剪的稀疏正则化项；Step 5, construct the sparse regularization term for spectral library clipping;

步骤6，建立光谱库裁剪和协同稀疏回归模型；Step 6: Establish spectral library clipping and collaborative sparse regression model;

步骤7，迭代求解光谱库裁剪和协同稀疏回归模型；Step 7, iteratively solve the spectral library clipping and collaborative sparse regression model;

步骤8，输出裁剪后的光谱库和端元丰度图。Step 8, output the cropped spectral library and endmember abundance map.

进一步的，步骤1输入光谱库和待解混高光谱数据的具体过程为：Further, the specific process of inputting the spectral library and the hyperspectral data to be unmixed in step 1 is as follows:

输入一个包含高光谱端元的光谱库

Figure GDA0001783923190000031
其中整数K>0表示每个光谱信号的波段数，整数N>0表示该光谱库所包含的端元光谱信号的数量；Enter a spectral library containing hyperspectral endmembers
Figure GDA0001783923190000031
The integer K>0 represents the number of bands of each spectral signal, and the integer N>0 represents the number of endmember spectral signals contained in the spectral library;

输入待解混的高光谱数据

Figure GDA0001783923190000032
其中，整数L>0表示高光谱数据的波段数，整数W>0和整数H>0分别表示高光谱数据空间维的宽度和高度。Input hyperspectral data to be unmixed
Figure GDA0001783923190000032
Among them, the integer L>0 represents the number of bands of the hyperspectral data, and the integer W>0 and the integer H>0 represent the width and height of the spatial dimension of the hyperspectral data, respectively.

进一步的，步骤2光谱库初始化与光谱数据矩阵构造的具体步骤为：Further, in step 2, the specific steps of spectral library initialization and spectral data matrix construction are:

步骤2-1，输入一个光谱库端元

Figure GDA0001783923190000033
波段数为K>0，波段为
Figure GDA0001783923190000034
而输入一个待解混高光谱信号
Figure GDA0001783923190000035
波段数为L>0，波段为 Band＝[b1,b2,...,bL]。根据待解混高光谱信号
Figure GDA0001783923190000036
的波段Band，保留光谱库端元
Figure GDA0001783923190000037
波段
Figure GDA0001783923190000038
中对应的波段，其余的波段均丢弃；经过这样的一步操作，光谱库中端元信号
Figure GDA0001783923190000039
的波段就和待解混高光谱信号
Figure GDA00017839231900000310
中的波段保持一致，得到剔除波段后的光谱库端元
Figure GDA00017839231900000311
波段数为L>0，波段为Band＝[b1,b2,...,bL]；Step 2-1, enter a spectral library endmember
Figure GDA0001783923190000033
The number of bands is K>0, and the band is
Figure GDA0001783923190000034
And input a hyperspectral signal to be unmixed
Figure GDA0001783923190000035
The number of bands is L>0, and the band is Band=[b 1 ,b 2 ,...,b L ]. According to the hyperspectral signal to be unmixed
Figure GDA0001783923190000036
Band, retains spectral library endmembers
Figure GDA0001783923190000037
band
Figure GDA0001783923190000038
The corresponding bands in the spectral library are discarded, and the rest of the bands are discarded; after such a one-step operation, the endmember signals in the spectral library
Figure GDA0001783923190000039
is the same as the hyperspectral signal to be unmixed
Figure GDA00017839231900000310
The bands in are consistent, and the spectral library endmembers after the bands are removed are obtained.
Figure GDA00017839231900000311
The number of bands is L>0, and the band is Band=[b 1 ,b 2 ,...,b L ];

步骤2-2，输入的光谱库为

Figure GDA00017839231900000312
波段数为K>0，包含的端元个数为N>0，对于光谱库
Figure GDA0001783923190000041
中的第i个端元信号
Figure GDA0001783923190000042
1≤i≤N，使用步骤2-1中的方法剔除其波段，得到剔除后的端元信号为
Figure GDA0001783923190000043
对光谱库
Figure GDA0001783923190000044
中的每个端元信号运用上述方法，即可得到剔除波段后的光谱库
Figure GDA0001783923190000045
Step 2-2, the input spectral library is
Figure GDA00017839231900000312
The number of bands is K>0, and the number of endmembers included is N>0. For the spectral library
Figure GDA0001783923190000041
The ith endmember signal in
Figure GDA0001783923190000042
1≤i≤N, use the method in step 2-1 to remove the band, and the removed endmember signal is
Figure GDA0001783923190000043
against spectral library
Figure GDA0001783923190000044
Using the above method for each endmember signal in
Figure GDA0001783923190000045

步骤2-3，初始输入的待解混高光谱数据为

Figure GDA0001783923190000046
将待解混高光谱数据
Figure GDA0001783923190000047
逐像素排列形成光谱数据矩阵
Figure GDA0001783923190000048
其中，整数L>0表示波段数，M＝W×H表示高光谱像元的个数，
Figure GDA0001783923190000049
表示一个高光谱像元，以光谱数据矩阵Y作为模型的一个输入。Step 2-3, the initial input hyperspectral data to be unmixed is
Figure GDA0001783923190000046
The hyperspectral data to be unmixed
Figure GDA0001783923190000047
Pixel-by-pixel arrangement to form a spectral data matrix
Figure GDA0001783923190000048
Among them, the integer L>0 represents the number of bands, M=W×H represents the number of hyperspectral pixels,
Figure GDA0001783923190000049
Represents a hyperspectral pixel, with the spectral data matrix Y as an input to the model.

进一步的，步骤3中构建待解混高光谱数据与光谱库矩阵之间的拟合误差保真项：Further, in step 3, a fitting error fidelity term between the hyperspectral data to be unmixed and the spectral library matrix is constructed:

Figure GDA00017839231900000410
Figure GDA00017839231900000410

式中

Figure GDA00017839231900000411
表示待解混的高光谱数据矩阵，
Figure GDA00017839231900000412
表示待裁剪的光谱库，
Figure GDA00017839231900000413
是初始化的丰度系数矩阵，其中L>0表示波段数，M>0 表示高光谱数据矩阵中光谱信号的个数，N>0表示光谱库所包含的端元光谱的数量。in the formula
Figure GDA00017839231900000411
represents the hyperspectral data matrix to be unmixed,
Figure GDA00017839231900000412
represents the spectral library to be tailored,
Figure GDA00017839231900000413
is the initialized abundance coefficient matrix, where L>0 represents the number of bands, M>0 represents the number of spectral signals in the hyperspectral data matrix, and N>0 represents the number of endmember spectra contained in the spectral library.

进一步的，步骤4中构建丰度矩阵协同稀疏性约束项具体为：Further, in step 4, the collaborative sparsity constraint term of constructing the abundance matrix is as follows:

Figure GDA00017839231900000414
Figure GDA00017839231900000414

式中

Figure GDA00017839231900000415
表示丰度系数矩阵，||·||2,1表示l2,1范数，Ui,j表示待解混高光谱数据矩阵
Figure GDA00017839231900000416
中第j个像元所对应的光谱库
Figure GDA00017839231900000417
中第i个端元的丰度系数， 1≤j≤M，1≤i≤N。in the formula
Figure GDA00017839231900000415
represents the abundance coefficient matrix, ||·|| 2,1 represents the l 2,1 norm, U i,j represents the hyperspectral data matrix to be unmixed
Figure GDA00017839231900000416
The spectral library corresponding to the jth pixel in
Figure GDA00017839231900000417
The abundance coefficient of the ith endmember in , 1≤j≤M, 1≤i≤N.

进一步的，步骤5中构建光谱库裁剪的稀疏正则化项具体为：Further, the sparse regularization term for constructing the spectral library clipping in step 5 is specifically:

Figure GDA00017839231900000418
Figure GDA00017839231900000418

式中

Figure GDA00017839231900000419
表示初始化后的光谱库，
Figure GDA00017839231900000420
表示待裁剪的光谱库，
Figure GDA00017839231900000421
表示
Figure GDA00017839231900000422
中的第p行，第q列所对应的值，
Figure GDA00017839231900000423
表示D中第p行，第q列所对应的值， 1≤p≤L，1≤q≤N，||·||1,1表示l1,1范数。in the formula
Figure GDA00017839231900000419
represents the initialized spectral library,
Figure GDA00017839231900000420
represents the spectral library to be tailored,
Figure GDA00017839231900000421
express
Figure GDA00017839231900000422
In the pth row, the value corresponding to the qth column,
Figure GDA00017839231900000423
Represents the value corresponding to the pth row and the qth column in D, 1≤p≤L, 1≤q≤N, ||·|| 1,1 represents the l 1,1 norm.

进一步的，步骤6中建立光谱库裁剪和协同稀疏回归模型具体为：Further, the spectral library tailoring and collaborative sparse regression model established in step 6 are as follows:

Figure GDA0001783923190000051
Figure GDA0001783923190000051

式中参数λ>0,β>0，

Figure GDA0001783923190000052
表示待解混的高光谱数据矩阵，
Figure GDA0001783923190000053
表示初始化后的光谱库，
Figure GDA0001783923190000054
表示待裁剪的光谱库，
Figure GDA0001783923190000055
表示丰度系数矩阵，||·||F表示矩阵的F范数，||·||2,1表示l2,1范数，||·||1,1表示l1,1范数。where the parameters λ>0, β>0,
Figure GDA0001783923190000052
represents the hyperspectral data matrix to be unmixed,
Figure GDA0001783923190000053
represents the initialized spectral library,
Figure GDA0001783923190000054
represents the spectral library to be tailored,
Figure GDA0001783923190000055
represents the abundance coefficient matrix, ||·|| F represents the F-norm of the matrix, ||·|| 2,1 represents the l 2,1 norm, ||·|| 1,1 represents the l 1,1 norm .

进一步的，步骤7中迭代光谱库裁剪和协同稀疏回归模型的具体步骤为：Further, the specific steps of iterative spectral library clipping and collaborative sparse regression model in step 7 are:

(1)首先输入初始化后的光谱库

Figure GDA0001783923190000056
和待解混的高光谱数据矩阵Y＝[y1,y2,...,yM]∈RL×M；(1) First input the initialized spectral library
Figure GDA0001783923190000056
and the hyperspectral data matrix to be unmixed Y=[y 1 , y 2 ,...,y M ]∈R L×M ;

(2)初始化丰度系数矩阵

Figure GDA0001783923190000057
以及待裁剪的光谱库
Figure GDA0001783923190000058
(2) Initialize the abundance coefficient matrix
Figure GDA0001783923190000057
and the spectral library to be tailored
Figure GDA0001783923190000058

(3)引入一个松弛变量

Figure GDA0001783923190000059
将上述光谱库裁剪和协同稀疏回归模型变成：(3) Introduce a slack variable
Figure GDA0001783923190000059
The above spectral library clipping and co-sparse regression model becomes:

Figure GDA00017839231900000510
Figure GDA00017839231900000510

式中的参数α>0,λ>0,β>0，对该模型引入等式约束U＝V1,U＝V2，然后求解其增广拉格朗日函数最小化模型：The parameters α>0, λ>0, β>0 in the formula, introduce equality constraints U=V 1 , U=V 2 to the model, and then solve its augmented Lagrangian function minimization model:

Figure GDA00017839231900000511
Figure GDA00017839231900000511

式中参数α>0,λ>0,β>0,μ>0，变量

Figure GDA00017839231900000512
使用交换方向乘子法求解后可以得到丰度系数矩阵
Figure GDA00017839231900000513
为：In the formula, the parameters α>0, λ>0, β>0, μ>0, the variable
Figure GDA00017839231900000512
The abundance coefficient matrix can be obtained after solving using the exchange direction multiplier method
Figure GDA00017839231900000513
for:

Figure GDA00017839231900000514
Figure GDA00017839231900000514

其中U(k+1)表示第k+1次迭代过程中的丰度矩阵U，V1 (k)表示第k次迭代过程中的变量V1，D1 (k)表示第k次迭代过程中的拉格朗日乘子变量D1，V2 (k)表示第k次迭代过程中的变量V2，D2 (k)表示第k次迭代过程中的拉格朗日乘子变量D2，求得待裁剪光谱库

Figure GDA0001783923190000061
为：where U (k+1) represents the abundance matrix U in the k+1th iteration, V 1 (k) represents the variable V 1 in the k-th iteration, and D 1 (k) represents the k-th iteration The Lagrangian multiplier variables in D 1 , V 2 (k) represent the variable V 2 in the k-th iteration process, and D 2 (k) represent the Lagrangian multiplier variable D in the k-th iteration process 2. Obtain the spectral library to be tailored
Figure GDA0001783923190000061
for:

Figure GDA0001783923190000062
Figure GDA0001783923190000062

其中，

Figure GDA0001783923190000063
表示第k+1次迭代过程中待裁剪光谱库
Figure GDA0001783923190000064
的第i列，D(:,i)表示初始化光谱库D中的第i列，H(:,i) (k+1)表示第k+1次迭代过程中松弛变量H的第i列，soft函数是软阈值收缩函数。in,
Figure GDA0001783923190000063
Indicates the spectral library to be tailored during the k+1 iteration
Figure GDA0001783923190000064
The i-th column of , D (:,i) represents the i-th column in the initialized spectral library D, H (:,i) (k+1) represents the i-th column of the slack variable H during the k+1-th iteration, The soft function is a soft threshold shrink function.

进一步的，步骤8输出裁剪后的光谱库和端元丰度图的具体过程为：Further, the specific process of outputting the cropped spectral library and endmember abundance map in step 8 is as follows:

(1)输出裁剪后的光谱库

Figure GDA0001783923190000065
(1) Output the cropped spectral library
Figure GDA0001783923190000065

(2)将丰度系数矩阵

Figure GDA0001783923190000066
重组为三维数据
Figure GDA0001783923190000067
然后输出
Figure GDA0001783923190000068
其中H表示矿物丰度图的高度，W表示矿物丰度图的宽度，N表示裁剪后光谱库
Figure GDA0001783923190000069
中的端元个数。(2) The abundance coefficient matrix
Figure GDA0001783923190000066
Reorganized into 3D data
Figure GDA0001783923190000067
then output
Figure GDA0001783923190000068
where H is the height of the mineral abundance map, W is the width of the mineral abundance map, and N is the spectral library after clipping
Figure GDA0001783923190000069
The number of endmembers in .

下面结合具体实施例对本发明进行详细说明。The present invention will be described in detail below with reference to specific embodiments.

实施例Example

结合图1，一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，包括以下步骤：Combined with Figure 1, a hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression includes the following steps:

步骤1，输入光谱库和待解混高光谱数据：输入一个包含高光谱端元的光谱库

Figure GDA00017839231900000610
一个光谱库指的是包含许多纯净物质(也称为端元)光谱信号的集合，其中K表示每个光谱信号的波段数，N表示该光谱库所包含的端元光谱信号的数量。在本实施例使用USGS光谱库作为输入光谱库。在USGS光谱库中，K＝224，N＝342。输入待解混的高光谱数据
Figure GDA00017839231900000611
其中L表示高光谱数据的波段数，W和H分别表示高光谱数据空间维的宽度和高度。以Cuprite数据集作为待解混的高光谱数据，在Cuprite数据集中L＝188，W＝191，H＝250。Step 1, input spectral library and hyperspectral data to be unmixed: input a spectral library containing hyperspectral endmembers
Figure GDA00017839231900000610
A spectral library refers to a collection of spectral signals containing many pure substances (also called endmembers), where K represents the number of bands for each spectral signal, and N represents the number of endmember spectral signals contained in the spectral library. The USGS spectral library is used as the input spectral library in this example. In the USGS spectral library, K=224, N=342. Input hyperspectral data to be unmixed
Figure GDA00017839231900000611
where L represents the number of bands of the hyperspectral data, and W and H represent the width and height of the spatial dimension of the hyperspectral data, respectively. Taking the Cuprite dataset as the hyperspectral data to be unmixed, L=188, W=191, and H=250 in the Cuprite dataset.

步骤2，光谱库初始化与光谱数据矩阵构造：Step 2, spectral library initialization and spectral data matrix construction:

(1)输入一个光谱库端元

Figure GDA00017839231900000612
波段数为K>0，波段为
Figure GDA00017839231900000613
而输入一个待解混高光谱信号
Figure GDA0001783923190000071
波段数为L>0，波段为Band＝[b1,b2,...,bL]。根据待解混高光谱信号
Figure GDA0001783923190000072
的波段Band，保留光谱库端元
Figure GDA0001783923190000073
波段
Figure GDA0001783923190000074
中对应的波段，其余的波段均丢弃。USGS光谱中的端元信号波段范围为1～224，总计224个波段，而在Cuprite数据集中总计有118个波段。在Cuprite数据集中由于水分吸收和信噪比较低等原因，1-2，105-115，150-170，223-224等波段的值都被剔除，所以USGS光谱库端元波段
Figure GDA0001783923190000075
中也要将对应的1-2，105-115，150-170，223-224等波段剔除。经过这样的一步操作，光谱库USGS中端元信号的波段就和高光谱数据Cuprite中的波段保持一致。(1) Input a spectral library endmember
Figure GDA00017839231900000612
The number of bands is K>0, and the band is
Figure GDA00017839231900000613
And input a hyperspectral signal to be unmixed
Figure GDA0001783923190000071
The number of bands is L>0, and the band is Band=[b 1 ,b 2 ,...,b L ]. According to the hyperspectral signal to be unmixed
Figure GDA0001783923190000072
Band, retains spectral library endmembers
Figure GDA0001783923190000073
band
Figure GDA0001783923190000074
The corresponding band in the , and the rest of the bands are discarded. The endmember signal bands in the USGS spectrum range from 1 to 224, with a total of 224 bands, while there are a total of 118 bands in the Cuprite dataset. In the Cuprite dataset, due to the low moisture absorption and signal-to-noise ratio, the values in the bands 1-2, 105-115, 150-170, 223-224 are all excluded, so the endmember bands of the USGS spectral library are
Figure GDA0001783923190000075
The corresponding bands such as 1-2, 105-115, 150-170, 223-224 should also be eliminated. After such a one-step operation, the band of the endmember signal in the spectral library USGS is consistent with the band in the hyperspectral data Cuprite.

(2)输入的光谱库为

Figure GDA0001783923190000076
波段数为K>0，包含的端元个数为N>0，对于光谱库
Figure GDA0001783923190000077
中的第i个端元信号
Figure GDA0001783923190000078
1≤i≤N，使用(1)中的方法剔除其波段，得到剔除后的端元信号为
Figure GDA0001783923190000079
对光谱库
Figure GDA00017839231900000710
中的每个端元信号运用上述方法，即可得到剔除波段后(即初始化)的光谱库
Figure GDA00017839231900000711
经过这一步即可得到剔除波段后USGS光谱库。(2) The input spectral library is
Figure GDA0001783923190000076
The number of bands is K>0, and the number of endmembers included is N>0. For the spectral library
Figure GDA0001783923190000077
The ith endmember signal in
Figure GDA0001783923190000078
1≤i≤N, use the method in (1) to remove its band, and the removed endmember signal is
Figure GDA0001783923190000079
against spectral library
Figure GDA00017839231900000710
Using the above method for each endmember signal in , the spectral library after the band is eliminated (that is, initialized)
Figure GDA00017839231900000711
After this step, the USGS spectral library after the band is eliminated can be obtained.

(3)初始输入的待解混高光谱数据为

Figure GDA00017839231900000712
其中，整数L>0表示待解混高光谱数据的波段数，整数W>0和整数H>0分别表示待解混高光谱数据空间维的宽度和高度。将待解混高光谱数据
Figure GDA00017839231900000713
逐像素排列形成光谱数据矩阵
Figure GDA00017839231900000714
其中，整数L>0表示波段数，M＝W×H表示高光谱像元的个数，
Figure GDA00017839231900000715
表示一个高光谱像元，以光谱数据矩阵Y作为模型的一个输入。Cuprite数据集
Figure GDA00017839231900000716
中，L＝188，W＝191，H＝250，将Cuprite数据集逐像素排列后，得到光谱数据矩阵
Figure GDA00017839231900000717
L＝188,M＝191×250＝47750。(3) The initial input hyperspectral data to be unmixed is
Figure GDA00017839231900000712
The integer L>0 represents the number of bands of the hyperspectral data to be unmixed, and the integer W>0 and the integer H>0 represent the width and height of the spatial dimension of the hyperspectral data to be unmixed, respectively. The hyperspectral data to be unmixed
Figure GDA00017839231900000713
Pixel-by-pixel arrangement to form a spectral data matrix
Figure GDA00017839231900000714
Among them, the integer L>0 represents the number of bands, M=W×H represents the number of hyperspectral pixels,
Figure GDA00017839231900000715
Represents a hyperspectral pixel, with the spectral data matrix Y as an input to the model. Cuprite dataset
Figure GDA00017839231900000716
, L=188, W=191, H=250, after arranging the Cuprite data set pixel by pixel, the spectral data matrix is obtained
Figure GDA00017839231900000717
L=188, M=191×250=47750.

步骤3，构建光谱拟合误差保真项：Step 3, construct the spectral fitting error fidelity term:

Figure GDA00017839231900000718
Figure GDA00017839231900000718

式中

Figure GDA00017839231900000719
表示待解混的高光谱数据矩阵，
Figure GDA00017839231900000720
表示待裁剪的光谱库，
Figure GDA00017839231900000721
是初始化的丰度系数矩阵，其中L>0表示波段数，M>0 表示高光谱数据矩阵中光谱信号的个数，N>0表示光谱库所包含的端元光谱的数量。in the formula
Figure GDA00017839231900000719
represents the hyperspectral data matrix to be unmixed,
Figure GDA00017839231900000720
represents the spectral library to be tailored,
Figure GDA00017839231900000721
is the initialized abundance coefficient matrix, where L>0 represents the number of bands, M>0 represents the number of spectral signals in the hyperspectral data matrix, and N>0 represents the number of endmember spectra contained in the spectral library.

步骤4，构建丰度矩阵协同稀疏性约束项：Step 4, construct the abundance matrix synergistic sparsity constraint:

Figure GDA0001783923190000081
Figure GDA0001783923190000081

式中

Figure GDA0001783923190000082
表示丰度系数矩阵，||·||2,1表示l2,1范数，Ui,j表示待解混高光谱数据矩阵
Figure GDA0001783923190000083
中第j个像元所对应光谱库
Figure GDA0001783923190000084
中第i个端元的丰度系数，1≤j≤M，1≤i≤N，N>0表示光谱库中端元光谱的个数，M>0表示高光谱数据矩阵中像元的个数。in the formula
Figure GDA0001783923190000082
represents the abundance coefficient matrix, ||·|| 2,1 represents the l 2,1 norm, U i,j represents the hyperspectral data matrix to be unmixed
Figure GDA0001783923190000083
The spectral library corresponding to the jth pixel in
Figure GDA0001783923190000084
The abundance coefficient of the i-th endmember in , 1≤j≤M, 1≤i≤N, N>0 indicates the number of endmember spectra in the spectral library, M>0 indicates the number of pixels in the hyperspectral data matrix number.

步骤5，建立光谱库裁剪的稀疏正则化项：Step 5, establish a sparse regularization term for spectral library clipping:

Figure GDA0001783923190000085
Figure GDA0001783923190000085

式中

Figure GDA0001783923190000086
表示初始化后的光谱库，
Figure GDA0001783923190000087
表示待裁剪的光谱库，
Figure GDA0001783923190000088
表示
Figure GDA0001783923190000089
中的第p行，第q列所对应的值，
Figure GDA00017839231900000810
表示D中第p行，第q列所对应的值， 1≤p≤L，1≤q≤N，||·||1,1表示l1,1范数。其中L>0表示每个光谱信号的波段数，N>0 表示光谱库中端元光谱的个数。in the formula
Figure GDA0001783923190000086
represents the initialized spectral library,
Figure GDA0001783923190000087
represents the spectral library to be tailored,
Figure GDA0001783923190000088
express
Figure GDA0001783923190000089
In the pth row, the value corresponding to the qth column,
Figure GDA00017839231900000810
Represents the value corresponding to the pth row and the qth column in D, 1≤p≤L, 1≤q≤N, ||·|| 1,1 represents the l 1,1 norm. Where L>0 represents the number of bands of each spectral signal, and N>0 represents the number of endmember spectra in the spectral library.

步骤6，求解光谱库裁剪和协同稀疏回归模型：Step 6, solve the spectral library clipping and collaborative sparse regression model:

Figure GDA00017839231900000811
Figure GDA00017839231900000811

式中参数λ>0,β>0，

Figure GDA00017839231900000812
表示待解混的高光谱数据矩阵，
Figure GDA00017839231900000813
表示初始化后的光谱库，
Figure GDA00017839231900000814
表示待裁剪的光谱库，
Figure GDA00017839231900000815
表示丰度系数矩阵，||·||F表示矩阵的F范数，||·||2,1表示l2,1范数，||·||1,1表示l1,1范数。where the parameters λ>0, β>0,
Figure GDA00017839231900000812
represents the hyperspectral data matrix to be unmixed,
Figure GDA00017839231900000813
represents the initialized spectral library,
Figure GDA00017839231900000814
represents the spectral library to be tailored,
Figure GDA00017839231900000815
represents the abundance coefficient matrix, ||·|| F represents the F-norm of the matrix, ||·|| 2,1 represents the l 2,1 norm, ||·|| 1,1 represents the l 1,1 norm .

步骤7，求解光谱库裁剪和协同稀疏回归模型：Step 7, solve the spectral library clipping and collaborative sparse regression model:

(1)首先输入初始化后的光谱库

Figure GDA00017839231900000816
和待解混的高光谱数据矩阵Y＝[y1,y2,...,yM]∈RL×M，其中L>0表示每个光谱信号的波段数，N>0表示光谱库中端元光谱的个数，M>0表示待解混高光谱数据矩阵中像元的个数。(1) First input the initialized spectral library
Figure GDA00017839231900000816
and the hyperspectral data matrix to be unmixed Y=[y 1 , y 2 ,...,y M ]∈R L×M , where L>0 represents the number of bands of each spectral signal, and N>0 represents the spectral library The number of mid-endmember spectra, where M>0 represents the number of pixels in the hyperspectral data matrix to be unmixed.

(2)初始化丰度系数矩阵

Figure GDA00017839231900000817
以及待裁剪的光谱库
Figure GDA00017839231900000818
(2) Initialize the abundance coefficient matrix
Figure GDA00017839231900000817
and the spectral library to be tailored
Figure GDA00017839231900000818

(3)引入一个松弛变量

Figure GDA0001783923190000091
将上述光谱库裁剪和协同稀疏回归模型变成：(3) Introduce a slack variable
Figure GDA0001783923190000091
The above spectral library clipping and co-sparse regression model becomes:

Figure GDA0001783923190000092
Figure GDA0001783923190000092

式中的参数α>0,λ>0,β>0。对该模型引入等式约束U＝V1,U＝V2，然后求解其增广拉格朗日函数最小化模型：The parameters in the formula α>0, λ>0, β>0. Introduce equality constraints U=V 1 , U=V 2 to the model, and then solve its augmented Lagrangian function minimization model:

Figure GDA0001783923190000093
Figure GDA0001783923190000093

式中参数α>0,λ>0,β>0,μ>0，变量

Figure GDA0001783923190000094
使用交换方向乘子法求解后可以得到丰度系数矩阵
Figure GDA0001783923190000095
为：In the formula, the parameters α>0, λ>0, β>0, μ>0, the variable
Figure GDA0001783923190000094
The abundance coefficient matrix can be obtained after solving using the exchange direction multiplier method
Figure GDA0001783923190000095
for:

Figure GDA0001783923190000096
Figure GDA0001783923190000096

其中U(k+1)表示第k+1次迭代过程中的丰度矩阵U，V1 (k)表示第k次迭代过程中的变量V1，D1 (k)表示第k次迭代过程中的拉格朗日乘子变量D1，V2 (k)表示第k次迭代过程中的变量V2，D2 (k)表示第k次迭代过程中的拉格朗日乘子变量D2。求得待裁剪光谱库

Figure GDA0001783923190000097
为：where U (k+1) represents the abundance matrix U in the k+1th iteration, V 1 (k) represents the variable V 1 in the k-th iteration, and D 1 (k) represents the k-th iteration The Lagrangian multiplier variables in D 1 , V 2 (k) represent the variable V 2 in the k-th iteration process, and D 2 (k) represent the Lagrangian multiplier variable D in the k-th iteration process 2 . Obtain the spectral library to be tailored
Figure GDA0001783923190000097
for:

Figure GDA0001783923190000098
Figure GDA0001783923190000098

其中

Figure GDA0001783923190000099
表示第k+1次迭代过程中待裁剪光谱库
Figure GDA00017839231900000910
的第i列，D(:,i)表示初始化光谱库D中的第i列，H(:,i) (k+1)表示第k+1次迭代过程中松弛变量H的第i列，soft函数是著名的软阈值收缩函数。in
Figure GDA0001783923190000099
Indicates the spectral library to be tailored during the k+1 iteration
Figure GDA00017839231900000910
The i-th column of , D (:,i) represents the i-th column in the initialized spectral library D, H (:,i) (k+1) represents the i-th column of the slack variable H during the k+1-th iteration, The soft function is the well-known soft threshold shrink function.

步骤8，输出裁剪后的光谱库和端元丰度图：Step 8, output the cropped spectral library and endmember abundance map:

(1)输出裁剪后的光谱库

Figure GDA00017839231900000911
(1) Output the cropped spectral library
Figure GDA00017839231900000911

(2)将丰度系数矩阵

Figure GDA00017839231900000912
重组为三维数据
Figure GDA00017839231900000913
然后输出
Figure GDA00017839231900000914
其中H表示矿物丰度图的高度，W表示矿物丰度图的宽度，N表示裁剪后光谱库
Figure GDA0001783923190000101
中的端元个数。(2) The abundance coefficient matrix
Figure GDA00017839231900000912
Reorganized into 3D data
Figure GDA00017839231900000913
then output
Figure GDA00017839231900000914
where H is the height of the mineral abundance map, W is the width of the mineral abundance map, and N is the spectral library after clipping
Figure GDA0001783923190000101
The number of endmembers in .

本实施例仿真实验采用的真实高光谱数据集是航空可见光/红外成像仪(VisibleInfraRed Imaging Spectrometer，AVIRIS)拍摄的Cuprite数据集。实验中使用的数据是Cuprite数据的一部分，大小为250×191个像元，包含400到2500纳米之间的224个波段，光谱分辨率为10纳米。在进行本次实验前，波段1-2，105-115，150-170，223-224 由于水汽吸收以及噪声较大的原因被去除，剩余188个可用波段。仿真实验均在Windows 7操作系统下采用MATLAB R2015b完成。The real hyperspectral data set used in the simulation experiment in this embodiment is the Cuprite data set captured by an aerial visible light/infrared imager (Visible Infra Red Imaging Spectrometer, AVIRIS). The data used in the experiments are part of Cuprite data, 250 × 191 pixels in size, containing 224 bands between 400 and 2500 nanometers, with a spectral resolution of 10 nanometers. Before this experiment, bands 1-2, 105-115, 150-170, 223-224 were removed due to water vapor absorption and high noise, leaving 188 available bands. The simulation experiments are all completed under the Windows 7 operating system using MATLAB R2015b.

本发明采用真实高光谱数据集检验算法的解混效果。为测试本发明算法的性能，将提出的基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法(PSCSR)与目前国际上流行的解混算法对比。对比方法包括：Tetracorder软件内置的解混方法，基于变量分裂和增广拉格朗日的稀疏解混(SUnSAL)[Bioucas-Dias J M,Figueiredo M A T.Alternatingdirection algorithms for constrained sparse regression:Application tohyperspectral unmixing[C].Hyperspectral Image and Signal Processing:Evolutionin Remote Sensing (WHISPERS),2010 2nd Workshop on.IEEE,2010:1-4.]，基于变量分裂和增广拉格朗日的协同稀疏解混(CLSUnSAL)[Iordache M D,Bioucas-Dias J M,PlazaA.Collaborative sparse regression for hyperspectral unmixing[J].IEEETransactions on geoscience and remote sensing,2014,52(1):341-354]和基于裁剪字典的非凸鼓励稀疏回归框架(DANSER) [Fu X,Ma W K,Bioucas-Dias J M,etal.Semiblind hyperspectral unmixing in the presence of spectral librarymismatches[J].IEEE Transactions on Geoscience and Remote Sensing, 2016,54(9):5171-5184.]The present invention adopts the real hyperspectral data set to check the unmixing effect of the algorithm. In order to test the performance of the algorithm of the present invention, the proposed hyperspectral data unmixing method (PSCSR) based on spectral library clipping and cooperative sparse regression is compared with the currently popular unmixing algorithms in the world. The comparison methods include: unmixing methods built into Tetracorder software, sparse unmixing based on variable splitting and augmented Lagrangian (SUnSAL) [Bioucas-Dias J M, Figueiredo M A T. Alternating direction algorithms for constrained sparse regression: Application tohyperspectral unmixing[ C]. Hyperspectral Image and Signal Processing: Evolution in Remote Sensing (WHISPERS), 2010 2nd Workshop on. IEEE, 2010: 1-4.], Cooperative Sparse Unmixing Based on Variable Splitting and Augmented Lagrangian (CLSUnSAL)[ Iordache M D, Bioucas-Dias J M, PlazaA. Collaborative sparse regression for hyperspectral unmixing [J]. IEEE Transactions on geoscience and remote sensing, 2014, 52(1): 341-354] and a non-convex encouraged sparse regression framework based on clipping dictionary ( DANSER) [Fu X,Ma W K,Bioucas-Dias J M,etal.Semiblind hyperspectral unmixing in the presence of spectral librarymismatches[J].IEEE Transactions on Geoscience and Remote Sensing, 2016,54(9):5171-5184.]

图2(a)～图2(o)为Cuprite数据集在不同解混算法下的三种矿物明矾石、水铵长石、玉髓的丰度图。可以看到本发明提出的解混算法得到的三种矿物的丰度图与使用Tetracorder软件得到的参考图基本一致，而且估计出的玉髓丰度图比其他三种解混算法(SUnSAL,CLSUnSAL,DANSER)更接近参考图，这证明本发明的解混算法能对真实的地物分布提供一个准确的估计。图3(a)～图3(c)为明矾石、水铵长石、玉髓三种矿物在初始光谱库和裁剪后光谱库中的光谱曲线对比。由图3(a)可以看出明矾石光谱曲线与裁剪过后的明矾石光谱曲线的不同之处大多发生在吸收峰，也就是说同一种物质的光谱变化具有一定稀疏性。Figures 2(a) to 2(o) are the abundance maps of the three minerals alumite, feldspar, and chalcedony in the Cuprite dataset under different unmixing algorithms. It can be seen that the abundance map of the three minerals obtained by the unmixing algorithm proposed in the present invention is basically consistent with the reference map obtained by using the Tetracorder software, and the estimated chalcedony abundance map is better than the other three unmixing algorithms (SUnSAL, CLSUnSAL). , DANSER) is closer to the reference map, which proves that the unmixing algorithm of the present invention can provide an accurate estimate for the real distribution of ground objects. Figures 3(a) to 3(c) show the comparison of spectral curves of three minerals, alumite, feldspar, and chalcedony in the initial spectral library and the clipped spectral library. It can be seen from Figure 3(a) that the difference between the spectral curve of alumite and the spectral curve of the cut alumite mostly occurs in the absorption peak, that is to say, the spectral changes of the same substance have a certain sparseness.

Claims (3)

1.一种基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，其特征在于，包含如下步骤：1. a hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression, is characterized in that, comprises the following steps: 步骤1，输入光谱库和待解混高光谱数据；具体过程为：Step 1, input the spectral library and the hyperspectral data to be unmixed; the specific process is: 输入一个包含高光谱端元的光谱库
Figure FDA0003288269660000011
其中整数K＞0表示每个光谱信号的波段数，整数N>0表示该光谱库所包含的端元光谱的数量；Enter a spectral library containing hyperspectral endmembers
Figure FDA0003288269660000011
The integer K>0 represents the number of bands of each spectral signal, and the integer N>0 represents the number of endmember spectra contained in the spectral library; 输入待解混的高光谱数据
Figure FDA0003288269660000012
其中，整数L>0表示高光谱数据的波段数，整数W>0和整数H＞0分别表示高光谱数据空间维的宽度和高度；Input hyperspectral data to be unmixed
Figure FDA0003288269660000012
Among them, the integer L>0 represents the number of bands of the hyperspectral data, and the integer W>0 and the integer H>0 represent the width and height of the spatial dimension of the hyperspectral data, respectively; 步骤2，光谱库初始化与光谱数据矩阵构造；具体步骤为：Step 2, spectral library initialization and spectral data matrix construction; the specific steps are: 步骤2-1，输入一个光谱库端元
Figure FDA0003288269660000013
波段数为K＞0，波段为
Figure FDA0003288269660000014
而输入一个待解混高光谱信号
Figure FDA0003288269660000015
波段数为L>0，波段为Band＝[b1,b2,...,bL]；根据待解混高光谱信号
Figure FDA0003288269660000016
的波段Band，保留光谱库端元
Figure FDA0003288269660000017
波段
Figure FDA0003288269660000018
中对应的波段，其余的波段均丢弃；经过这样的一步操作，光谱库中端元信号
Figure FDA0003288269660000019
的波段和待解混高光谱信号
Figure FDA00032882696600000110
中的波段保持一致，得到剔除波段后的光谱库端元
Figure FDA00032882696600000111
波段数为L>0，波段为Band＝[b1,b2,...,bL]；Step 2-1, enter a spectral library endmember
Figure FDA0003288269660000013
The number of bands is K>0, and the band is
Figure FDA0003288269660000014
And input a hyperspectral signal to be unmixed
Figure FDA0003288269660000015
The number of bands is L>0, and the band is Band=[b 1 ,b 2 ,...,b L ]; according to the hyperspectral signal to be unmixed
Figure FDA0003288269660000016
Band, retains spectral library endmembers
Figure FDA0003288269660000017
band
Figure FDA0003288269660000018
The corresponding bands in the spectral library are discarded, and the rest of the bands are discarded; after such a one-step operation, the endmember signals in the spectral library
Figure FDA0003288269660000019
and the hyperspectral signal to be unmixed
Figure FDA00032882696600000110
The bands in are consistent, and the spectral library endmembers after the bands are removed are obtained.
Figure FDA00032882696600000111
The number of bands is L>0, and the band is Band=[b 1 ,b 2 ,...,b L ]; 步骤2-2，输入的光谱库为
Figure FDA00032882696600000112
波段数为K＞0，包含的端元个数为N＞0，对于光谱库
Figure FDA00032882696600000113
中的第i个端元信号
Figure FDA00032882696600000114
使用步骤2-1中的方法剔除其波段，得到剔除后的端元信号为
Figure FDA00032882696600000115
对光谱库
Figure FDA00032882696600000116
中的每个端元信号运用上述方法，即可得到初始化后的光谱库
Figure FDA00032882696600000117
Step 2-2, the input spectral library is
Figure FDA00032882696600000112
The number of bands is K>0, and the number of endmembers included is N>0. For the spectral library
Figure FDA00032882696600000113
The ith endmember signal in
Figure FDA00032882696600000114
Use the method in step 2-1 to remove the band, and the removed endmember signal is
Figure FDA00032882696600000115
against spectral library
Figure FDA00032882696600000116
Using the above method for each endmember signal in , the initialized spectral library can be obtained
Figure FDA00032882696600000117
 步骤2-3，初始输入的待解混高光谱数据为
Figure FDA00032882696600000118
将待解混高光谱数据
Figure FDA00032882696600000119
逐像素排列形成光谱数据矩阵
Figure FDA00032882696600000120
其中，整数L>0表示波段数，M＝W×H表示高光谱像元的个数，
Figure FDA00032882696600000121
表示一个高光谱像元，以光谱数据矩阵Y作为模型的一个输入；Step 2-3, the initial input hyperspectral data to be unmixed is
Figure FDA00032882696600000118
The hyperspectral data to be unmixed
Figure FDA00032882696600000119
Pixel-by-pixel arrangement to form a spectral data matrix
Figure FDA00032882696600000120
Among them, the integer L>0 represents the number of bands, M=W×H represents the number of hyperspectral pixels,
Figure FDA00032882696600000121
Represents a hyperspectral pixel, using the spectral data matrix Y as an input to the model; 步骤3，构建待解混高光谱数据与光谱库矩阵之间的拟合误差保真项：Step 3, construct the fitting error fidelity term between the hyperspectral data to be unmixed and the spectral library matrix:
Figure FDA0003288269660000021
Figure FDA0003288269660000021
 式中
Figure FDA0003288269660000022
表示待解混的高光谱数据矩阵，
Figure FDA0003288269660000023
表示待裁剪的光谱库，
Figure FDA0003288269660000024
是丰度系数矩阵，其中L>0表示波段数，M＞0表示高光谱数据矩阵中光谱信号的个数；in the formula
Figure FDA0003288269660000022
represents the hyperspectral data matrix to be unmixed,
Figure FDA0003288269660000023
represents the spectral library to be tailored,
Figure FDA0003288269660000024
is the abundance coefficient matrix, where L>0 represents the number of bands, and M>0 represents the number of spectral signals in the hyperspectral data matrix; 步骤4，构建丰度矩阵协同稀疏性约束项；具体为：Step 4: Construct the synergistic sparsity constraint term of the abundance matrix; specifically:
Figure FDA0003288269660000025
Figure FDA0003288269660000025
 式中
Figure FDA0003288269660000026
表示丰度系数矩阵，||·||2,1表示l2,1范数，Ui,j表示待解混高光谱数据矩阵
Figure FDA0003288269660000027
中第j个像元所对应的光谱库
Figure FDA0003288269660000028
中第i个端元的丰度系数，1≤j≤M，1≤i≤N；in the formula
Figure FDA0003288269660000026
represents the abundance coefficient matrix, ||·|| 2,1 represents the l 2,1 norm, U i,j represents the hyperspectral data matrix to be unmixed
Figure FDA0003288269660000027
The spectral library corresponding to the jth pixel in
Figure FDA0003288269660000028
The abundance coefficient of the ith endmember in , 1≤j≤M, 1≤i≤N; 步骤5，构建光谱库裁剪的稀疏正则化项；具体为：Step 5, construct a sparse regularization term for spectral library clipping; specifically:
Figure FDA0003288269660000029
Figure FDA0003288269660000029
 式中
Figure FDA00032882696600000210
表示初始化后的光谱库，
Figure FDA00032882696600000211
表示待裁剪的光谱库，
Figure FDA00032882696600000212
表示
Figure FDA00032882696600000213
中的第p行，第q列所对应的值，
Figure FDA00032882696600000214
表示D中第p行，第q列所对应的值，1≤p≤L，1≤q≤N，||·||1,1表示l1,1范数；in the formula
Figure FDA00032882696600000210
represents the initialized spectral library,
Figure FDA00032882696600000211
represents the spectral library to be tailored,
Figure FDA00032882696600000212
express
Figure FDA00032882696600000213
In the pth row, the value corresponding to the qth column,
Figure FDA00032882696600000214
Represents the value corresponding to the pth row and the qth column in D, 1≤p≤L, 1≤q≤N, ||·|| 1,1 represents the l 1,1 norm; 步骤6，建立光谱库裁剪和协同稀疏回归模型；具体为：Step 6: Establish spectral library clipping and collaborative sparse regression model; specifically:
Figure FDA00032882696600000215
Figure FDA00032882696600000215
 式中参数λ＞0,β＞0，
Figure FDA00032882696600000216
表示待解混的高光谱数据矩阵，
Figure FDA00032882696600000217
表示初始化后的光谱库，
Figure FDA00032882696600000218
表示待裁剪的光谱库，
Figure FDA00032882696600000219
表示丰度系数矩阵，||·||F表示矩阵的F范数，||·||2,1表示l2,1范数，||·||1,1表示l1,1范数；In the formula, the parameters λ>0, β>0,
Figure FDA00032882696600000216
represents the hyperspectral data matrix to be unmixed,
Figure FDA00032882696600000217
represents the initialized spectral library,
Figure FDA00032882696600000218
represents the spectral library to be tailored,
Figure FDA00032882696600000219
represents the abundance coefficient matrix, ||·|| F represents the F-norm of the matrix, ||·|| 2,1 represents the l 2,1 norm, ||·|| 1,1 represents the l 1,1 norm ; 步骤7，迭代求解光谱库裁剪和协同稀疏回归模型；Step 7, iteratively solve the spectral library clipping and collaborative sparse regression model; 步骤8，输出裁剪后的光谱库和端元丰度图。Step 8, output the cropped spectral library and endmember abundance map. 2.根据权利要求1所述的基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，其特征在于，步骤7中迭代光谱库裁剪和协同稀疏回归模型的具体步骤为：2. the hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression according to claim 1, is characterized in that, in step 7, the concrete steps of iterative spectral library clipping and collaborative sparse regression model are: (1)输入初始化后的光谱库
Figure FDA0003288269660000031
和待解混的高光谱数据矩阵Y＝[y1,y2,...,yM]∈RL×M；(1) Input the initialized spectral library
Figure FDA0003288269660000031
and the hyperspectral data matrix to be unmixed Y=[y 1 , y 2 ,...,y M ]∈R L×M ; (2)初始化丰度系数矩阵
Figure FDA0003288269660000032
以及待裁剪的光谱库
Figure FDA0003288269660000033
(2) Initialize the abundance coefficient matrix
Figure FDA0003288269660000032
and the spectral library to be tailored
Figure FDA0003288269660000033
 (3)引入一个松弛变量
Figure FDA0003288269660000034
将上述光谱库裁剪和协同稀疏回归模型变成：(3) Introduce a slack variable
Figure FDA0003288269660000034
The above spectral library clipping and co-sparse regression model becomes:
Figure FDA0003288269660000035
Figure FDA0003288269660000035
 式中的参数α＞0,λ＞0,β＞0，对该模型引入等式约束U＝V1,U＝V2，然后求解其增广拉格朗日函数最小化模型：In the formula, the parameters α>0, λ>0, β>0, introduce equality constraints U=V 1 , U=V 2 to the model, and then solve its augmented Lagrangian function minimization model:
Figure FDA0003288269660000036
Figure FDA0003288269660000036
 式中参数μ＞0，变量
Figure FDA0003288269660000037
使用交换方向乘子法求解后可以得到丰度系数矩阵
Figure FDA0003288269660000038
为：where the parameter μ＞0, the variable
Figure FDA0003288269660000037
The abundance coefficient matrix can be obtained after solving using the exchange direction multiplier method
Figure FDA0003288269660000038
for:
Figure FDA0003288269660000039
Figure FDA0003288269660000039
 其中U(k+1)表示第k+1次迭代过程中的丰度矩阵U，V1 (k)表示第k次迭代过程中的变量V1，D1 (k)表示第k次迭代过程中的拉格朗日乘子变量D1，V2 (k)表示第k次迭代过程中的变量V2，D2 (k)表示第k次迭代过程中的拉格朗日乘子变量D2，求得待裁剪光谱库
Figure FDA00032882696600000310
为：where U (k+1) represents the abundance matrix U in the k+1th iteration, V 1 (k) represents the variable V 1 in the k-th iteration, and D 1 (k) represents the k-th iteration The Lagrangian multiplier variables in D 1 , V 2 (k) represent the variable V 2 in the k-th iteration process, and D 2 (k) represent the Lagrangian multiplier variable D in the k-th iteration process 2. Obtain the spectral library to be tailored
Figure FDA00032882696600000310
for:
Figure FDA00032882696600000311
Figure FDA00032882696600000311
 其中，
Figure FDA00032882696600000312
表示第k+1次迭代过程中待裁剪光谱库
Figure FDA00032882696600000313
的第i列，D(:,i)表示初始化光谱库D中的第i列，H(:,i) (k+1)表示第k+1次迭代过程中松弛变量H的第i列，soft函数是软阈值收缩函数。in,
Figure FDA00032882696600000312
Indicates the spectral library to be tailored during the k+1 iteration
Figure FDA00032882696600000313
The i-th column of , D (:,i) represents the i-th column in the initialized spectral library D, H (:,i) (k+1) represents the i-th column of the slack variable H during the k+1-th iteration, The soft function is a soft threshold shrink function. 3.根据权利要求2所述的基于光谱库裁剪与协同稀疏回归的高光谱数据解混方法，其特征在于，步骤8输出裁剪后的光谱库和端元丰度图的具体过程为：3. the hyperspectral data unmixing method based on spectral library clipping and collaborative sparse regression according to claim 2, is characterized in that, the concrete process of step 8 output clipping spectral library and endmember abundance map is: (1)输出裁剪后的光谱库
Figure FDA0003288269660000041
(1) Output the cropped spectral library
Figure FDA0003288269660000041
 (2)使用经裁剪的光谱库计算得到结果丰度系数矩阵
Figure FDA0003288269660000042
将其重组为三维数据
Figure FDA0003288269660000043
然后输出
Figure FDA0003288269660000044
其中H表示矿物丰度图的高度，W表示矿物丰度图的宽度，
Figure FDA0003288269660000045
表示裁剪后光谱库
Figure FDA0003288269660000046
中的端元个数。(2) Calculate the resulting abundance coefficient matrix using the tailored spectral library
Figure FDA0003288269660000042
Restructure it into 3D data
Figure FDA0003288269660000043
then output
Figure FDA0003288269660000044
where H represents the height of the mineral abundance map, W represents the width of the mineral abundance map,
Figure FDA0003288269660000045
Represents the clipped spectral library
Figure FDA0003288269660000046
The number of endmembers in .
CN201810016715.4A 2018-01-08 2018-01-08 Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression Active CN108734672B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810016715.4A CN108734672B (en) 2018-01-08 2018-01-08 Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810016715.4A CN108734672B (en) 2018-01-08 2018-01-08 Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression

Publications (2)

Publication Number Publication Date
CN108734672A CN108734672A (en) 2018-11-02
CN108734672B true CN108734672B (en) 2022-02-18

Family

ID=63940371

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810016715.4A Active CN108734672B (en) 2018-01-08 2018-01-08 Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression

Country Status (1)

Country Link
CN (1) CN108734672B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109580495B (en) * 2018-11-28 2021-08-24 重庆工商大学 Demixing device and method based on hyperspectral image
CN111260576B (en) * 2020-01-14 2022-07-05 哈尔滨工业大学 Hyperspectral unmixing algorithm based on de-noising three-dimensional convolution self-coding network
CN113094645B (en) * 2021-03-11 2024-07-30 南京理工大学 Hyperspectral data unmixing method based on morphological component constraint optimization

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104463223A (en) * 2014-12-22 2015-03-25 西安电子科技大学 Hyperspectral image group sparse demixing method based on empty spectral information abundance restraint
CN105825227A (en) * 2016-03-11 2016-08-03 南京航空航天大学 Hyperspectral image sparseness demixing method based on MFOCUSS and low-rank expression
US20160371563A1 (en) * 2015-06-22 2016-12-22 The Johns Hopkins University System and method for structured low-rank matrix factorization: optimality, algorithm, and applications to image processing

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104463223A (en) * 2014-12-22 2015-03-25 西安电子科技大学 Hyperspectral image group sparse demixing method based on empty spectral information abundance restraint
US20160371563A1 (en) * 2015-06-22 2016-12-22 The Johns Hopkins University System and method for structured low-rank matrix factorization: optimality, algorithm, and applications to image processing
CN105825227A (en) * 2016-03-11 2016-08-03 南京航空航天大学 Hyperspectral image sparseness demixing method based on MFOCUSS and low-rank expression

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Joint Local Abundance Sparse Unmixing for Hyperspectral Images;Mia Rizkinia 等;《Remote Sensing》;20171127;第1-22页 *
Semiblind Hyperspectral Unmixing in the Presence of Spectral Library Mismatches;Xiao Fu 等;《IEEE Transactions on Geoscience and Remote Sensing》;20160510;第5171-5184页 *

Also Published As

Publication number Publication date
CN108734672A (en) 2018-11-02

Similar Documents

Publication Publication Date Title
CN111260576B (en) Hyperspectral unmixing algorithm based on de-noising three-dimensional convolution self-coding network
Zhong et al. Iterative support vector machine for hyperspectral image classification
CN110084159B (en) Hyperspectral image classification method based on joint multi-level spatial spectral information CNN
CN109241843B (en) Space-spectrum combined multi-constraint optimization non-negative matrix unmixing method
CN104978573B (en) A kind of non-negative matrix factorization method applied to Hyperspectral imagery processing
CN103413292B (en) Based on the hyperspectral image nonlinear abundance estimation method of constraint least square
CN109376753B (en) Probability calculation method for three-dimensional spatial spectrum space dimension pixel generic
CN104182978B (en) A kind of high spectrum image object detection method based on empty spectrum nuclear sparse expression
CN113094645B (en) Hyperspectral data unmixing method based on morphological component constraint optimization
CN107292258B (en) High-spectral image low-rank representation clustering method based on bilateral weighted modulation and filtering
CN104463223B (en) Hyperspectral image group sparse unmixing method based on space spectrum information abundance constraint
CN108427934A (en) A kind of Hyperspectral imaging mixed pixel decomposition method
CN108734672B (en) Hyperspectral data unmixing method based on spectral library cutting and collaborative sparse regression
CN105184302B (en) A kind of high optical spectrum image end member extraction method
CN104268561B (en) High spectrum image solution mixing method based on structure priori low-rank representation
CN105825227B (en) A sparse unmixing method for hyperspectral images based on MFOCUSS and low-rank representation
CN116403046A (en) Hyperspectral image classification device and method
CN102540271A (en) Semi-supervised hyperspectral sub-pixel target detection method based on enhanced constraint sparse regression method
CN104794457A (en) Hyperspectral image target detection method based on sparse error matrix
CN108256557B (en) Hyperspectral image classification method combining deep learning and neighborhood integration
CN105957112A (en) Hyper-spectral sub pixel detection method based on fast UNCLS
CN108414468A (en) Infrared spectrum wave band feature Enhancement Method based on wavelet transformation and nonlinear transformation
CN103065310B (en) Based on the high spectrum image edge information extracting method of three-dimensional light spectral corner statistics
CN109522918B (en) Hyperspectral image feature extraction method based on improved local singular spectrum analysis
CN105427351B (en) Compression of hyperspectral images cognitive method based on manifold structure sparse prior

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant