Uniform Random Permutations and the Fisher-Yates Shuffle

A shuffle must produce each of the n! orderings with probability exactly 1/n!. This is a much stronger requirement than “looks random”, and most short implementations of a shuffle fail it. The failures are not subtle in the sense of being hard to detect: the incorrect algorithm below visits only a fraction of the permutations, and the fraction shrinks as n grows. The correct algorithm is barely longer, so the value in understanding the distinction is mostly in knowing which mutations of the correct loop are safe.

Two loops that look equivalent

The correct loop walks a position from the end down to the second element and swaps the element at that position with an element drawn uniformly from positions 0..i inclusive. After the loop, positions 0..n-1 hold a uniformly random permutation.

import random

def shuffle(a):
    for i in range(len(a) - 1, 0, -1):
        j = random.randrange(i + 1)   # 0 <= j <= i
        a[i], a[j] = a[j], a[i]

The common incorrect variant draws j from the entire array, 0..n-1, at every step.

def shuffle_biased(a):
    n = len(a)
    for i in range(n):
        j = random.randrange(n)       # 0 <= j <= n-1
        a[i], a[j] = a[j], a[i]

Both perform n - 1 or n swaps and both leave every element somewhere in the array. The first produces each permutation with probability 1/n!; the second produces each with probability n^{-n} times an integer count, and for n = 3 that count is 4 for some permutations and 5 for others, out of 3^3 = 27 equally likely swap sequences. Three of the six permutations get probability 4/27 and three get 5/27. For n = 4 the counts range from 8 to 15 out of 4^4 = 256 sequences. The bias does not vanish as n grows; the ratio between the most and least likely permutation grows.

The reason is a counting argument. The correct loop makes exactly one choice per position, with i + 1 options at step i, so the number of possible choice sequences is n × (n-1) × ... × 2 = n!, and the mapping from choice sequence to permutation is a bijection. The incorrect loop makes n choices with n options each, giving n^n sequences, and n^n is not a multiple of n! for any n ≥ 3. Any algorithm whose randomness has a sample space that is not a multiple of n! cannot assign equal probability to all permutations, no matter how the choices are mapped.

Why the shrinking bound is the whole trick

It helps to view the loop as building the permutation from the back. At step i, positions i+1 .. n-1 are already fixed. The element placed at position i is chosen uniformly from the i + 1 elements not yet fixed, which is exactly the conditional distribution required for a uniform permutation: given any assignment of the last n-1-i positions, every remaining element is equally likely to occupy position i.

The invariant is that after processing position i, the suffix is a uniformly random ordered selection of n - i elements from the input, and the prefix 0..i-1 holds the rest in arbitrary order. The step that establishes the invariant for i-1 chooses one of the i prefix elements uniformly and moves it to position i-1; since the prefix is an arbitrary arrangement of the remaining elements, each is equally likely.

This is also why the ascending form of the correct loop, drawing j from i..n-1 and swapping with position i, is equally correct. The direction is a convention. What matters is that the draw at step i ranges over exactly the positions that are still undecided, and that the element moved into the decided position is drawn from that set.

The inside-out variant for streams

When the input arrives one element at a time and its length is unknown, the array cannot be indexed by position. The inside-out algorithm builds the permutation as it reads, keeping the partial permutation in an array that grows. For each new element at index i, it draws j uniformly from 0..i; if j < i it moves the element currently at j to position i and stores the new element at j, otherwise it stores the new element at i.

def shuffle_stream(iterable):
    out = []
    for i, x in enumerate(iterable):
        j = random.randrange(i + 1)
        if j < i:
            out.append(out[j])
            out[j] = x
        else:
            out.append(x)
    return out

After i + 1 elements have been read, out holds a uniformly random permutation of them, and the argument is the same suffix invariant run in reverse: the new element is inserted at a uniformly chosen position among the i + 1 available, and the element displaced to the end takes the position the new element vacated. The output array is exactly the input length, so the memory cost is one element per input element and no second pass is needed.

The variant matters for reservoir sampling, where the same displacement trick keeps a uniform sample of fixed size from a stream of unknown length. It also matters for correctness under concurrent iteration: the stream is consumed once, in order, and no element is revisited.

The permutation as a mapping and its cycle structure

A shuffle can be represented as a function from positions to positions. Composing that function with itself repeatedly partitions the n positions into cycles, and the distribution of cycle lengths is a property of the uniform distribution rather than of the algorithm. The expected number of cycles in a uniform permutation of n elements is the n-th harmonic number, approximately ln n + 0.577, and it is obtained by writing the count as a sum of indicators and applying linearity, the same decomposition used throughout Randomized Recurrences and the Indicator Variable Method. The expected length of the cycle containing a given element is (n + 1)/2, and the probability that a uniform permutation has no fixed point is about 1/e.

These numbers are useful as a test. In a uniform permutation a given position holds its original element with probability exactly 1/n, and the expected number of fixed points over the whole permutation is exactly 1 for every n. The probability that at least one fixed point exists tends to 1 - 1/e, about 0.632, but the expectation stays at 1 because the count is usually 0 or 1 and rarely larger. An algorithm with a systematic tendency to leave elements near their original position shows an expected count well above 1, and a few thousand trials are enough to see the difference.

Cycle structure also explains why some shuffle applications care about more than uniformity. A permutation used as a key schedule or as a substitution table should not have short cycles, because a short cycle limits the number of distinct states reachable. A uniform permutation has short cycles with the probabilities the cycle-length distribution predicts, so uniformity does not by itself rule out a fixed point or a 2-cycle. If those are unacceptable, the permutation must be sampled from a conditioned distribution, which requires either rejection or a different construction.

Cost and the draws it consumes

The loop performs n - 1 swaps and consumes n - 1 bounded random draws, each with a different bound. The draw at step i needs a uniform integer in [0, i], so the word width required shrinks as the loop proceeds; a 32-bit draw suffices until i exceeds 2^32. The total number of random bits consumed is sum of ceil(log2(i+1)) for i from 1 to n-1, which is approximately n log2 n - 1.44n.

That figure is the reason a shuffle is not free at large n: for a million elements it consumes on the order of twenty million bits, and each bounded draw that uses rejection may consume more than one word. Where the generator is a system call or a hardware device with a per-call cost, the number of calls dominates. Batching draws into a buffer, or using a generator that fills a block at once, changes the constant substantially without changing the algorithm.

The swaps themselves are n - 1 exchanges, each touching two positions, for 2(n-1) memory accesses. The index j is drawn from a range that shrinks but is not confined to a locality window, so each access is likely to land on a cache line that is not resident. With 64-byte lines and 8-byte elements, one pass over the array is n/8 lines; a shuffle makes 2(n-1) accesses over those lines, so the miss count is close to the access count rather than to the line count. Blocked variants that restrict swaps to a window improve locality but are no longer exactly uniform, because positions outside the window can never be reached from inside it.

The draws themselves come from the primitive in Unbiased Random Integers from a Finite Random Source, and the off-by-one in randrange(i + 1) versus randrange(i) is the most common defect in hand-written shuffles. A shuffle that draws from 0..i-1 at step i never leaves position i holding the element that started there, so the output is always a derangement: no element stays in its original position. The sample space has size (n-1)!, while the number of derangements of n elements is the nearest integer to n!/e. For n = 4 those are 6 and 9, so the procedure cannot reach every derangement at all, let alone with equal probability. The output looks shuffled, which is why the defect survives casual inspection.