[R] Confidence intervals with nls()
Ben Bolker
bolker at ufl.edu
Sat Aug 2 19:04:42 CEST 2008
My best guess is that the code you were copying
was from an objective (fitted) function that computed
the gradient itself, so the "predict" result had
a "gradient" attribute. Here's an alternative that
(I think) works OK, but it comes with the usual
warranty (none).
## simulate "data"
O.age = rep(1:12,each=5)
k = 0.2
Linf = 700
t0 = 0
O.length = round(rnorm(length(O.age),
Linf*(1-exp(-k*O.age-t0)),sd=5))
OS = data.frame(O.age,O.length)
plot(O.length~O.age, data=OS)
## fit model
Oto = nls(O.length~Linf*(1-exp(-k*(O.age-t0))), data=OS,
start=list(Linf=1000, k=0.1, t0=0.1), trace=TRUE)
mod <- seq(0, 12)
mod=seq(1,12, by=0.01)
lines(mod, predict(Oto, list(O.age = mod)))
## use emdbook::deltavar to use delta method
## to compute standard errors
library(emdbook)
se.fit = sqrt(deltavar(
Linf*(1-exp(-k*mod-t0)),
meanval=coef(Oto),
Sigma=vcov(Oto)))
plot(O.length~O.age, data=OS)
matlines(mod, predict(Oto,list(O.age=mod))+
outer(se.fit,qnorm(c(.5, .025,.975))),type="l",
lty=1,col=c("black","gray","gray"))
I notice that I have messed up the VB formula
here -- it should be -k*(mod-t0) -- but you can
correct that
More information about the R-help
mailing list