46Randomization & the probabilistic mind

In Chapter 45 we turned strings and number theory into the machinery of crypto and quant. Every algorithm we've built this volume has made its choices deterministically: the same pivot, the same hash, the same next node, every single time. Now we do the one thing that sounds like surrender but turns out to be a superpower — we flip a coin. Here's the plan. We'll armor a sort so no attacker can find its worst case, compute a number like π by throwing darts at a wall, price a financial option by simulating a million possible futures, and draw a perfectly fair sample from a stream too big to hold. One question drives the whole chapter: how can deliberately adding randomness make an algorithm more reliable, not less? By the end you'll reach for a coin on reflex — when being predictable is the weakness, when the exact answer is out of reach, when you need fairness over data you can't even count.

01Every fixed plan has a nemesis

Let's start with a word we've leaned on all volume without saying it out loud: deterministic. A deterministic algorithm always makes the same choice on the same input. Feed it the same list twice and it does the exact same thing, step for step. Usually that's a virtue — it's what makes code testable and predictable. But watch what it costs us. A deterministic algorithm has a fixed worst case: one specific arrangement of input that drives it to its slowest behavior, every time. And if someone can read your code — it's open-source, or a public API, or a well-known library — they can construct that exact arrangement and hand it to you on purpose.

Quicksort is the classic victim. Tell it "always pick the first element as the pivot," then hand it an already-sorted list. Watch it unravel. Every partition peels off just one element, the recursion sinks n levels deep, and the tidy O(n log n) collapses into O(n²). (You stepped quicksort's partition in trace T28 — here we attack it.) That's not a thought experiment. It's a denial-of-service. Put a naive quicksort behind a web endpoint, and an attacker who sends sorted data can pin your CPU flat with a tiny request.

The escape is not a cleverer fixed rule. The adversary just studies the new rule and defeats that one too. The escape is to stop being predictable. If the algorithm flips a coin to make its choices, no single input is reliably its worst case — because the worst case now depends on coins the attacker cannot see.

So how much does a single coin flip actually buy? Let's stop asserting and start counting →

02Randomness is armor — the randomized pivot

Count the fixed-pivot disaster first. Partitioning a block of m elements costs m−1 comparisons. If every pivot is the extreme value — exactly what a sorted list guarantees — the block sizes march down n−1, n−2, …, 1, so the total is (n−1)+(n−2)+…+1 = n(n−1)/2. That's the textbook O(n²), and it isn't approximate. I ran a comparison-counting quicksort on a sorted list of n = 2000, and the fixed-pivot version reported 1,999,000 comparisons — exactly n(n−1)/2, to the comparison.

Now make the pivot random. Each pivot is equally likely to be any rank, so partitions come out balanced on average, and the expected number of comparisons works out to about 2n ln n (≈ 1.39·n·log₂n) — squarely the O(n log n) class. Same sorted input, random pivot: about 24,600 comparisons, averaged over twenty runs. That's roughly 80× fewer — on the very input that murdered the deterministic version.

The deep move: we didn't change the input, we changed who chooses. The chance that a random shuffle lands on any one specific bad arrangement is 1/n!. For just 20 items that's 1 in 2.4×10¹⁸; for a shuffled deck of 52 cards it's 1 in 8×10⁶⁷ — more arrangements than there are atoms in the galaxy. The worst case didn't vanish. It just became something you will not witness before the heat death of the universe.

def quicksort_count(a, pivot="first"):
    a = list(a); comps = 0; stack = [(0, len(a)-1)]
    while stack:
        lo, hi = stack.pop()
        if lo >= hi: continue
        if pivot == "random":                 # ← the one line that armors it
            p = random.randint(lo, hi); a[p], a[hi] = a[hi], a[p]
        pv = a[hi]; i = lo
        for j in range(lo, hi):
            comps += 1                         # count every key comparison
            if a[j] < pv:
                a[i], a[j] = a[j], a[i]; i += 1
        a[i], a[hi] = a[hi], a[i]
        stack.append((lo, i-1)); stack.append((i+1, hi))
    return comps
# sorted input, n=2000:  fixed → 1,999,000   random → ~24,600  (≈80× fewer)

Line by line: we sort a copy using an explicit stack instead of recursion, so the sorted worst case can't overflow Python's call stack. The if pivot == "random" line swaps a randomly chosen element into the pivot slot — that's the entire defense. The inner for loop is the partition, and comps += 1 tallies each key comparison, so we measure cost directly instead of trusting a formula. Everything else is ordinary Lomuto partitioning (stepped in T28).

Armor is one use of a coin. The other is stranger: using randomness to compute a number you can't otherwise reach →

03Compute by throwing darts — Monte Carlo

Now turn the idea inside out. Instead of using randomness to protect an exact algorithm, use it to estimate an answer that's hard to compute directly. The technique is called Monte Carlo, after the casino, and its motto is simple: to learn about a system, sample random outcomes from it and count.

The cleanest demonstration computes π. Draw a 1×1 square, and inside it the quarter of a circle of radius 1 centered at one corner. The square's area is 1; the quarter-circle's area is (1/4)·π·r² = π/4. Throw darts uniformly at random into the square. On average, the fraction that land inside the quarter circle must equal the ratio of the areas — π/4. So π ≈ 4 × (darts inside) / (total darts). A dart at (x, y) is inside exactly when x² + y² ≤ 1: one multiply, one compare, no trigonometry, and the symbol π appears nowhere in the code. You measure it by counting.

def estimate_pi(n):
    inside = 0
    for _ in range(n):
        x, y = random.random(), random.random()   # a dart in the unit square
        if x*x + y*y <= 1.0:                        # inside the quarter circle?
            inside += 1
    return 4 * inside / n                           # fraction inside ≈ π/4

#      n         estimate      error
#     100        3.120000     0.021593
#    1000        3.128000     0.013593
#   10000        3.126000     0.015593   ← note: WORSE than n=1000!
#  100000        3.137280     0.004313
# 1000000        3.140592     0.001001
# 10000000       3.142310     0.000717

Four lines do it. random.random() draws a coordinate in [0,1), so (x, y) is a uniform dart. The one comparison x*x + y*y <= 1.0 asks "is it inside?", inside tallies the hits, and 4 * inside / n turns the ratio into π. The table is a real run, and look closely: at n = 10,000 the estimate was actually worse than at n = 1,000. That's not a bug. A single random run can get lucky or unlucky. What improves reliably is never one run — it's the typical error, and pinning that down is the next section.

04The √N law — and pricing the future

Here's the law that governs every Monte Carlo estimate, and it's worth burning into memory: the typical error shrinks like 1/√N, where N is the number of samples. I measured it directly. Run the π estimator 200 times at N = 100, and the estimates scatter with a standard deviation of about 0.165. Quadruple to N = 400 and the scatter halves to 0.081. Quadruple again to 1600 and it's 0.039; once more to 6400 and it's 0.020. Every ×4 in samples halves the error — the fingerprint of 1/√N.

Flip that around and it stings. To win one more correct decimal digit — 10× smaller error — you need 100× more samples. Monte Carlo hands you the first digits cheaply and charges brutally for the last ones. Why 1/√N? Each sample is an independent draw, and averaging N independent draws is exactly the setup where the spread of the average narrows as 1/√N. It's the same statistics that make a poll of 1,000 people accurate to a few percent, not a few thousandths. This isn't special to π. It's the signature of averaging independent randomness.

And this is exactly why Wall Street runs on Monte Carlo. The fair price of a financial option is an expected value — the average payoff across every possible future path of the underlying stock, discounted back to today. For the simplest options there's a closed-form answer (the Black–Scholes formula). For the messy real ones — path-dependent, many assets, exotic payoffs — the integral has no formula at all. So you don't solve it; you sample it: simulate thousands of random futures, compute the payoff in each, and average. I priced a European call both ways. Black–Scholes says 10.4506. Monte Carlo, averaging simulated futures, gives 11.02 at 1,000 paths, 10.49 at 100,000 paths, and 10.4458 at a million paths — the √N law crawling toward the true price. Risk desks simulate millions of paths every night to answer "how much could we lose tomorrow?" That's expected-value thinking made computational: when you cannot reason about the average, sample it.

def mc_call(n, S0, K, r, sigma, T):
    drift = (r - 0.5*sigma**2)*T
    vol   = sigma*(T**0.5)
    total = 0.0
    for _ in range(n):
        z  = random.gauss(0, 1)               # one random shock
        ST = S0*math.exp(drift + vol*z)        # one simulated future price
        total += max(ST - K, 0.0)              # the call's payoff on that path
    return math.exp(-r*T) * total/n            # discounted average payoff
# Black–Scholes = 10.4506;  mc_call(1_000_000) → 10.4458

The loop is the whole idea. random.gauss(0, 1) draws one random "shock", S0*exp(drift + vol·z) turns it into one possible price at expiry, and max(ST − K, 0) is what the call pays on that future. Average over a million futures, discount to today, and you've evaluated an integral no one wrote down — by sampling it.

Quicksort's coins and π's darts gamble in opposite directions. Naming that difference is the whole taxonomy of randomized algorithms →

05Two ways to gamble — Las Vegas vs Monte Carlo

Look at what the coins actually put at risk. Randomized quicksort always returns a correctly sorted list — the coins only affect how long it takes. A bad streak makes it slower, never wrong. The π estimator is the mirror image. It always finishes in exactly N steps, but the answer is only approximate; a bad streak makes it wrong-ish, never slow. Those are the two families of randomized algorithm, and knowing which one you're holding tells you exactly what can go wrong.

There's a beautiful bridge between them. If you can check an answer cheaply, you can convert Monte Carlo into Las Vegas. Guess with a fast randomized method, then verify the guess; if it's wrong, guess again. The verify step turns "probably right" into "definitely right, eventually" — you've traded a small error probability for a small bit of extra runtime.

06A fair sample from an endless stream — reservoir sampling

The problem sounds impossible. You're watching a firehose — every tweet, every log line, every sensor reading — and you want to keep a uniformly random sample of, say, one item. When the stream finally ends, every item that ever flowed past should have had an equal chance of being the one you kept. The catch: you don't know how long the stream is (it could be billions of items, or never end), you get exactly one pass, and you have O(1) memory — room for the one item, nothing more.

The solution is one line of genius. Reservoir sampling: keep the first item. When the i-th item arrives, replace whatever you're holding with it with probability exactly 1/i. That's the entire algorithm — and notice that n, the stream's length, never appears.

def reservoir_one(stream):
    kept = None
    for i, item in enumerate(stream, start=1):   # i = 1, 2, 3, ...
        if random.random() < 1.0 / i:            # newest barges in w.p. 1/i
            kept = item
    return kept
# length-10 stream, 1,000,000 trials — how often each item was the survivor:
#   item 0: 10.00%   item 1: 10.04%   item 2:  9.98%   item 3:  9.99%
#   item 4:  9.97%   item 5: 10.00%   item 6: 10.06%   item 7:  9.98%  ...
# every item ≈ 1/10, dead fair — with room for exactly one item in memory

The loop reads one item at a time. At position i, random.random() < 1.0/i is true with probability 1/i. So item 1 is kept for sure, item 2 has a 1-in-2 shot at bumping it, item 3 a 1-in-3 shot, and so on. The comment is a real million-trial run: every one of the ten items ended up the survivor almost exactly a tenth of the time, deviating by at most 0.06 percentage points. Perfectly fair, in one pass, in constant memory.

Here's the magic, worked out. Consider item i. It gets kept at step i with probability 1/i. To still be the survivor at the end, it must then not be replaced by items i+1, i+2, …, n. Item j replaces with probability 1/j, so it fails to replace with probability (j−1)/j. Multiply the whole chain of survival:

(1/i) · (i/(i+1)) · ((i+1)/(i+2)) · … · ((n−1)/n) = 1/n

and it telescopes: every numerator cancels the denominator before it, leaving exactly 1/n. Every item, whether it arrived first or billionth, ends with probability 1/n. The math is fair by construction, and n is never needed.

Step back and see the shape of this whole chapter. Randomness bought three different things: armor (no adversary can aim at a random pivot), computation (sample the darts, the futures, the light rays, and average your way to an answer no formula gives), and fairness under ignorance (a uniform sample without ever knowing how much there is). The transferable move — the one the 1% run on reflex — is this: when being predictable is your weakness, add a coin; when the exact answer is out of reach, estimate it by trying many random cases and averaging; and when you want a property like fairness, work backwards from the property to the rule that produces it. Stanislaw Ulam invented Monte Carlo not by out-computing the neutron equations but by asking a humbler question. Recovering from illness and playing solitaire, he wondered about his odds of winning, and realized that instead of grinding through the combinatorics he could just deal many random hands and count. That question — "what if I just try it thousands of times and count?" — is available to you on every hard problem you'll ever meet.

Six chapters of patterns; one volume of thinking. Next we assemble them into a single operating system for the algorithmic mind — how to look at any new problem and know, in seconds, which technique it wants. Chapter 47 is the difference between the 99% and the 1%, and it closes the volume →