advertisement

Proceedings of the Twenty-Fourth International Joint Conference on Artificial Intelligence (IJCAI 2015) Semi-Orthogonal Multilinear PCA with Relaxed Start Qiquan Shi and Haiping Lu Department of Computer Science Hong Kong Baptist University, Hong Kong, China csqqshi@comp.hkbu.edu.hk, haiping@hkbu.edu.hk Abstract 2007]. However, PCA on multi-dimensional data requires reshaping tensors into vectors first. This vectorization often leads to breaking of original data structures, more complex model with lots of parameters, and high computational and memory demands [Lu et al., 2013]. Many researchers address this problem via multilinear extensions of PCA to deal with tensors directly, and there are two main approaches. One approach is based on Tensor-to-Tensor Projection (TTP) that learns low-dimensional tensors from highdimensional tensors. The two-dimensional PCA (2DPCA) [Yang et al., 2004] is probably the first PCA extension to deal with images without vectorization. The generalized low rank approximation of matrices (GLRAM) [Ye, 2005] and the generalized PCA (GPCA) [Ye et al., 2004] further generalize 2DPCA from single-sided projections to two-sided projections via reconstruction error minimization and variance maximization, respectively. Concurrent subspace analysis (CSA) [Xu et al., 2005] and multilinear PCA (MPCA) [Lu et al., 2008] extend GLRAM and GPCA to general higher-order tensors, respectively. Another approach is based on Tensor-to-Vector Projection (TVP) that learns low-dimensional vectors from highdimensional tensors in a successive way. The tensor rankone decomposition (TROD) [Shashua and Levin, 2001] minimizes reconstruction error via (greedy) successive residue calculation. The uncorrelated multilinear PCA (UMPCA) [Lu et al., 2009] maximizes variance with the zero-correlation constraint, following the successive derivation of PCA. However, the number of features that can be extracted by UMPCA is upper-bounded by the lowest mode dimension. For example, for a tensor of size 300×200×3, UMPCA can only extract three features, which have very limited usage. Orthogonality constraint is popular in feature extraction [Hua et al., 2007; Kokiopoulou and Saad, 2007; Gao et al., 2013], tensor decomposition [Kolda, 2001], and low-rank tensor approximation [Edelman et al., 1998; Wang et al., 2015]. PCA also obtains orthogonal projections, and the TTP-based PCA methods produce orthogonal projection vectors in each mode. However, none of the existing TVP-based PCA methods derive orthogonal projections. Our study found that it is indeed ineffective to impose full orthogonality in all the modes for TVP-based PCA, due to low captured variance and limited number of extracted features. In this paper, we present a new TVP-based multilin- Principal component analysis (PCA) is an unsupervised method for learning low-dimensional features with orthogonal projections. Multilinear PCA methods extend PCA to deal with multidimensional data (tensors) directly via tensor-to-tensor projection or tensor-to-vector projection (TVP). However, under the TVP setting, it is difficult to develop an effective multilinear PCA method with the orthogonality constraint. This paper tackles this problem by proposing a novel Semi-Orthogonal Multilinear PCA (SO-MPCA) approach. SO-MPCA learns low-dimensional features directly from tensors via TVP by imposing the orthogonality constraint in only one mode. This formulation results in more captured variance and more learned features than full orthogonality. For better generalization, we further introduce a relaxed start (RS) strategy to get SO-MPCA-RS by fixing the starting projection vectors, which increases the bias and reduces the variance of the learning model. Experiments on both face (2D) and gait (3D) data demonstrate that SO-MPCA-RS outperforms other competing algorithms on the whole, and the relaxed start strategy is also effective for other TVP-based PCA methods. Introduction Principal component analysis (PCA) is a classical unsupervised dimensionality reduction method [Jolliffe, 2002]. It transforms input data into a new feature space of lower dimension via orthogonal projections, while keeping most variance of the original data. PCA is widely used in areas such as data compression [Kusner et al., 2014], computer vision [Ke and Sukthankar, 2004], and pattern recognition [Anaraki and Hughes, 2014; Deng et al., 2014]. Many real-world data are multi-dimensional, in the form of tensors rather than vectors [Kolda and Bader, 2009]. The number of dimensions of a tensor is the order and each dimension is a mode of it. For example, gray images are second-order tensors (matrices) and video sequences are third-order tensors [Lu et al., 2013]. Tensor data are also common in applications such as data center monitoring, social network analysis, and network forensics [Faloutsos et al., 3805 where the pth component of y is obtained from the pth EMP as: ear PCA algorithm, Semi-Orthogonal Multilinear PCA (SOMPCA) with Relaxed Start, or SO-MPCA-RS, to be detailed in Sec. 3. There are two main contributions: T T ) (n) yp = y(p) = X ×1 u(1) · · · ×N u(N = X ×N p p n=1 {up }. (4) • We propose a novel SO-MPCA approach to maximize the captured variance via TVP with orthogonality constraint in only one mode, which is called semiorthogonality according to [Wang et al., 2015]. The semi-orthogonality results in more captured variance and more learned features than full-orthogonality. For the same tensor of size 300×200×3 discussed earlier, SO-MPCA can extract 300 features while fullorthogonal multilinear PCA can only extract three features (similar to UMPCA). SO-MPCA with Relaxed Start This section presents the proposed SO-MPCA-RS by first formulating the SO-MPCA problem, then deriving the solutions with a successive and conditional approach, and finally introducing the relaxed start strategy for better generalization. Formulation of Semi-Orthogonal MPCA We define the SO-MPCA problem with orthogonality constraint in only one mode, i.e., semi-orthogonality [Wang et al., 2015], as follows: The SO-MPCA problem: A set of M tensor data samples {X1 , X2 , · · · , XM } are available for training. Each 2 ×···×IN sample Xm ∈ RI1 ×I can be viewed a point in N N the tensor space RI1 N RI2 · · · RIN , where In is the nmode dimension and denotes the Kronecker product. SO(n) MPCA considers a TVP, which consists of P EMPs {up ∈ In ×1 P R , n = 1, · · · , N }p=1 , that projects the input tensor N N I2 R · · · RIN into a vector subspace RP , i.e., space RI1 • We introduce a Relaxed Start (RS) strategy to get SOMPCA-RS by fixing the starting projection vectors for better generalization [Abu-Mostafa et al., 2012]. This strategy constrains the hypothesis space to a smaller set, leading to increased bias and reduced variance of the learning model. The experimental results in Sec. 4 show that SO-MPCA-RS outperforms other competing PCAbased methods on the whole. In addition, this is a new strategy for tensor-based algorithms and we show its effectiveness for other TVP-based PCA methods. In the following, we cover the necessary background first. (n) P ym = Xm ×N n=1 {up , n = 1, · · · , N }p=1 Background for m = 1, · · · , M . The objective is to find a TVP to maximize the variance of the projected samples in each projection direction, subject to the orthogonality constraint in only one mode, denoted as the ν-mode. The variance is measured by the total scatter Sp defined as: Notations and basic operations: We follow the notations in [Lathauwer et al., 2000] to denote vectors by lowercase boldface letters, e.g., x; matrices by uppercase boldface letters, e.g., X; and tensors by calligraphic letters, e.g., X . We denote their elements with indices in parentheses, and indices by lowercase letters spanning the range from 1 to the uppercase letter of the index, e.g., n = 1, · · · , N . An N thorder tensor A ∈ RI1 ×···×IN is addressed by N indices {in }. Each in addresses the n-mode of A. The n-mode product of an N th-order tensor A by a vector u ∈ RIn , denoted by B = A ×n uT , is a tensor with entries: B(i1 , · · · , in−1 , 1, in+1 , · · · , iN ) = X (5) Sp = M X (ymp − ȳp )2 , (6) m=1 P (n) 1 where ymp = Xm ×N n=1 {up }, and ȳp = M m ymp . In other words, the objective of SO-MPCA is to obtain the P EMPs, with the pth EMP determined as: A(i1 , · · · , iN ) · u(in ). in {u(n) p , n = 1, · · · , N } (1) = arg max M X (ymp − y p )2 , (7) m=1 Tensor-to-vector projection: Elementary multilinear projections (EMPs) are the building blocks of a TVP. We denote an EMP as {u(1) , u(2) , · · · , u(N ) }, consisting of one unit projection vector in each mode, i.e., k u(n) k= 1 for n = 1, · · · , N , where k · k is the Euclidean norm for vectors. It projects a tensor X ∈ RI1 ×I2 ×···×IN to a scalar y through the N unit projection vectors as [Lu et al., 2013]: T T T y = X ×1 u(1) ×2 u(2) · · · ×N u(N ) . T s.t. up(n) u(n) =1 p T u(ν) p =0 for n = 1, · · · , N and for p > 1 and q = 1, · · · , p − 1, (9) (8) where the orthogonality constraint (9) is imposed only in the ν-mode and there is no such constraint for the other modes (n = 1, · · · , N , n 6= ν). The normalization constraint (8) is imposed for all modes. Bound on the number of features: Based on the proof of Corollary 1 in [Lu et al., 2009], we can derive that the number of features P that can be extracted by SO-MPCA is upper-bounded by the ν-mode dimension Iν : P ≤ Iν . Since we can choose any n as ν, we have the upper bound of P as P ≤ maxn In (i.e., the highest mode dimension). Selection of mode ν: Although we are free to choose any mode n as ν to impose the orthogonality constraint (9), it is (2) The TVP of a tensor X to a vector y ∈ RP consists of P (1) (N ) EMPs {up , · · · , up }, p = 1, · · · , P , which can be written (n) (n) P concisely as {up , n = 1, · · · , N }P p=1 or {up }p=1 : (n) P y = X ×N n=1 {up , n = 1, · · · , N }p=1 , u(ν) q (3) 3806 (n) often good to have more features in practice. Thus, in this paper, we choose the mode with the highest dimension as ν: ν = arg max In , captured. Then the total scatter matrix S̃p (n) {ỹmp , m = 1, · · · , M } becomes: (10) n S̃(n) p such that P = maxnP In = Iν . On the other hand, we can also obtain a total of n In features by running SO-MPCA N times with ν = 1, · · · , N . In this paper, we only focus on SO-MPCA with ν determined by (10). Semi-orthogonality vs. full-orthogonality: If we impose the orthogonality constraint (9) in all modes, we can get Full-Orthogonal Multilinear PCA (FO-MPCA). However, our study found that FO-MPCA is not effective primarily due to two reasons: = (n) ¯p where ỹ = 1 M s.t. T u(ν) u(ν) p p = 1 and = 0, q = 1, · · · p − 1. (ν) (ν) = λu(ν) Γ(ν) p S̃p up p , (14) where, p−1 X Γ(ν) p = [ I In − T (ν) u(ν) ], q uq (15) q=1 and IIn is an identity matrix of size In × In . Proof. First, we use Lagrange multipliers to transform the problem (13) to include all the constraints as: Conditional subproblem: In order to obtain the pth EMP (n) {up , n = 1, · · · , N }, we need to determine N vectors. We follow the approach of alternating least squares [Harshman, 1970]. Thus, we can only obtain locally optimal solutions as in many other tensor-based methods. For the pth EMP, (n) the parameters of the n-mode projection vector up are estimated one mode by one mode separately conditioned on the projection vectors in all the other modes. Assuming the pth projection vectors in all but n-mode are given, we project the input tensor samples in these (N − 1) modes to obtain the partial multilinear projections as in [Lu et al., 2013]: Lν T T (ν) (ν) = u(ν) S̃(ν) u(ν) p p up − λ(up p − 1) p−1 X − T µq u(ν) u(ν) p q , (16) q=1 where λ and {µq , q = 1, · · · , p−1} are Lagrange multipliers. (ν) Then we set the partial derivative of Lν with respect to up to zero: ∂Lν (ν) ∂up T ×n−1 u(n−1) p T (13) Theorem 1. The solution to the problem (13) is the (unitlength) eigenvector corresponding to the largest eigenvalue of the following eigenvalue problem: Step p (p = 2, · · · , P ): (n) Determine the pth EMP {up , n = 1, · · · , N } by maximizing Sp with the constraints (8) and (9). ) ×n+1 u(n+1) · · · ×N u(N , p p T u(ν) u(ν) p q We solve this problem by the following theorem: (n) ··· (n) ỹmp . T Step 1 (p = 1): Determine the first EMP {u1 , n = 1, · · · , N } by maximizing S1 with the constraint (8). T m=1 (ν) (ν) u(ν) S̃(ν) p = arg max up p up To solve the SO-MPCA problem, we follow the successive derivation in [Jolliffe, 2002; Lu et al., 2009] to determine EMPs one by one in P steps: = PM (n) Successive Derivation of SO-MPCA Xm ×1 u(1) p (12) For p = 1 (step 1), the solution for u1 , where n = (n) 1, · · · , N , is obtained as the unit eigenvector of S̃1 associated with the largest eigenvalue. For p ≥ 2, we need to deal with the ν-mode and other modes differently. For modes other than ν, the solution for (n) up , where n = 1, · · · , N , n 6= ν, is obtained as the unit (n) eigenvector of S̃p associated with the largest eigenvalue. Constrained optimization for ν-mode and p ≥ 2: When (ν) p ≥ 2, we need to determine up by solving the following constrained optimization problem: • Similar to UMPCA, the number of features that can be extracted by FO-MPCA is upper-bounded by the lowest mode dimension minn In , which can be quite limited. For instance, FO-MPCA can extract only three features for a tensor of size 300×200×3 while SO-MPCA can extract 300 features by choosing ν = 1 for the same tensor. This can be observed in Fig. 1 as well. (n) ỹm p (n) (n) ¯ p(n) )(ỹm ¯ p(n) )T , (ỹm − ỹ − ỹ p p m=1 • Due to the heavy constraints, the variance captured by FO-MPCA is quite low, even lower than UMPCA. In contrast, SO-MPCA can capture more variance than both FO-MPCA and UMPCA. This is illustrated in Fig. 1 in Sec. 4. T M X corresponding to (ν) (ν) = 2S̃(ν) p up − 2λup − Premultiplying (17) by up get (n) where ỹmp ∈ RIn . This conditional subproblem then be(n) comes to determine up that projects the vector samples (n) {ỹmp , m = 1, · · · , M } onto a line to maximize the variance µq u(ν) q = 0. , the third term vanishes and we T T (ν) (ν) 2u(ν) S̃(ν) u(ν) p p up − 2λup p =0 ⇒ 3807 (17) q=1 (ν) T (11) p−1 X T (ν) λ = u(ν) S̃(ν) p p up , (18) Relaxed Start for Better Generalization Algorithm 1 Semi-Orthogonal Multilinear PCA with Relaxed Start (SO-MPCA-RS) When we use SO-MPCA features for classification, we find the performance is limited. Therefore, we further introduce a simple relaxed start (RS) strategy to get SO-MPCA-RS by (n) fixing the first EMP {u1 , n = 1, · · · , N } (the starting projection vectors), without variance maximization. In this pa(n) per, we set this starting EMP u1 (for n = 1, · · · , N ) to the normalized uniform vector 1/ k 1 k for simplicity. This idea is motivated by the theoretical studies in Chapter 4 of [Abu-Mostafa et al., 2012] showing that constraining a learning model could lead to better generalization. By fixing the first EMP as simple vectors, the following EMPs have less freedom due to the imposed semi-orthogonality, which increases the bias and reduces the variance of the learning model. Thus, the SO-MPCA-RS model has a smaller hypothesis set than the SO-MPCA model. The two algorithms differ only in how to determine the first (starting) EMP though the following EMPs will all be different due to their dependency on the first EMP. This relaxed start strategy is not specific to SO-MPCA but generally applicable to any TVP-based subspace learning algorithm. We run controlled experiments in Sec. 4 to show that it can improve the performance of not only SO-MPCA but also TROD and UMPCA. Algorithm 1 summarizes the SO-MPCA-RS algorithm.1 The SO-MPCA algorithm can be obtained from Algorithm 1 by removing line 3, changing p = 2 in line 4 to p = 1 and (ν) setting Γ1 (p = 1) in line 10 to an identity matrix. 1: Input: A set of tensor samples {Xm ∈ RI1 ×···×IN , m = 1, · · · , M }, and the maximum number of iterations K. 2: Set ν = arg maxn In . (n) 3: Set the first EMP: u1 = 1/ k 1 k for n = 1, · · · , N . 4: for p = 2 to P do (n) 5: Initialize up = 1/ k 1 k for n = 1, · · · , N . 6: for k = 1 to K do 7: for n = 1 to N do (n) 8: Calculate the partial multilinear projection {ỹmp } for m = 1, · · · , M according to (11). 9: if n == ν then (ν) (ν) 10: Calculate Γp and S̃p according to (15) and (12), (ν) respectively. Then, set up to the eigenvector of (ν) (ν) Γp S̃p associated with the largest eigenvalue. 11: else (n) (n) 12: Calculate S̃p by (12). Set up to the eigenvector (n) of S̃p associated with the largest eigenvalue. 13: end if 14: end for 15: end for 16: end for (n) 17: Output The TVP {up , n = 1, · · · , N }P p=1 . which indicates that λ is exactly the criterion to be maximized, with the orthogonality constraint. Next, a set of (p − 1) equations are obtained by premulti(ν) T plying (17) by uq T (ν) (ν) 2u(ν) S̃p up q − , q = 1, · · · , p − 1, respectively, T 2λuq(ν) up(ν) − p−1 X Experiments T µs uq(ν) us(ν) = 0. (19) This section evaluates the proposed methods on both secondorder and third-order tensor data in terms of recognition rate, the number of extracted features, captured variance, and convergence. In addition, we also study the effectiveness of the relaxed start strategy on other TVP-based PCA algorithms. Data:2 For second-order tensors, we use the same subset of the FERET database [Phillips et al., 2000] as in [Lu et al., 2009], with 721 face images from 70 subjects. Each face image is normalized to 80 × 60 graylevel pixels. For thirdorder tensors, we use a subset of the USF HumanID “Gait Challenge” database [Sarkar et al., 2005]. We use the same gallery set (731 samples from 71 subjects) and probe A (727 samples from 71 subjects) as in [Lu et al., 2009], and we also test probe B (423 samples from 41 subjects) and probe C (420 samples from 41 subjects). Each gait sample is a (binary) silhouette sequence of size of 32×22×10. Experiment setup: In face recognition experiments, we randomly select L = 1, 2, 3, 4, 5, 6, 7 samples from each subject as the training data and use the rest for testing. We repeat such random splits (repetitions) ten times and report the mean correct recognition rates. In gait recognition experiments, we follow the standard setting and use the gallery set as the training data and probes A, B, and C as the test data (so there is s=1 The second term vanishes and the summand in the third term is non-zero only for s = q. Thus, we get T T (ν) (ν) (ν) 2u(ν) S̃(ν) S̃(ν) q p up − µq = 0 ⇒ µq = 2uq p up . (20) Substituting (20) into (17), we get (ν) 2S̃(ν) p up − 2λu(ν) p − p−1 X T (ν) (ν) u(ν) S̃(ν) q · 2uq p up = 0 q=1 ⇒ λu(ν) p = (ν) S̃(ν) p up − p−1 X T (ν) (ν) u(ν) S̃(ν) (21) q uq p up q=1 ⇒ λu(ν) p = [ I In − p−1 X T (ν) (ν) u(ν) ] S̃(ν) q uq p up . (22) q=1 Using the definition in (15), (22) can be rewritten as: (ν) (ν) Γ(ν) = λu(ν) p S̃p up p . (23) Since λ is the criterion to be maximized, this maximization is (ν) (ν) achieved by setting up to the (unit) eigenvector of Γ(ν) p S̃p associated with its corresponding largest eigenvalue . 1 Matlab code is available at: http://www.comp.hkbu.edu.hk/ ∼haiping/codedata.html 2 Both face and gait data are downloaded from: http://www.dsp. utoronto.ca/∼haiping/MSL.html 3808 Table 1: Face recognition rates in percentage (mean ± std) by the nearest neighbor classifier on the FERET subset. The top two results are highlighted with bold fonts and ‘-’ indicates that no enough features can be extracted. L 1 2 3 4 5 6 7 P 1 5 10 20 50 80 1 5 10 20 50 80 1 5 10 20 50 80 1 5 10 20 50 80 1 5 10 20 50 80 1 5 10 20 50 80 1 5 10 20 50 80 PCA 2.60±0.66 15.12±1.31 22.69±2.21 27.62±2.57 31.38±2.58 – 2.69±0.73 20.17±1.25 32.03±2.49 39.07±1.87 43.86±2.53 45.28±2.39 2.72±0.45 23.89±1.64 37.20±1.91 46.05±2.11 51.35±2.53 52.66±2.67 2.68±0.85 25.26±1.77 41.54±2.02 49.34±1.69 56.85±2.09 58.16±2.46 2.91±0.91 28.95±2.07 47.06±1.54 55.66±1.94 63.91±1.71 64.61±1.67 2.86±0.89 30.30±2.17 48.97±2.96 58.57±2.84 66.88±2.31 68.31±2.33 2.68±1.19 29.52±1.38 51.21±2.11 59.57±2.64 68.10±2.21 69.70±2.84 CSA 3.87±1.02 11.90±1.25 21.35±2.76 26.16±2.38 31.37±1.92 31.95±1.84 3.36±0.54 15.15±1.03 29.45±1.95 37.07±2.27 44.61±2.34 45.82±2.76 4.07±0.80 16.58±0.95 36.05±1.50 43.87±1.98 51.60±2.50 52.84±2.70 3.92±1.04 18.93±1.29 40.39±2.36 49.39±2.35 57.51±2.98 59.05±2.65 4.37±1.06 20.75±1.77 45.77±2.17 56.01±2.19 64.58±2.13 65.58±1.92 3.89±0.67 21.89±2.04 49.14±2.57 58.97±2.60 67.84±2.48 69.44±2.49 4.55±1.07 22.51±1.29 49.39±3.18 60.91±2.75 69.35±1.89 70.39±1.76 MPCA 2.52±0.76 16.65±1.82 22.24±1.94 27.16±1.47 31.29±1.71 32.17±2.09 2.63±0.55 21.53±0.90 28.04±1.69 38.86±2.12 44.54±2.74 46.02±2.67 2.25±0.44 25.95±1.26 34.91±2.40 45.48±2.46 52.00±2.70 53.31±2.59 2.22±0.62 28.71±1.91 39.43±2.05 50.18±3.03 57.48±2.72 58.91±2.50 2.72±0.88 32.99±2.47 43.29±3.07 56.79±2.14 64.37±2.27 65.85±2.02 2.49±1.01 33.42±2.52 45.65±2.85 59.73±2.96 67.84±2.63 69.70±2.35 2.12±1.07 34.72±2.86 46.19±2.43 61.69±2.57 69.26±2.22 70.65±1.97 TROD 2.70±0.44 15.88±1.20 21.52±2.83 26.30±2.49 29.63±2.21 31.14±2.42 2.65±0.77 21.34±1.56 30.05±1.70 36.87±2.37 42.67±2.16 44.46±2.53 2.95±0.61 24.48±1.79 34.83±2.96 43.37±2.51 48.92±2.69 51.17±2.65 3.11±0.83 27.37±2.34 38.82±3.91 47.57±2.70 54.58±2.56 57.30±2.46 2.59±0.64 31.75±2.79 43.80±3.51 54.47±1.66 61.54±2.75 64.02±2.40 2.86±1.01 33.59±2.63 45.88±2.97 57.11±3.22 64.05±2.95 66.78±2.94 2.81±0.82 32.03±2.46 46.10±2.60 57.58±2.75 65.54±2.79 67.97±2.45 UMPCA 5.98±2.65 23.23±4.49 31.83±5.17 35.94±5.65 36.14±5.73 – 7.28±2.44 26.90±5.23 40.17±6.76 44.51±6.54 45.47±6.65 – 7.42±1.17 33.86±3.65 49.39±2.83 55.83±3.32 56.42±3.11 – 7.96±2.15 37.66±4.88 55.10±4.55 62.40±4.49 63.13±4.15 – 7.41±2.14 40.78±5.82 60.49±6.37 66.90±6.23 67.71±6.31 – 9.07±0.83 42.52±4.99 63.16±5.32 70.73±5.39 72.09±5.18 – 11.39±1.85 44.98±5.32 65.67±5.82 73.16±4.28 74.11±4.52 – no random splits/repetitions), and report the rank 1 and rank 5 recognition rates [Sarkar et al., 2005]. Algorithms and their settings: We first evaluate SOMPCA and SO-MPCA-RS against five existing PCA-based methods: PCA [Jolliffe, 2002], CSA [Xu et al., 2005], MPCA [Lu et al., 2008], TROD [Shashua and Levin, 2001], and UMPCA [Lu et al., 2009].3 CSA and MPCA produce tensorial features so they need to be vectorized. MPCA uses the full projection. For TROD and UMPCA, we use the uniform initialization [Lu et al., 2009]. For SO-MPCA and SOMPCA-RS, we set the selected mode ν = 1 for the maximum number of features. For iterative algorithms, we set the number of iterations to 20. All features are sorted according to the scatters (captured variance) in descending order for classification. We use the Nearest Neighbor Classifier with the Euclidean distance measure to classify the top P features. We test up to P = 80 features in face recognition and up to P = 32 features in gait recognition. The performance of FOMPCA is much worse than SO-MPCA so it is not included in the comparisons (except variance study) to save space. Face recognition results: Table 1 shows the face recognition results for P = 1, 5, 10, 20, 50, 80 and L = 1, 2, 3, 4, 5, 6, 7, including both the mean and the standard deviation (std) over ten repetitions. We highlight the top two re- SO-MPCA 2.73 ±0.69 20.06±2.34 28.77±2.72 31.94±2.95 32.33±2.78 32.26±2.71 2.69±0.46 24.34±1.59 36.94±2.58 41.70±2.48 42.24±2.39 42.20±2.39 2.56±0.56 28.30±2.05 43.31±1.51 49.49±2.16 49.80±2.29 49.77±2.36 3.13±0.90 29.61±2.16 47.10±2.88 53.85±2.89 54.56±3.14 54.47±3.17 3.07±0.60 34.20±2.67 53.45±2.75 61.40±2.43 62.26±2.88 62.18±2.83 2.56±0.74 35.18±1.32 56.15±2.36 64.72±3.11 65.32±2.86 65.08±2.73 2.38±1.33 35.58±2.04 56.84±2.07 65.37±2.24 65.37±2.02 65.37±2.02 SO-MPCA-RS 6.85±1.44 27.27±2.36 36.34±3.56 40.32±3.40 40.48±3.09 40.41±3.09 7.97±1.10 33.82±1.45 46.63±1.95 52.19±2.11 52.22±1.73 52.19±1.80 7.55±1.17 36.93±1.76 54.38±3.09 61.25±2.85 61.08±2.83 60.98±2.83 8.34±1.37 40.25±1.52 59.30±2.49 66.60±3.07 67.03±2.86 66.89±2.90 8.38±0.97 42.35±3.04 63.23±3.37 69.97±2.55 70.70±2.47 70.70±2.45 8.97±0.97 43.32±1.82 65.75±2.76 74.39±2.79 74.75±2.60 74.72±2.56 10.91±1.37 45.89±2.34 67.53±2.34 74.89±1.97 75.24±2.12 75.19±2.18 TROD-RS 2.63±0.82 16.68±1.50 21.86±3.03 26.51±1.97 29.80±1.48 31.21±1.94 2.81±0.79 21.62±1.93 30.83±1.93 37.64±2.19 42.99±2.28 44.78±2.48 2.74±0.81 26.14±2.17 34.64±3.01 42.99±2.65 49.47±2.76 51.45±2.63 3.08±0.66 28.25±2.82 39.25±2.95 47.73±2.96 54.42±2.89 57.48±2.46 2.96±0.75 33.45±1.41 44.69±2.71 54.64±1.96 61.51±1.92 64.02±2.03 2.56±0.99 35.18±2.52 47.24±2.55 57.81±2.18 65.12±2.92 68.01±2.90 1.99±1.19 35.02±1.82 48.44±3.20 58.31±2.62 66.06±3.10 68.14±2.96 UMPCA-RS 6.04±2.00 24.78±4.76 33.16±5.36 36.65±5.46 37.05±5.46 – 6.82±1.56 29.19±3.55 44.01±3.22 48.67±3.57 49.07±3.61 – 7.34±1.08 32.68±4.67 50.55±4.81 55.89±4.16 56.38±4.62 – 7.03±1.68 38.84±2.97 57.30±4.86 64.08±4.94 64.85±5.01 – 7.20±1.35 41.67±1.95 62.88±2.59 70.19±3.25 70.81±3.08 – 7.21±1.35 44.65±4.14 66.51±3.10 74.52±3.10 75.12±2.79 – 8.57±1.68 45.93±1.93 66.93±2.01 75.80±2.24 76.88±1.72 – sults in each row in bold fonts for easy comparison. Only SOMPCA-RS consistently achieves the top 2 results in all cases. Compared with existing methods (PCA, CSA, MPCA, TROD and UMPCA), SO-MPCA-RS outperforms the best performing existing algorithm (UMPCA) by 3.79% on average. Furthermore, for larger L = 5, 6, 7, SO-MPCA-RS outperforms the other five methods at least by 2.26% on average. For smaller L = 1, 2, 3, SO-MPCA-RS achieves a greater improvement of at least 5.28% over existing methods, indicating that SO-MPCA-RS is more superior in dealing with the small sample size (overfitting) problem. Gait recognition results: Similarly, the gait recognition results are reported in Table 2 with the top two results highlighted. Again, only SO-MPCA-RS consistently achieves the top 2 results in all cases. In rank 1 rate, the best performing existing algorithm is PCA, which outperforms SO-MPCARS by 3.73% on average. While in rank 5 rate, SO-MPCARS outperforms the best performing existing algorithm (still PCA) by 6.76% on average. Number of features: In the tables, we use ‘-’ to indicate that there are not enough features. For PCA, there are at most 69 features for face data when L = 1 since there are only 70 samples for training. UMPCA can only extract 60 or 10 features for face and gait data, respectively. In contrast, SOMPCA and SO-MPCA-RS (with ν = 1) can learn 80 features for face data and 32 features for gait data. 3 For second-order tensors, CSA and MPCA are equivalent to GLRAM [Ye, 2005] and GPCA [Ye et al., 2004], respectively. 3809 Table 2: Rank 1 and rank 5 gait recognition rates in percentage (mean ± std ) by the nearest neighbor classifier on the USF subset. The top two results are highlighted with bold fonts and ‘-’ indicates that no enough features can be extracted. P robe A 1 B C A 5 B C PCA 30.99 52.11 67.61 71.83 26.83 48.78 65.85 68.29 12.20 29.27 34.15 46.34 57.75 80.28 87.32 87.32 48.78 73.17 78.05 78.05 51.22 53.66 65.85 65.85 PCA UMPCA FO−MPCA SO−MPCA SO−MPCA−RS 8 10 Variance Captured P 5 10 20 32 5 10 20 32 5 10 20 32 5 10 20 32 5 10 20 32 5 10 20 32 6 10 4 10 0 20 40 60 CSA 22.54 43.66 57.75 59.15 17.07 39.02 53.66 60.98 9.76 14.63 31.71 34.15 56.34 74.65 80.28 81.69 48.78 70.73 75.61 73.17 34.15 46.34 56.10 60.98 MPCA 32.39 49.30 60.56 61.97 24.39 46.34 58.54 58.54 14.63 19.51 29.27 29.27 66.20 77.46 81.69 83.10 53.66 70.73 78.05 78.05 41.46 43.90 60.98 58.54 TROD 28.17 42.25 53.52 60.56 19.51 39.02 53.66 65.85 4.88 14.63 21.95 31.71 54.93 77.46 76.06 78.87 48.78 60.98 75.61 80.49 36.59 48.78 48.78 56.10 UMPCA 39.44 57.75 – – 26.83 46.34 – – 24.39 29.27 – – 73.24 83.10 – – 58.54 65.85 – – 41.46 43.90 – – 6 10 4 10 80 0 20 Feature Index 40 60 80 Feature Index (a) Sorted variance (b) Unsorted variance Figure 1: The captured variance on face data with L = 1. 8 7 x 10 2.6 1.231 1.2308 1.2306 1.2304 1.2302 p=2 2 4 6 8 10 12 14 16 18 20 Number of Iterations (a) p = 2 Variance Captured Variance Captured 1.2312 x 10 2.59 2.58 2.57 2.56 p=5 2.55 2.54 2 4 6 SO-MPCA 30.99 49.30 54.93 54.93 29.27 53.66 58.54 60.98 12.20 29.27 31.71 31.71 67.61 77.46 80.28 80.28 63.41 70.73 70.73 70.73 43.90 48.78 48.78 48.78 SO-MPCA-RS 40.85 59.15 67.61 69.01 41.46 63.41 65.85 68.29 14.63 29.27 34.15 39.02 84.51 88.73 92.96 92.96 68.29 75.61 80.49 80.49 56.10 70.73 78.05 78.05 TROD-RS 18.31 33.80 53.52 64.79 19.51 36.59 43.90 63.41 7.32 19.51 24.39 39.02 42.25 59.15 77.46 81.69 53.66 75.61 78.05 80.49 29.27 43.90 46.34 56.10 UMPCA-RS 39.44 63.38 – – 39.02 51.22 – – 14.63 21.95 – – 73.24 85.92 – – 58.54 68.29 – – 36.59 56.10 – – The variance captured by the first (fixed) EMP of SO-MPCARS is much less than other EMPs, which is not surprising since the variance is not maximized. Convergence: We demonstrate the convergence of SOMPCA-RS in Fig. 2 for face data with L = 1. We can see that SO-MPCA-RS converges in just a few iterations. SO-MPCA has a similar convergence rate. Effectiveness of relaxed start: To evaluate the proposed relaxed start strategy, we apply it to two other TVPbased methods, TROD and UMPCA, getting TROD-RS and UMPCA-RS, respectively. We summarize their performance in the last two columns of Tables 1 and 2 for face and gait recognition experiments, respectively. Both tables show that relaxed start can help both TROD and UMPCA to achieve better recognition rates. From Table 1, TROD-RS improves over TROD by 0.32% and UMPCARS improves over UMPCA by 2.15% on average. From Table 2, TROD-RS achieves 3.03 % improvement over TROD and UMPCA-RS achieves 1.07% improvement over UMPCA for rank 1 rate. For rank 5 rate, TROD-RS improves 0.94% over TROD, and UMPCA-RS improves 5.28% over UMPCA. The relaxed start is most effective for our SO-MPCA. SO-MPCARS has an improvement of 9.97% on face data over SOMPCA. On gait data, SO-MPCA-RS outperforms SO-MPCA by 9.56% in rank 1 rate and 17.26% in rank 5 rate on average. In addition, SO-MPCA-RS has better face recognition performance than TROD-RS and UMPCA-RS with an improvement of 8.08% and 1.64%, respectively. On gait data, SOMPCA-RS improves the rank 1 recognition rate by 3.03% over TROD-RS and 13.25% over UMPCA-RS on average, and SO-MPCA-RS improves the rank 5 recognition rate by 11.07% and 13.73% over TROD-RS and UMPCA-RS. These controlled experiments show the effectiveness of relaxed start on SO-MPCA and other TVP-based multilinear PCA methods (UMPCA and TROD). One possible explanation is that RS increases the bias and reduces the variance of PCA UMPCA FO−MPCA SO−MPCA SO−MPCA−RS 8 10 Variance Captured Rank 8 10 12 14 16 18 20 Number of Iterations (b) p = 5 Figure 2: Illustration of the SO-MPCA-RS algorithm’s convergence performance on the face data with L = 1. Feature variance: We illustrate the variance captured by PCA, UMPCA, FO-MPCA, SO-MPCA, and SO-MPCA-RS in Fig. 1 for face data with L = 1 (not all methods are shown for clarity). Figure 1(a) shows the sorted variance. It is clear that semi-orthogonality captures more variance than full-orthogonality, as we discussed in Sec. 3.1. Moreover, both SO-MPCA and SO-MPCA-RS can capture more variance than UMPCA, but less than PCA (and also CSA, MPCA, and TROD, which are not shown). Though capturing less variance, SO-MPCA-RS achieves better overall classification performance than other PCA-based methods, with results consistently in the top two in all experiments. We also show the unsorted captured variance in Fig. 1(b). 3810 [Jolliffe, 2002] I. T. Jolliffe. Principal Component Analysis. Springer Series in Statistics, second edition, 2002. [Ke and Sukthankar, 2004] Y. Ke and R. Sukthankar. PCA-SIFT: A more distinctive representation for local image descriptors. In Proc. IEEE Int. Conf. on Computer Vision and Pattern Recognition, volume II, pages 506–513, 2004. [Kokiopoulou and Saad, 2007] E. Kokiopoulou and Y. Saad. Orthogonal neighborhood preserving projections: A projectionbased dimensionality reduction technique. IEEE Trans. Pattern. Anal. Mach. Intell., 29(12):2143–2156, 2007. [Kolda and Bader, 2009] T. G. Kolda and B. W. Bader. Tensor decompositions and applications. SIAM Rev., 51(3):455–500, 2009. [Kolda, 2001] T.G. Kolda. Orthogonal tensor decompositions. SIAM J. Matrix Anal. Appl., 23(1):243–255, 2001. [Kusner et al., 2014] M. Kusner, S. Tyree, K.Q. Weinberger, and K. Agrawal. Stochastic neighbor compression. In Proc. 31st Int. Conf. Machine Learning, pages 622–630, 2014. [Lathauwer et al., 2000] L. De Lathauwer, B. De Moor, and J. Vandewalle. On the best rank-1 and rank-(R1 , R2 , ..., RN ) approximation of higher-order tensors. SIAM J. Matrix Anal. Appl., 21(4):1324–1342, 2000. [Lu et al., 2008] H. Lu, K. N. Plataniotis, and A. N. Venetsanopoulos. MPCA: Multilinear principal component analysis of tensor objects. IEEE Trans. Neural Networks, 19(1):18–39, 2008. [Lu et al., 2009] H. Lu, K. N. Plataniotis, and A. N. Venetsanopoulos. Uncorrelated multilinear principal component analysis for unsupervised multilinear subspace learning. IEEE Trans. Neural Networks, 20(11):1820–1836, 2009. [Lu et al., 2013] H. Lu, K. N. Plataniotis, and A. N. Venetsanopoulos. Multilinear Subspace Learning: Dimensionality Reduction of Multidimensional Data. CRC Press, 2013. [Phillips et al., 2000] P. J. Phillips, H. Moon, S. A. Rizvi, and P. Rauss. The FERET evaluation method for face recognition algorithms. IEEE Trans. Pattern. Anal. Mach. Intell., 22(10):1090– 1104, 2000. [Sarkar et al., 2005] S. Sarkar, P.J. Phillips, Z. Liu, I. R. Vega, P. Grother, and K. W.Bowyer. The human ID gait challenge problem: Data sets, performance, and analysis. IEEE Trans. Pattern. Anal. Mach. Intell., 27(2):162–177, 2005. [Shashua and Levin, 2001] A. Shashua and A. Levin. Linear image coding for regression and classification using the tensor-rank principle. In Proc. IEEE Int. Conf. on Computer Vision and Pattern Recognition, volume I, pages 42–49, 2001. [Wang et al., 2015] L. Wang, M.T. Chu, and B. Yu Wang. SIAM J. Matrix Anal. Appl., 36(1):1–19, 2015. [Xu et al., 2005] D. Xu, S. Yan, L. Zhang, H.-J. Zhang, Z. Liu, and H.-Y. Shum;. Concurrent subspaces analysis. In Proc. IEEE Int. Conf. on Computer Vision and Pattern Recognition, volume II, pages 203–208, 2005. [Yang et al., 2004] J. Yang, D. Zhang, A. F Frangi, and J. Yang. Two-dimensional PCA: a new approach to appearance-based face representation and recognition. IEEE Trans. Pattern. Anal. Mach. Intell., 26(1):131–137, 2004. [Ye et al., 2004] J. Ye, R. Janardan, and Q. Li. GPCA: An efficient dimension reduction scheme for image compression and retrieval. In Proc. ACM SIGKDD Int. Conf. on Knowledge Discovery and Data Mining, pages 354–363, 2004. [Ye, 2005] J. Ye. Generalized low rank approximations of matrices. Machine Learning, 61(1-3):167–191, 2005. the learning model, while further investigation is needed. Conclusion This paper proposes a novel multilinear PCA algorithm under the TVP setting, named as semi-orthogonal multilinear PCA with relaxed start (SO-MPCA-RS). The proposed SO-MPCA approach learns features directly from tensors via TVP to maximize the captured variance with the orthogonality constraint imposed in only one mode. This semi-orthogonality can capture more variance and learn more features than fullorthogonality. Furthermore, the introduced relaxed start strategy can achieve better generalization by fixing the starting projection vectors to uniform vectors to increase the bias and reduce the variance of the learning model. Experiments on face (2D data) and gait (3D data) recognition show that SOMPCA-RS achieves the best overall performance compared with competing algorithms. In addition, relaxed start is also effective for other TVP-based PCA methods. In this paper, we studied semi-orthogonality in only one mode. A possible future work is to learn SO-MPCA-RS features from each mode separately and then do a feature/scorelevel fusion. Acknowledgments This research was supported by Research Grants Council of the Hong Kong Special Administrative Region (Grant 22200014 and the Hong Kong PhD Fellowship Scheme). References [Abu-Mostafa et al., 2012] Y. S. Abu-Mostafa, M. Magdon-Ismail, and H.-T. Lin. Learning from Data, chapter 4, pages 119–166. AMLBook, 2012. [Anaraki and Hughes, 2014] F.P. Anaraki and S. Hughes. Memory and computation efficient PCA via very sparse random projections. In Proc. 31st Int. Conf. on Machine Learning, pages 1341– 1349, 2014. [Deng et al., 2014] W. Deng, J. Hu, J. Lu, and J. Guo. Transforminvariant PCA: A unified approach to fully automatic face alignment, representation, and recognition. IEEE Trans. Pattern. Anal. Mach. Intell., 36(6):1275–1284, 2014. [Edelman et al., 1998] A. Edelman, T.A. Arias, and S.T. Smith. The geometry of algorithms with orthogonality constraints. SIAM J. Matrix Anal. Appl., 20(2):303–353, 1998. [Faloutsos et al., 2007] C. Faloutsos, T. G. Kolda, and J. Sun. Mining large time-evolving data using matrix and tensor tools. In Int. Conf. on Data Mining Tutorial, 2007. [Gao et al., 2013] Q. Gao, J. Ma, H. Zhang, X. Gao, and Y. Liu. Stable orthogonal local discriminant embedding for linear dimensionality reduction. IEEE Trans. Image Processing, 22(7):2521– 2531, 2013. [Harshman, 1970] R. A. Harshman. Foundations of the parafac procedure: Models and conditions for an “explanatory” multi-modal factor analysis. UCLA Working Papers in Phonetics, 16:1–84, 1970. [Hua et al., 2007] G. Hua, P.A. Viola, and S.M. Drucker. Face recognition using discriminatively trained orthogonal rank one tensor projections. In Proc. IEEE Int. Conf. on Computer Vision and Pattern Recognition, pages 1–8, 2007. 3811