# Algorithm in RandomChoice[{w1,w2,…}->{e1,e2,…},n]?

I wonder if someone happens to know the actual algorithm implemented in the particular overload of RandomChoice that handles explicitly presented discrete distributions, namely

RandomChoice[{w1, w2, ...}->{e1, e1, ...}]
gives a pseudorandom choice weighted by the w1.


It's blazing fast, as illustrated in my test notebook here: https://www.dropbox.com/s/rwm80fut60v926b/FastNonUniformPseudoRandoms002.cdf

(btw, such things are really useful in, say, particle filters).

I know of at least four algorithms for solving this problem. A greedy-space array of outcomes has O(1) time but O(S) space where S is the sum of the frequency counts (and only works when the weights are integers or can be scaled up to integers). Two algorithms explicitly invert the CDF, either by linear search -- O(N) time, O(N) space, or binary search -- O(log N) time, O(N) space. The best I know is Walker's Alias Method http://en.wikipedia.org/wiki/Alias_method

which is O(N) space, O(1) time. I implemented this in my notebook and it beats the pants off the others (as expected) but is still 30 times slower than the built-in (though still just O(1)).

I wonder if the built-in is just an optimized Walker's or is it something better?

-
Regarding the performance: it's probably implemented in C, not in Mathematica, which might be the reason why it's faster than your Mathematica-based implementation. – Szabolcs Mar 4 '14 at 2:53
You might find (if you've not already read it) wolfram.com/learningcenter/tutorialcollection/… interesting. Slightly dated but relevant, covers how MMA generally uses inverse CDF for such things until it deems that expensive, then switching to table-lookup or direct generation. As Szabolcs said, low-level internal code in C will always beat an interpreted interpretation... – ciao Mar 4 '14 at 3:05
@rasher +1 on your comment for "interpreted interpretation" ... clearly clear! – Dr. belisarius Mar 4 '14 at 3:17
@rasher Very interesting, it actually answered a related question I had today. Interesting enough that it's worth putting in a more prominent answer IMO. – Szabolcs Mar 4 '14 at 3:23