r - R 中的快速方差分析计算

Question

我有一个具有以下尺寸的数据框:

dim(b)  
[1]    974 433685

这些列代表我想要对其运行方差分析的变量(即，我想要运行 433,685 个方差分析)。样本大小为 974。最后一列是“组”变量。

我想出了 3 种不同的方法，但由于测试次数太多，所有方法都太慢了。

首先，让我们生成一个小的练习数据集来玩:

dat = as.data.frame(matrix(runif(10000*500), ncol = 10000, nrow = 500))
dat$group = rep(letters[1:10], 5000)

方法一(基于'sapply'):

system.time(sapply(dat[,-length(dat)], function(x) aov(x~group, data=dat) ))

   user  system elapsed 
 143.76    0.33  151.79 

方法 2(基于“parallel”包中的“mclapply”):

library(parallel)
options(mc.cores=3)
system.time(mclapply(dat[,-length(dat)], function(x) aov(x~group, data=dat) ))

   user  system elapsed 
 141.76    0.21  142.58 

方法 3(基于 'cbind'-ing the LHS ):

formula = as.formula( paste0("cbind(", paste(names(dat)[-length(dat)],collapse=","), ")~group") ) 
system.time(aov(formula, data=dat))

  user  system elapsed 
  10.00    0.22   10.25 

在实践数据集中，方法 3 明显胜出。然而，当我对我的实际数据执行此操作时，使用方法 3 计算仅 10 个(共 433,685 个)列需要这么长时间:

   user  system elapsed
119.028   5.430 124.414

不确定为什么我的实际数据需要更长的时间。我可以访问具有超过 16 个内核和 72GB RAM 的 Linux 集群。

有什么方法可以更快地计算出来吗？

Answer 1

使用相同的 design matrix 同时拟合许多一般线性模型(如方差分析) , Bioconductor/R limma package提供非常快速的 lmFit() 函数。这是如何使用 limma 拟合方差分析模型:

library(limma)

# generate some data 
# (same dimensions as in your question)
nrows <- 1e4
ncols <- 5e2
nlevels <- 10
dat <- matrix(
  runif(nrows * ncols), 
  nrow = nrows, 
  ncol = ncols
)
group <- factor(rep(
  letters[1:nlevels], 
  ncols / nlevels
))

# construct the design matrix
# (same as implicitly used in your question)
dmat <- model.matrix(~ group)
# fit the ANOVA model
fit <- lmFit(dat, dmat)

在我的笔记本电脑上，它在 0.4 - 0.45 秒内完成，数据维度与您问题中的数据相同。