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"); 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
______________________________________________
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.
