Размножаем реальность bootsrapping два метода
Бутстреп-процедура (или bootstrap) была предложена (Efron, 1979) как некоторое обобщение алгоритма «складного ножа», чтобы не уменьшать каждый раз число элементов по сравнению с исходной совокупностью. По одной из версий слово «bootstrap» означает кожаную полоску в виде петли, прикрепляемую к заднику походного ботинка для облегчения его натягивания на ногу. Благодаря этому термину появилась английская поговорка 1930-х г.: «Lift oneself by the bootstrap», которую можно трактовать как «Пробить себе дорогу благодаря собственным усилиям» (или подобно барону Мюнхгаузену вытянуть себя из болота за шнурки от ботинок). Основная идея бутстрепа по Б. Эфрону (1988) состоит в том, чтобы методом статистических испытаний Монте-Карло многократно извлекать повторные выборки из эмпирического распределения. А именно: берется конечная совокупность из n членов исходной выборки x1, x2, …, xn-1, xn , откуда на каждом шаге из n последовательных итераций с помощью датчика случайных чисел, равномерно распределенных на интервале [1, n], «вытягивается» произвольный элемент xk, который снова «возвращается» в исходную выборку (т. е. может быть извлечен повторно). Например, при n = 8 одна из таких комбинаций имеет вид x4, x2, x8, x2, x1, x2, x4, x5, т. е. отдельные элементы могут повторяться. Этим способом можно сформировать любое, сколь угодно большое число бутстреп-выборок.
В языке R бутстреп-процедуру можно реализовать, написав собственно свою функцию (которая и будет извлекать повторные выборки из исходных данных эмпирического распределения, создавая новые наборы данных), либо воспользоваться уже готовой библиотекой boot. Давайте рассмотрим описанное выше на конкретном примере и из одной выборки (x1) получим две бутстреп выборки, которые и сравним между собой.
library(boot)
set.seed(455)
options(digits = 10)
x1 = rnorm(50, 200, 5) # выборка исходных данных
Rboot = 50 # число повторов бутстрепинга
# второй способ применим функцию boot
xb1 <- boot(x1, function(u,i) (u[i]), R = Rboot)
X2 = as.vector(xb1$t)
summary(x1)
summary(X1)
summary(X2)
> summary(x1) Min. 1st Qu. Median Mean 3rd Qu. Max. 186.9704 199.0556 201.8258 202.0853 205.3851 213.4328 > summary(X1) Min. 1st Qu. Median Mean 3rd Qu. Max. 186.9704 199.3623 201.8307 202.2443 205.6955 213.4328 > summary(X2) Min. 1st Qu. Median Mean 3rd Qu. Max. 186.9704 199.3623 201.8307 202.1562 205.6955 213.4328 |
Видно, что описательная статистика между вектором исходных х1 данных и двух векторов (X1,X2) полученных в результате бутстрепинга практически не отличаются. Проведем тест сравнения двух групп (X1,X2) полученных из бутстрепинга между собой, приняв нулевую гипотезу (группы статистически не различаются) на уровне значимости 0.95
t.test(X1,X2)
> t.test(X1,X2) Welch Two Sample t-test data: X1 and X2 t = 0.6279096, df = 4997.9039, p-value = 0.5300918 alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval: -0.1869445647 0.3631272587 sample estimates: mean of x mean of y 202.2443309 202.1562396 |
Значение p-value = 0.5300918 на высшем уровне показывает что нулевая гипотеза остается в силе: статистически значимых различий в группах (X1,X2) полученных в результате бутстрепинга реализованного двумя различными способами нет. Так что можно применять на практике оба способа бутстреп-процедуры: в виде функции или в виде библиотеки boot. Полный скрипт :
library(boot)
set.seed(455)
options(digits = 10)
x1 = rnorm(50, 200, 5)
Rboot = 50
# первый способ напишем свою функцию
Xboot = function(data){
n = length(data)
boot = sample(n, replace=TRUE)
data_boot = (data[boot])
return(data_boot)
}
X1 = c()
X1 = as.vector(replicate(Rboot, (Xboot(data = x1))))
# второй способ применим функцию boot
xb1 <- boot(x1, function(u,i) (u[i]), R = Rboot)
X2 = as.vector(xb1$t)
summary(x1)
summary(X1)
summary(X2)
t.test(X1,X2)