Menu
▾
▴
#33 Improved algorithm for non-linear fitting
I think that interactions aren't being made correctly or they aren't working.
The attached file contains a data set that I tried to fit using the formula a/(x^5*(exp(b/x)-1)). The initial values used were: a = 2e+24 and b = 45941 (very close to the expected values). The result obtained were:
a = 2e+24 +/- 0
b = 45936.6863 +/- 32.437943
Iterations = 0
There is no problem with the b coefficient. But the a coefficient wasn't changed and the interactions number were zero.
This problem occurs with Linux and Windows versions.
In addition, this doesn't occur with 0.1.3 or oldest versions.
You've got a somewhat problematic example there - the derivative of the formula with respect to a is very small (roughly 1e-23 to 1e-25). I think this explains why fitting tends to leave "a" untouched - changing it simply has a very small effect on the residuals.
Due to the above, what exactly happens depends sensitively on the way the derivative is computed numerically. I've even produced different results with SciDAVis 0.1.3 running on different systems. Also, LabPlot simply fails to produce a fit with the initial values you gave.
I've verified the fitting code against NIST StRD[1] data, so I'm pretty certain there's nothing fundamental wrong with it.
One way to improve fitting accuracy is to compute derivatives analytically and write a fit plugin. While there are probably ways to improve our current approach to numerical differentiation, I don't think a generic, easy-to-use solution will ever provide the same accuracy as code written specifically for the problem at hand.
[1] http://www.itl.nist.gov/div898/strd/index.html
Following up on a thread on scidavis-contributors[1], I'm reopening this as a feature request. While I do not agree that it's an actual bug in SciDAVis (non-linear fitting cannot be guaranteed to find the global chi^2 minimum, nor can any convergence criterion be guaranteed to make sense in every case), it seems that Grace[2] finds a good solution for a greater range of initial parameters. It would be desirable to improve non-linear fitting in SciDAVis by supporting the same algorithm as Grace, and possibly by giving the user more choice for convergence criteria.
Note that with SciDAVis 0.1.4 and beyond, numerical differentation uses a routine from GSL instead of muParser's built-in diffirentiation. This was necessary for the fixed version of constant parameters, and of course it's necessary if we ever want to support using Python for fit formulas. I think this change explains why older versions of SciDAVis give somewhat different results. Still, SciDAVis 0.1.3 also does not work reliably for all sensible starting values, so it's not just a problem with the above-mentioned bugfix.
Studying the source code of Grace, I've found that it uses an automatic C translation of the Fortran subroutine LMDIF from MINPACK[3]. This routine implements a modified Levenberg-Marquardt algorithm (which is broadly the same as the one we're currently using from GSL, although the "modification" might make a difference). It calls the subroutine FDJAC2 (also from MINPACK) for numerical differentiation. I guess the simplest solution would be to add an option for fitting using exactly those routines, as in Grace.
[1] http://sourceforge.net/mailarchive/forum.php?thread_name=973128.62836.qm@web65616.mail.ac4.yahoo.com&forum_name=scidavis-contributors
[2] http://plasma-gate.weizmann.ac.il/Grace/
[3] http://www.netlib.org/minpack/