使用并行编程在 R 中编程 50000 个回归

Question

我有以下作业问题，我已经完成但似乎需要很长时间才能完成：

假设Y、X1、···、X1000均为均值为0、标准差为1的正态随机变量，且相互独立。生成 Y、X1、···、X1000 的 30 个样本。现在重复以下 50000 次：从 X1 中随机选取十个变量，. . ., X1000, 运行 Y 对这十个变量的线性回归并记录 R2。计算 50000 个 R2 的最大值。

这是我的代码，它适用于 8000 次回归（在我的 macbook pro 的每个核心上 1000 次回归），但似乎无法在每个核心上完成 6250 次回归（总共 50000 次回归）。这是我的代码：

    library(snow)
    cl <- makeCluster(8, type = "SOCK")
    invisible(clusterEvalQ(cl, reg_cluster <- function(rep, samples, n) {
     X <- list()
      R <- rep(0, rep)
      for (k in 1:rep) {
        Y <- rnorm(samples)
        for (j in 1:n) {
          X[[j]] <- rnorm(samples)
    }
   X_1 <- sample(X, 10, replace = FALSE)
   X_1_unlist <- unlist(X_1)
   X.1 <- matrix(X_1_unlist[1:30], ncol = 1)
   X.2 <- matrix(X_1_unlist[31:60], ncol = 1)
   X.3 <- matrix(X_1_unlist[61:90], ncol = 1)
   X.4 <- matrix(X_1_unlist[91:120], ncol = 1)
   X.5 <- matrix(X_1_unlist[121:150], ncol = 1)
   X.6 <- matrix(X_1_unlist[151:180], ncol = 1)
   X.7 <- matrix(X_1_unlist[181:210], ncol = 1)
   X.8 <- matrix(X_1_unlist[211:240], ncol = 1)
   X.9 <- matrix(X_1_unlist[241:270], ncol = 1)
   X.10 <- matrix(X_1_unlist[271:300], ncol = 1)
   X_data <- cbind(X.1, X.2, X.3, X.4, X.5, X.6, X.7, X.8, X.9, X.10)
   X_data <- as.data.frame (X_data)
   names(X_data) <- c("X1", "X2", "X3", "X4", "X5", "X6", "X7", "X8", "X9", "X10")
   attach(X_data)
   reg <- lm(Y ~ X1 + X2 + X3 + X4 + X5 + X6 + X7 + X8 + X9 + X10)
   R[k] <- summary(reg)$r.squared
  }
return(max(R))
}))
results <- clusterEvalQ(cl, reg_cluster(1000, 30, 1000))
results <-clusterEvalQ(cl, reg_cluster(6250, 30, 1000))
stopCluster(cl)
max_results <- c(results[[1]], results[[2]], results[[3]], results[[4]],
                 results[[5]], results[[6]], results[[7]], results[[8]])
max(max_results)

这里还有一点需要注意。每次我运行一个新的回归，Y 和所有的 X 都会重新生成。没有随机变量从一个回归转移到下一个。

所以我的问题是，我怎样才能使这个运行更快？

另外，谁能告诉我为什么它在 12 分钟后完成了 8000 次回归，但在 2.5 小时后仍然没有完成 50000 次回归？

编辑：以下程序已被教授确认：

1) 为Y, X1, ..., X1000 中的每一个生成30 个随机标准正态变量。 Y 总共有 30 个随机正态变量，所有 X 总共有 30 x 1,000 = 30,000 个随机正态变量（每个 30 个）

2) 随机 select X 的 1000 个选择中的十个（例如 X726、X325、X722、X410、X46、X635、X822、X518、X773、X187）

3) 运行使用 R 中的 lm 函数对 Y ~ 10 X 进行线性回归。Y 将有 30 个观测值，而每个 X 也将有 30 个观测值。本质上，我们会尝试拟合 Y = B0 + B1 * X1 + B2 * X2 + ... + B10 * X10，其中每个 X 代表第 2 部分中随机 selected 之一。

4) 将R2值记录在向量中

5) 重复步骤 1-4 50,000 次

6) 找出记录的 50,000 个中的最大 R2

Answer 1

这是一个似乎可以解决您的问题的替代代码。

ns <- 30

rvals <- replicate(50000, {
  y <- rnorm(ns)
  xvals <- replicate(1000, rnorm(ns))
  selecteds <- xvals[,sample(1:1000, 10)]
  df <- data.frame(y = ys, selecteds)
  summary(lm(paste("y ~", paste0("X", 1:10, collapse = "+")), data = df))$r.squared
})

我对集群不是很有经验，但以下是您的代码可能太慢的几个原因：

您嵌套了 foor 循环来创建 X，我使用了 replicate，这可能比使用列表快一些。
您生成的是一个空列表，X，这非常糟糕。 (勾选The R inferno - Circle 2)
您取消列出了几个列表元素只是为了使它们成为 1 列矩阵，然后将它们全部绑定并最后命名列。尽管这些步骤似乎是必要的，但我认为一次一个地执行这些步骤可能很慢。例如，colnames 会自动设置为 X1:X10.
没有必要使用 attach 并且可能会减慢速度。
如果您 open/close 集群太多，会消耗大量处理，并且会比非并行更慢。不过好像不是这样。

最后一点，请确保我和你做的一样，因为这个问题对我来说仍然有点令人困惑。

使用并行编程在 R 中编程 50000 个回归

Programming 50000 regressions in R using parallel programming

parallel-processing

performance

r