Changes

set FIDSCWD=$PWDset DIRS=`cat fid_files.ls| sed 's/\/fid//g'`set DIRN=`dirname cd ${FIDSDIRS[1]}`cd $DIRN

# Remove '''sleep 5''' from the script.# Click 'Save script' to make '''fid.com''' file, and 'Quit', and run the script.

== Spectral processing ==This step can be done by following wiki page [[Spectral_processing]].<br>Start '''nmrDraw''' by command nmrDraw == Convert and spectral processing all ==Now it is time we want to convert all the spectra.<br>You should have a '''fid from varian format to NMRPipe with .com''' and '''nmrproc.com''' in the first FID folder.<br>We now copy these script [[Tutorial_for_Relaxation_dispersion_analysis_cpmg_fixed_time_recorded_on_varian_as_fid_interleaved_scripts#CPMG_2_convert_and_process.sh | CPMG_2_convert_and_process.sh]] into all of the experimental folders, and execute them.

CPMG_2_convert_and_processcd $CWD set FIDS=`cat fid_files.shls`set DIRN1=`dirname $PWD/${FIDS[1]}` foreach I (`seq 2 ${#FIDS}`)set FID=${FIDS[$I]}; set DIRN=`dirname $FID`cd $DIRNecho $DIRNcp -f $DIRN1/fid.com .cp -f $DIRN1/nmrproc.com ../fid.com./nmrproc.comcd ..end</source> == Convert NMRPipe to Sparky ==Next we also want to convert them to SPARKY format.<source lang="bash">set FTS=`ls -v -d -1 */*.ft2` foreach FT ($FTS) set DNAME=`dirname $FT` set BNAME=`basename $FT` set FNAME=`echo $BNAME | cut -d'.' -f1` echo $FT $DNAME $BNAME $FNAME pipe2ucsf $FT ${DNAME}/${FNAME}.ucsfend</source> = Working with peaks = == Check the peak list matches ==Check that your peak list matches your spectrum.<br>Read the section in [[Tutorial_for_Relaxation_dispersion_analysis_cpmg_fixed_time_recorded_on_varian_as_fid_interleaved#Check_the_peak_list_matches | Check the peak list matches]]. <source lang="bash">set DIRS=`cat fid_files.ls | sed 's/\/fid//g'`sparky ${DIRS[1]}/test.ucsf

The final peak list is expected to be in:<source lang="bash">/peak_lists/peaks_corr_final.list</source>And have been saved by SPARKY, so it is in this format [[SPARKY_list]]. == Extract Check for peak movement ==Your should check, that the peaks do not move at the different experiments. Try opening some random spectra settings , and overlay them in SPARKY.<br>Read the section in [[Tutorial_for_Relaxation_dispersion_analysis_cpmg_fixed_time_recorded_on_varian_as_fid_interleaved#Check_for_peak_movement | Check for peak movement]]. <source lang="bash">sparky ${DIRS[1]}/test.ucsf ${DIRS[10]}/test.ucsf ${DIRS[25]}/test.ucsf ${DIRS[50]}/test.ucsf</source> == Measuring peak heights ==We will use the program [[NMRPipe_seriesTab | NMRPipe seriesTab]] to measure the intensities. '''seriesTab''' needs a input file, where the ppm values from Varian procpar a [[SPARKY_list | SPARKY list]] has been converted to spectral points.<br>The spectral points value depends on the spectral processing parameters. === Generate spectral point file ===Now we want to make Create a settings file we can read in relaxwith spectral point information with script [[Tutorial_for_Relaxation_dispersion_analysis_cpmg_fixed_time_recorded_on_varian_as_fid_interleaved_scripts#stPeakList.pl | stPeakList.pl]] . <source lang="bash">stPeakList.pl ${DIRS[1]}/test.ft2 ../peak_lists/peaks_corr_final.list > peaks_list.tabcat peaks_list.tab</source>

set FIDS=`cat fid_filesecho "test.ft2" > ft2_file.ls`set OUT=$PWD</exp_parameters.txtsource>

echo === Measure the height or sum in a spectral point box ===<source lang="# DIRN I deltadof2 dpwr2slock ncyc ss sfrqbash" > $OUTmkdir peak_lists foreach I line ("`seq 1 ${#FIDS}tail -n+2 exp_parameters.txt`") set FIDargv=( ${FIDS[$I]}; line ) set DIRN=`dirname $FID`1cd set I=$DIRN2 set deltadof2=`awk '/^deltadof2 /{f=1;next}f{print $2;exit}' procpar`3 set dpwr2slock=`awk '/^dpwr2slock /{f=1;next}f{print $2;exit}' procpar`4 set ncyc=`awk '/^ncyc /{f$5 set trim=1;next}f{print $2;exit}' procpar`6 set ss=`awk '/^ss /{f=1;next}f{print $2;exit}' procpar`7 set sfrq=`awk '/^sfrq /{f=1;next}f{print $2;exit}' procpar`8 echo "$DIRN I set FNAME=${I }_${deltadof2 }_${dpwr2slock }_${ncyc } cd $ss DIRN seriesTab -in ../peaks_list.tab -out ${FNAME}_max_standard.ser -list ../ft2_file.ls -max seriesTab -in ../peaks_list.tab -out ${FNAME}_max_dx1_dy1.ser -list ../ft2_file.ls -max -dx 1 -dy 1 seriesTab -in ../peaks_list.tab -out $sfrq" >> {FNAME}_sum_dx1_dy1.ser -list ../ft2_file.ls -sum -dx 1 -dy 1 cp $OUT{FNAME}_max_standard.ser ../peak_lists cd ..

You can check, if you have repetitions of experiments, by sorting the === See unique parameters, and see if they are dublicated.<br>We do this, by numerical sort columns 3,4 and 5 with the values for "deltadof2, dpwr2slock, ncyc".===

tail -n +2 exp_parameters.txt | awk '{print $3}' | sort -b -k 3k1,3n 1n | uniqtail -k n +2 exp_parameters.txt | awk '{print $4}' | sort -k1,4n 1n | uniqtail -k n +2 exp_parameters.txt | awk '{print $5}' | sort -k1,5n 1n | uniqtail -n +2 exp_parameters.txt | awk '{print $36}' | sort -k1, 1n | uniqtail -n +2 exp_parameters.txt | awk '{print $47}' | sort -k1, 1n | uniqtail -n +2 exp_parameters.txt | awk '{print $8}' | sort -k1,1n | uniq</source> == Scripts ===== 1_setup_r1rho.py ===This a script file to be able to call the setup. file: '''1_setup_r1rho.py'''.<source lang="Python"># Python module imports.from os import getcwd, sep # relax module imports.from data_store import Relax_data_store; ds = Relax_data_store() ############################################# Setup# The pipe names.if not (hasattr(ds, 'pipe_name') and hasattr(ds, 'pipe_bundle') and hasattr(ds, 'pipe_type')): # Set pipe name, bundle and type. ds.pipe_name = 'base pipe' ds.pipe_bundle = 'relax_disp' ds.pipe_type = 'relax_disp' # The data pathif not hasattr(ds, 'data_path'): ds.data_path = getcwd() ############################################ Start setup# Create the data pipe.pipe.create(pipe_name=ds.pipe_name, bundle=ds.pipe_bundle, pipe_type=ds.pipe_type) # Read the spins.spectrum.read_spins(file='1_0_46_0_max_standard.ser', dir=ds.data_path+sep+'peak_lists') # Name the isotope for field strength scaling.spin.isotope(isotope='15N') # Load the experiments settings file.expfile = open(ds.data_path+sep+'exp_parameters_sort.txt', 'r')expfileslines = expfile.readlines()expfile.close() # In MHzyOBS = 81.050# In ppmyCAR = 118.078centerPPM_N15 = yCAR ## Read the chemical shift data.chemical_shift.read(file='1_0_46_0_max_standard.ser', dir=ds.data_path+sep+'peak_lists') ## 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. spectrum.read_intensities(file=FNAME, dir=ds.data_path+sep+'peak_lists', spectrum_id=sp_id, int_method='height')  # Set the peak intensity errors, as defined as the baseplane RMSD. spectrum.baseplane_rmsd(error=apod_rmsd, spectrum_id=sp_id)  # Set the relaxation dispersion experiment type. relax_disp.exp_type(spectrum_id=sp_id, exp_type='R1rho')  # Set The spin-lock field strength, nu1, in Hz 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. relax_disp.spin_lock_offset(spectrum_id=sp_id, offset=omega_rf_ppm)  # Set the relaxation times (in s). relax_disp.relax_time(spectrum_id=sp_id, time=time_sl)  # Set the spectrometer frequency. spectrometer.frequency(id=sp_id, frq=set_sfrq, units='MHz')  # Read the R1 data# We do not read the R1 data, but rather with R1.# relax_data.read(ri_id='R1', ri_type='R1', frq=cdp.spectrometer_frq_list[0], file='R1_fitted_values.txt', dir=data_path, mol_name_col=1, res_num_col=2, res_name_col=3, spin_num_col=4, spin_name_col=5, data_col=6, error_col=7)</source> === 2_pre_run_r2eff.py ===This a script file to run the R2eff values only, with a high number of Monte Carlo simulations. file: '''2_pre_run_r2eff.py'''.<source lang="Python"># Python module imports.from os import getcwd, sepimport re # relax module imports.from auto_analyses.relax_disp import Relax_dispfrom data_store import Relax_data_store; ds = Relax_data_store()from specific_analyses.relax_disp.variables import MODEL_R2EFF  ############################################# Setup# The data pathif not hasattr(ds, 'data_path'): ds.data_path = getcwd() # The models to analyse.if not hasattr(ds, 'models'): ds.models = [MODEL_R2EFF] # The number of increments per parameter, to split up the search interval in grid search.if not hasattr(ds, 'grid_inc'): ds.grid_inc = 21 # The number of Monte-Carlo simulations, for the error analysis in the 'R2eff' model when exponential curves are fitted.# For estimating the error of the fitted R2eff values,# a high number should be provided. Later the high quality R2eff values will be read for subsequent model analyses.if not hasattr(ds, 'exp_mc_sim_num'): ds.exp_mc_sim_num = 2000 # The result directory.if not hasattr(ds, 'results_dir'): ds.results_dir = getcwd() + sep + 'results_R2eff' ## The optimisation function tolerance.## This is set to the standard value, and should not be changed.#if not hasattr(ds, 'opt_func_tol'):# ds.opt_func_tol = 1e-25#Relax_disp.opt_func_tol = ds.opt_func_tol #if not hasattr(ds, 'opt_max_iterations'):# ds.opt_max_iterations = int(1e7)#Relax_disp.opt_max_iterations = ds.opt_max_iteration  ############################################ Run script with setup.script(file='1_setup_r1rho.py', dir=ds.data_path) # To speed up the analysis, only select a few spins.deselect.all() # Load the experiments settings file.residues = open(ds.data_path+sep+'global_fit_residues.txt', 'r')residueslines = residues.readlines()residues.close() # Split the line string into number and text.r = re.compile("([a-zA-Z]+)([0-9]+)([a-zA-Z]+)(-)([a-zA-Z]+)") for i, line in enumerate(residueslines): if line[0] == "#": continue else: re_split = r.match(line) #print re_split.groups() resn = re_split.group(1) resi = int(re_split.group(2)) isotope = re_split.group(3)  select.spin(spin_id=':%i@%s'%(resi, isotope), change_all=False) # Run the analysis.Relax_disp(pipe_name=ds.pipe_name, pipe_bundle=ds.pipe_bundle, results_dir=ds.results_dir, models=ds.models, grid_inc=ds.grid_inc, exp_mc_sim_num=ds.exp_mc_sim_num)</source> === 3_analyse_models.py ===This a script file to analyse the models. file: '''3_analyse_models.py'''.<source lang="Python"># Python module imports.from os import getcwd, sepimport re # relax module imports.from auto_analyses.relax_disp import Relax_dispfrom data_store import Relax_data_store; ds = Relax_data_store()from specific_analyses.relax_disp.variables import MODEL_R2EFF, MODEL_NOREX_R1RHO, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05 ############################################# Setup# The pipe names.if not (hasattr(ds, 'pipe_name') and hasattr(ds, 'pipe_bundle') and hasattr(ds, 'pipe_type')): # Set pipe name, bundle and type. ds.pipe_name = 'base pipe' ds.pipe_bundle = 'relax_disp' ds.pipe_type = 'relax_disp' # The data pathif not hasattr(ds, 'data_path'): ds.data_path = getcwd() # The models to analyse.if not hasattr(ds, 'models'): #ds.models = [MODEL_NOREX_R1RHO, MODEL_MP05, MODEL_DPL94, MODEL_TP02, MODEL_TAP03] ds.models = [MODEL_NOREX_R1RHO, MODEL_DPL94] # The number of increments per parameter, to split up the search interval in grid search.if not hasattr(ds, 'grid_inc'): ds.grid_inc = 10 # The number of Monte-Carlo simulations for estimating the error of the parameters of the fitted models.if not hasattr(ds, 'mc_sim_num'): ds.mc_sim_num = 10 # The model selection technique. Either: 'AIC', 'AICc', 'BIC'if not hasattr(ds, 'modsel'): ds.modsel = 'AIC' # The previous result directory with R2eff values.if not hasattr(ds, 'pre_run_dir'): ds.pre_run_dir = getcwd() + sep + 'results_R2eff' + sep + 'R2eff' # The result directory.if not hasattr(ds, 'results_dir'): ds.results_dir = getcwd() + sep + 'results_models' ## The optimisation function tolerance.## This is set to the standard value, and should not be changed.#if not hasattr(ds, 'opt_func_tol'):# ds.opt_func_tol = 1e-25#Relax_disp.opt_func_tol = ds.opt_func_tol #if not hasattr(ds, 'opt_max_iterations'):# ds.opt_max_iterations = int(1e7)#Relax_disp.opt_max_iterations = ds.opt_max_iteration ########################################## Create the data pipe.pipe.create(pipe_name=ds.pipe_name, bundle=ds.pipe_bundle, pipe_type=ds.pipe_type) # Load the previous results into the base pipe.results.read(file='results', dir=ds.pre_run_dir) # If R1 is not measured, then do R1 fitting.r1_fit=True # Run the analysis.Relax_disp(pipe_name=ds.pipe_name, pipe_bundle=ds.pipe_bundle, results_dir=ds.results_dir, models=ds.models, grid_inc=ds.grid_inc, mc_sim_num=ds.mc_sim_num, modsel=ds.modsel, r1_fit=r1_fit)</source> === 4_inspect_results.py ===This a script file to inspect results in relax. file: '''4_inspect_results.py'''.<source lang="Python"># Python module imports.from os import getcwd, sepimport re # relax module imports.from pipe_control.mol_res_spin import generate_spin_string, return_spin, spin_loopfrom specific_analyses.relax_disp.data import generate_r20_key, loop_exp_frqfrom specific_analyses.relax_disp.variables import MODEL_R2EFF, MODEL_NOREX_R1RHO, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05 ############################################# Setupresults_dir = getcwd() + sep + 'results_models' # Load the previous statestate.load(state='final_state.bz2', dir=results_dir) # Display all pipespipe.display() # Define models which have been analysed.#MODELS = [MODEL_NOREX_R1RHO MODEL_MP05, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05]MODELS = [MODEL_NOREX_R1RHO, MODEL_DPL94] # Print results for each model.print("\n################")print("Printing results")print("################\n") # Store all the pipe names.pipes = [] for model in MODELS: # Skip R2eff model. if model == MODEL_R2EFF: continue  # Switch to pipe. pipe_name = '%s - relax_disp' % (model) pipes.append(pipe_name) pipe.switch(pipe_name=pipe_name) print("\nModel: %s" % (model))  # Loop over the spins. for cur_spin, mol_name, resi, resn, spin_id in spin_loop(full_info=True, return_id=True, skip_desel=True): # Generate spin string. spin_string = generate_spin_string(spin=cur_spin, mol_name=mol_name, res_num=resi, res_name=resn)  # Loop over the parameters. print("\nOptimised parameters for spin: %s" % (spin_string)) for param in cur_spin.params + ['chi2']: # Get the value. if param in ['r1_fit', 'r2']: for exp_type, frq, ei, mi in loop_exp_frq(return_indices=True): # Generate the R20 key. r20_key = generate_r20_key(exp_type=exp_type, frq=frq)  # Get the value. value = getattr(cur_spin, param)[r20_key]  # Print value. print("%-10s %-6s %-6s %3.8f" % ("Parameter:", param, "Value:", value))  # For all other parameters. else: # Get the value. value = getattr(cur_spin, param)  # Print value. print("%-10s %-6s %-6s %3.8f" % ("Parameter:", param, "Value:", value))  # Print the final pipe.pipe.switch(pipe_name='%s - relax_disp' % ('final'))print("\nFinal pipe") # Loop over the spins.for cur_spin, mol_name, resi, resn, spin_id in spin_loop(full_info=True, return_id=True, skip_desel=True): # Generate spin string. spin_string = generate_spin_string(spin=cur_spin, mol_name=mol_name, res_num=resi, res_name=resn)  # Loop over the parameters. print("\nOptimised model for spin: %s" % (spin_string)) param = 'model'  # Get the value. value = getattr(cur_spin, param) print("%-10s %-6s %-6s %6s" % ("Parameter:", param, "Value:", value))  # Print the model selectionprint("Printing the model selection")model_selection(method='AIC', modsel_pipe='test', pipes=pipes)pipe.display()</source> === 5_clustered_analyses.py ===This a script file to do a clustered analysis. file: '''5_clustered_analyses.py'''.<source lang="Python"># Python module imports.from os import getcwd, sepimport re # relax module imports.from auto_analyses.relax_disp import Relax_dispfrom data_store import Relax_data_store; ds = Relax_data_store()from pipe_control.mol_res_spin import spin_loopfrom specific_analyses.relax_disp.variables import MODEL_R2EFF, MODEL_NOREX_R1RHO, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05 ############################################# Setup# The pipe names.if not (hasattr(ds, 'pipe_name') and hasattr(ds, 'pipe_bundle') and hasattr(ds, 'pipe_type') and hasattr(ds, 'pipe_bundle_cluster')): # Set pipe name, bundle and type. ds.pipe_name = 'base pipe' ds.pipe_bundle = 'relax_disp' ds.pipe_type = 'relax_disp' ds.pipe_bundle_cluster = 'cluster' # The data pathif not hasattr(ds, 'data_path'): ds.data_path = getcwd() # The models to analyse.if not hasattr(ds, 'models'): #ds.models = [MODEL_NOREX_R1RHO, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05] ds.models = [MODEL_DPL94] # The number of increments per parameter, to split up the search interval in grid search.# This is not used, when pointing to a previous result directory.# Then an average of the previous values will be used.if not hasattr(ds, 'grid_inc'): ds.grid_inc = 10 # The number of Monte-Carlo simulations for estimating the error of the parameters of the fitted models.if not hasattr(ds, 'mc_sim_num'): ds.mc_sim_num = 10 # The model selection technique. Either: 'AIC', 'AICc', 'BIC'if not hasattr(ds, 'modsel'): ds.modsel = 'AIC' # The previous result directory with R2eff values.if not hasattr(ds, 'pre_run_dir'): ds.pre_run_dir = getcwd() + sep + 'results_models' + sep + ds.models[0] # The result directory.if not hasattr(ds, 'results_dir'): ds.results_dir = getcwd() + sep + 'results_clustering' ## The optimisation function tolerance.## This is set to the standard value, and should not be changed.#if not hasattr(ds, 'opt_func_tol'):# ds.opt_func_tol = 1e-25#Relax_disp.opt_func_tol = ds.opt_func_tol #if not hasattr(ds, 'opt_max_iterations'):# ds.opt_max_iterations = int(1e7)#Relax_disp.opt_max_iterations = ds.opt_max_iteration ########################################## Create the data pipe.ini_pipe_name = '%s - %s' % (ds.models[0], ds.pipe_bundle)pipe.create(pipe_name=ini_pipe_name, bundle=ds.pipe_bundle, pipe_type=ds.pipe_type) # Load the previous results into the base pipe.results.read(file='results', dir=ds.pre_run_dir) # Create a new pipe, where the clustering analysis will happen.# We will copy the pipe to get all information.pipe.copy(pipe_from=ini_pipe_name, pipe_to=ds.pipe_name, bundle_to=ds.pipe_bundle_cluster)pipe.switch(ds.pipe_name) pipe.display() # Now cluster spins.#relax_disp.cluster('model_cluster', ":1-100")for cur_spin, mol_name, resi, resn, spin_id in spin_loop(full_info=True, return_id=True, skip_desel=True): # Here one could write some advanced selecting rules. relax_disp.cluster('model_cluster', spin_id) # See the clustering in the current data pipe "cdp".for key, value in cdp.clustering.iteritems(): print key, value # Print parameter kex before copying.for cur_spin, mol_name, resi, resn, spin_id in spin_loop(full_info=True, return_id=True, skip_desel=True): print(cur_spin.kex) ## Make advanced parameter copy.# It is more advanced than the value.copy user function, in that clustering is taken into account. # When the destination data pipe has spin clusters defined, then the new parameter values, when required, will be taken as the median value.relax_disp.parameter_copy(pipe_from=ini_pipe_name, pipe_to=ds.pipe_name) # Print parameter kex after copying.for cur_spin, mol_name, resi, resn, spin_id in spin_loop(full_info=True, return_id=True, skip_desel=True): print(cur_spin.kex) pipe.display() # Run the analysis.Relax_disp(pipe_name=ds.pipe_name, pipe_bundle=ds.pipe_bundle_cluster, results_dir=ds.results_dir, models=ds.models, grid_inc=ds.grid_inc, mc_sim_num=ds.mc_sim_num, modsel=ds.modsel)