Why would I bother bootstrapping instead of just using a standard formula?

In a perfect world, your data follows a neat bell curve and you can just use a Z-score or a T-test to find your confidence intervals. But in the real world—especially in software engineering data—things are messy. Take bug-fix durations, for example. Most bugs are fixed in a few hours, but a handful of "nightmare" bugs take three weeks. That's a heavily skewed distribution.

If you want the standard error of the mean, there's a formula for that. But what if you want the standard error of the median? There isn't a simple, one-size-fits-all formula for the median that works across all distribution types. That's where bootstrapping comes in. Instead of relying on a theoretical formula, we treat our sample as a miniature population and resample from it thousands of times. We're essentially simulating the process of taking new samples from the real world.

How do I actually implement a bootstrap loop in R?

You could use a for loop, but in R, replicate() is your best friend here. It's cleaner and tells anyone reading your code exactly what's happening: you're repeating an operation a specific number of times.

Let's say we have a small vector of days it took to close 10 tickets. I'll show you how to bootstrap the median of this set:

# Our skewed sample of bug-fix days
fix_times <- c(1, 2, 1, 1, 3, 2, 1, 12, 15, 2)

# We want to find the median 10,000 times
boot_medians <- replicate(10000, {
  # resample with replacement
  resample <- sample(fix_times, size = length(fix_times), replace = TRUE)
  median(resample)
})

# Let's see what the distribution of our bootstrapped medians looks like
hist(boot_medians, main = "Distribution of Bootstrapped Medians", xlab = "Median Days")

The key here is replace = TRUE. If you don't resample with replacement, you're just shuffling the same numbers, and your median will be the same every single time. That's a common mistake I see early on.

How do I turn these thousands of samples into a confidence interval?

Once you have your vector of bootstrapped statistics (like boot_medians above), you don't need any complex calculus to find the confidence interval. You just use the percentile method. If you want a 95% confidence interval, you just find the 2.5th percentile and the 97.5th percentile of your bootstrapped distribution.

I usually use the quantile() function for this. It's direct and hard to mess up:

# Calculate the 95% confidence interval
ci <- quantile(boot_medians, probs = c(0.025, 0.975))
print(ci)

I'll be honest: it feels like "cheating" the first time you do it. You're just picking the edges of a list of numbers. But that's the beauty of bootstrapping—it lets the data speak for itself without forcing it into a Gaussian mold that doesn't actually fit.

📋 Practical Task

Calculating Confidence Intervals for API Latency

You've been given a dataset representing the response times (in milliseconds) of a specific API endpoint over 20 requests. The data is highly irregular, and your manager wants a 95% confidence interval for the mean response time using the bootstrap method.

Dataset: latency_ms <- c(120, 135, 110, 150, 450, 125, 130, 115, 140, 110, 120, 130, 600, 110, 125, 130, 140, 115, 120, 130)

Your Task:

1. Use replicate() to create a bootstrap distribution of the mean of latency_ms with 10,000 iterations.
2. Ensure you are resampling with replacement and that each resample is the same size as the original dataset.
3. Calculate and print the 95% confidence interval using the quantile() function.