Changes
Jump to navigation
Jump to search
← Older edit
Tutorial for Relaxation dispersion analysis r1rho fixed time recorded on varian as sequential spectra
12,846 bytes added
,
17:13, 6 November 2015
Using the {{caution}} template.
{{caution|This tutorial is incomplete.}}
= Intro =
This tutorial is not complete.<br>
This tutorial is based on the analysis of R1rho data, analysed in a master thesis.
== 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.
# 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, sep
import re
# relax module imports.
from auto_analyses.relax_disp import Relax_disp
from data_store import Relax_data_store; ds = Relax_data_store()
from specific_analyses.relax_disp.variables import MODEL_R2EFF
#########################################
#### Setup
# The data path
if 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, sep
import re
# relax module imports.
from auto_analyses.relax_disp import Relax_disp
from 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 path
if 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, sep
import re
# relax module imports.
from pipe_control.mol_res_spin import generate_spin_string, return_spin, spin_loop
from specific_analyses.relax_disp.data import generate_r20_key, loop_exp_frq
from specific_analyses.relax_disp.variables import MODEL_R2EFF, MODEL_NOREX_R1RHO, MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05
#########################################
#### Setup
results_dir = getcwd() + sep + 'results_models'
# Load the previous state
state.load(state='final_state.bz2', dir=results_dir)
# Display all pipes
pipe.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 selection
print("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, sep
import re
# relax module imports.
from auto_analyses.relax_disp import Relax_disp
from data_store import Relax_data_store; ds = Relax_data_store()
from pipe_control.mol_res_spin import spin_loop
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') 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 path
if 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)
</source>
= See also =
[[Category:Tutorials]]
Bugman
Trusted,
Bureaucrats
4,223
edits
Navigation menu
Personal tools
Create account
Log in
Namespaces
Page
Discussion
Variants
Views
Read
View source
View history
More
Search
Navigation
Main page
Homepage
Installation guides
Community portal
Recent changes
Random page
Help
Tools
Upload file
Special pages
Printable version