Доверительные интервалы для летальной дозы (LD) для логистической регрессии в R
Я хочу найти смертельную дозу (LD50
) с доверительным интервалом вR
, Другие линии программного обеспечения Minitab, SPSS, SAS предоставляют три разные версии таких доверительных интервалов. Я не мог найти такие интервалы в любом пакете вR
(Я также использовалfindFn
функция отsos
пакет).
Как я могу найти такие интервалы? Я кодировал для одного типа интервалов на основе метода Дельта (как не уверен в его правильности), но хотел бы использовать любую установленную функцию изR
пакет. Спасибо
MWE:
dose <- c(10.2, 7.7, 5.1, 3.8, 2.6, 0)
total <- c(50, 49, 46, 48, 50, 49)
affected <- c(44, 42, 24, 16, 6, 0)
finney71 <- data.frame(dose, total, affected)
fm1 <- glm(cbind(affected, total-affected) ~ log(dose),
family=binomial(link = logit), data=finney71[finney71$dose != 0, ])
summary(fm1)$coef
Estimate Std. Error z value Pr(>|z|)
(Intercept) -4.886912 0.6429272 -7.601035 2.937717e-14
log(dose) 3.103545 0.3877178 8.004650 1.198070e-15
library(MASS)
xp <- dose.p(fm1, p=c(0.50, 0.90, 0.95)) # from MASS
xp.ci <- xp + attr(xp, "SE") %*% matrix(qnorm(1 - 0.05/2)*c(-1,1), nrow=1)
zp.est <- exp(cbind(xp, attr(xp, "SE"), xp.ci[,1], xp.ci[,2]))
dimnames(zp.est)[[2]] <- c("LD", "SE", "LCL","UCL")
zp.est
LD SE LCL UCL
p = 0.50: 4.828918 1.053044 4.363708 5.343724
p = 0.90: 9.802082 1.104050 8.073495 11.900771
p = 0.95: 12.470382 1.133880 9.748334 15.952512