基于RPCA模型的P范数优化算法

第37卷第4期 温 州 大 学 学 报（自 然 科 学 版） 2016年11月 V ol 37, No 4 Journal of Wenzhou University (Natural Sciences) Nov, 2016基于RPCA 模型的P 范数优化算法刘 园，王 迪(温州大学数学与信息科学学院，浙江温州 325035)摘 要：在图像修复和视频处理中，低秩矩阵恢复有着非常广泛的应用．RPCA 模型是低秩矩阵恢复的经典模型，其基本思想是将一个数值矩阵分解为一个低秩矩阵与一个稀疏矩阵和的形式再进行求解．然而，RPCA 问题是NP 难的，一个通用的处理方式就是将RPCA 模型中矩阵的秩函数和0L 范数分别松弛为矩阵的核范数和1L 范数，从而将其近似转化为凸优化问题来求解，但这种由凸优化近似方法得到的解在相对较弱的非相干性条件下会使原始问题的解退化．针对这个问题，本文首先提出一种更接近于原始问题的非凸近似模型，即用矩阵的Schatten -p 范数和p L 范数（01p <<）分别代替矩阵的秩函数和0L 范数，然后针对提出的非凸近似模型，进一步给出有效的优化算法，最后，在人工数据集和真实图像数据集上进行实验，结果表明，所提出的模型是有效的． 关键词：低秩矩阵恢复；RPCA 模型；Schatten -p 范数；p L 范数中图分类号：TP391 文献标志码：A 文章编号：1674-3563(2016)04-0025-08DOI ：10.3875/j.issn.1674-3563.2016.04.005 本文的PDF 文件可以从 获得近些年来，随着互联网和信息科学的快速发展，人们可以获得大量的可视化数据，从这些可视化数据中，研究者们可以提取出与图像和视频处理等领域相关的有用信息．然而，在一些图像和视频的处理过程中，很多时候所给的图像与视频数据都伴随着光照强度变化、部分遮挡以及损坏等一系列噪音，能够从这些带有噪音污染的数据中准确恢复出原有的图像数据，是一个非常值得研究的问题．最近几年的研究表明，一些可视化数据可以表示为低秩部分与稀疏部分的和的形式，例如，对于一张带有光照阴影的人脸图片来说，低秩部分和稀疏部分分别捕获了人脸和光照阴影；在监控视频剪辑中，低秩部分和稀疏部分分别对应着背景和移动目标．最近几年，这些研究在图像和视频处理等领域中得到了高度的关注．1 背景知识假定一个观测数据矩阵n m R D ⨯∈，目标是恢复出它的低秩部分n m R A ⨯∈和稀疏部分E ∈m n R ⨯，使得E A D +=，也即下列鲁棒主成分分析（Robust Principle Component Analysis, RPCA ）优化问题：0,min rank()||||A EA E λ+，..s t D A E =+， （1）其中rank()A 表示矩阵A 的秩，0||||E 表示矩阵E 的0L 范数，即矩阵E 中非零元素的个数，λ是收稿日期：2015-09-18作者简介：刘园(1989- )，男，安徽六安人，硕士研究生，研究方向：应用分析与最优化理论温州大学学报（自然科学版）(2016)第37卷第4期26一个正的权衡参数．由于（1）中的目标函数是非连续且非凸函数，且优化问题（1）是NP-难问题，因此不易求解．一个常用的方法就是将秩函数松弛为凸的核范数，将0L 范数松弛为1L 范数，从而将问题（1）松弛为以下凸优化问题：1,min ||||||||A EA E λ*+，..s t D A E =+， （2）即主成分追踪（Principle Component Pursuit ，PCP ）问题[1]，其中*||||A 表示矩阵A 的核范数，即矩阵A 的奇异值之和，1||||E 表示矩阵E 的1L 范数，即矩阵E 所有元素的绝对值之和．文献[2]证明了在满足较强的非相干性条件下，通过求解凸优化问题（2）能获得非凸优化问题（1）的最优解．进一步，文献[3]给出了求解问题（2）的交替方向算法（Alternating Direction Method ，ADM ）．然而，当这些假设条件不满足时，凸优化问题（2）的解就会偏离原始问题（1）的最优解．相比较而言，p L 范数（01p <<）比1L 范数更接近于0L 范数，并且已有的相关文献[4]也通过理论分析和大量实验证明了p 范数在相对较弱的非干性条件下就能保证得到原始问题（1）的最优解．更重要的是，它们指出虽然通过求解p L 范数非凸优化问题得到的解是局部最优解，但此局部最优解要比1L 范数凸优化问题的全局最优解的效果更好，即运用p L 范数可成功恢复出具有更低信噪比的信号来[5]．因此，在这篇文章中我们运用非凸的替代函数，即将0L 范数松弛为p L 范数，将秩函数松弛为Schatten -p 范数，提出了一个新的矩阵恢复模型．大量实验表明本文算法要好于PCP 优化算法．本文的贡献可总结为以下几点：1）为低秩矩阵恢复问题提出了一个Schatten -p 范数与p L 范数联合正则化的主成分追踪模型，此模型为非凸模型，它比凸的PCP 模型更接近于原RPCA 模型；2）针对新的非凸模型，基于ADM 方法，提出了一个快速求解的算法； 3）大量的人工数据和真实可视化数据上的实验验证了本文算法的有效性．2 L p 优化算法容易看出，p L 范数（01p <<）比1L 范数更接近于0L 范数，特别地，当0→p 时，p L 范数就退化成了0L 范数．同理，Schatten -p 范数比核范数更接近于秩函数．因此，我们考虑以下基于Schatten -p 范数与p L 范数联合正则化的主成分追踪模型：,,min ||||||||p p p p p A EA E λ+，..s t D A E =+， （3） 其中n m RD ⨯∈（不妨假设）是一个观测数据矩阵，p A ||||表示矩阵A 的Schatten -p 范数，即p mi p i p A A /11))((||||∑==σ，)(A i σ为矩阵A 的第i 个奇异值，p p E ,||||表示矩阵E 的p L 范数，即p m i nj p ij p p e E /111,)||(||||∑∑===．有很多非凸惩罚函数可以用到非凸主成分追踪模型中，这里用p L 范数和p Schatten -范数是因为它们与1L 范数和Schatten -p 范数具有类似的性质．下面运用ADM 算法解决上述非凸优化问题（3）．首先引入拉格朗日乘子矩阵Y ，并给出问题（3）的拉格朗日函数：2,||||2,||||||||),,(F p p p p p E A D E A D E A Y E A L --+--++=μλ， （4）刘园等：基于RPCA 模型的P 范数优化算法27其中0>μ表示给定的参数，F A ||||表示矩阵A 的Frobenius 范数，即A A A F ,||||=．为了优化上述拉格朗日函数，下面给出基于交替方向的迭代求解算法：()(1)()()(1)(1)()(1)()(1)(1)arg min (,,) arg min (,,)k k k A k k k Ek k k k A L A E Y E L A E Y Y Y D A E μ++++++⎧=⎪⎪=⎨⎪=+--⎪⎩． （6） 具体地，对于子优化问题（5），根据文献[6]有：(1)()()arg min (,,,)k k k AA L A E Y μ+==()()()2arg min ||||,||||2p k k k p F A YD AE D A E μ+--+--= ()21arg min ||||||||2k p F p G A A η-+(diag())T UT V η=∑ ， 其中μη1=，μ)()()(k k k Y E D G +-=，U ~和V ~分别是)(k G 的左奇异特征向量矩阵和右奇异特征向量矩阵，12diag(,,,,0,,0)r λλλ∑=⋅⋅⋅⋅⋅⋅为)(k G 的奇异值矩阵，其中),,2,1(r i i ⋅⋅⋅=λ为G G T 的特征值，T k V U G~~)(∑=，ηT 为软阈值算子，即0||(){0,sgn()}||sgn()||a a a aif z h T z z if z h z if z h δββ*<⎧⎪==⎨⎪>⎩， 其中pa p --=21)]1(2[δβ，1-+=p a a a p h βδβ，|)|,(*z a ββ∈是以下非线性方程10q q z βλββ-+=>的较大的根，可以通过迭代算子(1)()()k k βρβ+=，1()q z q ρβλβ-=-求解．类似，对于子优化问题（6）有：(1)(1)()arg min (,,,)k k k EE L A E Y μ++==()(1)(1)2,arg min ||||,||||2p k k k p pF E E Y D A E D A E μλ+++--+--=()2,1arg min (||||||||)2k p F p p E Z E E μδ-+2(),,1arg min ()||2k p ij ijij E i j i jz e e δ=-+=∑∑ ()()2(),1arg min [()||]2k p k ij ij ij ij E i jz e e T z δδ-+=∑，其中μλδ=，μ)()1()(k k k Y A D Z +-=+．本文方法的具体步骤如算法1所述．在算法1中，μ首先被初始化为一个较小的数0μ，然（5）温州大学学报（自然科学版）(2016)第37卷第4期28后每一次迭代都几何式增大，直到它达到一个事先设定的较大阈值μ．这种连续性处理技术能够很大程度上加快算法的收敛速度[7]．算法1：求解问题（3）的基于交替方向法的迭代算法 输入：观测数据矩阵n m R D ⨯∈，正则参数λ． 初始化：0A ；0E ；00>μ；1>ρ；0=k ． 当不收敛时作以下迭代步骤：1）固定其它变量为最新值，由（5）式更新变量A ． 2）固定其它变量为最新值，由（6）式更新变量E ． 3）固定其它变量为最新值，更新拉格朗日乘子()k Y ；令1k k μρμ+=．输出：),(k k E A ．3 实验结果实验分三个部分对本文模型与PCP 模型进行对比．首先，由于人工数据能够提供真实无噪音的数据（Ground truth ），因此可以对矩阵恢复模型的精度在量上进行对比；其次，给出它们在单张图片恢复上的对比；最后在人脸识别数据集YaleB 上恢复带有椒盐噪声的人脸图像，并用恢复后的数据进行分类，比较分类效果． 3.1 人工数据首先生成一个秩为r 的n m ⨯数据矩阵作为Ground truth ，记为0A ，然后生成一个稀疏矩阵，记为0E ，其中0E 中的非零元素的位置是随机的，非零元素均匀分布在区间[-500,500]，最后令00E A D +=，并记A ˆ和E ˆ分别为算法计算出来的数据恢复矩阵和稀疏噪音矩阵．在以下实验中，取m =100，n =100，r =10，取相对误差2020||||/||ˆ||A A A Err -=作为算法精度的度量标准．为了消除生成数据的随机性因素对算法误差的影响，对每一组参数都做20次实验，取误差的均值来度量模型的精度．首先，通过实验来观察参数p 的变化对实验效果的影响情况．固定噪音数据所占的比例为ratio =0.3，即稀疏噪音矩阵非零元素的个数为ratio m n ⨯⨯，取参数(p 0.1p ≤≤0.9)．实验结果如图1所示，红色线表示本文所提出的模型在不同参数p 下的相对误差，蓝色线表示PCP 模型的相对误差（由于PCP 模型中没有参数p ，因此它的相对误差不变）．从图1可以看出，p 值在（0,1）范围内变化时，PCP 模型的相对误差在0.15左右，而本文所给出模型的相对误差在0.005之下．更重要的是，本文模型的相对误差基本不受参数p 的变化影响．鉴于所提出的模型对于参数p 的鲁棒性，在以下的实验中固定p =0.5．进一步，通过实验给出两种方法在噪音数据比例ratio(0.05ratio ≤≤0.5)下的相对误差．结果如图2所示，两模型的相对误差都随着噪音比例的增大而逐渐增大．噪音比例在0.05 – 0.2之间时，两种方法的相对误差差别不大，但当噪音比例在0.3 – 0.5之间时，PCP 模型的相对误差急剧增加，而本文模型随着噪音比例的增加，其相对误差却一直趋于平稳状态，并且误差一直保持在0.001以下．表1给出了两种模型相对误差的数值结果，从表中可以明显看出，在整个噪音比例范围内，本文p L 模型的相对误差都远远小于PCP 模型的相对误差．刘园等：基于RPCA 模型的P 范数优化算法29表1 PCP 模型与本文模型的相对误差比较相对误差相对误差噪音比例 PCP 模型 本文模型 / 10-5噪音比例 PCP 模型 本文模型 / 10-50.05 0.000 3 1.32 0.30 0.180 0 5.68 0.10 0.001 2 1.92 0.35 0.422 0 8.87 0.15 0.006 6 2.78 0.40 0.729 0 9.72 0.20 0.022 0 4.16 0.45 1.204 7 13.70 0.250.063 05.010.502.289 116.093.2 单张图片数据下面将给出两种模型在单张图片恢复上的效果比较．首先选取无明显噪声污染的图片作为原始观测图像，然后在这些图像上加入噪音比例ratio =0.3的“椒盐噪声”作为待恢复的图片，最后通过不同的模型将噪音滤除掉，恢复原有的图像．在这里，选取不同风格的四幅图像作为原始观测图像，其中包括灰度图片和RGB 彩色图片，如图3（a ）所示．对RGB 彩色图片，我们用模型分别恢复出R 通道、G 通道、B 通道的图像，然后再将它们合并起来，即为恢复出的彩色图片．从图3可以看出，通过PCP 模型恢复出来的图像中还含有很多噪音点，看起来比较模糊，与原始数据有一定的差距，而用本文模型恢复的图像更接近于原始观测图像，看起来较为清晰．因此，本文所提出的模型在单张图片恢复方面要优于PCP 模型． 3.3 人脸识别在YaleB 数据集[8]上进行实验．数据集包含38个人的2 414幅人脸照片，其中每个人又包含了60多幅不同光照、不同表情下的照片，每张照片是32 × 32像素的灰度图像．实验设计如下：1）对于每个人，随机选取30张人脸照片作为训练集样本，剩余的作为测试样本；2）在每幅测试样本图片上加入噪音比例为ratio(0.05ratio ≤≤0.5)的椒盐噪声，并将其排列为1 024维的列向量；3）将训练样本和加噪声后的测试样本按列排成一个数据矩阵，记为D ； 4）用PCP 模型和本文模型分别对D 进行矩阵恢复； 5）用K 近邻（KNN ）方法对恢复后的测试样本进行分类．图1 不同p 值下本文算法与PCP相对误差的比较图2 不同噪音比例中本文算法与PCP相对误差的比较温州大学学报（自然科学版）(2016)第37卷第4期30图4给出了5个人的部分人脸图像，每一行是一个人脸不同表情的图像，（a ）列为原始人脸数据图像，伴随有严重的光照阴影，（b ）列为增加噪音比例ratio = 0.3的“椒盐噪音”后的人脸图像，（c ）列是用本文模型恢复出来的人脸图像，（d ）列是用PCP 模型恢复出来的人脸图像．(a) 原始观测图像(b) 加了椒盐噪音后的图像(c) 用本文模型恢复出来的图像(d) 用PCP 模型恢复出来的图像图3 单张图片数据集上的实验结果(a) 原始数据图像(b) 加了椒盐噪音后的图像(c) 用本文模型恢复出来的图像(d) 用PCP 模型恢复出来的图像图4 人脸数据集上的实验结果刘园等：基于RPCA 模型的P 范数优化算法31从图4中可以看出，PCP 模型仅仅能够去除部分光照阴影，本文的模型能够去除几乎所有的光照阴影，更为重要的是，PCP 模型恢复出来的图像仍然伴随着许多的噪音点，使得恢复的图像显得模糊不清，而本文模型能够较好地去除原始图像中的椒盐噪音，效果更接近于真实人脸图像．图5给出了PCP 模型与本文模型在人脸识别正确率方面的实验结果．可以看到，随着噪音比例的增加，直接用“椒盐噪声”遮挡的图像进行分类所得到的正确率急剧下降，用通过PCP 模型处理后的图像进行分类所得的正确率从0.55逐渐下滑到0.42左右，而用本文模型处理后的图像进行分类的识别正确率从0.6平缓下滑到0.55左右，远远高于实用PCP 模型的正确率．表2给出了20次实验的识别正确率的均值以及方差，可以看出，在噪音比例范围内，采用本文模型得出的图像识别正确率远远高于采用PCP 模型的．4 结 论本文通过对RPCA 问题的研究，将其中的秩范数与0L 范数替换成非凸的Schatten -p 范数和p L 范数（01p <<），为低秩矩阵恢复提出一种更接近于RPCA 问题的非凸近似模型，并针对此模型，给出了一种有效的求解算法．在人工数据集和真实数据集上的实验结果验证了本文所给模型的有效性．参考文献[1] Ganesh A, Ma Y , Rao S, et al. Robust principal component analysis: exact recovery of corrupted low-rank matricesvia convex optimization [J]. Advances in Neural Information Processing Systems, 2009, 87(4): 2080-2088.[2] Candes E J, Romberg J, Tao T. Robust uncertainty principles: exact signal reconstruction from highly incompletefrequency information [J]. Ranaon on Nformaon Hory, 2006, 52(2): 489-509.[3] Lin Z, Chen M, Ma Y . The augmented lagrange multiplier method for exact recovery of corrupted low-rank matrices表2 PCP 模型与本文模型对于人脸识别的正确率均值及方差比较噪音图像RPCA 算法本文算法噪音比例 识别正确率 方 差 识别正确率 方 差 识别正确率 方 差 0.05 0.532 5 0.005 7 0.542 9 0.002 7 0.600 0 0.003 3 0.10 0.457 0 0.012 1 0.532 7 0.001 5 0.595 8 0.002 0 0.15 0.358 9 0.011 6 0.523 4 0.002 0 0.592 3 0.003 1 0.20 0.296 4 0.013 7 0.504 4 0.002 7 0.588 7 0.002 5 0.25 0.233 8 0.005 9 0.493 4 0.004 4 0.579 9 0.003 1 0.30 0.175 4 0.009 7 0.478 0 0.002 4 0.569 7 0.003 3 0.35 0.135 5 0.005 4 0.459 7 0.005 0 0.557 3 0.005 0 0.400.100 20.005 90.423 90.004 80.544 70.004 2图5 PCP 模型与本文模型对人脸识别的正确率对比温州大学学报（自然科学版）(2016)第37卷第4期32[J]. arXiv preprint arXiv, 2010, 23(5): 109-115.[4] Chen X, Xu F, Ye Y . Lower bound theory of nonzero entries in solutions of l2-lp minimization [J]. Siam Journal onScientific Computing, 2010, 32(5): 2832-2852.[5] Candès E J, Wakin M B, Boyd S P . Enhancing sparsity by reweighted l 1 minimization [J]. Journal of Fourier Analysisand Applications, 2008, 14(5-6): 877-905.[6] Marjanovic G , Solo V . On l q optimization and matrix completion [J]. IEEE Transactions on Signal Processing, 2012,60(11): 5714-5724.[7] Toh K C, Y un S. An accelerated proximal gradient algorithm for nuclear norm regularized least squares problems [J].Pacific J Optimization, 2010, 6(3): 615-640.[8] Georghiades A S, Belhumeur P N, Kriegman D J. From few to many: illumination cone models for face recognitionunder variable lighting and pose [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2001, 23(6): 643-660.The Study of P -norm Optimization AlgorithmBased on RPCA ModelLIU Y uan, W ANG Di(College of Mathematics and Information Science, Wenzhou University, Wenzhou, China 325035)Abstract: The model of low-rank matrix recovery has been widely used in image in painting and video processing field. RPCA model is a classical model of low rank matrix recovery. Its basic idea is to decompose a numerical matrix into the form of a low-rank matrix and a sparse matrix sum before it is solved. However, RPCA problem is NP-hard. A commonly-used processing mode is to loose the matrix rank function and the 0L norm from the RPCA model respectively into the matrix nuclear norm and the 1L norm. Thus, theproblem is solved through the way to transform its approximation into the convex optimization problem. However, the solution of the convex optimization approximation approach degrades the solution of original problem under the relatively weak incoherent condition. This paper firstly proposes a non-convex approximation model toward the problem, namely, applying the matrix Schatten -p norm and p L norm (01p <<) to replace the matrix rank function and 0L norm. And then the further effective optimization algorithm is given based on the non-convex approximation model. Lastly, the experiment based on the synthetic datasets and the real image datasets is made. It turns out that the model proposed is valid and effective.Key words: Low-rank Matrix Recovery; RPCA Model; Schatten -p -norm; p L -norm(编辑：王一芳)。