相关系数的排列

Posted 2023-03-12

【中文标题】相关系数的排列

【英文标题】：Permutations of correlation coefficients

【发布时间】：2013-09-29 11:24:56

【问题描述】：

我的问题是关于。

        A<-data.frame(A1=c(1,2,3,4,5),B1=c(6,7,8,9,10),C1=c(11,12,13,14,15 ))

        B<-data.frame(A2=c(6,7,7,10,11),B2=c(2,1,3,8,11),C2=c(1,5,16,7,8))

          cor(A,B)

          #           A2        B2       C2
          # A1 0.9481224 0.9190183 0.459588
          # B1 0.9481224 0.9190183 0.459588
          # C1 0.9481224 0.9190183 0.459588

我获得了这种相关性，然后想执行置换测试以检查相关性是否仍然成立。

我做了如下排列：

              A<-as.vector(t(A))
              B<-as.vector(t(B))

     corperm <- function(A,B,1000) 
     # n is the number of permutations
     # x and y are the vectors to correlate
    obs = abs(cor(A,B))
    tmp = sapply(1:n,function(z) abs(cor(sample(A,replace=TRUE),B)))
   return(1-sum(obs>tmp)/n)
     

结果是

         [1] 0.645

并使用“cor.test”

cor.test(A,B)

Pearson's product-moment correlation

data:  A and B
t = 0.4753, df = 13, p-value = 0.6425
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0.4089539  0.6026075
sample estimates:
cor 
0.1306868 

如何绘制图表或直方图以显示实际相关性和置换数据中的置换相关值???

【参考方案1】：

首先，你不可能完全按照这种方式来做......

> corperm = function(A,B,1000) 
Error: unexpected numeric constant in "corperm = function(A,B,1000"

第三个参数没有名字，但它应该有一个！也许你的意思是

> corperm <- function(A, B, n=1000) 
# etc

然后你需要考虑你想要实现什么。最初，您有两个包含 3 个变量的数据集，然后将它们折叠成两个向量并计算置换向量之间的相关性。为什么有意义？置换数据集的结构应与原始数据集相同。

obs = abs(cor(A,B))
tmp = sapply(1:n,function(z) abs(cor(sample(A,replace=TRUE),B)))
return(1-sum(obs>tmp)/n)

为什么在这里使用 replace=TRUE？如果您想要引导 CI-s，这是有道理的，但是 (a) 最好使用专用功能，然后例如从包引导引导，并且 (B) 您需要对 B 执行相同操作，即样本（B，替换=真）。

对于置换测试，您无需替换即可采样，无论您对 A 和 B 进行还是仅对 A 进行都没有区别。

以及如何获得直方图？好吧， hist(tmp) 会为您绘制置换值的直方图，而 obs 是观察到的相关性的绝对值。

HTHAB

（编辑）

corperm <- function(x, y, N=1000, plot=FALSE)
    reps <- replicate(N, cor(sample(x), y))
    obs <- cor(x,y)
    p <- mean(reps > obs) # shortcut for sum(reps > obs)/N
    if(plot)
        hist(reps)
        abline(v=obs, col="red")
        
     p
     

现在您可以在一对变量上使用它：

corperm(A[,1], B[,1])

要将其应用于所有对，请使用 for 或 mapply。 for 更容易理解，所以我不会坚持使用 mapply 来获取所有可能的配对。

res <- matrix(NA, nrow=NCOL(A), ncol=NCOL(B))
for(iii in 1:3) for(jjj in 1:3) res[iii,jjj] <- corperm(A[,iii], B[,jjj], plot=FALSE)
rownames(res)<-names(A)
colnames(res) <- names(B)
print(res)

要制作所有直方图，请在上面使用 plot=TRUE。

【讨论】：

@lebatsnok ：您能否建议如何排列这些相关性并获得直方图以显示相关性仍然成立。我按照教程完成了上述操作。

嗨，保罗。你能说出你的确切研究问题吗？您是否有兴趣测试数据框 A 和 B 中每对变量之间相关性的重要性？或者是别的什么？或者您所说的“相关性仍然成立”是什么意思？

@lebastsnok ：是的，测试两个数据帧中变量对之间相关性的显着性。

一种方法可能是首先为一对变量编写一个执行此操作的函数，然后在您需要的所有变量对上使用它。 - 见上面的编辑

在 corperm 函数中，未正确计算 p 值。而不是p <- mean(obs > reps)，它应该是p <- mean(reps > obs)。

【参考方案2】：

我认为对两个变体的相关性分析做置换检验意义不大，因为cor.test()function提供的“p.value”与置换检验的效果相同。