As John said, sums of exponentials are hard. One thing that often helps a lot is to use the partially linear structure: given a and b, you've got a linear model to compute A and B. Now that you're down to two nonlinear parameters, you can draw a contour plot of nearby values to see how much of a mess you're dealing with.

Duncan Murdoch

On 07/05/2020 9:12 a.m., PIKAL Petr wrote:
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.

dput(temp)
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")


______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.


______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
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.

Reply via email to