CN104376533A - 一种基于正则化主成分追踪的图像去噪方法 - Google Patents

Abstract

本发明公开了一种基于正则化主成分追踪的图像去噪方法，通过构建正则化主成分追踪模型，并采用迭代求解方法求解出图像包含的低秩分量和稀疏分量；其中稀疏分量即对应图像中的噪声；低秩分量对应去噪后的图像。本发明解决了当前基于鲁棒主成分分析的图像去噪方法去噪效果差的问题，能有效地提高图像中噪声的去除效果。

Description

一种基于正则化主成分追踪的图像去噪方法

技术领域

本发明涉及图像处理技术领域，具体是一种基于正则化主成分追踪的图像去噪方法。 

背景技术

21世纪是数字化时代，数字光电成像设备已经广泛应用于人们日常生活中。每天人们通过数码相机以及智能手机中的照相机所获取的图像不计其数。然而人们所获得的图像经常会受到元器件噪声、热噪声、光照变化等噪声的干扰，图像的质量会不可避免地受到影响。为了获得更清晰的图像和更佳的视觉体验，就需要对受到噪声污染的图像进行去噪声处理。 

当前的图像去噪方法中，为了获得最佳的信号表示和传输效果而利用图像的统计属性来建立图像模型受到了越来越多研究者的关注。因为根据噪声的统计特征，在自然环境下拍摄的图片其中所包含的噪声多数是稀疏的，而且图像本身其结构往往又是低秩的。比如，一张受到光线污染的图片中，低秩分量对应人脸部分，而稀疏分量则往往对应图片中阴影或者脸上的反光，即所说的噪声。根据以上所提的两点特性，图像去噪问题就被转化成了低秩矩阵和稀疏矩阵恢复的数学问题。目前通常使用鲁棒主成分分析对这一问题进行求解。鲁棒主成分分析用一个核范数对秩函数进行逼近，并且用L₁范数对噪声项进行约束。尽管鲁棒主成分分析可以获得收敛解而且可以有效计算，但是鲁棒主成分分析所用的核范数和L₁范数并不能很好地逼近原有的秩函数和L₀范数，直接导致了鲁棒主成分分析获得的解偏离了问题的真实解，进一步影响了图像去噪的效果，影响了用户的视觉体验。 

发明内容

本发明的目的是提供一种基于正则化主成分追踪的图像去噪方法，以解决现有技术基于鲁棒主成分分析的图像去噪方法去噪效果差的问题。 

为了达到上述目的，本发明所采用的技术方案为： 

一种基于正则化主成分追踪的图像去噪方法，其特征在于：构建正则化主成 分追踪模型，并采用迭代求解方法求解出图像包含的低秩分量和稀疏分量，其中低秩分量对应去噪后的图像，稀疏分量对应图像中的噪声；具体步骤如下： 

(1)、输入包含噪声的图像，以黑白图像为例，提取所述图像的灰度矩阵X∈R^m×n； 

(2)、建立如公式(1)所示的初始正则化主成分追踪模型的优化函数，求解X中的低秩分量L∈R^m×n和稀疏分量S∈R^m×n： 

$\underset{L,S}{\min }{\mathrm{\λ}}_1\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_i^p\left(L\right)+{\mathrm{\λ}}_2\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{\left|S_{\mathrm{ij}}\right|}^q,s.t.X=L+S---\left(1\right),$

公式(1)中，σ_i(L)表示矩阵L的第i个奇异值，i＝1,2,L,m；0≤p,q≤1，λ₁,λ₂>0； 

(3)、将公式(1)所示的优化函数放松到如公式(2)所示的优化函数进行求解： 

$\underset{L,S}{\min }{\mathrm{\λ}}_1\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{\σ}}_i\left(L\right)+\mathrm{\ϵ})}^p+{\mathrm{\λ}}_2\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{|S_{\mathrm{ij}}+\mathrm{\ϵ}|}^q+\frac{1}{2}{\left|\right|L+S+X\left|\right|}_F^2---\left(2\right),$

公式(2)中，参数0<ε≤1； 

(4)、利用迭代方法求解式(2)中的矩阵L和S，包括以下步骤： 

(4.1)、初始化迭代次数变量k＝0；利用matlab随机产生两个m×n维的随机矩阵分别作为低秩分量L和稀疏分量S的初始矩阵L⁰和S⁰； 

(4.2、)利用公式(3)更新L^k为L^k+1： 

$L^{k+1}=f_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(H^k\right)={\mathrm{US}}_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)V^T---\left(3\right),$

公式(3)中：权值向量参数μ₁>1；矩阵矩阵U、Σ和V^T是按照H^k＝UΣV^T分解形式对矩阵H^k进行奇异值分解所得，软阈值收缩算子如公式(4)所示： 

$S_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)=\mathrm{Diag}\left\{{({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}\right\}---\left(4\right),$

公式(4)中： ${({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}=\max ({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1,0);$ 其中Diag{ }表示对角化；权值向量中各元素取值为： 

$w_i^k=\frac{p}{{({\mathrm{\&sigma;}}_i\left(L^k\right)+\mathrm{\&epsiv;})}^{1-p}},i=1,L,m---\left(5\right),$

公式(5)中：ε是一个大于0的参数； 

(4.3)、利用公式(6)更新S^k为S^k+1： 

公式(6)中：权值矩阵参数μ₂>0.5；矩阵  是一个阈值收缩算子；利用公式(7)计算矩阵中第i行第j列元素的值： 

公式(7)中：表示矩阵T^k的第i行第j列的元素，1≤i≤m，1≤j≤n；表示矩阵M^k的第i行第j列的元素；利用公式(8)计算

$M_{\mathrm{ij}}^k=\frac{q}{{\left(\right|S_{\mathrm{ij}}^k|+\mathrm{\&epsiv;})}^{1-q}}---\left(8\right),$

(4.4)、令k＝k+1，即使迭代次数变量以1为单位进行递增；如果迭代次数变量k<20，则返回步骤(4.2)，继续进行迭代计算；如果k＝20，则迭代结束，并将此时的L²⁰和S²⁰作为式(1)的解L和S； 

至此，包含噪声的图像矩阵X中的低秩分量L和稀疏分量S就已经被求解出，低秩分量L即对应去除噪声后的图像，稀疏分量S即对应图像中的噪声部分。 

与已有技术相比，本发明有益效果体现在： 

(1)本发明提出了一种基于正则化主成分追踪的图像去噪方法，引入L_p,q范数，建立正则化主成分追踪模型来近似鲁棒主成分分析问题，克服了鲁棒主成分分析方法所求出的解偏离真实解的缺点，有效地恢复出了图像的低秩分量和稀疏分量，提高了图像去噪的效果。 

(2)本发明首次将L_p,q范数用来解决鲁棒主成分分析问题，并且提出了一种迭 代加权求解算法用来求解所建立的正则化主成分追踪模型，本发明方法在多次迭代之后，所求出的解能够很好地比逼近真实解。 

附图说明

图1为本发明的图像去噪方法工作流程图； 

图2(a)为不包含噪声的原始图像； 

图2(b)为加入噪声的图像； 

图2(c)为应用对比方法1进行去噪后的图像； 

图2(d)为应用对比方法2进行去噪后的图像； 

图2(e)为应用本发明进行去噪后的图像； 

具体实施方式

本实施例中，一种基于正则化主成分追踪的图像去噪方法可用于对黑白图像以及彩色图像的噪声进行去除，其特征是构建正则化主成分追踪模型，并采用迭代求解方法求解出图像包含的低秩分量和稀疏分量；所述稀疏分量即对应图像中的噪声；低秩分量对应去噪后的图像。 

如图1所示，本发明图像去噪方法具体实施步骤如下： 

步骤1、输入包含噪声的图像，以黑白图像为例，提取所述图像的灰度矩阵X∈R^m×n； 

对于彩色图像，彩色图像包含红绿蓝三个颜色通道，对每个通道的图像矩阵都提取低秩分量和稀疏分量，然后将三个通道的低秩分量合并起来，即构成去噪后的彩色图像。 

为了测试正则化主成分追踪的去噪效果，从Google图片搜索引擎中下载了50幅清晰的图片(认为是不含噪声的图片)，每幅图片都归一化为300×300大小，即m＝n＝300。对每幅上50％的像素叠加高斯噪声，该噪声满足正态分布N(0,0.2²)。 

步骤2、将图像去噪问题转化为低秩矩阵和稀疏矩阵恢复问题，可以用式(1)表示： 

$\underset{L,S}{\min }\mathrm{rank}\left(L\right)+\mathrm{\&lambda;}{\left|\right|S\left|\right|}_0,s.t.X=L+S---\left(1\right)$

式(1)中：rank(L)是秩函数，表示低秩分量L的秩；||g||₀为0范数，表示非零 元素的个数，即||S||₀表示稀疏分量S的稀疏度。 

为了求解式(1)所示的低秩矩阵和稀疏矩阵恢复问题，建立了如式(2)所示的初始正则化主成分追踪模型，求解X中的低秩分量L∈R^m×n和稀疏分量S∈R^m×n： 

$\underset{L,S}{\min }{\mathrm{\&lambda;}}_1\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i^p\left(L\right)+{\mathrm{\&lambda;}}_2\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{\left|S_{\mathrm{ij}}\right|}^q,s.t.X=L+S---\left(2\right)$

式(2)中，σ_i(L)表示矩阵L的第i个奇异值，i＝1,2,L,m；0≤p,q≤1，λ₁,λ₂>0； 

在测试中，参数设置为p＝q＝0.5，λ₁＝1，λ₂＝0.1； 

步骤3、将式(2)所示的优化函数放松到如式(3)所示的优化函数进行求解： 

$\underset{L,S}{\min }{\mathrm{\&lambda;}}_1\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{({\mathrm{\&sigma;}}_i\left(L\right)+\mathrm{\&epsiv;})}^p+{\mathrm{\&lambda;}}_2\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{|S_{\mathrm{ij}}+\mathrm{\&epsiv;}|}^q+\frac{1}{2}{\left|\right|L+S+X\left|\right|}_F^2---\left(3\right)$

式(3)中，参数0<ε≤1；测试中设置为ε＝0.001；为了便于求解，式(3)将式(2)中的放松到(σ_i(L)+ε)^p，将|S_ij|^q放松到|S_ij+ε|^q。 

式(3)是一个非平滑的问题，包含两个变量。本发明中提出一个迭代加权算法专门用于求解式(3)。具体如步骤4.1到步骤4.4所示。 

步骤4、利用迭代加权方法求解式(3)中的矩阵L和S，在每一次的迭代计算中交替更新低秩分量L和稀疏分量S，在迭代次数达到一个较大的数值时，函数会收敛到一个很小的数值，此时迭代就可以结束了，就将此时得到的低秩分量矩阵和稀疏分量矩阵作为式(3)的解。 

步骤4.1、初始化迭代次数变量k＝0；迭代次数k在每一次迭代后递增。利用matlab随机产生两个m×n维的随机矩阵分别作为低秩分量L和稀疏分量S的初始矩阵L⁰和S⁰； 

步骤4.2、首先在第k次迭代时，先固定稀疏分量矩阵S＝S^k，则式(3)即可简化为一个变量的问题，如式(4)所示： 

$L^{k+1}=\arg \underset{L}{\min }{\mathrm{\&lambda;}}_1\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{({\mathrm{\&sigma;}}_i\left(L\right)+\mathrm{\&epsiv;})}^p+\frac{1}{2}{\left|\right|L+S^k-X\left|\right|}_F^2---\left(4\right)$

利用式(5)求解式(4)，获得L^k+1： 

$L^{k+1}=f_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(H^k\right)={\mathrm{US}}_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)V^T---\left(5\right)$

式(5)中：权值向量参数μ₁＝2.1；矩阵矩阵U、Σ和V^T是按照H^k＝UΣV^T分解形式对矩阵H^k进行奇异值分解所得。软阈值收缩算子如式(6)所示： 

$S_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)=\mathrm{Diag}\left\{{({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}\right\}---\left(6\right)$

式(6)中： ${({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}=\max ({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1,0);$ 其中Diag{ }表示对角化；权值向量中各元素取值为： 

$w_i^k=\frac{p}{{({\mathrm{\&sigma;}}_i\left(L^k\right)+\mathrm{\&epsiv;})}^{1-p}},i=1,L,m---\left(7\right)$

步骤4.3、然后固定低秩变量矩阵L＝L^k，则式(3)转化为如式(8)所示的单变量的问题： 

$S^{k+1}=\arg \underset{S}{\min }{\mathrm{\&lambda;}}_2\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{\left(\right|S_{\mathrm{ij}}|+\mathrm{\&epsiv;})}^q+\frac{1}{2}{\left|\right|L^k+S-X\left|\right|}_F^2---\left(8\right)$

利用式(9)计算式(8)，获得S^k+1： 

式(9)中：权值矩阵参数μ₂＝2.1；矩阵 是一个阈值收缩算子；利用式(10)计算矩阵中第i行第j列元素的值： 

式(10)中：表示矩阵T^k的第i行第j列的元素，1≤i≤m，1≤j≤n；表示矩阵M^k的第i行第j列的元素；利用式(8)计算

$M_{\mathrm{ij}}^k=\frac{q}{{\left(\right|S_{\mathrm{ij}}^k|+\mathrm{\&epsiv;})}^{1-q}}---\left(8\right)$

步骤4.4：k＝k+1；如果迭代次数变量k<20，则返回步骤4.2，继续进行迭代计算；如果k＝20，则迭代结束，并将此时的L²⁰和S²⁰作为式(1)的解L和S； 

至此，包含噪声的图像矩阵X中的低秩分量L和稀疏分量S就已经被求解出，低秩分量L即对应去除噪声后的图像，稀疏分量S即对应图像中的噪声部分。 

通过图2(a)、图2(b)、图2(c)、图2(d)和图2(e)几幅实验图像的对比，可以看出本发明方法的去噪效果和明显优势。图2(a)为直接从Google图片搜索引擎上下载的图片，不包含噪声；图2(b)是在图2(a)的基础上叠加了高斯噪声；图2(c)和图2(d)为两种基于鲁棒主成分分析的图像去噪方法的去噪效果；图2(e)是本发明的基于正则化主成分追踪的图像去噪方法的去噪效果。很明显，本发明的方法的去噪效果优于两种基于鲁棒主成分分析的图像去噪方法，是最接近原始图像的，能够去除大部分噪声。 

本发明主要提出了一种图像去噪的方法，在市场化应用时，可以考虑将其移植手机的应用程序中，每次手机拍完照片，自动检测检测照片质量，当发现照片存在“噪声”时，自动执行去噪程序。本发明的图像去噪方法尤其对自然界常见的“高斯噪声”和“椒盐噪声”两种图像噪声具有突出的效果。目前在全球范围内，截止到2014年，以智能手机和平板电脑为代表的便携式互联网设备保守估计已达到20亿部，这些智能设备都拍照程序，如果每一个设备中都加入本发明的图像去噪程序，市场规模将会很大。 

以上，仅为本发明较佳的一种实施方式，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明的技术方案及其发明构思加以等同替换或相关参数改变，都应涵盖在本发明的保护范围之内。 

Claims (1)

1.一种基于正则化主成分追踪的图像去噪方法，其特征在于：构建正则化主成分追踪模型，并采用迭代求解方法求解出图像包含的低秩分量和稀疏分量，其中低秩分量对应去噪后的图像，稀疏分量对应图像中的噪声；具体步骤如下：

(1)、输入包含噪声的图像，以黑白图像为例，提取所述图像的灰度矩阵X∈R^m×n；

(2)、建立如公式(1)所示的初始正则化主成分追踪模型的优化函数，求解X中的低秩分量L∈R^m×n和稀疏分量S∈R^m×n：

$\underset{L,S}{\min }{\mathrm{\&lambda;}}_1\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\mathrm{\&sigma;}}_i^p\left(L\right)+{\mathrm{\&lambda;}}_2\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{\left|S_{\mathrm{ij}}\right|}^q,s.t.X=L+S---\left(1\right),$

公式(1)中，σ_i(L)表示矩阵L的第i个奇异值，i＝1,2,L,m；0≤p,q≤1，λ₁,λ₂>0；

(3)、将公式(1)所示的优化函数放松到如公式(2)所示的优化函数进行求解：

$\underset{L,S}{\min }{\mathrm{\&lambda;}}_1\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{({\mathrm{\&sigma;}}_i\left(L\right)+\mathrm{\&epsiv;})}^p+{\mathrm{\&lambda;}}_2\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{|S_{\mathrm{ij}}+\mathrm{\&epsiv;}|}^q+\frac{1}{2}{\left|\right|L+S-X\left|\right|}_F^2---\left(2\right),$

公式(2)中，参数0<ε≤1；

(4)、利用迭代方法求解式(2)中的矩阵L和S，包括以下步骤：

(4.1)、初始化迭代次数变量k＝0；利用matlab随机产生两个m×n维的随机矩阵分别作为低秩分量L和稀疏分量S的初始矩阵L⁰和S⁰；

(4.2、)利用公式(3)更新L^k为L^k+1：

$L^{k+1}=f_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(H^k\right)={\mathrm{US}}_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)V^T---\left(3\right),$

公式(3)中，权值向量参数μ₁>1；矩阵矩阵U、Σ和V^T是按照H^k＝UΣV^T分解形式对矩阵H^k进行奇异值分解所得，软阈值收缩算子如公式(4)所示：

$S_{{\mathrm{\&lambda;}}_1{w}^k/{\mathrm{\&mu;}}_1}\left(\mathrm{\&Sigma;}\right)=\mathrm{Diag}\left\{{({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_i{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}\right\}---\left(4\right),$

公式(4)中： ${({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1)}_{+}=\max ({\mathrm{\&Sigma;}}_{\mathrm{ii}}-{\mathrm{\&lambda;}}_1{w}_i^k/{\mathrm{\&mu;}}_1,0);$ 其中Diag{}表示对角化；权值向量中各元素取值为：

$w_i^k=\frac{p}{{({\mathrm{\&sigma;}}_i\left(L^k\right)+\mathrm{\&epsiv;})}^{1-p}},i=1,L,m---\left(5\right),$

公式(5)中：ε是一个大于0的参数；

(4.3)、利用公式(6)更新S^k为S^k+1：

$S^{k+1}=S_{{\mathrm{\&lambda;}}_2{M}^k/{\mathrm{\&mu;}}_2}^{\%}\left(T^k\right)---\left(6\right),$

公式(6)中：权值矩阵参数μ₂>0.5；矩阵 是一个阈值收缩算子；利用公式(7)计算矩阵中第i行第j列元素的值：

${\left[S_{{\mathrm{\&lambda;}}_2{M}^k/{\mathrm{\&mu;}}_2}^{\%}\left(T^k\right)\right]}_{\mathrm{ij}}=\left\{\begin{array}{cc}{T}_{\mathrm{ij}}^k-{\mathrm{\&lambda;}}_2{M}_{\mathrm{ij}}^k/{\mathrm{\&mu;}}_2& \mathrm{if}{T}_{\mathrm{ij}}^k>{\mathrm{\&lambda;}}_2{M}_{\mathrm{ij}}^k/{\mathrm{\&mu;}}_2\\ T_{\mathrm{ij}}^k+{\mathrm{\&lambda;}}_2{M}_{\mathrm{ij}}^k/{\mathrm{\&mu;}}_2& \mathrm{if}{T}_{\mathrm{ij}}^k<-{\mathrm{\&lambda;}}_2{M}_{\mathrm{ij}}^k/{\mathrm{\&mu;}}_2\\ 0& \mathrm{other}\end{array}\right.---\left(7\right),$

公式(7)中：表示矩阵T^k的第i行第j列的元素，1≤i≤m，1≤j≤n；表示矩阵M^k的第i行第j列的元素；利用公式(8)计算

$M_{\mathrm{ij}}^k=\frac{q}{{\left(\right|S_{\mathrm{ij}}^k|+\mathrm{\&epsiv;})}^{1-q}}---\left(8\right),$

(4.4)、令k＝k+1，即使迭代次数变量以1为单位进行递增；如果迭代次数变量k<20，则返回步骤(4.2)，继续进行迭代计算；如果k＝20，则迭代结束，并将此时的L²⁰和S²⁰作为式(1)的解L和S；

至此，包含噪声的图像矩阵X中的低秩分量L和稀疏分量S就已经被求解出，低秩分量L即对应去除噪声后的图像，稀疏分量S即对应图像中的噪声部分。