1. 引言

图像去噪是图像处理邻域的一个经典病态反问题，其目的是从受噪声干扰的退化图像中恢复出原始的真实图像。图像去噪在去除噪声的同时也需要保留图像信息，因此如何既能够有效去除图像中的噪声，又能尽可能多地保留图像细节信息一直是热门的科研问题。

在现有的图像去噪问题研究中，图像的非局部自相似(Nonlocal Self-Similarity, NSS)先验 [1] 的算法得到了广泛的研究与应用。利用图像中包含大量自相似结构的性质，文献 [2] 提出了一种非局部均值(Nonlocal Means, NLM)滤波算法，该算法在有效减少高斯噪声的同时较好地保护了原图像中的局部纹理特征。NLM算法具有良好的去噪效果，吸引了众多学者，大量基于NLM的改进算法也随之被提出。比如对NSS先验的进一步研究激发了基于字典学习的稀疏性先验。常用的基于字典学习的稀疏性先验算法为KSVD字典学习去噪算法 [3] ，利用图像的稀疏性和内部自相似性先验来避免伪影。基于稀疏性先验的方法将图像当成向量来处理，忽视了图像自身块与块之间存在的结构相似特性。为了更好地利用图像块的NSS特性，通常将相似度较高的图像块组合成图像组矩阵。如文献 [4] 提出的非局部中心化稀疏表示(Nonlocally Centralized Sparse Representation, NCSR)方法的图像复原模型，文献 [5] 在非局部相似块分组技术的基础上，提出了一种三维块匹配(Block Matching 3D, BM3D)协同滤波算法理论。该算法是首个将图像的非局部自相似性和稀疏表示结合起来的去噪算法，虽然其去噪效果极好，但是在边缘纹理等方面表现并不是很好，甚至出现了伪纹理现象。

利用图像自身的非局部自相似性进行图像去噪的另一个研究思路是低秩矩阵恢复(Low-Rank Matrix Recovery, LRMR)算法，而低秩矩阵恢复算法可以分成两类：低秩矩阵分解 [6] 和秩最小化 [7] 。其中，低秩矩阵分解方法是以SVD分解为核心的，通过限制奇异值的数量或大小来约束矩阵结构的稀疏性，从而给出图像块矩阵的低秩逼近 [8] 。由于秩最小化方法是一个非凸优化问题，较难求解，一般采用核范数最小化凸近似代替秩最小化实现低秩矩阵恢复。文献 [9] 提出了核范数最小化(Nuclear Norm Minimization, NNM)算法，利用核范数近似秩函数来构建低秩去噪模型，通过求解奇异值软阈值来恢复低秩矩阵实现图像去噪。然而，该算法对奇异值同时均等地减去一个固定的阈值，忽略了不同的奇异值具有不同物理意义，造成了数据偏差。为了解决此问题，文献 [10] 提出了加权核范数最小化(Weighted Nuclear Norm Minimization, WNNM)算法，权重的分配随奇异值的不同而改变，能更好地恢复低秩矩阵用于改善去噪性能。由于容易求解和去噪效果好的特点，许多基于低秩矩阵恢复的图像去噪算法都以WNNM为基础进一步展开研究。文献 [11] 论证了非凸松弛对于核范数的优点，提出加权Schatten-p范数最小化(WSNM)模型，该模型更加灵活，也能更好地逼近原始的低秩矩阵逼近问题，获得了比NNM更好的低秩矩阵估计结果。文献 [12] 将WNNM模型中的核范数替换为伽马范数，提出基于伽马范数最小化的图像去噪算法。针对WNNM去噪算法存在的问题，文献 [13] 提出快速加权核范数最小化算法。文献 [9] 利用相似图像块张量表示低秩特性，基于张量核范数最小化估算低秩张量实现图像去噪，得到较好的去噪性能。文献 [14] 在WNNM的基础上结合彩色图像三通道特点，提出了多通道加权核范数最小化(Multi-Channel Weighted Nuclear Norm Minimization, MCWNNM)的彩色图像去噪算法。

虽然WNNM算法基于奇异值的物理特性对奇异值赋予不同的权值，但权值是由与低秩矩阵的秩相关的经验参数所确定的。选择一个合适的经验参数适应现实场景的情况是十分困难的，一旦参数选择不当，会影响图像的去噪效果。针对此问题，文章提出一种自适应加权低秩矩阵恢复的图像去噪算法，通过自适应加权策略来准确有效地恢复出原始图像的低秩矩阵部分，达到去噪目的。与其它方法不同的是，文章通过改进的Gerschgorin圆盘估计方法来估计低秩矩阵的秩，利用估计的秩设计自适应权重对奇异值不同程度的惩罚，能够更加准确有效地恢复原始图像中的低秩结构。在大量的图像去噪实验中，通过与主流算法进行主观与客观评价，文章算法均能够表现出优秀的去噪性能。

2. 加权核范数最小化模型

加权核范数最小化是一种约束矩阵奇异值稀疏性的低秩算法，在核范数最小化的理论基础上增加了核范数权重，有效利用了自然图像的先验信息并提高了核范数的灵活性。作为NNM的推广，加权核范数最小化(WNNM)模型被描述为：

$\hat{X}=\underset{X}{\arg \min }\frac{1}{{\sigma }_n^2}{\Vert Y-X\Vert }_F^2+{\Vert X\Vert }_{W,*}$ (1)

${\Vert X\Vert }_{W,*}={\displaystyle \underset{i=1}{\sum }{w}_i}{\sigma }_i\left(X\right)$ (2)

式中 ${\Vert X\Vert }_{W,*}$ 是矩阵X的加权核范数， $$ 为矩阵X中第i个奇异值对应的权重向量。文献 [15] 给出了式(1)在 ${\sigma }_n\left(X\right)\le {\sigma }_{n-1}\left(X\right)\le \cdots \le {\sigma }_0\left(X\right)$ 且 $0\le w_1\le \cdots \le w_n$ 时，上述加权核范数最小化模型的全局化最优解：

$\hat{X}=U{S}_w[\sum ]V^{\text{T}}$ (3)

式中， $S_w[\sum ]$ 是权值向量w的广义软阈值算子，也可表示为：

$S_w[\sum ]=\max \left({\sigma }_i\left(X\right)-w_i,0\right)$ (4)

3. 自适应加权核范数最小化恢复模型

3.1. 提出的模型

WNNM模型因其优良的性能在近年来被广泛的应用，但是这一类方法中权重取决于经验上选择的一个正则化参数，并与低秩矩阵的秩相关，对于不同的任务通常需要反复调节参数。针对这一问题，文章提出一种自适应加权核范数最小化图像去噪模型，能够在通过观测数据本身自适应地加权，并能够准确有效恢复出低秩矩阵。

如前所述，在低秩矩阵恢复问题中，低秩矩阵的秩是一个非常重要的参数。在某些场景中，该参数是已知的，但是绝大部分的场景中，该参数是未知的。为了解决低秩矩阵的秩未知这一问题，文章引入Gerschgorin盘估计 [15] 的思想来估计低秩矩阵的秩。根据文献 [16] 可知，在噪声环境中传感器阵列接收到的信号M同样可以表示为低秩信源信号矩阵L和稀疏噪声信号矩阵S的和。秩为r的矩阵M的协方差矩阵MR可以定义为：

$R_M=M{M}^{\text{T}}$ (5)

对 $R_M$ 进行特征值分解得：

$R_M=U_{R_M}{\sum }_{R_M}{U}_{R_M}^H$ (6)

式中 $U_{R_M}=\left[u_1,u_2,\cdots ,u_n\right]$ 是特征向量构成的矩阵， ${\sum }_{R_M}=diag\left({\lambda }_1,{\lambda }_2,\cdots ,{\lambda }_r\right)$ 是特征值构成的对角阵。在无噪环境中 $R_M$ 秩为r，但在真实环境中由于噪声影响， $R_M$ 秩为n $\left(n\gg r\right)$ 。为准确估计低秩矩阵的秩r，引入Gerschgorin盘理论。首先，划分协方差矩阵：

$R_{\ddot{M}}=\left[\begin{array}{cc}{R}_{M1}& R\\ R^H& R_{nn}\end{array}\right]$ (7)

式中 $R_{M1}\in R^{\left(n-1\right)\times \left(n-1\right)}$ 是通过删去 $R_M$ 的最后一行和最后一列得到的。对矩阵 $R_{M1}$ 进行特征值分解可得：

$R_{M1}=U_{R_{M1}}{\sum }_{R_{M1}}{U}_{R_{M1}}^{\text{T}}$ (8)

$U_{R_{M1}}=\left[{q^{\prime }}_1,{q^{\prime }}_2,\cdots ,{q^{\prime }}_{n-1}\right]$ 为矩阵 $R_{M1}$ 特征向量矩阵， ${\sum }_{R_{M1}}=diag\left({{\lambda }^{\prime }}_1,{{\lambda }^{\prime }}_2,\cdots ,{{\lambda }^{\prime }}_{n-1}\right)$ 矩阵 $R_{M1}$ 特征值矩阵。

定义一个酉变化矩阵 $U\in R^{n\times n}\left(U{U}^H=I\right)$ 如下所示：

$U=\left(\begin{array}{cc}{U}_{R_{M1}}& 0\\ 0^T& 1\end{array}\right)$ (9)

则协方差矩阵经过酉变化之后为：

$R_T=U^H{R}_{M}U=\left(\begin{array}{cc}{\sum }_{R_{M1}}& U_{R_{M1}}^{H}R\\ R^H{U}_{R_{M1}}& R_{nn}\end{array}\right)=\left(\begin{array}{cccccc}{{\lambda }^{\prime }}_1& 0& 0& \cdots & 0& {\rho }_1\\ 0& {{\lambda }^{\prime }}_2& 0& \cdots & 0& {\rho }_2\\ 0& 0& {{\lambda }^{\prime }}_3& \cdots & 0& {\rho }_3\\ ⋮& ⋮& ⋮& \ddots & ⋮& ⋮\\ 0& 0& 0& \cdots & {{\lambda }^{\prime }}_{n-1}& {\rho }_{n-1}\\ {\rho }_1^{*}& {\rho }_2^{*}& {\rho }_3^{*}& \cdots & {\rho }_{n-1}^{*}& R_{nn}\end{array}\right)$ (10)

式中， ${\rho }_i={q^{\prime }}_i^{H}R$

为了进一步提高秩估计的准确性，文章提出一种新方法来缩小Gerschgorin Disk的有效半径。因此，构造一个对角矩阵D如下：

$D=diag\left({{\lambda }^{\prime }}_1,{{\lambda }^{\prime }}_2,\cdots ,{{\lambda }^{\prime }}_{n-1},{{\lambda }^{\prime }}_n\right)$ (11)

${{\lambda }^{\prime }}_n=\sqrt{{\displaystyle \underset{1}{\overset{n-1}{\sum }}{{\lambda }^{\prime }}_i^2}}$ (12)

定义一个新的转换矩阵 $R_{TD}$ 为：

$R_{TD}=D{R}_T{D}^{-1}=\left(\begin{array}{cccccc}{{\lambda }^{\prime }}_1& 0& 0& \cdots & 0& \frac{{{\lambda }^{\prime }}_1}{{{\lambda }^{\prime }}_n}{\rho }_1\\ 0& {{\lambda }^{\prime }}_2& 0& \cdots & 0& \frac{{{\lambda }^{\prime }}_2}{{{\lambda }^{\prime }}_n}{\rho }_2\\ 0& 0& {{\lambda }^{\prime }}_3& \cdots & 0& \frac{{{\lambda }^{\prime }}_3}{{{\lambda }^{\prime }}_n}{\rho }_3\\ ⋮& ⋮& ⋮& \ddots & ⋮& ⋮\\ 0& 0& 0& \cdots & {{\lambda }^{\prime }}_{n-1}& \frac{{{\lambda }^{\prime }}_{n-1}}{{{\lambda }^{\prime }}_n}{\rho }_{n-1}\\ \frac{{{\lambda }^{\prime }}_1}{{{\lambda }^{\prime }}_n}{\rho }_1^{*}& \frac{{{\lambda }^{\prime }}_2}{{{\lambda }^{\prime }}_n}{\rho }_2^{*}& \frac{{{\lambda }^{\prime }}_3}{{{\lambda }^{\prime }}_n}{\rho }_3^{*}& \cdots & \frac{{{\lambda }^{\prime }}_{n-1}}{{{\lambda }^{\prime }}_n}{\rho }_{n-1}^{*}& R_{nn}\end{array}\right)$ (13)

其中， $R_{TD}$ 和 $R_T$ 是相似矩阵，这两个矩阵具有相同的特征值。通过公式(13)可知，经过转换后的Gerschgorin Disk其圆心位置并没有改变，但是Gerschgorin Disk的半径受到不同程度的压缩。稀疏部分的Gerschgorin Disk半径被压缩的程度比低秩部分Gerschgorin Disk被压缩的程度要大。这意味着，低秩特征值定位集与稀疏特征值定位集的间距变得更大，两个定位集的边界变得更为清晰。由此，提出了新的启发式决策规则：

$RANK\left(k\right)=\frac{1}{\sqrt{{\displaystyle \underset{1}{\overset{n-1}{\sum }}{{\lambda }^{\prime }}_i^2}}}\left[\left|{{\lambda }^{\prime }}_k\right|r_k-\frac{D^{\left(k\right)}\left(M\right)}{n-1}{\displaystyle \underset{i=1}{\overset{n-1}{\sum }}\left|{{\lambda }^{\prime }}_i\right|r_i}\right]$ (14)

式中， $k=1,2,\cdots ,n-2$ ，当公式在k处首次达到负值时，低秩矩阵的秩为 $k-1$ 。调整因子 $0<D^{\left(k\right)}\left(M\right)<1$ 如下所示：

$D^{\left(k\right)}\left(M\right)=\frac{2\left|{{\lambda }^{\prime }}_{k+1}\right|}{\sqrt{{\displaystyle \underset{1}{\overset{n-1}{\sum }}{{\lambda }^{\prime }}_i^2}}}$ (15)

根据公式(15)定义可知， $D^{\left(k\right)}\left(M\right)$ 的值无需人为调整，完全由矩阵 $R_{TD}$ 自身所定义。因此，文章改进的Gerschgorin Disk估计法使其可以自动识别低秩矩阵的秩。

通过上述所提出的改进秩估计算法得到低秩矩阵的秩，结合加权核范数最小化算法，文章提出了一种基于秩估计的自适应加权核范数最小化恢复模型，即：

$\hat{X}=\underset{X}{\arg \min }\frac{1}{{\sigma }_n^2}{\Vert Y-X\Vert }_F^2+\mu {\Vert X\Vert }_{W_X,*}$ (16)

式中 ${\Vert X\Vert }_{W_X,*}={\displaystyle \underset{i}{\sum }{\omega }_{X,i}}{\lambda }_i\left(X\right)$ ， $W_X=diag\left({\left\{{\omega }_{X,i}\right\}}_{1\le i\le \min \left(m,n\right)}\right)$ 是低秩矩阵X的权重。

对于噪声图像Y，先确定它的某个图像块，然后在尺寸为 $X\times X$ 的搜索窗口中搜索它的相似块，把它们聚集为一个矩阵，记为 $Y_j$ 。对于每组相似块矩阵 $Y_j$ 求解则有如下优化问题：

${\hat{X}}_j=\underset{X_j}{\arg \min }\frac{1}{{\sigma }_n^2}{\Vert Y_j-X_j\Vert }_F^2+\mu {\Vert X_j\Vert }_{W_X,*}$ (17)

对于最小化问题(15)，通过广义软阈值算法可以得出闭式解为：

${\sigma }_i\left({\hat{X}}_j\right)=\max \left({\sigma }_i\left(Y_j\right)-\mu {\omega }_{X,i},0\right)$ (18)

再结合给出的权重，可得到最后的奇异值闭式解为：

${\sigma }_i\left({\hat{X}}_j\right)=\{\begin{array}{l}{\sigma }_i\left(Y_j\right)-\frac{{\left({\sigma }_{r+1}\left(Y_j\right)\right)}^2}{{\sigma }_i\left(Y_j\right)},i\le r\\ 0,i>r\end{array}$ (19)

式(17)提供了估计出低秩矩阵 ${\hat{X}}_j$ 的解，最后通过组合所有组矩阵 ${\left\{{\hat{X}}_j\right\}}_i^n$ 来得到去噪后图像X。在实际去噪过程中，可以通过进行多次迭代上述去噪过程来获取更清晰的图像。

3.2. 去噪算法流程

算法：自适应加权核范数最小化图像去噪

输入：噪声图像Y

输出：去噪图像 $X^k$

初始化： ${\hat{X}}_0=Y,Y_0=X$

步骤1：设置k = 1，最大迭代次数s

步骤2：迭代正则化 $Y^k=X^{k-1}+\epsilon \left(Y-X^{k-1}\right)$

步骤3：通过式(14)估计r

步骤4：通过式(19)估计 $X_j$

步骤5：将 $X_j$ 聚合形成去噪图像 $X^k$

4. 实验结果与分析

实验的硬件设备为：处理器Intel Core(TM)i5-4210U CPU@1.70 GHz 2.40 GHz内存为4.00 GB；系统为64位Windows 10操作系统；运行环境为Matlab R2021a。

为了验证文章提出的自适应加权低秩矩阵恢复算法在图像去噪任务上的有效性，文章选用标准图像库与Berkeley数据集 [17] 中的图片进行实验测试，并使用峰值信噪比(PSNR)和结构相似性(SSIM)两个客观指标与经典算法KSVD [3] 、NCSR [4] 、BM3D [5] 和WNNM [10] 进行定量评价与分析。实验参数设置如下，搜索窗口为30 × 30。对于噪声标准差为 ${\sigma }_n<30,30\le {\sigma }_n<50,50\le {\sigma }_n<75,{\sigma }_n\ge 75$ 的含噪图像，图像块尺寸分别设置为6 × 6，7 × 7，8 × 8和9 × 9搜索相似块个数m分别为60、80、110和130，算法迭代次数分别为7、9、11、13。

4.1. 主观视觉评价

为了验证文章提出方法的有效性，选用标准数据库的测试图像如图1所示，分别对所有图像添加均值为0、标准差 ${\sigma }_n$ 为10、30、50和100的高斯噪声来生成噪声测试图像。图2为 $\sigma =10$ 时多种算法的去噪效果比较图；图3为 $\sigma =30$ 时多种算法的去噪效果比较图；图4为 $\sigma =50$ 时多种算法的去噪效果比较图；图5为 $\sigma =100$ 时多种算法的去噪效果比较图。

通过对以上各种方法实验结果的对比，可以看出，KSVD算法得到的图像质量较差，BM3D算法产生了一些人工瑕疵，并光滑掉了一些细节，NCSR和WNNM算法都在不同程度上损失了一些图像的细节信息。文章算法在有效抑制图像噪声的同时，较好地保护了图像的边缘结构和细节信息。

为更好地进行主观视觉评价，在图6和图7中，展示了Berkeley数据集House图像和Monarch图像在噪声标准差50和100的情况下不同算法的去噪结果局部细节放大图。

从图6 House图像的细节对比图片来看，KSVD、NCSR算法产生了不同程度的伪影，BM3D、WNNM算法产生了模糊的纹理结构。而本文算法得到的去噪图片整体观感较为清晰且对局部结构的处理也较为平滑。图7是Monarch图像在 $\sigma =100$ 时几种算法去噪结果的局部翅膀结构细节放大图。在高噪声的情况下，KSVD算法去噪效果最差，去噪后的细节图像模糊不清。BM3D算法在去噪后的图像中出现噪声白斑。WNNM和NCSR算法在条纹处出现大量伪影。相比之下文章算法仍然有较为清晰的翅膀细节结构，且伪影更少。

4.2. 客观评价指标

为客观地验证文章算法的有效性，采用峰值信噪比和结构相似度两个指标进行评价。峰值信噪比越大，结构相似度越接近1，说明去噪图像与原始图像越接近，去噪效果越好；反之效果越差。表1给出了噪声水平为10、30、50、100时8幅图像去噪恢复后的峰值性噪比(PSNR)和结构相似度(SSIM)的结果。从表1可以看出，KSVD算法去噪效果欠佳；NCSR算法和BM3D算法去噪性能相对一般，随着噪声水平的增加，去噪性能明显下降；WNNM算法的去噪性能相比其它算法都有较好的效果，但是WNNM算法需要通过调整经验参数来适用于不同图像。文章算法在大多数情况下评价指标取得最优值，由此可知，文章所提算法相比其他算法，在不同噪声水平情况下都有较好的去噪性能。

5. 结论

针对加权核范数最小化图像去噪对低秩矩阵秩约束不强导致图像去噪效果不佳的问题，文章提出一种基于自适应加权核范数最小化的图像去噪算法。该算法通过引入改进的Gercshgorin Disk理论思想来估计低秩矩阵的秩，再通过观察矩阵的奇异值自适应地更新加权核范数的权重，从而更有效地恢复出低秩矩阵。在含噪图像上的实验结果表明，与现有的经典去噪算法KSVD、NCSR、BM3D、WNNM等相比，文章算法在不同噪声强度下具有更好的去噪性能，能在有效去除噪声的情况下，更好地保留原图像细节信息。

基金项目

国家自然科学基金(61861040)、甘肃省科学院应用研究与开发(2018JK-02)、甘肃省重点研发计划(20YF8GA125)、兰州市科技计划项目(2018-4-35)。

NOTES

^*Email: 1424923385@qq.com

参考文献