用apply函数替换在R中的dfs列表上迭代运行向后逐步回归的for循环，以减少计算时间

Question

R脚本和文件夹的简化版本，其中包含多个csv文件格式的数据集，都可以在我的GitHub存储库中找到。

在我的名为"LASSO代码"的脚本中，在将一个充满N个csv文件格式的数据集的文件夹加载到R中并将它们全部分配给一个名为"数据集"的列表后，我运行以下代码来拟合N个LASSO回归，每个数据集一个:

set.seed(11)     # to ensure replicability
LASSO_fits <- lapply(dfs, function(i) 
enet(x = as.matrix(select(i, starts_with("X"))), 
y = i$Y, lambda = 0, normalize = FALSE))

现在，我想为向后消除逐步回归复制相同的过程，我们将通过使用stats库中的step()函数来保持简单，使用另一个apply函数，而不必使用循环。问题是这样的，我唯一知道如何做到这一点是通过初始化或准备它之前运行它首先建立:

set.seed(100)      # for reproducibility
full_fits <- vector("list", length = length(dfs))
Backward_Stepwise_fits <- vector("list", length = length(dfs))

然后才拟合所有的Backward_Stepwise_fits，但我不知道如何把full_fits和Backward_Stepwise_fits都放入同一个apply函数，我能想到的唯一方法是使用for循环并将它们堆叠在一起，但这将是非常计算效率低下。我要运行这两种方法的数据集N的数量是26万!

我写了一个for循环，它实际上可以运行，但它花了12个多小时才在58,500个数据集上运行完，这是不可接受的慢。我使用的代码如下:

set.seed(100)      # for reproducibility
for(i in seq_along(dfs)) {
full_fits[[i]] <- lm(formula = Y ~ ., data = dfs[[i]])
Backward_Stepwise_fits[[i]] <- step(object = full_fits[[i]], 
scope = formula(full_fits[[i]]),
direction = 'backward', trace = 0) }

我尝试了以下操作，但在控制台中得到相应的错误信息:

> full_model_fits <- lapply(dfs, function(i)
+   lm(formula = Y ~ ., data = dfs))
Error in terms.formula(formula, data = data) : 
duplicated name 'X1' in data frame using '.'

Answer 1

有没有想过将整个过程并行化?

首先，您可以更简洁地定义代码。

system.time(
res <- lapply(lst, (X) {
full <- lm(Y ~ ., X)
back <- step(full, scope=formula(full), dir='back', trace=FALSE)
})
)
#  user  system elapsed 
# 3.895   0.008   3.897 
system.time(
res1 <- lapply(lst, (X) step(lm(Y ~ ., X), dir='back', trace=FALSE))
)
#  user  system elapsed 
# 3.820   0.016   3.833 
stopifnot(all.equal(res, res1))

结果相等，但没有时差。

现在，使用parallel::parLapply.

library(parallel)
CL <- makeCluster(detectCores() - 1L)
system.time(
res2 <- parLapply(CL, lst, (X) step(lm(Y ~ ., X), dir='back', trace=FALSE))
)
#  user  system elapsed 
# 0.075   0.032   0.861 
stopCluster(CL)
stopifnot(all.equal(res, res2))

在这台机器上大约快4.5倍。

如果我们现在想使用res2作为scope=运行正向逐步选择，，我们需要parallel::clusterMap,parLapply的多元变量:

# CL <- makeCluster(detectCores() - 1L)
res3 <- clusterMap(CL, (X, Y) step(lm(Y ~ 1, X), scope=formula(Y), dir='forw', trace=FALSE),
lst, res2)
# stopCluster(CL)

注意:这与使用您在注释中显示的for循环产生相同的系数:

stopifnot(all.equal(lapply(FS_fits, coef), unname(lapply(res3, coef))))

---

错误duplicated name 'X1' in data frame using '.'的意思是，在一些数据集中有两个列命名为"X1"。下面是找到它们的方法:

names(lst$dat6)[9] <- 'X1'  ## producing duplicated column X1 for demo 
sapply(lst, (x) anyDuplicated(names(x)))
# dat1  dat2  dat3  dat4  dat5  dat6  dat7  dat8  dat9 dat10 dat11 
# 0     0     0     0     0     9     0     0     0     0     0 
# ...

结果显示，在数据集dat6中，9列是(第一个)重复。其他的都是干净的。

---

数据:

set.seed(42)
n <- 50
lst <- replicate(n, {dat <- data.frame(matrix(rnorm(500*30), 500, 30))
cbind(Y=rowSums(as.matrix(dat)%*%rnorm(ncol(dat))) + rnorm(nrow(dat)), dat)}, simplify=FALSE) |> 
setNames(paste0('dat', seq_len(n)))