Небольшое предисловие¶
В погоне за сложными моделями машинного обучения, которые могут делать различные крутые штуки, многие аналитики незаслуженно забывают про статистику — фундамент, на котором стоит весь Data Science.
Статистика делает тебя сильнее как аналитика, так как она даёт уверенность в результатах, учит критическому мышлению и помогает объяснять выводы бизнесу или научному сообществу. Даже самая продвинутая модель будет бесполезна, если мы неверно интерпретируем её предсказания.
Так что давай вместе укреплять базу — ведь именно с неё начинается настоящий профессионализм в Data Science.
# Bot check
# HW_ID: phds_sem6
# Бот проверит этот ID и предупредит, если случайно сдать что-то не то.
# Status: not final
# Перед отправкой в финальном решении удали "not" в строчке выше.
# Так бот проверит, что ты отправляешь финальную версию, а не промежуточную.
# Никакие значения в этой ячейке не влияют на факт сдачи работы.
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as sps
import seaborn as sns
from matplotlib.lines import Line2D
sns.set(font_scale=1.5)
1. Несмещенность оценок¶
1. Пусть $X_1, ..., X_n$ — выборка из распределения $U[0, \theta]$. Рассмотрим оценки $X_{(n)}, \frac{n+1}{n}X_{(n)}, 2\overline{X}$ параметра $\theta$.
Давайте попробуем определить какие из этих оценок являются несмещенными, не используя теоретических выкладок.
Для начала вспомним,в чем смысл несмещенности
Ответ:
Теперь, чтобы решить задачу, визуализируем свойство немещенности.
Реализуем выборку из равномерного распределения размером sample_size=100.
sample_size = 100
X = sps.uniform.rvs(size=sample_size)
В этом случае мы сгенерировали выборку из распределения $U[0, 1]$.
Реализуйте три функции, каждая из которых на вход берет несколько выборок, а на выход выдает массив оценок
(первая функция - $X_{(n)}$, вторая функция - $\frac{n+1}{n}X_{(n)}$, третья функция - $2\overline{X}$) для каждой выборки.
def estimate_X_n(X):
"""
Принимает на вход массив размером (n_samples, sample_size), выдает массив оценок размера (n_samples,)
"""
return <...>
def estimate_X_n_corrected(X):
"""
Принимает на вход массив размером (n_samples, sample_size), выдает массив оценок размера (n_samples,)
"""
return <...>
def estimate_2_mean(X):
"""
Принимает на вход массив размером (n_samples, sample_size), выдает массив оценок размера (n_samples,)
"""
return <...>
Проверим, что ваши функции реализованы корректно (ячейка не должна выдавать ошибок).
X = np.array([[1, 3, 3902, 6], [2, 5, 69751, 89]])
assert(sum(estimate_X_n(X) != np.array(( 3902, 69751))) == 0)
assert(sum(estimate_X_n_corrected(X) != np.array((4877.5 ,87188.75))) == 0)
assert(sum(estimate_2_mean(X) != np.array([ 1956. , 34923.5])) == 0)
Зададим список оценок и разные параметры для отрисовки графика
estimators = [
(estimate_X_n, r'$X_{(n)}$', 'blue', 0.0),
(estimate_X_n_corrected, r'$\frac{n+1}{n}X_{(n)}$', 'purple', 0.1),
(estimate_2_mean, r'$2\overline{X}$', 'green', 0.2)
]
Мы хотим понять, являются ли оценки параметра $\theta$ смещенными, для этого нам нужно провести множество экспериментов (сгенерировать выборку много раз, так как при каждой генерации получаются разные числа).
Напишем функцию построения таких графиков в общем виде: на вход функция может получать любые распределения (distributions) и любые оценки (estimators)
def est_plot(distribution, estimators, sample_size, sample_count):
'''
Построение графика разброса реализаций оценок и их средних значений.
distribution -- распределение формата scipy.stats
estimators -- список оценок и параметров для отрисовки графиков
sample_size -- размер выборок
sample_count -- количество генерируемых выборок
'''
X = distribution.rvs(size=(sample_count, sample_size))
plt.figure(figsize=(15, 0.7*len(estimators)))
for estimator, name, color, y in estimators:
E = estimator(X)
plt.scatter(E, np.zeros(sample_count) + y, alpha=0.1,
s=100, color=color, label=name)
plt.scatter(E.mean(), y, marker='*', s=300,
color='w', edgecolors='black')
plt.vlines(1, -1, 1, color='r')
plt.title('Размер выборки = %d' % sample_size)
plt.yticks([])
plt.legend()
plt.xlim((0.8, 1.2))
plt.ylim((-0.1, 0.1 * len(estimators)))
Постройте три графика, с помощью функции est_plot. В каждом эксперименте будем генерировать 500 выборок размера (10, 100, 500).
sample_size_list = (10, 100, 500) # размеры выборок
sample_count = 500 # количество экспериментов
for sample_size in sample_size_list:
est_plot(sps.uniform, estimators, sample_size, sample_count)
Что можете сказать о поведении оценкок при различных размерах выборок?
Ответ:
2. Изучим поведение среднего оценок из первого пункта при росте размера $n$ выборки. Для вычисления зависимости нужно один раз сгенерировать выборки из п. 1.1 достаточно большого размера и посчитать оценки по префиксам, используя функции из numpy. Какие из оценок являются асимптотически несмещёнными (т.е. $\forall \theta \in \Theta\colon \mathsf{E}_\theta \widehat{\theta} \to \theta$ при $n\to +\infty$)?
def mean_plot(distribution, estimators, n_grid, sample_count):
'''
distribution -- распределение формата scipy.stats
estimators -- список оценок и параметров для отрисовки графиков
n_grid -- массив значений размера выборки
sample_count -- количество генерируемых выборок
'''
n_grid = n_grid.astype(int)
X = distribution.rvs(size = (sample_count, len(n_grid)))
plt.figure(figsize=(15, 3*len(estimators)))
plt.hlines(1, 0, len(n_grid), color='r', label=r'$\theta$')
for estimator, name, color, y in estimators:
E = np.array([np.mean(estimator(X[0:len(n_grid), 0:(n_grid[i]-1)])) for i in range(len(n_grid))])
plt.plot(n_grid, E, color= color, label= name)
plt.xlabel('размер выборки')
plt.legend()
plt.title('зависимость среднего оценок от размера выборки')
plt.xlim((10, len(n_grid)))
n_grid = np.linspace(10, 300, 100)
mean_plot(sps.uniform, estimators, n_grid, sample_count)
Сделайте выводы
Вывод:
Почему в лабораторных работах часто используют именно скорректированные оценки?
Ответ:
2. Состоятельность оценок¶
Пусть $X_1, ..., X_n$ — выборка из распределения $Exp(\theta)$. Как известно из теории, оценка $\widehat{\theta} = 1/\overline{X}$ является состоятельной и асимтотически нормальной оценкой параметра $\theta$ с асимптотической дисперсией $\theta^2$. В этой задаче вам необходимо визуализировать данные свойства.
Зададим параметры эксперимента
theta = 2 # истинное значение параметра
sample_size = 300 # размер выборок
sample_count = 500 # количество выборок
n_range = (np.arange(sample_size) + 1) # размеры подвыборок
Сгенерируем множество выборок (количество выборок - sample_count, размер каждой выборки - sample_size)
# генерируем множество выборок,
# параметр theta является обратным к параметру масштаба
samples = <...>
Вычислим оценки в зависимости от размера выборки, усредняем по последней оси
estimation = <...>
Визуализируем для каждой выборки все оценки в зависимости от размера выборки.
Что можно сказать о состоятельности данной оценки?
def est_plot(estimation, sample_count, sample_size, left=None, right=None, xlim = sample_size ):
'''
estimation -- массив оценок от размера выборки
sample_count -- количество генерируемых выборок
sample_size -- размер каждой выборки
left, rigth -- границы доверительного интеравала, нужны будут далее
xlim = область по x
'''
plt.figure(figsize=(15, 7))
for i in range(sample_count):
plt.plot(np.arange(sample_size) + 1, estimation[i], color='blue', alpha=0.05)
if type(left) and type(right) is np.ndarray:
plt.plot(np.arange(sample_size) + 1, left, color='black')
plt.plot(np.arange(sample_size) + 1, right, color='black')
labels = [r'$\hat{\theta}$', r'$\theta$',
'Границы доверительного интервала']
handels = [Line2D([0], [0], color='blue', lw=2),
Line2D([0], [0], color='red', lw=2),
Line2D([0], [0], color='black', lw=2)]
else:
labels = [r'$\hat{\theta}$', r'$\theta$']
handels = [Line2D([0], [0], color='blue', lw=2),
Line2D([0], [0], color='red', lw=2),]
plt.hlines(theta, 0 ,sample_size, color='red')
plt.title('Поведение оценки для разных реализаций')
plt.xlabel('Размер выборки')
plt.ylim((0, 5))
plt.legend(handels, labels)
plt.xlim((0, xlim));
est_plot(estimation, sample_count, sample_size)
Ответ: