[R] NLS equation self starting non linear

Gabor Grothendieck ggrothendieck at gmail.com
Fri Sep 3 06:32:27 CEST 2010


On Thu, Sep 2, 2010 at 7:39 PM, Marlin Keith Cox <marlinkcox at gmail.com> wrote:
> This data are kilojoules of energy that are consumed in starving fish over a
> time period (Days).  The KJ reach a lower asymptote and level off and I
> would like to use a non-linear plot to show this leveling off.  The data are
> noisy and the sample sizes not the largest.  I have tried selfstarting
> weibull curves and tried the following, both end with errors.
>
> Days<-c(12, 12, 12, 12, 22, 22, 22, 22, 32, 32, 32, 32, 37, 38, 38, 39, 39,
> 39, 39, 39, 39, 39,  1,  1,  1,  1)
> Joules<-c(8.782010,  8.540524,  8.507526, 11.296904,  4.232690, 13.026588,
> 10.213342,  4.771482,  4.560359,  6.146684, 9.651727,  8.064329,  9.419335,
> 7.129264,  6.079012,  7.095888, 17.996794,  7.028287,  8.028352,  5.497564,
> 11.607090,  9.987215, 11.065307, 21.433094, 20.366385, 22.475717)
> X11()
> par(cex=1.6)
> plot(Joules~Days,xlab="Days after fasting was initiated",ylab="Mean energy
> per fish (KJ)")
> model<-nls(joules~a+b*exp(-c*Days),start=list(a=8,b=9,c=-.229),
> control=list(minFactor=1e-12),trace=TRUE)
> summary(model)

Note that Joules is defined above but joules is used as well.  We have
assumed they are the same.

Also note that the formula is linear in log(joules) if a=0 so lets
assume a=0 and use lm to solve.  Also switch to the "plinear"
algorithm which only requires starting values for the nonlinear
parameters -- in this case only c (which we get from the lm).

> joules <- Joules
> coef(lm(log(joules) ~ Days))
(Intercept)        Days
 2.61442015 -0.01614845
>
> # so use 0.016 as starting value of c
> fm <- nls(joules~cbind(1, exp(-c*Days)), start = list(c = 0.016), alg = "plinear"); fm
Nonlinear regression model
  model:  joules ~ cbind(1, exp(-c * Days))
   data:  parent.frame()
      c   .lin1   .lin2
 0.2314  8.3630 13.2039
 residual sum-of-squares: 290.3

-- 
Statistics & Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com



More information about the R-help mailing list