Esiste una funzione R che applica una funzione a ciascuna coppia di colonne?

Spesso ho bisogno di applicare una funzione a ciascuna coppia di colonne in un dataframe / matrice e restituire i risultati in una matrice. Ora scrivo sempre un ciclo per farlo. Ad esempio, per creare una matrice contenente i p-value delle correlazioni che scrivo:

df <- data.frame(x=rnorm(100),y=rnorm(100),z=rnorm(100)) n <- ncol(df) foo <- matrix(0,n,n) for ( i in 1:n) { for (j in i:n) { foo[i,j] <- cor.test(df[,i],df[,j])$p.value } } foo[lower.tri(foo)] <- t(foo)[lower.tri(foo)] foo [,1] [,2] [,3] [1,] 0.0000000 0.7215071 0.5651266 [2,] 0.7215071 0.0000000 0.9019746 [3,] 0.5651266 0.9019746 0.0000000 

che funziona, ma è piuttosto lento per matrici molto grandi. Posso scrivere una funzione per questo in R (senza preoccuparmi di ridurre il tempo a metà assumendo un risultato simmetrico come sopra):

 Papply <- function(x,fun) { n <- ncol(x) foo <- matrix(0,n,n) for ( i in 1:n) { for (j in 1:n) { foo[i,j] <- fun(x[,i],x[,j]) } } return(foo) } 

O una funzione con Rcpp:

 library("Rcpp") library("inline") src <- ' NumericMatrix x(xR); Function f(fun); NumericMatrix y(x.ncol(),x.ncol()); for (int i = 0; i < x.ncol(); i++) { for (int j = 0; j < x.ncol(); j++) { y(i,j) = as(f(wrap(x(_,i)),wrap(x(_,j)))); } } return wrap(y); ' Papply2 <- cxxfunction(signature(xR="numeric",fun="function"),src,plugin="Rcpp") 

Ma entrambi sono piuttosto lenti anche su un piccolo set di dati di 100 variabili (ho pensato che la funzione Rcpp sarebbe stata più veloce, ma suppongo che la conversione tra R e C ++ richieda tutto il tempo):

 > system.time(Papply(matrix(rnorm(100*300),300,100),function(x,y)cor.test(x,y)$p.value)) user system elapsed 3.73 0.00 3.73 > system.time(Papply2(matrix(rnorm(100*300),300,100),function(x,y)cor.test(x,y)$p.value)) user system elapsed 3.71 0.02 3.75 

Quindi la mia domanda è:

  1. A causa della semplicità di queste funzioni presumo che questo sia già da qualche parte in R. Esiste una funzione apply o plyr che fa questo? L’ho cercato ma non sono stato in grado di trovarlo.
  2. Se è così, è più veloce?

Non sarebbe più veloce, ma puoi usare l’ outer per semplificare il codice. Richiede una funzione vettoriale, quindi qui ho usato Vectorize per realizzare una versione vettoriale della funzione per ottenere la correlazione tra due colonne.

 df < - data.frame(x=rnorm(100),y=rnorm(100),z=rnorm(100)) n <- ncol(df) corpij <- function(i,j,data) {cor.test(data[,i],data[,j])$p.value} corp <- Vectorize(corpij, vectorize.args=list("i","j")) outer(1:n,1:n,corp,data=df) 

Non sono sicuro se questo risolva il problema in modo corretto, ma dai un’occhiata al pacchetto psych William Revelle. corr.test restituisce l’elenco delle matrici con i coefficienti di correlazione, # of obs, statistica t-test e p-value. So che lo uso sempre (e AFAICS sei anche uno psicologo, quindi può soddisfare anche le tue esigenze). Scrivere loop non è il modo più elegante per farlo.

 library(psych) corr.test(mtcars) ( k < - corr.test(mtcars[1:5]) ) Call:corr.test(x = mtcars[1:5]) Correlation matrix mpg cyl disp hp drat mpg 1.00 -0.85 -0.85 -0.78 0.68 cyl -0.85 1.00 0.90 0.83 -0.70 disp -0.85 0.90 1.00 0.79 -0.71 hp -0.78 0.83 0.79 1.00 -0.45 drat 0.68 -0.70 -0.71 -0.45 1.00 Sample Size mpg cyl disp hp drat mpg 32 32 32 32 32 cyl 32 32 32 32 32 disp 32 32 32 32 32 hp 32 32 32 32 32 drat 32 32 32 32 32 Probability value mpg cyl disp hp drat mpg 0 0 0 0.00 0.00 cyl 0 0 0 0.00 0.00 disp 0 0 0 0.00 0.00 hp 0 0 0 0.00 0.01 drat 0 0 0 0.01 0.00 str(k) List of 5 $ r : num [1:5, 1:5] 1 -0.852 -0.848 -0.776 0.681 ... ..- attr(*, "dimnames")=List of 2 .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... $ n : num [1:5, 1:5] 32 32 32 32 32 32 32 32 32 32 ... ..- attr(*, "dimnames")=List of 2 .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... $ t : num [1:5, 1:5] Inf -8.92 -8.75 -6.74 5.1 ... ..- attr(*, "dimnames")=List of 2 .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... $ p : num [1:5, 1:5] 0.00 6.11e-10 9.38e-10 1.79e-07 1.78e-05 ... ..- attr(*, "dimnames")=List of 2 .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... .. ..$ : chr [1:5] "mpg" "cyl" "disp" "hp" ... $ Call: language corr.test(x = mtcars[1:5]) - attr(*, "class")= chr [1:2] "psych" "corr.test" 

Il 92% delle volte viene speso in cor.test.default e le routine chiamate così è inutile cercare di ottenere risultati più rapidi semplicemente riscrivendo Papply (oltre ai risparmi calcolando solo quelli sopra o sotto la diagonale assumendo che la tua funzione sia simmetrica in y ).

 > M < - matrix(rnorm(100*300),300,100) > Rprof(); junk < - Papply(M,function(x,y) cor.test( x, y)$p.value); Rprof(NULL) > summaryRprof() $by.self self.time self.pct total.time total.pct cor.test.default 4.36 29.54 13.56 91.87 # ... snip ... 

È ansible utilizzare mapply , ma poiché le altre risposte indicano che è improbabile che sia molto più veloce poiché la maggior parte delle volte viene utilizzata da cor.test .

 matrix(mapply(function(x,y) cor.test(df[,x],df[,y])$p.value,rep(1:3,3),sort(rep(1:3,3))),nrow=3,ncol=3) 

È ansible ridurre la quantità di lavoro mapply utilizzando l’ipotesi di simmetria e annotando la diagonale zero, ad es

 v < - mapply(function(x,y) cor.test(df[,x],df[,y])$p.value,rep(1:2,2:1),rev(rep(3:2,2:1))) m <- matrix(0,nrow=3,ncol=3) m[lower.tri(m)] <- v m[upper.tri(m)] <- v