Community
cancel
Showing results for 
Search instead for 
Did you mean: 
junuylia
Beginner
96 Views

How do I use MKL to sample from discrete probability distribution

For example,

I want to sample from a distribution with

x=1, 2, 3
p(x)=1/2, 1/3, 1/6

how do I do it with mkl library?

Thank you very much.

Jun
0 Kudos
4 Replies
Andrey_N_Intel
Employee
96 Views

Hello Jun,

You may want to follow the algorithm described below.

1. Split the interval I=[0,1) into three sub-intervals:
I1 = [0, 1/2)
I2 = [1/2, 5/6), 5/6=1/2+1/3
I3 = [5/6,1)

2. Generate uniformly distributed on the interval I random numbers u(1),...,u(n) using VSL RNG routine
v[d|s]RngUniform.

3. For each number u(i) find the interval Ij which contains u(i) and set y(i) to j

4. Output of the algorithm is sequence y(1),...,y(n)

Please, let me know if this answers your question.

Best,
Andrey

junuylia
Beginner
96 Views

Thank you Andrey.

I guess it will work. I'll make some follow up notes after I try this algorithm. Thanks!

Jun
mecej4
Black Belt
96 Views

No need for guesses! For a continuous random variable X with range (a,b) with density function f(x), define the cumulative probability function F(x) = \int_a^x f(t) dt, where \int_a^x stands for "integral from t=a to t=x". Then, the random variable U defined by U = F(X) / F(b) is uniformly distributed with range (0,1).

For discrete random distributions, replace integrals by sums in the definition of F(x). The algorithm that Andrey described for you essentially solves F(X) = U, where U is obtained from a uniform random number generator. The solution X is the result that you sought.

To supplement the mathematical argument that I stated, you can easily run a Monte-Carlo simulation, generating a few million random numbers with your distribution. Count the numbers of times that 1, 2 and 3 occur in the sequence, and divide by the number of random numbers generated. The resulting fractions should be close to 1/2, 1/3 and 1/6.
junuylia
Beginner
96 Views

I just saw this, which explains a lot. I've used the previous method, and it works.

This method is new for me, and I'll try it later. Thanks a lot!
Reply