Бутстреп-методи для стандартного відхилення

Бутстреп: статистична революція комп’ютерної ери

Бутстреп-ресемплінг — це потужний статистичний метод, який дозволяє оцінити вибірковий розподіл будь-якої статистики шляхом багаторазового повторного відбору зі спостережуваних даних. Запропонований Бредлі Ефроном у 1979 році, він здійснив революцію в статистичному висновуванні, даючи змогу аналізувати складні статистики без опори на математичні формули чи припущення щодо розподілу.

Ключова ідея бутстрепу елегантно проста: ваша вибірка є найкращою оцінкою генеральної сукупності. Повторно відбираючи зі своєї вибірки (з поверненням), ви імітуєте те, що було б, якби можна було багаторазово вибирати з генеральної сукупності. Цей підхід особливо цінний для стандартного відхилення, оскільки традиційні формули довірчих інтервалів припускають нормальність — припущення, яке на практиці часто не виконується.

Бутстреп став незамінним у сучасній науці про дані, оскільки працює з будь-якою статистикою (медіана, кореляція, коефіцієнти регресії, ваги нейронної мережі) і не висуває жодних припущень щодо розподілу ваших даних.

Навіщо бутстреп для стандартного відхилення?

Традиційні довірчі інтервали для стандартного відхилення передбачають, що дані походять із нормального розподілу. Коли це припущення порушується (а це трапляється часто), такі інтервали можуть бути вкрай неточними. Бутстреп надає альтернативу, вільну від розподільних припущень.

Основні переваги бутстрепу для стандартного відхилення:

• Без припущень щодо розподілу: однаково добре працює з нормальними, асиметричними та важкохвостими даними
• Ефективність на малих вибірках: часто точніший за параметричні методи при n < 30
• Робота зі складними статистиками: той самий підхід працює для тримованого СВ, MAD або власних мір варіабельності
• Візуальне розуміння: бутстреп-розподіл показує, що відбувається, а не лише кінцеві числа

Процедура бутстрепу

Алгоритм бутстрепу надзвичайно простий. З вашої початкової вибірки обсягом n спостережень:

Скільки бутстреп-вибірок? B = 1 000 часто достатньо для грубих оцінок і перевірки гіпотез. Для довірчих інтервалів B = 10 000 забезпечує стабільні перцентилі. Для якісних BCa-інтервалів для публікацій рекомендується B = 15 000+.

Методи побудови довірчих інтервалів бутстрепу

Існує кілька методів побудови довірчих інтервалів із бутстреп-вибірок, кожен зі своїми компромісами:

1. Перцентильний метод (найпростіший)

Найінтуїтивніший підхід: просто візьміть перцентилі бутстреп-розподілу.

Для 10 000 бутстреп-вибірок це 250-те та 9 750-те впорядковані значення. Простий, але може бути зміщеним, коли бутстреп-розподіл асиметричний.

2. Базовий (стрижневий) бутстреп

Використовує зв’язок між вибірковою статистикою та бутстреп-статистиками:

Де θ̂ — стандартне відхилення вихідної вибірки. Цей метод “відзеркалює” перцентильний інтервал відносно вибіркової оцінки.

3. BCa (із корекцією зміщення та прискоренням)

Золотий стандарт точності. BCa коригує як зміщення в бутстреп-розподілі, так і прискорення (те, як стандартна похибка змінюється зі значенням параметра). Складніший в обчисленні, але забезпечує інтервали з точністю другого порядку.

Розібраний приклад: ненормальні дані

Розглянемо 15 вимірювань часу реакції (у мс): 245, 312, 287, 456, 234, 298, 267, 523, 289, 301, 278, 645, 256, 289, 312. Ці дані мають правосторонню асиметрію (деякі дуже повільні відповіді).

Бутстреп-ДІ є асиметричним (ширший з правого боку), що відображає правосторонню асиметрію даних. ДІ за хі-квадратом не вловлює цю асиметрію.

Реалізація на Python

Повна реалізація бутстрепу з кількома методами побудови довірчих інтервалів:

import numpy as np
from scipy import stats

def bootstrap_sd_ci(data, n_bootstrap=10000, ci=0.95, method='percentile'):
    """
    Bootstrap confidence interval for standard deviation.

    Parameters:
    -----------
    data : array-like - Original sample
    n_bootstrap : int - Number of bootstrap samples
    ci : float - Confidence level (e.g., 0.95)
    method : str - 'percentile', 'basic', or 'bca'

    Returns:
    --------
    tuple : (lower_bound, upper_bound, bootstrap_sds)
    """
    data = np.array(data)
    n = len(data)
    original_sd = np.std(data, ddof=1)

    # Generate bootstrap samples and calculate SDs
    bootstrap_sds = np.array([
        np.std(np.random.choice(data, size=n, replace=True), ddof=1)
        for _ in range(n_bootstrap)
    ])

    alpha = 1 - ci

    if method == 'percentile':
        lower = np.percentile(bootstrap_sds, 100 * alpha/2)
        upper = np.percentile(bootstrap_sds, 100 * (1 - alpha/2))

    elif method == 'basic':
        lower = 2*original_sd - np.percentile(bootstrap_sds, 100*(1-alpha/2))
        upper = 2*original_sd - np.percentile(bootstrap_sds, 100*alpha/2)

    elif method == 'bca':
        # Bias correction
        prop_less = np.mean(bootstrap_sds < original_sd)
        z0 = stats.norm.ppf(prop_less)

        # Acceleration (jackknife estimate)
        jackknife_sds = np.array([
            np.std(np.delete(data, i), ddof=1) for i in range(n)
        ])
        jack_mean = jackknife_sds.mean()
        a = np.sum((jack_mean - jackknife_sds)**3) / \
            (6 * np.sum((jack_mean - jackknife_sds)**2)**1.5)

        # Adjusted percentiles
        z_alpha = stats.norm.ppf([alpha/2, 1-alpha/2])
        adj_percentiles = stats.norm.cdf(
            z0 + (z0 + z_alpha) / (1 - a*(z0 + z_alpha))
        ) * 100
        lower = np.percentile(bootstrap_sds, adj_percentiles[0])
        upper = np.percentile(bootstrap_sds, adj_percentiles[1])

    return lower, upper, bootstrap_sds

# Example usage
response_times = [245, 312, 287, 456, 234, 298, 267, 523, 289, 301, 278, 645, 256, 289, 312]

for method in ['percentile', 'basic', 'bca']:
    lower, upper, _ = bootstrap_sd_ci(response_times, method=method)
    print(f"{method.upper():12s} 95% CI: [{lower:.1f}, {upper:.1f}]")