非凸稀疏低秩约束的高光谱解混方法

孔繁锵，卞陈鼎，李云松，郭文骏

(1. 南京航空航天大学 航天学院，江苏 南京 210016; 2. 西安电子科技大学 综合业务网理论及关键技术国家重点实验室，陕西 西安 710071)



非凸稀疏低秩约束的高光谱解混方法

孔繁锵1，卞陈鼎1，李云松2，郭文骏1

(1. 南京航空航天大学 航天学院，江苏 南京 210016; 2. 西安电子科技大学 综合业务网理论及关键技术国家重点实验室，陕西 西安 710071)

针对高光谱混合像元的丰度矩阵具有行稀疏特性，提出一种非凸稀疏低秩约束的高光谱解混方法．首先，建立高光谱图像非凸稀疏低秩约束模型，将丰度系数矩阵的非凸p范数作为稀疏约束，并将丰度系数矩阵奇异值的非凸p范数作为低秩约束；其次，构建联合低秩性先验与稀疏性先验的非凸极小化模型，并提出求解的增广拉格朗日交替极小化算法，将复合正则化问题分解成多个单一正则化问题，交替迭代求解．实验仿真结果表明，该算法比贪婪算法和凸优化算法能获得更高的解混精度，并且适用于信噪比较高的高光谱数据．

图像处理;稀疏解混;稀疏表示;低秩;凸优化

高光谱遥感图像是由几十乃至数百个连续波段图像组成的三维图像数据，有丰富的地物空间特性和光谱特性，可定性、定量地分析和识别被测对象，在环境监测、土地调查、军事探测以及海洋遥感等领域有着广泛的应用．但受光学器件性能的限制和地物复杂多样性的影响，使得高光谱遥感图像单个像元中往往包含不止一种地物覆盖类型，从而形成混合像元．混合像元在很大程度上影响了地物的识别和分类精度，因此，如何有效地解决混合像元分解问题，已成为高光谱遥感图像定量分析的关键．

在基于线性混合模型的高光谱混合像元分解中，一般采用基于几何的方法和基于统计的方法等[1]．近几年来，随着压缩传感和稀疏表示理论的快速发展，学者们将稀疏性约束加入到高光谱解混模型中，形成了第3种方法——稀疏解混方法[2]．稀疏解混方法主要有凸优化算法和贪婪算法等，在凸优化算法方面有SUnSAL(Sparse Unmixing by variable Splitting and Augmented Lagrangian)[3]、CLSUnSAL[4]和SUnSAL-TV[5]等，这些方法用L1范数来表示信号的稀疏性，能够进行高效的求解，然而数据在低噪声情况下，这些方法对噪声是不鲁棒的．在贪婪算法方面，其代表性算法有匹配追踪算法(Matching Pursuit，MP)和正交基追踪算法(Orthogonal Matching Pursuit，OMP)等[2]，然而匹配追踪算法及正交基追踪算法容易陷入局部最优问题．低秩表示方法[6-7]能在全局性的低秩约束下求解系数矩阵，来减小局部最优问题的影响．文献[8]中提出一种低秩约束的高光谱数据解混方法，对丰度系数矩阵进行低秩约束，其低秩约束解混模型为凸优化模型，能够进行高效的求解．但近几年来，非凸优化理论[7，9-10]证明在某些条件下相较于凸优化方法，非凸优化方法可进一步表征系数矩阵的稀疏性，在稀疏信号重建中对噪声具有更强的鲁棒性．

基于以上分析，为了提高高光谱解混的精度，笔者提出了一种新的高光谱图像解混方法，建立联合低秩性先验与稀疏性先验的非凸极小化模型，并使用增广拉格朗日乘子法对解混模型进行快速求解，提高求解精度．

1 联合稀疏解混模型

线性高光谱混合模型假设混合像元为光谱特性稳定的纯像元(端元)的线性组合．假设观测值包含L个波段，那么线性光谱混合模型可以表示为

由于高光谱图像中包含的端元数一般远小于光谱库中的端元数，而且高光谱数据中相邻像素所含的物质相似，因此丰度系数矩阵中非零元素所在位置具有结构化特征，即行稀疏特性．由此文献[4]中提出了联合稀疏解混模型:

其中，Y∈RL×K，为高光谱图像数据，包含L个波段，K个像元；X∈ RM×K，为丰度系数矩阵；N为噪声矩阵．那么，联合稀疏解混问题可以表示为如下优化问题:

其中，Xrow-0代表丰度矩阵X中非零元素的行数，·F表示Frobenius范数．

2 非凸稀疏低秩约束的高光谱解混方法

2.1 非凸稀疏低秩约束的高光谱解混模型

近几年来非凸正则化理论[7，9-10]说明: 在一定条件下，非凸优化在理论和实际上优于凸优化，能进一步表征矩阵或向量的稀疏性，在信号重构过程中具有更强的抗噪性能．因此，针对式(3)的线性高光谱混合像元模型，提出一种非凸稀疏低秩约束的高光谱解混模型:

其中，非负约束项lR+(·)定义为

∀i，j ．

通过重写最优化问题式(5)，利用拉格朗日乘子法将约束项作为二次罚函数加入目标函数:

其中，μ为惩罚因子，V1=X，V2=X．式(6)将无约束问题转化为多项无约束最优化问题，当惩罚因子选择适当时，能够较好地逼近最优解．首先求解参数X，然后利用X依次对每个单一正则项进行求解更新．假设当前为第 k+1 次迭代，则对X的更新是求解最优化问题:

然后，利用X(k+1)，将变量V1和V2的最优化问题分别表示为

由此，将求解式(6)的复合正则化问题分解为式(7)～(9)各个独立的、求解相对简单的单一正则化问题．

2.2 算法实现与步骤

f(σ)=f(σ(k))+f(σ(k))，σ-σ(k)=p

针对式(9)，可利用硬阈值法进行迭代更新:

对于D1和D2的迭代更新，由下式可得

3 实验结果与分析

实验1 模拟数据1是由50×50个混合像元构成的，分别包括5，7，9，11，13，15端元的6组高光谱模拟数据，丰度系数服从Dirichlet随机分布，每个混合像元满足非负性与和为一性．在生成的模拟图像中加入零均值高斯白噪声，信噪比分别为 20 dB，25 dB，30 dB，35 dB，40 dB，45 dB 和 50 dB．模拟数据实验从不同端元数解混性能和抗噪声性能两个方面分别对以上算法的解混结果进行对比分析．

不同端元数时解混性能的结果如图1所示，对 30 dB 高斯噪声下不同端元数时模拟数据解混性能进行比较．在不同端元数的情况下，笔者提出的算法和SOMP算法性能最好，SUnSAL-TV算法的性能要高于CLSUnSAL算法的，但都低于笔者提出的算法的．笔者提出的算法性能要高于SOMP算法的，在 30 dB 高斯噪声下不同端元数时平均提高了 4.042 dB．抗噪声性能的结果如图2所示，是在信噪比为 20～ 50 dB 的高斯噪声下，9端元模拟数据1解混结果比较．与其他算法相比可以看出，笔者提出的算法的信号重构误差值在大部分情况下都是最高的．与SOMP算法相比，笔者提出的算法的信号重构误差值均高于SOMP算法的，在 20～ 50 dB 高斯噪声下平均提高了 5.35 dB．因此，笔者提出的算法将联合稀疏特性和低秩特性用于稀疏解混中，相比CLSUnSAL和SOMP算法能有效地提高解混精度．

图1 30dB高斯噪声下模拟数据1解混结果比较图2 9端元时的模拟数据1解混结果比较

图4 不同信噪比高斯噪声的模拟数据2解混结果比较

实验2 模拟数据2采用文献[6]中的模拟数据实验2中的数据，由 100× 100个混合像元构成，并由USGS光谱库中随机选取的9个光谱曲线线性组合而成，其中3个原始丰度图像如图3(a)中原始丰度图像所示．在生成模拟图像后，加入零均值高斯白噪声，信噪比分别为 20 dB，25 dB，30 dB，35 dB，40 dB，45 dB 和 50 dB．

图4是在信噪比为20～50 dB的高斯噪声下，模拟数据2解混结果比较．可看出，随着信噪比的提高，各种算法的信号重构误差逐步提高．与其他算法相比，CLSUnSAL的信号重构误差最低，笔者提出的算法和SOMP算法的信号重构误差值最高．与SOMP算法相比，笔者提出的算法的信号重构误差值均高于SOMP算法的，在 20～ 50 dB 噪声下平均提高了 1.45 dB．图3为 30 dB 高斯噪声下各种算法的丰度重建图像比较，从图中可以看出，笔者提出的算法与SUnSAL-TV算法相比，其丰度重建图像包含的噪声点更多，但保留了更多的边缘区域，更接近真实丰度图像．与SOMP和CLSUnSAL算法相比，笔者提出的算法的丰度重建图像噪声点较少，并更多地保留了边缘区域，具有更好的视觉效果．

4 总 结

针对高光谱数据矩阵的秩与所观测场景的地物类别之间存在近似相等关系，即高光谱丰度系数矩阵存在低秩属性，笔者提出一种非凸稀疏低秩约束的高光谱解混方法．该方法将非凸稀疏约束和非凸低秩约束融入稀疏解混模型中，并采用增广拉格朗日交替极小化方法对解混模型进行快速求解．使用两组高光谱模拟数据进行仿真实验，用信号重构误差评价指标做比较，笔者提出的算法在性能和鲁棒性方面均优于贪婪算法和凸优化算法的．

[1] BIOUCAS-DIAS J M, PLAZA A, DOBIGEON N, et al. Hyperspectral Unmixing Overview: Geometrical, Statistical, and Sparse Regression-based Approaches[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2012, 5(2): 354-379.

[2]IORDACHE M D, BIOUCAS-DIASJ M, PLAZA A. Sparse Unmixing of Hyperspectral Data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2011, 49(6): 2014-2039.

[3]AFONSO M V, BIOUCAS-DIAS J M, FIGUEIREDO M A. An Augmented Lagrangian Approach to the Constrained Optimization Formulation of Imaging Inverse Problems[J]. IEEE Transactions on Image Processing, 2011, 20(3): 681-695.

[4]IORDACHE M D, BIOUCAS-DIAS J M, PLAZA A. Collaborative Sparse Regression for Hyperspectral Unmixing[J]. IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(1): 341-354.

[5]IORDACHE M D, BIOUCAS-DIAS J M, PLAZA A. Total Variation Spatial Regularization for Sparse Hyperspectral Unmixing[J]. IEEE Transactions on Geoscience and Remote Sensing, 2012, 50(11): 4484-4502.

[6]LIU G C, LIN Z C, YAN S C, et al. Robust Recovery of Subspace Structures by Low-rank Representation[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2013, 35(1): 171-184.

[7]MAJUMDAR A, WARD R K, ABOULNASR T. Non-convex Algorithm for Sparse and Low-rank Recovery: Application to Dynamic MRI Reconstruction[J]. Magnetic Resonance Imaging, 2013, 31(3): 448-455.

[8]QU Q, NASRABADI N M, TRAN T D. Abundance Estimation for Bilinear Mixture Models via Joint Sparse and Low-rank Representation[J]. IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(7): 4404-4423.

[9]张文娟, 冯象初. 非凸低秩稀疏约束的图像超像素分割方法[J]. 西安电子科技大学学报, 2013, 40(5): 86-91.

ZHANG Wenjuan, FENG Xiangchu. Image Super-pixels Segmentation Method Based on the Non-convex Low-rank and Sparse Constraint[J]. Journal of Xidian University, 2013, 40(5): 86-91.

[10]CHARTRAND R, STANEVA V. Restricted Isometry Properties and Nonconvex Compressive Sensing[J]. Inverse Problems, 2008, 24(3): 1-14.

[11]DONG W S, SHI G M, LI X, et al. Compressive Sensing via Nonlocal Low-rank Regularization[J]. IEEE IEEE Transactions on Image Processing, 2014, 23(8): 3618-3632.

[12]SHI Z W, TANG W, DUREN Z, et al. Subspace Matching Pursuit for Sparse Unmixing of Hyperspectral Data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(6): 3256-3274.

[13]CLARK R N, SWAYZE G A, WISE R A, et al. USGS Digital Spectral Library Splib06a [EB/OL]. [2014-12-28]. http://speclab.cr.usgs.gov/spectral.lib06.

(编辑：郭 华)

Hyperspectral unmixing method based on the non-convex sparse and low-rank constraints

KONGFanqiang1，BIANChending1，LIYunsong2，GUOWenjun1

(1. College of Astronautics, Nanjing Univ. of Aeronautics and Astronautics, Nanjing 210016, China; 2. State Key Lab. of Integrated Service Networks, Xidian Univ., Xi’an 710071, China)

Aiming at the row sparse feature of the abundance matrix of hyperspectral mixed pixels, a hyperspectral unmixing method based on the non-convex sparse and low-rank constraints is presented. The non-convex sparse representation and non-convex low-rank models are first constructed, which take the non-convex p-norm of the abundance matrix as the sparse constraint and the non-convex p-norm of the singular values of the abundance matrix as the low-rank constraint. Then the low-rank prior and the sparse prior are jointly utilized to construct a non-convex minimization model. An augmented lagrange alternating minimization method is proposed to solve the unmixing model, the compound regularization problem is decomposed into multiple single regularization problems solved by the variable separation method. Experimental results demonstrate that the proposed method outperforms the greedy algorithm and the convex algorithms with a better spectral unmixing accuracy, and is suitable for high signal-to-noise ratio hyperspectral data.

image processing; sparse unmixing; sparse representation; low-rank; convex optimization

2015-11-09

时间：2016-04-01

国家自然科学基金资助项目(61401200，61201365)

孔繁锵(1980-)，男，讲师， E-mail: kongfq@nuaa.edu.cn．

http://www.cnki.net/kcms/detail/61.1076.tn.20160401.1622.040.html

10.3969/j.issn.1001-2400.2016.06.020

TP751.1

A

1001-2400(2016)06-0116-06