Доверительные интервалы для прогнозов логистической регрессии
в R предсказать.lm вычисляет Прогнозы на основе результатов линейной регрессии, а также предлагает вычислить доверительные интервалы для этих прогнозов. Согласно руководству, эти интервалы основаны на дисперсии ошибок подгонки, но не на интервалах ошибок коэффициента.
с другой стороны предсказать.glm, который вычисляет Прогнозы на основе логистической и пуассоновской регрессии (среди нескольких других), не имеет возможности для доверительных интервалов. И я даже трудно представить себе, как такие доверительные интервалы могут быть вычислены, чтобы обеспечить значимое понимание для Пуассона и логистической регрессии.
есть ли случаи, в которых имеет смысл предоставлять доверительные интервалы для таких прогнозов? Как их можно интерпретировать? И каковы предположения в этих случаях?
1 ответ:
обычным способом является вычисление доверительного интервала в масштабе линейного предиктора, где все будет более нормальным (гауссовым), а затем применить обратную функцию связи для отображения доверительного интервала от масштаба линейного предиктора до масштаба ответа.
для этого вам нужно две вещи;
- вызов
predict()
Сtype = "link"
и- вызов
predict()
Сse.fit = TRUE
.первый производит Прогнозы на масштаб линейного предиктора, второй возвращает стандартные ошибки предсказаний. В псевдокоде
## foo <- mtcars[,c("mpg","vs")]; names(foo) <- c("x","y") ## Working example data mod <- glm(y ~ x, data = foo, family = binomial) preddata <- with(foo, data.frame(x = seq(min(x), max(x), length = 100))) preds <- predict(mod, newdata = preddata, type = "link", se.fit = TRUE)
preds
это тогда список с компонентамиfit
иse.fit
.доверительный интервал на линейном предикторе тогда
critval <- 1.96 ## approx 95% CI upr <- preds$fit + (critval * preds$se.fit) lwr <- preds$fit - (critval * preds$se.fit) fit <- preds$fit
critval
выбирается из t или z (нормальное) распределение по мере необходимости (я точно забыл, что использовать для какого типа GLM и какие свойства) с необходимым покрытием. Элемент1.96
- это значение гауссовского распределения, дающее 95% покрытия:> qnorm(0.975) ## 0.975 as this is upper tail, 2.5% also in lower tail [1] 1.959964
теперь
fit
,upr
иlwr
нам нужно применить к ним обратную функцию связи.fit2 <- mod$family$linkinv(fit) upr2 <- mod$family$linkinv(upr) lwr2 <- mod$family$linkinv(lwr)
теперь вы можете построить все три и данные.
preddata$lwr <- lwr2 preddata$upr <- upr2 ggplot(data=foo, mapping=aes(x=x,y=y)) + geom_point() + stat_smooth(method="glm", method.args=list(family=binomial)) + geom_line(data=preddata, mapping=aes(x=x, y=upr), col="red") + geom_line(data=preddata, mapping=aes(x=x, y=lwr), col="red")