Why is there such a big difference between the fitting results of 'exp2' and the same form of custom function?

10 views (last 30 days)
Here I get the 'test.txt' data, and want to have a good fitting result.
I found that using
[exp_lm,gof_lm] = fit(x,y,'exp2',Algorithm='Levenberg-Marquardt')
to fit the data will get a obviously better result than
[exp_lm,gof_lm] = fit(x,y,ft,Algorithm='Levenberg-Marquardt')
where ft was been set as
ft = a*exp(b*x) + c*exp(d*x)
which is the same val of 'exp2'.
May I ask what is the reason for this, and what do I need to do to make the custom function obtain the same fitting result?
and
  1 Comment
Walter Roberson
Walter Roberson on 31 Mar 2023
As best I can tell, some of the fit types must be custom-implemented using other algorithms rather than some kind of least-squared fit. However, I do not presently know which fit types or what the actual algorithms are.

Sign in to comment.

Accepted Answer

dpb
dpb on 1 Apr 2023
Edited: dpb on 1 Apr 2023
The biggest difference is probably found in the documentation under the 'Nonlinear Least-Squares Options" parameter 'StartPoint' which says
"If no start points (the default value of an empty vector) are passed to the fit function, starting points for some library models are determined heuristically. For rational and Weibull models, and all custom nonlinear models, the toolbox selects default initial values for coefficients uniformly at random from the interval (0,1). As a result, multiple fits using the same data and model might lead to different fitted coefficients. To avoid this, specify initial values for coefficients with a fitoptionsobject or a vector value for the StartPoint value."
A random value between [0 1] for the exponents for the above dataset is about as poor a choice as one could make; in fact, it's clear that the fit would probably be pretty good as Y0*[1-model] instead. Using a heuristic approach to roughly estimate some starting points would undoubtedly do much better; in particular ensuring the exponents are <0 to start with.
>> ft = fittype('a*exp(b*x)+c*exp(d*x)');
>> A0=[max(xy(:,2))-min(xy(:,2))]/exp(-xy(1,1))
A0 =
810007.04
>> [f,~,o]=fit(xy(:,1),xy(:,2),ft,'Algorithm','Levenberg-Marquardt','StartPoint',[A0 -1 A0/10 -1/10])
f =
General model:
f(x) = a*exp(b*x)+c*exp(d*x)
Coefficients (with 95% confidence bounds):
a = 1.11e+06 (1.061e+06, 1.159e+06)
b = -1.307 (-1.34, -1.274)
c = 3.036e+04 (2.799e+04, 3.272e+04)
d = -0.1388 (-0.1476, -0.1301)
o =
struct with fields:
numobs: 49.00
numparam: 4.00
residuals: [49×1 double]
Jacobian: [49×4 double]
exitflag: 3.00
iterations: 8.00
funcCount: 45.00
stepsize: 49.95
cgiterations: []
firstorderopt: 35569.25
algorithm: 'levenberg-marquardt'
message: 'Success, but fitting stopped because change in residuals less than tolerance (TolFun).'
>> plot(xy(:,1),xy(:,2),'x')
>> hold on
>> plot(f)
>>
does in fact, almost identically reproduce the builtin. The amplitude is based on the difference between Ymax and Yinf; the two exponentials are set to be <0; I didn't try to get fancy and figure out heuristics for the second exponential but just took a rough clue from the original solution; if you didn't have that, it might take some more work to get a good starting point; I also haven't taken the time to see how far off can be and still succeed. I did run into some cases where the internal algorithm crashed with an intermediary evaluation that overflowed so have to take some care to avoid that issue or revise guesses if run into the problem. You can bound coefficients, but not with the Marquardt algorithm.
The above produces
Seems like I recall another very similar Q? a while back, @Walter Roberson, but didn't find it in a quick search.
  3 Comments
dpb
dpb on 1 Apr 2023
Edited: dpb on 2 Apr 2023
No problem; glad to help; I remembered this peculiarity from some time previous. It seems a very poor way to have designed the toolbox; one has some sympathy for the designers/implementors for the custom model in that it's a tough nut to crack when the user could write anything to have some way to make a more informed guess without some additional user help/input.
Given the importance of the feature, however, it would seem perhaps the note should be at the top level of the documentation instead of only buried down where it is...
I didn't try to see whether for this case only ensuring the exponents were negative might be enough to keep the convergence from losing its way or not...that's an interesting question.
>> [f,~,o]=fit(xy(:,1),xy(:,2),ft,'Algorithm','Levenberg-Marquardt','StartPoint',[rand -1 rand -1])
Error using fit>iFit (line 348)
Inf computed by model function, fitting cannot continue.
Try using or tightening upper and lower bounds on coefficients.
Error in fit (line 116)
[fitobj, goodness, output, convmsg] = iFit( xdatain, ydatain, fittypeobj, ...
>>
Well, not so good, indeed...the first attempt does the thing of overflowing when tries to calculate the model and quits. So needs better starting point than just that hint.
John D'Errico
John D'Errico on 1 Apr 2023
Exactly so. If I knew the model as a sum of two exponential terms, I could arguably choose a scheme to find starting values that would be intelligently chosen, at least some of the time. Even then, I can give you data where no starting values in a finite number of digits will be viable.
But a black box model, unless I was willing to try to write code to parse the model, there is no way to choose intelligent starting values. And then you get complete crap at least some of the time. The burden of course is to force the user to provide intelligent starting values. They would/should learn that quickly. But that is a literal truism for any nonlinear model - crappy starting values will often result in crap for a fit.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!