3

我尝试比较多达数千个估计的 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()

任何帮助表示赞赏!

4

2 回答 2

1

这可以像并行化随机森林或引导程序一样并行化。您只需在每个工作人员上执行顺序代码,但每个工作人员使用较少数量的迭代。这比将 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,]

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

于 2013-10-30T22:35:02.557 回答
0

嗯,功能就在那里。我不确定每次迭代会返回什么。

也许试试这个?

myFunc <- function(xx, N) {
  pos <- order(rbeta(N, shape1=df$alpha, shape2=df$beta) )
  vec[pos] + 1:N
}

使用 doParallel 将允许您添加结果:

require(doParallel)
registerDoParallel(cores=4)
foreach(i=1:NDRAWS, .combine='+') %dopar% myFunc(i, N)
于 2013-10-29T20:08:34.280 回答