Avoiding Optimal Mean Robust PCA/2DPCA with Non-Greedy l1-Norm Maximization

1 Introduction

High-dimensional data are frequently generated in many scientific domains, such as image processing, visual description, remote sensing, time series prediction and gene expression. However, it is usually computationally expensive to handle high-dimensional data due to the curse of dimensionality [Parsons et al., 2004; Chang et al., 2015; Nie and Huang, 2016; Nie et al., 2010]. Therefore, dimension reduction techniques are typically used to extract meaningful features from high-dimensional data without degrading performance. Among these methods, principal component analysis (PCA) learns a set of projections that constitute a low-dimensional linear subspace. It has been widely used in many applications for its simplicity and effectiveness [Jolliffe, 2002].

⇤Corresponding author. Email: feipingnie@gmail.com

Typically, standard PCA is based on mean square error, and thus it is disproportionately affected by the presence of outliers which occur often in high-dimensional data. For this issue, multiple robust PCA methods have been proposed to enhance the robustness of PCA by replacing the `2-norm with `1-norm distance [Wright et al., 2009; Torre and Black, 2001; Jolliffe, 2002; Ke and Kanade, 2005; Chang et al., 2016]. However, `1-norm based robust PCA methods usually perform worse due to their lack of rotational invariance and expensive computation. To solve this problem, a rotational invariant `1-norm based robust PCA, namely R1-PCA, has been proposed to soften the contributions from outliers by re-weighting each data point iteratively [Ding et al., 2006]. This method was further extended to its 2D version in [Huang and Ding, 2008; Yang et al., 2004]. However, the R1PCA/2DPCA models are solved with subspace iteration algorithm, which requires a lot of time to achieve convergence [Kwak, 2008]. Kwak proposed an intuitive method to ensure both the robustness and rotational invariance of PCA by maximizing the `1-norm of variance with a greedy algorithm [Kwak, 2008]. The corresponding 2D and supervised versions can be found in [Liu et al., 2010; Li et al., 2010; Pang et al., 2010]. However, the greedy algorithm optimizes the projection directions one by one, which makes it easy to get stuck in a local solution. For this issue, Nie et al. [Nie et al., 2011] exploited an efficient non-greedy optimization algorithm to optimize all projection directions simultaneously for the `1-norm maximization problem; Kwak [Kwak, 2014] extended the non-greedy algorithm to an `p-norm based maximization problem. The corresponding robust 2DPCA with non-greedy algorithm can be found in [Wang et al., 2015].

However, the `1-norm based robust PCA mentioned above usually uses the mean of the data as the optimal mean [Nie et al., 2014]; moreover, it assumes that the data are already centered, i.e., the average of the data is zero. However, this assumption is unreasonable for the following three reasons: (1) It's hard to ensure the zero mean in real-world applications; (2) The outliers in high-dimensional data often make the data mean biased, which degrades the robustness of PCA [He et al., 2011]. (3) It ignores the data mean calculation and incorrectly uses the average of the data as the optimal mean for `1-norm based robust PCA. However, the average of the data is the optimal mean for conventional PCA based on a Euclidean distance.

There are relatively few studies which take these important issues into consideration. To the best of our knowledge, He et al. [He et al., 2011] proposed a robust PCA based on n maximum correntropy criterion and handled non-centered data with an estimation of optimal mean; Nie et al. [Nie et al., 2014] introduced a mean variable and exploited a novel robust PCA objective with optimal mean. Nevertheless, both of these methods integrate the mean calculation into the optimization objective and lead to expensive computation.

In this paper, we develop an equivalent reformulation of the maximization of `2 variances for conventional PCA, such that the optimal projection directions are learned via maximizing the sum of projected difference between each pair of instances instead of the difference between each instance and the mean of the data. Based on this reformulation, we propose a new robust PCA by maximizing the sum of projected differences between each pair of instances based on the `1norm distance. This method automatically avoids calculating the `1-norm based optimal mean and makes the assumption of centered data unnecessary. An efficient non-greedy method is further exploited to maximize the objective with fast convergence in practical application. Intuitively, we also extend the proposed robust PCA to its 2D version for image recognition. It is noteworthy that the proposed algorithms keep linear computation complexity with respect to the number of data points in practical application.

2 Principal Component Analysis Review

Suppose the given data matrix is X = [x1↪ x2↪ · · · ↪ xn] 2 Rd⇥n , where each instance xi is represented by a vector with d-dimensionality; n refers to the number of instances. Conventional PCA learns a transformation to map high dimensional data to low dimensional representations. Specifically, let W = [w1↪w2↪ · · · ↪wm] 2 Rd⇥m be a semi-orthogonal transformation matrix, the idea of traditional PCA is formulated by minimizing the reconstruction error based on `2norm distance in the original high-dimensional space, i.e.

where m is the mean of the data. By setting the derivative of objective function (1) with respect to m to zero, we obtain the optimal mean of the data as m = ¯x = 1n Pni=1 xi. With some evident transformations, optimization problem (1)

can be reformulated as the maximization of covariance in the projected space, i.e.,

Generally, the mean ¯x is supposed to be zero; otherwise, it is subtracted from each instance of optimization problem (2). In such a way, the orthogonal transformation matrix W can be solved by maximizing the following optimization problem

where the instances xi (i = 1↪ 2↪ · · · ↪ n) are centered. PCA has been widely applied in many applications for its efficiency and simplicity However, the high computational complexity and the outlier sensitivity induced by `2-norm make it hard to apply for a large scale data with high dimensionality [Nie et al., 2011]. As a result, robust PCA is proposed by directly substitute `1-norm for `2-norm maximization in optimization problem (3)[Kwak, 2008; Galpin and Hawkins, 1987; Nie et al., 2011], i.e.,

Nevertheless, the existing robust PCA methods and its various neglect the optimal mean calculation and incorrectly utilize ¯x as the optimal mean. However, ¯x is definitely the optimal mean in terms of `2-norm distance, rather than the `1norm used in the objective functions of RPCA. Nie et al. [Nie et al., 2014] consider this issue and propose a new robust PCA by integrating the optimization of optimal mean into the dimension reduction objective; however, it leads to expensive computation.

3 The ProposedMethodology

In this section, we consider a general case that the mean of the data is not zero, and propose a novel robust PCA based on `1-norm distance. This method automatically avoids calculating the optimal mean with `1-norm distance and makes the assumption of centered data unnecessary. For a better representation, we first introduce the following theorem. It reformulates the objective of conventional PCA as maximizing the sum of the projected difference between each pair of instances.

Theorem 1. Let X = [x1↪ x2↪ · · · ↪ xn] 2 Rd⇥n be the data matrix and ¯x = 1n Pni=1 xi be the mean of X. The solution W of conventional PCA which minimizes the reconstruction error based on Euclidean distance, i.e.,

is also the solution of the following optimization problem

Proof. With fixed W satisfying W >W = I, we set the derivative of objective function (5) with respect to variable m to zero. Then we have m = ¯x. We substitute m = ¯x into the objective function (5) and reformulate it as the maximization of covariance in the projected space, i.e.,

We go further to substitute ¯x = 1n Pni=1 xi into the objective function above and achieve the following equation

On the other hand, it is evident to reformulate the objective function of optimization problem (6) through equivalent transformation as

According to Eq. (8) and Eq. (9), the proof is completed.

Based on Theorem 1, we equivalently reformulate the `2norm based PCA as the following optimization problem,

As opposed to the conventional PCA formulated by optimization problem (5) or (2), the alternative formulation (10) estimate the transformation matrix with the calculation of optimal mean avoided automatically. However, considering the sensitivity of `2-norm distance to outliers, in this paper, we propose a novel robust PCA based on `1-norm by solving the following optimization problem

Note that existing robust PCA which directly replaces the `2 norm in optimization problem (3) with `1-norm and incorrectly employs ¯x as the optimal mean of `1-norm based robust PCA. The proposed objective (11) automatically avoids calculating the `1-norm based optimal mean and makes the assumption on centered data unnecessary.

3.1 Optimization procedure

In this section, we use an efficient iterative re-weighted algorithm [Nie et al., 2011] to solve the non-smooth optimization problem (11). This algorithm is first proposed to solve general optimization problem,

Algorithm 1 Non-greedy `1-norm maximization.

Initialize: z(1) 2 C, k = 1.

• 1: while not converge do
• 2: v(k) i = sgn(gi(z(k) )) for each i;
• 5: end while
• 3: z(k+1) = argmaxz2C f(z) + Pi v(k) i gi(z); 4: k = k +1;

Output: z(k) .

where z 2 C is an arbitrary constraint; f and gi are arbitrary functions defined on C for each i. Let vi = sgn(g(z)) with element-wise sign function sgn(·), then the objective function L(z) can be reformulated as

As a result, the general optimization problem (12) can be solved with an non-greedy re-weighted algorithm which is described in Algorithm 1.

We follow Algorithm 1 and exploit an non-greedy method to solve the proposed optimization problem (11). The key step lies in addressing the following optimization problem

where vij = sgn((W (k) )>(xi - xj)) is a vector with mdimensionality. Let R = Pij (xi - xj)v> ij 2 Rd⇥m, then we have Pi↪j v> ijW >(xi - xj) = Tr(W >R). As a result, optimization problem (14) can be rewritten as

We solve optimization problem (15) based on the following Theorem 2.

Theorem 2. Suppose the SVD of R is R = P⇤Q>, where P 2 Rd⇥d , ⇤ 2 Rd⇥m and Q 2 Rm⇥m. The solution of optimization problem (15) is derived as W = P [I; 0]Q>.

Proof. Based on the SVD of R, we have

where = Q>W >P; λkk and kk represent the (k↪ k)-th element of matrices ⇤ and , respectively. Recall that the constraint W >W = I, so we have > = I and kk  1, where I is an m by m identity matrix. Combining with the fact λkk ≥ 0 since λkk is singular value of R, we arrive at the following inequality

where the equality holds when kk = 1 (1  k  m). As a result, the objective function (15) reaches its maximum when = [I↪ 0]. Recall that = Q>W >P, the optimal solution to optimization problem (15) is W = P >Q> = P[I; 0]Q>. The proof is completed.

Algorithm 2 Robust PCA with non-greedy `1-norm maximization

Input: data set ¶xi 2 Rd : i = 1↪ 2↪ · · · ↪ n♦, m. Initialize: W (1) 2 Rd⇥m s.t. (W (1) )>W (1) = I, t = 1.

• 1: while not converge do
• 2: vi↪j = sgn((W (t) )>(xi -xj)) (8i < j);
• 3: R = Pi↪j (xi -xj)v> ij ; 4: Calculate the SVD of R as R = P⇤Q>, then W (t+1) = PQ>;
• 5: t = t + 1;

6: end while

Output: W t 2 Rd⇥m.

In summary, we describe the non-greedy `1-norm maximization algorithm in Algorithm 2 for optimization problem (14). Theoretical analysis in [Nie et al., 2011] ensures that the proposed non-greedy algorithm will convergence and usually obtain a local maximum solution within ten iterations in practical application.

3.2 Complexity discussion

In this section, we analyze the computational complexity of the proposed Algorithm 2. Given F = [f1↪ f2↪ · · · ↪ fn] = W >X 2 Rm⇥n with computational complexity O(ndm), we have vij = sgn(fi - fj) 2 Rm for i↪ j = 1↪ 2↪ · · · ↪ n. It seems like the computational cost of vij (i↪ j = 1↪ 2↪ · · · ↪ n) is O(mn2 ). In fact, it can be avoided with some strategies. Note that the computation of R in Algorithm 2 only depends on m-dimensional vectors vi· := Pj v> ij and v·j := Pi vij due to the following equivalent transformation:

Consequently, the computational cost of R is O(nmd) when vi· and v·i are given. Recall that vij = sgn(fi - fj), thus we have vi· = Pj vij = Pj (fi - fj). As a result, the k-th entry of vectors vi·↪ v·i 2 Rm can be obtained efficiently by sorting the k-th entries of f1↪ f2↪ · · · ↪ fn with computational complexity O(n log(n)). Based on this ranking, the entire computational cost over vi· and v·i (i = 1↪ 2↪ · · · ↪ n) is O(nmlog(n)). As a result, the computational complexity of Algorithm 2 is O(nm(log(n) + d)t). However, since log(n) is usually much smaller than d in practical application, the computational complexity of Algorithm 2 reduces to O(nmdt). Therefore, the proposed robust PCA does not require an additional computational cost in contrast to the state of the art robust PCA in [Kwak, 2008; Nie et al., 2011].

4 Extensions to 2DVersion of Robust PCA

To keep the structural information of two dimensional (2D) image matrix, 2DPCA is proposed to construct an covariance matrix using the original 2D image matrices directly. In this section, we extend the proposed robust PCA to its 2D version and develop a novel robust 2DPCA with the calculation of the optimal mean avoided automatically.

We suppose that the data is denoted by X = [X1↪X2↪ · · · ↪Xn] 2 Rc⇥d⇥n , where each component Xi 2 Rc⇥d (i = 1↪ 2↪ · · · ↪ n) refers to a image matrix; n is the number of data. Let U = [u1↪ u2↪ · · · ↪ um] 2 Rd⇥m be the collection of projection directions uk 2 Rd for k = 1↪ 2↪ · · · ↪m, then conventional 2DPCA learns transformation matrix U by maximizing the following optimization problem with Frobenius norm,

where M = 1n Pi Xi is the Frobenius norm based optimal mean which is usually assumed to be zero in conventional 2DPCA. However, as mentioned in the precious section, the mean of the data is not always zero in practical application.

Following Theorem 1, we rewrite the objective of optimization problem (17) and learn the directions via maximizing the sum of the projected difference between each pair of images, i.e.,

Considering the poor robustness of Euclidean distance, we go further to replace the `F -norm used in optimization problem (1) with `1-norm and propose an alternative formulation of robust 2DPCA asmax

where xik↪ xjk 2 Rd (k = 1↪ 2↪ · · · ↪ c) are the transpose of the k-th row of image matrices Xi and Xj, respectively. We can also solve the optimization problem (20) through Algorithm 1. Thus, the key step lies in addressing the following optimization problem

where vijk = sgn (U >(xik -xjk)) is an m-dimensional vector. Denote S = Pi↪j Pck=1 (xik - xjk)v> ijk , the optimization problem (21) can be rewritten as

Assume the SVD of S is S = ADB>, we have the optimal solution to problem (22) is U = A[I; 0]B> according to Theorem 2. We summarize the robust 2DPCA with nongreedy `1-norm maximization in Algorithm 3. Similar analysis in [Nie et al., 2011; Wang et al., 2015] ensures the convergence of the proposed Algorithm 3. Since d ≫ log(n) and n ≫ m in practical applications, we derive the computational complexity of Algorithm 3 as O(nmcdt) according to similar strategies used for Algorithm 2.

Algorithm 3 Robust 2DPCA with non-greedy `1-norm maximization

Input: data set ¶Xi 2 Rc⇥d : i = 1↪ 2↪ · · · ↪ n♦, m. Initialize: U (1) 2 Rd⇥m s.t. (U (1) )>U (1) = I, t = 1.

• 1: while not converge do
• 3: S = Pij Pck=1 (xik -xjk)v> ijk 2 Rd⇥m; 4: Calculate the SVD of S as S = ADB>, then
• 2: vijk = sgn ((U (t) )>(xik -xjk)) 2 Rm (8i < j↪ 8k);

structed image [Nie et al., 2014], i.e., 1n Pni=1 kxri - xoi k2 where n is the number of images, xri denotes the reconstructed image and xoi is the original image without occlusion.

We compare the reconstruction error of the proposed RPCA-AOM with robust PCA with non-greedy `1-norm maximization (RPCA) [Nie et al., 2011] and optimal mean robust PCA (RPCA-OM) [Nie et al., 2014] in Table 1. Specifically, the reconstruction errors with respect to nine different reduced dimensions from 10 to 50 are reported over 5 benchmark data sets, including the Japanese Female Facial Expression Database (JAFFE) [Dailey et al., 2010], UMIST face data set [Wechsler et al., 2012], the ORL database of faces [Cai et al., 2007], Columbia Object Image Library-20 (COIL-20) data set [Nene et al., 1996] and the USPS handwritten digit database [Liu et al., 2003]. All of the image data sets are downloaded from different web sites.

From Table 1, we can observe that: (1) The proposed RPCA-AOM algorithm performs better than RPCA over all data sets, except for slightly worse performance on a few projected dimensions. Note that RPCA also employs the nongreedy `1-norm maximization algorithm. However, the proposed RPCA-AOM requires no assumption on the zero-mean of the data, while the RPCA method definitely depends on this assumption. (2) The proposed RPCA-AOMmethod performs better than RPCA-OM, except for projected dimension 15 in ORL data set and projected dimension 35 in USPS data set. Note that both RPCA-OM and PRCA-AOM take the incorrect optimal mean into consideration. However, RPCAOM integrates the optimization of the optimal mean into

• 5: t = t + 1;
• 6: end while

Output: U t 2 Rd⇥m.

5 Experimental Analysis

In this section, we conduct thorough experimental evaluations of the proposed Robust PCA and Robust 2DPCA with Avoiding Optimal Mean, abbreviated as RPCA-AOM and 2DRPCA-AOM, respectively.

5.1 Reconstruction error comparison for RPCA

Regarding the experiments on reconstruction with different robust PCA methods, we normalize each initial feature of image into [0↪ 1] and randomly select 20% images to be occluded with randomly place of 1◁4 size for fair comparison. The evaluation metric is defined as the average reconstruction error between an original unoccluded image and the reconthe procedure of dimension reduction, which leads to expensive computational cost. Instead, the proposed PRCA-AOM avoids the calculation of the optimal mean automatically.

5.2 Recognition comparison for robust 2DPCA

Regarding face recognition task, we design a series of experiments to evaluate the performance of different robust 2DPCA methods over the ORL database and UMIST data set. These methods includes robust 2DPCA with greedy algorithm (2DRPCA-G) [Li et al., 2010], robust 2DPCA with non-greedy algorithm (2DRPCA-NG) [Wang et al., 2015] and the proposed 2DRPCA-AOM. Note that both 2DRPCAG and 2DPCA-NG depend on the assumption of zero-mean of the data and incorrectly employ the mean of the data as the optimal mean of `1-norm based 2DPCA.

The ORL data set consists of 400 face images of 40 objects, and each object contains ten images. The UMIST data sets consists of 575 face images of 20 objects, and each object contains a varying number of images ranging from 48 to 19. We randomly select half of the images from each object to form the training set and retain the rest as testing set for both of the data sets. To illustrate the robustness of 2DRPCA-AOM, we corrupt a varying percentage of training images with outliers and recognize testing face images in the reduced space with the nearest neighbor (NN) classifier.

With 0, 20, 40 and 60 percentage of images corrupted in training set, Figure 1 and Figure 2 demonstrate the comparisons of recognition accuracy with respect to different methods over ORL data set and UMIST data set, respectively. It indicates that 2DRPCA-G and 2DRPCA-NG achieve worse performance than the proposed 2DRPCA-AOM, especially when the percentage of corrupted training image becomes larger. As a result, the proposed method shows better robustness to outliers. Note that 2DRPCA-NG and 2DRPCA-AOM are based on non-greedy algorithm. They perform much better than 2DRPCA which is based on the greedy algorithm. This illustrates the efficiency and superiority of non-greedy algorithm for `1-maximization.

6 Conclusion

In this paper, we propose a novel robust PCA to learn the projection directions by maximizing the `1-norm based projected difference between each pair of instances, instead of the difference between each instance and the mean of data. This method automatically avoids calculating the optimal mean based on the `1-norm distance and makes any assumption of centered data unnecessary. To solve the proposed non-smooth objective, a non-greedy algorithm is exploited with fast convergence and low computational cost. In a natural way, we extend the proposed method to its 2D version and study the corresponding robust 2DPCA for image recognition. Extensive experimental results illustrate the effectiveness and superiority of the proposed robust PCA and robust 2DPCA on image reconstruction and recognition.

Acknowledgments

This research was funded by the National Science Foundation of China under Grant Nos. 61502377, 61532015, 91418205 and 61532004, the National Science Foundation (NSF) under grant No. IIS-1251187, the U.S. Army Research Office

(W911NF-13-1-0277), China Postdoctoral Science Foundation under Grant No. 2015M582662 and Ministry of Education Innovation Research Team under Grant No. IRT13035.