Pseudorandom Number Generation

Question

Why do I need pseudorandom number generators?

Answer

Pseudorandom number generators play a central role in a computer simulation, a flexible and powerful tool that can be used in different ways, e.g., to get a better understanding of real-world data, to appreciate the functioning of complex systems, or for scenario analysis. The students who want to familiarize themselves with the role of computer simulation in understanding economics and social formations may want to refer to i) Axelrod, R. The complexity of cooperation, Princeton University Press, 1997; ii) Schelling, T. C. Micromotives and macrobehavior. WW Norton & Company, 1978; iii) Epstein, J. M. & Axtell, R. Growing Artificial Societies—Social Science from the Bottom Up. Artif Life 3, 237–242, 1997.

Question

What is the gamut of pseudorandom numbers available in NumPy?

Answer

NumPy has three families of pseudorandom number generators:

Simple random data
Permutations
Distributions

The tables below provide a summary of the available generators.

Simple Random Data

Routine	Synopsis
`integers(low[, high, size, dtype, endpoint])`	Return random integers from low (inclusive) to high (exclusive), or if endpoint=True, low (inclusive) to high (inclusive)
`random([size, dtype, out])`	Return random floats in the half-open interval [0.0, 1.0)
`choice(a[, size, replace, p, axis, shuffle])`	Generates a random sample from a given array
`bytes(length)`	Return random bytes

Note: the statements included in the ‘Routine’ column assume numpy.random.Generator has been imported.

Permutations

Routine	Synopsis
`shuffle(x[, axis])`	Modify an array or sequence in place by shuffling its contents
`permutation(x[, axis])`	Randomly permute a sequence, or return a permuted range
`permuted(x[, axis, out])`	Randomly permute x along axis

Note: the statements included in the ‘Routine’ column assume numpy.random.Generator has been imported.

Distributions

Routine	Synopsis
`beta(a, b[, size])`	Draw samples from a Beta distribution
`binomial(n, p[, size])`	Draw samples from a binomial distribution
`chisquare(df[, size])`	Draw samples from a chi-square distribution
`dirichlet(alpha[, size])`	Draw samples from the Dirichlet distribution
`exponential([scale, size])`	Draw samples from an exponential distribution
`f(dfnum, dfden[, size])`	Draw samples from an F distribution
`gamma(shape[, scale, size])`	Draw samples from a Gamma distribution
`geometric(p[, size])`	Draw samples from the geometric distribution
`gumbel([loc, scale, size])`	Draw samples from a Gumbel distribution
`hypergeometric(ngood, nbad, nsample[, size])`	Draw samples from a Hypergeometric distribution
`laplace([loc, scale, size])`	Draw samples from the Laplace or double exponential distribution with specified location (or mean) and scale (decay)
`logistic([loc, scale, size])`	Draw samples from a logistic distribution
`lognormal([mean, sigma, size])`	Draw samples from a log-normal distribution
`logseries(p[, size])`	Draw samples from a logarithmic series distribution
`multinomial(n, pvals[, size])`	Draw samples from a multinomial distribution
`multivariate_hypergeometric(colors, nsample)`	Generate variates from a multivariate hypergeometric distribution
`multivariate_normal(mean, cov[, size, ..])`	Draw random samples from a multivariate normal distribution
`negative_binomial(n, p[, size])`	Draw samples from a negative binomial distribution
`noncentral_chisquare(df, nonc[, size])`	Draw samples from a noncentral chi-square distribution
`noncentral_f(dfnum, dfden, nonc[, size])`	Draw samples from the noncentral F distribution
`normal([loc, scale, size])`	Draw random samples from a normal (Gaussian) distribution
`pareto(a[, size])`	Draw samples from a Pareto II or Lomax distribution with a specified shape
`poisson([lam, size])`	Draw samples from a Poisson distribution
`power(a[, size])`	Draws samples in [0, 1] from a power distribution with positive exponent a - 1
`rayleigh([scale, size])`	Draw samples from a Rayleigh distribution
`standard_cauchy([size])`	Draw samples from a standard Cauchy distribution with mode = 0
`standard_exponential([size, dtype, method, out])`	Draw samples from the standard exponential distribution
`standard_gamma(shape[, size, dtype, out])`	Draw samples from a standard Gamma distribution
`standard_normal([size, dtype, out])`	Draw samples from a standard Normal distribution (mean=0, stdev=1)
`standard_t(df[, size])`	Draw samples from a standard Student’s t distribution with df degrees of freedom
`triangular(left, mode, right[, size])`	Draw samples from the triangular distribution over the interval [left, right]
`uniform([low, high, size])`	Draw samples from a uniform distribution
`vonmises(mu, kappa[, size])`	Draw samples from a von Mises distribution
`wald(mean, scale[, size])`	Draw samples from a Wald, or inverse Gaussian, distribution
`weibull(a[, size])`	Draw samples from a Weibull distribution
`zipf(a[, size])`	Draw samples from a Zipf distribution

Note: the statements included in the ‘Routine’ column assume numpy.random.Generator has been imported.

Permutation in NumPy

Question

Permutation in NumPy

Answer

Permuting means changing the order of the elements in an array. Array permutations can be done ‘in-place’ (see the code snippet below, line 13), so the original array is modified, or a copy is created (see line 17). Note the outcome of the np.random functions are not reproducible. Simply, running the same NumPy generator n times might return n different outcomes. To ensure the reproducibility of NumPy code, we must initialize a random generator instance, as shown in the lower section of the code snippet below. Specifically, in line 23, we assign the object rgn to the outcome of .default_rng. The parameter passed to the generator, known as ‘seed,’ is an arbitrary number object. If you reproduce lines 23-25, you will get the outcome I got, displayed in line 26.

# import numpy with the socially accepted alias 'np'
>>> import numpy as np

# shuffling/permuting an existing array
# --+ the array
>>> A = np.arange(10)
>>> A
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
# --+ .shuffle produces in-place changes
>>> np.random.shuffle(A)
>>> A
array([6, 9, 3, 8, 7, 2, 4, 1, 0, 5])
# --+ .permutation creates a copy
>>> np.random.permutation(A)
array([7, 6, 5, 2, 1, 8, 0, 4, 9, 3])
>>> A
array([6, 9, 3, 8, 7, 2, 4, 1, 0, 5])

# shuffling an existing array and ensuring reproducibility
# --+ the array
>>> B = [0, 1, 2, 3]
# --+ initialize a random generator instance
>>> rng = np.random.default_rng(12345)
>>> rng.shuffle(B)
>>> B
[2, 0, 3, 1]

Sampling

Question

Sampling what?

Answer

Random sampling is the process of selecting a random subset of elements from an array. NumPy allows us:

to sample the elements belonging to an existing array
to sample the elements included in a certain interval
to sample from a theoretical distribution

Lines 9 and 12 included in the code snippet below show how to sample a given number of elements (see parameter size) from an existing array. The code in line 12 differs from the code in line 9 because of the optional parameter replace that is set to False (default is True). In so doing, an element can be sampled once and once only.

In lines 22 and 27, we create two arrays of shape (10000,) from the random normal (see .random.normal) and Poisson (see .random.poisson) distribution respectively. In the interest of redundancy, .random.normal takes three mandatory arguments: the mean of the distribution (loc), the standard deviation of the distribution (scale), and the size of the array (size); .random.poisson takes two mandatory arguments: the expected value/variance of the distribution (lam) and the size of the distribution (size). In lines 32 - 50, we use Matplotlib to visualize the distribution of the two arrays.

# import numpy with the socially accepted alias 'np'
>>> import numpy as np

# the array
>>> A = np.arange(10)

# random sampling from an existing array
# --+ sampling with replacement (items can be drawn multiple times)
>>> np.random.choice(A, size=7)
array([3, 5, 7, 3, 1, 5, 2])
# --+ sampling without replacement (items can be drawn only once)
>>> np.random.choice(A, size=7, replace=False)
array([4, 2, 6, 3, 1, 0, 5])

# sampling from a range
# --+ legacy top-level function
>>> np.random.randint(0, 10, size=4)
array([1, 4, 1, 9])
# --+ or, on a generator instance --- 'integers' is a Generator method,
#     there is no np.random.integers
>>> rng = np.random.default_rng(42)
>>> rng.integers(0, 10, size=4)
array([0, 7, 6, 4])

# sampling from a theoretical distribution
# --+ 10,000 items from a normal distribution with a mean of 10 and
# standard deviation of 10
>>> N = np.random.normal(loc=10, scale=10, size=10000)
>>> N
array([-4.32898522,  7.35694695, 19.76873692, ..., 10.91842644,
       -2.23323711,  3.83215207])
# --+ 10,000 items from the Poisson distribution with lambda=10
>>> P = np.random.poisson(lam=10, size=10000)
>>> P
array([ 6, 13, 13, ...,  6,  8, 13])
# --+ visualize the two arrays
# ----+ Normal distribution
fig = plt.figure(figsize=(3, 3))
ax = fig.add_subplot(111)
ax.hist(N, color='blue', bins=50)
ax.set_ylabel('Count')
ax.set_xlabel('Value')
plt.title('Normal distribution data')
plt.grid(True,  ls="--")
plt.show()
# ----+ Poisson distribution
from collections import Counter
P_FR = Counter(P)
fig = plt.figure(figsize=(3, 3))
ax = fig.add_subplot(111)
ax.scatter(P_FR.keys(), P_FR.values(), color='blue')
ax.set_ylabel('Count')
ax.set_xlabel('Value')
plt.title('Poisson distribution data')
plt.grid(True, ls="--")
plt.show()