Sparse PCA: reproduction of the synthetic example

The paper:

Hui Zou, Trevor Hastie, and Robert tibshirani,

Sparse principal component analysis,

Journal of computational and graphical statistics, 15 (2): 265-286,200 6.

 

Reproduction of the synthetic example in section 5.2 using R programming:

 1 library(elasticnet) 2  3 ## sample version of SPCA 4 n = 1000 5 v1 = rnorm(n,0,sqrt(290)) 6 v2 = rnorm(n,0,sqrt(300)) 7 v3 = -.3*v1 + 0.925*v2 + rnorm(n) 8 x1 = v1 + rnorm(n) 9 x2 = v1 + rnorm(n)10 x3 = v1 + rnorm(n)11 x4 = v1 + rnorm(n)12 13 x5 = v2 + rnorm(n)14 x6 = v2 + rnorm(n)15 x7 = v2 + rnorm(n)16 x8 = v2 + rnorm(n)17 18 x9 = v3 + rnorm(n)19 x10 = v3 + rnorm(n)20 21 x = cbind(x1,x2,x3,x4,x5,x6,x7,x8,x9,x10)22 x.cov = t(x) %*% x/n; head(x.cov)23 a = spca(x, 2, type=‘predictor‘, sparse=‘varnum‘, para=c(4,4), lambda=0)24 a
25 ## population version of SPCA26 g1 = matrix(290, 4, 4)27 diag(g1) = 29128 29 g2 = matrix(300, 4, 4)30 diag(g2) = 30131 32 g3 = matrix(283.7875, 2, 2)33 diag(g3) = diag(g3)+134 35 36 g1g3 = matrix(-87, 4, 2)37 g2g3 = matrix(277.5, 4, 2)38 39 # construct the exact covariance matrix40 x.cov = matrix(0, 10, 10)41 x.cov[1:4,1:4] = g142 x.cov[5:8,5:8] = g243 x.cov[9:10,9:10] = g344 x.cov[1:4,9:10] = g1g345 x.cov[9:10,1:4] = t(g1g3)46 x.cov[5:8,9:10] = g2g347 x.cov[9:10,5:8] = t(g2g3)48   49 50 b = spca(x.cov, 2, type=‘Gram‘, sparse=‘varnum‘, para=c(4,4), lambda=0)51 b

The results of the population version using exact covariance matrix are exactly as in the paper:

> bCall:spca(x = x.cov, K = 2, para = c(4, 4), type = "Gram", sparse = "varnum",     lambda = 0)2 sparse PCs Pct. of exp. var. : 40.9 39.5 Num. of non-zero loadings :  4 4 Sparse loadings       PC1 PC2 [1,] 0.0 0.5 [2,] 0.0 0.5 [3,] 0.0 0.5 [4,] 0.0 0.5 [5,] 0.5 0.0 [6,] 0.5 0.0 [7,] 0.5 0.0 [8,] 0.5 0.0 [9,] 0.0 0.0[10,] 0.0 0.0

But the sample version may randomly vary a little.

> aCall:spca(x = x, K = 2, para = c(4, 4), type = "predictor", sparse = "varnum",     lambda = 0)2 sparse PCs Pct. of exp. var. : 37.9 37.6 Num. of non-zero loadings :  4 4 Sparse loadings        PC1    PC2x1   0.000 -0.303x2   0.000 -0.533x3   0.000 -0.576x4   0.000 -0.540x5  -0.492  0.000x6  -0.287  0.000x7  -0.481  0.000x8  -0.666  0.000x9   0.000  0.000x10  0.000  0.000

 

Having fun learning sparse PCA!

Sparse PCA: reproduction of the synthetic example