advertisement

Principal Component Analysis on non-Gaussian Dependent Data Fang Han Johns Hopkins University, 615 N.Wolfe Street, Baltimore, MD 21205 USA Han Liu Princeton University, 98 Charlton Street, Princeton, NJ 08544 USA Abstract In this paper, we analyze the performance of a semiparametric principal component analysis named Copula Component Analysis (COCA) (Han & Liu, 2012) when the data are dependent. The semiparametric model assumes that, after unspecified marginally monotone transformations, the distributions are multivariate Gaussian. We study the scenario where the observations are drawn from non-i.i.d. processes (m-dependency or a more general φ-mixing case). We show that COCA can allow weak dependence. In particular, we provide the generalization bounds of convergence for both support recovery and parameter estimation of COCA for the dependent data. We provide explicit sufficient conditions on the degree of dependence, under which the parametric rate can be maintained. To our knowledge, this is the first work analyzing the theoretical performance of PCA for the dependent data in high dimensional settings. Our results strictly generalize the analysis in Han & Liu (2012) and the techniques we used have the separate interest for analyzing a variety of other multivariate statistical methods. 1. Introduction In much of studies on Principal Component Analysis (PCA) it is assumed that the n observations x1 , . . . , xn of a random vector X ∈ Rd are independent. Moreover, in high dimensions, it is commonly assumed that X follows a multivariate Gaussian or sub-Gaussian distribution such that the estimators are consistent. Proceedings of the 30 th International Conference on Machine Learning, Atlanta, Georgia, USA, 2013. JMLR: W&CP volume 28. Copyright 2013 by the author(s). fhan@jhsph.edu hanliu@princeton.edu In this paper we focus on a semiparametric method built on the nonparanormal model. A continuous random vector X := (X1 , . . . , Xd )T follows a nonparanormal distribution if there exists a set of univariate monotone functions f := {fj }dj=1 such that f (X) := (f1 (X1 ), . . . , fd (Xd ))T follows a Gaussian distribution. In this paper we show that the proposed method can loosen both the data independence and Gaussian/sub-Gaussian assumptions. Let Σ be the covariance matrix of X. PCA aims at recovering the top m leading vectors u1 , . . . , um of Σ. The usual procedures are to estimate the top m leadb1 , . . . , u b m of the sample covariance ing eigenvectors u matrix S. However, there are two drawbacks: (i) when d > n, Johnstone & Lu (2009) show that PCA is inconb 1 be the leading sistent. More specifically, let u1 and u eigenvectors of Σ and S. For two vectors v1 , v2 ∈ Rd , we denote the angle between v1 and v2 by ∠(v1 , v2 ). b 1 ) does not Johnstone & Lu (2009) prove that ∠(u1 , u converge to 0. (ii) The performances of the estimators rely on independence of the n observations. Their performance is unknown if dependence is present among x1 , . . . , xn . A remedy for the inconsistency problem when d > n is to assume that u1 = (u11 , . . . , u1d )T is sparse, i.e., card(supp(u1 )) := card({j : u1j 6= 0}) = s < n, where card(·) represents the cardinality of a given set. Different sparse PCA procedures have been developed to exploit the sparsity structure: greedy algorithms (d’Aspremont et al., 2008), lasso-type methods including SCoTLASS (Jolliffe et al., 2003), SPCA (Zou et al., 2006) and sPCA-rSVD (Shen & Huang, 2008), a number of power methods (Journée et al., 2010; Yuan & Zhang, 2011; Ma, 2011), the biconvex algorithm PMD (Witten et al., 2009) and the semidefinite relaxation DSPCA (d’Aspremont et al., 2004). For data dependence, in low dimensions, Skinner et al. (1986) study the behavior of PCA when the selection PCA for the Dependent Data of a sample of observations depends on a vector of latent covariates as, for example, in stratified sampling. Their analysis is based on the normality assumption and that the knowledge of the survey design is known. There are even fewer literatures in high dimensions for the dependent data. Loh & Wainwright (2011) study the high dimensional regression for Gaussian random vectors following a stationary vector regressive process. Very recently, Fan et al. (2012) analyze the penalized least square estimators, taking a weakly dependence structure, called α-mixing, of the noisy term into consideration. There are several drawbacks of PCA (or sparse PCA): (i) It is not scale-invariant, i.e., changing the measurement scale of variables makes the estimates different; (ii) Most estimating procedures require the data to be either Gaussian or sub-Gaussian so that the sample covariance matrix S converges to Σ in a fast rate; (iii) It cannot handle dependent data. Compared with PCA and sparse PCA, Han & Liu (2012) exploit a nonparametric Kendall’s tau based regularization procedure, named Copula Component Analysis (COCA), for parameter estimation. They show that COCA is scale-invariant, able to deal with continuous data with arbitrary margins. In this paper, we further generalize their results, showing that COCA can allow weak data dependence. In particular, we provide the generalization bounds of convergence for both support recovery and parameter estimation for dependent data using our method. We provide explicit sufficient conditions on the the degree of dependence, under which the same parametric rate can be achieved. To our knowledge, this is the first work of analyzing the theoretical performance of PCA for dependent data in high dimensional settings. The rest of the paper is organized as follows. In the next section, we briefly review the nonparanormal family (Liu et al., 2009; 2012) and the data dependence structure. In Section 3, we introduce the models and rank-based estimators proposed by Han & Liu (2012). We provide our main theoretical analysis of the rankbased estimators for the dependent data in Section 4. In Section 5, we employ the on synthetic data to illustrate the robustness of COCA to data dependence. 2. Background We start with notations: Let M = [Mjk ] ∈ Rd×d and v = (v1 , ..., vd )T ∈ Rd . Let v’s subvector with entries indexed by I be denoted by vI . Let M’s submatrix with rows indexed by I and columns indexed by J be denoted by MIJ . Let MI∗ and M∗J be the submatrix of M with rows in I, and the submatrix of M with columns in J. For 0 < q < ∞, we define the `q and `∞ vector norms as ||v||q := X d q |vi | 1/q and ||v||∞ := max |vi |, 1≤i≤d i=1 and we define ||v||0 := card(supp(v)) · ||v||2 . We define the matrix `max norm as the elementwise maximum value: ||M||max := max{|Mij |}. Let Λj (M) be the j-th largest eigenvalue of M. In particular, Λmin (M) := Λd (M) and Λmax (M) := Λ1 (M) are the smallest and largest eigenvalues of M. The vectorized matrix of M, denoted by vec(M), is defined as: vec(M) := (MT∗1 , . . . , MT∗d )T . Let Sd−1 := {v ∈ Rd : ||v||2 = 1} be the d-dimensional `2 sphere. For any two vectors a, b ∈ Rd and any two squared matrices A, B ∈ Rd×d , denote the inner product of a and b, A and B by ha, bi := aT b and hA, Bi := Tr(AT B). d The sign = denotes that the two sides of the equality have the same distributions. For a sequence of random vector {Xt }∞ t=−∞ and two integers L < U , we denote XLU := {Xt }U t=L . The notation P(·) represents the probability if it is a set inside the brackets and the law (distribution) if it is a random vector inside the brackets. 2.1. The Nonparanormal We first provide the definition of the nonparanormal following Liu et al. (2012). Definition 2.1 (The nonparanormal). Let f = {fj }dj=1 be a set of monotone univariate functions. We say that a d dimensional random variable X = (X1 , . . . , Xd )T follows a nonparanormal distribution N P Nd (Σ, f ), if f (X):=(f1 (X1 ), . . . , fd (Xd ))T ∼Nd (0, Σ), diag(Σ)=1. We call Σ the latent correlation matrix. We next proceed to the invariance property of the rank-based estimator Kendall’s tau in the nonparanormal family. Let X = (X1 , . . . , Xd )T be a random ej and X ek be two independent copies of vector. Let X Xj and Xk . The population version of the Kendall’s tau statistic is: ej ), sign(Xk − X ek ) . τ (Xj , Xk ) := Corr sign(Xj − X PCA for the Dependent Data Let x1 , . . . , xn ∈ Rd be n observed data points. The sample version Kendall’s tau statistic is defined as: τbjk (x1 , . . . , xn ) := 2 n(n − 1) X sign(xij − xi0 j ) 1≤i<i0 ≤n · (xik − xi0 k ), (2.1) which is monotone transformation-invariant correlation between the empirical realizations of two random variables Xj and Xk . For x1 , . . . , xn independent, it is easy to verify that Eb τjk (x1 , . . . , xn ) = τ (Xj , Xk ). b = [R bjk ] ∈ Rd×d with We denote R bjk = sin( π τbjk (x1 , . . . , xn )) R 2 the Kendall’s tau matrix. Another interpretation of the Kendall’s tau statistic is that it is an association measure based on the notion of concordance. We call two pairs of real numbers (s, t) and (u, v) concordant if (s − t)(u − v) > 0 and disconcordant if (s − t)(u − v) < 0. Kruskal (1958) show that ej )(Xk − X ek ) > 0 τ (Xj , Xk ) =P (Xj − X ej )(Xk − X ek ) < 0 . (2.2) − P (Xj − X The following theorem, coming from Kruskal (1958), states the invariance property of the relationship between the population Kendall’s tau statistic τ (Xj , Xk ) and the latent correlation coefficient Σjk in the nonparanormal family. Theorem 2.2. Let X := (X1 , . . . , Xd )T ∼ N P Nd (Σ, f ). We denote τ (Xj , Xk ) to be the population Kendall’s tau statistic between Xj and Xk . Then Σjk = sin π2 τ (Xj , Xk ) . Proof. To prove Theorem 2.2, we actually have ej )(Xk − X ek ) > 0 τ (Xj , Xk ) = P (Xj − X ej )(Xk − X ek ) < 0 − P (Xj − X ej ))(fk (Xk ) − fk (X ek )) > 0 = P (fj (Xj ) − fj (X ej ))(fk (Xk ) − fk (X ek )) < 0 − P (fj (Xj ) − fj (X = 2 arcsin(Σjk ). π The last equality is coming from Kruskal (1958)’s result for Gaussian distribution. 2.2. Mixing Conditions In this section we provide definitions of several models of non-independent data. In particular, we will introduce the notions of strong mixing conditions: φmixing and η-mixing. These will be utilized later in analyzing the performance of the proposed method for the dependent data. We first introduce the stationary m-dependence sequences as follows. Definition 2.3 (m-dependence). A stationary sequence X1 , . . . , Xn is said to be m-dependence if and d only if (i) Xi = X for i ∈ {1, . . . , n} and some random vector X ∈ Rd ; (ii) For any s, t ∈ {1, . . . , n}, Xs is independent of Xt whenever |s − t| > m. m-dependency is of particular interest in several fields. For example, in genetics and epidemiology, data might contain samples of families and there are correlation between members of the same family. Moreover, mdependence can be seen as a simplified version of a time series, where Xs and Xt will often be highly dependent if |s − t| is small, with decreasing dependence as |s − t| increases. Next we proceed to some more general weak dependency conditions. In particular, we build the dependence structure on the mixing sequences. To this end, we first introduce the φ measure of dependence as follows. Definition 2.4 (φ measure of dependence). Let (Ω, F, P) be the probability space and A, B ⊂ F be two σ-fields. We define the φ measures of dependence as: φ(A, B) := sup |P(B | A) − P(B)|. (2.3) A∈A,B∈B,P(A)>0 We then describe the strong mixing conditions. Let X = {Xt }∞ t=−∞ be a sequence of random vectors. For −∞ ≤ L ≤ U ≤ ∞, define the σ-field FLU to be FLU := σ(Xi , L ≤ i ≤ U ). For two probability measures µ1 , µ2 on the measurable space (Ω, F), the total variation distance between µ1 and µ2 is defined as: ||µ1 − µ2 ||T V := sup |µ1 (A) − µ2 (A)|. (2.4) A∈F With the above notations, the φ and η dependence coefficients are defined as: Definition 2.5. Let X = {Xi }∞ i=0 be a sequence of random vectors defined in the probability space (Ω, F, P) and Xij be the subsequence. For 0 ≤ L ≤ U ≤ ∞, remind that FLU := σ(Xi , L ≤ i ≤ U ). The φ PCA for the Dependent Data and η dependence coefficients are defined as: ∞ ); φ(m) := sup φ(F0j , Fj+m j∈Z ηij := ess sup ||P(Xjn | X1i−1 = y, Xi = x) y,x,x0 − P(Xjn | X1i−1 = y, Xi = x0 )||T V . The following lemma connects ηij to φ(m). This is coming from Kontorovich & Ramanan (2008) for discrete case and a stronger version applicable to continuous case can be traced back to Samson (2000). This lemma will play an important role later in analyzing our proposed method. For self-containedness, we include a proof here. 3.2. Models One of the intuition of PCA is coming from the Gaussian distribution. In geometric, the principal components define the major axes of the contours of constant probability for the multivariate Gaussian (Jolliffe, 2005). However, such a nice interpretation does not exist anymore when the distributions are away from the Gaussian. Balasubramanian & Schwartz (2002) construct examples such that PCA loses in the sense of preserving the structure of the data to the most. However, under the nonparanormal model and considering the monotone transformation f as a type of data contamination, the geometric intuition of PCA comes back. Proof. By definition, we have In particular, for a positive definite matrix Σ with diag(Σ) = 1, let λ1 ≥ λ2 ≥, . . . , ≥ λd be eigenvalues of Σ and θ1 , . . . , θd be the corresponding eigenvectors. For 0 ≤ q ≤ 1, the `q ball Bq (Rq ) is defined as: ||P(Xjn |X1i−1=y,Xi = x)−P(Xjn |X1i−1=y,Xi = x0 )||T V when q = 0, B0 (R0 ):={v ∈ Rd : card(supp(v)) ≤ R0 }; ≤ ||P(Xjn | X1i−1 = y, Xi = x) − P(Xjn )||T V when 0 < q ≤ 1, Bq (Rq ) := {v ∈ Rd : ||v||qq ≤ Rq }. Lemma 2.6 (Samson (2000)). With the above notations, we have ηij ≤ 2φ(j − i). + ||P(Xjn | X1i−1 = y, Xi = x0 ) − P(Xjn )||T V ≤ 2φ(j − i). Therefore, by the continuity of the probability, we have ηij ≤ 2φ(j − i). This completes the proof. 3. Methods In this section, we briefly review the statistical models of Copula Component Analysis (COCA) proposed by Han & Liu (2012). The aim of COCA is to recover the principal components of the latent Gaussian distribution in the nonparanormal family. 3.1. Scale-invariant PCA PCA is not scale invariant, meaning that variables measured in different scales will result in different estimators (Jolliffe, 2005). To attack this problem, PCA conducted on the sample correlation matrix S0 instead of the sample covariance matrix S is commonly used. We call the procedure of conducting PCA on S0 the scale-invariant PCA. It is realized that a large portion of works claiming doing PCA are actually doing the scale-invariant PCA. It is under debate whether PCA or the scale-invariant PCA are preferred in different circumstances and we refer to Jolliffe (2005) for more discussions on it. In the population level, the scaleinvariant PCA aims at recovering the leading eigenvectors of the correlation matrix, which has the same sparsity pattern as the leading eigenvectors of the covariance matrix. Accordingly, the model M(q, Rq , Σ, f ) is considered: M(q, Rq , Σ, f ) = X : X ∼ N P Nd (Σ, f ), θ1 ∈ Sd−1 ∩ Bq (Rq ) . (3.1) The `q ball induces a (weak) sparsity pattern when 0 ≤ q ≤ 1 and has been analyzed in linear regression (Raskutti et al., 2011) and sparse PCA (Vu & Lei, 2012). Moreover, the data are assumed to come from a nonparanormal distribution, which is a strict extension to the Gaussian distribution. Inspired by the model M(q, Rq , Σ, f ), we consider the following global estimator θe1 , which maximizes the following equation with the constraint that θe1 ∈ Bq (Rq ) for some 0 ≤ q ≤ 1: b θe1 = argmax v T Rv, v∈Rd subject to v ∈ Sd−1 ∩ Bq (Rq ). (3.2) b is the estimated Kendall’s tau matrix. The corHere R responding estimator θe1 can be considered as a nonlinear dimensional reduction procedure and has the potential to gain more flexibility compared with PCA, as shown in the analysis of Han & Liu (2012). 4. Theoretical Properties In this section we provide the theoretical properties of the proposed COCA estimator θe1 as obtained in Equation (3.2) for the dependent data. To our knowledge, PCA for the Dependent Data this is the first work analyzing the theoretical performance of PCA for the dependent data in high dimensions. To provide some insights, we first deliver the rate of θe1 converging to θ1 when the data points are independent from each other. This theorem is coming from Han & Liu (2012). Theorem 4.1 (Independence). Let θe1 be the global optimum in Equation (3.2) and X ∈ M(q, Rq , Σ, f ). Let x1 , . . . , xn be an independent sequence of realizab := [sin( π τbjk (x1 , . . . , xn ))] be defined tions of X and R 2 d−1 as in Equation (2.1). For any two vectors p v1 ∈ S d−1 T and v2 ∈ S , let | sin ∠(v1 , v2 )| = 1 − (v1 v2 )2 . Then we have, with probability at least 1 − 1/d2 , d stationary sequence with Xt = X. Let the global optimum in Equation (3.2), b := [sin( π τbjk (X1 , . . . , Xn ))] is defined R 2 Equation (2.1). Let the parameter θe1 be where as in γm,n := 2(m + 1)2 (1 − m/n) + m(m + 1)(2m + 1)/(3n) represent the effect of dependence on the rate of convergence. Then we have, for any n ≥ 4m2 /(γm,n log d), with probability at least 1 − 1/d2 2 sin ∠(θe1 , θ1 ) ≤ γq Rq2 4π 2 γm,n log d · (λ1 − λ2 )2 n 2−q 2 , (4.2) , (4.1) √ where γq = 2 · I(q = 1) + 4 · I(q = 0) + (1 + 3)2 · I(0 < q < 1) and λj = Λj (Σ) for j = 1, 2. √ where γq = 2 · I(q = 1) + 4 · I(q = 0) + (1 + 3)2 · I(0 < q < 1) and λj = Λj (Σ) for j = 1, 2. Proof. The key is to estimate the convergence rate of b to Σ in the m-dependence setting. The detailed R proof is presented in Han & Liu (2013). Proof. The key idea of the proof is to utilize the `max b to Σ. See Han & Liu norm convergence result of R (2012) for a detailed proof. Remark 4.4. It can be observed in Equation (4.2) that θe1 converges to θ1 in a rate related to both (n, d) and m. Generally, when R q and λ1 , λ2 do not scale sin2∠(θe1 , θ1 ) ≤ γq Rq2 2π 2 log d · (λ1 − λ2 )2 n 2−q 2 2 It can be observed that the convergence rate of θe1 to θ1 will be faster when θ1 lies in a more sparse ball. It makes sense because the effect of “the curse of dimensionality” will be decreasing when the parameters are more and more sparse. Generally, whenRq and λ1 , λ2 with (n, d), the rate is OP ( m is fixed, the rate is optimal. log d 1−q/2 ) n . When m Using Theorem 4.3, we provide the support recovery result using a similar technique as the proof of Corollary 4.2. do not scale with (n, d), the rate is OP ( logn d )1−q/2 , which is the parametric rate Vu & Lei (2012) obtain. In the following, we show that Theorem 4.1 can be applied to derive a support recovery result. Corollary 4.5 (m-dependency). With the settings and notations in Theorem 4.3 held, let Corollary 4.2 (Independence). With the settings and notations in Theorem 4.1 held, let . If we further have b := supp(θe1 ) Θ := supp(θ1 ) and Θ 4R0 πγm,n min |θ1j | ≥ j∈Θ λ1 − λ2 b := supp(θe1 ). Θ := supp(θ1 ) and Θ If we further have r √ 2 2R0 π log d , min |θ1j | ≥ j∈Θ λ1 − λ2 n b = Θ) ≥ 1 − d−2 . then we have P(Θ We then generalize Theorem 4.1 to the nonindependent cases. Here the notions of m-dependence and φ-mixing sequences as defined in Section 2.2 are exploited. We first provide an upper bound for the estimator θe1 when the data points form an mdependence sequence. Theorem 4.3 (m-dependence). Let X ∈ M(q, Rq , Σ, f ) and {Xt }nt=1 be a m-dependence r log d , n b = Θ) ≥ 1 − d−2 . then we have P(Θ We next proceed to bound the angle between θe1 and θ1 in a more general setting of data dependence. Here the dependence is quantified by φ measure as defined in Definitions 2.4 and 2.5. Theorem 4.6 (φ-mixing). Let X ∈ M(q, Rq , Σ, f ) d and {Xt }nt=1 be a stationary sequence with Xt = X. We assume that {Xt }nt=1 satisfies that for any j 6= k ∈ {1, . . . , d} and m ∈ N, sup φ σ((X1 ){j,k} , . . . , (Xi ){j,k} ), i∈Z σ((Xi+m ){j,k} , . . . , (Xn ){j,k} ) ≤ φ(m). PCA for the Dependent Data Here {φ(i)}n−1 i=1 is a sequence of positive numbers. Let θe1 be the global optimum in Equation (3.2), where b := [sin( π τbjk (X1 , . . . , Xn ))] is defined as in EquaR 2 tion (2.1). Let !2 n−1 X γφ := 1 + 2 φ(i) (4.3) i=1 represent the effect of dependence on the rate of convergence. Then we have, supposing i2 hP n−1 n−i φ(i) i=1 n−1 n≥ , γφ log d with probability at least 1 − 1/d2 , sin2∠(θe1 , θ1 ) ≤ γq Rq2 log d 16π 2 γφ · (λ1 − λ2 )2 n 2−q 2 , (4.4) √ where γq = 2 · I(q = 1) + 4 · I(q = 0) + (1 + 3)2 · I(0 < q < 1) and λj = Λj (Σ) for j = 1, 2. Proof. The main proof is to quantify the bias term, and then utilize Lemma 2.6 to build a bridge between the results in Section 5.1 in Kontorovich (2007) and the desired concentration inequality we need. Detailed proof is presented in Han & Liu (2013). Remark 4.7. Sequences satisfying m-dependence is a subset of φ-mixing sequences. However, Theorem 4.3 provides a faster convergence rate than the result in Theorem 4.6. Moreover, the proof techniques in Theorem 4.3 has interesting points itself. It can be observed that when φ(i) decreases fast (i.e., weakdependence), the rate is near-optimal. For example, supposing that Rq and λ1 , λ2 do not scale (n, d), with 1−q/2 log d −2 when φ(i) = O(i ), the rate is OP , n −1 which When is optimal. φ(i) = O(i ), the rate is 1−q/2 log2 n log d OP . n Again, a support recovery result can be provided using a similar technique as the proof of Corollaries 4.2 and 4.5. Corollary 4.8 (φ-mixing). With the settings and notations in Theorem 4.3 held, let b := supp(θe1 ) Θ := supp(θ1 ) and Θ . If we further have 8R0 πγφ min |θ1j | ≥ j∈Θ λ1 − λ2 r b = Θ) ≥ 1 − d−2 . then we have P(Θ log d , n 5. Experiments In this section we investigate the robustness of COCA to data dependence on the synthetic data. We use the truncated power method proposed by Yuan & Zhang (2011) to approximate the global estimator θe1 obtained in Equation (3.2). Two procedures are considered here: Pearson: the classic sparse PCA using the Pearson sample correlation matrix; Kendall: the proposed rank-based scale-invariant PCA method using the Kendall’s tau correlation matrix. In the simulation study we sample n data points x1 , . . . , xn from a certain random vector X ∈ Rd with some type of data dependence. Here we set d = 100. We follow a similar generating scheme as in Han & Liu (2012). A covariance matrix Σ is firstly synthesized through the eigenvalue decomposition, where the first two eigenvalues are given and the corresponding eigenvectors are pre-specified to be sparse. In detail, we let ω1 ≥ ω2 ≥ . . . ≥ ωd be the d eigenvalues of Σ and u1 , u2 , . . . , ud be the corresponding eigenvectors. Suppose that the first two dominant eigenvectors of Σ, u1 and u2 , are sparse in the sense that only the first s = 10 entries of u1 and the second s = 10 entries of u2 nonzero, i.e., 1 1 √ √ 1 ≤ j ≤ 10 11 ≤ j ≤ 20 10 10 and u2j = u1j = 0 otherwise 0 otherwise and ω1 = 5, ω2 = 2, ω3 = . . . = ωd = 1. The remaining eigenvectors are chosen arbitrarily. The correlation matrix Σ0 is accordingly generated and the leading eigenvector of Σ0 is sparse. We aim at recovering the leading eigenvector θ1 . To sample data from the nonparanormal, we also need the transformation functions: f = {fj }dj=1 . Here the following transformation function is considered: There exist five univariate monotone functions h1 , h2 , . . . , h5 : R → R and h = {h1 , h2 , h3 , h4 , h5 , h1 , h2 , h3 , h4 , h5 , . . .}, where h1 (x) := x, h2 (x) := sign(x)|x|1/2 , h3 (x) := x3 , h4 (x) := Φ(x), h5 (x) := exp(x). Here Φ is defined to be the cumulative distribution functions of the standard Gaussian. We then generate n = 100 data points y1 , . . . , yn such that yi ∼ Nd (0, Σ) where Σ is defined as above. To evaluate the robustness of different methods for dependent data, we PCA for the Dependent Data ρ = 0.2 ρ = 0.4 ρ = 0.6 0.4 0.6 0.8 1.0 0.8 0.0 0.2 0.4 0.6 0.8 0.0 1.0 0.0 0.4 0.6 NPN NPN NPN 0.6 0.8 1.0 1.0 0.8 0.2 0.4 TPR 0.6 0.8 0.0 Pearson Kendall 0.0 0.2 FPR 0.4 0.6 0.8 FPR 1.0 Pearson Kendall 0.0 0.2 0.4 TPR 0.6 0.8 0.6 TPR 0.4 0.2 0.4 0.8 1.0 FPR 0.0 0.2 0.2 FPR Pearson Kendall 0.0 Pearson Kendall FPR 1.0 0.2 1.0 0.0 Pearson Kendall 0.0 0.0 Pearson Kendall 0.2 0.4 TPR 0.6 0.8 0.6 TPR 0.4 0.2 0.2 0.4 TPR 0.6 0.8 1.0 Normal 1.0 Normal 1.0 Normal 0.0 0.2 0.4 0.6 0.8 1.0 FPR Figure 1. ROC curves for the Gaussian (Scheme 1) and non-Gaussian (Scheme 2) data (above and below) using the truncated power algorithm are presented. Here the data dependence degrees are at different levels (ρ = 0.2, 0.4, 0.6). n is fixed to be 100 and d = 100. suppose that z1 , . . . , zn follow a stationary vector autoregressive process as defined in Loh & Wainwright (2011). In detail, we assume that z1 = y1 and for some real number 0 ≤ ρ ≤ 1 p zt+1 = ρ · zt + 1 − ρ2 yt+1 , for t = 1, . . . , n − 1. Here we have that zi ∼ Nd (0, Σ) forms a dependent random sequence. Finally, we have the data points x1 , . . . , xn : [Scheme 1] {xi }ni=1 = {zi }ni=1 , with xi ∼ Nd (0, Σ); [Scheme 2] {xi }ni=1 = {h(zi )}ni=1 where h := {h1 , h2 , h3 , h4 , h5 , . . .}, with xi follows a non-Gaussian nonparanormal distribution. The final data matrix we obtained is X = (x1 , . . . , xn )T ∈ Rn×d . The truncated power algorithm is then employed on X to computer the estimated leading eigenvector θe1 . To evaluate the empirical variable selection property of different methods, we define S := {1 ≤ j ≤ d : θ1j = 6 0}, Sbδ := {1 ≤ j ≤ d : θe1j 6= 0}, (5.1) (5.2) to be the support sets of the true leading eigenvector θ1 and the estimated leading eigenvector θe1 using the tuning parameter δ. In this way, the False Positive Number (FPN) and False Negative Number (FNN) of δ are defined as: FPN(δ) FNN(δ) := the number of features in Sbδ not in S, := the number of features in S not in Sbδ . Then we can further define the False Positive Rate(FPR) and False Negative Rate (FNR) corresponding to the tuning parameter δ to be FPR(δ) := FPN(δ)/(d − s), FNR(δ) := FNN(δ)/s. Under the Scheme 1 and Scheme 2 with different levels of dependence (ρ = 0, 0.2, 0.4, 0.6, 0.8), we repeatedly generate the data matrix X for 1,000 times and compute the averaged False Positive Rates and False Negative Rates using a path of tuning parameters δ. The feature selection performances of different methods are then evaluated by plotting (FPR(δ), 1−FNR(δ)). The corresponding ROC curves are presented in Figure 1. There are several observations we can see from Figure 1: (i) With the increase of the data dependence level, both methods’ performance decreases. (ii) Compared with the Gaussian case (Scheme 1), the difference between Pearson and Kendall are larger when the data PCA for the Dependent Data Table 1. Quantitative comparison on the dataset under the generating Scheme 1 and Scheme 2. The means of the oracle false positive and false negative rates with their standard deviations in parentheses are presented. Here n = 100, d = 100 and the dependence degree ρ is increasing from 0 to 0.8. Gaussian(Scheme 1) Pearson non-Gaussian (Scheme 2) Kendall Pearson Kendall ρ FPR(%) FNR FPR FNR FPR FNR FPR FNR 0.0 0.2 0.4 0.6 0.8 1.1(0.80) 1.8(0.89) 4.4(1.30) 10.3(2.19) 20.2(4.30) 0.0(0.00) 0.0(0.14) 0.2(0.47) 2.8(1.77) 20.8(5.87) 1.2(0.52) 1.3(0.69) 2.7(1.06) 8.0(1.92) 18.8(4.17) 0.0(0.00) 0.0(0.00) 0.0(0.20) 1.7(1.39) 18.7(5.52) 17.2(3.60) 17.4(3.42) 18.5(3.56) 20.8(4.40) 24.4(4.93) 7.3(3.72) 7.6(3.66) 10.5(4.52) 16.9(5.70) 27.7(6.38) 1.2(0.52) 1.3(0.69) 2.7(1.06) 8.0(1.92) 18.8(4.17) 0.0(0.00) 0.0(0.00) 0.0(0.20) 1.7(1.39) 18.7(5.52) are generated from Scheme 2. This coincides with the observations in Han & Liu (2012). (iii) When the data dependence degree ρ increases, Kendall performs better than Pearson in both the Gaussian and Nonparanormal cases, meaning that Kendall is more robust to the data dependence than Pearson. To explore the empirical performances of difference methods using the truncated power method more, we define an oracle tuning parameter δ ∗ to be the δ with the lowest FPR(δ) + FNR(δ): δ ∗ := argmin( FPR(δ) + FNR(δ) ). (5.3) δ In this way, an estimator θe1 using the oracle tuning parameter δ ∗ can be calculated and we compute the oracle false positive and false negative rates as: FPR∗ = FPR(δ ∗ ) and FNR∗ = FNR(δ ∗ ). (5.4) We present the means and standard deviations of (FPR∗ , FNR∗ ) in Table 1. There are several observations we can see from Table 1: (i) When ρ is increasing, both methods’ oracle positive and negative rates decrease. (ii) In the perfect Gaussian case (Scheme 1) where the data points are independent of each other (ρ = 0), there is no statistically significant difference between Kendall and Pearson. (iii) There exist statistically significant differences between Kendall and Pearson in Scheme 2, no matter how large the degree of data dependence (ρ) is. (iv) There is a statistically significant difference between Pearson and Kendall for the Gaussian case when ρ = 0.4, and Kendall performs constantly better than Pearson when ρ > 0. In all, Kendall is more robust to the data dependence than Pearson. 6. Conclusion In this paper we analyze both theoretical and empirical performance of a newly proposed high dimensional semiparametric principal component analysis, named Copula Component Analysis (COCA), when the data are dependent. We provide explicit upper bounds of convergence for COCA estimators when the observations are drawn from several different types of noni.i.d. processes. Our results show that COCA can allow weak dependence. To our knowledge, this is the first work analyzing the theoretical performance of PCA for the dependent data in high dimensions. Our result strictly generalize the analysis in Han & Liu (2012) and the techniques we used have the separate interest for analyzing a variety of other multivariate statistical methods. 7. Acknowledgement This research was supported by NSF award IIS1116730. PCA for the Dependent Data References Balasubramanian, M. and Schwartz, E.L. The isomap algorithm and topological stability. Science, 295 (5552):7–7, 2002. d’Aspremont, A., El Ghaoui, L., Jordan, M.I., and Lanckriet, G.R.G. A direct formulation for sparse PCA using semidefinite programming. Computer Science Division, University of California, 2004. d’Aspremont, A., Bach, F., and Ghaoui, L.E. Optimal solutions for sparse principal component analysis. The Journal of Machine Learning Research, 9:1269– 1294, 2008. Fan, J., Qi, L., and Tong, X. Penalized least squares estimation with weakly dependent data. 2012. Han, F. and Liu, H. Semiparametric principal component analysis. In Advances in Neural Information Processing Systems 25, pp. 171–179, 2012. Han, F. and Liu, H. Principal component analysis on high dimensional complex and noisy data. Technical Report, 2013. Johnstone, I.M. and Lu, A.Y. On consistency and sparsity for principal components analysis in high dimensions. Journal of the American Statistical Association, 104(486):682–693, 2009. Liu, H., Lafferty, J., and Wasserman, L. The nonparanormal: Semiparametric estimation of high dimensional undirected graphs. The Journal of Machine Learning Research, 10:2295–2328, 2009. Liu, H., Han, F., Yuan, M., Lafferty, J., and Wasserman, L. High dimensional semiparametric gaussian copula graphical models. Annals of Statistics, 2012. Loh, P.L. and Wainwright, M.J. High-dimensional regression with noisy and missing data: Provable guarantees with non-convexity. Arxiv preprint arXiv:1109.3714, 2011. Ma, Z. Sparse principal component analysis and iterative thresholding. Arxiv preprint arXiv:1112.2432, 2011. McDiarmid, C. On the method of bounded differences. Surveys in combinatorics, 141(1):148–188, 1989. Raskutti, G., Wainwright, M.J., and Yu, B. Minimax rates of estimation for high-dimensional linear regression over ell q-balls. Information Theory, IEEE Transactions on Information Theory, 57(10):6976– 6994, 2011. Samson, P.M. Concentration of measure inequalities for markov chains and phi-mixing processes. The Annals of Probability, 28(1):416–461, 2000. Jolliffe, I. Principal component analysis. Wiley Online Library, 2005. Shen, H. and Huang, J.Z. Sparse principal component analysis via regularized low rank matrix approximation. Journal of multivariate analysis, 99(6):1015– 1034, 2008. Jolliffe, I.T., Trendafilov, N.T., and Uddin, M. A modified principal component technique based on the lasso. Journal of Computational and Graphical Statistics, 12(3):531–547, 2003. Skinner, CJ, Holmes, DJ, and Smith, TMF. The effect of sample design on principal component analysis. Journal of the American Statistical Association, 81 (395):789–798, 1986. Journée, M., Nesterov, Y., Richtárik, P., and Sepulchre, R. Generalized power method for sparse principal component analysis. The Journal of Machine Learning Research, 11:517–553, 2010. Vu, V.Q. and Lei, J. Minimax rates of estimation for sparse pca in high dimensions. Arxiv preprint arXiv:1202.0786, 2012. Kontorovich, L. Measure concentration of strongly mixing processes with applications. PhD thesis, Carnegie Mellon University, 2007. Kontorovich, L.A. and Ramanan, K. Concentration inequalities for dependent random variables via the martingale method. The Annals of Probability, 36 (6):2126–2158, 2008. Kruskal, W.H. Ordinal measures of association. Journal of the American Statistical Association, 53(284): 814–861, 1958. Witten, D.M., Tibshirani, R., and Hastie, T. A penalized matrix decomposition, with applications to sparse principal components and canonical correlation analysis. Biostatistics, 10(3):515–534, 2009. Yuan, X.T. and Zhang, T. Truncated power method for sparse eigenvalue problems. Arxiv preprint arXiv:1112.2679, 2011. Zou, H., Hastie, T., and Tibshirani, R. Sparse principal component analysis. Journal of computational and graphical statistics, 15(2):265–286, 2006.