Generate Random(a, b) making calls to Random(0, 1)

Question

There is known Random(0,1) function, it is a uniformed random function, which means, it will give 0 or 1, with probability 50%. Implement Random(a, b) that only makes calls to Random(0,1)

What I though so far is, put the range a-b in a 0 based array, then I have index 0, 1, 2...b-a.

then call the RANDOM(0,1) b-a times, sum the results as generated idx. and return the element.

However since there is no answer in the book, I don't know if this way is correct or the best. How to prove that the probability of returning each element is exactly same and is 1/(b-a+1) ?

And what is the right/better way to do this?

Answer 1

If you have a RNG that returns {0, 1} with equal probability, you can easily create a RNG that returns numbers {0, 2^n} with equal probability.

To do this you just use your original RNG n times and get a binary number like 0010110111. Each of the numbers are (from 0 to 2^n) are equally likely.

Now it is easy to get a RNG from a to b, where b - a = 2^n. You just create a previous RNG and add a to it.

Now the last question is what should you do if b-a is not 2^n?

---

Good thing that you have to do almost nothing. Relying on rejection sampling technique. It tells you that if you have a big set and have a RNG over that set and need to select an element from a subset of this set, you can just keep selecting an element from a bigger set and discarding them till they exist in your subset.

So all you do, is find b-a and find the first n such that b-a <= 2^n. Then using rejection sampling till you picked an element smaller b-a. Than you just add a.

Answer 2

Divide and conquer could help us in generating a random number in range [a,b] using random(0,1). The idea is

1. if a is equal to b, then random number is a
2. Find mid of the range [a,b]
3. Generate random(0,1)
4. If above is 0, return a random number in range [a,mid] using recursion
5. else return a random number in range [mid+1, b] using recursion

The working 'C' code is as follows.

int random(int a, int b)
{ 
    if(a == b)
        return a;

    int c = RANDOM(0,1); // Returns 0 or 1 with probability 0.5
    int mid = a + (b-a)/2;

    if(c == 0)
       return random(a, mid);
    else
       return random(mid + 1, b);
}

Answer 3

If your RANDOM(0, 1) returns either 0 or 1, each with probability 0.5 then you can generate bits until you have enough to represent the number (b-a+1) in binary. This gives you a random number in a slightly too large range: you can test and repeat if it fails. Something like this (in Python).

def rand_pow2(bit_count):
    """Return a random number with the given number of bits."""
    result = 0
    for i in xrange(bit_count):
        result = 2 * result + RANDOM(0, 1)
    return result

def random_range(a, b):
    """Return a random integer in the closed interval [a, b]."""
    bit_count = math.ceil(math.log2(b - a + 1))
    while True:
        r = rand_pow2(bit_count)
        if a + r <= b:
            return a + r

Answer 4

Your inclination to put the range from 0 to a-b first is correct. However, you cannot do it as you stated. This question asks exactly how to do that, and the answer utilizes unique factorization. Write m=a-b in base 2, keeping track of the largest needed exponent, say e. Then, find the biggest multiple of m that is smaller than 2^e, call it k. Finally, generate e numbers with RANDOM(0,1), take them as the base 2 expansion of some number x, if x < k*m, return x, otherwise try again. The program looks something like this (simple case when m<2^2):

int RANDOM(0,m) {

    // find largest power of n needed to write m in base 2
    int e=0;
    while (m > 2^e) {
        ++e;
    }

    // find largest multiple of m less than 2^e
    int k=1;
    while (k*m < 2^2) {
        ++k
    }
    --k; // we went one too far

    while (1) {
        // generate a random number in base 2
        int x = 0;
        for (int i=0; i<e; ++i) {
            x = x*2 + RANDOM(0,1); 
        }
        // if x isn't too large, return it x modulo m
        if (x < m*k) 
            return (x % m);
    }
}

Now you can simply add a to the result to get uniformly distributed numbers between a and b.

Answer 5

When you sum random numbers, the result is not longer evenly distributed - it looks like a Gaussian function. Look up "law of large numbers" or read any probability book / article. Just like flipping coins 100 times is highly highly unlikely to give 100 heads. It's likely to give close to 50 heads and 50 tails.