r - R中使用降雪的并行蒙特卡罗模拟

Question

我尝试比较多达数千个估计的 Beta 分布。每个 beta 分布都以两个形状参数 alpha 和 beta 为特征。我现在为每个分布抽取 100,000 个样本。作为最终结果，我想获得每个抽样中概率最高的分布顺序。我的第一种方法是使用 lapply 生成 N * NDRAWS 数值矩阵，当 N 超过 10,000 时，它消耗了太多内存。 (10,000 * 100,000 * 8 字节)

所以我决定使用顺序方法对每个抽奖进行排序，然后对所有抽奖的顺序求和并获得最终顺序，如下例所示:

set.seed(12345)
N=100
NDRAWS=100000
df <- data.frame(alpha=sample(1:20, N, replace=T), beta=sample(1:200, N, replace=T))  

vec    <- vector(mode = "integer", length = N )

for(i in 1:NDRAWS){
  # order probabilities after a single draw for every theta
  pos <- order(rbeta(N, shape1=df$alpha, shape2=df$beta) )

  # sum up winning positions for every theta
  vec[pos] <- vec[pos] + 1:N
}

# order thetas
ord <- order(-vec)

df[ord,]

这只消耗 N * 4 字节的内存，因为没有巨大的矩阵，只有一个长度为 N 的向量。我现在的问题是，如何使用降雪(或任何其他多核包)通过采取利用我的 4 个 CPU 核心，而不是只使用一个核心？？？

# parallelize using snowfall pckg
library(snowfall)
sfInit( parallel=TRUE, cpus=4, type="SOCK")
sfLapply( 1:NDRAWS, function(x) ?????? )
sfStop()

感谢任何帮助!

Answer 1

这可以用与并行化随机森林或 Bootstrap 相同的方式进行并行化。您只需对每个工作人员执行顺序代码，但每个工作人员使用较少的迭代次数。这比将 for 循环的每次迭代拆分为一个单独的并行任务要高效得多。

这是您的完整示例，已转换为使用带有 doParallel 后端的 foreach 包:

set.seed(12345)
N=100
NDRAWS=100000
df <- data.frame(alpha=sample(1:20, N, replace=T),
                 beta=sample(1:200, N, replace=T))
library(doParallel)
nworkers <- detectCores()
cl <- makePSOCKcluster(nworkers)
clusterSetRNGStream(cl, c(1,2,3,4,5,6,7))
registerDoParallel(cl)

vec <- foreach(ndraws=rep(ceiling(NDRAWS/nworkers), nworkers),
               .combine='+') %dopar% {
  v <- integer(N)
  for(i in 1:ndraws) {
    pos <- order(rbeta(N, shape1=df$alpha, shape2=df$beta) )
    v[pos] <- v[pos] + 1:N
  }
  v
}
ord <- order(-vec)
df[ord,]

请注意，这给出了与顺序版本不同的结果，因为 worker 生成了不同的随机数。我使用了并行包提供的并行随机数支持，因为这是很好的做法。