Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

user defined discrete random distribution?

55 views
Skip to first unread message

Studer Christoph

unread,
May 4, 2004, 11:43:35 AM5/4/04
to
is there any function in matlab, where i can define for example a
vector p=( 2 4 1 3 2 ) with the discrete probability density
p/sum(p), such that the result of the function returns a random
number (index of vector p is the random result) according to the
defined distribution? this would be very usefull!

thanks in advance

christoph

Ken Davis

unread,
May 4, 2004, 11:49:36 AM5/4/04
to
"Studer Christoph" <m...@rune.ch> wrote in message
news:eedb4...@webx.raydaftYaTP...

Without thinking through the details, it seems like you could use hist or
histc to select what to return after letting X, the bin boundaries be
defined by cumsum(p/sum(p)) and feeding it uniform variates generated by
rand.


Tom Lane

unread,
May 4, 2004, 1:11:45 PM5/4/04
to
"Studer Christoph" <m...@rune.ch> wrote in message
news:eedb4...@webx.raydaftYaTP...

Christoph, if you have the Statistics Toolbox, try "help randsample" and
look at the example at the end. This example returns random values from the
letters 'ACGT', but you could choose instead to return random values from
1:4.

-- Tom


Peter Boettcher

unread,
May 5, 2004, 10:24:57 AM5/5/04
to
"Studer Christoph" <m...@rune.ch> writes:

So the usual way to do this is to construct a mapping function from
a uniform rv to the rv with the distribution you want. This mapping
function is the inverse of the CDF of the distribution.

The CDF for discrete rvs is the cumsum of the PDF:

my_cdf = cumsum(p)

% my_cdf = [2 6 7 10 12]

Now if you take a bunch of uniform samples from [0 to 12), and take
the index of the next greatest element of my_cdf, you'll have your
distribution. You could do it manually (using find and a loop), but a
look-up table is better:

clear q;
q(my_cdf) = 1;
my_lut = cumsum([1 q]);

Now samples may be generated using:

drawn_samples = my_lut(floor(my_cdf(end)*rand(1,N)) + 1);

--
Peter Boettcher <boet...@ll.mit.edu>
MIT Lincoln Laboratory
MATLAB FAQ: http://www.mit.edu/~pwb/cssm/

SeBy

unread,
May 5, 2004, 11:13:27 AM5/5/04
to

A fast one line way :

[toto , draw_samples] = histc(rand(1 , N), [0 cumsum(p)./sum(p)])

Sébastien

0 new messages