Hi Troels, Why not call this user function error_analysis.covariance_matrix? By using the specific analysis API, you can easily generalise this. It requires only a few steps, mainly just shifting your existing code around:
1) Rename the user function. 2) Have a backend in the new module pipe_control.error_analysis (I may then merge the pipe_control.monte_carlo module into it for a general module for all error estimate techniques). A good example of how to get the specific analysis API and use it is the Monte Carlo error_analysis() function (http://www.nmr-relax.com/api/3.2/pipe_control.monte_carlo-pysrc.html#error_analysis). 3) Add a new API method to specific_analyses.api_base which, like most of the others, raises a RelaxImplementError (e.g. covariance_matrix()). This allows each analysis type to implement it as desired, and allows the user function to fail cleanly when used on an analysis that doesn't support it. 4) Shift the covariance matrix calculation code of multifit_covar() into the relax library, maybe in the lib.statistics module. Btw, the docstring of this function needs a bit of work to allow the API documentation to be compiled. 5) Create the dispersion specific covariance_matrix() method which obtains the Jacobian and calls the relax library function. This would just be the estimate_r2eff_err() function converted into a method. This is how you add functionality to this API, and it would be a useful learning exercise. If you know how to control this API, it will be very powerful when you are adding new features to relax. In this case, most of your code is unrelated to the dispersion analysis and it is a perfect example of such an independent and powerful feature. It would be a rather simple exercise, changing the new system tests and just shifting the current functions/methods to fit into this infrastructure. Note that the feature of estimating the initial exponential curve parameters by linearisation is also partly analysis independent as it is of benefit for the relaxation curve-fitting analysis type as well. So it could be spun out as a separate user function (maybe exponential_curve.linear_estimate), use the specific API, and putting the specific_analyses/relax_disp/estimate_r2eff.py Exp.estimate_x0_exp() method into the relax library. This too would be a simple code shifting exercise. What do you think? The specific analysis API is one of the most important design concepts in relax, and being able to control and extend it will allow new ideas to be added much more quickly. Regards, Edward On 27 August 2014 20:06, <tlin...@nmr-relax.com> wrote: > Author: tlinnet > Date: Wed Aug 27 20:06:23 2014 > New Revision: 25347 > > URL: http://svn.gna.org/viewcvs/relax?rev=25347&view=rev > Log: > Added front-end to the new user function relax_disp.r2eff_err_estimate(), > which will estimate the R2eff errors > from a pipe and spins with optimised values of R2eff and i0. > > The co-variance matrix can be calculated from the optimised parameters, and > the Jacobian. > > Big care should be taken not to directly trust these results, since the > errors are quite different compared to the Monte-Carlo simulations. > > This implementation, will reach the exact same error estimation as > scipy.optimize.leastsq. > > But with much better control over the data, and insight into the calculations. > > task #7822(https://gna.org/task/index.php?7822): Implement user function to > estimate R2eff and associated errors for exponential curve fitting. > > Modified: > trunk/user_functions/relax_disp.py > > Modified: trunk/user_functions/relax_disp.py > URL: > http://svn.gna.org/viewcvs/relax/trunk/user_functions/relax_disp.py?rev=25347&r1=25346&r2=25347&view=diff > ============================================================================== > --- trunk/user_functions/relax_disp.py (original) > +++ trunk/user_functions/relax_disp.py Wed Aug 27 20:06:23 2014 > @@ -40,6 +40,7 @@ > from pipe_control.mol_res_spin import get_spin_ids > from specific_analyses.relax_disp.catia import catia_execute, catia_input > from specific_analyses.relax_disp.cpmgfit import cpmgfit_execute, > cpmgfit_input > +from specific_analyses.relax_disp.estimate_r2eff import estimate_r2eff_err > from specific_analyses.relax_disp.data import cpmg_setup, insignificance, > plot_disp_curves, plot_exp_curves, r2eff_read, r2eff_read_spin, relax_time, > set_exp_type, r20_from_min_r2eff, spin_lock_field, spin_lock_offset, > write_disp_curves > from specific_analyses.relax_disp.data import INTERPOLATE_DISP, > INTERPOLATE_OFFSET, X_AXIS_DISP, X_AXIS_W_EFF, X_AXIS_THETA, Y_AXIS_R2_R1RHO, > Y_AXIS_R2_EFF > from specific_analyses.relax_disp.nessy import nessy_input > @@ -627,6 +628,50 @@ > uf.gui_icon = "oxygen.status.object-locked" > uf.wizard_size = (800, 500) > uf.wizard_image = ANALYSIS_IMAGE_PATH + 'relax_disp_200x200.png' > + > + > +# The relax_disp.r2eff_err_estimate user function. > +uf = uf_info.add_uf('relax_disp.r2eff_err_estimate') > +uf.title = "Estimate R2eff errors by the Jacobian matrix." > +uf.title_short = "Estimate R2eff errors." > +uf.add_keyarg( > + name = "spin_id", > + py_type = "str", > + arg_type = "spin ID", > + desc_short = "spin ID to restrict value setting to", > + desc = "The spin ID string to restrict value setting to.", > + can_be_none = True > +) > +uf.add_keyarg( > + name = "epsrel", > + py_type = "float", > + default = 0.0, > + desc_short = "parameter to remove linear-dependent columns.", > + desc = "The parameter to remove linear-dependent columns when J is rank > deficient.", > + can_be_none = False > +) > +uf.add_keyarg( > + name = "verbosity", > + default = 1, > + py_type = "int", > + desc_short = "amount of information to print.", > + desc = "The higher the value, the greater the verbosity.", > + can_be_none = False > +) > +# Description. > +uf.desc.append(Desc_container()) > +uf.desc[-1].add_paragraph("This is a new experimental feature from version > 3.3, and should only be tried out with big care.") > +uf.desc[-1].add_paragraph("This will estimate R2eff errors by using the > exponential decay Jacobian matrix 'J' to compute the covariance matrix of the > best-fit parameters.") > +uf.desc[-1].add_paragraph("This can be an huge time saving step, when > performing model fitting in R1rho. Errors of R2eff values, are normally > estimated by time-consuming Monte-Carlo simulations.") > +uf.desc[-1].add_paragraph("This method is inspired from the GNU Scientific > Library (GSL).") > +uf.desc[-1].add_paragraph("The covariance matrix is given by: covar = Qxx = > (J^T J)^{-1}.") > +uf.desc[-1].add_paragraph("Qxx is computed by QR decomposition, Qxx=QR, > Qxx^-1=R^-1 Q^T. The columns of R which satisfy: |R_{kk}| <= epsrel |R_{11}| > are considered linearly-dependent and are excluded from the covariance matrix > (the corresponding rows and columns of the covariance matrix are set to > zero).") > +uf.desc[-1].add_paragraph("The parameter 'epsrel' is used to remove > linear-dependent columns when J is rank deficient.") > +uf.backend = estimate_r2eff_err > +uf.menu_text = "&r2eff_err_estimate" > +uf.gui_icon = "relax.relax_fit" > +uf.wizard_size = (800, 800) > +uf.wizard_image = ANALYSIS_IMAGE_PATH + sep + 'blank_150x150.png' > > > # The relax_disp.r2eff_read user function. > > > _______________________________________________ > relax (http://www.nmr-relax.com) > > This is the relax-commits mailing list > relax-comm...@gna.org > > 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@gna.org 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
