Доверительный интервал для линии регрессии mle2

Question

Я выполняю нелинейную регрессию с mle2 в R, и я хочу создать 95% -ную точную полосу доверия вокруг кривой наилучшего соответствия. Ниже приведен упрощенный пример того, что я пытаюсь сделать, когда я пытаюсь использовать прогноз для объекта mle2 fit. Любые предложения о том, как это сделать?

library(bbmle) 

# Fabricated data 
e.u <- function(x, k) { exp(-k * x) } 
n <- 40 
t.bio <- 1:n 
bio <- 10*e.u(t.bio,log(2)/10) + rnorm(n,0,sqrt(e.u(t.bio,log(2)/10))) 

#Use mle2 to estimate the parameters 
intake.guess <- 10 
rc.guess <- 0.07 

n.log.like <- function(intake,k) { 
    sum.y <- 0 
    for (i in 1:length(bio)) { 
    x <- intake * e.u(t.bio[i],k) 
    y <- bio[i] 
    sum.y <- sum.y + log(dnorm(y,x,0.1*sqrt(y))) } 
    return(-sum.y) 
} 

b <- mle2(n.log.like, 
     start=list(
      intake=intake.guess, 
      k=rc.guess), 
     data=list(
      t.bio=t.bio, 
      bio=bio), 
     method="Nelder-Mead", 
     skip.hessian=FALSE) 

intake <- coef(summary(b))[1,1]   
rc <- coef(summary(b))[2,1]   
summary(b) 

#Scatter plot 
bio.p <- numeric(n) 
x <- 1:n 
for (i in 1:n) { bio.p[i] <- intake * e.u(x[i],rc) } 
plot(x,bio.p,type="l",log="xy",main="", 
    xlab="Days After Intake",ylab="Excretion") 
points(t.bio,bio) 

#I want to generate a confidence interval on the regression line 
bio.hat <- predict(b) 

Answer 1

Вы можете попробовать этот ресурс:

https://stat.ethz.ch/R-manual/R-patched/library/stats/html/predict.lm.html

Существует se.fit параметр, передаваемый с predict(). Если TRUE, он рассчитает SE. Я не пробовал.

Другой способ - использовать library(ggplot2), который с одним из параметров позволяет получить уровень доверия, наложенный на ваш график автоматически. Пример:

c <- ggplot(mtcars, aes(qsec, wt)) 
c + stat_smooth(se = TRUE) + geom_point() 

Это будет помещать полосу вокруг участка, чтобы указать доверительный интервал. По умолчанию значение se установлено в TRUE.

Ссылка: http://svitsrv25.epfl.ch/R-doc/library/ggplot2/html/stat_smooth.html

В-третьих, вы можете рассчитать доверительный интервал и накладываться на участок, установив пар (новый = TRUE). Это может быть громоздким, но все же осуществимым.