Question

Chaque fois que je lance des simulations à grande échelle de monte carlo dans S-Plus, je finis toujours par me laisser pousser la barbe en l’attendant.

Quelles sont les meilleures astuces pour exécuter des simulations de monte carlo en R? De bons exemples de processus en cours de manière distribuée?

Était-ce utile?

La solution

  • L'utilisation de plusieurs cœurs / machines devrait être simple si vous utilisez uniquement des réplications indépendantes parallèles , tout en étant conscient des faiblesses courantes des générateurs de nombres aléatoires (par exemple, si vous utilisez l'heure actuelle comme valeur de départ). , engendrer de nombreux processus avec un RNG pour chacun pourrait produire des nombres aléatoires corrélés, ce qui conduit à des résultats non valides - voir par exemple cet article )

  • Vous pouvez utiliser la réduction de la variance pour réduire le nombre de réplications requises , c’est-à-dire réduire la taille de l’échantillon requis. Des techniques plus avancées de réduction de la variance peuvent être trouvées dans de nombreux manuels, par exemple. celui-ci .

Autres conseils

Préalignez vos vecteurs!

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

L’échantillonnage de l’hypercube latin s’applique facilement et a une influence majeure sur les résultats. Fondamentalement, vous prenez un échantillon d'hypercube latin à partir d'une distribution uniforme (par exemple, en utilisant random LHS () dans le paquetage lhs) et vous le transformez en votre distribution souhaitée en utilisant, par exemple, qnorm (uniformsample).

Je sais que ce fil est vraiment vieux, mais si quelqu'un y tombe par hasard et cherche une méthode encore plus rapide, je pense que ce qui suit fonctionne:

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
Licencié sous: CC-BY-SA avec attribution
Non affilié à StackOverflow
scroll top