赵 辉，刘衍舟，黄 橙，王天龙

(1.重庆邮电大学 通信与信息工程学院，重庆 400065；2.重庆邮电大学 图像与通信信号处理实验室，重庆 400065)

0 引 言

压缩感知(compressed sensing，CS)[1]理论区别于传统Nyquist采样定理，可以以低采样速率对信号同时进行采样与压缩，是一种更为简洁、有效的采样技术。CS理论由稀疏表示、观测矩阵的设计以及重构算法3部分组成，本文着重研究重构算法。最近研究表明，基于图像非局部自相似性(nonlocal self-similarity，NSS)的重构模型能有效地保留细节，并在图像重构中显示出巨大的优势。NSS描述了在图像中全局定位的高层次模式(例如边缘和纹理)的重复性，非局部低秩(nonlocal low-rank，NLR)[2]模型得益于此特性，并广泛应用于图像处理[2-4]。

NLR模型的核心内容主要包括两方面，一是低秩矩阵构造。传统低秩矩阵构造利用基于欧氏距离的相似块匹配方法，在固定窗口内搜索固定数目的相似块，该相似块匹配方法适合变化比较平滑的图像块，并且忽略了图像结构的多样性以及差异性，导致重构效果不佳。另一个核心内容为低秩矩阵恢复，低秩矩阵恢复可借助矩阵秩最优化问题求解，而矩阵秩最优化问题中的秩函数是非凸非连续的，为NP-hand问题，直接对NP-hand问题求解难度较大，为此寻求一种原始秩函数替代方案成为了研究的重点。

针对以上问题，本文提出了一种自适应非局部低秩的图像压缩感知重构算法，算法通过一种基于欧氏距离和结构相似性[5](structural similarity，SSIM)相似块匹配方法，自适应地选择相似块数目和搜索窗口尺寸，并用加权Schatten p-范数[6]作为秩函数求解矩阵秩最小化问题，从而提高重构算法重构性能。

1 CS理论及NLR-CS模型

CS理论指出，在某些情况下，通过求解优化问题，可以精确的从一组被噪声破坏的测量数据中恢复出原始信号。设维度为N的原始信号u∈RN×1， 通过以下变换获得原始信号u的一些M(M≪N) 线性非自适应测量y=Φu+v， 其中，y是观测向量即观测值，Φ∈RM×N是已知的观测矩阵，v∈RM是观测过程中可能引入的噪声(如加性高斯白噪声)。观测矩阵Φ在观测过程中可以当作一种退化操作，所以信号重构从本质上来说，属于一种病态的逆问题。然而，如果原始信号在稀疏变换域上是足够稀疏(或可压缩)的，且观测矩阵与稀疏矩阵高度不相关，即满足限制等距性质(RIP)的条件，则可以从欠采样量中精准地重构出原始信号。为获得对噪声具有更为稳定性的估计，通常需要依赖于原始信号的先验信息，然后，可以将原始信号的估计表述为一个基于正则化的优化问题，其中，正则项中包含了图像的某些先验信息。求解此问题的一般优化模型可以表述为

(1)

Gi={ic|d(ui,uic)

(2)

式中：Gi表示样本块ui的所有相似块在图像中的位置， d(·) 表示欧氏距离，T为相似度量阈值，uic表示样本块ui的第c个相似块[7]。找到ui的所有相似块后，再将相似块按相似度从高到低的顺序依次按列展开聚合成数据矩阵Pi=[ui,1,ui,2,…,ui,c]，Pi∈Rn×m。 由于这些图像块间是相似的，因此ui,1,ui,2,…,ui,c各元素值基本一致，所以数据矩阵Pi具有低秩特性。实际应用中，Pi通常会受到噪声的影响，因此可以将数据矩阵Pi建模为Pi=Li+Si， 其中Li为低秩矩阵，Si为高斯噪声矩阵，因此，NLR-CS[2]算法模型可以表示为

(3)

(4)

式中：Pi=Riu=[Ri1u,Ri2u,…Ricu] 表示相似图像块组成的矩阵， rank(·) 为矩阵的秩，后一项为低秩矩阵正则项。

2 基于欧氏距离和SSIM的相似块匹配方法

经典的相似块匹配方法，如式(2)所示，即欧氏距离匹配法，其相似性与欧氏距离成反比，度量依据为两点间的灰度差，由于图像结构复杂多变，致使该方法误判率较高。而图像相连像素点间灰度值其相关性较高，导致人眼不可避免的存在一些误差，而对于图像结构复杂的边缘等区域，则显得十分敏感。SSIM作为评价图像质量的指标，相比于PSNR、MSE、SAD等经典图像质量评价指标而言，SSIM更加符合人眼视觉系统[7]，SSIM的数学定义可参见文献[5]。

SSIM从不同的方面比较图像的失真，即其相似性是通过多种影响因素进行综合评估，其中，图像的结构特征占主导，并且SSIM能够有效的从图像结构特征对两幅图像进行相似性判断，因此，本文提出一种图像相似性度量方案，即SSIM联合欧氏距离。其主要思想是：首先，基于欧氏距离法在指定的区域内选出样本块ui最相似的图像块2×c个(其中，c为需要的样本块个数)构成图像块集Ci； 接着，获取当前样本块ui与图像块集Ci中各图像块的SSIM值；然后，从图像块集Ci中选c个图像块构成图像块集Cui， 以SSIM值大小为选取标准；最后，将Cui中各图像块按列展开，以图像块的SSIM值大小为顺序进行排列，然后构造成数据矩阵Pi。

3 自适应相似块数目

将相似块组聚合就可以构造成低秩矩阵，而相似块数目或多或少均会影响低秩矩阵特性，从而影响重构效果。为此，选择合适的相似块数目对于重构效果至关重要，而本文基于图像的不同结构区域，通过自适应策略选取相似块数目。该策略主要利用图像块的块结构稀疏度S(ρ)∈(0,1][8]， 根据块结构稀疏度值的大小来判断图像块所在区域的结构特征，当图像块位于图像边缘等复杂结构区域时，其值偏大；当处于图像平滑区域时，其值较小。经过实验仿真分析，相似块数目设置如下

(5)

式中：λ为经验值，一般可设置为0.55，kmax和kmin为所有点的块结构稀疏度值中最大值与最小值[7]。

4 自适应搜索窗口大小

相似块搜索方法主要分为两种，全局搜索[8]和局部搜索[2]，前者复杂度较高且耗时；而后者则利用图像块与其邻域相似性，在指定图像区域内进行搜索，降低较低复杂度的同时提高了搜索效率。但由于图像结构特征导致不同位置的图像块与其邻域相似性不一致，结构复杂区域其邻域一致性较低，如边缘等；结构简单区域，如平滑区域，其邻域一致性较高。为此，本文提出一种搜索窗口策略，即利用图像自身结构特征的自适应匹配搜索窗口大小。该策略利用图像块结构稀疏度对图像结构区域进行判定，结构复杂区域块结构稀疏度值较大，结构简单区域其稀疏度值较小，即图像结构复杂度与块结构稀疏度之间存在正比关系。自适应搜索窗口匹配策略表示如下

(6)

式中：Win=win×win为窗口大小。

5 低秩矩阵恢复

核范数已证明是原始秩函数的最紧凸性松弛,然而，由于核范数极小化对所有奇异值给予相同程度的处理，所以它并不是最适合替代原始秩函数的方法。为弥补核范数对原始秩函数逼近的不足，相继提出了一系列原始秩函数的替代方案。文献[9]提出一种全新的矩阵范数用来替换传统的核范数，即截断核范数(truncated nuclear norm，TRNN)。Zha等[10]利用奇异值的特性，对每个奇异值使用不同权重处理，提出了著名的加权核范数(weighted nuclear norm minimization，WNNM)。在TRNN基础上，定义一种最小min(m,n)-r奇异值的p(0≤p≤1) 次幂之和的方案，即截断Schatten p-范数[11]。接着，Zhang等[12]对每个奇异值分配不同的幂，然后给与不同权重的处理，即Schatten p-范数的方法。Xie等[6]提出了加权Schatten p-范数,该范数结合了Schatten p-范数和加权核范数的优点，取得了最好的效果。为此，本文将原始秩函数替换为加权Schatten p-范数。如下所示

(7)

式中： 0

6 重构算法模型

本文重构算法模型可以表示为

(8)

对于式(8)，可通过引入辅助变量转换为其等价带约束的形式，如下所示

(9)

接着，写出式(9)的增广拉格朗日函数形式，表示为

(10)

式中：a为拉格朗日常数，x为辅助变量，α为正常数。原始问题(8)的最优解即为式(10)的鞍点，因此，可利用下式对式(8)涉及的问题进行迭代求解

(11)

a(t+1)=a(t)-(u(t+1)-x(t+1))

(12)

其中，t为迭代次数。最后使用ADMM方将原问题转换成求解Li，x和u这3个子问题，具体求解过程如下

6.1 x-子问题求解

首先，固定Li和u， 由于式(13)本质上是一个严格的凸二次函数的极小化问题，即

(13)

为得到式(13)的封闭解，对式(13)进行一阶求导，令其导数为零，可得

(14)

6.2 Li-子问题求解

固定u和x，Li的优化求解可以表示为

(15)

式中：Xi,s=Ri,sx。 式(15)可通过计算加权Schatten p-范数重写如下

(16)

就奇异值特性而言，较大的奇异值往往占据着更大的权重，主要信息基本上集中在大奇异值上，即大奇异值更为重要。文献[6]给出了加权Schatten p-范数的权重更新策略，如下所示

(17)

(18)

(19)

忽略等式(18)中的约束。式(16)可重新表示为

(20)

其中

(21)

且式中的h(Li) 函数为利普希茨(Lipschiz)连续的，对于这类函数可以利用文献[14]中的方法进行处理，因此，可以将式(20)的解按下式进行处理

(22)

Li求解过程总结如算法1所示。

算法1：L(t)的求解

j=j+1；

Endfor

6.3 u-子问题求解

首先，固定Li和x，u的优化问题可以表示为

(23)

该式也为一个严格的凸二次优化问题，令其导数为零，可得其封闭解，具体如下

我国资源禀赋是富煤、贫油、少气，目前已探明的常规天然气储量不够丰富，人均水平更低；页岩气等非常规天然气储量虽然丰富，但目前开发程度还有限，产量还较低，短期内在天然气产量中占比不会太高。近年来，我国天然气进口比重逐步提升，2016年对外依存度为 34%，2017年提高到接近38%，能源供应安全存在一定风险。此外，我国天然气需求的季节性不平衡较高，特别是在冬季民生用气和发电用气之间存在矛盾，这些都是制约我国天然气发电产业发展的基础因素。

(24)

式中：I为单位矩阵。由于式(24)中包含了矩阵求逆的运算，计算复杂度较高，而图像包含了大量数据，直接求解的话，对内存和时间有较大的需求，因此，本文采用共轭梯度法[3]求解。

综上所述，本文所提图像重构算法具体过程如算法2所示。

算法2： 自适应NLR的图像压缩感知重构算法

输入： 输入原始图像u采样得到的CS观测值y， 参数T，T0，K，p；

初始化：u0=ΦTy，λ，η，a；

Outerloop:fort=1∶T

(1) 将图像划分为重叠尺寸相同的样本图像块ui；

(2) 获取样本图像块ui的块结构稀疏度值；

(3) 由式(5)选出合适的相似块数目；

(4) 由式(6)选出合适的搜索窗口大小；

(5) 对每个样本图像块ui使用本文所提的匹配方法搜索相似块， 然后利用图像块集构成矩阵低秩Pi

(6)Ift≤T0,Wi=[1,1,…,1]Τ;Endif

(7)InnerloopFork=1,2,…,K

1)通过式(14)更新x(t+1)

4) 通过式(24)更新u(t+1)；

Endfor

(8) 根据式(12)更新a；

Endfor

7 实验结果与分析

所有仿真实验均在MATLAB 2016a上进行的，计算机运行内存为6 GB，Intel(R) Core(TM) i5-4590，中央处理器(CPU) 3.00 GHz，并装载windows10 64位操作系统。为了客观评价重构算法的性能，选用TVNLR[16]、NLR-CS[2]和GSR-NCR[17]作为对比算法，峰值信噪比(PSNR/dB)、SSIM和完成一次图像重构所需要的运行时间(Time/min)作为图像重构质量客观评价指标，选用Barbara、Boats、Cameraman、Lena、Monarch、Parrots这6幅结构丰富且彼此间差异较大的图像作为仿真图像。

7.1 参数设置

本文进行图像CS重构时，观测矩阵为部分傅里叶矩阵，相关参数设置为：划分的样本图像块尺寸均为6×6；每隔5个像素点选取一个样本块，即样本块间的重叠间距为5个像素点；相似块的数目按照第3小节所提自适应策略选取，结构复杂区域的相似块数目为20个，平滑区域的样本相似块数目为45个；搜索窗口以当前样本块为中心自适应的选择大小，结构复杂区域搜索窗口大小可以设置为30×30尺寸，即图像的边缘和纹理等区域，平滑区域搜索窗口大小设置为15×15；总的迭代次数为270次。

7.2 相似块匹配方法的影响

只改变文献[2]中算法的相似块匹配方法，不区分图像结构区域，即将相似块数目和搜索窗口大小均固定，分别为45个和20×20的尺寸。在采样率为0.20时进行实验仿真，将本文所提相似块匹配方法与传统匹配方法进行对比，仿真结果见表1。

表1 不同匹配方法重构图像PSNR/dB和SSIM比较

从表1中数据可知，本文所提相似块匹配方法均能获得最高的PSNR和SSIM，与对比方法相比，其效果更佳。因此，有效地验证了本文所提相似块匹配方法有效性。

7.3 相似块数目的影响

文献[2]中，选取的样本相似块的数目为45个，并没有考虑图像的结构区域，本文基于图像结构区域自适应选取相似块数目，为使复杂结构区域样本图像块聚合而成的低秩矩阵能够获得更强的低秩特性，因此适当减少复杂结构区域相似块的数目，这种策略有效地保证了图像的结构区域。基于以上考虑，平滑区域相似块数目仍固定为45个，以结构复杂的Barbara图片为仿真图片，采样率为0.20时，仿真结果见表2。

表2 结构区域不同相似块数目重构图像PSNR/dB比较

从表2中数据可知，复杂结构区域样本块数目为20个时，能够获得最佳重构效果。因此，当样本块位于图像复杂的结构区域时，其相似块数目取20个，当样本块位于图像平滑区域时，样本块个数仍固定为45个。

7.4 搜索窗口大小的影响

当样本块位于图像的边缘、纹理等结构复杂的区域时，其邻域一致性较低，以当前样本块为中心，在指定的局部搜索区域内匹配相似块，为了尽可能搜索到足够多的图像块，应该扩大搜索区域；当样本块邻域一致性较高时，即位于图像的平滑区域，为了减小计算的复杂度，可以适当减小搜索区域的尺寸。传统搜索窗口为固定大小，且在局部区域搜索，如文献[2]，其搜索窗口大小为20×20。为了找出合适的平滑区域搜索窗口大小，可先固定复杂结构样本块的搜索窗口大小，因此，选取20×20大小的局部窗口进行仿真。其仿真结果见表3。

表3 平滑区域不同大小搜索窗口重构图像的PSNR/dB和Time/min比较

从表3中数据可知，随着搜索窗口尺寸的改变，PSNR值的变换波动较小，基本上处于一个平稳的状态，即搜索窗口大小的改变对图像重构性能影响可忽略不计。但从重构时间而已，随着搜索窗口的增大，重构需要的时间越来越多，当搜索窗口不大于15×15时，改变搜索窗口耗时变化相对较小。因此，基于重构性能与重构时间考虑，平滑区域搜索窗口的尺寸确定为15×15。接着，固定平滑区域样本块搜索窗口大小后，再对复杂结构区域样本块搜索窗口大小进行调整，其仿真结果见表4。

表4 复杂结构区域不同大小搜索窗口重构图像的PSNR/dB比较

从表4中可知，窗口大小为30×30时，PSNR可以取得最高值，即位于复杂结构区域的样本图像块，适当增大其搜索窗口的大小有利于重构效果。为此，本文在不同的结构区域，相似块搜索区域大小也不一样，当样本块位于结构复杂的区域时，搜索窗口可固定为30×30；当样本块位于平滑区域时，搜索窗口尺寸可规定为15×15。

7.5 算法重构性能评估

对无噪声的自然图像进行CS重构，图1为Barbara测试图像在采样率为0.05时，各重构算法重构图像的主观视觉对比图。该采样率下获得的图像信息较少，GSR-NCR算法和TV算法重构出的图像边缘等细节模糊较严重，重构效果较差，利用图像NSS的NLR-CS重构算法，在低采样率时，能够较好地重构出原始图像，但图像的边缘等细节仍然不能取得较好的效果，而本文所改进的重构算法，在低秩矩阵的构造中，充分利用了图像自身结构特征，以及最优的替代秩函数求解，取得了最佳重构效果，最接近真实图像，视觉效果优于另外几种算法。从图中局部放大区域可以看出，相对于其它几种算法重构出来的图像，本文改进后的算法有更好的细节保护性能。为了更客观评价各重构算法的性能，表5分别给出了在采样率在0.05、0.10和0.20条件下，自然图像无噪声重构结果的客观评价指标对比。从表5中仿真数据可以看出，与其它几种对比算法相比，本文提出重构算法在客观评价指标PSNR和SSIM上，均取得了最大值，就Barbara测试图像而言，在采样率为0.05时，较GSR-NCR算法、TVNLR算法和NLR-CS法PSNR增益分别高达8.83 dB、5.84 dB和0.54 dB。

图1 Barbara标准测试图像在采样率为0.05时各算法重构主观视觉对比

表5 各重构算法在不同采样率下重构图像的PSNR/dB和SSIM比较

通过以上仿真实验分析，有效地验证了本文所提自适应重构算法性能的优异性。

8 结束语

为了弥补基于NSS的非局部低秩模型在图像重构过程中忽略图像自身的结构特征，使得重构图像的效果不佳的问题，本文提出了一种图像重构算法，即自适应非局部低秩重构算法。首先，计算样本图像块的块结构稀疏度值，根据块结构稀疏度值来判定样本块所属的区域结构特征，根据不同的结构区域自适应选取合适的搜索窗口尺寸和相似块的数目；接着，联合SSIM和欧氏距离，提出一种更为有效的相似块匹配方案，即先利用欧氏距离在指定窗口区域内选出一定数目的预选相似块集，再计算预选相似块的SSIM值，根据SSIM值大小，选出所需要的相似块；最后，针对低秩矩阵最优化问题，本文引入加权Schatten p-范数作为秩函数的逼近。实验结果表明，提出的相似块匹配方法优于传统的相似块方法，且改进后的算法重构性能优于对比算法。本文中所提出的相似块组策略不仅可以应用在低秩矩阵构造中，同样也可以应用在组稀疏表示中相似块组的构造。由于本文算法只利用了一种图像先验信息，其适应性不强，后续的研究将联合图像其它先验信息构造联合正则项约束条件，提高重构算法的适应性，实现更好的重构效果。