Dear all
I started to use nlxb instead of nls to get rid of singular gradient error.
I try to fit double exponential function to my data, but results I obtain
are strongly dependent on starting values.
tsmes ~ A*exp(a*plast) + B* exp(b*plast)
Changing b from 0.1 to 0.01 gives me completely different results. I usually
check result by a plot but could the result be inspected if it achieved good
result without plotting?
Or is there any way how to perform such task?
Cheers
Petr
Below is working example.
temp <- structure(list(tsmes = c(31, 32, 32, 32, 32, 32, 32, 32, 33,
34, 35, 35, 36, 36, 36, 37, 38, 39, 40, 40, 40, 40, 40, 41, 43,
44, 44, 44, 46, 47, 47, 47, 47, 48, 49, 51, 51, 51, 52, 53, 54,
54, 55, 57, 57, 57, 59, 59, 60, 62, 63, 64, 65, 66, 66, 67, 67,
68, 70, 72, 74, 76, 78, 81, 84, 85, 86, 88, 90, 91, 92, 94, 96,
97, 99, 100, 102, 104, 106, 109, 112, 115, 119, 123, 127, 133,
141, 153, 163, 171), plast = c(50, 51, 52, 52, 53, 53, 53, 54,
55, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 64, 64, 65, 65, 66,
66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 75, 76, 76, 77, 77, 78,
78, 79, 80, 81, 82, 83, 84, 85, 85, 86, 86, 87, 88, 88, 89, 90,
91, 91, 93, 93, 94, 95, 96, 96, 97, 98, 98, 99, 100, 100, 101,
102, 103, 103, 104, 105, 106, 107, 107, 108, 109, 110, 111, 112,
112, 113, 113, 114, 115, 116)), row.names = 2411:2500, class = "data.frame")
library(nlsr)
fit <- nlxb(tsmes ~ A*exp(a*plast) + B* exp(b*plast), data=temp,
start=list(A=1, B=15, a=0.025, b=0.01))
coef(fit)
A B a b
3.945167e-13 9.772749e+00 2.802274e-01 2.179781e-02
plot(temp$plast, temp$tsmes, ylim=c(0,200))
lines(temp$plast, predict(fit, newdata=temp), col="pink", lwd=3)
ccc <- coef(fit)
lines(0:120,ccc[1]*exp(ccc[3]*(0:120)))
lines(0:120,ccc[2]*exp(ccc[4]*(0:120)), lty=3, lwd=2)
# wrong fit with slightly different b
fit <- nlxb(tsmes ~ A*exp(a*plast) + B* exp(b*plast), data=temp,
start=list(A=1, B=15, a=0.025, b=0.1))
coef(fit)
A B a b
2911.6448377 6.8320597 -49.1373979 0.0261391
lines(temp$plast, predict(fit, newdata=temp), col="red", lwd=3)
ccc <- coef(fit)
lines(0:120,ccc[1]*exp(ccc[3]*(0:120)), col="red")
lines(0:120,ccc[2]*exp(ccc[4]*(0:120)), lty=3, lwd=2, col="red")