是否有一个R函数将函数应用于每对列？

Question

我经常需要对数据框/矩阵中的每对列应用函数，并将结果以矩阵形式返回。现在我总是写一个循环来做到这一点。例如，为了使含我写相关的p值的矩阵：

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 

其作品，但对于非常大的矩阵相当缓慢。

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) 
} 

或用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<double>(f(wrap(x(_,i)),wrap(x(_,j)))); 
    } 
} 
return wrap(y); 
' 

Papply2 <- cxxfunction(signature(xR="numeric",fun="function"),src,plugin="Rcpp") 

但两者都相当我可以在R（不与假设如上对称的结果切削时间缩短了一半打扰）写一个函数为这个减缓甚至在100个变量的一个非常小的数据集（我认为RCPP功能会更快，但我猜R和C之间的转换++所有的时间采取它的通行费）：

> 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 

所以我的问题是：

1. 由于这些函数的简单性，我认为这已经在R的某个地方了。是否有应用程序或plyr函数执行此操作？我一直在寻找它，但一直没能找到它。
2. 如果是这样，它是否更快？

Answer 1

它不会更快，但您可以使用outer来简化代码。它确实需要一个矢量化函数，所以在这里我使用Vectorize来创建函数的矢量化版本以获得两列之间的相关性。

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) 

Answer 2

我不确定这是否以正确的方式解决您的问题，但看看William Revelle的psych包。 corr.test返回具有相关系数，obs数，t检验统计量和p值的矩阵列表。我知道我一直都在使用它（而AFAICS你也是一名心理学家，所以它也可以满足你的需求）。编写循环并不是这样做的最优雅的方式。

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" 

Answer 3

的时间

92％被消耗在cor.test.default和例程调用所以它没有希望通过简单地重写Papply（除储蓄从计算只有那些高于或低于对角线假设以获得更快的结果您函数在x和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 ... 

Answer 4

您可以使用mapply，但其他的答案陈述其不太可能更快，因为大多数的时间是由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) 

你可以通过使用对称的假设，并指出零对角线减少工作mapply做多少，例如

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