Qual è il modo più veloce per calcolare i primi due componenti principali in R?

Question

Sto usando princomp in R per eseguire PCA. La mia matrice di dati è enorme (10K x 10K con ogni valore fino a 4 punti decimali). Sono necessarie ~ 3,5 ore e ~ 6,5 GB di memoria fisica su un processore Xeon 2,27 GHz.

Poiché desidero solo i primi due componenti, c'è un modo più veloce per farlo?

Aggiornamento:

Oltre alla velocità, C'è un modo efficiente della memoria per fare questo?

Sono necessarie ~ 2 ore e ~ 6.3 GB di memoria fisica per il calcolo dei primi due componenti utilizzando svd(,2,).

Answer 1

A volte si accede alle cosiddette scomposizioni "economiche" che consentono di limitare il numero di autovalori/autovettori. Sembra che eigen() e prcomp() non lo offrono, ma svd() consente di specificare il numero massimo da calcolare.

Su piccoli matrici, i guadagni sembrano modesti:

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> 

ma il fattore tre rispetto al princomp() può essere vale la pena ricostruire princomp() da svd() come svd() consente di fermare dopo due valori.

Answer 2

È possibile scrivere autonomamente la funzione e arrestarsi su 2 componenti. Non è troppo difficile. Ce l'ho in giro da qualche parte, se lo trovo lo posterò.

Answer 3

è possibile utilizzare l'approccio di rete neurale per trovare il componente principale. descrizione di base è data qui .. http://www.heikohoffmann.de/htmlthesis/node26.html

primo principale componente, y = w1 * x1 + w2 * x2 e secondo componente ortogonale può essere calcolato come q = w2 * x1-x2 w1 *.

Answer 4

Il power method potrebbe essere quello che vuoi. Se lo si codifica in R, che non è affatto difficile, penso che si possa scoprire che non è più veloce dell'approccio SVD suggerito in un'altra risposta, che fa uso di routine compilate LAPACK.

Answer 5

Il pacchetto 'svd' fornisce le routine per SVD/eigendecomposition troncati tramite l'algoritmo di Lanczos. Puoi usarlo per calcolare solo i primi due componenti principali.

Qui ho:

> 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

ho provato implementazione del pacchetto di pcaMethods dell'algoritmo NIPALS. Di default calcola i primi 2 componenti principali. Risulta essere più lento rispetto agli altri metodi suggeriti.

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