R Academy · Урок

Бутстреп-ресэмплинг в R

Реализуйте доверительные интервалы методом бутстрепа с помощью sample() с возвращением

Урок 4 из 413 шагов

«Бутстреп-ресэмплинг в R» — бесплатный урок R Academy на CoddyKit. Это урок 4 из 4. Ты можешь прочитать весь урок бесплатно ниже — а потом практиковать его прямо в браузере с встроенным редактором кода и ИИ-репетитором 24/7. Это часть пути обучения R Academy, и твой прогресс синхронизируется между веб-версией и приложением CoddyKit. Курс R Academy содержит 4 уроков всего.

Что такое бутстрэп

Бутстрэп — это метод повторной выборки, который оценивает выборочное распределение статистики, выполняя выборку с возвращением из наблюдаемых данных. Он не требует предположений о распределении: сами данные служат распределением.

set.seed(42)
# Original sample of 20 observations
original_data <- c(14, 18, 11, 13, 6, 8, 2, 17,
                   20, 10, 15, 12, 9, 16, 7, 3,
                   19, 5, 4, 1)
n <- length(original_data)

# One bootstrap resample: sample WITH replacement
bootstrap_sample1 <- sample(original_data, n, replace = TRUE)
print(bootstrap_sample1)
# Contains repeats! Some original values are missing.

mean(original_data)    # original mean
mean(bootstrap_sample1)  # bootstrap mean (different)

Бутстрэп-распределение среднего

Сгенерируйте множество бутстрэп-выборок, вычислите статистику для каждой и соберите результаты. Это бутстрэп-распределение приближает выборочное распределение статистики.

set.seed(42)
x <- c(14, 18, 11, 13, 6, 8, 2, 17, 20, 10,
        15, 12, 9, 16, 7, 3, 19, 5, 4, 1)
n <- length(x)
B <- 5000  # number of bootstrap replicates

# Bootstrap distribution of the mean
boot_means <- replicate(B, {
  x_star <- sample(x, n, replace = TRUE)
  mean(x_star)
})

cat('Original mean:', mean(x), '\n')
cat('Bootstrap mean:', mean(boot_means), '\n')
cat('Bootstrap SE:', sd(boot_means), '\n')
cat('Theoretical SE:', sd(x)/sqrt(n))

Бутстрэп для любой статистики

Бутстрэп работает с любой статистикой — медианой, корреляцией, коэффициентами регрессии, усечённым средним. Просто замените mean() на интересующую Вас статистику в цикле повторений.

set.seed(42)
x <- c(2, 4, 6, 8, 10, 12, 100, 1, 3, 5)
# Notice the outlier 100 - makes the mean unreliable
B <- 5000; n <- length(x)

# Bootstrap distribution of the median
boot_medians <- replicate(B, {
  median(sample(x, n, replace = TRUE))
})

# Bootstrap distribution of the trimmed mean (10%)
boot_tmeans <- replicate(B, {
  mean(sample(x, n, replace = TRUE), trim = 0.1)
})

cat('Sample median:', median(x), '\n')
cat('Bootstrap SE of median:', sd(boot_medians), '\n')
cat('Bootstrap SE of trimmed mean:', sd(boot_tmeans))

Бутстрэп-доверительные интервалы: метод процентилей

Простейший бутстрэп-ДИ: возьмите квантили бутстрэп-распределения. Для 95%-го ДИ используются 2,5-й и 97,5-й процентили. Метод быстрый и интуитивно понятный, но для асимметричных статистик покрытие может быть неточным.

set.seed(42)
x <- rnorm(30, mean = 5, sd = 2)
B <- 5000; n <- length(x)

# Bootstrap CI for the mean
boot_means <- replicate(B, mean(sample(x, n, replace = TRUE)))

# Percentile CI
ci_percentile <- quantile(boot_means, c(0.025, 0.975))
cat('Percentile CI:', round(ci_percentile, 3), '\n')

# For comparison: t-interval (parametric)
t_ci <- mean(x) + qt(c(0.025, 0.975), df = n-1) * sd(x)/sqrt(n)
cat('t-interval CI:', round(t_ci, 3), '\n')

cat('True mean: 5')

Стандартная ошибка по бутстрэпу

Бутстрэп-оценка стандартной ошибки — это просто sd(boot_statistics). Она оценивает, насколько изменялась бы статистика в повторных выборках из генеральной совокупности, не требуя предположений о распределении.

set.seed(42)
# Bootstrap SE for the ratio of two means
group_a <- rnorm(20, mean = 10, sd = 2)
group_b <- rnorm(20, mean = 8, sd = 2)
all_data <- data.frame(
  value = c(group_a, group_b),
  group = rep(c('A', 'B'), each = 20)
)

B <- 4000; n <- nrow(all_data)
boot_ratio <- replicate(B, {
  idx <- sample(1:n, n, replace = TRUE)
  d <- all_data[idx, ]
  mean(d$value[d$group == 'A']) /
    mean(d$value[d$group == 'B'])
})

cat('Observed ratio:', mean(group_a)/mean(group_b), '\n')
cat('Bootstrap SE:', round(sd(boot_ratio), 4), '\n')
cat('95% CI:', round(quantile(boot_ratio, c(.025,.975)), 3))

Доверительный интервал BCa

ДИ BCa (с поправками на смещение и ускорение) точнее процентильного ДИ, особенно для асимметричных или смещённых статистик. Он учитывает смещение и асимметрию бутстрэп-распределения.

# BCa CI implementation
bca_ci <- function(data, stat_fn, B = 2000, alpha = 0.05) {
  n <- length(data)
  theta_hat <- stat_fn(data)
  # Bootstrap replicates
  boot_vals <- replicate(B, stat_fn(sample(data, n, replace = TRUE)))
  # Bias correction z0
  z0 <- qnorm(mean(boot_vals < theta_hat))
  # Acceleration a (jackknife)
  jk <- sapply(seq_len(n), function(i) stat_fn(data[-i]))
  num <- sum((mean(jk) - jk)^3)
  den <- 6 * (sum((mean(jk) - jk)^2))^(3/2)
  a <- num / den
  # Adjusted quantiles
  z_alpha <- qnorm(c(alpha/2, 1 - alpha/2))
  p_adj <- pnorm(z0 + (z0 + z_alpha) / (1 - a * (z0 + z_alpha)))
  quantile(boot_vals, p_adj)
}

set.seed(42)
x <- rexp(25, rate = 0.5)  # skewed distribution
bca_ci(x, median)

Бутстрэп для регрессии

Стройте бутстрэп-ДИ для регрессии, выполняя повторную выборку строк таблицы данных. Это учитывает ненормальные ошибки и позволяет делать корректный статистический вывод без предположения о нормальности.

set.seed(42)
n <- 40
x <- runif(n, 0, 10)
y <- 2 + 1.5 * x + rnorm(n, sd = 2)
df <- data.frame(x = x, y = y)
B <- 3000

# Bootstrap the slope coefficient
boot_slopes <- replicate(B, {
  idx <- sample(1:n, n, replace = TRUE)
  d_boot <- df[idx, ]
  coef(lm(y ~ x, data = d_boot))[2]
})

cat('OLS slope:', round(coef(lm(y~x, data=df))[2], 3), '\n')
cat('Bootstrap SE of slope:', round(sd(boot_slopes), 4), '\n')
cat('Bootstrap 95% CI:', round(quantile(boot_slopes, c(.025,.975)), 3))

Число бутстрэп-повторений

Сколько бутстрэп-повторений B Вам нужно? Для стандартных ошибок: B=200–500. Для процентильных ДИ: B=1000–2000. Для ДИ BCa: B=2000–5000. Большее число обычно лучше, но после 10 000 прирост качества уменьшается.

set.seed(42)
x <- rnorm(50, mean = 3, sd = 1)
B_values <- c(100, 500, 1000, 2000, 5000, 10000)

# See how CI width varies with B
results <- sapply(B_values, function(B) {
  boot_m <- replicate(B, mean(sample(x, length(x), replace = TRUE)))
  q <- quantile(boot_m, c(0.025, 0.975))
  q[2] - q[1]  # CI width
})

result_df <- data.frame(B = B_values, ci_width = round(results, 4))
print(result_df)
# CI width stabilizes as B increases

Парный бутстрэп

В задачах с двумя выборками повторно выбирайте пары (x, y) вместе, а не x и y независимо. Это сохраняет структуру зависимости внутри пар и даёт корректные ДИ для корреляции и парных разностей.

set.seed(42)
# Paired data: before/after treatment
before <- c(10, 12, 8, 15, 11, 9, 14, 13)
after  <- c(12, 14, 9, 16, 11, 10, 15, 12)
n <- length(before)
B <- 4000

# Bootstrap paired mean difference
boot_diff <- replicate(B, {
  idx <- sample(1:n, n, replace = TRUE)
  mean(after[idx] - before[idx])
})

observed_diff <- mean(after - before)
cat('Observed mean diff:', observed_diff, '\n')
cat('Bootstrap SE:', round(sd(boot_diff), 4), '\n')
cat('95% CI:', round(quantile(boot_diff, c(0.025, 0.975)), 3))

Бутстрэп-тест гипотезы

Бутстрэп можно использовать и для проверки гипотез. Чтобы проверить H₀: θ=θ₀, получите распределение при нулевой гипотезе, сдвигая или перецентрируя бутстрэп-выборки, а затем вычислите p-значение как долю значений, превышающих наблюдаемую статистику.

set.seed(42)
# Test H0: mean = 0 vs H1: mean != 0
x <- rnorm(30, mean = 0.5, sd = 2)
theta0 <- 0  # null hypothesis
observed_t <- (mean(x) - theta0) / (sd(x) / sqrt(length(x)))

B <- 5000; n <- length(x)
# Shift x to have mean = theta0 under H0
x_centered <- x - mean(x) + theta0

boot_t <- replicate(B, {
  x_star <- sample(x_centered, n, replace = TRUE)
  (mean(x_star) - theta0) / (sd(x_star) / sqrt(n))
})

p_value <- mean(abs(boot_t) >= abs(observed_t))
cat('Observed t:', round(observed_t, 3), '\n')
cat('Bootstrap p-value:', round(p_value, 4))

Обзор пакета boot

Пакет boot предоставляет boot(data, statistic, R) с расширенными возможностями: ДИ BCa с помощью boot.ci(), параллельные вычисления, стратифицированную выборку. Функция statistic получает данные и вектор индексов.

# library(boot)  # uncomment to use

# The boot() statistic function signature:
# statistic(data, indices) -> scalar or vector

# Example (conceptual - requires boot package):
# mean_stat <- function(data, indices) {
#   mean(data[indices])
# }
# results <- boot(data = x, statistic = mean_stat, R = 5000)
# boot.ci(results, type = c('perc', 'bca'))  # CIs
# plot(results)  # histogram of bootstrap distribution

# Manual equivalent:
set.seed(42)
x <- rnorm(30, mean = 5)
B <- 2000
boot_vals <- replicate(B, mean(sample(x, length(x), replace = TRUE)))
quantile(boot_vals, c(0.025, 0.975))

Быстрая проверка

Проверьте своё понимание повторной выборки методом бутстрэпа в R.

Итоги: повторная выборка методом бутстрэпа

Главные выводы: бутстрэп выполняет выборку данных с возвращением с помощью sample(x, n, replace=TRUE). Используйте replicate(B, ...), чтобы получить B бутстрэп-статистик. sd(boot_stats) возвращает бутстрэп-оценку стандартной ошибки. Процентильный ДИ: quantile(boot_stats, c(0.025, 0.975)). Для ДИ используйте B≥2000. Пакет boot добавляет ДИ BCa и параллельные вычисления. Бутстрэп работает с любой статистикой без предположений о распределении.

set.seed(42)
x <- c(3, 5, 7, 2, 9, 4, 6, 8, 1, 10)
n <- length(x); B <- 3000

# Generic bootstrap template:
boot_stat <- replicate(B, {
  x_star <- sample(x, n, replace = TRUE)
  median(x_star)  # replace with any statistic
})

# Summary
cat('Original median:', median(x), '\n')
cat('Bootstrap SE:', round(sd(boot_stat), 3), '\n')
ci <- quantile(boot_stat, c(0.025, 0.975))
cat('95% Percentile CI: [', ci[1], ',', ci[2], ']')
Можно начать бесплатно

Изучай R с ИИ-репетитором — бесплатно

Пиши и запускай код прямо в браузере, получай мгновенную помощь от ИИ-репетитора 24/7 и продолжи учиться на сайте или в приложении.

Курсы
43
Уроки
159

Часто задаваемые вопросы

Урок «Бутстреп-ресэмплинг в R» бесплатный?

Да — полный текст урока «Бутстреп-ресэмплинг в R» бесплатно доступен здесь в веб-версии. Чтобы практиковать его интерактивно (встроенный редактор кода и ИИ-репетитор 24/7) и разблокировать остальной курс R Academy, подпишись на CoddyKit PRO. Курс R Academy содержит 4 уроков всего.

Чему я научусь в уроке «Бутстреп-ресэмплинг в R»?

Реализуйте доверительные интервалы методом бутстрепа с помощью sample() с возвращением Ты практикуешь R Academy с помощью реального кода, который запускаешь прямо в браузере, и ИИ-репетитор 24/7 отвечает на твои вопросы во время урока.

Нужен ли мне опыт, чтобы начать R Academy?

Предыдущий опыт не требуется. R Academy на CoddyKit структурирован для всех уровней — от новичков до продвинутых, поэтому ты можешь начать отсюда или с самого начала и учиться в своем темпе. Это урок 4 из 4.

Сколько времени занимает урок «Бутстреп-ресэмплинг в R»?

Большинство уроков CoddyKit занимают около 5–10 минут. Каждый из них компактный и интерактивный, поэтому ты постоянно делаешь прогресс и продолжаешь с того же места в веб-версии и приложении.

Можно ли писать и запускать код в этом уроке R Academy?

Да. Каждый урок R Academy включает встроенный редактор кода, поэтому ты пишешь и запускаешь реальный код прямо в браузере и получаешь моментальную обратную связь от AI — локальная установка не требуется.

Все уроки этого курса

  1. set.seed() и воспроизводимость
  2. Генерация случайных распределений
  3. Основы метода Монте-Карло
  4. Бутстреп-ресэмплинг в R
← Назад к R Academy