卷积鲁棒主成分分析

王 心，朱浩华，刘光灿

（南京信息工程大学自动化学院，南京 210044）

（*通信作者电子邮箱xinwang@nuist.edu.cn）

0 引言

随着大数据技术的发展，大量的数据资源不断涌现，但由于目前数据的获取没有约束，所以获得的可观测数据带有大量噪声，例如数据遭到严重腐蚀或者数据含有离群点。总的来说，这些误差会显著降低数据样本的代表性，从而严重扭曲数据分析，因此如何恢复目标数据成为了一个十分重要的问题，也带来了很大挑战。Candès 等［1］和Zhang 等［2］提出了鲁棒主成分分析（Robust Principal Component Analysis，RPCA）算法，它又称作主成分追求（Principal Component Pursuit，PCP）。假设通过可观察的数据矩阵M获得先验的相关信息，然后根据先验信息建立合适的模型，同时设计模型的优化算法，最后求解得到恢复的数据矩阵。RPCA 采用低秩模型，将数据矩阵分解为1 个低秩矩阵和1 个稀疏矩阵，即M=L+E，那么利用RPCA，在满足一些条件时就可以精确地恢复出L和E。广义鲁棒主成分分析（Generalized Robust Principal Component Analysis，GRPCA）［3］是RPCA 的拓展。当数据被两种以上的混合噪声污染时，RPCA处理结果并不理想，而GRPCA通过最小化核范数、l1范数和l21范数的组合问题，可以分离出被混合噪声污染的低秩矩阵。RPCA 和GRPCA 的应用广泛［3-4］，但是RPCA 和GRPCA 算法在解决数据结构非低秩性时效果会很差，即无法分离出1 个低秩矩阵L。为了处理L不满足低秩的问题，目前已有的解决方法是假设L在特征映射后是低秩的，这意味着L在某些特征空间是潜在低秩的，而在原始空间中本身的秩比较高甚至是满秩的。相关的代表算法有核主成分分析（Kernel Principal Component Analysis，KPCA）。KPCA 是传统主成分分析（Principal Component Analysis，PCA）的扩展，它主要寻求核空间中数据点的近似低秩。与传统PCA 相似的是，即使在映射后，KPCA对离群点（outliers）也非常敏感，如果数据存在大量噪声或者离群点，KPCA 处理问题的效果会非常差，因此一些鲁棒的核低秩算法进一步被提出。如文献［5-6］等的研究中就提供了一种处理子空间聚类的核低秩方法，并且验证了核低秩近似确实有利于非线性数据的聚类。尽管这些方法在聚类和线性低秩恢复任务上取得了很大的成功，但是它们都无法恢复原始空间中的非线性或超低维数据，所以本质上不能直接应用于非线性数据恢复问题。目前Xie等［7-8］提出了一种更加鲁棒的核低秩算法来处理非低秩数据恢复问题，该算法基于RPCA，考虑当L从原始空间映射到可再生核希尔伯特空间H中能得到Φ(L)，此时非线性的观测点可以被认为是线性的，因此Φ(L)此时是低秩的。噪声矩阵E是逐列稀疏的，因此采用l21范数作为E的约束。在求解这个模型时，采用了核函数，例如凸的多项式核函数以及非凸的高斯核，采用核函数的目的主要是为了方便计算，降低计算难度。Xie 等［7-8］的算法相比之前的KPCA 算法结果更鲁棒，但是算法的计算复杂度高，计算过程中损耗较大，且仅解决了矩阵恢复问题，对于视频恢复即张量处理没有提及。本文提出了一种新的处理非线性数据的算法——卷积鲁棒主成分分析（Convolution Robust Principal Component Analysis，CRPCA），不仅可以处理矩阵恢复问题，也可以应用于视频恢复，这是其他算法无法比拟的。CRPCA 保留了RPCA 本身的优良性质，同时将RPCA 与卷积矩阵（Convolution Matrix，CM）结合，不需要将非线性的数据映射到某个特征空间，而是利用卷积矩阵的低秩性，即A(L)是低秩的，同时考虑噪声E是稀疏的，然后使用乘子交替方向法（Alternating Direction Method of Multipliers，ADMM）［9-12］来求解目标函数，最后将CRPCA 应用于去噪实验，结果表明，相对于RPCA、GRPCA、KRPCA 算法，CRPCA 算法的去噪效果更好、复杂度更低、计算损耗更少。

1 主要符号及公式

大写字母表示1 阶或更高阶的张量，包括向量、矩阵和高阶张量。在大多数情况下，本文考虑矩阵问题（即2 阶张量），对于矩阵X，X：，j表示矩阵X的第j列，Xi，：表示矩阵X的第i行，矩阵的奇异值分解为X=UΣVT，‖X‖*=Σiσi(X)表示矩阵的核范数，也就是矩阵X的奇异值之和。‖X‖1=Σij‖Xij‖表示X的l1范数，表示X的l21范数。两个其中XT是矩阵的转置，tr(·)是矩阵的迹，I定义为单位矩阵。手写体字体，例如A 表示线性算子，I 表示单位算子，对于希尔伯特空间之间的线性算子L：H1→H2，它的Hermitian 伴随（或共轭）算子定义为L*，并且有以下定义：矩阵的欧氏内积为

2 鲁棒主成分分析

鲁棒主成分分析（RPCA）很早就被提出，但至今依然受大家关注。在许多实际应用中，假设所给定的很大的数据矩阵M，它可以分解成：

其中：L是低秩的，E是稀疏的。此时并不知道L的列空间和行空间，甚至不知道它们的维数。同样，也不知道E中非零项的位置以及个数。为了恢复矩阵当矩阵M的低秩结构L，本文考虑当E中元素服从独立同分布的高斯分布，那么可以使用经典的PCA来求解这个问题，即求解以下最优化问题：

对矩阵M进行奇异值分解（Singular Value Decomposition，SVD）可以得到上述优化问题的最优解。但是当M有大量噪声时，传统PCA方法无法解决这个问题。此时优化问题为：

该优化问题是NP 难问题，无法在多项式时间内求解，所以对此优化问题的目标函数进行凸松弛。矩阵的核范数是秩函数的凸替代，矩阵的l1范数是l0范数的凸包，故将上述问题转化成以下凸优化问题：

Candès 等经过研究证明了在满足一定的条件下，矩阵的低秩部分L和稀疏部分E可以被准确地恢复出来。

3 卷积鲁棒主成分分析

卷积鲁棒主成分分析（CRPCA）是在鲁棒主成分分析（RPCA）的基础上加入卷积矩阵，解决了RPCA 无法求解数据不满足低秩结构的问题，在处理非线性数据上有很好的效果。

3.1 卷积矩阵

在信号处理中，（离散）卷积是最基础、最重要的概念。它的定义虽然大多是唯一的，但是也有很多种不同的解释，这取决于使用的边界条件。本文考虑的是循环卷积，即考虑循环边界条件的卷积［13-14］。循环卷积过程就是将X∈Rm和K∈Rk(k≤m)转换成X*K∈Rm：

其中：*表示卷积算子，假设[X]i-s=[X]i-s+m(i≤s)，这就是循环边界条件。在本文一直假设k≤m，并称K为卷积核。普遍来说，卷积算子是线性的并且可以转换成矩阵相乘的形式：

其中：Ak(X)是卷积矩阵，下标k表示卷积矩阵总是与核大小k相关。根据循环卷积，向量X=[x1，x2，…，xm]的卷积矩阵是大小为m×k的截断循环矩阵：

换句话说，Ak(X)的第j列就等于Sj-1(X)，S是循环移位运算符：该移位算子可以由Matlab 中circshit 函数实现，在k=m的特殊情况下，卷积矩阵Ak(X)是m×m的循环方阵。vec(·)是将X*K向量化，因此，可以得到：

3.2 基本性质

考虑最简单的例子，向量X∈Rm，那么它的卷积矩阵为Ak(X)∈Rm×k，则：

3.3 潜在低秩结构恢复

当n=1并且L为向量时，卷积矩阵Ak(L)的第j列其实是L中元素循环移位到j-1 个位置。准确地说，当L具有基本连续性，且移动程度相对较小时，循环移动前后的信号大多是低秩的，所以低秩的卷积矩阵将存在。为了处理数据潜在低秩问题，解决以下凸优化问题：

这个优化问题是凸的，可以被不同的方法解决。为了计算的有效性，在本文中采用乘子方向交替法（ADMM）。首先将式（9）转换成下列等式问题：

其增广拉格朗日函数为：

上述的问题是无约束的。通过固定其他变量，并分别对变量R、L、E求导，然后更新拉格朗日乘子Y1和Y2，μ＞0为惩罚项参数。本文算法概述了ADMM 求解凸优化问题（9）。在更新变量时利用了奇异值阈值（Singular Value Thresholding，SVT）算子和收缩算子（Shrinkage Operator，SO）［15-16］求解。

本文算法的实现过程如下所示：

1）输入数据M和参数λ。

2）对相关参数进行初始化R=L=0，E=0，Y1=Y2=0，μ=10-6，μmax=1010，ρ=1.05，ε=10-6，iter=0，Itermax=103。

3）固定其他变量，更新R：

4）固定其他变量，更新E：

5）固定其他变量，更新L：

6）更新拉格朗日乘子：

7）更新参数：

8）如未达到收敛条件

4 实验结果及分析

实验测试了卷积鲁棒主成分分析用于图像去噪、视频去噪的能力，分别在合成数据、MNIST、COIL-2 数据集以及DeepVideoDeblurring 视频数据集上进行测试，显示了本算法对复杂数据处理的高效性。

4.1 实验设置

1）基线（Baseline）：本文比较了几种方法，包括鲁棒主成分分析（RPCA）、核鲁棒主成分分析（KRPCA）、广义鲁棒主成分分析（GRPCA）来评估模型的性能。

2）评估指标（Evaluation Metric）：峰值信噪比（Peak Signal-to-Noise Ratio，PSNR），单位为dB。

其中：X0、Xrec∈Rm×n分别为原始数据和恢复数据。

4.2 算法复杂度分析

假设输入矩阵为M∈Rm×n，在RPCA 计算过程中单步迭代计算复杂度为O(m2n)，并且需要O(1/ε)的迭代次数以达到准确度ε。KRPCA 中计算奇异值分解以及n的立方根的总计算复杂度为O(n3+rn2)。在求解GRPCA 算法中，每一次循环中，奇异值分解的运算量为O((m+n)3)，矩阵求逆的运算量为O((m+n)3)。此外GRPCA 算法循环内的多次乘法加法和截断阈值的运算量为O(5n(m+n)+14n)。因此，GRPCA 算法的计算复杂度为O((m+n)3+I()5n(m+n)+14n)，其中I为循环的次数。

4.3 数据去噪

1）本文先对合成数据进行实验，图1 展示了向量恢复的结果，即对一维数组处理结果。实验中随机选取了100 个点，形成f(x)=sinx离散函数；再选取20 个点作为离群点。在这个实验中，恢复出的数据L十分接近真实值truth，验证了本文算法的准确性和有效性。

图1 对向量的恢复效果Fig.1 Effect of vector recovery

2）表1 以及图2 使用MNIST 数据集和COIL-20 数据集，黑白图片作为二维数组，可视为矩阵。该部分验证了算法对图片的恢复能力，即恢复潜在矩阵低秩的高效性。可以看出本文算法效果相较于其他三种算法，结果更好，恢复的准确度更高。

图2 矩阵恢复的效果Fig.2 Effect of matrix recovery

表1 不同方法在MNIST和COIL-20数据集上峰值信噪比的对比 单位：dBTab.1 Comparison of PSNR of different methods on MNIST and COIL-20 datasets unit：dB

3）张量是多维数组，视频比图片多了时间信息，因此视频可视作张量。由于摄像头的晃动，导致拍摄的视频背景变化，因此无法将背景看作低秩的，无法将前景从背景中分离出来。该部分利用DeepVideoDeblurring 视频序列，同时加入噪声作为输入。而RPCA、KRPCA、GRPCA 无法对视频（即张量）整体进行操作，图3展示了本文算法对视频序列的处理。

图3 视频恢复的效果（按顺序截取视频的前12帧）Fig.3 Effect of video recovery（capturing 12 frames of the video in sequence）

5 结语

本文介绍了一种新的去噪算法卷积鲁棒主成分分析（CRPCA）。在解决数据潜在低秩的情况下，卷积鲁棒主成分分析利用数据矩阵将数据转换为低秩，并与鲁棒主成分分析（RPCA）结合，同时采用乘子交替方向法（ADMM）进行计算。相较于传统的鲁棒主成分分析（RPCA）和核鲁棒主成分分析（KRPCA）以及广义鲁棒主成分分析（GRPCA），在解决非线性数据及去噪问题上，由表和图可以看出，本文算法效果更好，恢复出的数据更接近原始数据，在对向量、矩阵去噪处理上均取得了不错的效果；同时，解决了其余算法无法处理的视频问题，即张量去噪问题。除了去噪问题外，本文算法也可用于图像的去模糊问题。