Random Number Generation and Sampling Methods
Peter Occil - 22/Jun/2020
Peter Occil - 22/Jun/2020
[SHOWTOGROUPS=4,20]
Has many ways applications can sample from a random number generator, with pseudocode in many cases.
This page discusses many ways applications can generate and sample random content using an underlying random number generator (RNG), often with pseudocode. Those methods include ways to generate uniform random numbers from an underlying RNG; randomized content and conditions; and non-uniform random numbers, including weighted choice, the Poisson distribution, and other probability distributions. Sample Python code that implements many of the methods in this document is available.
This page discusses many ways applications can generate and sample random content using an underlying random number generator (RNG), often with pseudocode. Those methods include—
All the random number methods presented on this page are ultimately based on an underlying RNG; however, the methods make no assumptions on that RNG's implementation (e.g., whether that RNG uses only its input and its state to produce numbers) or on that RNG's statistical quality or predictability.
In general, this document does not cover:
In this document, a random number generator (RNG) means software and/or hardware that seeks to generate numbers with the property that each possible outcome is as likely as any other without influence by anything else(1).
The pseudocode conventions apply to this document.
This document uses the following notation for intervals:
Uniform Random Integers
This section describes how an underlying RNG can be used to generate independent uniform random integers. This section describes four methods: RNDINT, RNDINTEXC, RNDINTRANGE, RNDINTEXCRANGE. Of these, RNDINT, described next, can serve as the basis for the remaining methods.
RNDINT: Random Integers in [0, N]
In this document, RNDINT(maxInclusive) is the core method for using an underlying RNG to generate independent uniform random integers in the interval [0, maxInclusive].(3). The following are some ways to implement RNDINT:
[/SHOWTOGROUPS]
Has many ways applications can sample from a random number generator, with pseudocode in many cases.
This page discusses many ways applications can generate and sample random content using an underlying random number generator (RNG), often with pseudocode. Those methods include ways to generate uniform random numbers from an underlying RNG; randomized content and conditions; and non-uniform random numbers, including weighted choice, the Poisson distribution, and other probability distributions. Sample Python code that implements many of the methods in this document is available.
This page discusses many ways applications can generate and sample random content using an underlying random number generator (RNG), often with pseudocode. Those methods include—
- ways to generate uniform random numbers from an underlying RNG (such as the core method, RNDINT(N)),
- ways to generate randomized content and conditions, such as true/false conditions, shuffling, and sampling unique items from a list, and
- generating non-uniform random numbers, including weighted choice, the Poisson distribution, and other probability distributions.
All the random number methods presented on this page are ultimately based on an underlying RNG; however, the methods make no assumptions on that RNG's implementation (e.g., whether that RNG uses only its input and its state to produce numbers) or on that RNG's statistical quality or predictability.
In general, this document does not cover:
- How to choose an underlying RNG for a particular application, including in terms of security, performance, and quality. I have written more on RNG recommendations in another document.
- How to generate random security parameters such as encryption keys.
- Randomness extraction (also known as unbiasing, deskewing, or whitening), such as hash functions and von Neumann unbiasing.
- Variance reduction techniques such as importance sampling or common random numbers.
In this document, a random number generator (RNG) means software and/or hardware that seeks to generate numbers with the property that each possible outcome is as likely as any other without influence by anything else(1).
The pseudocode conventions apply to this document.
This document uses the following notation for intervals:
- [a, b) means "a or greater, but less than b".
- (a, b) means "greater than a, but less than b".
- (a, b] means "greater than a and less than or equal to b".
- [a, b] means "a or greater and b or less".
- log1p(x) is equivalent to ln(1 + x) and is a robust alternative where x is a floating-point number (Pedersen 2018)(2).
- MakeRatio(n, d) creates a rational number with the given numerator n and denominator d.
- Sum(list) calculates the sum of the numbers in the given list. Note, however, that summing floating-point numbers naïvely can introduce round-off errors. Kahan summation along with parallel reductions can be a more robust way than the naïve approach to compute the sum of three or more floating-point numbers.
Uniform Random Integers
This section describes how an underlying RNG can be used to generate independent uniform random integers. This section describes four methods: RNDINT, RNDINTEXC, RNDINTRANGE, RNDINTEXCRANGE. Of these, RNDINT, described next, can serve as the basis for the remaining methods.
RNDINT: Random Integers in [0, N]
In this document, RNDINT(maxInclusive) is the core method for using an underlying RNG to generate independent uniform random integers in the interval [0, maxInclusive].(3). The following are some ways to implement RNDINT:
- Rejection sampling, which roughly means: sample in a bigger range until a sampled number fits the smaller range. This method is unbiased but has a variable running time.
- Reduction method. Generate bignumber, an N-bit random integer with many more bits than maxInclusive + 1 has, then find—
- rem(bignumber, maxInclusive + 1) (modulo reduction), or
- (bignumber * (maxInclusive + 1)) >> N (see (Lemire 2016)(4)).
- Either method's running time is theoretically constant, but can introduce a so-called modulo bias (some numbers are slightly more likely to be chosen than others), which, however, gets smaller the more bits bignumber has.
| If the underlying RNG produces: | Then RNG() is: | And MODULUS is: |
|---|---|---|
| 32-bit non-negative integers. | The underlying RNG. | 232. |
| 64-bit non-negative integers. | The underlying RNG. | 264. |
| Integers in the interval [0, n). | The underlying RNG. | n. |
| Numbers in the interval [0, 1) known to be evenly spaced by a number p (e.g., dSFMT). | The underlying RNG, except with its outputs multiplied by p. | 1/p. |
| Numbers in the interval [0, 1), where numbers in [0, 0.5) and those in [0.5, 1) are known to occur with equal probability (e.g., Java's Math.random()). | An RNG that outputs 0 if the underlying RNG outputs a number less than 0.5, or 1 otherwise. | 2. |
| Numbers not specified above. | A new RNG formed by writing the underlying RNG's outputs to a stream of memory units (such as 8-bit bytes) and using a randomness extraction technique to transform that stream to n-bit integers. | 2n. |
Code:
METHOD RndIntHelperNonPowerOfTwo(maxInclusive)
if maxInclusive <= MODULUS - 1:
// NOTE: If the programming language implements
// division with two integers by discarding the result's
// fractional part, the division can be used as is without
// using a "floor" function.
nPlusOne = maxInclusive + 1
maxexc = floor((MODULUS - 1) / nPlusOne) * nPlusOne
while true
ret = RNG()
if ret < nPlusOne: return ret
if ret < maxexc: return rem(ret, nPlusOne)
end
else
cx = floor(maxInclusive / MODULUS) + 1
while true
ret = cx * RNG()
// NOTE: The addition operation below should
// check for integer overflow and should reject the
// number if overflow would result.
ret = ret + RNDINT(cx - 1)
if ret <= maxInclusive: return ret
end
end
END METHOD
METHOD RndIntHelperPowerOfTwo(maxInclusive)
// NOTE: Finds the number of bits minus 1 needed
// to represent MODULUS (in other words, the number
// of random bits returned by RNG() ). This will
// be a constant here, though.
modBits = ln(MODULUS)/ln(2)
// Fast dice roller algorithm.
x = 1
y = 0
nextBit = modBits
rngv = 0
while true
if nextBit >= modBits
nextBit = 0
rngv = RNG()
end
x = x * 2
y = y * 2 + rem(rngv, 2)
rngv = floor(rngv / 2)
nextBit = nextBit + 1
if x > maxInclusive
if y <= maxInclusive: return y
x = x - maxInclusive - 1
y = y - maxInclusive - 1
end
end
END METHOD
METHOD RNDINT(maxInclusive)
// maxInclusive must be 0 or greater
if maxInclusive < 0: return error
if maxInclusive == 0: return 0
if maxInclusive == MODULUS - 1: return RNG()
// NOTE: Finds the number of bits minus 1 needed
// to represent MODULUS (if it's a power of 2).
// This will be a constant here, though. If modBits
// is an integer, MODULUS is a power of 2, which
// is checked below.
modBits=ln(MODULUS)/ln(2)
if floor(modBits) == modBits
return RndIntHelperPowerOfTwo(maxInclusive)
else
return RndIntHelperNonPowerOfTwo(maxInclusive)
end
END METHOD