We should then definitely have this tested in the test suite! I just had a look at it and will commit the fix very soon. The problem is that I have recently introduced the spin-lock (or hard pulse for the CPMG) offset dimension into all of the dispersion data structures used in the target functions. But the grid_search_setup() method was not updated correctly for this. Note the loop_exp_frq_point() function call, the fix is to replace it with loop_exp_frq_offset_point() to loop over all of these! Cheers, Edward On 10 December 2013 11:43, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:
Hi Edward. I can say that running this: relax ./test_suite/shared_data/dispersion/r1rho_on_res_m61/r2eff_calc.py Give the same error. File "/Users/tlinnet/software/relax_trunk/pipe_control/minimise.py", line 152, in grid_search grid_search(lower=lower, upper=upper, inc=inc, constraints=constraints, verbosity=verbosity) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 1093, in grid_search self.minimise(min_algor='grid', lower=lower, upper=upper, inc=inc, constraints=constraints, verbosity=verbosity, sim_index=sim_index) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 1142, in minimise self._minimise_r2eff(min_algor=min_algor, min_options=min_options, func_tol=func_tol, grad_tol=grad_tol, max_iterations=max_iterations, constraints=constraints, scaling=scaling, verbosity=verbosity, sim_index=sim_index, lower=lower, upper=upper, inc=inc) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 357, in _minimise_r2eff grid_size, inc_new, lower_new, upper_new = grid_search_setup(spins=[spin], spin_ids=[spin_id], param_vector=param_vector, lower=lower, upper=upper, inc=inc, scaling_matrix=scaling_matrix) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/optimisation.py", line 282, in grid_search_setup lower_new.append(lower[i] / scaling_matrix[i, i]) IndexError: list index out of range 2013/12/10 Edward d'Auvergne <edward@xxxxxxxxxxxxx>Hi, I saw that, but I'm not sure why it is saying that. Maybe you should add some debugging printouts in the grid_search_setup() function and see what is happening? It won't be because of the missing R1 data. I would suggest that if you would like to develop the code for optimising R1 as a parameter of the R1rho models (i.e. all of them), that we should create a special branch for that. I know exactly what needs to be done and I can direct you in the right direction. But it needs to be developed in a branch as it will be quite disruptive. Regards, Edward On 10 December 2013 09:56, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:Hi Edward. I now receive this error: relax -s Relax_disp.test_r1rho_kjaergaard But I guess this is expected, since I do not load any R1 data in. ? Best Troels --------- relax> relax_disp.select_model(model='R2eff') R2eff/R1rho value and error determination. Optimisation ============ Nesting and model equivalence checks ------------------------------------ No model nesting or model equivalence detected. relax> grid_search(lower=None, upper=None, inc=4, constraints=True, verbosity=1) Traceback (most recent call last): File "/Users/tlinnet/software/relax_trunk/test_suite/system_tests/relax_disp.py", line 2591, in test_r1rho_kjaergaard relax_disp.Relax_disp(pipe_name=pipe_name, pipe_bundle=pipe_bundle, results_dir=ds.tmpdir, models=MODELS, grid_inc=GRID_INC, mc_sim_num=MC_NUM, modsel=MODSEL) File "/Users/tlinnet/software/relax_trunk/auto_analyses/relax_disp.py", line 116, in __init__ self.run() File "/Users/tlinnet/software/relax_trunk/auto_analyses/relax_disp.py", line 451, in run self.optimise(model=model) File "/Users/tlinnet/software/relax_trunk/auto_analyses/relax_disp.py", line 350, in optimise self.interpreter.grid_search(inc=self.grid_inc) File "/Users/tlinnet/software/relax_trunk/prompt/uf_objects.py", line 221, in __call__ self._backend(*new_args, **uf_kargs) File "/Users/tlinnet/software/relax_trunk/pipe_control/minimise.py", line 152, in grid_search grid_search(lower=lower, upper=upper, inc=inc, constraints=constraints, verbosity=verbosity) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 1093, in grid_search self.minimise(min_algor='grid', lower=lower, upper=upper, inc=inc, constraints=constraints, verbosity=verbosity, sim_index=sim_index) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 1142, in minimise self._minimise_r2eff(min_algor=min_algor, min_options=min_options, func_tol=func_tol, grad_tol=grad_tol, max_iterations=max_iterations, constraints=constraints, scaling=scaling, verbosity=verbosity, sim_index=sim_index, lower=lower, upper=upper, inc=inc) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/api.py", line 357, in _minimise_r2eff grid_size, inc_new, lower_new, upper_new = grid_search_setup(spins=[spin], spin_ids=[spin_id], param_vector=param_vector, lower=lower, upper=upper, inc=inc, scaling_matrix=scaling_matrix) File "/Users/tlinnet/software/relax_trunk/specific_analyses/relax_disp/optimisation.py", line 282, in grid_search_setup lower_new.append(lower[i] / scaling_matrix[i, i]) IndexError: list index out of range 2013/12/10 Edward d'Auvergne <edward@xxxxxxxxxxxxx>Ah, ok, that would explain the strangeness. Then the simple solution would be to check that the correct data pipe type is being used at the start of both these user functions! Then a user should never reach this failure point. Regards, Edward On 10 December 2013 09:34, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:Uf... There was the error: # Set the relaxation times (in s). self.interpreter.relax_fit.relax_time(time=time_sl, spectrum_id=sp_id) was used instead of: # Set the relaxation times (in s). self.interpreter.relax_disp.relax_time(spectrum_id=sp_id, time=time_sl) 2013/12/10 Edward d'Auvergne <edward@xxxxxxxxxxxxx>Hi, That's exactly the point of failure. Maybe it would be useful to create one or two simple unit tests for the get_curve_type() function, using some of the data in the current data pipe where your script is failing. The problem is that for the experiment ID being sent into it, it can only find one matching relaxation time. Are you sure that all data is loaded correctly? Regards, Edward On 10 December 2013 09:19, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:Hi Edward. The function: has_exponential_exp_type() is returning False, as the get_curve_type(id) is returning: fixed time So that test look good enough for me? Is it the: fixed time which is a problem? Best Troels 2013/12/10 Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx>Hi Edward. I will try to look into it today. Best Troels 2013/12/9 Edward d'Auvergne <edward@xxxxxxxxxxxxx>Hi Troels, Note that I will soon release relax 3.1.1 with all of the recent fixes (possibly tomorrow). This includes the 'NS MMQ 3-site', 'NS MMQ 3-site linear', 'NS R1rho 3-site', and 'NS R1rho 3-site linear' dispersion models which are now implemented, tested, and find similar results to Dmitry Korzhnev's cpmg_fit software. If the test is not functional by then, I will disable it for the release. It would be very useful to have fixed though as it looks like variable relaxation time data (SQ CPMG, MMQ CPMG, or R1rho) will all currently fail in the auto-analysis. I've had a look and have another hint for you - the problem is that has_exponential_exp_type() function is not returning the correct answer ;) Regards, Edward On 9 December 2013 18:12, Edward d'Auvergne <edward@xxxxxxxxxxxxx> wrote:Hi, This is a strange failure! The auto-analysis should not be running the calc user function for this. That is the problem, you cannot run the calc user function for non-constant relaxation time experiments. The error message is the standard one to tell the user that. The grid_search and minimise user functions should be used instead (maybe the error message can be modified to clarify this and include this info). But it is the auto-analysis that is running this. So the problem is there. relax can of course handle variable relaxation times for any dispersion data type. Regards, Edward On 9 December 2013 18:02, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:Hi Edward. The R1rho data is not constant time. I have only included the models: MODELS = ['R2eff', 'No Rex', 'DPL94'] I have set it up, to find a solution for analysing R1rho data, where R1 data has not been acquired, but for different It actually also fails at the moment, and will probably do for some time. -------- relax -s Relax_disp.test_r1rho_kjaergaar relax> calc(verbosity=1) Traceback (most recent call last): File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/test_suite/system_tests/relax_disp.py", line 2581, in test_r1rho_kjaergaard relax_disp.Relax_disp(pipe_name=pipe_name, pipe_bundle=pipe_bundle, results_dir=ds.tmpdir, models=MODELS, grid_inc=GRID_INC, mc_sim_num=MC_NUM, modsel=MODSEL) File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/auto_analyses/relax_disp.py", line 116, in __init__ self.run() File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/auto_analyses/relax_disp.py", line 447, in run self.interpreter.calc() File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/prompt/uf_objects.py", line 221, in __call__ self._backend(*new_args, **uf_kargs) File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/pipe_control/minimise.py", line 86, in calc calculate(verbosity=verbosity) File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/specific_analyses/relax_disp/api.py", line 717, in calculate self._calculate_r2eff() File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/specific_analyses/relax_disp/api.py", line 182, in _calculate_r2eff check_exp_type_fixed_time() File "/sbinlab2/tlinnet/software/NMR-relax/relax_trunk/specific_analyses/relax_disp/checks.py", line 112, in check_exp_type_fixed_time raise RelaxError("The experiment '%s' is not of the fixed relaxation time period data type." % exp_type) RelaxError: RelaxError: The experiment 'R1rho' is not of the fixed relaxation time period data type. -------------- Is the R1rho analysis only implemented for fixed time periods? Best Troels 2013/12/9 Edward d'Auvergne <edward@xxxxxxxxxxxxx>Hi Troels, When looking at this data and analysis, remember that I have not implemented Dmitry Korzhnev's "correction" for constant time R1rho data. I don't know if that was used in the original publication for your data. More details are given in the 'To do' section of the manual (I only recently added this info). I also don't know what the rest of the field think of his correction and how it applies to later models from the Palmer group. Regards, Edward On 9 December 2013 17:49, <tlinnet@xxxxxxxxxxxxx> wrote:Author: tlinnet Date: Mon Dec 9 17:49:49 2013 New Revision: 21920 URL: http://svn.gna.org/viewcvs/relax?rev=21920&view=rev Log: Added system test for the analysis of optimisation of the Kjaergaard et al., 2013 Off-resonance R1rho relaxation dispersion experiments using the 'DPL' model. Work in progress for Support Request #3083, (https://gna.org/support/index.php?3083) - Addition of Data-set for R1rho analysis. Modified: trunk/test_suite/system_tests/relax_disp.py Modified: trunk/test_suite/system_tests/relax_disp.py URL: http://svn.gna.org/viewcvs/relax/trunk/test_suite/system_tests/relax_disp.py?rev=21920&r1=21919&r2=21920&view=diff ============================================================================== --- trunk/test_suite/system_tests/relax_disp.py (original) +++ trunk/test_suite/system_tests/relax_disp.py Mon Dec 9 17:49:49 2013 @@ -2450,6 +2450,137 @@ self.assertAlmostEqual(cdp.mol[0].res[0].spin[0].chi2, 0.030959849811015544, 3) + def test_r1rho_kjaergaard(self): + """Optimisation of the Kjaergaard et al., 2013 Off-resonance R1rho relaxation dispersion experiments using the 'DPL' model. + + This uses the data from Kjaergaard's paper at U{DOI: 10.1021/bi4001062<http://dx.doi.org/10.1021/bi4001062>}. + + """ + + # The path to the data files. + data_path = status.install_path + sep+'test_suite'+sep+'shared_data'+sep+'dispersion'+sep+'Kjaergaard_et_al_2013' + + # Set pipe name, bundle and type. + pipe_name = 'base pipe' + pipe_bundle = 'relax_disp' + pipe_type= 'relax_disp' + + # Create the data pipe. + self.interpreter.pipe.create(pipe_name=pipe_name, bundle=pipe_bundle, pipe_type=pipe_type) + + # Read the spins. + self.interpreter.spectrum.read_spins(file='1_0_46_0_max_standard.ser', dir=data_path+sep+'peak_lists') + + # Test some of the sequence. + self.assertEqual(len(cdp.mol), 1) + self.assertEqual(cdp.mol[0].name, None) + self.assertEqual(len(cdp.mol[0].res), 48) + + # Name the isotope for field strength scaling. + self.interpreter.spin.isotope(isotope='15N') + + # Set number of experiments to be used. + NR_exp = -1 + + # Load the experiments settings file. + expfile = open(data_path+sep+'exp_parameters_sort.txt','r') + expfileslines = expfile.readlines()[:NR_exp] + expfile.close() + + # In MHz + yOBS = 81.050 + # In ppm + yCAR = 118.078 + centerPPM_N15 = yCAR + + ## Read the chemical shift data. + self.interpreter.chemical_shift.read(file='1_0_46_0_max_standard.ser', dir=data_path+sep+'peak_lists') + + # Test the chemical shift data. + cs = [122.223, 122.162, 114.250, 125.852, 118.626, 117.449, 119.999, 122.610, 118.602, 118.291, 115.393, + 121.288, 117.448, 116.378, 116.316, 117.263, 122.211, 118.748, 118.103, 119.421, 119.317, 119.386, 117.279, + 122.103, 120.038, 116.698, 111.811, 118.639, 118.285, 121.318, 117.770, 119.948, 119.759, 118.314, 118.160, + 121.442, 118.714, 113.080, 125.706, 119.183, 120.966, 122.361, 126.675, 117.069, 120.875, 109.372, 119.811, 126.048] + + i = 0 + for spin, spin_id in spin_loop(return_id=True): + print spin.name, spin.num, spin_id, spin.chemical_shift, cs[i] + # Check the chemical shift. + self.assertEqual(spin.chemical_shift, cs[i]) + + # Increment the index. + i += 1 + + # The lock power to field, has been found in an calibration experiment. + spin_lock_field_strengths_Hz = {'35': 431.0, '39': 651.2, '41': 800.5, '43': 984.0, '46': 1341.11, '48': 1648.5} + + # Apply spectra settings. + for i in range(len(expfileslines)): + line=expfileslines[i] + if line[0] == "#": + continue + else: + # DIRN I deltadof2 dpwr2slock ncyc trim ss sfrq + DIRN = line.split()[0] + I = int(line.split()[1]) + deltadof2 = line.split()[2] + dpwr2slock = line.split()[3] + ncyc = int(line.split()[4]) + trim = float(line.split()[5]) + ss = int(line.split()[6]) + set_sfrq = float(line.split()[7]) + apod_rmsd = float(line.split()[8]) + spin_lock_field_strength = spin_lock_field_strengths_Hz[dpwr2slock] + + # Calculate spin_lock time + time_sl = 2*ncyc*trim + + # Define file name for peak list. + FNAME = "%s_%s_%s_%s_max_standard.ser"%(I, deltadof2, dpwr2slock, ncyc) + sp_id = "%s_%s_%s_%s"%(I, deltadof2, dpwr2slock, ncyc) + + # Load the peak intensities. + self.interpreter.spectrum.read_intensities(file=FNAME, dir=data_path+sep+'peak_lists', spectrum_id=sp_id, int_method='height') + + # Set the peak intensity errors, as defined as the baseplane RMSD. + self.interpreter.spectrum.baseplane_rmsd(error=apod_rmsd, spectrum_id=sp_id) + + # Set the relaxation dispersion experiment type. + self.interpreter.relax_disp.exp_type(spectrum_id=sp_id, exp_type='R1rho') + + # Set The spin-lock field strength, nu1, in Hz + self.interpreter.relax_disp.spin_lock_field(spectrum_id=sp_id, field=spin_lock_field_strength) + + # Calculating the spin-lock offset in ppm, from offsets values provided in Hz. + frq_N15_Hz = yOBS * 1E6 + offset_ppm_N15 = float(deltadof2) / frq_N15_Hz * 1E6 + omega_rf_ppm = centerPPM_N15 + offset_ppm_N15 + + # Set The spin-lock offset, omega_rf, in ppm. + self.interpreter.relax_disp.spin_lock_offset(spectrum_id=sp_id, offset=omega_rf_ppm) + + # Set the relaxation times (in s). + self.interpreter.relax_fit.relax_time(time=time_sl, spectrum_id=sp_id) + + # Set the spectrometer frequency. + self.interpreter.spectrometer.frequency(id=sp_id, frq=set_sfrq, units='MHz') + + # The dispersion models. + MODELS = ['R2eff', 'No Rex', 'DPL94'] + + # The grid search size (the number of increments per dimension). + GRID_INC = 4 + + # The number of Monte Carlo simulations to be used for error analysis at the end of the analysis. + MC_NUM = 3 + + # Model selection technique. + MODSEL = 'AIC' + + # Run the analysis. + relax_disp.Relax_disp(pipe_name=pipe_name, pipe_bundle=pipe_bundle, results_dir=ds.tmpdir, models=MODELS, grid_inc=GRID_INC, mc_sim_num=MC_NUM, modsel=MODSEL) + + def test_r2eff_read(self): """Test the operation of the relax_disp.r2eff_read user function.""" _______________________________________________ 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