Package specific_analyses :: Package model_free :: Module api
[hide private]
[frames] | no frames]

Source Code for Module specific_analyses.model_free.api

   1  ############################################################################### 
   2  #                                                                             # 
   3  # Copyright (C) 2007-2014 Edward d'Auvergne                                   # 
   4  # Copyright (C) 2007 Gary S Thompson (https://gna.org/users/varioustoxins)    # 
   5  #                                                                             # 
   6  # This file is part of the program relax (http://www.nmr-relax.com).          # 
   7  #                                                                             # 
   8  # This program is free software: you can redistribute it and/or modify        # 
   9  # it under the terms of the GNU General Public License as published by        # 
  10  # the Free Software Foundation, either version 3 of the License, or           # 
  11  # (at your option) any later version.                                         # 
  12  #                                                                             # 
  13  # This program is distributed in the hope that it will be useful,             # 
  14  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
  15  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
  16  # GNU General Public License for more details.                                # 
  17  #                                                                             # 
  18  # You should have received a copy of the GNU General Public License           # 
  19  # along with this program.  If not, see <http://www.gnu.org/licenses/>.       # 
  20  #                                                                             # 
  21  ############################################################################### 
  22   
  23  # Module docstring. 
  24  """The Lipari-Szabo model-free analysis API object.""" 
  25   
  26   
  27  # Python module imports. 
  28  import bmrblib 
  29  from copy import deepcopy 
  30  from math import pi 
  31  from minfx.grid import grid_split 
  32  from numpy import array, dot, float64, int32, zeros 
  33  from numpy.linalg import inv 
  34  from re import match, search 
  35  import string 
  36  from types import MethodType 
  37  from warnings import warn 
  38   
  39  # relax module imports. 
  40  from lib.arg_check import is_num_list, is_str_list 
  41  from lib.errors import RelaxError, RelaxFault, RelaxNoModelError, RelaxNoSequenceError, RelaxNoTensorError 
  42  from lib.float import isInf 
  43  from lib.physical_constants import h_bar, mu0, return_gyromagnetic_ratio 
  44  from lib.warnings import RelaxDeselectWarning, RelaxWarning 
  45  from multi import Processor_box 
  46  from pipe_control import diffusion_tensor, interatomic, mol_res_spin, pipes, relax_data, sequence 
  47  from pipe_control.bmrb import list_sample_conditions 
  48  from pipe_control.exp_info import bmrb_write_citations, bmrb_write_methods, bmrb_write_software 
  49  from pipe_control.interatomic import return_interatom_list 
  50  from pipe_control.mol_res_spin import count_spins, exists_mol_res_spin_data, find_index, get_molecule_names, return_spin, return_spin_from_index, return_spin_indices, spin_loop 
  51  from specific_analyses.api_base import API_base 
  52  from specific_analyses.api_common import API_common 
  53  from specific_analyses.model_free.bmrb import sf_csa_read, sf_model_free_read, to_bmrb_model 
  54  from specific_analyses.model_free.data import compare_objects 
  55  from specific_analyses.model_free.molmol import Molmol 
  56  from specific_analyses.model_free.model import determine_model_type 
  57  from specific_analyses.model_free.parameters import are_mf_params_set, assemble_param_names, assemble_param_vector, assemble_scaling_matrix, linear_constraints 
  58  from specific_analyses.model_free.optimisation import MF_grid_command, MF_memo, MF_minimise_command, grid_search_config, minimise_data_setup, relax_data_opt_structs, reset_min_stats 
  59  from specific_analyses.model_free.parameter_object import Model_free_params 
  60  from specific_analyses.model_free.pymol import Pymol 
  61  from target_functions.mf import Mf 
  62   
  63   
64 -class Model_free(API_base, API_common):
65 """Parent class containing all the model-free specific functions.""" 66 67 # Class variable for storing the class instance (for the singleton design pattern). 68 instance = None 69
70 - def __init__(self):
71 """Initialise the class by placing API_common methods into the API.""" 72 73 # Place methods into the API. 74 self.base_data_loop = self._base_data_loop_spin 75 self.return_error = self._return_error_relax_data 76 self.return_value = self._return_value_general 77 self.sim_pack_data = self._sim_pack_relax_data 78 79 # Initialise the macro classes. 80 self._molmol_macros = Molmol() 81 self._pymol_macros = Pymol() 82 83 # Alias the macro creation methods. 84 self.pymol_macro = self._pymol_macros.create_macro 85 self.molmol_macro = self._molmol_macros.create_macro 86 87 # Place a copy of the parameter list object in the instance namespace. 88 self._PARAMS = Model_free_params()
89 90
91 - def back_calc_ri(self, spin_index=None, ri_id=None, ri_type=None, frq=None):
92 """Back-calculation of relaxation data from the model-free parameter values. 93 94 @keyword spin_index: The global spin index. 95 @type spin_index: int 96 @keyword ri_id: The relaxation data ID string. 97 @type ri_id: str 98 @keyword ri_type: The relaxation data type. 99 @type ri_type: str 100 @keyword frq: The field strength. 101 @type frq: float 102 @return: The back calculated relaxation data value corresponding to the index. 103 @rtype: float 104 """ 105 106 # Get the spin container. 107 spin, spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 108 109 # Missing interatomic vectors. 110 if hasattr(cdp, 'diff_tensor') and (cdp.diff_tensor.type == 'spheroid' or cdp.diff_tensor.type == 'ellipsoid'): 111 interatoms = interatomic.return_interatom_list(spin_id) 112 for i in range(len(interatoms)): 113 # No dipolar relaxation mechanism. 114 if not interatoms[i].dipole_pair: 115 continue 116 117 # Check the unit vectors. 118 if not hasattr(interatoms[i], 'vector') or interatoms[i].vector == None: 119 warn(RelaxDeselectWarning(spin_id, 'missing structural data')) 120 return 121 122 # Execute the over-fit deselection. 123 self.overfit_deselect(data_check=False, verbose=False) 124 125 # Get the relaxation value from the minimise function. 126 value = self.minimise(min_algor='back_calc', min_options=(spin_index, ri_id, ri_type, frq)) 127 128 # Return the relaxation value. 129 return value
130 131
132 - def bmrb_read(self, file_path, version=None, sample_conditions=None):
133 """Read the model-free results from a BMRB NMR-STAR v3.1 formatted file. 134 135 @param file_path: The full file path. 136 @type file_path: str 137 @keyword version: The BMRB version to force the reading. 138 @type version: None or str 139 @keyword sample_conditions: The sample condition label to read. Only one sample condition can be read per data pipe. 140 @type sample_conditions: None or str 141 """ 142 143 # Initialise the NMR-STAR data object. 144 star = bmrblib.create_nmr_star('relax_model_free_results', file_path, version) 145 146 # Read the contents of the STAR formatted file. 147 star.read() 148 149 # The sample conditions. 150 sample_conds = list_sample_conditions(star) 151 if sample_conditions and sample_conditions not in sample_conds: 152 raise RelaxError("The sample conditions label '%s' does not correspond to any of the labels in the file: %s" % (sample_conditions, sample_conds)) 153 if not sample_conditions and len(sample_conds) > 1: 154 raise RelaxError("Only one of the sample conditions in %s can be loaded per relax data pipe." % sample_conds) 155 156 # The diffusion tensor. 157 diffusion_tensor.bmrb_read(star) 158 159 # Generate the molecule and residue containers from the entity records. 160 mol_res_spin.bmrb_read(star) 161 162 # Read the relaxation data saveframes. 163 relax_data.bmrb_read(star, sample_conditions=sample_conditions) 164 165 # Read the model-free data saveframes. 166 sf_model_free_read(star, sample_conditions=sample_conditions) 167 168 # Read the CSA data saveframes. 169 sf_csa_read(star)
170 171
172 - def bmrb_write(self, file_path, version=None):
173 """Write the model-free results to a BMRB NMR-STAR v3.1 formatted file. 174 175 @param file_path: The full file path. 176 @type file_path: str 177 @keyword version: The BMRB NMR-STAR dictionary format to output to. 178 @type version: str 179 """ 180 181 # Alias the current data pipe. 182 cdp = pipes.get_pipe() 183 184 # Initialise the NMR-STAR data object. 185 star = bmrblib.create_nmr_star('relax_model_free_results', file_path, version) 186 187 # Global minimisation stats. 188 global_chi2 = None 189 if hasattr(cdp, 'chi2'): 190 global_chi2 = cdp.chi2 191 192 # Rex frq. 193 rex_frq = cdp.spectrometer_frq[cdp.ri_ids[0]] 194 195 # Initialise the spin specific data lists. 196 mol_name_list = [] 197 res_num_list = [] 198 res_name_list = [] 199 atom_name_list = [] 200 201 csa_list = [] 202 r_list = [] 203 isotope_list = [] 204 element_list = [] 205 206 local_tm_list = [] 207 s2_list = [] 208 s2f_list = [] 209 s2s_list = [] 210 te_list = [] 211 tf_list = [] 212 ts_list = [] 213 rex_list = [] 214 215 local_tm_err_list = [] 216 s2_err_list = [] 217 s2f_err_list = [] 218 s2s_err_list = [] 219 te_err_list = [] 220 tf_err_list = [] 221 ts_err_list = [] 222 rex_err_list = [] 223 224 chi2_list = [] 225 model_list = [] 226 227 # Store the spin specific data in lists for later use. 228 for spin, mol_name, res_num, res_name, spin_id in spin_loop(full_info=True, return_id=True): 229 # Skip the protons. 230 if spin.name == 'H' or (hasattr(spin, 'element') and spin.element == 'H'): 231 warn(RelaxWarning("Skipping the proton spin '%s'." % spin_id)) 232 continue 233 234 # Check the data for None (not allowed in BMRB!). 235 if res_num == None: 236 raise RelaxError("For the BMRB, the residue of spin '%s' must be numbered." % spin_id) 237 if res_name == None: 238 raise RelaxError("For the BMRB, the residue of spin '%s' must be named." % spin_id) 239 if spin.name == None: 240 raise RelaxError("For the BMRB, the spin '%s' must be named." % spin_id) 241 if not hasattr(spin, 'isotope') or spin.isotope == None: 242 raise RelaxError("For the BMRB, the spin isotope type of '%s' must be specified." % spin_id) 243 if not hasattr(spin, 'element') or spin.element == None: 244 raise RelaxError("For the BMRB, the spin element type of '%s' must be specified. Please use the spin user function for setting the element type." % spin_id) 245 246 # The molecule/residue/spin info. 247 mol_name_list.append(mol_name) 248 res_num_list.append(res_num) 249 res_name_list.append(res_name) 250 atom_name_list.append(spin.name) 251 252 # CSA values. 253 if hasattr(spin, 'csa'): 254 csa_list.append(spin.csa * 1e6) # In ppm. 255 else: 256 csa_list.append(None) 257 258 # Interatomic distances. 259 interatoms = return_interatom_list(spin_id) 260 for i in range(len(interatoms)): 261 # No relaxation mechanism. 262 if not interatoms[i].dipole_pair: 263 continue 264 265 # Add the interatomic distance. 266 if hasattr(interatoms[i], 'r'): 267 r_list.append(interatoms[i].r) 268 else: 269 r_list.append(None) 270 271 # Stop adding. 272 break 273 274 # The nuclear isotope. 275 if hasattr(spin, 'isotope'): 276 isotope_list.append(int(spin.isotope.strip(string.ascii_letters))) 277 else: 278 isotope_list.append(None) 279 280 # The element. 281 if hasattr(spin, 'element'): 282 element_list.append(spin.element) 283 else: 284 element_list.append(None) 285 286 # Diffusion tensor. 287 if hasattr(spin, 'local_tm'): 288 local_tm_list.append(spin.local_tm) 289 else: 290 local_tm_list.append(None) 291 if hasattr(spin, 'local_tm_err'): 292 local_tm_err_list.append(spin.local_tm_err) 293 else: 294 local_tm_err_list.append(None) 295 296 # Model-free parameter values. 297 s2_list.append(spin.s2) 298 s2f_list.append(spin.s2f) 299 s2s_list.append(spin.s2s) 300 te_list.append(spin.te) 301 tf_list.append(spin.tf) 302 ts_list.append(spin.ts) 303 if spin.rex == None: 304 rex_list.append(None) 305 else: 306 rex_list.append(spin.rex / (2.0*pi*rex_frq**2)) 307 308 # Model-free parameter errors. 309 params = ['s2', 's2f', 's2s', 'te', 'tf', 'ts', 'rex'] 310 for param in params: 311 # The error list. 312 err_list = locals()[param+'_err_list'] 313 314 # Append the error. 315 if hasattr(spin, param+'_err'): 316 # The value. 317 val = getattr(spin, param+'_err') 318 319 # Scaling. 320 if param == 'rex' and val != None: 321 val = val / (2.0*pi*rex_frq**2) 322 323 # Append. 324 err_list.append(val) 325 326 # Or None. 327 else: 328 err_list.append(None) 329 330 331 # Opt stats. 332 chi2_list.append(spin.chi2) 333 334 # Model-free model. 335 model_list.append(to_bmrb_model(spin.model)) 336 337 # Convert the molecule names into the entity IDs. 338 entity_ids = zeros(len(mol_name_list), int32) 339 mol_names = get_molecule_names() 340 for i in range(len(mol_name_list)): 341 for j in range(len(mol_names)): 342 if mol_name_list[i] == mol_names[j]: 343 entity_ids[i] = j+1 344 345 346 # Create Supergroup 2 : The citations. 347 ###################################### 348 349 # Generate the citations saveframe. 350 bmrb_write_citations(star) 351 352 353 # Create Supergroup 3 : The molecular assembly saveframes. 354 ########################################################## 355 356 # Generate the entity saveframe. 357 mol_res_spin.bmrb_write_entity(star) 358 359 360 # Create Supergroup 4: The experimental descriptions saveframes. 361 ################################################################# 362 363 # Generate the method saveframes. 364 bmrb_write_methods(star) 365 366 # Generate the software saveframe. 367 software_ids, software_labels = bmrb_write_software(star) 368 369 370 # Create Supergroup 5 : The NMR parameters saveframes. 371 ###################################################### 372 373 # Generate the CSA saveframe. 374 star.chem_shift_anisotropy.add(entity_ids=entity_ids, res_nums=res_num_list, res_names=res_name_list, atom_names=atom_name_list, atom_types=element_list, isotope=isotope_list, csa=csa_list) 375 376 377 # Create Supergroup 6 : The kinetic data saveframes. 378 #################################################### 379 380 # Generate the relaxation data saveframes. 381 relax_data.bmrb_write(star) 382 383 384 # Create Supergroup 7 : The thermodynamics saveframes. 385 ###################################################### 386 387 # Get the relax software id. 388 for i in range(len(software_ids)): 389 if software_labels[i] == 'relax': 390 software_id = software_ids[i] 391 392 # Generate the model-free data saveframe. 393 star.model_free.add(global_chi2=global_chi2, software_ids=[software_id], software_labels=['relax'], entity_ids=entity_ids, res_nums=res_num_list, res_names=res_name_list, atom_names=atom_name_list, atom_types=element_list, isotope=isotope_list, bond_length=r_list, local_tc=local_tm_list, s2=s2_list, s2f=s2f_list, s2s=s2s_list, te=te_list, tf=tf_list, ts=ts_list, rex=rex_list, local_tc_err=local_tm_err_list, s2_err=s2_err_list, s2f_err=s2f_err_list, s2s_err=s2s_err_list, te_err=te_err_list, tf_err=tf_err_list, ts_err=ts_err_list, rex_err=rex_err_list, rex_frq=rex_frq, chi2=chi2_list, model_fit=model_list) 394 395 396 # Create Supergroup 8 : The structure determination saveframes. 397 ############################################################### 398 399 # Generate the diffusion tensor saveframes. 400 diffusion_tensor.bmrb_write(star) 401 402 403 # Write the contents to the STAR formatted file. 404 star.write()
405 406
407 - def calculate(self, spin_id=None, verbosity=1, sim_index=None):
408 """Calculation of the model-free chi-squared value. 409 410 @keyword spin_id: The spin identification string. 411 @type spin_id: str 412 @keyword verbosity: The amount of information to print. The higher the value, the greater the verbosity. 413 @type verbosity: int 414 @keyword sim_index: The optional MC simulation index. 415 @type sim_index: int 416 """ 417 418 # Test if sequence data is loaded. 419 if not exists_mol_res_spin_data(): 420 raise RelaxNoSequenceError 421 422 # Determine the model type. 423 model_type = determine_model_type() 424 425 # Test if diffusion tensor data exists. 426 if model_type != 'local_tm' and not diffusion_tensor.diff_data_exists(): 427 raise RelaxNoTensorError('diffusion') 428 429 # Test if the PDB file has been loaded. 430 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere' and not hasattr(cdp, 'structure'): 431 raise RelaxNoPdbError 432 433 # Loop over the spins. 434 for spin, id in spin_loop(spin_id, return_id=True): 435 # Skip deselected spins. 436 if not spin.select: 437 continue 438 439 # Test if the model-free model has been setup. 440 if not spin.model: 441 raise RelaxNoModelError 442 443 # Test if the nuclear isotope type has been set. 444 if not hasattr(spin, 'isotope'): 445 raise RelaxSpinTypeError 446 447 # Test if the model-free parameter values exist. 448 unset_param = are_mf_params_set(spin) 449 if unset_param != None: 450 raise RelaxNoValueError(unset_param) 451 452 # Test if the CSA value has been set. 453 if not hasattr(spin, 'csa') or spin.csa == None: 454 raise RelaxNoValueError("CSA") 455 456 # Test the interatomic data. 457 interatoms = return_interatom_list(spin_id) 458 for interatom in interatoms: 459 # No relaxation mechanism. 460 if not interatom.dipole_pair: 461 continue 462 463 # Test if unit vectors exist. 464 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere' and not hasattr(interatom, 'vector'): 465 raise RelaxNoVectorsError 466 467 # Test if multiple unit vectors exist. 468 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere' and hasattr(interatom, 'vector') and is_num_list(interatom.vector[0], raise_error=False): 469 raise RelaxMultiVectorError 470 471 # The interacting spin. 472 if spin_id != interatom.spin_id1: 473 spin_id2 = interatom.spin_id1 474 else: 475 spin_id2 = interatom.spin_id2 476 spin2 = return_spin(spin_id2) 477 478 # Test if the nuclear isotope type has been set. 479 if not hasattr(spin2, 'isotope'): 480 raise RelaxSpinTypeError 481 482 # Test if the interatomic distance has been set. 483 if not hasattr(interatom, 'r') or interatom.r == None: 484 raise RelaxNoValueError("interatomic distance", spin_id=spin_id, spin_id2=spin_id2) 485 486 # Skip spins where there is no data or errors. 487 if not hasattr(spin, 'ri_data') or not hasattr(spin, 'ri_data_err'): 488 continue 489 490 # Make sure that the errors are strictly positive numbers. 491 for ri_id in cdp.ri_ids: 492 # Alias. 493 err = spin.ri_data_err[ri_id] 494 495 # Checks. 496 if err != None and err == 0.0: 497 raise RelaxError("Zero error for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (id, ri_id)) 498 elif err != None and err < 0.0: 499 raise RelaxError("Negative error of %s for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (err, id, ri_id)) 500 501 # Create the initial parameter vector. 502 param_vector = assemble_param_vector(spin=spin, sim_index=sim_index) 503 504 # The relaxation data optimisation structures. 505 data = relax_data_opt_structs(spin, sim_index=sim_index) 506 507 # The spin data. 508 ri_data = [array(data[0])] 509 ri_data_err = [array(data[1])] 510 num_frq = [data[2]] 511 num_ri = [data[3]] 512 ri_labels = [data[4]] 513 frq = [data[5]] 514 remap_table = [data[6]] 515 noe_r1_table = [data[7]] 516 gx = [return_gyromagnetic_ratio(spin.isotope)] 517 if sim_index == None: 518 csa = [spin.csa] 519 else: 520 csa = [spin.csa_sim[sim_index]] 521 522 # The interatomic data. 523 interatoms = return_interatom_list(spin_id) 524 for i in range(len(interatoms)): 525 # No relaxation mechanism. 526 if not interatoms[i].dipole_pair: 527 continue 528 529 # The surrounding spins. 530 if spin_id != interatoms[i].spin_id1: 531 spin_id2 = interatoms[i].spin_id1 532 else: 533 spin_id2 = interatoms[i].spin_id2 534 spin2 = return_spin(spin_id2) 535 536 # The data. 537 if sim_index == None: 538 r = [interatoms[i].r] 539 else: 540 r = [interatoms[i].r_sim[sim_index]] 541 542 # Vectors. 543 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere': 544 xh_unit_vectors = [interatoms[i].vector] 545 else: 546 xh_unit_vectors = [None] 547 548 # Gyromagnetic ratios. 549 gh = [return_gyromagnetic_ratio(spin2.isotope)] 550 551 # Count the number of model-free parameters for the residue index. 552 num_params = [len(spin.params)] 553 554 # Repackage the parameter values as a local model (ignore if the diffusion tensor is not fixed). 555 param_values = [assemble_param_vector(model_type='mf')] 556 557 # Package the diffusion tensor parameters. 558 if model_type == 'local_tm': 559 diff_params = [spin.local_tm] 560 diff_type = 'sphere' 561 else: 562 # Diff type. 563 diff_type = cdp.diff_tensor.type 564 565 # Spherical diffusion. 566 if diff_type == 'sphere': 567 diff_params = [cdp.diff_tensor.tm] 568 569 # Spheroidal diffusion. 570 elif diff_type == 'spheroid': 571 diff_params = [cdp.diff_tensor.tm, cdp.diff_tensor.Da, cdp.diff_tensor.theta, cdp.diff_tensor.phi] 572 573 # Ellipsoidal diffusion. 574 elif diff_type == 'ellipsoid': 575 diff_params = [cdp.diff_tensor.tm, cdp.diff_tensor.Da, cdp.diff_tensor.Dr, cdp.diff_tensor.alpha, cdp.diff_tensor.beta, cdp.diff_tensor.gamma] 576 577 # Initialise the model-free function. 578 mf = Mf(init_params=param_vector, model_type='mf', diff_type=diff_type, diff_params=diff_params, num_spins=1, equations=[spin.equation], param_types=[spin.params], param_values=param_values, relax_data=ri_data, errors=ri_data_err, bond_length=r, csa=csa, num_frq=num_frq, frq=frq, num_ri=num_ri, remap_table=remap_table, noe_r1_table=noe_r1_table, ri_labels=ri_labels, gx=gx, gh=gh, h_bar=h_bar, mu0=mu0, num_params=num_params, vectors=xh_unit_vectors) 579 580 # Chi-squared calculation. 581 try: 582 chi2 = mf.func(param_vector) 583 except OverflowError: 584 chi2 = 1e200 585 586 # Global chi-squared value. 587 if model_type == 'all' or model_type == 'diff': 588 cdp.chi2 = cdp.chi2 + chi2 589 else: 590 spin.chi2 = chi2
591 592
593 - def create_mc_data(self, data_id=None):
594 """Create the Monte Carlo Ri data. 595 596 @keyword data_id: The spin identification string, as yielded by the base_data_loop() generator method. 597 @type data_id: str 598 @return: The Monte Carlo simulation data. 599 @rtype: list of floats 600 """ 601 602 # Initialise the MC data structure. 603 mc_data = [] 604 605 # Get the spin container and global spin index. 606 spin = return_spin(data_id) 607 global_index = find_index(data_id) 608 609 # Skip deselected spins. 610 if not spin.select: 611 return 612 613 # Test if the model is set, returning None to signal that the spin should be skipped. 614 if not hasattr(spin, 'model') or not spin.model: 615 return None 616 617 # Loop over the relaxation data. 618 for ri_id in cdp.ri_ids: 619 # Back calculate the value. 620 value = self.back_calc_ri(spin_index=global_index, ri_id=ri_id, ri_type=cdp.ri_type[ri_id], frq=cdp.spectrometer_frq[ri_id]) 621 622 # Append the value. 623 mc_data.append(value) 624 625 # Return the data. 626 return mc_data
627 628
629 - def data_init(self, data_cont, sim=False):
630 """Initialise the spin specific data structures. 631 632 @param data_cont: The spin data container. 633 @type data_cont: SpinContainer instance 634 @keyword sim: The Monte Carlo simulation flag, which if true will initialise the simulation data structure. 635 @type sim: bool 636 """ 637 638 # Loop over the data structure names. 639 for name in self._PARAMS.loop(scope='spin'): 640 # Blacklisted data structures. 641 if name in ['ri_data', 'ri_data_bc', 'ri_data_err']: 642 continue 643 644 # Data structures which are initially empty arrays. 645 list_data = [ 'params' ] 646 if name in list_data: 647 init_data = [] 648 649 # Set everything else initially to None or False. 650 init_data = None 651 if self._PARAMS.type(name) == bool: 652 init_data = False 653 if name == 'select': 654 init_data = True 655 656 # If the name is not in 'data_cont', add it. 657 if not hasattr(data_cont, name): 658 setattr(data_cont, name, init_data)
659 660
661 - def deselect(self, model_info, sim_index=None):
662 """Deselect models or simulations. 663 664 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 665 @type model_info: int 666 @keyword sim_index: The optional Monte Carlo simulation index. If None, then models will be deselected, otherwise the given simulation will. 667 @type sim_index: None or int 668 """ 669 670 # Determine the model type. 671 model_type = determine_model_type() 672 673 # Local models. 674 if model_type == 'mf' or model_type == 'local_tm': 675 # Get the spin. 676 spin = return_spin_from_index(model_info) 677 678 # Spin deselection. 679 if sim_index == None: 680 spin.select = False 681 682 # Simulation deselection. 683 else: 684 spin.select_sim[sim_index] = False 685 686 # Global models. 687 else: 688 # Global model deselection. 689 if sim_index == None: 690 raise RelaxError("Cannot deselect the global model.") 691 692 # Simulation deselection. 693 else: 694 # Deselect. 695 cdp.select_sim[sim_index] = False
696 697
698 - def duplicate_data(self, pipe_from=None, pipe_to=None, model_info=None, global_stats=False, verbose=True):
699 """Duplicate the data specific to a single model-free model. 700 701 @keyword pipe_from: The data pipe to copy the data from. 702 @type pipe_from: str 703 @keyword pipe_to: The data pipe to copy the data to. 704 @type pipe_to: str 705 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 706 @type model_info: int 707 @keyword global_stats: The global statistics flag. 708 @type global_stats: bool 709 @keyword verbose: A flag which if True will cause info about each spin to be printed out as the sequence is generated. 710 @type verbose: bool 711 """ 712 713 # Arg tests. 714 if model_info == None: 715 raise RelaxError("The model_info argument cannot be None.") 716 717 # First create the pipe_to data pipe, if it doesn't exist, but don't switch to it. 718 if not pipes.has_pipe(pipe_to): 719 pipes.create(pipe_to, pipe_type='mf', switch=False) 720 721 # Get the data pipes. 722 dp_from = pipes.get_pipe(pipe_from) 723 dp_to = pipes.get_pipe(pipe_to) 724 725 # Duplicate all non-sequence specific data. 726 for data_name in dir(dp_from): 727 # Skip the container objects. 728 if data_name in ['diff_tensor', 'mol', 'interatomic', 'structure', 'exp_info']: 729 continue 730 731 # Skip special objects. 732 if search('^_', data_name) or data_name in list(dp_from.__class__.__dict__.keys()): 733 continue 734 735 # Get the original object. 736 data_from = getattr(dp_from, data_name) 737 738 # The data already exists. 739 if hasattr(dp_to, data_name): 740 # Get the object in the target pipe. 741 data_to = getattr(dp_to, data_name) 742 743 # The data must match! 744 if data_from != data_to: 745 raise RelaxError("The object " + repr(data_name) + " is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 746 747 # Skip the data. 748 continue 749 750 # Duplicate the data. 751 setattr(dp_to, data_name, deepcopy(data_from)) 752 753 # Diffusion tensor comparison. 754 if hasattr(dp_from, 'diff_tensor'): 755 # Duplicate the tensor if it doesn't exist. 756 if not hasattr(dp_to, 'diff_tensor'): 757 setattr(dp_to, 'diff_tensor', deepcopy(dp_from.diff_tensor)) 758 759 # Otherwise compare the objects inside the container. 760 else: 761 # Loop over the modifiable objects. 762 for data_name in dp_from.diff_tensor._mod_attr: 763 # Get the original object. 764 data_from = None 765 if hasattr(dp_from.diff_tensor, data_name): 766 data_from = getattr(dp_from.diff_tensor, data_name) 767 768 # Get the target object. 769 if data_from and not hasattr(dp_to.diff_tensor, data_name): 770 raise RelaxError("The diffusion tensor object " + repr(data_name) + " of the " + repr(pipe_from) + " data pipe is not located in the " + repr(pipe_to) + " data pipe.") 771 elif data_from: 772 data_to = getattr(dp_to.diff_tensor, data_name) 773 else: 774 continue 775 776 # The data must match! 777 if data_from != data_to: 778 raise RelaxError("The object " + repr(data_name) + " is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 779 780 # Structure comparison. 781 if hasattr(dp_from, 'structure'): 782 # Duplicate the tensor if it doesn't exist. 783 if not hasattr(dp_to, 'structure'): 784 setattr(dp_to, 'structure', deepcopy(dp_from.structure)) 785 786 # Otherwise compare the objects inside the container. 787 else: 788 # Modifiable object checks. 789 compare_objects(dp_from.structure, dp_to.structure, pipe_from, pipe_to) 790 791 # Tests for the model and molecule containers. 792 if len(dp_from.structure.structural_data) != len(dp_from.structure.structural_data): 793 raise RelaxError("The number of structural models is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 794 795 # Loop over the models. 796 for i in range(len(dp_from.structure.structural_data)): 797 # Alias. 798 model_from = dp_from.structure.structural_data[i] 799 model_to = dp_to.structure.structural_data[i] 800 801 # Model numbers. 802 if model_from.num != model_to.num: 803 raise RelaxError("The structure models are not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 804 805 # Molecule number. 806 if len(model_from.mol) != len(model_to.mol): 807 raise RelaxError("The number of molecules is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 808 809 # Loop over the models. 810 for mol_index in range(len(model_from.mol)): 811 # Modifiable object checks. 812 compare_objects(model_from.mol[mol_index], model_to.mol[mol_index], pipe_from, pipe_to) 813 814 # No sequence data, so skip the rest. 815 if dp_from.mol.is_empty(): 816 return 817 818 # Duplicate the sequence data if it doesn't exist. 819 if dp_to.mol.is_empty(): 820 sequence.copy(pipe_from=pipe_from, pipe_to=pipe_to, preserve_select=True, verbose=verbose) 821 822 # Duplicate the interatomic data if it doesn't exist. 823 if dp_to.interatomic.is_empty(): 824 interatomic.copy(pipe_from=pipe_from, pipe_to=pipe_to, verbose=verbose) 825 826 # Determine the model type of the original data pipe. 827 pipes.switch(pipe_from) 828 model_type = determine_model_type() 829 830 # Sequence specific data. 831 spin, spin_id = return_spin_from_index(model_info, pipe=pipe_from, return_spin_id=True) 832 if model_type == 'mf' or (model_type == 'local_tm' and not global_stats): 833 # Duplicate the spin specific data. 834 for name in dir(spin): 835 # Skip special objects. 836 if search('^__', name): 837 continue 838 839 # Get the object. 840 obj = getattr(spin, name) 841 842 # Skip methods. 843 if isinstance(obj, MethodType): 844 continue 845 846 # Duplicate the object. 847 new_obj = deepcopy(getattr(spin, name)) 848 setattr(dp_to.mol[spin._mol_index].res[spin._res_index].spin[spin._spin_index], name, new_obj) 849 850 # Duplicate the relaxation active spins which have not been copied yet. 851 interatoms = interatomic.return_interatom_list(spin_id) 852 for interatom in interatoms: 853 # No relaxation mechanism. 854 if not interatom.dipole_pair: 855 continue 856 857 # The interacting spin. 858 if spin_id != interatom.spin_id1: 859 spin_id2 = interatom.spin_id1 860 else: 861 spin_id2 = interatom.spin_id2 862 spin2 = return_spin(spin_id2) 863 864 # Duplicate the spin specific data. 865 mol_index, res_index, spin_index = return_spin_indices(spin_id2) 866 dp_to.mol[mol_index].res[res_index].spin[spin_index] = deepcopy(spin2) 867 868 # Other data types. 869 else: 870 # Duplicate all the spin specific data. 871 dp_to.mol = deepcopy(dp_from.mol)
872 873
874 - def eliminate(self, name, value, model_info, args, sim=None):
875 """Model-free model elimination, parameter by parameter. 876 877 @param name: The parameter name. 878 @type name: str 879 @param value: The parameter value. 880 @type value: float 881 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 882 @type model_info: int 883 @param args: The c1 and c2 elimination constant overrides. 884 @type args: None or tuple of float 885 @keyword sim: The Monte Carlo simulation index. 886 @type sim: int 887 @return: True if the model is to be eliminated, False otherwise. 888 @rtype: bool 889 """ 890 891 # Default values. 892 c1 = 50.0 * 1e-9 893 c2 = 1.5 894 895 # Depack the arguments. 896 if args != None: 897 c1, c2 = args 898 899 # Determine the model type. 900 model_type = determine_model_type() 901 902 # Can't handle this one yet! 903 if model_type != 'mf' and model_type != 'local_tm': 904 raise RelaxError("Elimination of the global model is not yet supported.") 905 906 # Get the spin and it's id string. 907 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 908 909 # Get the tm value. 910 if model_type == 'local_tm': 911 tm = spin.local_tm 912 else: 913 tm = cdp.diff_tensor.tm 914 915 # No tm value set, so skip the tests (no elimination). 916 if tm == None: 917 return False 918 919 # Local tm. 920 if name == 'local_tm' and value >= c1: 921 if sim == None: 922 print("Data pipe '%s': The local tm parameter of %.5g is greater than %.5g, eliminating spin system '%s'." % (pipes.cdp_name(), value, c1, spin_id)) 923 else: 924 print("Data pipe '%s': The local tm parameter of %.5g is greater than %.5g, eliminating simulation %i of spin system '%s'." % (pipes.cdp_name(), value, c1, sim, spin_id)) 925 return True 926 927 # Internal correlation times. 928 if match('t[efs]', name) and value >= c2 * tm: 929 if sim == None: 930 print("Data pipe '%s': The %s value of %.5g is greater than %.5g, eliminating spin system '%s'." % (pipes.cdp_name(), name, value, c2*tm, spin_id)) 931 else: 932 print("Data pipe '%s': The %s value of %.5g is greater than %.5g, eliminating simulation %i of spin system '%s'." % (pipes.cdp_name(), name, value, c2*tm, sim, spin_id)) 933 return True 934 935 # Accept model. 936 return False
937 938
939 - def get_param_names(self, model_info=None):
940 """Return a vector of parameter names. 941 942 @keyword model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 943 @type model_info: int 944 @return: The vector of parameter names. 945 @rtype: list of str 946 """ 947 948 # Determine the model type. 949 model_type = determine_model_type() 950 951 # Get the spin ids. 952 if model_type == 'mf' or model_type == 'local_tm': 953 # Get the spin and it's id string. 954 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 955 else: 956 spin_id = None 957 958 # Assemble and return the parameter names. 959 return assemble_param_names(model_type, spin_id=spin_id)
960 961
962 - def get_param_values(self, model_info=None, sim_index=None):
963 """Return a vector of parameter values. 964 965 @keyword model_info: The model index from model_info(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 966 @type model_info: int 967 @keyword sim_index: The Monte Carlo simulation index. 968 @type sim_index: int 969 @return: The vector of parameter values. 970 @rtype: list of str 971 """ 972 973 # Test if the model-free models have been set up. 974 for spin in spin_loop(): 975 # Skip deselected spins. 976 if not spin.select: 977 continue 978 979 # Not setup. 980 if not spin.model: 981 raise RelaxNoModelError 982 983 # Determine the model type. 984 model_type = determine_model_type() 985 986 # Set the spin container (to None if the model is global). 987 if model_type == 'mf' or model_type == 'local_tm': 988 spin = return_spin_from_index(model_info) 989 else: 990 spin = None 991 992 # Assemble the parameter values and return them. 993 return assemble_param_vector(spin=spin, sim_index=sim_index, model_type=model_type)
994 995
996 - def grid_search(self, lower=None, upper=None, inc=None, constraints=True, verbosity=1, sim_index=None):
997 """The model-free grid search function. 998 999 @keyword lower: The lower bounds of the grid search which must be equal to the 1000 number of parameters in the model. 1001 @type lower: array of numbers 1002 @keyword upper: The upper bounds of the grid search which must be equal to the 1003 number of parameters in the model. 1004 @type upper: array of numbers 1005 @keyword inc: The increments for each dimension of the space for the grid search. 1006 The number of elements in the array must equal to the number of 1007 parameters in the model. 1008 @type inc: array of int 1009 @keyword constraints: If True, constraints are applied during the grid search (eliminating 1010 parts of the grid). If False, no constraints are used. 1011 @type constraints: bool 1012 @keyword verbosity: A flag specifying the amount of information to print. The higher 1013 the value, the greater the verbosity. 1014 @type verbosity: int 1015 @keyword sim_index: The index of the simulation to apply the grid search to. If None, 1016 the normal model is optimised. 1017 @type sim_index: int 1018 """ 1019 1020 # Minimisation. 1021 self.minimise(min_algor='grid', lower=lower, upper=upper, inc=inc, constraints=constraints, verbosity=verbosity, sim_index=sim_index)
1022 1023
1024 - def map_bounds(self, param, spin_id=None):
1025 """Create bounds for the OpenDX mapping function. 1026 1027 @param param: The name of the parameter to return the lower and upper bounds of. 1028 @type param: str 1029 @param spin_id: The spin identification string. 1030 @type spin_id: str 1031 @return: The upper and lower bounds of the parameter. 1032 @rtype: list of float 1033 """ 1034 1035 # Diffusion tensor bounds. 1036 if self._PARAMS.scope(param) == 'global': 1037 return diffusion_tensor.map_bounds(param) 1038 1039 # Get the spin. 1040 spin = return_spin(spin_id) 1041 1042 # {s2, s2f, s2s}. 1043 if search('^s2', param): 1044 return [0.0, 1.0] 1045 1046 # {local tm, te, tf, ts}. 1047 elif search('^t', param) or param == 'local_tm': 1048 return [0.0, 1e-8] 1049 1050 # Rex. 1051 elif param == 'rex': 1052 return [0.0, 30.0 / (2.0 * pi * cdp.spectrometer_frq[cdp.ri_ids[0]])**2] 1053 1054 # Interatomic distances. 1055 elif param == 'r': 1056 return [1.0 * 1e-10, 1.1 * 1e-10] 1057 1058 # CSA. 1059 elif param == 'csa': 1060 return [-100 * 1e-6, -300 * 1e-6]
1061 1062
1063 - def minimise(self, min_algor=None, min_options=None, func_tol=None, grad_tol=None, max_iterations=None, constraints=False, scaling=True, verbosity=0, sim_index=None, lower=None, upper=None, inc=None):
1064 """Model-free minimisation function. 1065 1066 Three categories of models exist for which the approach to minimisation is different. These 1067 are: 1068 1069 Single spin optimisations: The 'mf' and 'local_tm' model types which are the 1070 model-free parameters for single spins, optionally with a local tm parameter. These 1071 models have no global parameters. 1072 1073 Diffusion tensor optimisations: The 'diff' diffusion tensor model type. No spin 1074 specific parameters exist. 1075 1076 Optimisation of everything: The 'all' model type consisting of all model-free and all 1077 diffusion tensor parameters. 1078 1079 1080 @keyword min_algor: The minimisation algorithm to use. 1081 @type min_algor: str 1082 @keyword min_options: An array of options to be used by the minimisation algorithm. 1083 @type min_options: array of str 1084 @keyword func_tol: The function tolerance which, when reached, terminates optimisation. 1085 Setting this to None turns of the check. 1086 @type func_tol: None or float 1087 @keyword grad_tol: The gradient tolerance which, when reached, terminates optimisation. 1088 Setting this to None turns of the check. 1089 @type grad_tol: None or float 1090 @keyword max_iterations: The maximum number of iterations for the algorithm. 1091 @type max_iterations: int 1092 @keyword constraints: If True, constraints are used during optimisation. 1093 @type constraints: bool 1094 @keyword scaling: If True, diagonal scaling is enabled during optimisation to allow 1095 the problem to be better conditioned. 1096 @type scaling: bool 1097 @keyword verbosity: The amount of information to print. The higher the value, the 1098 greater the verbosity. 1099 @type verbosity: int 1100 @keyword sim_index: The index of the simulation to optimise. This should be None if 1101 normal optimisation is desired. 1102 @type sim_index: None or int 1103 @keyword lower: The lower bounds of the grid search which must be equal to the 1104 number of parameters in the model. This optional argument is only 1105 used when doing a grid search. 1106 @type lower: array of numbers 1107 @keyword upper: The upper bounds of the grid search which must be equal to the 1108 number of parameters in the model. This optional argument is only 1109 used when doing a grid search. 1110 @type upper: array of numbers 1111 @keyword inc: The increments for each dimension of the space for the grid search. 1112 The number of elements in the array must equal to the number of 1113 parameters in the model. This argument is only used when doing a 1114 grid search. 1115 @type inc: array of int 1116 """ 1117 1118 # Test if sequence data is loaded. 1119 if not exists_mol_res_spin_data(): 1120 raise RelaxNoSequenceError 1121 1122 # Test if the model-free model has been setup, and that the nuclear isotope types have been set. 1123 for spin in spin_loop(): 1124 # Skip deselected spins. 1125 if not spin.select: 1126 continue 1127 1128 # Not setup. 1129 if not spin.model: 1130 raise RelaxNoModelError 1131 1132 # Test if the nuclear isotope type has been set. 1133 if not hasattr(spin, 'isotope'): 1134 raise RelaxSpinTypeError 1135 1136 # Reset the minimisation statistics. 1137 if sim_index == None and min_algor != 'back_calc': 1138 reset_min_stats() 1139 1140 # Containers for the model-free data and optimisation parameters. 1141 data_store = Data_container() 1142 opt_params = Data_container() 1143 1144 # Add the imported parameters. 1145 data_store.h_bar = h_bar 1146 data_store.mu0 = mu0 1147 opt_params.min_algor = min_algor 1148 opt_params.min_options = min_options 1149 opt_params.func_tol = func_tol 1150 opt_params.grad_tol = grad_tol 1151 opt_params.max_iterations = max_iterations 1152 1153 # Add the keyword args. 1154 opt_params.verbosity = verbosity 1155 1156 # Determine the model type. 1157 data_store.model_type = determine_model_type() 1158 if not data_store.model_type: 1159 return 1160 1161 # Model type for the back-calculate function. 1162 if min_algor == 'back_calc' and data_store.model_type != 'local_tm': 1163 data_store.model_type = 'mf' 1164 1165 # Test if diffusion tensor data exists. 1166 if data_store.model_type != 'local_tm' and not diffusion_tensor.diff_data_exists(): 1167 raise RelaxNoTensorError('diffusion') 1168 1169 # Tests for the PDB file and unit vectors. 1170 if data_store.model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere': 1171 # Test if the structure file has been loaded. 1172 if not hasattr(cdp, 'structure'): 1173 raise RelaxNoPdbError 1174 1175 # Test if unit vectors exist. 1176 for spin, spin_id in spin_loop(return_id=True): 1177 # Skip deselected spins. 1178 if not spin.select: 1179 continue 1180 1181 # Get the interatomic data container. 1182 interatoms = return_interatom_list(spin_id) 1183 1184 # Unit vectors. 1185 for i in range(len(interatoms)): 1186 # No relaxation mechanism. 1187 if not interatoms[i].dipole_pair: 1188 continue 1189 1190 # Check for the vectors. 1191 if not hasattr(interatoms[i], 'vector'): 1192 raise RelaxNoVectorsError 1193 1194 # Test if the model-free parameter values are set for minimising diffusion tensor parameters by themselves. 1195 if data_store.model_type == 'diff': 1196 # Loop over the sequence. 1197 for spin in spin_loop(): 1198 unset_param = are_mf_params_set(spin) 1199 if unset_param != None: 1200 raise RelaxNoValueError(unset_param) 1201 1202 # Print out. 1203 if verbosity >= 1: 1204 if data_store.model_type == 'mf': 1205 print("Only the model-free parameters for single spins will be used.") 1206 elif data_store.model_type == 'local_mf': 1207 print("Only a local tm value together with the model-free parameters for single spins will be used.") 1208 elif data_store.model_type == 'diff': 1209 print("Only diffusion tensor parameters will be used.") 1210 elif data_store.model_type == 'all': 1211 print("The diffusion tensor parameters together with the model-free parameters for all spins will be used.") 1212 1213 # Test if the CSA and interatomic distances have been set. 1214 for spin, spin_id in spin_loop(return_id=True): 1215 # Skip deselected spins. 1216 if not spin.select: 1217 continue 1218 1219 # CSA value. 1220 if not hasattr(spin, 'csa') or spin.csa == None: 1221 raise RelaxNoValueError("CSA") 1222 1223 # Get the interatomic data container. 1224 interatoms = return_interatom_list(spin_id) 1225 1226 # Interatomic distances. 1227 count = 0 1228 for i in range(len(interatoms)): 1229 # No relaxation mechanism. 1230 if not interatoms[i].dipole_pair: 1231 continue 1232 1233 # Check for the distances. 1234 if not hasattr(interatoms[i], 'r') or interatoms[i].r == None: 1235 raise RelaxNoValueError("interatomic distance", spin_id=interatoms[i].spin_id1, spin_id2=interatoms[i].spin_id2) 1236 1237 # Count the number of interactions. 1238 count += 1 1239 1240 # Too many interactions. 1241 if count > 1: 1242 raise RelaxError("The spin '%s' has %s dipolar relaxation interactions defined, but only a maximum of one is currently supported." % (spin_id, count)) 1243 1244 # Number of spins, minimisation instances, and data sets for each model type. 1245 if data_store.model_type == 'mf' or data_store.model_type == 'local_tm': 1246 num_instances = count_spins(skip_desel=False) 1247 num_data_sets = 1 1248 data_store.num_spins = 1 1249 elif data_store.model_type == 'diff' or data_store.model_type == 'all': 1250 num_instances = 1 1251 num_data_sets = count_spins(skip_desel=False) 1252 data_store.num_spins = count_spins() 1253 1254 # Number of spins, minimisation instances, and data sets for the back-calculate function. 1255 if min_algor == 'back_calc': 1256 num_instances = 1 1257 num_data_sets = 0 1258 data_store.num_spins = 1 1259 1260 # Get the Processor box singleton (it contains the Processor instance) and alias the Processor. 1261 processor_box = Processor_box() 1262 processor = processor_box.processor 1263 1264 # Loop over the minimisation instances. 1265 ####################################### 1266 1267 for i in range(num_instances): 1268 # Get the spin container if required. 1269 if data_store.model_type == 'diff' or data_store.model_type == 'all': 1270 spin_index = None 1271 spin, data_store.spin_id = None, None 1272 elif min_algor == 'back_calc': 1273 spin_index = opt_params.min_options[0] 1274 spin, data_store.spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 1275 else: 1276 spin_index = i 1277 spin, data_store.spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 1278 1279 # Individual spin stuff. 1280 if spin and (data_store.model_type == 'mf' or data_store.model_type == 'local_tm') and not min_algor == 'back_calc': 1281 # Skip deselected spins. 1282 if not spin.select: 1283 continue 1284 1285 # Skip spins missing relaxation data or errors. 1286 if not hasattr(spin, 'ri_data') or not hasattr(spin, 'ri_data_err'): 1287 continue 1288 1289 # Skip spins missing the dipolar interaction. 1290 if spin and (data_store.model_type == 'mf' or data_store.model_type == 'local_tm'): 1291 interatoms = return_interatom_list(data_store.spin_id) 1292 if not len(interatoms): 1293 continue 1294 1295 # Parameter vector and diagonal scaling. 1296 if min_algor == 'back_calc': 1297 # Create the initial parameter vector. 1298 opt_params.param_vector = assemble_param_vector(spin=spin, model_type=data_store.model_type) 1299 1300 # Diagonal scaling. 1301 data_store.scaling_matrix = None 1302 1303 else: 1304 # Create the initial parameter vector. 1305 opt_params.param_vector = assemble_param_vector(spin=spin, sim_index=sim_index) 1306 1307 # The number of parameters. 1308 num_params = len(opt_params.param_vector) 1309 1310 # Diagonal scaling. 1311 data_store.scaling_matrix = assemble_scaling_matrix(num_params, model_type=data_store.model_type, spin=spin, scaling=scaling) 1312 if len(data_store.scaling_matrix): 1313 opt_params.param_vector = dot(inv(data_store.scaling_matrix), opt_params.param_vector) 1314 1315 # Configure the grid search. 1316 opt_params.inc, opt_params.lower, opt_params.upper = None, None, None 1317 if match('^[Gg]rid', min_algor): 1318 opt_params.inc, opt_params.lower, opt_params.upper = grid_search_config(num_params, spin=spin, lower=lower, upper=upper, inc=inc, scaling_matrix=data_store.scaling_matrix) 1319 1320 # Scaling of values for the set function. 1321 if match('^[Ss]et', min_algor): 1322 opt_params.min_options = dot(inv(data_store.scaling_matrix), opt_params.min_options) 1323 1324 # Linear constraints. 1325 if constraints: 1326 opt_params.A, opt_params.b = linear_constraints(num_params, model_type=data_store.model_type, spin=spin, scaling_matrix=data_store.scaling_matrix) 1327 else: 1328 opt_params.A, opt_params.b = None, None 1329 1330 # Get the data for minimisation. 1331 minimise_data_setup(data_store, min_algor, num_data_sets, opt_params.min_options, spin=spin, sim_index=sim_index) 1332 1333 # Setup the minimisation algorithm when constraints are present. 1334 if constraints and not match('^[Gg]rid', min_algor): 1335 algor = opt_params.min_options[0] 1336 else: 1337 algor = min_algor 1338 1339 # Initialise the function to minimise (for back-calculation and LM minimisation). 1340 if min_algor == 'back_calc' or match('[Ll][Mm]$', algor) or match('[Ll]evenburg-[Mm]arquardt$', algor): 1341 mf = Mf(init_params=opt_params.param_vector, model_type=data_store.model_type, diff_type=data_store.diff_type, diff_params=data_store.diff_params, scaling_matrix=data_store.scaling_matrix, num_spins=data_store.num_spins, equations=data_store.equations, param_types=data_store.param_types, param_values=data_store.param_values, relax_data=data_store.ri_data, errors=data_store.ri_data_err, bond_length=data_store.r, csa=data_store.csa, num_frq=data_store.num_frq, frq=data_store.frq, num_ri=data_store.num_ri, remap_table=data_store.remap_table, noe_r1_table=data_store.noe_r1_table, ri_labels=data_store.ri_types, gx=data_store.gx, gh=data_store.gh, h_bar=data_store.h_bar, mu0=data_store.mu0, num_params=data_store.num_params, vectors=data_store.xh_unit_vectors) 1342 1343 # Levenberg-Marquardt minimisation. 1344 if match('[Ll][Mm]$', algor) or match('[Ll]evenburg-[Mm]arquardt$', algor): 1345 # Total number of ri. 1346 number_ri = 0 1347 for k in range(len(ri_data_err)): 1348 number_ri = number_ri + len(ri_data_err[k]) 1349 1350 # Reconstruct the error data structure. 1351 lm_error = zeros(number_ri, float64) 1352 index = 0 1353 for k in range(len(ri_data_err)): 1354 lm_error[index:index+len(ri_data_err[k])] = ri_data_err[k] 1355 index = index + len(ri_data_err[k]) 1356 1357 opt_params.min_options = opt_params.min_options + (mf.lm_dri, lm_error) 1358 1359 # Back-calculation. 1360 if min_algor == 'back_calc': 1361 return mf.calc_ri() 1362 1363 # Parallelised grid search for the diffusion parameter space. 1364 if match('^[Gg]rid', min_algor) and data_store.model_type == 'diff': 1365 # Print out. 1366 print("Parallelised diffusion tensor grid search.") 1367 1368 # Loop over each grid subdivision. 1369 for subdivision in grid_split(divisions=processor.processor_size(), lower=opt_params.lower, upper=opt_params.upper, inc=opt_params.inc): 1370 # Set the points. 1371 opt_params.subdivision = subdivision 1372 1373 # Grid search initialisation. 1374 command = MF_grid_command() 1375 1376 # Pass in the data and optimisation parameters. 1377 command.store_data(deepcopy(data_store), deepcopy(opt_params)) 1378 1379 # Set up the model-free memo and add it to the processor queue. 1380 memo = MF_memo(model_free=self, model_type=data_store.model_type, spin=spin, sim_index=sim_index, scaling=scaling, scaling_matrix=data_store.scaling_matrix) 1381 processor.add_to_queue(command, memo) 1382 1383 # Execute the queued elements. 1384 processor.run_queue() 1385 1386 # Exit this method. 1387 return 1388 1389 # Normal grid search (command initialisation). 1390 if search('^[Gg]rid', min_algor): 1391 command = MF_grid_command() 1392 1393 # Minimisation of all other model types (command initialisation). 1394 else: 1395 command = MF_minimise_command() 1396 1397 # Pass in the data and optimisation parameters. 1398 command.store_data(deepcopy(data_store), deepcopy(opt_params)) 1399 1400 # Set up the model-free memo and add it to the processor queue. 1401 memo = MF_memo(model_free=self, model_type=data_store.model_type, spin=spin, sim_index=sim_index, scaling=scaling, scaling_matrix=data_store.scaling_matrix) 1402 processor.add_to_queue(command, memo) 1403 1404 # Execute the queued elements. 1405 processor.run_queue()
1406 1407
1408 - def model_desc(self, model_info):
1409 """Return a description of the model. 1410 1411 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1412 @type model_info: int 1413 @return: The model description. 1414 @rtype: str 1415 """ 1416 1417 # Determine the model type. 1418 model_type = determine_model_type() 1419 1420 # Global models. 1421 if model_type == 'all': 1422 return "Global model - all diffusion tensor parameters and spin specific model-free parameters." 1423 elif model_type == 'diff': 1424 return "Diffusion tensor model." 1425 1426 # Spin specific model. 1427 else: 1428 # Get the spin container. 1429 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 1430 1431 # Return the description. 1432 return "Model-free model of spin '%s'." % spin_id
1433 1434
1435 - def model_loop(self):
1436 """Generator method for looping over the models (global or local). 1437 1438 If the model type is 'all' or 'diff', then this yields the single value of zero. Otherwise 1439 the global spin index is yielded. 1440 1441 1442 @return: The model index. This is zero for the global models or equal to the global spin 1443 index (which covers the molecule, residue, and spin indices). 1444 @rtype: int 1445 """ 1446 1447 # Determine the model type. 1448 model_type = determine_model_type() 1449 1450 # Global model. 1451 if model_type == 'all' or model_type == 'diff': 1452 yield 0 1453 1454 # Spin specific models. 1455 else: 1456 # Loop over the spins. 1457 global_index = -1 1458 for spin in spin_loop(): 1459 # Increment the global spin index. 1460 global_index = global_index + 1 1461 1462 # Yield the spin index. 1463 yield global_index
1464 1465
1466 - def model_statistics(self, model_info=None, spin_id=None, global_stats=None):
1467 """Return the k, n, and chi2 model statistics. 1468 1469 k - number of parameters. 1470 n - number of data points. 1471 chi2 - the chi-squared value. 1472 1473 1474 @keyword model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1475 @type model_info: int 1476 @keyword spin_id: The spin identification string. Either this or the instance keyword argument must be supplied. 1477 @type spin_id: None or str 1478 @keyword global_stats: A parameter which determines if global or local statistics are returned. If None, then the appropriateness of global or local statistics is automatically determined. 1479 @type global_stats: None or bool 1480 @return: The optimisation statistics, in tuple format, of the number of parameters (k), the number of data points (n), and the chi-squared value (chi2). 1481 @rtype: tuple of (int, int, float) 1482 """ 1483 1484 # Bad argument combination. 1485 if model_info == None and spin_id == None: 1486 raise RelaxError("Either the model_info or spin_id argument must be supplied.") 1487 elif model_info != None and spin_id != None: 1488 raise RelaxError("The model_info arg " + repr(model_info) + " and spin_id arg " + repr(spin_id) + " clash. Only one should be supplied.") 1489 1490 # Determine the model type. 1491 model_type = determine_model_type() 1492 1493 # Determine if local or global statistics will be returned. 1494 if global_stats == None: 1495 if model_type in ['mf', 'local_tm']: 1496 global_stats = False 1497 else: 1498 global_stats = True 1499 1500 # Statistics for a single residue. 1501 if not global_stats: 1502 # Get the SpinContainer. 1503 if spin_id: 1504 spin = return_spin(spin_id) 1505 else: 1506 spin = return_spin_from_index(model_info) 1507 1508 # Skip deselected residues. 1509 if not spin.select: 1510 return None, None, None 1511 1512 # Missing data sets. 1513 if not hasattr(spin, 'ri_data'): 1514 return None, None, None 1515 1516 # Count the number of parameters. 1517 param_vector = assemble_param_vector(spin=spin) 1518 k = len(param_vector) 1519 1520 # Count the number of data points. 1521 n = len(spin.ri_data) 1522 1523 # The chi2 value. 1524 chi2 = spin.chi2 1525 1526 # Global stats. 1527 elif global_stats: 1528 # Count the number of parameters. 1529 param_vector = assemble_param_vector() 1530 k = len(param_vector) 1531 1532 # Count the number of data points. 1533 n = 0 1534 chi2 = 0 1535 for spin in spin_loop(): 1536 # Skip deselected residues. 1537 if not spin.select: 1538 continue 1539 1540 # Skip residues with no relaxation data. 1541 if not hasattr(spin, 'ri_data') or not len(spin.ri_data): 1542 continue 1543 1544 n = n + len(spin.ri_data) 1545 1546 # Local tm models. 1547 if model_type == 'local_tm': 1548 chi2 = chi2 + spin.chi2 1549 1550 # The chi2 value. 1551 if model_type != 'local_tm': 1552 if not hasattr(cdp, 'chi2'): 1553 raise RelaxError("Global statistics are not available, most likely because the global model has not been optimised.") 1554 chi2 = cdp.chi2 1555 1556 # Return the data. 1557 return k, n, chi2
1558 1559
1560 - def model_type(self):
1561 """Return the type of the model, either being 'local' or 'global'. 1562 1563 @return: The model type, one of 'local' or 'global'. 1564 @rtype: str 1565 """ 1566 1567 # Determine the model type. 1568 model_type = determine_model_type() 1569 1570 # Global models. 1571 if model_type in ['all', 'diff']: 1572 return 'global' 1573 1574 # Local models. 1575 else: 1576 return 'local'
1577 1578
1579 - def num_instances(self):
1580 """Function for returning the number of instances. 1581 1582 @return: The number of instances used for optimisation. Either the number of spins if 1583 the local optimisations are setup ('mf' and 'local_tm'), or 1 for the global 1584 models. 1585 @rtype: int 1586 """ 1587 1588 # Test if sequence data exists. 1589 if not exists_mol_res_spin_data(): 1590 return 0 1591 1592 # Determine the model type. 1593 model_type = determine_model_type() 1594 1595 # Sequence specific data. 1596 if model_type == 'mf' or model_type == 'local_tm': 1597 return count_spins() 1598 1599 # Other data. 1600 elif model_type == 'diff' or model_type == 'all': 1601 return 1 1602 1603 # Should not be here. 1604 else: 1605 raise RelaxFault
1606 1607
1608 - def overfit_deselect(self, data_check=True, verbose=True):
1609 """Deselect spins which have insufficient data to support minimisation. 1610 1611 @keyword data_check: A flag to signal if the presence of base data is to be checked for. 1612 @type data_check: bool 1613 @keyword verbose: A flag which if True will allow printouts. 1614 @type verbose: bool 1615 """ 1616 1617 # Print out. 1618 if verbose: 1619 print("\nOver-fit spin deselection:") 1620 1621 # Test if sequence data exists. 1622 if not exists_mol_res_spin_data(): 1623 raise RelaxNoSequenceError 1624 1625 # Is structural data required? 1626 need_vect = False 1627 if hasattr(cdp, 'diff_tensor') and (cdp.diff_tensor.type == 'spheroid' or cdp.diff_tensor.type == 'ellipsoid'): 1628 need_vect = True 1629 1630 # Loop over the sequence. 1631 deselect_flag = False 1632 spin_count = 0 1633 for spin, spin_id in spin_loop(return_id=True): 1634 # Skip deselected spins. 1635 if not spin.select: 1636 continue 1637 1638 # The interatomic data. 1639 interatoms = interatomic.return_interatom_list(spin_id) 1640 1641 # Loop over the interatomic data. 1642 dipole_relax = False 1643 for i in range(len(interatoms)): 1644 # No dipolar relaxation mechanism. 1645 if not interatoms[i].dipole_pair: 1646 continue 1647 1648 # The surrounding spins. 1649 if spin_id != interatoms[i].spin_id1: 1650 spin_id2 = interatoms[i].spin_id1 1651 else: 1652 spin_id2 = interatoms[i].spin_id2 1653 spin2 = return_spin(spin_id2) 1654 1655 # Dipolar relaxation flag. 1656 dipole_relax = True 1657 1658 # No relaxation mechanism. 1659 if not dipole_relax or not hasattr(spin, 'csa') or spin.csa == None: 1660 warn(RelaxDeselectWarning(spin_id, 'an absence of relaxation mechanisms')) 1661 spin.select = False 1662 deselect_flag = True 1663 continue 1664 1665 # Data checks. 1666 if data_check: 1667 # The number of relaxation data points (and for infinite data). 1668 data_points = 0 1669 inf_data = False 1670 if hasattr(cdp, 'ri_ids') and hasattr(spin, 'ri_data'): 1671 for id in cdp.ri_ids: 1672 if id in spin.ri_data and spin.ri_data[id] != None: 1673 data_points += 1 1674 1675 # Infinite data! 1676 if isInf(spin.ri_data[id]): 1677 inf_data = True 1678 1679 # Infinite data. 1680 if inf_data: 1681 warn(RelaxDeselectWarning(spin_id, 'infinite relaxation data')) 1682 spin.select = False 1683 deselect_flag = True 1684 continue 1685 1686 # Relaxation data must exist! 1687 if not hasattr(spin, 'ri_data'): 1688 warn(RelaxDeselectWarning(spin_id, 'missing relaxation data')) 1689 spin.select = False 1690 deselect_flag = True 1691 continue 1692 1693 # Require 3 or more relaxation data points. 1694 elif data_points < 3: 1695 warn(RelaxDeselectWarning(spin_id, 'insufficient relaxation data, 3 or more data points are required')) 1696 spin.select = False 1697 deselect_flag = True 1698 continue 1699 1700 # Require at least as many data points as params to prevent over-fitting. 1701 elif hasattr(spin, 'params') and spin.params and len(spin.params) > data_points: 1702 warn(RelaxDeselectWarning(spin_id, 'over-fitting - more parameters than data points')) 1703 spin.select = False 1704 deselect_flag = True 1705 continue 1706 1707 # Test for structural data if required. 1708 for i in range(len(interatoms)): 1709 # No dipolar relaxation mechanism. 1710 if not interatoms[i].dipole_pair: 1711 continue 1712 1713 # Check the unit vectors. 1714 if need_vect: 1715 if not hasattr(interatoms[i], 'vector') or interatoms[i].vector == None: 1716 warn(RelaxDeselectWarning(spin_id, 'missing structural data')) 1717 spin.select = False 1718 deselect_flag = True 1719 continue 1720 1721 # Increment the spin number. 1722 spin_count += 1 1723 1724 # No spins selected, so fail hard to prevent the user from going any further. 1725 if spin_count == 0: 1726 warn(RelaxWarning("No spins are selected therefore the optimisation or calculation cannot proceed.")) 1727 1728 # Final printout. 1729 if verbose and not deselect_flag: 1730 print("No spins have been deselected.")
1731 1732
1733 - def set_error(self, model_info, index, error):
1734 """Set the parameter errors. 1735 1736 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1737 @type model_info: int 1738 @param index: The index of the parameter to set the errors for. 1739 @type index: int 1740 @param error: The error value. 1741 @type error: float 1742 """ 1743 1744 # Parameter increment counter. 1745 inc = 0 1746 1747 # Determine the model type. 1748 model_type = determine_model_type() 1749 1750 # Get the parameter object names. 1751 param_names = self.data_names(set='params', scope='spin') 1752 1753 1754 # Diffusion tensor parameter errors. 1755 #################################### 1756 1757 if model_type == 'diff' or model_type == 'all': 1758 # Spherical diffusion. 1759 if cdp.diff_tensor.type == 'sphere': 1760 # Return the parameter array. 1761 if index == 0: 1762 cdp.diff_tensor.set(param='tm', value=error, category='err') 1763 1764 # Increment. 1765 inc = inc + 1 1766 1767 # Spheroidal diffusion. 1768 elif cdp.diff_tensor.type == 'spheroid': 1769 # Return the parameter array. 1770 if index == 0: 1771 cdp.diff_tensor.set(param='tm', value=error, category='err') 1772 elif index == 1: 1773 cdp.diff_tensor.set(param='Da', value=error, category='err') 1774 elif index == 2: 1775 cdp.diff_tensor.set(param='theta', value=error, category='err') 1776 elif index == 3: 1777 cdp.diff_tensor.set(param='phi', value=error, category='err') 1778 1779 # Increment. 1780 inc = inc + 4 1781 1782 # Ellipsoidal diffusion. 1783 elif cdp.diff_tensor.type == 'ellipsoid': 1784 # Return the parameter array. 1785 if index == 0: 1786 cdp.diff_tensor.set(param='tm', value=error, category='err') 1787 elif index == 1: 1788 cdp.diff_tensor.set(param='Da', value=error, category='err') 1789 elif index == 2: 1790 cdp.diff_tensor.set(param='Dr', value=error, category='err') 1791 elif index == 3: 1792 cdp.diff_tensor.set(param='alpha', value=error, category='err') 1793 elif index == 4: 1794 cdp.diff_tensor.set(param='beta', value=error, category='err') 1795 elif index == 5: 1796 cdp.diff_tensor.set(param='gamma', value=error, category='err') 1797 1798 # Increment. 1799 inc = inc + 6 1800 1801 1802 # Model-free parameter errors for the model type 'all'. 1803 ####################################################### 1804 1805 if model_type == 'all': 1806 # Loop over the spins. 1807 for spin in spin_loop(): 1808 # Skip deselected spins. 1809 if not spin.select: 1810 continue 1811 1812 # Loop over the residue specific parameters. 1813 for param in param_names: 1814 # Return the parameter array. 1815 if index == inc: 1816 setattr(spin, param + "_err", error) 1817 1818 # Increment. 1819 inc = inc + 1 1820 1821 1822 # Model-free parameters for the model types 'mf' and 'local_tm'. 1823 ################################################################ 1824 1825 if model_type == 'mf' or model_type == 'local_tm': 1826 # Get the spin container. 1827 spin = return_spin_from_index(model_info) 1828 1829 # Skip deselected residues. 1830 if not spin.select: 1831 return 1832 1833 # Loop over the residue specific parameters. 1834 for param in param_names: 1835 # Return the parameter array. 1836 if index == inc: 1837 setattr(spin, param + "_err", error) 1838 1839 # Increment. 1840 inc = inc + 1
1841 1842
1843 - def set_param_values(self, param=None, value=None, index=None, spin_id=None, error=False, force=True):
1844 """Set the model-free parameter values. 1845 1846 @keyword param: The parameter name list. 1847 @type param: list of str 1848 @keyword value: The parameter value list. 1849 @type value: list 1850 @keyword index: The index for parameters which are of the list-type. This is unused. 1851 @type index: None or int 1852 @keyword spin_id: The spin identification string, only used for spin specific parameters. 1853 @type spin_id: None or str 1854 @keyword error: A flag which if True will allow the parameter errors to be set instead of the values. 1855 @type error: bool 1856 @keyword force: A flag which if True will cause current values to be overwritten. If False, a RelaxError will raised if the parameter value is already set. 1857 @type force: bool 1858 """ 1859 1860 # Checks. 1861 is_str_list(param, 'parameter name') 1862 1863 # Separate out the diffusion tensor parameters from the model-free parameters. 1864 diff_params = [] 1865 diff_vals = [] 1866 mf_params = [] 1867 mf_vals = [] 1868 for i in range(len(param)): 1869 # Diffusion tensor parameter. 1870 if self._PARAMS.scope(param[i]) == 'global': 1871 if error: 1872 diff_params.append(param[i] + '_err') 1873 else: 1874 diff_params.append(param[i]) 1875 diff_vals.append(value[i]) 1876 1877 # Model-free parameter. 1878 else: 1879 mf_params.append(param[i]) 1880 mf_vals.append(value[i]) 1881 1882 # Set the diffusion tensor parameters. 1883 if diff_params: 1884 diffusion_tensor.set(value=diff_vals, param=diff_params) 1885 1886 # Set the model-free parameters. 1887 for i in range(len(mf_params)): 1888 # Check if it is a model-free parameter. 1889 if mf_params[i] not in self.data_names(set='params', scope='spin') and mf_params[i] not in self.data_names(set='generic', scope='spin'): 1890 raise RelaxError("The parameter '%s' is unknown. It should be one of %s or %s" % (mf_params[i], self.data_names(set='params', scope='spin'), self.data_names(set='generic', scope='spin'))) 1891 1892 # The error object name. 1893 if error: 1894 mf_params[i] += '_err' 1895 1896 # Set the parameter. 1897 for spin in spin_loop(spin_id): 1898 setattr(spin, mf_params[i], mf_vals[i])
1899 1900
1901 - def set_selected_sim(self, model_info, select_sim):
1902 """Set all simulation selection flags. 1903 1904 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1905 @type model_info: int 1906 @param select_sim: The selection flags. 1907 @type select_sim: bool 1908 """ 1909 1910 # Determine the model type. 1911 model_type = determine_model_type() 1912 1913 # Global model. 1914 if model_type == 'all' or model_type == 'diff': 1915 cdp.select_sim = select_sim 1916 1917 # Spin specific model. 1918 else: 1919 # Get the spin container. 1920 spin = return_spin_from_index(model_info) 1921 1922 # Skip if deselected. 1923 if not spin.select: 1924 return 1925 1926 # Set the simulation flags. 1927 spin.select_sim = deepcopy(select_sim)
1928 1929
1930 - def set_update(self, param, spin):
1931 """Function to update the other model-free parameters. 1932 1933 @param param: The name of the parameter which has been changed. 1934 @type param: str 1935 @param spin: The SpinContainer object. 1936 @type spin: SpinContainer 1937 """ 1938 1939 # S2f parameter. 1940 if param == 's2f': 1941 # Update S2 if S2s exists. 1942 if hasattr(spin, 's2s') and spin.s2s != None: 1943 spin.s2 = spin.s2f * spin.s2s 1944 1945 1946 # S2s parameter. 1947 if param == 's2s': 1948 # Update S2 if S2f exists. 1949 if hasattr(spin, 's2f') and spin.s2f != None: 1950 spin.s2 = spin.s2f * spin.s2s
1951 1952
1953 - def sim_init_values(self):
1954 """Initialise the Monte Carlo parameter values.""" 1955 1956 # Determine the model type. 1957 model_type = determine_model_type() 1958 1959 # Get the parameter object names. 1960 param_names = self.data_names(set='params', scope='spin') 1961 1962 # Get the minimisation statistic object names. 1963 min_names = self.data_names(set='min', scope='spin') 1964 1965 # List of diffusion tensor parameters. 1966 if model_type == 'diff' or model_type == 'all': 1967 # Spherical diffusion. 1968 if cdp.diff_tensor.type == 'sphere': 1969 diff_params = ['tm'] 1970 1971 # Spheroidal diffusion. 1972 elif cdp.diff_tensor.type == 'spheroid': 1973 diff_params = ['tm', 'Da', 'theta', 'phi'] 1974 1975 # Ellipsoidal diffusion. 1976 elif cdp.diff_tensor.type == 'ellipsoid': 1977 diff_params = ['tm', 'Da', 'Dr', 'alpha', 'beta', 'gamma'] 1978 1979 1980 # Test if Monte Carlo parameter values have already been set. 1981 ############################################################# 1982 1983 # Diffusion tensor parameters and non spin specific minimisation statistics. 1984 if model_type == 'diff' or model_type == 'all': 1985 # Loop over the parameters. 1986 for object_name in diff_params: 1987 # Name for the simulation object. 1988 sim_object_name = object_name + '_sim' 1989 1990 # Test if the simulation object already exists. 1991 if hasattr(cdp.diff_tensor, sim_object_name): 1992 raise RelaxError("Monte Carlo parameter values have already been set.") 1993 1994 # Loop over the minimisation stats objects. 1995 for object_name in min_names: 1996 # Name for the simulation object. 1997 sim_object_name = object_name + '_sim' 1998 1999 # Test if the simulation object already exists. 2000 if hasattr(cdp, sim_object_name): 2001 raise RelaxError("Monte Carlo parameter values have already been set.") 2002 2003 # Spin specific parameters. 2004 if model_type != 'diff': 2005 for spin in spin_loop(): 2006 # Skip deselected spins. 2007 if not spin.select: 2008 continue 2009 2010 # Loop over all the parameter names. 2011 for object_name in param_names: 2012 # Name for the simulation object. 2013 sim_object_name = object_name + '_sim' 2014 2015 # Test if the simulation object already exists. 2016 if hasattr(spin, sim_object_name): 2017 raise RelaxError("Monte Carlo parameter values have already been set.") 2018 2019 2020 # Set the Monte Carlo parameter values. 2021 ####################################### 2022 2023 # Loop over the global minimisation stats objects. 2024 for object_name in min_names: 2025 # Skip non-existent objects. 2026 if not hasattr(cdp, object_name): 2027 continue 2028 2029 # Name for the simulation object. 2030 sim_object_name = object_name + '_sim' 2031 2032 # Create the simulation object. 2033 setattr(cdp, sim_object_name, []) 2034 2035 # Get the simulation object. 2036 sim_object = getattr(cdp, sim_object_name) 2037 2038 # Loop over the simulations. 2039 for j in range(cdp.sim_number): 2040 # Get the object. 2041 object = getattr(cdp, object_name) 2042 2043 # Copy and append the data. 2044 sim_object.append(deepcopy(object)) 2045 2046 # Diffusion tensor parameters and non spin specific minimisation statistics. 2047 if model_type == 'diff' or model_type == 'all': 2048 # Set up the number of simulations. 2049 cdp.diff_tensor.set_sim_num(cdp.sim_number) 2050 2051 # Loop over the parameters, setting the initial simulation values to those of the parameter value. 2052 for object_name in diff_params: 2053 for j in range(cdp.sim_number): 2054 cdp.diff_tensor.set(param=object_name, value=deepcopy(getattr(cdp.diff_tensor, object_name)), category='sim', sim_index=j) 2055 2056 # Spin specific parameters. 2057 if model_type != 'diff': 2058 for spin in spin_loop(): 2059 # Skip deselected spins. 2060 if not spin.select: 2061 continue 2062 2063 # Loop over all the data names. 2064 for object_name in param_names: 2065 # Name for the simulation object. 2066 sim_object_name = object_name + '_sim' 2067 2068 # Create the simulation object. 2069 setattr(spin, sim_object_name, []) 2070 2071 # Get the simulation object. 2072 sim_object = getattr(spin, sim_object_name) 2073 2074 # Loop over the simulations. 2075 for j in range(cdp.sim_number): 2076 # Copy and append the data. 2077 sim_object.append(deepcopy(getattr(spin, object_name))) 2078 2079 # Loop over all the minimisation object names. 2080 for object_name in min_names: 2081 # Name for the simulation object. 2082 sim_object_name = object_name + '_sim' 2083 2084 # Create the simulation object. 2085 setattr(spin, sim_object_name, []) 2086 2087 # Get the simulation object. 2088 sim_object = getattr(spin, sim_object_name) 2089 2090 # Loop over the simulations. 2091 for j in range(cdp.sim_number): 2092 # Copy and append the data. 2093 sim_object.append(deepcopy(getattr(spin, object_name)))
2094 2095
2096 - def sim_return_chi2(self, model_info, index=None):
2097 """Return the simulation chi-squared values. 2098 2099 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2100 @type model_info: int 2101 @keyword index: The optional simulation index. 2102 @type index: int 2103 @return: The list of simulation chi-squared values. If the index is supplied, only a single value will be returned. 2104 @rtype: list of float or float 2105 """ 2106 2107 # Determine the model type. 2108 model_type = determine_model_type() 2109 2110 # Single instance. 2111 if model_type == 'all' or model_type == 'diff': 2112 return cdp.chi2_sim 2113 2114 # Multiple instances. 2115 else: 2116 # Get the spin container. 2117 spin = return_spin_from_index(model_info) 2118 2119 # Return the list. 2120 return spin.chi2_sim
2121 2122
2123 - def sim_return_param(self, model_info, index):
2124 """Return the array of simulation parameter values. 2125 2126 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2127 @type model_info: int 2128 @param index: The index of the parameter to return the array of values for. 2129 @type index: int 2130 @return: The array of simulation parameter values. 2131 @rtype: list of float 2132 """ 2133 2134 # Parameter increment counter. 2135 inc = 0 2136 2137 # Determine the model type. 2138 model_type = determine_model_type() 2139 2140 # Get the parameter object names. 2141 param_names = self.data_names(set='params', scope='spin') 2142 2143 2144 # Diffusion tensor parameters. 2145 ############################## 2146 2147 if model_type == 'diff' or model_type == 'all': 2148 # Spherical diffusion. 2149 if cdp.diff_tensor.type == 'sphere': 2150 # Return the parameter array. 2151 if index == 0: 2152 return cdp.diff_tensor.tm_sim 2153 2154 # Increment. 2155 inc = inc + 1 2156 2157 # Spheroidal diffusion. 2158 elif cdp.diff_tensor.type == 'spheroid': 2159 # Return the parameter array. 2160 if index == 0: 2161 return cdp.diff_tensor.tm_sim 2162 elif index == 1: 2163 return cdp.diff_tensor.Da_sim 2164 elif index == 2: 2165 return cdp.diff_tensor.theta_sim 2166 elif index == 3: 2167 return cdp.diff_tensor.phi_sim 2168 2169 # Increment. 2170 inc = inc + 4 2171 2172 # Ellipsoidal diffusion. 2173 elif cdp.diff_tensor.type == 'ellipsoid': 2174 # Return the parameter array. 2175 if index == 0: 2176 return cdp.diff_tensor.tm_sim 2177 elif index == 1: 2178 return cdp.diff_tensor.Da_sim 2179 elif index == 2: 2180 return cdp.diff_tensor.Dr_sim 2181 elif index == 3: 2182 return cdp.diff_tensor.alpha_sim 2183 elif index == 4: 2184 return cdp.diff_tensor.beta_sim 2185 elif index == 5: 2186 return cdp.diff_tensor.gamma_sim 2187 2188 # Increment. 2189 inc = inc + 6 2190 2191 2192 # Model-free parameters for the model type 'all'. 2193 ################################################# 2194 2195 if model_type == 'all': 2196 # Loop over the spins. 2197 for spin in spin_loop(): 2198 # Skip deselected spins. 2199 if not spin.select: 2200 continue 2201 2202 # Loop over the spin specific parameters. 2203 for param in param_names: 2204 # Return the parameter array. 2205 if index == inc: 2206 return getattr(spin, param + "_sim") 2207 2208 # Increment. 2209 inc = inc + 1 2210 2211 2212 # Model-free parameters for the model types 'mf' and 'local_tm'. 2213 ################################################################ 2214 2215 if model_type == 'mf' or model_type == 'local_tm': 2216 # Get the spin container. 2217 spin = return_spin_from_index(model_info) 2218 2219 # Skip deselected spins. 2220 if not spin.select: 2221 return 2222 2223 # Loop over the spin specific parameters. 2224 for param in param_names: 2225 # Return the parameter array. 2226 if index == inc: 2227 return getattr(spin, param + "_sim") 2228 2229 # Increment. 2230 inc = inc + 1
2231 2232
2233 - def sim_return_selected(self, model_info):
2234 """Return the array of selected simulation flags for the spin. 2235 2236 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2237 @type model_info: int 2238 @return: The array of selected simulation flags. 2239 @rtype: list of int 2240 """ 2241 2242 # Determine the model type. 2243 model_type = determine_model_type() 2244 2245 # Single instance. 2246 if model_type == 'all' or model_type == 'diff': 2247 return cdp.select_sim 2248 2249 # Multiple instances. 2250 else: 2251 # Get the spin container. 2252 spin = return_spin_from_index(model_info) 2253 2254 # Skip if deselected. 2255 if not spin.select: 2256 return 2257 2258 # Return the list. 2259 return spin.select_sim
2260 2261
2262 - def skip_function(self, model_info):
2263 """Skip certain data. 2264 2265 @param model_info: The model index from model_loop(). This is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2266 @type model_info: int 2267 @return: True if the data should be skipped, False otherwise. 2268 @rtype: bool 2269 """ 2270 2271 # Determine the model type. 2272 model_type = determine_model_type() 2273 2274 # Sequence specific data. 2275 if (model_type == 'mf' or model_type == 'local_tm') and not return_spin_from_index(model_info).select: 2276 return True 2277 2278 # Don't skip. 2279 return False
2280 2281 2282
2283 -class Data_container:
2284 """Empty class to be used for data storage."""
2285