¿Cuál es la forma más rápida de calcular los dos primeros componentes principales en R?

Question

Estoy usando princomp en R para realizar PCA. Mi matriz de datos es enorme (10K x 10K con cada valor hasta 4 decimales). Toma ~ 3.5 horas y ~ 6.5 GB de memoria física en un procesador Xeon a 2.27 GHz.

Como solo quiero los dos primeros componentes, ¿hay una forma más rápida de hacerlo?

Actualización:

Además de la velocidad, ¿Hay una manera eficiente de la memoria para hacer esto?

Se necesitan ~ 2 horas y ~ 6.3 GB de memoria física para calcular los dos primeros componentes usando svd(,2,).

Answer 1

A veces tiene acceso a descomposiciones llamadas 'económicas' que le permiten limitar el número de autovalores/eigenvectores. Parece que eigen() y prcomp() no ofrecen esto, pero svd() le permite especificar el número máximo para calcular.

En matrices pequeñas, las ganancias parecen modestos:

R> set.seed(42); N <- 10; M <- matrix(rnorm(N*N), N, N) 
R> library(rbenchmark) 
R> benchmark(eigen(M), svd(M,2,0), prcomp(M), princomp(M), order="relative") 
      test replications elapsed relative user.self sys.self user.child 
2 svd(M, 2, 0)   100 0.021 1.00000  0.02  0   0 
3 prcomp(M)   100 0.043 2.04762  0.04  0   0 
1  eigen(M)   100 0.050 2.38095  0.05  0   0 
4 princomp(M)   100 0.065 3.09524  0.06  0   0 
R> 

pero el factor de tres con respecto a princomp() puede valer la pena reconstruir princomp() de svd() como svd() le permite parar después de dos valores.

Answer 2

Puede escribir la función usted mismo y detenerse en 2 componentes. No es muy difícil. Lo tengo por ahí, si lo encuentro, lo publicaré.

Answer 3

puede utilizar el enfoque de red neuronal para encontrar el componente principal. Descripción básica se da aquí .. http://www.heikohoffmann.de/htmlthesis/node26.html

primer componente principal, y = w1 * x1 + x2 w2 * y componentes ortogonales En segundo lugar se pueden calcular como q = w2 * x1-x2 * W1.

Answer 4

El power method podría ser lo que usted desea. Si lo codifica en R, que no es difícil en absoluto, creo que puede encontrar que no es más rápido que el enfoque SVD sugerido en otra respuesta, que hace uso de las rutinas compiladas de LAPACK.

Answer 5

El paquete 'svd' proporciona las rutinas para la descomposición SVD/eigen truncada a través del algoritmo Lanczos. Puede usarlo para calcular solo los dos primeros componentes principales.

Aquí me tienen:

> library(svd) 
> set.seed(42); N <- 1000; M <- matrix(rnorm(N*N), N, N) 
> system.time(svd(M, 2, 0)) 
    user system elapsed 
    7.355 0.069 7.501 
> system.time(princomp(M)) 
    user system elapsed 
    5.985 0.055 6.085 
> system.time(prcomp(M)) 
    user system elapsed 
    9.267 0.060 9.368 
> system.time(trlan.svd(M, neig = 2)) 
    user system elapsed 
    0.606 0.004 0.614 
> system.time(trlan.svd(M, neig = 20)) 
    user system elapsed 
    1.894 0.009 1.910 
> system.time(propack.svd(M, neig = 20)) 
    user system elapsed 
    1.072 0.011 1.087 

Answer 6

me trataron aplicación del pcaMethods del paquete del algoritmo NIPALS. Por defecto, calcula los primeros 2 componentes principales. Resulta ser más lento que los otros métodos sugeridos.

set.seed(42); N <- 10; M <- matrix(rnorm(N*N), N, N) 
library(pcaMethods) 
library(rbenchmark) 
m1 <- pca(M, method="nipals", nPcs=2) 
benchmark(pca(M, method="nipals"), 
      eigen(M), svd(M,2,0), prcomp(M), princomp(M), order="relative") 

         test replications elapsed relative user.self sys.self 
3    svd(M, 2, 0)   100 0.02  1.0  0.02  0 
2     eigen(M)   100 0.03  1.5  0.03  0 
4     prcomp(M)   100 0.03  1.5  0.03  0 
5    princomp(M)   100 0.05  2.5  0.05  0 
1 pca(M, method = "nipals")   100 0.23  11.5  0.24  0