The Law of Large Numbers
What the theorem actually says
Draw observations one at a time from a population with mean \(\mu\) and keep a running average \(\bar{X}_n = \frac{1}{n}\sum_{i=1}^n X_i\). The Law of Large Numbers (LLN) says this running average converges to \(\mu\) as \(n\) grows: collect enough data and your sample mean gets arbitrarily close to the truth, and stays there. \[\bar{X}_n \;\xrightarrow{\;p\;}\; \mu \qquad \text{as } n \to \infty.\]
All it asks of the population is a finite mean (\(E|X| < \infty\)). The shape does not matter — skewed, bimodal, discrete, all converge. This is the theorem that makes estimation possible at all: it is the guarantee that more data means a better estimate.
Weak vs. strong
There are two versions, differing only in how the average converges:
| Weak LLN | Strong LLN | |
|---|---|---|
| Statement | \(\bar{X}_n \xrightarrow{p} \mu\) | \(\bar{X}_n \xrightarrow{a.s.} \mu\) |
| In words | For any tolerance, the chance the average is off by more than that shrinks to zero | The sequence of averages settles on \(\mu\) with probability 1 and does not wander back out |
| Due to | Bernoulli, Chebyshev | Kolmogorov (1933) |
For applied work the distinction rarely bites — both say “the average lands on the truth.” The strong version just promises it lands once and stays, rather than merely being probably close at each large \(n\).
Watch it converge
Add observations one at a time and watch the running average settle onto the true mean (the horizontal line). This is the LLN in action — the cumulative mean stabilizes regardless of the population shape. Draw several paths at once and they all funnel to the same \(\mu\) by different routes.
#| standalone: true
#| viewerHeight: 620
library(shiny)
ui <- fluidPage(
tags$head(tags$style(HTML("
.stats-box {
background: #f0f4f8; border-radius: 6px; padding: 14px;
margin-top: 12px; font-size: 14px; line-height: 1.9;
}
.stats-box b { color: #2c3e50; }
"))),
sidebarLayout(
sidebarPanel(
width = 3,
selectInput("pop", "Population:",
choices = c("Uniform(0, 1)",
"Exponential(1)",
"Bimodal",
"Bernoulli(0.3)",
"Heavy-tailed (t, df=2)")),
sliderInput("n", "Number of observations:",
min = 50, max = 5000, value = 500, step = 50),
sliderInput("paths", "Number of paths:",
min = 1, max = 10, value = 3, step = 1),
uiOutput("results")
),
mainPanel(
width = 9,
plotOutput("running_avg", height = "500px")
)
)
)
server <- function(input, output, session) {
draw_pop <- function(n, pop) {
switch(pop,
"Uniform(0, 1)" = runif(n),
"Exponential(1)" = rexp(n, rate = 1),
"Bimodal" = {
k <- rbinom(n, 1, 0.5)
k * rnorm(n, -2, 0.6) + (1 - k) * rnorm(n, 2, 0.6)
},
"Bernoulli(0.3)" = rbinom(n, 1, 0.3),
"Heavy-tailed (t, df=2)" = rt(n, df = 2)
)
}
pop_mu <- function(pop) {
switch(pop,
"Uniform(0, 1)" = 0.5,
"Exponential(1)" = 1,
"Bimodal" = 0,
"Bernoulli(0.3)" = 0.3,
"Heavy-tailed (t, df=2)" = 0
)
}
dat <- reactive({
n <- input$n
paths <- input$paths
pop <- input$pop
mu <- pop_mu(pop)
all_paths <- lapply(seq_len(paths), function(i) {
x <- draw_pop(n, pop)
cumsum(x) / seq_along(x)
})
list(all_paths = all_paths, n = n, mu = mu, pop = pop, paths = paths)
})
output$running_avg <- renderPlot({
d <- dat()
par(mar = c(4.5, 4.5, 3, 1))
ylim <- range(unlist(d$all_paths))
ylim <- ylim + c(-1, 1) * 0.1 * diff(ylim)
cols <- c("#e74c3c", "#3498db", "#27ae60", "#9b59b6", "#e67e22",
"#1abc9c", "#34495e", "#f39c12", "#2ecc71", "#c0392b")
plot(NULL, xlim = c(1, d$n), ylim = ylim,
xlab = "Number of observations",
ylab = "Cumulative mean",
main = paste0("LLN: Running average converges to μ = ", d$mu))
abline(h = d$mu, lty = 2, lwd = 2.5, col = "#2c3e50")
for (i in seq_along(d$all_paths)) {
lines(seq_along(d$all_paths[[i]]), d$all_paths[[i]],
col = adjustcolor(cols[i], 0.7), lwd = 1.5)
}
legend("topright", bty = "n", cex = 0.9,
legend = c(paste0("True mean (μ = ", d$mu, ")"),
paste0(d$paths, " sample path(s)")),
col = c("#2c3e50", cols[1]),
lwd = c(2.5, 1.5), lty = c(2, 1))
})
output$results <- renderUI({
d <- dat()
final_means <- sapply(d$all_paths, function(p) p[length(p)])
tags$div(class = "stats-box",
HTML(paste0(
"<b>True mean:</b> ", d$mu, "<br>",
"<b>Final running avg(s):</b><br>",
paste(round(final_means, 4), collapse = ", "), "<br>",
"<b>Max deviation:</b> ",
round(max(abs(final_means - d$mu)), 4)
))
)
})
}
shinyApp(ui, server)What the LLN does not promise
- It says nothing about how uncertain your estimate is. That averages converge is the LLN; how they scatter around the truth — the bell curve, the standard error \(\sigma/\sqrt{n}\) — is the Central Limit Theorem. You need the LLN first: if the mean hasn’t stabilized, knowing its distribution is normal doesn’t help.
- It is not the gambler’s fallacy. A long run of heads is not “corrected” by future tails. The average dilutes early deviations across ever more observations; it does not reverse them.
- It needs a finite mean. The Cauchy distribution (a \(t\) with 1 degree of freedom) has no mean at all — its running average never settles, no matter how much data you collect. Try Heavy-tailed (t, df=2) above and then imagine going heavier: at df=1 the path never converges.
From here
- See the two theorems run against each other — consistency vs. asymptotic normality — on LLN vs CLT.
- The distribution the average settles into is the Central Limit Theorem.
Did you know?
- Jacob Bernoulli proved the first version of the Law of Large Numbers in Ars Conjectandi, published posthumously in 1713. He called it his “golden theorem” and spent over 20 years on the proof. The result seems obvious in hindsight — of course averages converge — but rigorously proving why required entirely new mathematical machinery.
- The Strong LLN (almost sure convergence) was proved by Kolmogorov in 1933, more than two centuries after Bernoulli’s weak version. It says convergence happens with probability 1, not just “usually.”
- The Cauchy distribution (\(t\) with \(\text{df}=1\)) has no finite mean, so the LLN does not apply. The running average never settles down, no matter how much data you collect — it keeps jumping around forever.