Worked solutions — Probability & Discrete Distributions

Exercise 1 — Boy-or-girl paradox

Sample space: {BB, BG, GB, GG}, each with probability 1/4.

(a) Condition on "at least one girl": 3 of 4 outcomes survive. Only GG = both girls. So $P = 1/3$.

(b) Condition on "older is a girl": 2 of 4 outcomes survive ({GB, GG}). Only GG. So $P = 1/2$.

(c) The two conditioning events are different. The second is strictly stronger (more specific), shrinking the sample space more.

set.seed(1)
fams <- matrix(sample(c("B","G"), 2*1e6, replace=TRUE), ncol=2)
mean(rowSums(fams=="G")==2) / mean(rowSums(fams=="G")>=1) # ~ 0.333
mean(fams[,1]=="G" & fams[,2]=="G") / mean(fams[,1]=="G") # ~ 0.500

Exercise 2 — Taxi / Bayes

Let $B$ = blue taxi (prior 0.01), $W_B$ = witness says blue.

(0.99*0.01) / (0.99*0.01 + 0.05*0.99) # 0.1667

Interpretation: only ~17% chance the taxi really was blue — base-rate dominance.

Exercise 3 — Tallest among neighbors

Define $X_i = 1$ if person $i$ is shorter than both neighbors. The indicator approach:

For any three distinct people, the smallest is in the "middle" with probability 1/3 — by symmetry of orderings, the middle is equally likely to be smallest, middle, or largest among the three. So $E[X] = n/3$.

n <- 50; reps <- 1e5
samples <- replicate(reps, {
 h <- sample(n)
 sum(h < c(h[n], h[-n]) & h < c(h[-1], h[1]))
})
mean(samples) # ~ 16.67

Exercise 4 — Cystic fibrosis (basic)

(a) $P(\text{unaffected child}) = 1 - 1/4 = 3/4$.

(b) $P(\text{four unaffected}) = (3/4)^4 = 81/256 \approx 0.3164$.

(3/4)^4 # 0.3164
dbinom(0, 4, 0.25) # P(0 affected) = same number

Exercise 5 — Poisson vs Negative Binomial

Set λ = 10, simulate 1000 draws per distribution:

library(dplyr); library(tidyr); library(ggplot2)
set.seed(2025); n <- 1000; lambda <- 10
sims <- tibble(
 pois = rpois(n, lambda),
 nb_k5 = rnbinom(n, mu=lambda, size=5),
 nb_k1 = rnbinom(n, mu=lambda, size=1),
 nb_kp5= rnbinom(n, mu=lambda, size=0.5))
sims %>% summarise(across(everything(),
 list(mean=mean, var=var, p0=~mean(.==0), pgt20=~mean(.>20))))

(d) Many-zero, occasional-large data is the NB k=0.5 pattern.

Exercise 6 — Marten survival

dbinom(0:5, 5, 0.7)
# 0.00243 0.02835 0.13230 0.30870 0.36015 0.16807

This is Bin(5, 0.7). Mode at $k=4$ (rounding $(n+1)p = 4.2$ down). Mean = 3.5, variance = $5 \cdot 0.7 \cdot 0.3 = 1.05$.

Exercise 7 — CWD surveillance

(a) $Y \sim \text{Bin}(24, 0.12)$.

dbinom(3, 24, 0.12) # 0.2387
pbinom(4, 24, 0.12) # 0.8471
dbinom(0, 24, 0.12) # 0.0465

(d) The 4.65% chance of zero detections despite 12% prevalence is the critical surveillance message: sampling too few animals fails to certify "free of disease" even when present. Real-world surveillance programs use this exact calculation to set minimum sample sizes for declarations of freedom from disease.

Exercise 8 — Elk multinomial

p <- c(0.20, 0.45, 0.25, 0.10)
dmultinom(c(18, 38, 18, 6), prob=p) # 0.00117
set.seed(1)
flights <- rmultinom(5000, size=80, prob=p)
hist(flights[3,], main="Calves observed", xlab="count")
mean(flights[3,]) # ~ 20 (= 80 * 0.25)
var(flights[3,]) # ~ 15 (= 80*0.25*0.75)

(c) Each marginal of a multinomial is binomial: calves ~ Bin(80, 0.25).

Exercise 9 — Likert simulation

p <- c(.30, .20, .15, .20, .15)
rmultinom(1, 80, p)
sims <- rmultinom(500, 80, p)
mean(sims[5,]); sd(sims[5,]) # ~ 12 and ~ 3.2

The standard deviation across replicate surveys (~3) is large relative to the point estimate (~12), illustrating the noise in a single n=80 survey.

Exercise 10 — Cystic fibrosis generalized

$K \sim \text{Bin}(n, 0.25)$.

(b) $P(K \ge 1) = 1 - (0.75)^5 = 1 - 0.2373 = 0.7627$.

(c) Smallest $n$ with $P(\ge 1) \ge 0.5$ is $n = 3$ (since $0.75^3 = 0.4219 < 0.5$ — meaning $P(\ge 1) = 1 - 0.4219 = 0.578$).

1 - dbinom(0, 5, 0.25) # 0.7627
which(1 - dbinom(0, 1:10, 0.25) >= 0.5)[1] # 3

Stretch 11 — Sum of Bernoullis = Binomial

(a) Linearity of expectation: $E[\sum Y_i] = \sum p = np$.

(b) Under independence, $\mathrm{Var}(\sum Y_i) = \sum \mathrm{Var}(Y_i) = np(1-p)$.

(c) Linearity of expectation does NOT require independence. Variance addition does. The covariance terms $2\sum_{i<j} \mathrm{Cov}(Y_i, Y_j)$ vanish only under independence. Positive correlation inflates variance — overdispersion.

Stretch 12 — Binomial → Poisson limit

Set $p = \lambda/n$. Then

The three pieces:

1. $\binom{n}{k}\left(\frac{\lambda}{n}\right)^k = \frac{n(n-1)\dots(n-k+1)}{k!\, n^k} \lambda^k \to \frac{\lambda^k}{k!}$ as $n \to \infty$.
2. $(1-\lambda/n)^n \to e^{-\lambda}$.
3. $(1-\lambda/n)^{-k} \to 1$.

Product: $\frac{\lambda^k e^{-\lambda}}{k!}$, the Poisson PMF.

Every case below holds $np = 4$ fixed while $n$ grows. The mean is pinned at 4 throughout; the variance is what moves:

# Exact variances: Binomial is np(1-p), Poisson is lambda
20*0.2*0.8 # 3.20
100*0.04*0.96 # 3.84
1000*0.004*0.996 # 3.984
# Pois(4) # 4.00

set.seed(1)
for(pp in list(c(20,0.2), c(100,0.04), c(1000,0.004))) {
 x <- rbinom(1e4, pp[1], pp[2])
 cat("n=", pp[1], " mean=", mean(x), " var=", var(x), "\n")
}
y <- rpois(1e4, 4); cat("Pois mean=", mean(y), " var=", var(y), "\n")
# n= 20   mean= 4.0015  var= 3.256023
# n= 100  mean= 4.0021  var= 3.825078
# n= 1000 mean= 3.9707  var= 3.917433
# Pois    mean= 4.0035  var= 3.993687

The Binomial variance $np(1-p) = \lambda(1-p)$ is always short of the Poisson variance $\lambda$ by exactly the factor $(1-p)$. That is the whole story: Bin(20, 0.2) has variance 3.20, a full 20% below 4, and a histogram visibly narrower than Pois(4). Push $p$ down to 0.04 and the shortfall is 4%; at $p = 0.004$ it is 0.4%. The Poisson is the limit reached only as $p \to 0$ with $\lambda = np$ held fixed — so the distributions are not interchangeable at $n = 20$, and become so only as $n$ grows.

The PMFs tell the same story. Comparing each Binomial to Pois(4) over $k = 0, \dots, 15$, the largest absolute difference in probability is:

k <- 0:15
max(abs(dbinom(k, 20, 0.2) - dpois(k, 4))) # 0.0228
max(abs(dbinom(k, 100, 0.04) - dpois(k, 4))) # 0.0040
max(abs(dbinom(k, 1000, 0.004) - dpois(k, 4))) # 0.00039

The error is $O(p) = O(\lambda/n)$, so asymptotically each tenfold increase in $n$ should buy a tenfold reduction in the gap. Watch the two steps separately: $0.02283 / 0.004022 = 5.7$, but $0.004022 / 0.0003918 = 10.3$. The second step delivers the promised factor of ten; the first only manages 5.7. The reason is that the $O(p)$ rate is an asymptotic statement, and $p = 0.2$ is not small — $n = 20$ is not yet in the regime where the rate applies. That is exactly why Bin(20, 0.2) is the case that shows you the approximation has to be earned.

Stretch 13 — Monty Hall (camera trap)

P(initial guess right) = 1/3 regardless of what is revealed. The technician's reveal of an empty camera moves the remaining 2/3 onto the unopened alternative. So switch.

The simulation has to actually model the technician and the switch — otherwise it just restates the analytic claim. Draw the truth and the guess, let the technician open a camera that is neither, then switch to whatever is left:

set.seed(1)
trials <- 1e5
truth <- sample(1:3, trials, replace=TRUE) # camera holding the marten
guess <- sample(1:3, trials, replace=TRUE) # your initial pick

# Technician opens a camera that is neither your pick nor the marten's
opened <- vapply(seq_len(trials), function(i) {
 choices <- setdiff(1:3, c(guess[i], truth[i]))
 if (length(choices) == 1) choices else sample(choices, 1)
}, numeric(1))

# You switch to the one remaining unopened camera
switched <- vapply(seq_len(trials), function(i) {
 setdiff(1:3, c(guess[i], opened[i]))
}, numeric(1))

mean(guess == truth) # stay wins: 0.33257
mean(switched == truth) # switch wins: 0.66743

Real output: staying wins 0.33257 of 100,000 trials, switching wins 0.66743 — matching 1/3 and 2/3.