Principal component analysis (PCA) is a fundamental technique of exploratory data analysis and has been widely used in many data mining tasks. Given a data matrix , the classic PCA seeks the best rank- approximation with complexity
. It is well known that PCA is sensitive to outliers. To combat this drawback, in the past decade, a number of approaches to robust PCA (RPCA) have been proposed, including alternating minimization, random sampling techniques [2, 3], multivariate trimming , and others [5, 6]. More recently, a new type of RPCA method has emerged and become popular [7, 8]. It assumes that can be separated into two parts, i.e., a low-rank and a sparse , by solving the following problem:
where is the (pseudo) norm which counts the number of nonzero elements of a matrix, and is a balancing parameter. Because minimizing the rank and the norm in (1) is generally NP-hard, in practice, (1) is often relaxed to the following convex optimization problem :
where is the nuclear norm of with denoting the -th largest singular value of , and is the norm. Theoretically, under some mild conditions, (2) can exactly separate with the true rank from . A number of algorithms have been developed to solve (2), including singular value thresholding (SVT) , accelerated proximal gradient (APG) , and two versions of augmented Lagrange multipliers (ALM) based approaches : exact ALM and inexact ALM (IALM). Among these algorithms, the ALM based are the state-of-the-art ones for solving (2), which need to compute SVDs of matrices per iteration. To improve efficiency, another ALM based algorithm adopts PROPACK package , which solves only partial SVDs instead of full SVDs. However, this is still computationally costly when and are both large. Despite the elegant theory of the convex RPCA in (2), it has two major drawbacks: 1) When the underlying matrix has no incoherence guarantee 
, or the data get grossly corrupted, the results may be far from the underlying true ones; 2) The nuclear norm may lead to biased estimation of the rank. To combat these drawbacks,  fixes the rank of as a hard constraint, while  uses a non-convex rank approximation to more accurately approximate the rank of , which needs to solve full SVDs. [14, 15] need only to solve partial SVDs, which significantly reduces the complexity compared to computation of full SVD; for example, AltProj has a complexity of . However, if is not known a priori, [14, 15] may fail to recover correctly.
To further reduce the complexity, enhance the scalability and alleviate the dependence on the knowledge of , in this paper, we propose a factorization-based model for RPCA, where can be decomposed as with , , , and . This model relaxes the requirement on a priori knowledge of the rank of , assuming only that it is upper bounded by . Scalable algorithms are developed to optimize our models efficiently. After acceptance of this paper, it came to our attention that  also proposed a factorization-based approach for such a problem. However, our model is distinct from . Our model has two variants, one of which uses an explicit rank constraint by a matrix factorization in the case that the ground truth rank is known, while the other uses non-convex rank approximation in the case that the ground truth rank is unknown. Both variants of our model differ starkly from , which only considers the second case and simply uses the nuclear norm.
We summarize the contributions of this paper as follows: We propose a factorization-based model for RPCA, allowing the recovery of the low-rank component with or without a priori knowledge of its true rank. In the absence of true rank, a non-convex rank approximation is adopted, which can approximate the rank of a matrix more accurately than the nuclear norm. Efficient ALM-type optimization algorithms are developed with scalability in both and . This is contrasted to AltProj whose cost is . This difference is important when and are large. Empirically, extensive experiments testify to the effectiveness of our model and algorithms both quantitatively and qualitatively in various applications.
Ii Related Work
where is defined to be the sum of
norms of column vectors of a matrix. When a matrix has large singular values, the nuclear norm may be far from an accurate approximation of the rank. Non-convex rank approximations have been considered in a number of applications, such as subspace clustering[19, 20]. A non-convex rank approximation has also been used in RPCA , which replaces the nuclear norm in (3) by a non-convex rank approximation , with , and being the -th largest singular value of . The above approaches usually need to solve SVDs. When the matrix involved is large, the computation of SVD, in general, is intensive. To reduce the complexity of RPCA, several approaches have been attempted. For example, AlgProj  uses non-convex alternating minimization techniques in RPCA and admits complexity.  uses a factorization approach:
which solves SVDs of thin matrices and hence admits scalability, where
is identity matrix with proper size.
Iii Fast Factorization-based RPCA
In this section, we formulate the Fast Factorization-based RPCA (FFP). We model the data as , where is a sparse component and is a low-rank approximation of with being the core matrix, , satisfying with being an identity matrix of a proper size. It is seen that the factorization provides a natural upper bound for the rank of the low-rank component of . The upper bound, which is , can be used to relax the stringent requirement on the knowledge of the true rank by AltProj algorithm. To capture the sparse structure of , we adopt norm to obtain sparsity. Thus, we consider the following objective function:
For many applications, the true rank of , which is , is known. In the presence of this prior knowledge, we let in 5. However, when this information on is not present, 5 may lack the desired capability of recovering with a (unknown) rank of while an arbitrary used. In this situation, we propose to marry the advantages of the convex and fixed-rank RPCA approaches, yielding the following optimization problem:
where is a balancing parameter. Even in the case that knowledge of the precise value of is not available, a proper can still be chosen because, in the worst case, we may let . Oftentimes, with domain information on the application at hand, an upper bound can be obtained. It has been shown that the nuclear norm can not approximate the true rank well if there are dominant singular values in a matrix, and non-convex rank approximations may help improve learning performance . Here, we adopt a log-determinant rank approximation , , to obtain the following RPCA model:
Due to the fact that
model 7 can be reduced to the following model:
respectively. The derivations for optimization are similar to those in . We summarize the optimization in Algorithm 1. Here, we define and to be the left and right singular vectors of a matrix from thin SVD, and define the operator , with
with , and .
Iv-a Complexity Analysis
Given that , both F-FFP and U-FFP have complexity .
To evaluate the proposed model and algorithms, we consider three important applications: foreground-background separation, shadow removal from face images, and anomaly detection. We compare our algorithms with the state-of-the-art methods, including IALM111http://perception.csl.illinois.edu/matrix-rank/sample_code.html#RPCA.  and AltProj222http://www.personal.psu.edu/nsa10/codes.html. , both of which make use of the PROPACK package  to solve SVDs for efficiency. All experiments in this section are conducted using Matlab on a dual-core Intel Xeon E3-1240 V2 3.40 GHz Linux Server with 8 GB memory. For purpose of reproductivity, we provide our codes on the website333https://www.researchgate.net/publication/308174615_codes_icdm2016.
V-a Foreground-background separation
Foreground-background separation is to detect moving objects or interesting activities in a scene, and remove background(s) from a video sequence. For this task, we use 15 datasets, as listed in the first column of Table I, among which the first 11 contain a single background while the rest 4 have 2 backgrounds444Datasets used in this subsection can be found at:
http://research.microsoft.com/en-us/um/people/jckrumm/wallflower/testimages.htm.. In this case, we have and for the first 11 and last 4 datasets, respectively. For each dataset, the data matrix is constructed by treating all vectorized frames as columns555For computational ease, down-sampling is performed on Camera Parameter, Highway, Office, Shopping Mall, Pedestrian, and Time of Day data sets.. In the following, we test F-FFP and U-FFP in two cases according to whether is known.
V-A1 Case 1 ( is known)
We set for F-FFP and AltProj. We terminate all methods after 200 iterations or when is reached. For IALM, we fix and for all these data sets for fast convergence and good visual quality. The balancing parameter is chosen as the theoretical one . For Altproj, the default parameters are used. For F-FFP, we use the same and as IALM. Without specification, the parameter settings remain the same throughout this paper.
The numerical results are reported in Table I. It is observed that IALM separates more sparsely but fails to recover with low rank, while F-FFP and AltProj recover properly. F-FFP generates more sparse than AltProj for most of the datasets. All these methods have competitive fitting error. However, it is possible for F-FFP to obtain more accurate fitting if the same iterations or time is provided as IALM or AltProj. Furthermore, F-FFP needs the least amount of time on all these datasets.
|Data||Method||Rank()||# of Iter.||# of SVDs||Time|
|Time Of Day||AltProj||2||0.8646||4.72e-4||44||46||61.63|
For IALM and AltProj, (partial) SVDs are for matrices. For F-FFP, SVDs are for matrices, which are computationally far less expensive than those required by IALM and AltProj.
V-A2 Case 2 ( is unknown)
is specified as a tight upper bound of based on domain knowledge on the video. In this test, we set on all datasets and compare U-FFP with AltProj and IALM666Here the results of IALM are not shown since they are the same as Case 1.. For U-FFP, is chosen from . We show the numerical results in Table II. It is observed that U-FFP is still able to recover with the true rank, whereas AltProj fails in this case. Besides, the time cost of U-FFP increases slightly by less than 1 second for most of the datasets while the time cost of AltProj increases by about 10-20 seconds.
We also show some video frames in Figure 1 and the visual results in Figure 2. It is observed that the backgrounds recovered by IALM still have shadows; AltProj separates foreground and background well with known, but the backgrounds are not clean when is unknown; both F-FFP and U-FFP can successfully separate foreground and background from the video.
|Data||Method||Rank()||# of Iter.||# of SVDs||Time|
|Time Of Day||AltProj||5||0.8651||4.61e-4||46||51||73.35|
For AltProj, (partial) SVDs are performed on matrices. For U-FFP, SVDs are for , , and matrices, which are computationally far less expensive than those required by AltProj.
V-B Shadow removal from face images
Learning face is an important topic in pattern recognition; however, shadows make it more challenging. There are often shadows on face images due to varying lighting conditions. Therefore, it is crucial to handle shadows, peculiarities and saturations on face images, so as to improve the learning capability on face image data. This can be handled with RPCA, because clean images reside in a low-rank subspace while the shadows correspond to sparse components. Following [13, 8], we use Extended Yale B (EYaleB) dataset  for this test. Out of 38 persons in this dataset, we choose the first two subjects. For each subject, there are 64 heavily corrupted images of size taken under varying lighting conditions. Each image is vectorized into a column of a data matrix. It is natural to assume that the face images of the same person reside in a rank one space, hence we have .
First, we consider the case where is known. Following [13, 8], IALM, AltProj, and F-FFP are applied to each subject and the quantitative and visual results are shown in Table III and Figure 3, respectively. It is observed from Figure 3 that AltProj and F-FFP can successfully remove shadows from face images. The majority of shadows can be removed by IALM, but some still remain. From Table III, we can see that AltProj and F-FFP can recover the low-rank component with exactly rank 1, while IALM recovers which has a higher rank. Then, we consider the case where is unknown. Following the setting in V-A, we set and apply AltProj and U-FFP to each subject. The quantitative and visual results are shown in Table IV and Figure 4, respectively. It is observed from Table IV that U-FFP has similar visual results as F-FFP while AltProj appears unable to remove shadows completely. Quantitatively as shown in Table IV, AltProj gives an that has a higher rank while U-FFP gives an that has the true rank. Besides, the proposed U-FFP is the fastest among all these methods.
|Data||Method||Rank()||# of Iter.||# of SVDs||Time|
|Data||Method||Rank()||# of Iter.||# of SVDs||Time|
V-C Anomaly Detection
Given a number of images from one subject, they form a low-dimensional subspace. Any image significantly different from the majority of the images can be regarded as an outlier; besides, fewer images from another subject can be regarded as outliers. Anomaly detection is to identify such kinds of outliers. USPS dataset contains 9,298 images of hand-written digits of size . Following , among these images, we select the first 190 images of ‘1’s and the last 10 of ‘7’s and construct a data matrix of size by regarding each vectorized image as a column. Since the number of ‘1’s is far greater than ‘7’s, the former is the dominant digit, while the images of the latter are outliers. The true rank of should be 1. Some examples of these selected images are shown in Figure 5. It is observed that besides ‘7’s, some ‘1’s are quite different from the majority. Therefore, anomaly detection, in this case, is to detect not only the ‘7’s, but also the anomaly of ‘1’s. After applying F-FFP, the columns in that correspond to anomalies contain relatively larger values. We use norm to measure every column of and show the values in Figure 6. These outliers can be identified by finding the columns with the highest bars. For ease of visualization, we vanish all those values that are smaller than 5 in Figure 6. The corresponding digits are shown in Figure 7 as outliers, which include all ‘7’s and several unusual ‘1’s.
To numerically illustrate the scalability of F-FFP and U-FFP, we test the time cost versus the values of and , respectively. To test how the computation time increases with , we uniformly sample a partition of the frames with the sampling rate in and keep all pixels for each frame. To test the relationship of time with respect to , we use all frames and down-sample each frame with different rates varying over , to keep the spatial information of each frame. We run F-FFP and U-FFP with 50 iterations in 10 repeated runs and we report the average time in Figure 8. It is observed that the time cost increases essentially linearly with both and for both F-FFP and U-FFP.
In this paper, we propose a new, factorization-based RPCA model. Non-convex rank approximation is used to minimize the rank when the ground truth is unknown. ALM-type optimization is developed for solving the model. Our model and algorithms admit scalability in both data dimension and sample size, suggesting a potential for real world applications. Extensive experiments testify to the effectiveness and scalability of the proposed model and algorithms both quantitatively and qualitatively.
Qiang Cheng is the corresponding author. This work is supported by National Science Foundation under grant IIS-1218712, National Natural Science Foundation of China, under grant 11241005, and Shanxi Scholarship Council of China 2015-093.
-  Q. Ke and T. Kanade, “Robust l norm factorization in the presence of outliers and missing data by alternative convex programming,” in Computer Vision and Pattern Recognition, 2005. CVPR 2005. IEEE Computer Society Conference on, vol. 1. IEEE, 2005, pp. 739–746.
-  R. Maronna, D. Martin, and V. Yohai, Robust statistics. John Wiley & Sons, Chichester. ISBN, 2006.
-  F. De La Torre and M. J. Black, “A framework for robust subspace learning,” International Journal of Computer Vision, vol. 54, no. 1-3, pp. 117–142, 2003.
R. Gnanadesikan and J. R. Kettenring, “Robust estimates, residuals, and outlier detection with multiresponse data,”Biometrics, pp. 81–124, 1972.
-  L. Xu and A. L. Yuille, “Robust principal component analysis by self-organizing rules based on statistical physics approach,” Neural Networks, IEEE Transactions on, vol. 6, no. 1, pp. 131–143, 1995.
-  C. Croux and G. Haesbroeck, “Principal component analysis based on robust estimators of the covariance or correlation matrix: influence functions and efficiencies,” Biometrika, vol. 87, no. 3, pp. 603–618, 2000.
-  J. Wright, A. Ganesh, S. Rao, Y. Peng, and Y. Ma, “Robust principal component analysis: Exact recovery of corrupted low-rank matrices via convex optimization,” in Advances in neural information processing systems, 2009, pp. 2080–2088.
-  E. J. Candès, X. Li, Y. Ma, and J. Wright, “Robust principal component analysis?” Journal of the ACM (JACM), vol. 58, no. 3, p. 11, 2011.
-  J.-F. Cai, E. J. Candès, and Z. Shen, “A singular value thresholding algorithm for matrix completion,” SIAM Journal on Optimization, vol. 20, no. 4, pp. 1956–1982, 2010.
-  K.-C. Toh and S. Yun, “An accelerated proximal gradient algorithm for nuclear norm regularized linear least squares problems,” Pacific Journal of Optimization, vol. 6, no. 615-640, p. 15, 2010.
-  Z. Lin, M. Chen, and Y. Ma, “The augmented lagrange multiplier method for exact recovery of corrupted low-rank matrices,” arXiv preprint arXiv:1009.5055, 2010.
-  X. Ding, L. He, and L. Carin, “Bayesian robust principal component analysis,” Image Processing, IEEE Transactions on, vol. 20, no. 12, pp. 3419–3430, 2011.
-  Z. Kang, C. Peng, and Q. Cheng, “Robust pca via nonconvex rank approximation,” in Data Mining (ICDM), 2015 IEEE International Conference on. IEEE, 2015, pp. 211–220.
-  W. K. Leow, Y. Cheng, L. Zhang, T. Sim, and L. Foo, “Background recovery by fixed-rank robust principal component analysis,” in Computer Analysis of Images and Patterns. Springer, 2013, pp. 54–61.
-  P. Netrapalli, U. Niranjan, S. Sanghavi, A. Anandkumar, and P. Jain, “Non-convex robust pca,” in Advances in Neural Information Processing Systems, 2014, pp. 1107–1115.
-  G. Liu and S. Yan, “Active subspace: Toward scalable low-rank learning,” Neural computation, vol. 24, no. 12, pp. 3371–3394, 2012.
-  H. Xu, C. Caramanis, and S. Sanghavi, “Robust pca via outlier pursuit,” in Advances in Neural Information Processing Systems, 2010, pp. 2496–2504.
-  M. McCoy, J. A. Tropp et al., “Two proposals for robust pca using semidefinite programming,” Electronic Journal of Statistics, vol. 5, pp. 1123–1160, 2011.
-  C. Peng, Z. Kang, H. Li, and Q. Cheng, “Subspace clustering using log-determinant rank approximation,” in Proceedings of the 21th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining. ACM, 2015, pp. 925–934.
C. Peng, Z. Kang, M. Yang, and Q. Cheng, “Feature selection embedded subspace clustering,”IEEE Signal Processing Letters, vol. 23, no. 7, pp. 1018–1022, July 2016.
-  R. Basri and D. W. Jacobs, “Lambertian reflectance and linear subspaces,” Pattern Analysis and Machine Intelligence, IEEE Transactions on, vol. 25, no. 2, pp. 218–233, 2003.
A. S. Georghiades, P. N. Belhumeur, and D. J. Kriegman, “From few to many: Illumination cone models for face recognition under variable lighting and pose,”Pattern Analysis and Machine Intelligence, IEEE Transactions on, vol. 23, no. 6, pp. 643–660, 2001.