Re: r25228 - in /trunk/test_suite/shared_data/curve_fitting/profiling: profiling_relax_fit.py relax_fit.py -- August 25, 2014

If you would like such a once off solution rather than a general
solution - for example the general error_analysis.covariance_matrix
solution may have been useful for the errors for the dispersion models
themselves for a quick estimate - then I would recommend using the
same code but using different code paths to isolate it.  I would
suggest putting everything into a
specific_analyses.relax_disp.estimate_r2eff module.  That includes the
backend specific_analyses.relax_disp.optimisation.estimate_r2eff()
function and the code in the target_functions package.  The code in
target_functions.relax_disp_curve_fit is a lot more than just a target
function, so it doesn't really belong in this package.

In the new module, you would first import the dep_check module.  Then
if dep_check.scipy_model == True, import the scipy function.  At the
start of the estimate_r2eff() function, you should add new RelaxError
if dep_check.scipy_model == False.  This allows relax to continue
operating without Scipy being installed.  You should also make sure
that any tests in the test suite can operate without Scipy being
present.  This is important as I know that the majority of relax users
do not have Scipy installed.

Regards,

Edward




On 25 August 2014 11:55, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:

Hi Edward.

You have already convinced me that there are other options.

What I try to do here, is the give the 'availability' to the user, to
try different methods.

Here one method is: scipy.optimize.leastsq

The user function detailed describes, that this is by 
scipy.optimize.leastsq.

Anything which comes from that, can be tracked back to this function.

It is the same arguments, that one cannot change the published "model
function" in dispersion, but they stay as they are published.

Best
Troels


2014-08-25 11:43 GMT+02:00 Edward d'Auvergne <edward@xxxxxxxxxxxxx>:

Hi,

The scipy.optimize.leastsq() function contained a fatal bug back in
2003.  I reported it and it was not fixed the last time I checked,
probably around 2008.  I wrote minfx to replace scipy.optimize due to
the fatal bugs present in all 3 optimisation algorithms available back
in 2003 (Levenberg-Marquardt, Nelder-Mead simplex, and one other).
Using MINPACK means that they have changed the code for the LM
algorithm, but due to the disregard the developers had for fixing
fatal bugs back in the early days, there is not much someone can say
to convince me that the current code is also not problematic.  Showing
that scipy.optimize.leastsq() finds the correct result does not
disprove this.  The fatal failures occur rarely and hence normally you
don't see a problem.  But when the failure occurs, in case you are
carefully checking you will never know.  This is standard for
optimisation bugs, for example the fatal bug in Modelfree4 that I
fixed and was released as version 4.20
(http://www.palmer.hs.columbia.edu/software/modelfree.html).  This bug
was not seen in 300+ publications from 1991 to 2006 (over 15 years),
despite it affecting all these published results!

Minfx is far more extensive and complete when compared to
scipy.optimize.  And it is tested in the relax test suite.  The only
feature scipy.optimize has that minfx does not is that it calls
scipy.integrate.quadpack for numerical gradient estimates.  But
deriving these yourself is basic.  For example using Quickmath, open
http://www.quickmath.com/webMathematica3/quickmath/calculus/differentiate/advanced.jsp.
Set the equation to:

I*exp(-R*t)

Set the variables to I first, then the differentiate button will tell
you the df/dI partial derivative.  Set the variable to R,
differentiate, then the result is the df/dR partial derivative.  The
vector [df/dI, df/dR] is the gradient.  The Hessian is simple, just
set the variables to [I, I], then [I, R], then [R, I], then [R, R].
Put the result into a matrix and you have the Hessian.  I can't
emphasis enough how basic these equations and solution are.

Anyway, the optimisation solution in relax is minfx, not
scipy.optimize.  All algorithms in scipy are available in minfx.  And
the scipy ones have never been faster than the minfx ones.  The
MINPACK LM algorithm might be slightly faster than the minfx LM
algorithm, but this is not the reason for the speed ups you see!

Regards,

Edward

On 25 August 2014 11:06, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> 
wrote:

Hi Edward.

This post continues from:
http://thread.gmane.org/gmane.science.nmr.relax.devel/6814

As I wrote in this post, the new user function:
"relax_disp.r2eff_estimate()"

gives freedom to the user to try other algorithms.
Here scipy.optimize.leastsq which should be (or is) a well known
option in the python community for minimisation.

Did you find errors in MINPACK or in scipy ?

Best
Troels


2014-08-25 10:13 GMT+02:00 Edward d'Auvergne <edward@xxxxxxxxxxxxx>:

Hi Troels,

You should note that the 20x speed up is due to the change in
optimisation algorithm rather than scipy vs. minfx vs. MINPACK.  The
reason I wrote minfx is that I started with scipy optimisation but
found all implemented algorithms contained fatal bugs.  This was not
fixed for years after I reported it, and I don't know if the code had
changed since 2003 when I looked into it.  Anyway a scipy optimisation
solution is incompatible with the minfx optimisation solution in
relax.  If you derive and code the gradients into relax, then you can
use the minfx LM algorithm as a solution.

Regards,

Edward







On 25 August 2014 01:08,  <tlinnet@xxxxxxxxxxxxx> wrote:

Author: tlinnet
Date: Mon Aug 25 01:08:44 2014
New Revision: 25228

URL: http://svn.gna.org/viewcvs/relax?rev=25228&view=rev
Log:
Further improved the profiling of relax curve fit.

This profiling shows, that Python code is about twice as slow as the C 
code implemented.

But it also shows that optimising with scipy.optimize.leastsq is 20 X 
faster.
It also gives reasonable error values.

Combining a function for a linear fit to guess the initial values, 
together
with scipy optimise, will be an extreme time win for estimating R2eff 
values fast.

A further test would be to use relax Monte-Carlo simulations for say 
1000-2000 iterations,
and compare to the errors extracted from estimated covariance.

Added:
    trunk/test_suite/shared_data/curve_fitting/profiling/relax_fit.py
Modified:
    
trunk/test_suite/shared_data/curve_fitting/profiling/profiling_relax_fit.py

[This mail would be too long, it was shortened to contain the URLs 
only.]

Modified: 
trunk/test_suite/shared_data/curve_fitting/profiling/profiling_relax_fit.py
URL: 
http://svn.gna.org/viewcvs/relax/trunk/test_suite/shared_data/curve_fitting/profiling/profiling_relax_fit.py?rev=25228&r1=25227&r2=25228&view=diff

Added: trunk/test_suite/shared_data/curve_fitting/profiling/relax_fit.py
URL: 
http://svn.gna.org/viewcvs/relax/trunk/test_suite/shared_data/curve_fitting/profiling/relax_fit.py?rev=25228&view=auto


_______________________________________________
relax (http://www.nmr-relax.com)

This is the relax-commits mailing list
relax-commits@xxxxxxx

To unsubscribe from this list, get a password
reminder, or change your subscription options,
visit the list information page at
https://mail.gna.org/listinfo/relax-commits


_______________________________________________
relax (http://www.nmr-relax.com)

This is the relax-devel mailing list
relax-devel@xxxxxxx

To unsubscribe from this list, get a password
reminder, or change your subscription options,
visit the list information page at
https://mail.gna.org/listinfo/relax-devel

Re: r25228 - in /trunk/test_suite/shared_data/curve_fitting/profiling: profiling_relax_fit.py relax_fit.py

Header

Content

Related Messages