Va de muestrear los números que tienen asignadas probabilidades
. Una manera muy impráctica (en R, basta usar
sample
) y nada intuitiva de hacerlo es recurriendo a la distribución de Gumbel:
library(evd) pes <- runif(5) pes <- pes / sum(pes) gammas <- log(pes) + 2 x <- rgumbel(length(pes)) muestra <- which.max(gammas + x)
O, en masa, aplicando
get_samples <- function(n){ replicate(n, { x <- rgumbel(length(pes)) which.max(gammas + x) }) }
El seudocódigo está extraído de la Wikipedia y el motivo por el que la cosa funciona en lugar de no funcionar, que es la parte bonita del asunto, está explicado aquí.