r - 如何在多列之间进行偏相关分析并通过多个协变量进行控制？

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