On 2010-09-02 22:32, Gabor Grothendieck wrote:
On Thu, Sep 2, 2010 at 7:39 PM, Marlin Keith Cox<marlink...@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"); fmNonlinear 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
Keith, You would get Gabor's solution from your nls() call if you did: 1. replace 'Joules' with 'joules' 2. replace 'c=-.229' with 'c=.229' in your start vector. You already have the 'minus' in your formula. I find it useful to add the curve you get with your start values to the scatterplot to see how reasonable the values are: plot(Joules~Days) curve(8 + 9 * exp(-.229 * x), add=TRUE) After you get a convergent solution, you can add a curve with the model values. -Peter Ehlers ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
