Выборки из распределений

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

R позволяет также гененировать случайные выборки из различных распределений заданного объема. Для этого используются функции типа rdist(), где вместо dist указывается необходимое распределение. Первый аргумент этой функции всегда объем выборки, а затем – параметры распределения. Доступные распределения можно найти здесь.

Бинарное распределение

Например, сгененируем выборку объема 8 из бинарного распределения с параметром \(p=\frac{1}{4}\). Для этого нужно попросить R сгенерировать выборку из биномиального распределения с параметрами \(n=1\) и \(p=\frac{1}{4}\) (так как бинарное распределение – это частный случай биномиального распределения).

rbinom(8, 1, 1/4)
## [1] 0 0 0 0 0 0 1 0

Каждый раз при реализации функции, генерирующей случайную выборку, у Вас будет появляться разная последовательность значений. Так и должно быть, ведь каждое из наблюдений выборки – случайная величина. В среднем, мы ожидаем, что четвертая часть выборки, то есть 2 наблюдения, должна содержать значения 1. Но при малельком размере выборки возможными также будут выборки, где 1 встречается чаще или реже.

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

Мы можем примерно оценить, как часто нам выпадет выборка, где содержится 1 единица. Это вероятность одного успеха в серии из 8 испытаний Бернулли с тем же самым параметром вероятности успеха:

dbinom(1, 8, 1/4)
## [1] 0.2669678
dbinom(0, 8, 1/4)
## [1] 0.1001129
dbinom(8, 8, 1/4)
## [1] 1.525879e-05

Вы видим, что примерно в 10 случаях из 100 выборка не будет содержать ни одной 1, а вероятность сплошных единиц очень мала – следовательно, такая выборка (одни единицы) неправдоподобна.

Примечание. Для того, чтобы воспроизвести сгенерированную случайную выборку при повторной реализации функции, нужно прописать перед ней функцию set.seed(any_number). Например:

set.seed(111)

Нормальное распределение

Задание 3 из семинара 10

to be updated

Гистограммы в R

Гистограмма – график, который является аналогом графика функции плотности для распределний. Все наблюдения в выборке упорядочиваются по возрастанию. Затем упорядоченная выборка начиная от некоторого стартового значения (иногда это минимальное значение в выборке, а иногда – некоторое значение, меньшее, чем минимум) разбивается на \(n\) равных интервалов. После этого подсчитывается количество наблюдений, которые попали в каждый интервал. Это количество откладывется на оси Y. По оси X откладываются значения наблюдений в выборке.

Вызовем встроенную в R базу данных, где содержатся измерения температуры тела двух бобров в течение дня. После реализации команды data() в памяти R (окошко Environment) появятся две базы данных – по одной на бобра. Описание базы данных можно найти здесь.

data(beavers)

Посмотрим на базу данных. Выведем первые 6 строчек базы данных beaver1 – это делается с помощью функции head().

head(beaver1)

Нам нужна колонка (по-другому колонки в базе данных называются переменными) temp. Чтобы обратиться к колонке внутри базы данных, нужно сначала прописать название базы данных, затем поставить знак $, а затем прописать название переменной.

beaver1$temp
##   [1] 36.33 36.34 36.35 36.42 36.55 36.69 36.71 36.75 36.81 36.88 36.89
##  [12] 36.91 36.85 36.89 36.89 36.67 36.50 36.74 36.77 36.76 36.78 36.82
##  [23] 36.89 36.99 36.92 36.99 36.89 36.94 36.92 36.97 36.91 36.79 36.77
##  [34] 36.69 36.62 36.54 36.55 36.67 36.69 36.62 36.64 36.59 36.65 36.75
##  [45] 36.80 36.81 36.87 36.87 36.89 36.94 36.98 36.95 37.00 37.07 37.05
##  [56] 37.00 36.95 37.00 36.94 36.88 36.93 36.98 36.97 36.85 36.92 36.99
##  [67] 37.01 37.10 37.09 37.02 36.96 36.84 36.87 36.85 36.85 36.87 36.89
##  [78] 36.86 36.91 37.53 37.23 37.20 37.25 37.20 37.21 37.24 37.10 37.20
##  [89] 37.18 36.93 36.83 36.93 36.83 36.80 36.75 36.71 36.73 36.75 36.72
## [100] 36.76 36.70 36.82 36.88 36.94 36.79 36.78 36.80 36.82 36.84 36.86
## [111] 36.88 36.93 36.97 37.15

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

hist(beaver1$temp)

Настроим гистограмму.

hist(beaver1$temp,
     xlab = "Темература тела бобра", # изменить подпись оси OX
     ylab = "Количество наблюдений", # изменить подпись оси OY
     main = "Моя первая сливовая гистограмма про бобров", # изменить название гистограммы
     col = "plum" # поменять цвет гистограммы
     )

Дополнительно: кривая нормального распределения

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

Прежде всего поменяем шкалу оси OY, указав аргумент freq = F. Тогда по оси Y будут отложены значения плотности: так, чтобы площадь под гистограммой была равна 1.

hist(beaver1$temp,
     xlab = "Темература тела бобра", # изменить подпись оси OX
     ylab = "Плотность", # изменить подпись оси OY
     main = "Моя первая сливовая гистограмма про бобров", # изменить название гистограммы
     col = "plum", # поменять цвет гистограммы
     freq = F # плотность
     )