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

Statistics of location / Центральные тенденции

Statistics of dispersion / Меры разброса

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

set.seed(14)
apples <- rnorm(8, 8, 1) #Возьмем восемь случайных
                          #значений из нормального распределения
apples <- round(apples,1) #Перезапишем вектор, округлив измерения
                          #до первого знака после запятой
apples

## [1]  7.3  9.7 10.1  9.5  8.0  9.2  7.9  9.1

Наши данные в исходном виде выглядят примерно так:

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

sort(apples)

## [1]  7.3  7.9  8.0  9.1  9.2  9.5  9.7 10.1

В ранжированном ряду медиана расположена так, что слева и справа от нее находится равное число измерений.

• Если n нечетное, то медиана = значение с индексом \(\frac{n+1}{2}\).
• Если n четное, то медиана посередине между \(\frac{n}{2}\) и \(\frac{n+1}{2}\) значениями.

sort(apples)

## [1]  7.3  7.9  8.0  9.1  9.2  9.5  9.7 10.1

Медиана находится в промежутке между значениями 9.1 и 9.2, т.е. 9.15

Проверим себя:

median(apples)

## [1] 9.15

Медиана обладает чрезвычайно приятным свойством – устойчивостью к выбросам.

Представим, что наше распределение пострадало от неаккуратности. Допустим сотрудник, которому мы поручили измерять яблоки, измерил также диаметр арбуза и записал этот результат вместе со всеми остальными.

apples2 <- c(apples, 68) #Создадим вектор с новым значением
sort(apples2)

## [1]  7.3  7.9  8.0  9.1  9.2  9.5  9.7 10.1 68.0

Медиана теперь равна 9.2

median(apples2)

## [1] 9.2

Давайте для сравнения посмотрим на среднюю.

mean(apples)

## [1] 8.85

mean(apples2)

## [1] 15.4

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

Квантили делят распределение на четыре равные части, каждая из которых включает по 25% значений.

I квантиль отсекает как раз 25%. II - 50%. Это ни что иное, как медиана. III квантиль отсекает 75% значений.

Определим их с помощью команды quantile

quantile(apples)

##    0%   25%   50%   75%  100% 
##  7.30  7.97  9.15  9.55 10.10

Нам показали так же минимальную и максимальную величины.

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

summary(apples)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    7.30    7.97    9.15    8.85    9.55   10.10

summary(apples2)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##     7.3     8.0     9.2    15.4     9.7    68.0

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

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

quantile(apples, probs = c(0.1, 0.99))

##   10%   99% 
##  7.72 10.07

В файле data\2.1_diatome_count.csv даны количества диатомовых водорослей в желудках этих животных. Прочитаем эти данные и посмотрим на их структуру:

diatoms <- read.table("data/2.1_diatome_count.csv", 
                      header = T, sep = "")
head(diatoms)

##          sp count
## 1 host_free    10
## 2 host_free     0
## 3 host_free     1
## 4 host_free     0
## 5 host_free     2
## 6 host_free     0

В таблице 2 переменные: вид (sp) и число водорослей в желудке (count). В переменной sp есть три варианта значений: "host_free" "host_symbiont" и "symbiont"

Ваша задача рассчитать 5-number summary для количества диатомовых в желудках хозяев и симбионтов.

summary(diatoms$count[diatoms$sp == "host_free"])

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##     0.0     1.0     4.0     4.7     8.0    12.0       1

summary(diatoms$count [diatoms$sp == "host_symbiont"])

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00    0.25    4.50    5.35    8.00   20.00

summary(diatoms$count [diatoms$sp == "symbiont"])

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##     0.0     0.0     1.0     1.5     2.0     7.0       4

Формат данных несколько сложен для человеческого глаза, зато очень подходит для ggplot.

library("ggplot2")
ggplot(data = diatoms, aes(y = count, x = sp)) + geom_boxplot()

Главный плюс, но так же и минус связки медиана + квантили это ее независимость от формы распределения.

Будь оно симметричным или с хвостом, 5-number summary опишет, а боксплот нарисует его без искажений. Однако, далеко на этих характеристиках не уедешь.

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

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

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

Если:

• данные в выборке распределены нормально
• мы можем узнать или оценить параметры этого распределения

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

Так что же это за чудодейственные параметры? Средняя и стандартное отклонение конечно!

Арифметическое среднее

Формула

\[\bar{x}=\frac{\sum{x_i}}{n}\]

Рассчитаем вручную и проверим:

mean_apple_diam <- sum(apples) / length(apples)
mean_apple_diam

## [1] 8.85

round(mean_apple_diam,3) == round(mean(apples),3)

## [1] TRUE

Девиата (отклонение)

Первым кирпичиком в построении параметрических мер разброса будет разность между значением вариаты (измерения) и средней:

\[x_i - \bar{x}\]

raw.deviates <- c(apples - mean(apples))
raw.deviates

## [1] -1.55  0.85  1.25  0.65 -0.85  0.35 -0.95  0.25

Как теперь из этого вектора значений получить одну числовую характеристику?

К сожалению мы не можем просто сложить все значения девиат и поделить их на объем выборки. Сумма девиат всегда будет равна нулю.

round(sum(raw.deviates))

## [1] 0

Путь 1: среднее отклонение = Average Deviation

Первый выход, который бросается в глаза, это использовать абсолютные значения

\[AD = \frac{\sum{\lvert x_i - \bar{x} \rvert}}{n}\]

sum(abs(raw.deviates)) / length(apples)

Путь 2: сумма квадратов = SS, Sum of Squares

Избавиться от знака девиаты можно не только с помощью модуля, но и возведя значение в квадрат.

\[SS = \sum{{(x_i - \bar{x})}^2}\]

sum((raw.deviates)^2)

Дисперсия = MS, Mean Square, Variance

Если мы теперь поделим сумму квадратов на объем выборки, то получим дисперсию.

\[variance=\frac{\sum{(x_i - \bar{x})^2}}{n}= \frac{SS}{n}\]

sum((raw.deviates)^2) / length(apples)

Среднеквадратичное/стандартное отклонение = Standard Deviation

Квадратный корень из дисперсии позволит вернуться к исходным единицам измерения

\[SD=\sqrt{\frac{\sum{(x_i - \bar{x})^2}}{n}} = \sqrt{variance}\]

sqrt(sum((raw.deviates)^2) / length(apples))

• только вместе,

• чувствительны к выбросам,

• плохо работают с несимметричными распределениями.

Из пяти положительных чисел создайте выборку со средней = 10 и медианой = 7

В выборке с медианой = 7 и n = 5, мы точно знаем: (а) одно из значений должно быть равно 7, (б) два значения должны быть меньше, и два - больше 7.

Создадим вектор, в котором одно значение задано, а три других просто придумаем:

example <- c(2, 5, 7, 10)

Среднее это сумма всех значений выборки, поделенная на ее объем. Умножив среднюю на 5 получим сумму всех значений.

Определим недостающее и проверим себя:

10 * 5 - sum(example)

## [1] 26

example <- c(2, 5, 7, 10, 26) #перезапишем вектор
mean(example)

## [1] 10

Статистики

• описывают выборки
• их обозначают латинскими буквами: \(М\), \(s^2\).

Параметры

• описывают генеральную совокупность
• их обозначают греческими буквами: \(\mu\), \(\sigma^2\)

С помощью выборочных статистик мы стремимся описать популяционные параметры.

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

Сумма "сырых" девиат всегда равна нулю, поэтому зная n-1 значений, мы без труда рассчитаем каким должно быть последнее. А значит, оно не может варьировать "свободно".

Дисперсия (выборочная оценка параметра)

\[\sigma^2=\frac{\sum{(x_i - \bar{x})^2}}{n-1}\]

Стандартное отклонение (выборочная оценка параметра)

\[\sigma=\sqrt{\frac{\sum{(x_i - \bar{x})^2}}{n-1}}\]

Давайте проверим как работают степени свободы на примере.

Создадим простейшую выборку со следующими признаками:

sim.mean <- 100
sim.sd <- 2
sim.n <- 11

Сначала рассчитаем дисперсию и сумму квадратов

sim.var <- sim.sd^2
sim.ss <- sim.var * (sim.n - 1)

SD = 2 это маленький разброс, так что надо выбирать значения недалеко от среднего. Например, 101 - очень удобное число. Заполним ими нашу выборку почти до конца:

sim.sample <- rep(101,10)

Итак, 10 значений мы взяли практически с потолка. Последнее, одиннадцатое, должно быть предсказано, иначе наша затея не сработает, и выборка с желаемыми признаками не получится.

Рассчитаем значение суммы квадратов для уже созданных нами "свободных" вариат:

free.ss <- (10) * ((101-100)^2)
free.ss

## [1] 10

Из общей суммы квадратов вычтем уже имеющуюся:

determined.ss <- sim.ss - free.ss

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

sqrt(determined.ss)

## [1] 5.48

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

last.value <- 100 - sqrt(determined.ss)
last.value

## [1] 94.5

Добавим последнее, "детерминированное", значение к остальным:

sim.sample <- c(rep(101,10), last.value)
sim.sample

##  [1] 101.0 101.0 101.0 101.0 101.0 101.0 101.0 101.0 101.0 101.0  94.5

Можем теперь проверить себя:

round(mean(sim.sample))

## [1] 100

round(sd(sim.sample))

## [1] 2

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

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

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

Давайте для наших экспериментов создадим генеральную совокупность с параметрами \(\mu = 50\) и \(\sigma = 7\).

Вот график ее частотного распределения.

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

Выборочные средние нормально распределены с параметрами \(\mu\) и \(SD _\bar{x}\)

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

\[SE _\bar{x} = \sigma / \sqrt{N}\]

Выборку должны характеризовать центральная тенденция (средняя) и мера разброса.

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

Поэтому лучше характеризовать данные в формате

\[Mean \pm SD (N)\]

а не

\[Mean \pm SE (N)\]

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