Какой самый лучший трюк, чтобы ускорить моделирование Монте-Карло? [закрытый]

Всякий раз, когда я запускаю крупномасштабные симуляции Монте-Карло в S-Plus, я всегда заканчиваю тем, что отращиваю бороду, пока жду его завершения.

Каковы лучшие приемы для запуска моделирования Монте-Карло в R? Есть ли хорошие примеры выполнения процессов распределенным способом?

4 5
distributed-computing

4 ответа:

  • Использование нескольких ядер / машин должно быть простым, если вы просто используете параллельные независимые репликации , но помните об общих недостатках генераторов случайных чисел (например, при использовании текущего времени в качестве затравки, порождение многих процессов с одним ГСЧ для каждого может привести к коррелированным случайным числам, что приводит к недействительным результатам - см., например, эту статью)

  • Возможно, вы захотите использовать уменьшение дисперсии для уменьшения числа требуемая репликация , то есть уменьшение размера требуемого образца. Более продвинутые методы уменьшения дисперсии можно найти во многих учебниках, например, в этом .

10

Заранее распределите векторы!

> 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 
>
5

Латинская выборка гиперкубов легко применяется и оказывает значительное влияние на результаты. В основном вы берете образец латинского гиперкуба из равномерного распределения (например, используя randomLHS () в пакете lhs) и преобразуете его в желаемое распределение, используя, например, qnorm (uniformsample).

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
3