3.3 Randomness, Sampling, and Reproducibility

Randomness powers simulations, experiments, train/test splits, and sampling. In software, NumPy normally produces pseudorandom values: a deterministic algorithm expands an internal state into a sequence that behaves statistically like random data. This is useful because a controlled sequence can be replayed.

Reproducible randomness is not the same as “always use seed 42.” It is a record of the generator, its initial state, and the operations that consumed that state.

Create an explicit generator

Modern NumPy code should create a Generator with np.random.default_rng:

import numpy as np

rng = np.random.default_rng(2026)

The integer 2026 is a seed, an input used to initialize generator state. Creating a second generator with the same seed and making the same calls in the same order reproduces the sequence:

rng_a = np.random.default_rng(2026)
rng_b = np.random.default_rng(2026)

first = rng_a.integers(0, 10, size=5)
replay = rng_b.integers(0, 10, size=5)

assert np.array_equal(first, replay)

An extra draw changes every later position:

rng_c = np.random.default_rng(2026)
rng_c.integers(0, 10)              # consumes one value
shifted = rng_c.integers(0, 10, size=5)

The seed is the same, but the call sequence is not. Refactoring code, adding a diagnostic sample, or changing an array's shape can therefore change later results.

Passing the generator into a function makes this dependency visible:

def simulate_waits(rng, size):
    return rng.normal(loc=38, scale=7, size=size)

rng = np.random.default_rng(2026)
waits = simulate_waits(rng, size=1_000)

This approach avoids hidden reliance on the legacy module-level random state. It also makes tests easy: a test can provide a generator with a known seed.

Match the distribution to the question

A generator offers many distributions. Their parameters express different assumptions:

rng = np.random.default_rng(2026)

uniform_0_to_1 = rng.random(5)
dice = rng.integers(1, 7, size=20)  # high=7 is excluded
waits = rng.normal(loc=38, scale=7, size=1_000)

random draws uniformly from the half-open interval $[0, 1)$. integers(1, 7) draws whole numbers from 1 through 6. normal uses loc as the mean and scale as the standard deviation.

Choosing a convenient distribution does not make it a valid model. A normal distribution permits negative values and is symmetric, so it may be poor for strongly right-skewed waiting times. Explore real data, state modeling assumptions, and perform sensitivity checks.

For several independent simulation streams, do not invent nearby seeds such as 100, 101, and 102 and assume independence. A SeedSequence can spawn reproducibly managed child states:

root = np.random.SeedSequence(2026)
children = root.spawn(3)
generators = [np.random.default_rng(child) for child in children]

Store the root seed and the mapping between child streams and tasks. Parallel execution may finish in a different order, but each task can retain its own reproducible stream.

Sample with or without replacement

Generator.choice draws items from a population:

order_ids = np.array(["A104", "A105", "A106", "A107", "A108"])

sample = rng.choice(order_ids, size=3, replace=False)
bootstrap = rng.choice(order_ids, size=8, replace=True)

Without replacement, an item can appear at most once and the sample size cannot exceed the population size. With replacement, every draw returns its item before the next draw, so duplicates are possible and larger samples are allowed. Bootstrap methods deliberately use replacement; an audit sample of unique orders often does not.

Probabilities can be supplied with p, but they must be non-negative and sum to one:

weights = np.array([0.10, 0.15, 0.25, 0.20, 0.30])
weighted = rng.choice(order_ids, size=3, replace=False, p=weights)

Weights change the design and therefore the interpretation of estimates. Preserve them with the sample, and use appropriate survey or importance-weighting methods when drawing population conclusions.

permutation returns a shuffled copy or a permutation of positions, while shuffle modifies an array in place:

shuffled_copy = rng.permutation(order_ids)

working = order_ids.copy()
rng.shuffle(working)

Prefer a copy when later code still needs original order. In-place mutation is safe only when ownership is clear.

Reproducible does not mean representative

Suppose a population contains eight North-zone records and four South-zone records. A reproducible simple random sample may contain no South records. The seed lets another analyst replay that outcome; it does not repair coverage bias or guarantee balance.

When important groups are small, stratified sampling can sample within each group:

zones = np.array(["North"] * 8 + ["South"] * 4)
north = np.flatnonzero(zones == "North")
south = np.flatnonzero(zones == "South")

selected = np.concatenate([
    rng.choice(north, size=3, replace=False),
    rng.choice(south, size=3, replace=False),
])

This design gives equal sample counts, not population-proportional counts. An unweighted mean of the six sampled records would overrepresent South relative to the original population. Sampling design and analysis must be planned together.

For an auditable random workflow, record:

• NumPy version and generator family when exact replay matters.

• Root seed or saved generator state.

• Population definition, exclusions, and ordering.

• Sample size, replacement policy, probabilities, and stratification.

• The code version and operation sequence that consumed random state.

Do not use ordinary pseudorandom generators for passwords, tokens, or cryptographic secrets. Those tasks require security-focused facilities such as Python's secrets module.

The final section of this chapter looks beneath array syntax. Understanding which operations share storage and how many full-size buffers a pipeline creates will help you write code that is both correct and efficient.