RandomSample

This class generates random numbers from a discrete distribution and is optimised for sampling a small population many times.

The discrete statistical distribution can only take discrete values. It has probability function

(1)

cumulative distribution function

(2)

and population mean

(3)

The algorithm used firstly calculates the CDF of the given probability array, and then each time a new random number is requested, it generates an uniform random deviate between 0 and 1 then finds its corresponding abscissa value, as suggested in the graph below. Afterwards it uses binary search to detect the interval in which that value is situated, which in return gives the actual discrete random number. For instance if the ordinate would be approximately equal to as in the graph, the corresponding abscissa would be around , which corresponds to an index in the values array that finally gives the random number.

The given example displays 20 random floating point numbers from a discrete distribution. It uses two different generators to achieve this. The first generator uses a particular value to initialize the seed, while the second one is using the system timer. Notice that it was necessary to divide the timer by the <em> MERSENNEDIV </em> value in order to keep the seed in the interval from 0 to 1. Since the seed of the first generator is never changed, the first 10 numbers will always remain the same. However since the second generator is initialized via the system timer, the next 10 numbers will obviously vary with each execution of the program, Also notice the declaration of the probability array P and value array V. They are used as constructor arguments when declaring the generators to specify the probabilities P[i] associated to each value V[i] in the values table. Always make sure that the components of P sum up to 1. If there is no values table, the application assumes the following array

V[n] = {0, 1, 2, ..., n - 1}

where n is the size of the probabilities array.

In the example output you will find 10 numbers corresponding to the output of the first generator.

Speed

The average running time for generating 100,000,000 random numbers using this class on a 750MHz microprocessor is 16 seconds.

#include <iostream>
#include <time.h>
#include <codecogs/stats/dists/discrete/discrete/randomsample.h>
using namespace std;

int main()
{
    double P[4] = {   0.4,   0.1,   0.0,   0.5 };
    double V[4] = { 1.234, 2.345, 3.456, 4.567 };

    Stats::Dists::Discrete::Discrete::RandomSample<double> A(4, P, V, false, 0.3416);
    Stats::Dists::Discrete::Discrete::RandomSample<double> B(4, P, false, time(0) / MERSENNEDIV);

    for (int i = 0; i < 10; ++i)
        cout << A.genReal() << " ";
    cout << endl;

    for (int i = 0; i < 10; ++i)
        cout << B.genReal() << " ";
    cout << endl;

    return 0;
}

4.567 1.234 4.567 1.234 4.567 4.567 1.234 2.345 1.234 1.234

References

• MathWorld, http://mathworld.wolfram.com/DiscreteDistribution.html
• The Newran03 random number generator library of Robert Davies, http://www.robertnz.net/nr03doc.htm