7

每当我在 S-Plus 中运行大规模的蒙特卡罗模拟时,我总是在等待它完成时长出胡须。

在 R 中运行蒙特卡罗模拟的最佳技巧是什么?任何以分布式方式运行进程的好例子?

4

4 回答 4

11
  • 如果您只是使用并行独立复制,则使用多个核心/机器应该很简单,但请注意随机数生成器的常见缺陷(例如,如果使用当前时间作为种子,则为每个进程生成一个 RNG 可能会产生相关随机数字,这会导致无效的结果- 参见例如本文

  • 您可能希望使用方差减少减少所需的重复次数,即缩小所需样本的大小。更高级的方差减少技术可以在许多教科书中找到,例如在这本教科书中。

于 2009-09-10T06:06:16.393 回答
5

预先分配你的载体!

> nsims <- 10000
> n <- 100
> 
> system.time({
     res <- NULL
     for (i in 1:nsims) {
         res <- c(res,mean(rnorm(n)))
     }
 })
   user  system elapsed 
  0.761   0.015   0.783 
> 
> system.time({
     res <- rep(NA, nsims)
     for (i in 1:nsims) {
         res[i] <- mean(rnorm(n))
     }
 })
   user  system elapsed 
  0.485   0.001   0.488 
> 
于 2009-09-10T18:58:25.933 回答
3

拉丁超立方抽样很容易应用,并且对结果有重大影响。基本上,您从均匀分布中获取一个拉丁超立方体样本(例如,使用包 lhs 中的 randomLHS())并使用例如 qnorm(uniformsample) 将其转换为您想要的分布。

于 2009-09-10T07:13:53.447 回答
3

我知道这个线程真的很旧,但如果有人偶然发现它并正在寻找一种更快的方法,我认为以下方法有效:

library(data.table)
library(microbenchmark)

nsims <- 10000
n <- 100

# Answer from @Eduardo_Leoni:
preallocate<-function(nsims, n) {
  res <- rep(NA, nsims)
  for (i in 1:nsims) {
    res[i] <- mean(rnorm(n))
  }
  return(res)
}

# Answer using data.table:
datatable<-function(nsims,n) {
  dt <- data.table(i=1:nsims)[,list(res=mean(rnorm(1:n))),by=i]
  return(dt)
}

# Timing benchmark:
microbenchmark(preallocate(nsims,n), datatable(nsims,n), times=100)
#Unit: milliseconds
#                  expr      min       lq   median       uq      max neval
# preallocate(nsims, n) 428.4022 432.3249 434.2910 436.4806 489.2061   100
#   datatable(nsims, n) 238.9006 242.3517 244.1229 246.5998 303.6133   100
于 2013-09-02T13:30:38.133 回答