Just a few lines of code if you are willing to use numpy. Uses fact that any prob. distr. can be sampled by computing the cumulative distribution, drawing a random number from 0 to 1, and finding the x-value where that number is attained on the cumulative distribution. The searchsorted(..) function performs this search.
Python, 20 lines
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
from numpy import * from numpy.random import random def toss(c, sh): """Return random samples of cumulative vector (1-D numpy array) c The samples are placed in an array of shape sh. Each element of array returned is an integer from 0 through n-1, where n is the length of c.""" y = random(sh) #return array of uniform numbers (from 0 to 1) of shape sh x = searchsorted(c,y) return x #sample usage: p=array((0.1, 0.2, 0.6, 0.1)) #vector of probabilities, normalized to 1 c=cumsum(p) #cumulative probability vector print toss(c, (10,)) #want 10 samples in 1-D array
Bad style. I generally find the
from numpy import *
to be bad style. In particular, it would be helpful in a recipe to show exactly where everything comes from. If you want to be lazy later and import things with the star notation, that is fine, but here it just obfuscates things.