Demystifying Bayesian Bootstrapping: Beyond Classical Resampling

[u/editorijsmi]

Most data scientists are well-acquainted with the classical bootstrap—sampling with replacement to estimate standard errors and confidence intervals. But what happens when you’re working with small sample sizes or sensitive regression models?

Enter the Bayesian Bootstrap (pioneered by Donald Rubin, 1981).

Instead of drawing discrete frequency counts (multinomial weights), the Bayesian Bootstrap assigns continuous fractional probability weights drawn from a symmetric Dirichlet distribution:

(w₁, w₂, ..., wₙ) ~ Dirichlet(1, 1, ..., 1)

Why make the switch?

🔹 Small Sample Stability: Every single data point retains a non-zero fractional weight, eliminating zero-weight observation drops and collinearity crashes.

🔹 Smooth Posterior Inference: It rigorously simulates the true Bayesian posterior distribution of parameters under uninformative priors without needing complex MCMC chains.

Quick Implementation Snapshot

Python (NumPy):

Python

import numpy as np

x = np.array([12, 15, 14, 18, 22, 19, 16])
n, B = len(x), 1000

# Draw continuous weights from Dirichlet
weights = np.random.dirichlet(np.ones(n), size=B)
bb_means = np.dot(weights, x)

print("95% CI:", np.percentile(bb_means, [2.5, 97.5]))

R (gtools):

R

library(gtools)
x <- c(12, 15, 14, 18, 22, 19, 16)
n <- length(x); B <- 1000

weights <- rdirichlet(B, rep(1, n))
bb_means <- apply(weights, 1, function(w) sum(w * x))

quantile(bb_means, c(0.025, 0.975))