Bootstrap для оценки дов.интервалов

Вопросы по статистическому пакету R. Не обязательно гео.
Ответить
marratt
Интересующийся
Сообщения: 28
Зарегистрирован: 12 май 2012, 07:33
Репутация: 0

Bootstrap для оценки дов.интервалов

Сообщение marratt » 18 окт 2018, 12:13

Здравствуйте,
Передо мной стоит задача определить доверительные интервалы для оценки отношения двух выборочных переменных: R = x/y,
где х - число животных в данном квадрате с признаком L;
y - число животных в данном квадрате.
Распределение сильно ассиметрично.
Решил использовать bootstrap для определения доверительных интервалов полученной оценки.
Для этого использую пакет 'boot'. Исходные данные в виде data frame - lv, котороя состоит из 2 столбцов - х и у, соответственно.
Код следующий:

Код: Выделить всё

lv.boot <- boot(lv, sum(lv[,1])/sum(lv[,2]), R = 1000, stype = "f")
boot.ci(lv.boot, conf = c(0.90, 0.95),
        type = c("norm", "basic", "perc", "bca"))
Однако на этапе команды boot выходит ошибка:
Error in statistic(data, original, ...) :
could not find function "statistic"

Пробовал разные варианты, писать предварительные функции, но ничего не помогает. Не могу понять, где у меня ошибка.
Как сделать так, чтобы он сделал ресэмплинг исходных данных, посчитал их отношение и определил доверительные интервалы?
Спасибо

gamm
Гуру
Сообщения: 3015
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 705
Ваше звание: программист
Откуда: Казань

Re: Bootstrap для оценки дов.интервалов

Сообщение gamm » 18 окт 2018, 14:32

То, что нужно считать, нужно оформить в виде функции (это про sum(...) - у вас выражение). См. пример в документации.

А что мешает просто использовать glm(...,family="binomial")?

Получите среднее и доверительный интервал в предположении, что все квадраты однородны. Либо нужно добавлять в модель факторы (свойства квадратов)

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

marratt
Интересующийся
Сообщения: 28
Зарегистрирован: 12 май 2012, 07:33
Репутация: 0

Re: Bootstrap для оценки дов.интервалов

Сообщение marratt » 19 окт 2018, 11:47

Здравствуйте, спасибо за советы,
gamm писал(а):
18 окт 2018, 14:32
То, что нужно считать, нужно оформить в виде функции (это про sum(...) - у вас выражение). См. пример в документации.
Да, я попробовал оформить в виде функции:

Код: Выделить всё

> f3 <- function(lv, id) {
+   sum(lv$l)/sum(lv$m)}
> boot.out <- boot(lv,f3,1000)
После этого процесс пошёл, но я что-то неверное указал в функции, потому что он мне выдал:

Код: Выделить всё

ORDINARY NONPARAMETRIC BOOTSTRAP
Call:
boot(data = lv, statistic = f3, R = 1000)
Bootstrap Statistics :
     original  bias    std. error
t1* 0.2272727       0           0
Если он и делал ресэмплинг, то все 1000 проб дали одно и то же отношение - что вряд ли.
Видимо, у меня в функции есть ошибка.
gamm писал(а):
18 окт 2018, 14:32
А что мешает просто использовать glm(...,family="binomial")?
Я не настолько искушён в регрессионном анализе, даже и не подумал его использовать для оценки интервалов. У меня просто данные не дихотомические. У меня есть квадраты, в которых было несколько животных M и несколько животных с признаком L. То есть я не могу разделить всю выборку на два класса. А в биномиальной регрессии, насколько я понимаю, используется дихотомия признаков.
Квадраты у меня неоднородны - выборка была стратифицированная (4 страты).
gamm писал(а):
18 окт 2018, 14:32
А если уж очень хочется бутстрапа, кто мешает сделать его руками в цикле, выбирая индексы квадратов с возвратом. Это просто, и все будет под полным контролем.
Скорее всего так и сделаю, это, конечно, не очень изящно и эргономично, зато, действительно, всё на виду, всё под контролем.
Спасибо

gamm
Гуру
Сообщения: 3015
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 705
Ваше звание: программист
Откуда: Казань

Re: Bootstrap для оценки дов.интервалов

Сообщение gamm » 19 окт 2018, 13:20

marratt писал(а):
19 окт 2018, 11:47
Я не настолько искушён в регрессионном анализе, даже и не подумал его использовать для оценки интервалов. У меня просто данные не дихотомические. У меня есть квадраты, в которых было несколько животных M и несколько животных с признаком L.
данные у вас как раз дихотомические, каждое из М животных либо имеет признак, либо нет :D
Это совершенно классическая биномиальная модель, данные для нее готовятся в виде двух столбцов: two-column matrix with the columns giving the numbers of successes and failures, т.е. в первом столбце число животных с признаком, во втором - без него

Код: Выделить всё

M<-c(10,5,7,9)
L<-c( 2,4,2,3)
Y<-cbind(L,M-L)
cur.lm<-glm(Y~1,family="binomial")
summary(cur.lm)
# prediction on logit-scale
fitted.prob<-predict(cur.lm,newdata=data.frame(1),se.fit=TRUE)
# get fitted proportion and 95% CI
fitted.prop<-1/(1+exp(-fitted.prob$fit))
fitted.prop.lo<-1/(1+exp(-(fitted.prob$fit-qnorm(0.975)*fitted.prob$se.fit)))
fitted.prop.hi<-1/(1+exp(-(fitted.prob$fit+qnorm(0.975)*fitted.prob$se.fit)))
c(fitted.prop.lo,fitted.prop,fitted.prop.hi)
#get proportion
sum(L)/sum(M)
А теперь по существу задачи. У вас данные наверняка пространственно коррелированы между соседними квадратами, и лучше бы сочинить модель с пространственной корреляцией, и посмотреть, не меняется ли доля по территории или в зависимости от разных факторов вроде ландшафта или растительности ....

marratt
Интересующийся
Сообщения: 28
Зарегистрирован: 12 май 2012, 07:33
Репутация: 0

Re: Bootstrap для оценки дов.интервалов

Сообщение marratt » 19 окт 2018, 16:26

Здравствуйте, большое спасибо за столь ценные советы и за код для glm.
Вы знаете, я когда планировал это исследование, я думал просто расчитать долю животных с признаком L, а затем доверительные интервалы по Clopper-Pearson.
Но здесь есть один нюанс, который меня смутил: моя выборка включает 123 квадрата. Но животных мы нашли всего 44.
Если я расчитываю обычную долю и Clopper-Pearson CI, то дов.интервалы будут считаться исходя из n=44.
Я почитал У.Кокрена (Методы выборочного исследования, 1976), у него есть описания оценки отношения и стандартной ошибки для него. При таком подходе n=123. Вот только мои данные не нормально распределены и ассимитричны, поэтому эту оценку стандартной ошибки я не могу использовать для обычного расчёта дов.интервалов. Вот и решил попробовать бутстреп.
Но ваш совет о пространственной корреляции всё меняет, действительно, у меня есть все основания предполагать, что распределение признака L имеет пространственную приуроченность (по сути L - это наличие определённой инфекции в животном). Буду читать литературу по пространственному анализу (если регрессионный анализ я хоть чуть-чуть да понимаю, то с пространственным анализом у меня всё плохо - ничего о нём не знаю).
Спасибо ещё раз!

gamm
Гуру
Сообщения: 3015
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 705
Ваше звание: программист
Откуда: Казань

Re: Bootstrap для оценки дов.интервалов

Сообщение gamm » 19 окт 2018, 17:10

marratt писал(а):
19 окт 2018, 16:26
Но ваш совет о пространственной корреляции всё меняет, действительно, у меня есть все основания предполагать, что распределение признака L имеет пространственную приуроченность (по сути L - это наличие определённой инфекции в животном).
если можно, опишите, как собирались данные, и про кого.

marratt
Интересующийся
Сообщения: 28
Зарегистрирован: 12 май 2012, 07:33
Репутация: 0

Re: Bootstrap для оценки дов.интервалов

Сообщение marratt » 20 окт 2018, 11:33

gamm писал(а):
19 окт 2018, 17:10
если можно, опишите, как собирались данные, и про кого.
Отправил личным сообщением

Ответить

Вернуться в «R»

Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей и 1 гость