[R] confidence interval for nls

[Henrik Andersson]

I tried the example and it works fine,

but why o why, do I not get any gradient from another prediction?


############################################
#example code
############################################
yran <-
   c(0.0118821538311157, 0.323340819569374, 0.525108551769669, 
0.648173528359086,
     0.674941464557777, 0.745619897023202, 0.675263947547285, 
0.782342006150897,
     0.804472377613314, 0.813517628865957, 0.83178876929294, 
0.8353033444236,
     0.782341888336996, 0.782751195243735, 0.951449319233294, 
0.753939045638585,
     0.866298051948752, 0.90725708524503, 1.04254520185355, 
0.941419545223617,
     0.900965057020577)

x <- seq(0,10,length=21)

##plot(yran~x)
mich <- function(x,K,rmax) rmax*x/(x+K)
mm.nls <- nls(yran~mich(x,K,rmax),start=list(K=5,rmax=3))

xx=seq(0,10,length=200)
predict(mm.nls,newdata=list(x = xx))  #produces no gradient

mm.nls$m$gradient()                    #this does, but only 
at the 						       #observed data points

Peter Dalgaard wrote:

> Henrik Andersson <h.andersson at nioo.knaw.nl> writes:
> 
> 
>>Do I have the right impression that it's currently not possible to
>>produce confidence intervals for the nls predictions using R?
>>
>>I had a course were we used SAS PROC nlin and there you could get
>>intervals for the parameters and the prediction but I do not have
>>access to SAS.
>>
>>Would it be difficult to implement, I tried to dig into the help pages
>>of nls, vcov and nlsModel but I could not really make sense out of
>>this?
> 
> 
> It's pretty trivial (if you stick with the linear approximation), once
> you realize that predict.nls actually returns the gradient in an
> attribute:
> 
> example(predict.nls)
> se.fit <- sqrt(apply(attr(predict(fm,list(Time = tt)),"gradient"),1, 
>                   function(x) sum(vcov(fm)*outer(x,x))))
> matplot(tt, predict(fm,list(Time = tt))+
>                outer(se.fit,qnorm(c(.5, .025,.975))),type="l")
> 
> points(demand ~ Time, data = BOD)
> 
> One slight issue is that it doesn't work if "newdata" is omitted, but
> then you can easily get the gradient from fm$m$gradient()
> 


-- 
---------------------------------------------
Henrik Andersson
Netherlands Institute of Ecology -
Centre for Estuarine and Marine Ecology
P.O. Box 140
4400 AC Yerseke
Phone: +31 113 577473
h.andersson at nioo.knaw.nl
http://www.nioo.knaw.nl/ppages/handersson