On 26/08/10 19:48, David Winsemius wrote:
On Aug 26, 2010, at 1:35 PM, Marlin Keith Cox wrote:
I need the parameters estimated for a non-linear equation, an example
of the
data is below.
# rm(list=ls()) I really wish people would add comments to destructive
pieces of code.
Time<-c( 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 4, 4, 4,
4, 4, 5, 5, 5, 5, 5, 8, 8, 8, 8, 8)
Level<-c( 100, 110, 90, 95, 87, 60, 65, 61, 55, 57, 40, 41, 50,
47,
44, 44, 42, 38, 40, 37, 37, 35, 40, 34, 32, 20, 22, 25, 27,
29)
plot(Time,Level,pch=16)
You did not say what sort of "non-linear equation" would best suit, nor
did you offer any background regarding the domain of study. There must
be many ways to do this. After looking at the data, a first pass looks
like this:
> lm(log(Level) ~Time )
Call:
lm(formula = log(Level) ~ Time)
Coefficients:
(Intercept) Time
4.4294 -0.1673
> exp(4.4294)
[1] 83.88107
> points(unique(Time), exp(4.4294 -unique(Time)*0.1673), col="red", pch=4)
Maybe a Weibull model would be more appropriate.
... and to continue David's analysis:
Time <- c(0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 4, 4, 4,
4, 4, 5, 5, 5, 5, 5, 8, 8, 8, 8, 8)
Level <- c(100, 110, 90, 95, 87, 60, 65, 61, 55, 57, 40, 41, 50,
47, 44, 44, 42, 38, 40, 37, 37, 35, 40, 34, 32, 20, 22, 25, 27, 29)
plot(Time, Level, pch = 16)
# First look at log-transformed Level (what David did)
plot(Time, log(Level), pch = 16)
lmod <- lm(log(Level) ~ Time)
summary(lmod)
abline(lmod)
# It is not that clear if log transformation stabilizes variance...
# Given these results, I would tend to use something like this:
plot(Time, Level, pch = 16)
nlmod <- nls(Level ~ exp(a * Time + b) + c,
start = list(a = coef(lmod)[2], b = coef(lmod)[1], c = 0))
summary(nlmod)
Time2 <- seq(0, 8, by = 0.1)
lines(Time2, predict(nlmod, newdata = list(Time = Time2)), col = "red")
You don't tell us enough information about what you study to help
further. Choosing a regression model should also be dependent upon a
priori knowledge about what you study. This looks like an exponential
decay until a non-null value (c) estimated here around 20.6. Do you know
it in advance? On the contrary, are you looking for this value? If yes,
it would be nice, perhaps, to measure 'Level's also for a little bit
longer 'Time's (until 10 or 12 here).
Hope this helps
Best,
Philippe
______________________________________________
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.
