Unit 17: Chance

A probability is a count divided by a count

Two dice have 36 outcomes, all equally likely. Six of them add up to 7.

So the chance of 7 is 6 out of 36, which is one in six. When the outcomes are equally likely, probability is nothing but counting — which is the whole of unit 15 again.

Python

from fractions import Fraction

outcomes = [(a, b) for a in range(1, 7) for b in range(1, 7)]
sevens = [p for p in outcomes if p[0] + p[1] == 7]

print(len(outcomes), len(sevens))
print(Fraction(len(sevens), len(outcomes)))

It prints

36 6
1/6

The whole shape of two dice

Seven is the commonest because it has the most ways to happen. Two and twelve have one each.

Adding a third die makes the middle even more crowded, which is the reason a total of 10 or 11 is a good bet in a lot of board games.

Python

counts = [0] * 13
for a in range(1, 7):
    for b in range(1, 7):
        counts[a + b] += 1

for total in range(2, 13):
    print(total, counts[total])

It prints

2 1
3 2
4 3
5 4
6 5
7 6
8 5
9 4
10 3
11 2
12 1

Keep it as a fraction while you can

Contest answers are usually asked for as a fraction, or as a fraction under a modulus using unit 14's inverse.

Decimals lose exactness and a judge comparing to six places will notice. Fraction is exact and slow; use it to check, and integers to submit.

Every outcome, weighted by its chance

The expected value is what you would average over a very long run. Add up each outcome times its probability.

One die averages 3.5 — a number it can never actually show, which is the first surprise about expectations.

Python

expected = 0
for face in range(1, 7):
    expected += face * (1 / 6)

print(round(expected, 4))

It prints

3.5

Expectations simply add up

The expected total of two dice is 7 — and you can get that by adding 3.5 and 3.5 rather than by working through all 36 outcomes.

That is linearity of expectation, and it holds even when the things are not independent. It is the most useful single fact in contest probability.

Python

slow = 0
for a in range(1, 7):
    for b in range(1, 7):
        slow += (a + b) / 36

print(round(slow, 4))
print(round(3.5 + 3.5, 4))

It prints

7.0
7.0

The trick it makes possible

"How many of the 36 outcomes have a double somewhere" is a counting job. "What is the expected number of dice showing a six" is not: each die contributes 1/6, so two dice give 1/3. No outcomes are listed at all.

Break the thing you are counting into a sum of tiny yes/no questions, and add up their chances.

Python

print(round(2 * (1 / 6), 4))

slow = 0
for a in range(1, 7):
    for b in range(1, 7):
        sixes = (1 if a == 6 else 0) + (1 if b == 6 else 0)
        slow += sixes / 36

print(round(slow, 4))

It prints

0.3333
0.3333

A row of five squares, and a coin

You start on square 2. Each turn you step left or right, half and half. Squares 0 and 4 are traps: land there and you stop.

Instead of following one journey, follow the probability: keep a number for each square saying how likely you are to be standing on it.

Each step spreads it out

Every square hands half its probability to each neighbour — except the traps, which keep theirs.

After one step you are on 1 or 3. After two, back on 2 half the time and caught a quarter of the time each side. This is a Markov chain.

Python

probs = [0.0] * 5
probs[2] = 1.0

for step in range(3):
    nxt = [0.0] * 5
    nxt[0] = probs[0]
    nxt[4] = probs[4]
    for i in range(1, 4):
        nxt[i - 1] += probs[i] / 2
        nxt[i + 1] += probs[i] / 2
    probs = nxt
    print([round(x, 4) for x in probs])

It prints

[0.0, 0.5, 0.0, 0.5, 0.0]
[0.25, 0.0, 0.5, 0.0, 0.25]
[0.25, 0.25, 0.0, 0.25, 0.25]

Run it long enough and it settles

After a hundred steps almost everything has fallen into a trap: half at each end, because you started exactly in the middle.

The probabilities always add to 1, which is the check to make after every step. If they drift, a rule is handing out chance it does not have.

Python

probs = [0.0] * 5
probs[2] = 1.0

for step in range(100):
    nxt = [0.0] * 5
    nxt[0] = probs[0]
    nxt[4] = probs[4]
    for i in range(1, 4):
        nxt[i - 1] += probs[i] / 2
        nxt[i + 1] += probs[i] / 2
    probs = nxt

print([round(x, 4) for x in probs])
print(round(sum(probs), 4))

It prints

[0.5, 0.0, 0.0, 0.0, 0.5]
1.0

And it is a matrix, if you want it to be

One step is a fixed rule applied to a list of numbers — which is exactly unit 16's matrix multiplication.

So "where am I after a thousand million steps" is a matrix power, in about thirty multiplications. That is the join between the last unit and this one.

Throwing until a six

How many throws, on average? The answer is a sum with no end: one throw times the chance it works, two throws times the chance the first missed and the second worked, and so on for ever.

Add two hundred terms and it has clearly stopped moving: 6.

Python

p = 1 / 6
total = 0.0
for k in range(1, 200):
    total += k * p * (1 - p) ** (k - 1)

print(round(total, 6))

It prints

6.0

The one-line way to get the same 6

Throw once. Either it is a six — done — or it is not, and you are back exactly where you started with the same expected wait ahead of you.

So E = 1 + (5/6) E, which rearranges to E = 6. A self-referential equation, and no series at all. Expected waits are almost always solved this way.

All six faces, not just one

Waiting for every face is the same idea six times. When you have seen k faces, the chance a throw shows a new one is (6 - k) / 6, so the wait for it is 6 / (6 - k).

Add those six waits: 1, then 1.2, then 1.5, 2, 3, 6. It comes to 14.7 — twice as long as most people guess.

Python

total = 0.0
for k in range(6):
    wait = 6 / (6 - k)
    total += wait
    print(k, round(wait, 4), round(total, 4))

It prints

0 1.0 1.0
1 1.2 2.2
2 1.5 3.7
3 2.0 5.7
4 3.0 8.7
5 6.0 14.7

And a check by playing it out

Ten thousand pretend throws of two dice, counting the sevens — using the seeded generator from the AI track rather than random, so the number below is the same every time anyone runs it.

0.1614 against a true 0.1667. Simulation is a way to check an answer, never to produce one: it is out by 3% here, and a judge wants six decimal places.

Python

seed = 12345
hits = 0
trials = 10000

for i in range(trials):
    seed = (seed * 75 + 74) % 65537
    a = seed % 6 + 1
    seed = (seed * 75 + 74) % 65537
    b = seed % 6 + 1
    if a + b == 7:
        hits += 1

print(hits, round(hits / trials, 4), round(1 / 6, 4))

It prints

1614 0.1614 0.1667