Среднее, медиана и стандартное отклонение в R

Среднее и медиана

Две самые используемые сводки во всей статистике — по одному вызову функции каждая:

mpg <- mtcars$mpg

print(mean(mpg))
print(median(mpg))

mean() складывает всё и делит на количество. median() сортирует значения и берёт среднее по положению (или усредняет два средних, когда количество чётное). Для 32 машин в mtcars средняя экономичность составляет около 20.1 мили на галлон, а медиана 19.2 — они близки, и это говорит, что данные не сильно скошены. Когда же они сильно расходятся, это тоже информация; к этому мы вернёмся в конце.

Одна вещь удивляет всех: если вектор содержит хотя бы одно NA, обе функции возвращают NA:

x <- c(4, 8, NA, 6)

print(mean(x))
print(mean(x, na.rm = TRUE))

Это сделано намеренно — R отказывается тихо делать вид, что пропущенного значения нет. na.rm = TRUE говорит: «считай по тем значениям, что есть». Почти каждая сводная функция в этой статье его принимает. Полная история о том, как распространяется NA, — в статье о пропущенных значениях.

Стандартное отклонение и дисперсия

sd() измеряет разброс: грубо говоря, насколько типичное значение отстоит от среднего, в тех же единицах, что и данные. var() — её квадрат:

mpg <- mtcars$mpg

print(sd(mpg))
print(var(mpg))
print(sd(mpg)^2)

Стандартное отклонение около 6 миль на галлон означает, что машины обычно укладываются примерно в 6 единиц от среднего значения 20.1. Поскольку sd() выражена в собственных единицах данных, именно её и сообщают; var() в основном всплывает внутри других формул.

Вот деталь, важная для учёбы: sd() и var() вычисляют выборочную статистику — они делят сумму квадратов отклонений на n − 1, а не на n:

x <- c(2, 4, 4, 4, 5, 5, 7, 9)
n <- length(x)

sample_var <- sum((x - mean(x))^2) / (n - 1)
population_var <- sum((x - mean(x))^2) / n

print(sample_var)      # what var(x) gives you
print(var(x))
print(population_var)  # divide by n instead

Почему n − 1? Потому что среднее вы оценили по тем же данным, и отклонения вокруг этого оценённого среднего систематически чуть меньше нужного; деление на n − 1 это корректирует. Поскольку ваши данные почти всегда являются выборкой из чего-то большего, вам нужна версия с n − 1. Если у вас действительно вся генеральная совокупность (каждый студент в классе, каждый товар в каталоге), домножьте: var(x) * (n - 1) / n.

Стандартная ошибка среднего

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

mpg <- mtcars$mpg

se <- sd(mpg) / sqrt(length(mpg))
print(se)

Стандартная ошибка уменьшается с ростом выборки — соберите вчетверо больше данных, и SE уменьшится вдвое, — потому что большая выборка точнее пришпиливает среднее. А sd с ростом выборки не уменьшается: машины настолько разнообразны, насколько они разнообразны, сколько бы вы их ни измерили. SE — строительный блок доверительных интервалов, где она и окупается.

summary() — обзор одним вызовом

summary() за один раз даёт пятичисловую сводку плюс среднее и работает на целых датафреймах:

print(summary(mtcars$mpg))

print(summary(mtcars[, c("mpg", "wt", "hp")]))

Вызванная на датафрейме, она сводит каждый столбец: числовые получают минимум, квартили, среднее и максимум, а факторы — количества по уровням. Это первое, что стоит запускать на любом только что загруженном наборе данных: невозможные значения (отрицательный возраст, максимум 9999) сразу бросаются в глаза.

Квантили, размах и IQR

quantile() обобщает медиану на любую точку разреза:

mpg <- mtcars$mpg

print(quantile(mpg))
print(quantile(mpg, probs = c(0.1, 0.5, 0.9)))

print(range(mpg))
print(IQR(mpg))
print(max(mpg) - min(mpg))

Без аргументов она возвращает минимум, квартили и максимум. Передайте probs = для конкретных точек разреза — 10-й и 90-й процентили выше очерчивают, где живёт основная масса данных. IQR() (расстояние между 25-м и 75-м процентилями) — мера разброса, которая, в отличие от sd(), не таскается за выбросами. range() возвращает минимум и максимум парой.

Мода: mode() в R делает НЕ это

Эта ловушка ловит всех ровно один раз. В R есть функция mode(), и к статистике она не имеет никакого отношения — она сообщает тип хранения объекта:

x <- c(3, 7, 7, 2, 7, 5, 3)

print(mode(x))   # "numeric" - not what you wanted

counts <- table(x)
print(counts)
print(names(which.max(counts)))

Идиома, которую стоит запомнить: table(x) считает, как часто встречается каждое значение, which.max() находит наибольшее количество, а names() вытаскивает само значение. Учтите, что оно возвращается символьной строкой (имена таблицы всегда таковы); оберните в as.numeric(), если нужно с ним считать. Если два значения делят первое место, which.max() молча вернёт только первое — проверяйте таблицу сами, когда ничьи правдоподобны.

Среднее против медианы: что сообщать

Среднее использует каждое значение, и это его сила и его слабость — одно экстремальное значение тянет его за собой. Медиану интересует только середина, поэтому выбросы её почти не трогают:

salaries <- c(42000, 45000, 47000, 48000, 51000, 54000, 58000)
print(mean(salaries))
print(median(salaries))

# The CEO joins the dataset:
salaries <- c(salaries, 900000)
print(mean(salaries))
print(median(salaries))

Одно добавленное значение толкает среднее с примерно 49 300 до более чем 155 000 — числа, которое не описывает никого в данных, — тогда как медиана сдвигается лишь с 48 000 до 49 500. Поэтому доходы, цены на жильё и длительность госпитализации сообщают медианами: скошенные данные с длинным хвостом делают среднее вводящим в заблуждение. Для примерно симметричных данных обе меры сходятся, и среднее вполне годится (и статистически эффективнее). Быстрая гистограмма подскажет, в какой вы ситуации, — а сравнение вашего среднего с медианой само по себе является однострочной проверкой на скошенность.

Что вы уносите с собой

• mean(x) и median(x) — добавляйте na.rm = TRUE, когда есть пропущенные значения.
• sd(x) и var(x) считают выборочную статистику (знаменатель n − 1) — а это то, что вам нужно.
• Стандартная ошибка это sd(x) / sqrt(length(x)) — она измеряет, насколько хорошо вы знаете среднее, а не насколько разбросаны данные.
• summary() на свежем датафрейме — самая быстрая проверка здравого смысла в R.
• Мода это names(which.max(table(x))) — mode() в R про типы хранения.
• Скошенные данные или выбросы: сообщайте медиану. Симметричные данные: среднее вполне годится.

Дальше: измерение того, как две переменные движутся вместе, — корреляция через cor() и cor.test().