1

我有一个数据框 M,我想计算 M 的列之间的所有成对相关性。我可以使用应用函数轻松完成此操作,例如

pvals = laply(M, function(x) llply(M, function(y) cor.test(x, y)$p.value))

但是,此解决方案所做的工作量是所需工作的 2 倍,因为 x 和 y 之间的相关性与 y 和 x 之间的相关性相同。

我正在寻找一种快速、简单的方法来计算唯一列对之间的所有相关性。我希望结果是一个 NxN 矩阵,其中 N=ncol(M)。我在 Stack Overflow 上搜索了很长时间,但找不到任何这样做的东西。谢谢!

4

2 回答 2

5

对于虹膜数据,您可以执行以下操作:

data(iris)
r <- cor(iris[1:4])

得到相关矩阵。

你可以看看cor.test实际上做了什么stats:::cor.test并找到这个......

    df <- n - 2L
    ESTIMATE <- c(cor = r)
    PARAMETER <- c(df = df)
    STATISTIC <- c(t = sqrt(df) * r/sqrt(1 - r^2))
    p <- pt(STATISTIC, df)

这都是矢量化的,所以你可以运行它。

对维基百科上的不同测试有很好的讨论:http ://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient

于 2013-10-18T22:45:10.500 回答
1

你可以使用combn

#Some data:
DF <- USJudgeRatings
#transform to matrix for better subset performance:
m <- as.matrix(DF)
#use combn and its `FUN` argument: 
res <- matrix(nrow=ncol(DF), ncol=ncol(DF))
res[lower.tri(res)] <- combn(seq_along(DF), 2, function(ind) cor.test(m[, ind[[1]]], m[, ind[[2]]])$p.value)
res[upper.tri(res)] <- t(res)[upper.tri(res)]
diag(res) <- 0

基准:

corpRoland <- function(DF) {
  m <- as.matrix(DF)
  res <- matrix(nrow=ncol(DF), ncol=ncol(DF))
  res[lower.tri(res)] <- combn(seq_along(DF), 2, function(ind) cor.test(m[, ind[[1]]], m[, ind[[2]]])$p.value)
  res[upper.tri(res)] <- t(res)[upper.tri(res)]
  diag(res) <- 0
  res}

corpNeal <- function(DF) {
  cors <- cor(DF)
  df <- nrow(DF)-2
  STATISTIC <- c(t = sqrt(df) * cors/sqrt(1 - cors^2))
  p <- pt(STATISTIC, df)
  matrix(2 * pmin(p, 1 - p),nrow=ncol(DF))}


library(microbenchmark)
DF <- as.data.frame(matrix(rnorm(1e3), ncol=10))
microbenchmark(corpRoland(DF), corpNeal(DF))
#Unit: microseconds
#           expr       min         lq    median       uq       max neval
# corpRoland(DF) 14021.003 14228.2040 14950.212 15157.27 17013.574   100
#   corpNeal(DF)   342.631   351.6775   373.636   385.34   467.773   100

DF <- as.data.frame(matrix(rnorm(1e4), ncol=100))
microbenchmark(corpRoland(DF), corpNeal(DF), times=10)
# Unit: milliseconds
#           expr         min          lq      median          uq         max neval
# corpRoland(DF) 1595.878487 1601.221980 1615.391891 1633.746678 1637.373231    10
#   corpNeal(DF)    8.359662    8.751755    9.021532    9.509576    9.753154    10

因此,您应该使用@NealFultz 的答案。

于 2013-10-18T22:42:10.907 回答