Введение

В данной практической работе рассматривается метод бутстрапа на языках R и Python.

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

В работе используется встроенный датасет iris. Он содержит 150 наблюдений и 5 столбцов:

длина чашелистика;
ширина чашелистика;
длина лепестка;
ширина лепестка;
вид ириса.

Перед запуском отчета необходимые пакеты нужно один раз установить в консоли RStudio:

install.packages("ggplot2")
install.packages("dplyr")
install.packages("boot")
install.packages("reticulate")

1. Выбор датасета и разведочный анализ данных

Выполнение задания

# Загружаем библиотеки.
# ggplot2 нужен для построения графиков.
# dplyr нужен для удобной работы с данными.

library(ggplot2)
library(dplyr)

# Загружаем встроенный датасет iris.

data(iris)

# Смотрим первые строки датасета.
# Это помогает понять, какие столбцы есть в таблице.

head(iris)

##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa

# Выводим краткую статистику по всем столбцам.
# Для числовых признаков отображаются минимум, максимум, среднее, медиана и квартили.
# Для Species отображается количество объектов каждого вида.

summary(iris)

##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    :50  
##  versicolor:50  
##  virginica :50  
##                 
##                 
##

# Проверяем размер датасета.

cat("Размер датасета:", nrow(iris), "строк,", ncol(iris), "столбцов\n")

## Размер датасета: 150 строк, 5 столбцов

# Проверяем наличие пропущенных значений.
# is.na() находит пропуски, sum() считает их количество.

cat("Пропущенные значения:", sum(is.na(iris)), "\n")

## Пропущенные значения: 0

Визуализация данных

# Boxplot показывает распределение длины лепестка по видам ириса.
# По нему удобно сравнить, какие виды имеют большие или меньшие лепестки.

ggplot(iris, aes(x = Species, y = Petal.Length, fill = Species)) +
  geom_boxplot(alpha = 0.7) +
  labs(
    title = "Длина лепестка по видам ириса",
    x = "Вид",
    y = "Длина лепестка"
  ) +
  theme_minimal() +
  theme(legend.position = "none")

# Гистограмма показывает распределение длины лепестка.
# Цвет используется для разделения видов ириса.

ggplot(iris, aes(x = Petal.Length, fill = Species)) +
  geom_histogram(bins = 20, alpha = 0.6, position = "identity") +
  labs(
    title = "Гистограмма длины лепестка",
    x = "Длина лепестка",
    y = "Частота"
  ) +
  theme_minimal()

# Диаграмма рассеяния показывает связь между длиной и шириной лепестка.
# Каждая точка — отдельный цветок.
# Цвет точки показывает вид ириса.

ggplot(iris, aes(x = Petal.Length, y = Petal.Width, color = Species)) +
  geom_point(size = 2, alpha = 0.8) +
  labs(
    title = "Связь между длиной и шириной лепестка",
    x = "Длина лепестка",
    y = "Ширина лепестка"
  ) +
  theme_minimal()

Вывод по заданию 1

Для выполнения практической работы был выбран датасет iris.

Датасет содержит 150 наблюдений и 5 столбцов. Пропущенных значений в данных нет, поэтому дополнительная очистка от NA не требуется.

По графику boxplot видно, что вид setosa заметно отличается от остальных видов по длине лепестка. Виды versicolor и virginica расположены ближе друг к другу, но у virginica лепестки в среднем длиннее.

Гистограмма подтверждает, что распределения длины лепестка у разных видов отличаются. Диаграмма рассеяния показывает положительную связь между длиной и шириной лепестка: при увеличении длины лепестка обычно увеличивается и его ширина.

2. Бутстрап на языке Python

Выполнение задания

В данной части выполняется бутстрап на Python.

Задача: сравнить 90-й перцентиль длины лепестка у двух видов ириса:

versicolor;
virginica.

Считается разница:

90-й перцентиль virginica - 90-й перцентиль versicolor

Если доверительный интервал не включает 0, значит различие можно считать статистически значимым.

# Подключаем reticulate.
# Этот пакет позволяет запускать Python-код внутри RMarkdown.

library(reticulate)

# Импортируем библиотеки Python.
# numpy нужен для вычислений и бутстрап-выборок.
# pandas нужен для работы с таблицей.
# matplotlib нужен для построения графиков.
# load_iris нужен для загрузки датасета iris из sklearn.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.datasets import load_iris

# Загружаем датасет iris из sklearn.

iris_data = load_iris()

# Создаем таблицу DataFrame.
# В iris_data.data лежат числовые признаки.
# В iris_data.feature_names лежат названия признаков.

df = pd.DataFrame(
    iris_data.data,
    columns = iris_data.feature_names
)

# Добавляем столбец species.
# 0 = setosa, 1 = versicolor, 2 = virginica.

df["species"] = iris_data.target

# Выводим первые строки таблицы.

print(df.head())

##    sepal length (cm)  sepal width (cm)  ...  petal width (cm)  species
## 0                5.1               3.5  ...               0.2        0
## 1                4.9               3.0  ...               0.2        0
## 2                4.7               3.2  ...               0.2        0
## 3                4.6               3.1  ...               0.2        0
## 4                5.0               3.6  ...               0.2        0
## 
## [5 rows x 5 columns]

# Выбираем длину лепестка для versicolor.
# species == 1 соответствует versicolor.

values_a = df[df["species"] == 1]["petal length (cm)"].values

# Выбираем длину лепестка для virginica.
# species == 2 соответствует virginica.

values_b = df[df["species"] == 2]["petal length (cm)"].values

print(f"Versicolor: {len(values_a)} наблюдений, среднее = {values_a.mean():.2f}")

## Versicolor: 50 наблюдений, среднее = 4.26

print(f"Virginica: {len(values_b)} наблюдений, среднее = {values_b.mean():.2f}")

## Virginica: 50 наблюдений, среднее = 5.55

# Фиксируем случайность, чтобы результат можно было повторить.

np.random.seed(123)

# B — количество бутстрап-итераций.
# Чем больше B, тем стабильнее оценка.

B = 10000

# alpha = 0.05 соответствует 95% доверительному интервалу.

alpha = 0.05

# Точечная оценка:
# разница 90-х перцентилей в исходных данных.

point_estimate = np.quantile(values_b, 0.9) - np.quantile(values_a, 0.9)

# Создаем бутстрап-выборки для versicolor.
# np.random.choice берет элементы с возвращением.
# Размер каждой выборки равен размеру исходной группы.

bootstrap_a = np.quantile(
    np.random.choice(values_a, (B, len(values_a)), replace = True),
    0.9,
    axis = 1
)

# Создаем бутстрап-выборки для virginica.

bootstrap_b = np.quantile(
    np.random.choice(values_b, (B, len(values_b)), replace = True),
    0.9,
    axis = 1
)

# Для каждой итерации считаем разницу 90-х перцентилей.

bootstrap_diff = bootstrap_b - bootstrap_a

# Строим 95% доверительный интервал.
# Берем 2.5% и 97.5% квантили бутстрап-распределения.

ci_low, ci_high = np.quantile(bootstrap_diff, [alpha / 2, 1 - alpha / 2])

# Проверяем, входит ли 0 в доверительный интервал.
# Если 0 не входит, различие статистически значимо.

has_effect = not (ci_low < 0 < ci_high)

print("Разница 90-х перцентилей:", round(point_estimate, 3))

## Разница 90-х перцентилей: 1.51

print("95% доверительный интервал:", round(ci_low, 3), "-", round(ci_high, 3))

## 95% доверительный интервал: 1.11 - 1.99

print("Различие статистически значимо:", has_effect)

## Различие статистически значимо: True

Визуализация бутстрапа на Python

# Строим гистограмму бутстрап-распределения.
# Она показывает, какие значения разницы чаще всего получались
# при многократном создании выборок.

plt.figure(figsize = (10, 6))

plt.hist(
    bootstrap_diff,
    bins = 50,
    color = "skyblue",
    edgecolor = "black",
    alpha = 0.7
)

# Красная линия — точечная оценка.

plt.axvline(
    point_estimate,
    color = "red",
    linewidth = 2,
    label = "Точечная оценка"
)

# Зеленые пунктирные линии — границы доверительного интервала.

plt.axvline(
    ci_low,
    color = "green",
    linestyle = "--",
    linewidth = 2,
    label = "95% ДИ"
)

plt.axvline(
    ci_high,
    color = "green",
    linestyle = "--",
    linewidth = 2
)

# Черная пунктирная линия — нулевое значение.
# Если оно не попадает внутрь интервала, различие значимо.

plt.axvline(
    0,
    color = "black",
    linestyle = ":",
    linewidth = 2,
    label = "0"
)

plt.title("Бутстрап-распределение разницы 90-х перцентилей")
plt.xlabel("Virginica - Versicolor")
plt.ylabel("Частота")
plt.legend()
plt.grid(alpha = 0.3)

plt.show()

Вывод по заданию 2

В задании был реализован бутстрап на языке Python.

Сравнивался 90-й перцентиль длины лепестка у видов versicolor и virginica. Сначала была рассчитана точечная оценка разницы 90-х перцентилей, затем были созданы бутстрап-выборки с возвращением.

По результатам бутстрапа был построен 95% доверительный интервал. Доверительный интервал не включает 0, поэтому можно сделать вывод, что различие между 90-ми перцентилями статистически значимо.

Это означает, что у вида virginica верхняя часть распределения длины лепестка заметно больше, чем у versicolor.

3. Бутстрап на языке R

Выполнение задания

В данной части выполняется бутстрап на R.

Задача: оценить устойчивость коэффициента детерминации R² для модели:

Petal.Length ~ Petal.Width

То есть проверяется, насколько хорошо ширина лепестка объясняет длину лепестка.

# Подключаем пакет boot.
# Он используется для выполнения бутстрапа на R.

library(boot)

# ggplot2 уже подключался выше, но можно подключить повторно для ясности.

library(ggplot2)

# Создаем функцию для бутстрапа.
# data — исходные данные.
# indices — индексы строк, которые попали в бутстрап-выборку.
#
# Внутри функции:
# 1. Создается бутстрап-выборка d.
# 2. Строится линейная модель Petal.Length ~ Petal.Width.
# 3. Возвращается R² этой модели.

rsq_function <- function(data, indices) {
  d <- data[indices, ]

  fit <- lm(Petal.Length ~ Petal.Width, data = d)

  summary(fit)$r.square
}

# Фиксируем случайность для воспроизводимости результата.

set.seed(42)

# Запускаем бутстрап.
# R = 2000 означает 2000 бутстрап-повторений.

reps <- boot(
  data = iris,
  statistic = rsq_function,
  R = 2000
)

# Выводим результат бутстрапа.

reps

## 
## ORDINARY NONPARAMETRIC BOOTSTRAP
## 
## 
## Call:
## boot(data = iris, statistic = rsq_function, R = 2000)
## 
## 
## Bootstrap Statistics :
##      original       bias    std. error
## t1* 0.9271098 0.0005393181  0.00994126

# Строим доверительный интервал методом BCa.
# BCa — это уточненный бутстрап-доверительный интервал.

ci_result <- boot.ci(
  reps,
  type = "bca"
)

ci_result

## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 2000 bootstrap replicates
## 
## CALL : 
## boot.ci(boot.out = reps, type = "bca")
## 
## Intervals : 
## Level       BCa          
## 95%   ( 0.9035,  0.9435 )  
## Calculations and Intervals on Original Scale

# Сохраняем точечную оценку R².
# reps$t0 — значение статистики на исходной выборке.

point_r2 <- reps$t0

# Извлекаем границы доверительного интервала.

ci_low_r <- ci_result$bca[4]
ci_high_r <- ci_result$bca[5]

cat("Оригинальное R²:", round(point_r2, 4), "\n")

## Оригинальное R²: 0.9271

cat("95% ДИ (BCa):", round(ci_low_r, 4), "-", round(ci_high_r, 4), "\n")

## 95% ДИ (BCa): 0.9035 - 0.9435

Визуализация бутстрапа на R

# Преобразуем бутстрап-оценки R² в data.frame для ggplot.

boot_df <- data.frame(
  r_squared = reps$t
)

# Строим гистограмму бутстрап-распределения R².

ggplot(boot_df, aes(x = r_squared)) +
  geom_histogram(
    bins = 50,
    fill = "skyblue",
    color = "black",
    alpha = 0.7
  ) +
  geom_vline(
    xintercept = point_r2,
    color = "red",
    linewidth = 1
  ) +
  geom_vline(
    xintercept = ci_low_r,
    color = "green",
    linetype = "dashed",
    linewidth = 1
  ) +
  geom_vline(
    xintercept = ci_high_r,
    color = "green",
    linetype = "dashed",
    linewidth = 1
  ) +
  labs(
    title = "Бутстрап-распределение R²",
    subtitle = paste0(
      "95% ДИ: [",
      round(ci_low_r, 3),
      "; ",
      round(ci_high_r, 3),
      "]"
    ),
    x = "R²",
    y = "Частота"
  ) +
  theme_minimal()

# Дополнительно построим исходные данные с линией линейной регрессии.
# Это показывает связь между шириной и длиной лепестка.

ggplot(iris, aes(x = Petal.Width, y = Petal.Length, color = Species)) +
  geom_point(alpha = 0.7, size = 2) +
  geom_smooth(
    method = "lm",
    se = TRUE,
    color = "black"
  ) +
  labs(
    title = "Регрессия: длина лепестка от ширины лепестка",
    subtitle = paste("R² =", round(point_r2, 3)),
    x = "Ширина лепестка",
    y = "Длина лепестка"
  ) +
  theme_minimal()

Вывод по заданию 3

В задании был реализован бутстрап на языке R.

В качестве статистики использовался коэффициент детерминации R² для линейной модели Petal.Length ~ Petal.Width. Этот коэффициент показывает, какая доля изменения длины лепестка объясняется шириной лепестка.

Полученное значение R² близко к 0.93. Это означает, что ширина лепестка объясняет около 93% вариации длины лепестка.

Доверительный интервал для R² получился узким и находится выше 0.9. Это говорит о том, что оценка достаточно стабильна, а связь между шириной и длиной лепестка сильная.

Общий вывод

В практической работе был рассмотрен метод бутстрапа на языках Python и R.

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

В Python был выполнен бутстрап для сравнения 90-х перцентилей длины лепестка у видов virginica и versicolor. Доверительный интервал не включил 0, поэтому различие можно считать статистически значимым.

В R был выполнен бутстрап для коэффициента детерминации R² линейной модели Petal.Length ~ Petal.Width. Результат показал, что связь между шириной и длиной лепестка сильная, а модель достаточно устойчива.

Таким образом, бутстрап позволяет оценивать устойчивость статистик и строить доверительные интервалы без строгих предположений о распределении данных.

Практическая работа №6: Бутстрап на R и Python

Папоян Гагик Карапетович 231-333

2026-05-30

Введение

1. Выбор датасета и разведочный анализ данных

Выполнение задания

Визуализация данных

Вывод по заданию 1

2. Бутстрап на языке Python

Выполнение задания

Визуализация бутстрапа на Python

Вывод по заданию 2

3. Бутстрап на языке R

Выполнение задания

Визуализация бутстрапа на R

Вывод по заданию 3

Общий вывод