The StartPoint of Curve Fitting of the Function f(x)=a*exp(-x/b)+c

13 views (last 30 days)
I am trying to fit the curve f(x)=a*exp(-x/b)+c , but the startpoints influences the results, besides, the result of b is the same with the initial value, how can I solve this?
x=[10,30,60,120,180,300,420,600,900,1200,1500,1800,2100,2400,2700,3000,3600];
y=[7.23E-08,6.83E-08,6.58E-08,5.86E-08,5.52E-08,4.54E-08,3.97E-08,3.43E-08,2.77E-08,2.17E-08,
1.86E-08,1.44E-08,1.19E-08,1.02E-08,9.4E-09,8.7E-09,6.2E-09];
f1=fittype('a*exp(-x/b)+c','independent','x','coefficients',{'a','b','c'});
opts=fitoptions('Method','NonlinearLeastSquares');
opts.StartPoint=[1e-8,1000,1e-8];%if I set the value as 1 1 1 the result is poor%
opts.Display='Off';
[result,gof]=fit(x,y,f1,opts);

Accepted Answer

Walter Roberson
Walter Roberson on 15 Feb 2022
If we try to solve symbolically, and use calculus
x=[10,30,60,120,180,300,420,600,900,1200,1500,1800,2100,2400,2700,3000,3600];
y=[7.23E-08,6.83E-08,6.58E-08,5.86E-08,5.52E-08,4.54E-08,3.97E-08,3.43E-08,2.77E-08,2.17E-08,1.86E-08,1.44E-08,1.19E-08,1.02E-08,9.4E-09,8.7E-09,6.2E-09];
syms a b c
model_y = a*exp(-x/b)+c
r1 = vpa(sum((y - model_y).^2))
dr1 = diff(r1, a)
partiala = solve(dr1, a)
r2 = subs(r1, a, partiala)
dr2 = diff(r2, c)
partialc = solve(dr2, c)
r3 = subs(r2, c, partialc); %very long expression!
dr3 = diff(r3, b);
partialb = vpasolve(dr3, b);
bestb = partialb
bestc = vpa(subs(partialc, b, bestb))
besta = vpa(subs(subs(partiala, c, bestc),b,bestb))
but in practice, the attempt to solve() dr2 for c takes a very long time.
But r1 is a function of two variables, so you could try to plot it and look for minima. Or you could plot dr2 and look for minima.
When I explore, I find that there is a range of b values near 0 for which the values are imaginary, but that the functions are well defined if you take b a bit more negative or a bit more positive. And when I explore dr2 (the derivative with respect to c), I see that it is asymptopic towards 0 on both sides away from b = 0.
What that implies is that the system does not have a finite global minima. That any numeric exploration of the minima is going to stop when it hits the limit of numeric accuracy.
With the derivative going to 0 at the infinite limit of b, then examining the model we can see that the implication is that you can get a better fit from using the model that a = 0, b (not 0, not -infinity, could be +infinity, but finite non-zero is okay), c = mean(y) . Or a = finite, b = +infinity, c = mean(y) . To phrase that another way, dropping the term a*exp(-x/b) appears to give you a better fit. Not a good fit, but better than you can get with any finite non-zero b.

More Answers (1)

Alan Weiss
Alan Weiss on 15 Feb 2022
I think that your data is very small, about order 1e-8 instead of order 1. Here is a little script I knocked together that seems to work reasonably:
xdata = [10,30,60,120,180,300,420,600,900,1200,1500,1800,2100,2400,2700,3000,3600];
ydata = [7.23E-08,6.83E-08,6.58E-08,5.86E-08,5.52E-08,4.54E-08,3.97E-08,3.43E-08,2.77E-08,2.17E-08,...
1.86E-08,1.44E-08,1.19E-08,1.02E-08,9.4E-09,8.7E-09,6.2E-09];
zdata = ydata*1e8; % rescale
fun = @(x)x(1)*exp(-xdata/x(2)) + x(3) - zdata;
sol = lsqnonlin(fun,[1,1000,1]) % Tries to make fun = 0
Local minimum found. Optimization completed because the size of the gradient is less than the value of the optimality tolerance.
sol = 1×3
6.1506 698.6373 0.8581
plot(xdata,zdata,'ro',xdata,fun(sol)+zdata,'kx')
Alan Weiss
MATLAB mathematical toolbox documentation

Products


Release

R2020b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by