如何进行多柱间的偏相关分析和多个协变量的控制？

Question

我希望对由多个协变量控制的多列之间进行偏相关分析，然后提取r和p值。我的真实数据有一些丢失的值，。

我发现这个答案(Pairwise partial correlation of a matrix, controlling by one variable)可能很有用，所以我将这个方法调整到我的代码中。因为我缺少的值，所以我不能使用ppcor::pcor.test() (被描述为“缺失值不允许”)来实现我的目标。

在这里，我使用内置数据集mtcars来显示我遇到的问题。

# load "ggm" packages to perform partial correlation analysis
library(ggm)

# subset mtcars dataset and make some datapoints as missing values
mydata <- cbind(mtcars[1:8])
mydata[4:10,3] <- rep(NA,7)
mydata[1:5,4] <- NA

# perform partial correlation analysis among the first 6 columns with the last two columns as covariates
sapply(1:(ncol(mydata)-2), function(x) sapply(1:(ncol(mydata)-2), function(y) {
  if (x == y) 1
  else ggm::pcor(c(mydata[,x], mydata[,y], mydata[,7], mydata[,8]),var(mydata))
}))

# error:
Error in S[u, u] : subscript out of bounds

在这个步骤中我得到了错误，既不能执行偏相关，也不能提取r和p值。

非常感谢你的帮助！

艾拉

Answer 1

您不需要在pcor函数中传递列值。您可以传递列号或列名。试着：

sapply(1:(ncol(mydata)-2), function(x) sapply(1:(ncol(mydata)-2), function(y) {
   if (x == y) 1
   else ggm::pcor(c(x, y,7,8),var(mydata))
}))

#           [,1]       [,2] [,3] [,4]       [,5]       [,6]
#[1,]  1.0000000 -0.7208025   NA   NA  0.5717984 -0.8260219
#[2,] -0.7208025  1.0000000   NA   NA -0.6969510  0.7414846
#[3,]         NA         NA    1   NA         NA         NA
#[4,]         NA         NA   NA    1         NA         NA
#[5,]  0.5717984 -0.6969510   NA   NA  1.0000000 -0.5510354
#[6,] -0.8260219  0.7414846   NA   NA -0.5510354  1.0000000