一种基于快速随机投影的矩阵填充方法

冯雅莉 孙为军

(广东工业大学 广东 广州 510006)

0 引 言

近年来，矩阵的使用越来越广泛和越来越受重视，而矩阵分析方法的重要性也显得越来越突出，例如矩阵分解、非负矩阵分解、低秩矩阵近似、矩阵填充等。而本文主要是关注其中的矩阵填充问题，在实际问题中，矩阵填充问题主要表现在图像和视频处理、推荐系统、文本分析和机器学习[1-2]等方面。为了方便对数据进行处理，一般我们会把数据转化成矩阵的表示形式，但这些数据转为矩阵后，经常会面临缺失，损毁和噪声污染等问题，所以如何恢复数据和对数据进行去噪是我们需要解决和研究的问题，而矩阵填充技术就是解决这类问题的其中一种有效的方法。矩阵填充技术在图像修复[3]、核磁共振图像分割[4]、EEG信号处理、推荐系统[5]等领域中被广泛应用。

在近二十年里，矩阵填充技术越来越受关注。它的主要目的是通过矩阵中已知的元素来恢复未知的元素，通常利用矩阵中行与行或者列与列之间的线性关系对未知元素进行估计和恢复，所以需要被填充的矩阵是低秩的，这是矩阵填充的重要前提之一。近年来，在矩阵填充领域中涌现了许多有效的算法，例如近邻梯度下降算法[6](proximal gradient descent,PGD)、分块坐标下降算法[7]、矩阵空间Bregman迭代算法[8]、交替方向乘子法[9](alternating direction method of multipliers,ADMM),还有一些随机优化方法，例如随机近邻梯度下降算法(stochastic proximal gradient descent,SPGD)。

在解决实际的矩阵填充的问题中，对矩阵进行降维能够大大减少矩阵填充在计算中的成本，常用的方法有主成分分析(principal component analysis,PCA)和奇异值分解(singular value decomposition,SVD)。大多矩阵填充方法都是利用传统的SVD方法对原始矩阵进行降维的，但是传统的SVD的操作计算成本很高，所以不适用于大规模大尺寸的矩阵中。本文针对这个问题，提出了一种快速的随机投影方法，同时利用人工仿真数据和实际图片像素数据进行了多个实验，验证了本算法的可行性。

1 基础知识

对于秩为r的矩阵X∈Rm×n，它的奇异值分解为：X=UΣVT，其中U∈Rm×r和V∈Rn×r满足UTU=VTV=I。(I为r阶单位矩阵)，Σ=diag(σ1,σ2,…,σr)且其对角线上的元素满足σ1≥σ2≥…≥σr>0，σi表示X的第i个奇异值。

X的Frobenious范数满足：

式中：〈X,X〉表示X与X的内积，trace(·)表示矩阵的迹算子。

2 快速矩阵填充方法

2.1 矩阵填充模型

秩是矩阵的一个重要的全局信息[9]，它能体现出矩阵的信息冗余度，而矩阵能够被充分地填充的前提条件是矩阵的信息是冗余的，其冗余度越低，填充难度就越高。所以矩阵填充的问题一般是最小化矩阵的秩的问题，其模型如下：

(1)

s.t.XΩ=MΩ

式中：X,M∈Rm×n，M可观察的元素的索引存放在Ω集合中，而索引不在Ω集合中的元素均为未知的。由于函数rank(X)是非凸的，在实际优化中求解非常困难。文献[10]提出对rank(X)的最小化问题可以转化成对X的核范数最小化问题，模型如下：

(2)

s.t.XΩ=MΩ

在过去使用核范数优化模型的矩阵填充方法很多，例如奇异值阈值算法(Singular Value Thresholding,SVT)[11], 加速近端梯度算法(Accelerated Proximal Gradient,APG)[12]，增广拉格朗日乘子法(Augmented Lagrange Multipliers,ALM)。

本文引入一个新的变量Z，并把式(2)改写成：

(3)

rank(Z)≤rank(X)

2.2 随机快速矩阵填充方法

解式(3)，一般采用SVD的方法，但是采用传统的SVD的方法，计算的时间和计算复杂度都非常大和高，这不利于对大尺寸的矩阵进行填充和处理，所以本文提出了一种快速的随机投影方法。

假设一个矩阵X∈Rm×n，Ω和ΩC分别表示X矩阵中的已知和未知元素的索引集合。为了解式(3)，我们可以采用交替最小二乘法对Z进行更新。

本文方法分为以下几个部分：

第一部分是对数据进行预处理，初始化带未知元素的矩阵X，把矩阵X的未知元素用独立同分布的均值为0方差为1的高斯随机分布赋值，如下所示：

第二部分，如果对X进行传统的SVD操作，则计算复杂度为Ο(min(m2n,mn2)),随着矩阵的尺寸增大，计算的复杂度呈指数增长，所以需要对矩阵X进行降维，方法如下：

C=X×H

(4)

式中：H∈Rn×r，且H是由独立同分布的高斯随机分布组成，至于r是如何选择的可以参考文献[13]。

第三部分，就是求出X的近似左奇异值矩阵：

U=C×V×Σ-1

(5)

式中：V是利用对CTC求特征向量获得，Σ是通过对CTC求特征值开平方根获得。

第四部分，对X进行矩阵填充，或者是低秩近似。因为U在正常情况下是C的一个正交矩阵，所以UUT=I，而C是X的一个近似矩阵，所以我们以如下方式更新X：

X≈U×UT×X

(6)

第五部分是更新停止的条件的设定，本文采用均方差(mean squared error,MSE)来作为停止条件之一，定义如下：

(7)

当MSE小于或者等于我们设定的阈值时，或者迭代次数达到了我们设定的最大迭代次数时，迭代停止。

本文的随机投影算法的详细步骤如算法1所示。

算法1基于随机投影的矩阵填充算法

输入：带有未知元素的X矩阵(X∈Rm×n)，估计的秩r，阈值t，最大迭代次数M

输出：填充后的矩阵X

初始化：X中的未知元素都用独立同分布的均值为0方差为1的高斯随机分布赋值

循环： 迭代次数n=1:M

步骤1C=X×H，其中H∈Rn×r

i.i.d.N(0,1)

步骤2 [V,d]=eig(CTC)，其中V是C的右奇异矩阵，d是CTC的特征值组

步骤3U=C×V×Σ-1，其中

Σ=diag(diag(d).^0.5)

步骤4X=U×UT×X

当MSE≤t或者n=M时，迭代结束。

2.3 复杂度分析

算法1中，步骤1的计算复杂度为Ο(mnr),步骤2计算复杂度为Ο(nr+r2)，步骤3 的计算复杂度为Ο(2nr2)，步骤4的计算复杂度为Ο(m2r+m2n)，所以本文算法迭代一次的计算复杂度大约为Ο((n+r)m2+mnr+2nr2+nr+r2)。

3 实验分析

将本文算法与SVT[11]、APG[12]和ALM[14]三个经典的算法进行比较，由于精确增广拉格朗日乘子法(Exact Augmented Lagrange Multipliers，EALM)的计算成本较高，所以本文仅对ALM算法中的非精确增广拉格朗日乘子法(Inexact Augmented Lagrange Multipliers，IALM)作比较。实验中，我们假设需要填充的矩阵M∈Rm×n，其秩为r=(r1,r2)，实验的收敛条件为[15]：

(8)

或者是迭代次数达到设置的最大迭代次数。式中Vi表示我们需要迭代更新的第i个变量，ε表示收敛或者迭代停止的阈值，一般设定为10-5，本文中的最大迭代次数设置为100。

3.1 仿真数据实验

本节包含两个部分：第一部分是需要被恢复的矩阵的秩不同的情况下，四个算法的恢复情况；第二部分是在不同采样率下四个算法的恢复情况。对矩阵恢复有重要影响的两个参数为恢复矩阵的秩和采样率。我们用均方差(MSE)和运行时间来作为对矩阵恢复效果的衡量标准。

第一部分，假设矩阵M∈R100×100，在其秩从2递增至20间设计了10组实验，采样率p=0.5。实验结果如图1所示。由图1(a)可知，当矩阵的秩在2～20之间时，FRPMC方法的对矩阵的恢复的效果仅次于ALM算法，但图1(b)中，FRPMC的恢复的速度比ALM快了一倍。在仿真实验中虽然APG和IALM的运行速度相对比较快，但是它们矩阵恢复的相对误差较高，即恢复的精度相对较低。

第二部分，假设矩阵M∈R100×100，矩阵的秩为5，在采样率分别为0.4、0.6、0.8的条件下，用FRPMC、APG、IALM分别对矩阵进行恢复，结果如图2所示。由图2(a)可知，FRPMC方法的精度相对另外的两种算法要高，而且随着采样率的增加，恢复的精度越高，恢复效果越好，即FRPMC在高采样率的情况下，恢复效果比较显著。

(a)

(b)图1 不同算法在矩阵的秩不同的情况下矩阵恢复的结果

(a)

3.2 真实数据实验

本节主要利用对实际的缺失的黑白图片的恢复来说明FRPMC能够实现图像恢复的功能。因为黑白图片的像素相当于一个二维的矩阵。实验将FRPMC算法与SVT、APG、EALM、IALM进行比较，并用峰值信噪比(Peak Signal-to-Noise Ratio,PSNR)[17]作为衡量图像恢复的效果的完整性的标准，定义如下：

(9)

图3中的6幅图是分别在不同随机采样下生成的，采样率从左至右分别是0.3、0.4、0.5、0.6、0.7、0.8。图4的6幅图则是图3中相对应的用FRPMC算法恢复后的图像。图4(a)是采样率为0.3的情况下恢复后的图像，从图中可以清晰地看到房子的轮廓窗、房子在水中的倒影，说明FPRMC方法能够有效恢复带有丢失数据的图像。

(a) (b) (c) (d) (e) (f)图4 经过FRPMC处理后与图3相对应的恢复后的图像

接下来，随机抽取来自Berkeley Segmentation数据集中的50幅图来进行实验。首先把这50幅彩图进行灰度化，使其变成黑白图片，再进行实验。

对这50幅处理后的黑白图片进行随机采样，采样率分别是0.3、0.4、0.5、0.6、0.7、0.8，分别用FRPMC、SVT、APG和IALM四种方法进行实验。实验假设最大迭代次数为100。实验结果如表1和表2所示。

表1 在不同采样率下四个算法的PSNR值

表2 在不同采样率下四个算法的运行时间

由表1可知，在图像的尺寸为512×512的情况下，即FPRMC算法对黑白图像恢复的PSNR是最高的，即效果最好的。由表2可知，FRPMC算法在实际图像恢复中的速度也是最快的。

图5为在采样率为0.7的实验中随机选取的5幅图的实验情况。表3是图5中对应的5幅图的恢复结果，分别是PSNR值和各个算法对应的运行时间。

(a) 原始图片 (b) 加噪后图片(P=0.7) (c) FRPMC算法 (d) SVT算法 (e) APG算法 (f) IALM算法图5 随机选取5幅图的实验情况

表3 图5中5幅图相对应的PSNR值和对应算法的运行时间

4 结 语

传统的矩阵填充方法很难避免标准的SVD的计算操作，而随着矩阵维度的增加，SVD操作的计算复杂度呈指数增加，不适用于大规模、高维和高阶的矩阵的处理。如何设计可扩展的快速的算法是矩阵填充算法研究的核心。为了解决这个问题，本文提出了一种快速的随机投影方法，相对比于一般传统的算法。FRPMC算法主要是利用随机投影的方法对经过随机初试化的矩阵(带有缺失值)进行降维，然后再对其进行特征值和特征向量的计算，利用SVD的计算模型来设计了一个矩阵的低秩近似的模型，对矩阵进行恢复，大大减少了计算的成本。

在今后的矩阵填充算法研究中，还将关注以下几个方面：(1) 在确保精度的情况下，提高矩阵填充的速度；(2) 设计适合不同应用场景的矩阵填充的模型；(3) 将矩阵填充扩展到高维数组中，例如张量填充问题上；(4) 研究可以利用已知元素来确定或者估计矩阵的秩的方法。