Probability Theory and Distributions

[→] The three axioms, and what each one is for

1. Non-negativity: P(A) ≥ 0 for every event A A probability is never negative. This is a definition, not something derived — it rules out nonsense results, and any calculation that produces a negative "probability" contains an error somewhere upstream.
2. Normalisation: P(Ω) = 1 Something in the sample space happens with certainty. This is what turns probability into a number between 0 and 1 rather than an arbitrary positive quantity — it fixes the scale.
3. Countable additivity: for disjoint events A₁, A₂, … — P(A₁ ∪ A₂ ∪ …) = P(A₁) + P(A₂) + … If two events cannot both happen, the probability that either happens is the sum. This is the axiom that makes a PMF sum to 1 (20.4) and is the engine behind almost every "add up the cases" calculation in this chapter.

[→] Independence, and why it's a qualification rather than a definition of "unrelated"

Events A and B are independent if knowing one occurred tells you nothing about the other: P(A | B) = P(A). Substitute this into the multiplication rule above:

P(A ∩ B) = P(A | B) · P(B)      <- multiplication rule, always true
            = P(A) · P(B)          <- substitute P(A|B) = P(A), the definition of independence

P(A ∩ B) = P(A)·P(B) is therefore not a separate fact to memorise — it is what independence means, run through the definition of conditional probability. This is the single most-used step in this chapter: every Binomial, Geometric and Negative Binomial derivation multiplies probabilities across trials, and doing so is only legal because the trials are assumed independent. Change that assumption and every formula in Part B stops applying.

[→] Linearity of expectation, and why it needs no independence at all

E[X + Y] = E[X] + E[Y], for any X and Y — correlated, dependent, anything. This single fact is the most-used tool in this entire chapter, because it is what lets you compute E[X] for a complicated X by writing it as a sum of simple pieces.

E[X+Y] = ∑∑ (x+y)·P(X=x, Y=y)               <- definition of expectation, applied to X+Y
       = ∑∑ x·P(X=x,Y=y) + ∑∑ y·P(X=x,Y=y)   <- distribute the sum; legal because
                                                    sums are finite/absolutely convergent
       = ∑x · [∑y P(X=x,Y=y)] + ∑y · [∑x P(X=x,Y=y)]  <- factor out x and y, regroup
       = ∑x · P(X=x) + ∑y · P(Y=y)      <- the bracketed sum is exactly the
                                                    marginal PMF of X (or Y) -- summing
                                                    out the other variable
       = E[X] + E[Y]                              <- definition of expectation, again

The step worth noticing is the third line: summing P(X=x, Y=y) over every value of y recovers P(X=x) by definition of a marginal distribution — it doesn't matter whether X and Y are independent, because you're just regrouping a sum, not multiplying anything. That's why linearity holds unconditionally, and it's why the Binomial's mean (20.7) is np in three lines despite the trials being correlated through the fixed total count.

[→] Deriving Var(X) = E[X²] − E[X]²

Var(X) = E[(X - E[X])²]                        <- definition
       = E[X² - 2X·E[X] + E[X]²]              <- expand the square; E[X] is a
                                                  constant number, not a random variable
       = E[X²] - 2E[X]·E[X] + E[X]²         <- linearity: E[aX] = a·E[X] for a
                                                  constant a, and E[constant] = that constant
       = E[X²] - 2E[X]² + E[X]²               <- combine 2E[X]·E[X] = 2E[X]²
       = E[X²] - E[X]²                        <- combine the last two terms

Two steps are doing real work here. Treating E[X] as a plain number (not a random quantity) is what makes "E[aX] = a·E[X]" apply to it in step three — that scalar rule is itself a one-line consequence of linearity, since a constant is just a degenerate random variable. This is why nearly every derivation below computes E[X²] as its hard step and then simply subtracts E[X]²: the identity converts "find the spread" into "find one more moment," which is usually an easier sum or integral than working with the definition directly.

[→] Deriving E[X] and Var(X)

E[X]  = 0·(1-p) + 1·p = p                    <- definition of expectation: sum over
                                                the two possible values of X

E[X²] = 0²·(1-p) + 1²·p = p                  <- same sum, but X² takes the same
                                                two values (0 and 1) as X does

Var(X) = E[X²] - E[X]² = p - p² = p(1-p)      <- the computational identity from 19.5

The reason E[X²] = E[X] here specifically is that 0²=0 and 1²=1 — X is one of the very few random variables where squaring doesn't change the value, so this shortcut is a property of the {0,1} support, not a general trick. E[X] = p, Var(X) = p(1−p). Variance is maximised at p = 0.5 (a fair coin is the most unpredictable single trial) and is exactly 0 at p = 0 or p = 1, where the outcome is certain and there's nothing left to vary.

[→] The mean, the fast way: linearity beats the PMF

You could compute E[X] by summing k·P(X=k) over the whole PMF — it works, but it's a genuinely awkward sum. Linearity of expectation (20.5) makes it three lines instead.

Let Xi = 1 if trial i succeeds, 0 otherwise    <- Xi is Bernoulli(p), one per trial
Then X = X1 + Xi + ... + Xn             <- total successes is the sum of indicators

E[X] = E[X1] + E[Xi] + ... + E[Xn]       <- linearity (19.5) -- holds even though the
                                                Xi are dependent through the fixed n
     = p + p + ... + p   (n terms)          <- each Xi is Bernoulli(p), so E[Xi] = p (19.6)
     = np

The subtlety worth sitting with: the X_i are not independent in the sense that knowing nine of the twenty trial outcomes changes nothing about your belief in the tenth (they genuinely are independent, since each trial is separate) — but even if they weren't, linearity wouldn't care, because it never used independence in its proof (20.5). This is exactly why the trick generalises: any count of successes across n trials has mean n·p(per trial success chance), whether or not the trials are independent, as long as each trial has the same success probability.

[→] The variance: this step needs independence, and here's exactly where

Variance does not add as cleanly as expectation in general. For two random variables, Var(X+Y) = Var(X) + Var(Y) + 2Cov(X,Y), and the covariance term vanishes only when X and Y are independent. This is the first derivation in the chapter where independence is doing indispensable work rather than just being assumed for convenience.

Var(X) = Var(X1 + Xi + ... + Xn)
       = Var(X1) + Var(Xi) + ... + Var(Xn)    <- Cov(Xi,Xj) = 0 for i≠j because
                                                the trials are INDEPENDENT -- this is
                                                the step that fails for correlated trials
       = p(1-p) + p(1-p) + ... + p(1-p)  (n terms)   <- each Xi is Bernoulli(p), Var = p(1-p) (19.6)
       = np(1-p)

E[X] = np, Var(X) = np(1−p). The contrast with the mean derivation is the entire point of showing both: expectation adds for free regardless of dependence, but variance only adds this cleanly when the trials are independent. If your twenty API calls share a resource and fail together under load, this formula for the variance is wrong even though the formula for the mean is still fine — a distinction that only becomes visible once you've seen where each proof actually uses (or doesn't use) independence.

[→] Confirming the PMF sums to 1, using the geometric series

Before trusting a mean or variance formula, it's worth checking the PMF is legitimate — axiom 2 (20.2) requires it to sum to 1. This uses the standard infinite geometric series result ∑_k=0^∞ r^k = 1/(1−r) for |r| < 1, which is where the distribution's name comes from.

∑k=1∞ (1-p)k-1·p  = p · ∑k=1∞ (1-p)k-1    <- factor out the constant p
                            = p · ∑j=0∞ (1-p)j      <- substitute j = k-1, so the sum
                                                        restarts at j=0
                            = p · [1 / (1-(1-p))]           <- geometric series formula, r = 1-p
                            = p · (1/p) = 1                 <- simplify: 1-(1-p) = p

That substitution j = k−1 is the recurring move in every derivation below: it turns a sum starting at an inconvenient index into the standard geometric series starting at 0, where the textbook formula applies directly.

[→] Deriving E[X] = 1/p, by differentiating the series

This is the classic trick for geometric-series sums involving an extra factor of k: differentiate the basic series with respect to r, because differentiation brings a factor of k down out of the exponent.

Start from:     ∑k=0∞ rk = 1/(1-r)              <- the known geometric series

Differentiate   ∑k=0∞ k·rk-1 = 1/(1-r)²      <- d/dr of both sides; term-by-term
both sides                                              differentiation is legal here because
w.r.t. r:                                               the series converges nicely for |r|<1

Now compute E[X]:
E[X] = ∑k=1∞ k·(1-p)k-1·p           <- definition of expectation
     = p · ∑k=1∞ k·(1-p)k-1          <- factor out p
     = p · [1/(1-(1-p))²]                    <- this sum is EXACTLY the differentiated
                                                series above, with r = 1-p
     = p · (1/p²) = 1/p                      <- simplify: 1-(1-p) = p

E[X] = 1/p. This should feel right before you even finish the algebra: if success happens 3% of the time, you expect to wait 1/0.03 ≈ 33 trials for one — rarer successes mean longer waits, and the relationship is a clean reciprocal.

[→] Deriving Var(X) = (1−p)/p²

Differentiating once gave ∑k·r^k−1. Differentiating again brings down a factor of k(k−1), which is exactly the factorial-moment shape used throughout this chapter for handling E[X²] cleanly.

Differentiate    ∑k=0∞ k(k-1)·rk-2 = 2/(1-r)³    <- differentiate the first derivative's
again:                                                 result again w.r.t. r

E[X(X-1)] = ∑k=1∞ k(k-1)·(1-p)k-1·p
          = p·(1-p) · ∑k=1∞ k(k-1)·(1-p)k-2   <- pull out one factor of (1-p) so the
                                                        remaining sum matches the 2nd-derivative form
          = p·(1-p) · [2/(1-(1-p))³]         <- substitute r = 1-p into the 2nd-derivative series
          = p·(1-p) · (2/p³) = 2(1-p)/p²

Now assemble the variance:
E[X²] = E[X(X-1)] + E[X]              <- algebra: X² = X(X-1) + X, so this is just
                                        linearity applied to that identity
         = 2(1-p)/p² + 1/p
         = [2(1-p) + p] / p²           <- common denominator
         = (2-p) / p²                  <- simplify 2(1-p)+p = 2-2p+p = 2-p

Var(X) = E[X²] - E[X]²                <- the identity from 19.5
       = (2-p)/p² - 1/p²
       = (1-p)/p²

The identity X² = X(X−1) + X is worth remembering on its own: it is what converts "the hard moment" E[X²] into "the factorial moment" E[X(X−1)], which is almost always easier to compute for a distribution defined by products of probabilities, because the (X−1) cancels an awkward term in the sum. The same identity reappears in the Poisson derivation (20.10).

[→] The mean and variance, without touching the PMF at all

Waiting for r successes is the same as waiting for the 1st success, then — once that resets — waiting for the next one, r times over. Each individual wait is Geometric(p), and the total trial count is their sum.

Let W1 = trials for the 1st success
    W2 = ADDITIONAL trials for the 2nd success, after the 1st
    ...
    Wr = additional trials for the r-th success, after the (r-1)-th

Each Wi is Geometric(p), and by the memorylessness of Geometric (19.8), each
new waiting period restarts fresh regardless of how the previous one ended --
so the Wi are INDEPENDENT and IDENTICALLY DISTRIBUTED.

X = W1 + W2 + ... + Wr

E[X] = E[W1] + E[W2] + ... + E[Wr]     <- linearity (19.5); doesn't even need the Wi
                                        to be independent, just identically distributed
     = (1/p) + (1/p) + ... + (1/p)  (r terms)   <- each Wi is Geometric(p), E[Wi]=1/p (19.8)
     = r/p

Var(X) = Var(W1) + Var(W2) + ... + Var(Wr)   <- variance of a sum of INDEPENDENT terms
                                            adds (19.7) -- this step needs the
                                            independence that memorylessness gave us
       = r · (1-p)/p²                          <- each Wi has Var = (1-p)/p² (19.8)

E[X] = r/p, Var(X) = r(1−p)/p². Set r=1 and both formulas collapse exactly to the Geometric's 1/p and (1−p)/p² — which is a useful sanity check on any formula you derive this way: it should reduce to the simpler distribution it generalises when you plug in the boundary case.

[→] Deriving the PMF as a limit of Binomial(n, λ/n)

Imagine chopping one second into n tiny slices, each slice an independent Bernoulli trial with a tiny success chance p = λ/n — small enough that at most one event can plausibly land in any one slice. As n → ∞ the slices become infinitesimally fine and this becomes an exact model of events landing at random instants. The number of "successes" is Binomial(n, λ/n), and its mean stays fixed at n·(λ/n) = λ no matter how fine the slicing gets — which is exactly the property you want for "a rate of λ per second."

P(X=k) = lim  C(n,k) · (λ/n)k · (1-λ/n)n-k
       n→∞

Split C(n,k) = n!/(k!(n-k)!) and regroup into three separate limits:

  n! / [(n-k)! · nk]  →  1              <- this ratio is n(n-1)...(n-k+1)/n^k,
                                             a product of k terms each → 1 as n→∞
                                             (e.g. (n-1)/n → 1, (n-2)/n → 1, ...)

  (1 - λ/n)n          →  e-λ       <- the standard calculus limit that DEFINES e:
                                             lim (1+x/n)^n = e^x, with x = -λ

  (1 - λ/n)-k         →  1              <- the exponent -k is FIXED while n→∞,
                                             so this factor just goes to 1^(-k) = 1

Multiply the surviving pieces together with the untouched λk/k!:

P(X=k) = 1 · λk/k! · e-λ · 1 = e-λλk/k!

The e^−λ term is not a special constant chosen to make the formula work — it's the ordinary calculus limit lim(1+x/n)ⁿ = e^x showing up because "many independent trials, each individually unlikely" is precisely the setup that definition describes. Anywhere you see e^−λ in a probability formula, a hidden "n → ∞, p → 0" limit is almost always the reason.

[→] Mean and variance: take the limit of np and np(1−p) directly

Rather than re-deriving from the Poisson PMF, take the Binomial's already-proven mean and variance (20.7) and push p = λ/n through the same limit.

E[X] = lim n·p = lim n·(λ/n) = λ                    <- the n's cancel exactly, for every n

Var(X) = lim n·p·(1-p) = lim n·(λ/n)·(1-λ/n)
       = lim λ·(1-λ/n) = λ·(1-0) = λ      <- as n→∞, λ/n → 0, so the
                                                second factor just goes to 1

E[X] = λ, Var(X) = λ. Mean equals variance — a fact unique enough among named distributions that it functions as a diagnostic in its own right (20.17): if a real dataset's sample variance is far from its sample mean, that alone is evidence the data isn't really Poisson, independent of any other check.

[→] The mean, by the same indicator trick as Binomial

Write X as a sum of n indicators exactly as in 20.7 — X_i = 1 if the i-th draw is a success. Linearity doesn't care that the draws are dependent (20.5), so the mean comes out just as easily as it did for Binomial.

By symmetry, EVERY individual draw -- 1st, 2nd, ..., n-th -- has the same marginal
probability of being a success: K/N. (This is a genuinely non-obvious fact: the 5th
card drawn from a shuffled deck is exactly as likely to be an ace as the 1st card is,
even though which cards remain by then depends on everything drawn before it.)

E[X] = E[X1] + E[X2] + ... + E[Xn]    <- linearity (19.5), regardless of dependence
     = (K/N) + (K/N) + ... + (K/N)  (n terms)
     = n·K/N

E[X] = nK/N — identical in shape to the Binomial mean with p = K/N, which makes sense: on average, sampling without replacement still draws successes in proportion to how common they are in the population.

[→] E[X] and Var(X), by direct integration

The continuous analogue of ∑k·p(k) is ∫x·f(x)dx — same idea, sum replaced by integral because x now varies continuously rather than jumping between discrete values.

E[X] = ∫ab x · 1/(b-a) dx
     = 1/(b-a) · [x²/2]ab            <- antiderivative of x is x²/2
     = 1/(b-a) · (b²-a²)/2
     = (b-a)(b+a) / [2(b-a)]              <- factor the difference of squares b²-a²
     = (a+b)/2                            <- the (b-a) terms cancel

E[X²] = ∫ab x² · 1/(b-a) dx
       = 1/(b-a) · [x³/3]ab
       = (b³-a³) / [3(b-a)]
       = (a²+ab+b²)/3            <- factor b³-a³ = (b-a)(b²+ab+a²), cancel (b-a)

Var(X) = E[X²] - E[X]²                <- the identity from 19.5, unchanged for continuous X
       = (a²+ab+b²)/3 - (a+b)²/4
       = (b-a)²/12                    <- combine over a common denominator of 12; the
                                       algebra is mechanical but the cancellation
                                       is real -- try it with a=0,b=1 to check: 1/12

E[X] = (a+b)/2 is just the interval's midpoint, exactly as intuition demands for a flat density. Var(X) = (b−a)²/12 is worth keeping as a reference point: it's the variance formula every other continuous distribution in this chapter gets compared against, since Uniform is "spread with no preference at all" — the baseline case.

[→] Deriving the PDF from "zero events so far," the Poisson connection

The cleanest route to the Exponential is through its CDF, using the Poisson distribution already derived in 20.10 rather than starting from the PDF as given.

Let X = waiting time for the first event of a rate-λ Poisson process.
"X > t" means EXACTLY: zero events occurred in the interval [0, t].

The number of events in an interval of length t is Poisson(λt) -- the rate
λ scaled up by however long you wait (this scaling is a property of a
Poisson process, taken as given here).

P(X > t) = P(zero Poisson(λt) events)
         = e-λt(λt)0/0!        <- Poisson PMF at k=0, from 19.10
         = e-λt

So the CDF is F(t) = P(X ≤ t) = 1 - e-λt

f(t) = F'(t) = λe-λt          <- differentiate the CDF to recover the PDF;
                                       this is the defining relationship from 19.4

This route matters more than it looks: it shows Exponential isn't an independent assumption bolted onto Poisson — it's a restatement of the same process, just asking "how long until the first event" instead of "how many events in a fixed window." Whatever assumptions justified Poisson (20.10's independence and no-clustering conditions) automatically justify Exponential too, and whatever breaks one breaks the other.

[→] E[X] = 1/λ, by integration by parts

E[X] = ∫0∞ x · λe-λx dx

Integrate by parts: ∫u dv = uv - ∫v du
  Let u = x            dv = λe-λx dx
      du = dx           v = -e-λx          <- antiderivative of λe^(-λx) is -e^(-λx)

E[X] = [-x·e-λx]0∞  +  ∫0∞ e-λx dx
     = (0 - 0)                          + [-1/λ · e-λx]0∞
     = 0 + (0 - (-1/λ))
     = 1/λ

The boundary term [-x·e^−λx] at ∞ needs a moment's justification: exponential decay beats polynomial growth, so x·e^−λx → 0 as x → ∞ even though x itself is growing without bound. This same fact makes every moment of the Exponential finite, unlike some heavier-tailed continuous distributions where it doesn't hold.

[→] The Central Limit Theorem, stated precisely

Take any distribution with finite mean μ and finite variance σ² — it can be wildly skewed, discrete, bounded, anything. Draw n independent samples from it and average them. The CLT says:

As n → ∞,   (X̄ - μ) / (σ/√n)   converges in distribution to Normal(0, 1)

Equivalently:  X̄  is approximately Normal(μ, σ²/n)  for large n

The remarkable part is what's absent from the statement: no assumption about the shape of the original distribution. A sum or average of many small independent contributions looks approximately Normal almost regardless of where those contributions came from — which is why measurement errors, aggregated user behaviour, and sums of many small effects all tend toward a bell curve even when no individual effect is remotely bell-shaped. This is also precisely the fact chapter 17 relied on for confidence intervals on a sample mean (17.7), used there without re-deriving it.

[→] The continuity correction, and why skipping it is a real error

Binomial is discrete; Normal is continuous. Approximating P(X=k) for the discrete variable by plugging k directly into the continuous Normal CDF quietly asks the wrong question — a single point has probability exactly 0 under any continuous distribution (20.4). The fix is to treat the discrete value k as if it occupies the continuous interval [k−0.5, k+0.5].

Binomial(20, 0.5): exact P(X ≤ 10) = 0.5881

Without correction:  P(Z < (10-10)/√5)         = P(Z<0)     = 0.5000   ← visibly wrong
With correction:      P(Z < (10.5-10)/√5)       = P(Z<0.224)  = 0.5885   ← matches closely

The uncorrected version is off by nearly 9 percentage points here, entirely from forgetting that "at most 10" for a discrete variable means "up to but not including 11.5" on the continuous number line the Normal lives on. Always add or subtract 0.5 at the boundary you're evaluating, in the direction that includes the boundary value itself.

[→] Ask these in order, on every problem

1. 1. Discrete count, or continuous measurement? Is the answer an integer (number of failures, number of clicks) or a real number on a continuum (time, distance, amount)? This alone splits the whole chapter into Part B and Part C.
2. 2. Is there a fixed number of trials, n? "Out of 20 calls" is a fixed n. "Keep going until something happens" has no fixed n — the trial count is itself the random quantity. This is the single biggest fork in the discrete side: fixed n points to Binomial or Hypergeometric; no fixed n points to Geometric, Negative Binomial, or Poisson.
3. 3. Are trials independent, with replacement (or from an effectively infinite population)? Drawing cards without putting them back, or auditing without replacement, breaks the constant-p assumption every Bernoulli-based distribution needs. This is the question that separates Hypergeometric from everything else.
4. 4. Counting successes, or waiting for one to happen? "How many defects in this batch" counts occurrences. "How long until a defect" waits for one. The same underlying process produces a counting distribution and a waiting-time distribution as two sides of the same coin — Poisson counts events in a window, Exponential times the gap until the next one (20.13).
5. 5. Waiting for exactly one occurrence, or for the r-th? "The first success" is Geometric or Exponential. "The fifth success" is Negative Binomial — the direct generalisation proved in 20.9 by summing r Geometrics.