30

Cor.test()接受向量xy作为参数,但我有一个完整的数据矩阵,我想成对测试。Cor()将此矩阵作为参数就好了,我希望找到一种方法来为cor.test().

其他人的常见建议似乎是使用cor.prob()

https://stat.ethz.ch/pipermail/r-help/2001-November/016201.html

但是这些 p 值与由cor.test()!!!生成的 p 值不同。Cor.test()似乎也比cor.prob().

有人有其他选择cor.prob()吗?如果解决方案涉及嵌套的 for 循环,那就这样吧(我足够新R,即使这对我来说也是个问题)。

4

5 回答 5

45

corr.test包中psych的目的是为了做到这一点:

library("psych")
data(sat.act)
corr.test(sat.act)

如评论中所述,要将基函数的p值复制到cor.test()整个矩阵上,则需要关闭p值的调整以进行多重比较(默认是使用 Holm 的调整方法):

 corr.test(sat.act, adjust = "none")

[但在解释这些结果时要小心!]

于 2012-10-28T20:38:58.187 回答
15

如果您严格追求矩阵格式的 pvalues,那么这是从文森特(链接cor.test)无耻地窃取的解决方案:

cor.test.p <- function(x){
    FUN <- function(x, y) cor.test(x, y)[["p.value"]]
    z <- outer(
      colnames(x), 
      colnames(x), 
      Vectorize(function(i,j) FUN(x[,i], x[,j]))
    )
    dimnames(z) <- list(colnames(x), colnames(x))
    z
}

cor.test.p(mtcars)

注意:Tommy 还提供了一个更快的解决方案,但不太容易实现。哦,没有 for 循环 :)

编辑v_outer我的包中有一个函数qdapTools可以让这个任务变得非常简单:

library(qdapTools)
(out <- v_outer(mtcars, function(x, y) cor.test(x, y)[["p.value"]]))
print(out, digits=4)  # for more digits
于 2012-10-28T19:45:26.807 回答
7

可能最简单的方法是使用rcorr()来自 Hmisc 的。它只需要一个矩阵,因此rcorr(as.matrix(x))如果您的数据位于 data.frame 中,请使用它。它将返回一个列表,其中包含:1) r 成对矩阵,2) 成对 n 矩阵,3) r 的 p 值矩阵。它会自动忽略丢失的数据。

理想情况下,这种函数也应该采用 data.frames 并输出符合“新统计”的置信区间。

于 2015-05-26T18:58:58.383 回答
5

公认的解决方案(心理包中的 corr.test 函数)有效,但对于大型矩阵来说非常慢。我正在使用与药物敏感性矩阵(~1,000 x ~500)相关的基因表达矩阵(~20,000 x ~1,000),我不得不停止它,因为它需要永远。

我从 psych 包中获取了一些代码,并直接使用了 cor() 函数,并得到了更好的结果:

# find (pairwise complete) correlation matrix between two matrices x and y
# compare to corr.test(x, y, adjust = "none")
n <- t(!is.na(x)) %*% (!is.na(y)) # same as count.pairwise(x,y) from psych package
r <- cor(x, y, use = "pairwise.complete.obs") # MUCH MUCH faster than corr.test()
cor2pvalue = function(r, n) {
  t <- (r*sqrt(n-2))/sqrt(1-r^2)
  p <- 2*(1 - pt(abs(t),(n-2)))
  se <- sqrt((1-r*r)/(n-2))
  out <- list(r, n, t, p, se)
  names(out) <- c("r", "n", "t", "p", "se")
  return(out)
}
# get a list with matrices of correlation, pvalues, standard error, etc.
result = cor2pvalue(r,n)

即使有两个 100 x 200 矩阵,差异也是惊人的。一两秒对 45 秒。

> system.time(test_func(x,y))
   user  system elapsed 
  0.308   2.452   0.130 
> system.time(corr.test(x,y, adjust = "none"))
   user  system elapsed 
 45.004   3.276  45.814 
于 2017-09-14T01:59:39.053 回答
-1

“公认的解决方案(corr.test心理包中的函数)有效,但对于大型矩阵来说非常慢。”

如果你使用ci=FALSE,那么速度会快很多。默认情况下,会找到置信区间。但是,这会导致速度略有下降。所以,对于rs, tsand ps, set ci=FALSE

于 2018-09-25T09:23:06.597 回答