Pedro Mardones <mardones.p <at> gmail.com> writes:
>
> Dear R-users;
>
> I'm working with a a dataset that was previously used to fit a
> nonlinear model of the form:
>
> Y ~ a * (1 + b * log(1 - c * X^d))
>
> The parameters published elsewhere are:
>
> a = 1.758863, b = .217217, c = .99031, and d = .054589
>
> However, there is no way I can replicate this result. I've tried
> several options (including SAS) w/o success.
>
This is going to be a very difficult model to fit.
## plot data and putative (old) fit
oldparms <- list(a = 1.758863, b = .217217, c = .99031, d = .054589)
plot(Y~X,data=test)
with(oldparms,curve(a * (1 + b * log(1 - c * x^d)),add=TRUE, col=2))
## compute RSS of old fit:
oldfit <- with(c(oldparms,as.list(test)),
a * (1 + b * log(1 - c * X^d)))
oldresid <- test$Y-oldfit
sum(oldresid^2) ## 0.263
I used your nlsLM() code with a slight modification:
library(minpack.lm)
fit2 <- nlsLM(Y ~ a * (1 + b * log(1 - c * X^d)), data = test,
start = list(a = 1.75, b = .22, c = .99, d = .05),
control = nls.lm.control(maxfev = 1000,maxiter=1000))
sum(residuals(fit2)^2) ## 0.241
with(as.list(coef(fit2)),
curve(a * (1 + b * log(1 - c * x^d)),add=TRUE, col=4))
## calculate correlation among parameters
cov2cor(vcov(fit2))
a b c d
a 1.0000000 -0.9999982 0.9999960 -0.9999991
b -0.9999982 1.0000000 -0.9999995 0.9999999
c 0.9999960 -0.9999995 1.0000000 -0.9999988
d -0.9999991 0.9999999 -0.9999988 1.0000000
My conclusions:
* nlsLM() finds a slightly better fit to the data than
the quoted parameters;
* the parameters are all _extremely_ strongly
correlated with each other. With a little bit more
algebra/calculus you could probably figure out why, by
Taylor expanding etc.: c is very close to 1, d is very close to zero ...
______________________________________________
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.
