基于低秩和稀疏模型的高光谱图像快速去噪方法

杨 垚 , 黄 聪, 王华军

(成都理工大学 a.地球物理学院,b.地球勘探与信息技术教育部重点实验室,成都 610059)

0 引言

高光谱遥感影像(HSI)是当下遥感领域的研究热点，其影像包含大量空间信息和光谱信息，众多研究领域已广泛使用HSI提取可用信息。然而，HSI在成像、传输等过程中无可避免的受到噪声的干扰，复杂的噪声类型(高斯噪声、脉冲噪声、死线等),降低了图像的分辨率并给科学研究增加了难度。应对复杂噪声的干扰，尽量提高信噪比且最大程度保留有用信息数据的去噪流程是必需的。各领域由于研究目的和方法的不同，影像的去噪速度和去噪效果等众多因素各有所需。

为了更加适用于各领域，众多HSI去噪算法相继被提出，HSI数据有多个频带，每个频带可看作一张灰度图像，2-D图像去噪的方式是最直接的(如BM3D、WNNM)，此类方法利用了频带信息，忽略了频带之间的强相关性，实际中不仅不能有效地恢复图像，而且还会有伪像和失真的现象[1-2]。考虑HSI多频带相关特性，众多基于谱间相关性的3-D去噪方法被提出。3-D去噪方法分为两类，①基于变换域，该方法一般基于主成分分析，离散变换或小波变换(如基于主成分分析、基于离散余弦变换(DCT)、基于小波变换的滤波[3-5])；②基于空间域的方法，该方法更多地利用了空间和频谱信息。稀疏先验模型使用字典学习来有效利用HSI的RAC(the redundancy and correlation，冗余与相关性)特点，因HSI干净数据可以由少量原子线性表示，而噪声数据不能(如将高度相关的频段聚集到一起充分利用全局和局部的RAC，再用频谱和空间联合稀疏表示和联合字典学习的去噪办法)[6]。全变分(TV)先验模型运用像素之间的平滑特征，方法在去噪和保护边缘的统一表现突出，但是由于HSI噪声类型复杂，一般与其他模型一起使用(如提出全变分正则化和低秩约束模型结合使用，低秩和稀疏矩阵的分解为TV正则化提供噪声分量，TV正则化后获得了增强的干净图像[7])。低秩先验模型从高维数据中发现和利用低维结构来去噪，是HSI去噪最有力的约束模型。基于干净HSI低秩先验特性和非高斯噪声稀疏性，建立了低秩矩阵恢复模型(low-rank matrix recovery，LRMR)，使用GoDec和增广拉格朗日乘子法(ALM)求解提出的恢复模型，最终得到低秩的矩阵(干净的HSI)[8-9]；LRMR模型没有考虑不同频带噪声强度不同的事实，进而提出一种噪声调整的迭代LRMA改进方法(NAILRMA)，引入随机奇异值分解来求解提出的迭代模型，通过每次迭代的反馈数据不断更新噪声强度来提升去噪能力[10]；构建了l0梯度正则化低秩张量(LETL0)模型，使用ALM求解模型，用此方法降噪后HSI图像边缘大部分保留[11]。

由于频谱信息冗余的充分使用，低秩先验模型在HSI去噪上大放光彩，然而在HSI空间上的相关性却没有充分利用。这里考虑了低秩模型的长短处，提出一种基于子空间低秩和稀疏矩阵分解的迭代算法SLRSF(subspace low-rank matrix and sparse matrix factorization)，将HSI每个频带矢量化处理并将三维张量重塑为二维矩阵，再将二维矩阵映射到低维频谱子空间，使用低秩和稀疏张量正则化的方法迭代处理特征图像得到干净特征图像，再还原得到干净图像。该方法充分利用了RAC，算法效率高，去噪能力优秀。

1 方法

高光谱影像噪声类型众多，常见且对数据干扰较大的噪声有高斯噪声、脉冲噪声、坏点或条纹、死线等。考虑这些噪声数据特点(稀疏和低秩)可将HSI分解为：

Y=X+S+N

(1)

其中:Y为观测数据；X为干净的数据；S为稀疏噪声(包含脉冲噪声、条纹、死线)；N为高斯噪声；Y、X、S、N∈Rh×p；m、n分别是HSI图像空间的长和宽；p是光谱波段数；h=m*n。

2.1 频谱低秩分解

HSI光谱波段数较多且高度相关，为此使用一个k维子空间E来线性表示X。

X=EZ

s.t.ETE=Ik

(2)

其中:k≪p；Ik∈Rk×k是单位矩阵；X∈Rh×p；E∈Rh×k，E由HySime[12]算法计算得出，Z被称为RCI(representation coefficients image)，它的行称为本征图像，Z∈Rk×p。进而HSI退化模型表示为式(3)。

Y=EZ+S+N

(3)

1.2 迭代模型

假设N是满足高斯分布且均值为0，HSI的去噪问题表示为：

(4)

(5)

式(4)中第一项控制数据保真度，第二项控制约束稀疏噪声，第三项控制干净图像的低秩特性。

式(4)的低秩正则化和稀疏正则化的拉格朗日模型为式(6)。

rank(X)≤r

(6)

由式(2)，RCI的去噪模型可以表示为式(7)。

rank(Z)≤r

(7)

1.3 优化及实现

对于式(7)，我们使用增广拉格朗日乘子法(augmented lagrange method)来解决它。

λ‖S2‖1+<Λ,Y2-Z-S2>+

(8)

Y2=ETY,S2=ETS

(9)

式中：Λ为拉格朗日乘数，式(8)的解法为对其中一个变量迭代时固定其余变量。

1.3.1Z的迭代

包含Z的子问题为：

(10)

参考文献[13]中第二节给出了求解式(10)的具体方法。

(11)

1.3.2S2的迭代

包含S2的子问题为：

(12)

运用下面的软阈值运算来解式(12)，其中T>0。

(13)

(14)

1.3.2Λ的迭代

Λ(k+1)=Λ(k)+μ(Y2-Z(k+1)-S(k+1))

(15)

基于SLRSF的低秩先验模型算法流程如下。

输入：加噪HSI数据Y∈Rm×n×p，参数lamda，迭代停止条件ε

输出：干净HSI数据X

1：Y=reshape(Y,m*n,p)’

2： 使用HySime从Y中计算子空间E和subdim

3：E=E[:,1:subdim]

4：Y2=E’ *Y;

5：Y2=reshape(Y2’,m,n,subdim)’

6：Y2_clean=RCI_LR(Y2,lamda)

6.1：Z=0，S2=0，∧=0，μ=0.01，μmax=106，ρ=1.5

6.2： 用式(11)、式(14)、和式(15)分别更新Z、S2、∧

6.3：μ=min(ρμ,μmax)

6.4： 判断，若不满足不等式则继续从5.2开始执行，若满足则迭代完成，若超过一定迭代次数则自动退出。

7：Y2_clean=reshape(Y2_clean’,m,n,subdim)’

8：X=E*Y2_clean

用于混合噪声的高光谱图像降噪。

2 算法测试与结果分析

为验证方法和代码的可行性与效果，分别对Indian pines、KSC、Pavia、Pavia80、Salinas、Salinas A、Cuprite数据集加噪。加噪类型分为：①标准差相同的高斯噪声；②标准差不同的高斯噪声；③标准差不同的高斯噪声，并随机加入脉冲噪声；④标准差不同的高斯噪声，随机加脉冲噪声和死线。

这里选取了ALM RPCA、LRMR、BM4D方法对比本文提出的方法。量化指标采用平均峰值信噪比(MPSNR)、平均结构相似度(MSSIM)和光谱差异指标(ERGAS)来测评这些方法。视觉上图像细节保存尺度和噪声消除能力大小也是一重要指标。

2.1 Indian pines

Indian pines数据集是由AVIRIS传感器在印第安纳州西北部采集的数据，波长范围为0.4*10-6m到2.5*10-6m，224个光谱带，每个光谱带145*145像素。经过去除吸水条带的数据处理后，还剩220个光谱。

对于混合噪声的去噪后图像的所有波段，BM4D的去噪效果明显不佳，高斯噪声与死线均无法去除，LRMR去死线的效果不太好，其余方法均能去除死线、高斯噪声和椒盐噪声的影响。图1(b)、图1 (d)、图1 (e)和图2(b)、图2 (e)的右下图区域恢复不是很清晰，提出的SLRSF基本能保留原图像的细节信息且去噪效果良好。

图1 Indian pines数据集第129波段死线与高斯噪声情况实验结果对比

图2 Indian pines数据集第29波段椒盐噪声与高斯噪声情况实验结果对比

2.2 KSC

Kennedy Space Center(KSC)是NASA AVIRIS仪器在弗罗里达州肯尼迪航天中心采集的224个10 nm带宽的高光谱数据，中心波长范围在400 nm到2 500 nm之间。去除吸水条带和低SNR频段后，剩下176个频段可用。

对于Indian pines数据集(表1)，只有高斯噪声的情况下SLRSF和BM4D的量化指标都表现不错，其中SLRSF综合表现更为优秀。当加入脉冲噪声或死线混合噪声时SLRSF各项量化指标都比其余三种方法好，且花费的时间远远低于它们。

表1 基于Indian pines数据集的算法比较

对于KSC数据集(表2)(数据量是Indian pines的12倍)，当数据集增大时，SLRSF算法运行所需内存约为其余方法的0.1倍，且运行时间远小于其余方法。经过对各频段的图像对比后发现，BM4D在纯高斯噪声情况下表现优秀，混合噪声的去噪效果不好；ALM RPCA在各种噪声下平均结构相似度很高；SLRSF在4种噪声情况下量化指标均在中等偏上水平，且在噪声类型越复杂的数据集下表现优于其余方法。

表2 基于KSC数据集的算法比较

当将算法运用于数据集Salinas(AVIRIS成像光谱仪对美国加利福尼亚州Salinas山谷所成的像)和Pavia University(德国机载反射光学光谱成像仪对意大利的帕维亚城所成的像)时，发现若一个原始波段被噪声干扰特别严重(基本无图像显示)而相邻波段有数据时，这个波段在去噪后会有不同程度的假图像生成，可能是过度利用了光谱的相关性。总的来说，对于不同数据集SLRSF均能快速有效地去除人工合成噪声，各光谱去噪视觉效果表现不错，各量化指标均表现突出。

3 结论

笔者基于低秩和稀疏矩阵分解方法，针对该类方法在空间相关性信息上利用的不足，提出一种基于子空间低维投影特征图像低秩和稀疏矩阵分解的迭代方法。通过噪声迭代的方式逐步消除高斯噪声和非高斯噪声的影响，得到干净图像。测试模拟数据集比较结果后得出，该方法不仅能快速有效去除高斯噪声、椒盐噪声和死线的影响，恢复出高质量的干净图像，并且在定量评价指标的表现上也非常优秀，最突出的是算法空间复杂度远低于其余方法，用于HSI图像即时处理和批量处理等的选择上有很大的优先性。