310-infer_stats.qmd

# Статистика вывода {#sec-infer_stats}

## Введение в статистику вывода {#sec-infer_intro}

**Статистика вывода *(inferential statistics)*** - это основной раздел статистики, связанный уже не с описанием и суммированием имеющихся данных (описательная статистика), а с попытками сделать вывод о **генеральной совокупности *(population)*** на основе имеющихся данных по **выборке *(sample)****.* Короче говоря, статистика вывода - это о том, как выйти за пределы наших данных. Это именно то, зачем мы проводим исследования - мы не можем собрать информацию обо всей генеральной совокупности, но по тому, что имеем (т.е. по нашей выборке), можем попытаться как-то оценить **параметры** распределения в генеральной совокупности.

Итак, еще раз: генеральной совокупности - параметры **([P]{.underline}opulation - [P]{.underline}arameters),** а у выборки - статистики **([S]{.underline}ample - [S]{.underline}tatistics).** Параметры обычно обозначаются греческими буквами: $\mu$, $\sigma$ (или большими латинскими: $M$, $S$ ), а статистики - соответствующими латинскими: $m$, $s$.

## Случайные величины и их распределения {#sec-distributions}

**Случайная величина (*random variable*)** -- это переменная значения которой в зависимости от случая принимают различные значения. Представьте себе машину, которая по требованию выдает какое-то случайное число (или даже несколько). Это и будет **случайная величина.** Какие это могут быть значения, какие-то конкретные или любые в каком-то диапазоне? Могут ли какие-то значения выпадать чаще других? Очевидно, случайные переменные могут быть разными в зависимости от закона, который стоит за тем, какие значения она может принимать и с какой частотой. Этот закон называется распределением вероятности.

Если случайная переменная может принимать только какие-то конкретные значения (например, 3, 10, 0.25 или число пи), то случайная переменная имеет **дискретное распределение**. Если же случайная переменная может принимать любые значения (в каком-то диапазоне или же вообще любые), то такая случайная переменная будет иметь **непрерывное распределение**. По большей части нас будут интересовать именно непрерывные. Хотя для полного их понимания нужен матан (да, именно тот, который *calculus*), непрерывные распределения довольно интуитивно понятны.

Может ли случайная переменная принимать значения с разной частотой? Если нет, то речь идет о **равномерном распределении** **(*uniform distribution*).** Пример случайной величины с дискретным равномерным распределением -- игральный кубик, с вероятностью $\frac{1}{6}$ возвращает одно из значений: 1, 2, 3, 4, 5 или 6. Однако равномерное распределение может быть как дискретным, так и непрерывным.

Все остальные распределения (кроме **равномерного распределения)** обозначают, что какие-то значения имеют больший шанс выпадения, чем другие. Разные распределения могут быть описаны с помощью **функций распределения**. В R есть свой набор из трех основных функций для каждого распределения и еще одна функция: генератор случайных чисел из выбранного распределения.

Вот эти функции:

-   `d*()` -- **функция вероятности** **(*probability mass function*)** для дискретных распределений и **функция плотности вероятности** для непрерывных распределений;
-   `p*()` -- **функция накопленной плотности распределения** **(*cumulative distribution function; cdf*);**
-   `q*()` -- **квантильная функция** **(*quantile function*),** или **обратная функция накопленной плотности распределения** **(*inverse cumulative distribution function*).**

Ну и функция `r*()` для создания случайных выборок из выбранного распределения.

Во всех четырех случаях, вместо звездочки нужно поставить название соответствующего распределения в R. Таким образом, для каждого распределения в R есть четыре функции (с похожими названиями, различающимися первой буквой), и таких семейств функций в базовом R (точнее, во встроенном пакете `{stats}`) очень много (посмотреть можно с помощью `?Distributions`). Вот лишь некоторые из них:

![](images/rdists.PNG){width="400"}

Возьмем, например, нормальное распределение, и разберем на его примере все эти функции: `dnorm()`, `pnorm()`, `qnorm()` и `rnorm()`.

Для того, чтобы описать нормальное распределение нам нужно всего два параметра - его среднее $\mu$ и стандартное отклонение $\sigma$. Если же $\mu = 0$, а $\sigma = 1$, то такое **нормальное распределение** называется **стандартным *(standard normal distribution).*** Заметьте, здесь мы говорим о *параметрах* распределения в генеральной совокупности, а не о *статистиках* в конкретной выборке, хотя называются они одинаково.

Я думаю, все видели, как выглядит нормальное распределение. Это та самая "колоколообразная кривая"[^310-infer_stats-1].

[^310-infer_stats-1]: Если распределение имеет колоколообразнуюю форму, то это еще не значит, что оно нормальное. Есть очень много похожих по форме распределений, которые, тем не менее, значительно отличаются от нормального, в том числе своими свойствами. Например, распределение Коши тоже внешне похоже на колокол, но у него нет среднего и стандартного отклонения!

![](images/normalstandard.jpg){width="400"}

То, что Вы видите на картинке - это так называемая **функция плотности вероятности *(probability density function).*** Вот ее формула для **нормального распределения *(normal distribution):*** $$P(x) = \frac{e^{-(x - \mu)^{2}/(2\sigma^{2}) }} {\sigma\sqrt{2\pi}}$$

Это аналог ***probability mass function*** для дискретных величин - распределение вероятности того, что случайная величина имеет данное значение. Довольно очевидно, да? Тогда почему аналогичная функция для непрерывных распределений называется по-другому?

Дело в том, что вероятность получения конкретного значения для непрерывных величин равна нулю. Это немного необычный момент, но тем не менее: если распределение непрерывно, то получить точно 0,000000000000... с бесконечным количеством нулей после запятой невозможно, даже если среднее распределения равно нулю. Тем не менее, мы можем посчитать вероятность того, что случайное значение окажется в определенном промежутке. **Функция плотности вероятности** как раз и позволяет это сделать. Чтобы посчитать эту вероятность, надо посчитать площадь соответствующего участка под функцией. Соответственно, площадь под всей функцией плотности вероятности равна 1.

<!--# Здесь картинку с интервалом! -->

Конечно, функцию плотности можно получить в R - это функции `d*()`, где `*` означает соответствующее распределение. Для нормального распределения это будет функция `dnorm()`.

Функция `dnorm()` имеет следующие основные параметры: `x =` - вектор принимаемых значений, среднее и стандартное отклонение распределения (`mean =` и `sd =` соответственно). В качестве параметров по умолчанию используются 0 для среднего и 1 для стандартного отклонения, то есть **стандартное нормальное распределение.** Для других распределений параметры будут отличаться (например, вместо `mean =` и `sd =` в t-распределении будет параметр `df =` -- "степени свободы").

Давайте посмотрим, как работает эта функция, визуализировав результат ее выполнения на векторе от -3 до 3 с небольшим шагом (0.1). С помощью базовой функции `plot()` мы построим диаграмму рассеяния, где по оси *x* используем исходный вектор с последовательностью от -3 до 3 с шагом 0.1, а по оси *y* -- результат применения функции `dnorm()` на исходном векторе. Каждая отдельная точка -- это значение **функции плотности вероятности** (координата *y)* для конкретного *x*. Поскольку таких точек много, и они находятся близко друг к другу, то фактически получается кривая линия.

```{r}
vec <- seq(-3,3, 0.1)
plot(vec,dnorm(vec))
```

Вот мы и получили то, что называем **стандартным нормальным распределением,** точнее, **функцию плотности вероятности** стандартного нормального распределения.

Давайте теперь возьмем другое нормальное распределение с другими параметрами. Например, среднее будет равно 100, а стандартное отклонение -- 15. Это нормы для шкалы IQ:

```{r}
iq <- seq(50,150, 0.1)
plot(iq, dnorm(iq, mean = 100, sd = 15))
```

::: callout-tip
## *Полезное:* шкала IQ

Шкала IQ --- это вообще очень удобная шкала. Поскольку в измерениях интеллекта нет никаких объективных метрик (сантиметров, градусов, килограммов), с которыми его можно было бы сравнить, то единственная возможность дать какую-то оценку величины интеллекта у человека --- сравнить его с интеллектом другого человека, например, по количеству решенных заданий за определенное время. Поэтому шкала IQ так сконструирована, чтобы среднее значение количество решенных заданий обозначалось за 100, а стандартное отклонение за 15. Поэтому по баллу IQ (если используется хороший тест, разумеется) можно понять процент людей, которых респондент опережает по интеллекту, по крайней мере, сравнивая с выборкой, на которой был стандартизирован тест. Примерно так же устроены и другие психологические шкалы.
:::

Следующая функция это **функция накопленной плотности распределения *(cumulative distribution function; cdf)***[^310-infer_stats-2]. Это функция очень важная, потому что именно на ней основано **тестирование уровня значимости нулевой гипотезы,** которым мы будем заниматься в дальнейшем.

[^310-infer_stats-2]: Иногда эта функция называется *интегральной функцией вероятности*, а иногда просто *функцией вероятности.* Последний вариант может запутать, потому что часто под функцией вероятности подразумевают функцию плотности вероятности.

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

Для этой функции распределения используются функции вида `p*()`, в частности, функция `pnorm()` для нормального распределения:

```{r}
plot(iq, pnorm(iq, mean = 100, sd = 15))
```

Какова вероятность того, что полученное случайное значение IQ будет меньше или равно 100?

```{r}
pnorm(100, mean = 100, sd = 15)
```

А меньше или равно 130?

```{r}
pnorm(130, mean = 100, sd = 15)
```

Следующая функция --- это **квантильная функция *(quantile function)****,* или **обратная функция накопленной плотности распределения *(inverse cumulative distribution function):***

```{r}
prob <- seq(0,1, 0.01)
plot(prob, qnorm(prob, mean = 100, sd = 15))
```

Обратная функция означает, что если мы применим сначала одну, а потом другую, то (если не берем особо крайних значений) на выходе получим исходные числа [^310-infer_stats-3]. **Квантильная функция** возвращает значение, которое по заданной вероятности случайная переменная не будет превышать. Поскольку квантильная функция -- это функция от вероятности, квантильная функция определена на отрезке от 0 до 1.

[^310-infer_stats-3]: Если мы возьмем достаточно большие значения, то этот трюк не сработает, что свзязано с тем, что числа дробные числа в компьютерах хранятся с ограниченной точностью.

```{r}
qnorm(pnorm(-4:4))
```

Ну и последняя важная функция - это `rnorm()` - просто генерирует выборку значений из данного распределения заданной длины `n =`:

```{r}
set.seed(42)
samp <- rnorm(100, mean = 100, sd = 15)
samp
```

> `set.seed()` - это функция, которая позволяет получить нам воспроизводимые результаты при использовании генератора случайных чисел. Короче говоря, если все мы поставим `set.seed(42)` [^310-infer_stats-4], то одна и та же строчка выдаст нам один и тот же результат на разных компьютерах. Как это вообще возможно, это же случайные числа? Дело в том, что... нет. Они "псевдо-случайные". На самом деле, используются определенные алгоритмы, чтобы генерировать числа, которые выглядят как случайные. Например, можно брать цифры после запятой в числе пи после, например, 42го знака. Реальные алгоритмы создания псевдо-случайных чисел, конечно, гораздо сложнее, но суть примерно такая. Насколько подобные числа действительно получаются случайными - отдельный сложный математический вопрос. Но для нас этого вполне достаточно.

[^310-infer_stats-4]: А почему именно 42? Ну, можно брать любое число, которое Вам нравится. А 42 - это [ответ на главный вопрос жизни, вселенной и всего такого](https://ru.wikipedia.org/wiki/%D0%9E%D1%82%D0%B2%D0%B5%D1%82_%D0%BD%D0%B0_%D0%B3%D0%BB%D0%B0%D0%B2%D0%BD%D1%8B%D0%B9_%D0%B2%D0%BE%D0%BF%D1%80%D0%BE%D1%81_%D0%B6%D0%B8%D0%B7%D0%BD%D0%B8,_%D0%B2%D1%81%D0%B5%D0%BB%D0%B5%D0%BD%D0%BD%D0%BE%D0%B9_%D0%B8_%D0%B2%D1%81%D0%B5%D0%B3%D0%BE_%D1%82%D0%B0%D0%BA%D0%BE%D0%B3%D0%BE)

```{r}
hist(samp, breaks = 10)
```

## Оценки {#sec-estimates}

Самая основа статистики вывода - это **оценки *(estimates)****.* Как мы уже знаем, у распределений есть определенные параметры, которые описывают данное распределение.

Например, для того, чтобы описать нормальное распределение нам нужно всего два параметра - его среднее $\mu$ и стандартное отклонение $\sigma$.

Наша цель - как-нибудь оценить эти параметры, потому что обычно мы их не знаем. Допустим, человеческий рост распределен нормально \[\^isnotnormal\]. Какой средний рост в популяции? Какое у него стандартное отклонение? Для этого мы используем разного рода оценки: точечные и интервальные. В качестве оценок параметров часто используются статистики по выборке. Вот здесь легко запутаться, поэтому покажу картинку:\
\[\^isnotnormal\]: На самом деле, "чистое" нормальное распределение, так же как и любое другое "чистое" распределение, в природе встречается достаточно редко. Но довольно многие процессы могут апроксимироваться нормальным распределением благодаря центральной предельной теореме @sec-clt.

## Точечные оценки {#sec-sample_estimate}

Представим, что только что сгенерированные данные с помощью `rnorm()` - это и есть наша выборка. Перед нами стоит задача оценить IQ популяции по этой выборке. Попробуем **оценить *(estimate)*** среднее в генеральной совокупности.

Здесь все просто и очень очевидно - самой лучшей оценкой среднего в генеральной совокупности будет среднее по выборке:

$$\hat{\mu}=\frac{\sum_{i=1}^n x_i}{n}$$

Оценки обозначаются значком домика: $\hat{}$. $\sum_{i=1}^n x_i$ означает сумму всех элементов $x$ от $1$ до $n$. Подобная запись с $\sum_{i=1}^n$ нам будет встречаться очень часто!

Посчитаем оценку среднего в генеральной совокупности с помощью уже знакомой нам функции для выборочного среднего:

```{r}
mean(samp)
```

Конечно, не совсем точно, но лучше оценки не придумаешь. Что есть, то есть. Чем больше выборка, тем в ближе оценка будет к популяционному среднему, т.е. тем выше **точность *(efficiency)*** оценки.

А что с оценкой стандартного отклонения? Давайте воспользуемся нашей предыдущей формулой.

$$\hat{\sigma}= \sqrt\frac{\sum\limits_{i=1}^{n} (x_{i} - \overline{x})^2} {n}$$

```{r}
sqrt(sum((samp - mean(samp))^2)/length(samp))
```

В данном случае мы несколько промахнулись вверх, но обычно оценка по этой формуле дает небольшое **смещение** ***(bias)*** в меньшую сторону. В отличие от выборочного среднего как оценки среднего в генеральной совокупности, использование выборочного стандартного отклонения по выборке приводит к смещенной оценке! Само по себе это кажется странным, и это нормально. Но этому есть всякие серьезные математические доказательства. Кроме того, есть множество всяких демонстраций, например, [от Академии Хана](https://www.khanacademy.org/math/ap-statistics/summarizing-quantitative-data-ap/more-standard-deviation/v/simulation-showing-bias-in-sample-variance).

::: callout-warning
## *Для продвинутых:* откуда смещение оценки?

Одно из простых объяснений такое: поскольку мы оцениваем стандартное отклонение на основе оценки среднего, среднеквадратичные расстояния счиnаются не от реального среднего в генеральной совокупности, а от среднего по выборке, которое немного смещено в ту или иную сторону. Представьте, что так получилось, что в нашей выборке оно оказалось сильно смещено и выборочное среднее получилось около 90 (такое бывает). Это значит, что получилось много довольно низких значений около этого 90, а высчитываться среднеквадратичные разницы будут не от среднего 100 (которое мы не знаем), а от этого 90. Поэтому стандартное отклонение получится сильно меньше, чем должно было получиться.
:::

Чтобы получить **несмещенную (unbiased) оценку** стандартного отклонения или дисперсии, то нам нужно делить не на $n$, а на $n-1$.

$$\hat{\sigma} = \sqrt\frac{\sum\limits_{i=1}^{n} (x_{i} - \overline{x})^2} {n-1}$$ Заметьте, что именно так делает стандартная функция `sd()`:

```{r}
sqrt(sum((samp - mean(samp))^2)/(length(samp) - 1))
sd(samp)
```

Это называется поправкой Бесселя, и она настолько распространена, что именно скорректированное стандартное отклонение обычно называют стандартным отклонением. Поэтому в дальнейшем будет использоваться именно эта формула или же просто функция `sd()`.

Таким образом, у оценки есть два критерия качества -- ее **точность** и ее **несмещенность.**

## Интервальные оценки {#sec-interval_estimate}

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

Самая распространенная интервальная оценка называется **доверительным интервалом *(confidence interval)****.* Цель доверительного интервала -- "покрыть" параметр генеральной совокупности с определенной степенью уверенности.

Например, 95% доверительный интервал или просто $CI95\%$ означает, что примерно в 95% случаев подсчитанный на выборках интервал будет ловить значение параметра в популяции. Например, построив $CI95\%$ по нашей сгенерированной выборке, мы хотим чтобы примерно в 95% случаев выборок, построенных таким способом, этот интервал ловил истинное среднее (в данном случае --- 100).

Для того, чтобы научиться строить такие интервалы, нам нужно разобраться с **выборочным распределение *(sampling distribution).***

## Выборочное распределение {#sec-sample_dist}

Чтобы разобраться с тем, что стоит за **доверительными интервалами** и **тестированием уровня значимости нулевой гипотезы** (короче говоря, с основными инструментами статистики вывода), нам надо разобраться с очень абстрактной концепцией -- выборочными распределениями.

Представьте себе, что мы бы выбрали не одну, а сразу много выборок, а потом у каждой выборки посчитали бы среднее. Мы бы получили новый вектор данных - средние выборок из одного распределения. Давайте это сделаем[^310-infer_stats-5]:

[^310-infer_stats-5]: Конечно, на практике никто не будет собирать тысячу выборок. Но нам важно понимать, чтобы было, если бы мы повторяли один и тот же эксперимент в одинаковых условиях, даже если это и невозможно сделать. Это позволит нам понять, какие выводы мы можем сделать на основе тех данных, что мы имеем на самом деле.

```{r}
sample_means <- replicate(1000, mean(rnorm(100, mean = 100, sd = 15)))
```

Каждая выборка состоит из сотни "испытуемых", всего таких выборок 1000. По каждой мы посчитали среднее. Как распределены эти средние?

```{r}
hist(sample_means, breaks = 30)
```

Вот это распределение и есть **выборочное распределение средних *(sampling distribution of means).*** Точнее, было бы, если бы мы взяли не 1000 выборок, а бесконечное количество выборок. Так что у нас всего лишь аппроксимация выборочного распределения средних, к сожалению.

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

Например, среднее выборочного распределения средних будет равно популяционному среднему:

$$
\mu_{\overline{x}} = \mu
$$

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

```{r}
mean(sample_means)
```

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

Теперь перейдем к стандартному отклонению выборочного распределения средних. Звучит очень громоздко, не правда ли? Поэтому стандартное отклонение выборочного распределения обычно называют **стандартной ошибкой *(standard error),*** а стандартное отклонение выборочного распределения средних называют **стандартной ошибкой среднего** ***(standard error of the mean; s.e.m.).***

Чему будет равна **стандартная ошибка среднего?** Давайте посмотрим на стандартное отклонение в нашей симуляции выборочного распределения средних:

```{r}
sd(sample_means)
```

Очень похоже на стандартное отклонение в шкале IQ (15), но в 10 раз меньше. Это не случайно: **стандартная ошибка среднего** равна стандартному отклонению в генеральной совокупности, деленному на корень из размера выборки.

$$\sigma_{\overline{x}}= \frac{\sigma} {\sqrt{n}}$$

Как раз это мы и получили: размер нашей выборки - 100, а корень из 100 равен 10:

```{r}
15/sqrt(length(samp))
```

![](images/310-infer_stats_sampling_dist.png)

Из этой формулы следует, что стандартная ошибка тем меньше, чем больше выборка. И это довольно логично: чем больше у нас размер выборок, тем меньше будет разброс их средних. При этом уменьшение разброса не будет линейным: чтобы уменьшить его в 2 раза, нужно увеличить выборку в 4 раза, а чтобы уменьшить его в 10 раз, нужно увеличить выборку в 100 раз.

Ну а если мы не знаем стандартного отклонения в генеральной совокупности (что обычно и бывает в жизни), то можем оценить стандартную ошибку среднего с помощью оценки стандартного отклонения, посчитанного на выборке:

$$\hat{\sigma_{\overline{x}}}= \frac{\hat{\sigma}}{\sqrt{n}}$$

```{r}
sd(samp)/sqrt(length(samp))
```

Именно **оценку** **стандартной ошибки среднего** обычно используют на графиках в научных статьях в качестве **усов** ***(error bars)***[^310-infer_stats-6]***.***

[^310-infer_stats-6]: Для барплотов или точечных диаграмм. Для ящиков с усами используются другие правила, см. @sec-base_vis_boxplot.

## Важность нормального распределения {#sec-nordm_dist_importance}

В математике все числа равны, но некоторые все-таки заметно равнее других. Например, есть 0 и 1 - и очевидно, что это очень важные числа. Даже гугл выдает больше страниц по числам 0 и 1, чем по другим числам, например, 8 и 23218974. Не только целые числа могут быть важными. Вот, например, число пи и число Эйлера: они очень часто встречаются в самых разных, подчас неожиданных местах. Например, в тех же функциях распределений. На эти числа завязано очень интересные и важные свойства, поэтому такая зацикленность на них неудивительна. Например, число пи связывает радиус круг с длиной его окружности и площадью, поэтому когда перед нами что-то круглое, то и пи появится с большой вероятностью.

В статистике тоже есть распределения, которые "главнее" других распределений. Есть, например, равномерное распределение. Например, равномерно распределены исходы одного броска кубика.

Но есть и распределение, которое смело можно называть королем всех распределений - это уже знакомое нам нормальное распределение.

Нормальное распределение - это что-то в духе единицы из мира распределений. Это не просто одно из распределений, а это фундаментально важная штуковина, которая обладает почти магической силой. И сила эта зовется "Центральная Пределельная Теорема".

## Центральная предельная теорема {#sec-clt}

Давайте теперь сделаем много выборок не из нормального, а из логнормального распределения. Что-то похожее представляет собой распределение времени реакции на многие задачи (особенно связанные с выбором и когда нужно подумать перед нажатием кнопки). Одна выборка будет распределена примерно так:

```{r}
hist(rlnorm(10000), breaks = 100)
```

Распределение сильно ассиметрично, но его форма примерно понятна.

А как будут распределены средние многих выборок, взятых из логнормального распределения? Хочется сказать, что так же, но нет!

```{r}
many_means <- replicate(1000, mean(rlnorm(10000)))
hist(many_means, breaks = 100)
```

Удивительно, но средние по выборкам из логнормального распределения будут выглядеть почти нормально! Более того, для других распределений это тоже будет верно: **согласно центральной предельной теореме *(ЦПТ, central limit theorem)*, какой бы ни была форма распределения в генеральной совокупности, выборочное распределение средних будет стремиться к нормальному. При этом чем больше размер выборки, тем ближе выборочное распределение средних будет к нормальному.** Это очень важная фишка, на которой основаны многие статистические тесты.

![](images/IllustrationCentralTheorem.png){width="400"}

Не верите? Попробуйте сами! Можете поиграться с разными распределениями с помощью кода (можете посмотреть другие распределения в хэлпе: `?Distributions` ). А еще можно потыкать интерактивную [Shiny-демонстрацию магии ЦПТ](https://gallery.shinyapps.io/CLT_mean/). Я очень рекомендую поиграться с ней. Попробуйте разные значения, посмотрите что будет.

Как только наиграетесь, то сразу станет понятно, почему именно нормальное распределение занимает такое важное место в статистике.

::: callout-note
## *Полезное:* ЦПТ в реальной жизни

Если выходить за рамки выборочного распределения средних, то центральная предельная теорема говорит нам о том, что сумма слабо зависящих случайных величин, имеющих примерно одинаковое влияние, имеет распределение близкое к нормальному. Например, рост является следствием большого количества генетических и средовых факторов, поэтому он распределен примерно нормально. Однако это нормальное распределение портят такие факторы как пол --- он вносит очень сильный вклад в рост, поэтому распределение становится более плоским, т.к. это смесь двух похожих на нормальное распределений. И даже внутри этих распределений есть более и менее влияющие факторы. Поэтому какую бы более узкую генеральную совокупность мы ни рассматривали, распределение никогда не будет идеально нормальным. За исключением, возможно, фундаментальных физических процессов.
:::

## Строим доверительный интервал {#sec-ci_build}

Теперь мы знаем достаточно, чтобы построить доверительный интервал своими руками на основе стандартной ошибки. Мы построим самый стандартный вариант - 95% доверительный интервал.

Давайте еще раз посмотрим на нормальное распределение.

![](images/normalstandard.jpg){width="400"}

Мы хотим поймать симметрично 95% от площади под кривой. Для этого нам нужно отбросить по 2.5% с обоих сторон. Эти 2.5% соответствуют примерно двум стандартным отклонениям от среднего. Если быть точнее, то 1.96. Если быть еще точнее:

```{r}
qnorm(0.975)
```

Почему 0.975? Потому что мы смотрим квантильную функцию по верхней границе: отсекаем правые 0.025:

```{r}
qnorm(1 - (1 - 0.95)/2)
```

Давайте сохраним это число. Назовем его `zcr`:

```{r}
zcr <- qnorm(1 - (1 - 0.95)/2)
```

Это количество стандартных отклонений от среднего в нормальном распределении, которое включает в себя ровно 95% площади нормального распределения. Теперь давайте посчитаем стандартную ошибку. Здесь мы знаем стандартное отклонение в генеральной совокупности (это 15), его поделим на корень из размера выборки:

```{r}
sem <- 15/sqrt(length(samp))
```

Чтобы посчитать нижнюю и верхнюю границы доверительного интервала, нам нужно вычесть и прибавить соответственно нужное количество стандартных ошибок:

```{r}
mean(samp) - sem*zcr #нижняя граница
mean(samp) + sem*zcr #верхняя граница
```

Давайте теперь нарисуем сотню доверительных интервалов с помощью `ggplot2`! Цветом обозначим интервалы, которые не поймали истинное значение параметра в центральной совокупности.

```{r}
library(tidyverse)
sample_size <- 100
set.seed(42)
ci_simulations <- tibble(
  m = replicate(sample_size, mean(rnorm(sample_size, mean = 100, sd = 15))),
  se = 15/sqrt(sample_size),
  lower = m - se*zcr,
  higher = m + se*zcr,
  parameter_inside = lower<100 & higher>100
)
```

```{r}
many_ci_gg <- ggplot(data = ci_simulations, aes(x = 1:sample_size,y = m)) +
  geom_pointrange(aes(ymin = lower,ymax = higher,colour = parameter_inside))+
  geom_hline(yintercept = 100)+
  coord_flip() +
  theme_minimal()
many_ci_gg
```

Примерно 5% не ловят 100 в 95% доверительный интервал! Примерно это и означает доверительный интервал: где-то в 95% он ловит параметр в генеральной совокупности, а в 5% - нет.

::: callout-important
## *Осторожно:* определение доверительного интервала

Понятие доверительного интервала вызывает кучу недопонимания и ошибок. Очень многие его интерпретируют, например, как интервал, включающий в себя 95% значений популяции, но это неправильно.
:::

Еще я советую посмотреть вот эту визуализацию: [Доверительные интервалы сыпятся как из мешка](http://rpsychologist.com/d3/CI/)

## Тестирование значимости нулевой гипотезы {#sec-nhst}

**Тестирование значимости нулевой гипотезы (null hypothesis significance testing)** --- это основной подход в статистике вывода. Вы про него точно слышали (хотя, возможно, не знали, что он так называется), потому что де-факто он является стандартом в психологии, биологии, медицине и многих других науках.

Мы сейчас детально его проведем на примере одного из самых простых статистических тестов - **z-тестов *(z-test).*** Однако та же самая логика стоит и за остальными статистическими тестами.

![](images/310-infer_stats_NHST.png)

**1. Формулирование нулевой и альтернативной гипотезы.**

Сначала мы задаем две гипотезы о параметрах распределения. Одна из них называется нулевой: она обычно включает положение о том, что различий или связи нет или что это различие/связь равно определенному числу. Если мы хотим применить тестирование значимости к нашей "выборке", то нулевую гипотезу можно будет сформулировать так: $$H_0: \mu = 100$$

Наша нулевая гипотеза заключается в том, что выборка была взята из распределения со средним 100.

Альтернативная или ненулевая гипотеза либо говорит о том, что среднее в генеральной совокупности на самом деле не равно какому-то конкретному числу (в нашем случае --- 100) или что две выборки взяты из групп с различным средним и т.п.

$$H_1: \mu \ne 100$$

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

Заметьте, мы формулируем гипотезу не про *статистики* в выборке, а про *параметры* в генеральной совокупности, поэтому пользуемся греческими (или большими латинскими) буквами.

Вся дальнейшая логика расчетов будет строиться именно на нулевой гипотезе: мы будем пытаться понять, насколько реалистичны наши результаты при верности нулевой гипотезы, которую мы заранее обозначили. Это похоже на "доказательство от обратного" в геометрии: мы исходим из того, что эффекта не существует и пытаемся прийти к противоречию с данными, чтобы отбросить эту нулевую гипотезу и принять альтернативную гипотезу.

**2. Подсчет тестовой статистики по выборке**

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

Если мы знаем стандартное отклонение в генеральной совокупности, то можем посчитать $z$-статистику по формуле:

$$z = \frac{\overline{x} - \mu} {\sigma / \sqrt{N}} $$

```{r}
m <- mean(samp)
sem <- 15/sqrt(length(samp))
z <- (m - 100)/sem
z
```

$z$-статистика --- это выборочное среднее, из которого вычтено среднее в генеральной совокупности согласно нашей гипотезе. Получившуюся разницу мы делим на стандартную ошибку.

**3. Расчет p-value**

И вот мы подобрались к самому важному этапу -- расчет **p-value**.

> **p-value** -- это вероятность получить такую и более отклоняющуюся тестовую статистику при условии верности нулевой гипотезы.

Очень важно понять, как именно рассчитывается *p-value,* потому что на этом основывается сама идея тестирования значимости нулевой гипотезы! Для этого нам нужно вернуться к идее выборочного распределения, только теперь уже не среднего, а $z$-статистики. Впрочем, выглядеть оно будет абсолютно так же -- нормально! Благодаря тому, что мы вычли среднее и поделили на стандартное отклонение, среднее этого нормального распределения будет равно 0, а стандартное отклонение -- 1. То есть перед нами снова **стандартное нормальное распределение.**

Благодаря ЦПТ, даже если распределение в генеральной совокупности несколько отличается от нормального, а выборка достаточно большая, то выборочное распределение $z$-статистик при верности нулевой гипотезы будет (примерно) нормальным.

Теперь нам нужно соотнести нашу $z$-статистику с теоретическим выборочным распределением $z$-статистик. Это позволит нам оценить вероятность получить такие и более отличающиеся результаты при допущении верности нулевой гипотезы, т.е. p-value.

Какая вероятность получить $z$-статистику `r z`, если на самом деле нулевая гипотеза верна? Это вопрос с подвохом: как мы выяснили ранее, для непрерывных распределений вероятность получить отдельное число равна 0. Но мы можем посчитать, какая вероятность получить такую же или большую $z$-статистику!

Графически эту вероятность можно представить как площадь под кривой функции плотности распределения от полученной $z$-статистики до плюс бесконечности или от минус бесконечности до полученной $z$-статистики, если $z$-статистика отрицательная.

![](images/explanationsPValue.png){width="400"}

Или же можем воспользоваться **функцией накопленной плотности распределения.** Для нормального распределения это можно сделать с помощью уже знакомой нам функции `pnorm()`:

```{r}
pnorm(z)
```

`pnorm()` считает от минус бесконечности до заданного числа, а нам нужно наоборот --- от заданного числа до плюс бесконечности, потому что $z$ отличается от 0 в большую сторону. Этого можно добиться вычетанием из 1[^310-infer_stats-7]:

[^310-infer_stats-7]: или же можно поставить параметру `lower.tail =` функции `pnorm()` значение `FALSE`.

```{r}
1 - pnorm(z)
```

Обычно это число еще и умножают на 2, потому что мы заранее не знаем, в какую сторону будет отклоняться среднее по нашей выборке. Вернемся к шагу 1: мы сформулировали ненулевую гипотезу таким образом, что среднее генеральной совокупности не равно 100:

$$H_1: \mu \ne 100$$

Это означает, что $H_1$ включает в себя как случаи, когда среднее отклоняется в большую сторону, так и случаи когда среднее отклоняется в меньшую сторону.

```{r}
p <- (1 - pnorm(z))*2
p
```

Вот эта число и есть *p-value ---* *вероятность получения такого же и более экстремального значения тестовой статистики при условии, что нулевая гипотеза верна.*

**4. Принятие решения о гипотезах**

Отлично, мы посчитали *p-value*. В данном конкретном случае он оказался равен `r p`. Это много или мало? Фактически это означает, что если нулевая гипотеза верна, то в большинстве случаев мы будем получать z-статистики больше нашей. Короче говоря, это вполне реалистичный случай, если нулевая гипотеза верна.

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

Поэтому если *p-value* большой, мы не можем сделать выводов про верность нулевой и альтернативной гипотезы. Мы можем лишь сказать, что у нас нет оснований отклонить нулевую гипотезу.

Если же *p-value* очень маленький, то здесь у нас появляется больше однозначности. Например, если *p-value* равен .02, мы можем сказать, что ситуация малореалистичная: такие и более сильные отклонения от среднего мы можем получить только раз в 50 случаев, если $H_0$ верна. Поэтому мы отклоняем $H_0$ и принимаем $H_1$.

Насколько маленьким должен быть *p-value*, чтобы отклонить нулевую гипотезу? **Критическое значение *p-value*, при котором отклоняют нулевую гипотезу, называется уровнем** $\alpha$. Это максимальный уровень ошибки, который мы допускаем в исследовании.

Так получилось исторически, что стандартный уровень $\alpha$ равен .05. Нужно помнить, что .05 --- это просто общепринятая условность, за этим числом не стоит никакого сакрального знания. Просто так получилось.

Очевидно, что такой статистический подход к принятию решений будет периодически приводить к нас ошибкам. Если $H_0$ на самом деле верна, а мы ее отвергли и приняли $H_1$, то это **ошибка первого рода** (*type I error*). Вероятность этой ошибки и есть наше критическое значение $\alpha$. Однако есть вероятность ошибиться и в другую сторону, т.е. ошибочно не отклонить $H_0$ --- это **ошибка второго рода** (*type II error*), эта вероятность обозначается буквой $\beta$.

| Принятое решение  Реальность | $H_0$ верна                  | $H_1$ верна                   |
|------------------------|------------------------|-------------------------|
| Не отклоняем $H_0$           | Верный пропуск               | Ошибка 2 рода (type II error) |
| Отклоняем $H_0$              | Ошибка 1 рода (type I error) | Верное попадание              |