Bootstrap is a resampling method where we repeatedly draw samples from the same dataset (with replacement).
set.seed(123)
n <- 100
x <- rnorm(n, mean = 10, sd = 2)
y <- 3 + 1.5 * x + rnorm(n, mean = 0, sd = 2)
data <- data.frame(x, y)
data[sample(1:n, 10), "x"] <- NA
head(data)
## x y
## 1 8.879049 14.89776
## 2 9.539645 17.82323
## 3 13.117417 22.18274
## 4 10.141017 17.51644
## 5 10.258575 16.48463
## 6 13.430130 23.05514
Here we try to generate random data that follows normal distribution and linear
Bootstrap is applied to clean data by repeatedly resampling and fitting a regression model. It estimates the variability of coefficients (especially the slope) without relying on distributional assumptions. The resulting confidence interval shows how stable the relationship between x and y is.
clean_data <- na.omit(data)
boot_regression <- function(data, indices) {
d <- data[indices, ]
model <- lm(y ~ x, data = d)
return(coef(model))
}
library(boot)
boot_result <- boot(
data = clean_data,
statistic = boot_regression,
R = 1000
)
boot_result
##
## ORDINARY NONPARAMETRIC BOOTSTRAP
##
##
## Call:
## boot(data = clean_data, statistic = boot_regression, R = 1000)
##
##
## Bootstrap Statistics :
## original bias std. error
## t1* 3.581084 0.06067069 1.1482885
## t2* 1.412127 -0.00547455 0.1074228
mean_x <- mean(data$x, na.rm = TRUE)
data$ximp <- ifelse(is.na(data$x), mean_x, data$x)
model_imp <- lm(y ~ ximp, data = data)
summary(model_imp)
##
## Call:
## lm(formula = y ~ ximp, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1153 -1.4394 -0.0902 1.2053 6.5280
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.6538 1.2332 2.963 0.00383 **
## ximp 1.4121 0.1191 11.854 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.109 on 98 degrees of freedom
## Multiple R-squared: 0.5891, Adjusted R-squared: 0.5849
## F-statistic: 140.5 on 1 and 98 DF, p-value: < 2.2e-16
This regression shows that ximp has a strong and statistically significant positive effect on y (slope ≈ 1.41, p < 0.001), explaining about 59% of the variation in y, with reasonably small residual error.
boot_imp <- function(data, indices) {
d <- data[indices, ]
model <- lm(y ~ ximp, data = d)
return(coef(model))
}
boot_result_imp <- boot(data = data, statistic = boot_imp, R = 1000)
boot.ci(boot_result_imp, type = "perc", index = 2)
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 1000 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot_result_imp, type = "perc", index = 2)
##
## Intervals :
## Level Percentile
## 95% ( 1.188, 1.603 )
## Calculations and Intervals on Original Scale
library(mice)
## Warning: package 'mice' was built under R version 4.5.3
##
## Attaching package: 'mice'
## The following object is masked from 'package:stats':
##
## filter
## The following objects are masked from 'package:base':
##
## cbind, rbind
imp <- mice(
data[, c("x", "y")],
m = 5,
method = 'pmm',
seed = 123
)
##
## iter imp variable
## 1 1 x
## 1 2 x
## 1 3 x
## 1 4 x
## 1 5 x
## 2 1 x
## 2 2 x
## 2 3 x
## 2 4 x
## 2 5 x
## 3 1 x
## 3 2 x
## 3 3 x
## 3 4 x
## 3 5 x
## 4 1 x
## 4 2 x
## 4 3 x
## 4 4 x
## 4 5 x
## 5 1 x
## 5 2 x
## 5 3 x
## 5 4 x
## 5 5 x
model_mi <- with(imp, lm(y ~ x))
summary(pool(model_mi))
## term estimate std.error statistic df p.value
## 1 (Intercept) 3.619991 1.1112706 3.257524 78.99385 1.657655e-03
## 2 x 1.408248 0.1068028 13.185496 78.10532 1.472407e-21
The iteration table shows that MICE is repeatedly imputing missing values of x across 5 datasets and multiple iterations to stabilize the estimates, while the final result indicates that x has a strong, significant positive effect on y (slope ≈ 1.41, p < 0.001) with reliable standard error (~0.107).
MICE provides the most reliable estimates because it preserves variability, while mean imputation tends to bias the results despite producing similar coefficients.