<jupyter_start><jupyter_text>Задание 4 Для нормального распределения $X_i \sim N(\theta, \sigma ^2)$, если взять сумму независимых $X_i$, то доверительный интервал для $\theta$ будет:$$P_{\theta}\left(\overline{x} - \frac{c_{1-\frac{\alpha}{2}} \sigma}{\sqrt{n}} < \theta < \overline{x} + \frac{c_{1-\frac{\alpha}{2}} \sigma}{\sqrt{n}} \right) = 1 - \alpha$$Для **параметрического бутстрэпа** - оценим $\hat{\theta}$ и $\hat{\sigma ^2}$ по выборке и построим дополнительные выборки из распределения $X_i \sim N(\hat{\theta}, \hat{\sigma ^2})$Для **не параметрического бутстрэпа** - из исходной выборки будем делать случайный выбор с возвратом.  сделаем пункты 1-4 на 1 выборке<jupyter_code>import numpy as np
from scipy import stats # чтобы считать квантили
def ci_theoretical(data, alpha=0.05):
    """дает кортеж из нижней и верхней границы доверительного интервала для среднего в нормальных x_i

    Args:
        data (array like): данные sample
        alpha (float, optional): уровень значимости. Defaults to 0.05.
    """
    std = np.std(data, ddof=1, dtype=np.float64)
    n = len(data)
    avg = np.mean(data)


    CI_Theoretical = (avg + stats.norm.ppf(alpha/2)*std /np.sqrt(n), avg + stats.norm.ppf(1-alpha/2)*std /np.sqrt(n))
    return CI_Theoretical


n = 50 # размер выборки 
alpha = 0.05 # параметр ошибки

np.random.seed(80085)
theta = np.random.uniform(10,20) # параметр тета
theor_sigma = np.random.uniform(1,2) # параметр тета
print((theta, theor_sigma))
sample1 = np.random.normal(theta,theor_sigma,size=n) 
sample1<jupyter_output><empty_output><jupyter_text>3.1) теоретический доверительный интервал<jupyter_code>ci_result_theor = ci_theoretical(sample1)
ci_result_theor<jupyter_output><empty_output><jupyter_text>3.2) параметрический бутстрэп доверительный интервал<jupyter_code>def ci_param_bootstrap(data, alpha=0.05, number_of_bootstrap_samples=10, size_of_bootstrap_samples=20 ):
    """параметрический бутстрэп

    Args:
        data (array like): данные для оценки среднего
        alpha (float, optional): увроень доверия. Defaults to 0.05.
        number_of_bootstrap_samples (int, optional): сколько сэмплов для бутстрэпа делать. Defaults to 10.
        size_of_bootstrap_samples (int, optional): сколько наблюдений делать в каждый сэмпл. Defaults to 20.
    """
    # Оцениваем неизвестный параметр theta 
    sample_mean = np.mean(data) 
    sample_std = np.std(data, ddof=1)
    # print(sample_mean, sample_std)

    # Генерируем выборку из распределения N(sample_mean, sigma)
    bootstrap_samples = np.random.normal(sample_mean,sample_std,size=[number_of_bootstrap_samples,size_of_bootstrap_samples]) 

    #  Считаем среднее для каждой выборки 
    bootstrap_estimates = np.apply_along_axis(np.mean, 1, bootstrap_samples)

    # Вычисляем параметрический бутстрэп доверительный интервал
    CI_Bootstrap_Parametric = (np.quantile(bootstrap_estimates,alpha/2), np.quantile(bootstrap_estimates,1-alpha/2))

    return(CI_Bootstrap_Parametric)
    
ci_result_parambs = ci_param_bootstrap(sample1)
ci_result_parambs<jupyter_output><empty_output><jupyter_text>3.3)  непараметрический бутстрэп интервал<jupyter_code>def ci_non_param_bootstrap(data, alpha=0.05, number_of_bootstrap_samples=10, size_of_bootstrap_samples=20 ):
    """непараметрический бутстрэп

    Args:
        data (array like): данные для оценки среднего
        alpha (float, optional): увроень доверия. Defaults to 0.05.
        number_of_bootstrap_samples (int, optional): сколько сэмплов для бутстрэпа делать. Defaults to 10.
        size_of_bootstrap_samples (int, optional): сколько наблюдений делать в каждый сэмпл. Defaults to 20.
    """

    # Генерируем выборку из исходного распределения
    bootstrap_samples = np.random.choice(data,size=[number_of_bootstrap_samples,size_of_bootstrap_samples]) 

    #  Считаем среднее для каждой выборки 
    bootstrap_estimates = np.apply_along_axis(np.mean, 1, bootstrap_samples)

    # Вычисляем параметрический бутстрэп доверительный интервал
    CI_Bootstrap_Parametric = (np.quantile(bootstrap_estimates,alpha/2), np.quantile(bootstrap_estimates,1-alpha/2))

    return(CI_Bootstrap_Parametric)

ci_result_non_parambs = ci_non_param_bootstrap(sample1)
ci_result_non_parambs<jupyter_output><empty_output><jupyter_text>сравним их длины<jupyter_code>print(f'Длина интервала для теоретического доверительного интервала {ci_result_theor[1] - ci_result_theor[0]:.4f}')
print(f'Длина интервала для параметрического бутстрэпа {ci_result_parambs[1] - ci_result_parambs[0]:.4f}')
print(f'Длина интервала для непараметрического бутстрэпа {ci_result_non_parambs[1] - ci_result_non_parambs[0]:.4f}')<jupyter_output>Длина интервала для теоретического доверительного интервала 0.5064
Длина интервала для параметрического бутстрэпа 0.6527
Длина интервала для непараметрического бутстрэпа 0.6477<jupyter_text>Для этого случайного зерна и параметров, мы получили то, что длина теоретического доверительного интервала самай коротка, потом идет непараметрический бутстрэп, потом параметрический бутстрэп (самый широкий интервал)  Сделаем это 10000 раз<jupyter_code>import pandas as pd 
from tqdm.notebook import trange
iterations = 10000

theoretical = np.zeros([iterations,2]) # здесь будем хранить результаты для теоретического доверительного интервала
parametric_bootstrap = np.zeros([iterations,2]) # здесь будем хранить результаты для параметрического бутстрэпа 
nonparametric_bootstrap = np.zeros([iterations,2]) # здесь будем хранить результаты для непараметрического бутстрэпа 

for i in trange(iterations, desc='calc CI'):
    sample4 = np.random.normal(theta,theor_sigma,size=n) 

    CI_Theoretical = ci_theoretical(sample4,alpha)
    theoretical[i,0] = (theta >= CI_Theoretical[0]) and (theta <= CI_Theoretical[1])
    theoretical[i,1] = CI_Theoretical[1] - CI_Theoretical[0]

    CI_param = ci_param_bootstrap(sample4,alpha)
    parametric_bootstrap[i,0] = (theta >= CI_param[0]) and (theta <= CI_param[1])
    parametric_bootstrap[i,1] = CI_param[1] - CI_param[0]

    CI_nonparam = ci_non_param_bootstrap(sample4,alpha)
    nonparametric_bootstrap[i,0] = (theta >= CI_nonparam[0]) and (theta <= CI_nonparam[1])
    nonparametric_bootstrap[i,1] = CI_nonparam[1] - CI_nonparam[0]




prop_correct = {}
average_length = {}

prop_correct['theoretical'], average_length['theoretical'] = np.mean(theoretical, axis=0)
prop_correct['param'], average_length['param'] = np.mean(parametric_bootstrap, axis=0)
prop_correct['non_param'], average_length['non_param'] = np.mean(nonparametric_bootstrap, axis=0)


df = pd.DataFrame.from_records([prop_correct, average_length], index = ['доля попадений в CI', 'длина CI средняя']).T
df<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1044620205939412, 'fraction_numerical': 0.022384718698701685, 'mean_word_length': 5.66865671641791, 'pattern_counts': {'":': 0, '<': 29, '<?xml version=': 0, '>': 27, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '649205', 'license': 'MIT', 'n_tokens_mistral': 2651, 'n_tokens_neox': 2699, 'n_words': 492, 'path': '03 Math/Stat hw2/hw2.ipynb', 'repo_name': 'ksetdekov/HSE_DS', 'size': '6701'}