specific_analyses.model

71 """Parent class containing all the model-free specific functions.""" 72 73 # Class variable for storing the class instance (for the singleton design pattern). 74 instance = None 75

76 - def __init__(self):

77 """Initialise the class by placing API_common methods into the API.""" 78 79 # Place methods into the API. 80 self.base_data_loop = self._base_data_loop_spin 81 self.return_error = self._return_error_relax_data 82 self.return_value = self._return_value_general 83 self.sim_pack_data = self._sim_pack_relax_data 84 85 # Initialise the macro classes. 86 self._molmol_macros = Molmol() 87 self._pymol_macros = Pymol() 88 89 # Alias the macro creation methods. 90 self.pymol_macro = self._pymol_macros.create_macro 91 self.molmol_macro = self._molmol_macros.create_macro 92 93 # Place a copy of the parameter list object in the instance namespace. 94 self._PARAMS = Model_free_params()

95 96

97 - def back_calc_ri(self, spin_index=None, ri_id=None, ri_type=None, frq=None):

98 """Back-calculation of relaxation data from the model-free parameter values. 99 100 @keyword spin_index: The global spin index. 101 @type spin_index: int 102 @keyword ri_id: The relaxation data ID string. 103 @type ri_id: str 104 @keyword ri_type: The relaxation data type. 105 @type ri_type: str 106 @keyword frq: The field strength. 107 @type frq: float 108 @return: The back calculated relaxation data value corresponding to the index. 109 @rtype: float 110 """ 111 112 # Get the spin container. 113 spin, spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 114 115 # Missing interatomic vectors. 116 if hasattr(cdp, 'diff_tensor') and (cdp.diff_tensor.type == 'spheroid' or cdp.diff_tensor.type == 'ellipsoid'): 117 interatoms = interatomic.return_interatom_list(spin_hash=spin._hash) 118 for i in range(len(interatoms)): 119 # No dipolar relaxation mechanism. 120 if not interatoms[i].dipole_pair: 121 continue 122 123 # Check the unit vectors. 124 if not hasattr(interatoms[i], 'vector') or interatoms[i].vector is None: 125 warn(RelaxDeselectWarning(spin_id, 'missing structural data')) 126 return 127 128 # Execute the over-fit deselection. 129 self.overfit_deselect(data_check=False, verbose=False) 130 131 # Get the relaxation value from the minimise function. 132 value = self.minimise(min_algor='back_calc', min_options=(spin_index, ri_id, ri_type, frq)) 133 134 # Return the relaxation value. 135 return value

136 137

138 - def bmrb_read(self, file_path, version=None, sample_conditions=None):

139 """Read the model-free results from a BMRB NMR-STAR v3.1 formatted file. 140 141 @param file_path: The full file path. 142 @type file_path: str 143 @keyword version: The BMRB version to force the reading. 144 @type version: None or str 145 @keyword sample_conditions: The sample condition label to read. Only one sample condition can be read per data pipe. 146 @type sample_conditions: None or str 147 """ 148 149 # Initialise the NMR-STAR data object. 150 star = bmrblib.create_nmr_star('relax_model_free_results', file_path, version) 151 152 # Read the contents of the STAR formatted file. 153 star.read() 154 155 # The sample conditions. 156 sample_conds = list_sample_conditions(star) 157 if sample_conditions and sample_conditions not in sample_conds: 158 raise RelaxError("The sample conditions label '%s' does not correspond to any of the labels in the file: %s" % (sample_conditions, sample_conds)) 159 if not sample_conditions and len(sample_conds) > 1: 160 raise RelaxError("Only one of the sample conditions in %s can be loaded per relax data pipe." % sample_conds) 161 162 # The diffusion tensor. 163 diffusion_tensor.bmrb_read(star) 164 165 # Generate the molecule and residue containers from the entity records. 166 mol_res_spin.bmrb_read(star) 167 168 # Read the relaxation data saveframes. 169 relax_data.bmrb_read(star, sample_conditions=sample_conditions) 170 171 # Read the model-free data saveframes. 172 sf_model_free_read(star, sample_conditions=sample_conditions) 173 174 # Read the CSA data saveframes. 175 sf_csa_read(star)

176 177

178 - def bmrb_write(self, file_path, version=None):

179 """Write the model-free results to a BMRB NMR-STAR v3.1 formatted file. 180 181 @param file_path: The full file path. 182 @type file_path: str 183 @keyword version: The BMRB NMR-STAR dictionary format to output to. 184 @type version: str 185 """ 186 187 # Checks. 188 check_spectrometer_setup(escalate=2) 189 190 # Alias the current data pipe. 191 cdp = pipes.get_pipe() 192 193 # Initialise the NMR-STAR data object. 194 star = bmrblib.create_nmr_star('relax_model_free_results', file_path, version) 195 196 # Global minimisation stats. 197 global_chi2 = None 198 if hasattr(cdp, 'chi2'): 199 global_chi2 = cdp.chi2 200 201 # Rex frq. 202 rex_frq = cdp.spectrometer_frq[cdp.ri_ids[0]] 203 204 # Initialise the spin specific data lists. 205 mol_name_list = [] 206 res_num_list = [] 207 res_name_list = [] 208 atom_name_list = [] 209 210 csa_list = [] 211 r_list = [] 212 isotope_list = [] 213 element_list = [] 214 215 local_tm_list = [] 216 s2_list = [] 217 s2f_list = [] 218 s2s_list = [] 219 te_list = [] 220 tf_list = [] 221 ts_list = [] 222 rex_list = [] 223 224 local_tm_err_list = [] 225 s2_err_list = [] 226 s2f_err_list = [] 227 s2s_err_list = [] 228 te_err_list = [] 229 tf_err_list = [] 230 ts_err_list = [] 231 rex_err_list = [] 232 233 chi2_list = [] 234 model_list = [] 235 236 # Store the spin specific data in lists for later use. 237 for spin, mol_name, res_num, res_name, spin_id in spin_loop(full_info=True, skip_desel=True, return_id=True): 238 # Skip the protons. 239 if spin.name == 'H' or (hasattr(spin, 'element') and spin.element == 'H'): 240 warn(RelaxWarning("Skipping the proton spin '%s'." % spin_id)) 241 continue 242 243 # No model setup. 244 if not hasattr(spin, 'model'): 245 warn(RelaxWarning("Skipping the spin '%s' as no model-free model has been set up." % spin_id)) 246 continue 247 248 # Check the data for None (not allowed in BMRB!). 249 if res_num == None: 250 raise RelaxError("For the BMRB, the residue of spin '%s' must be numbered." % spin_id) 251 if res_name == None: 252 raise RelaxError("For the BMRB, the residue of spin '%s' must be named." % spin_id) 253 if spin.name == None: 254 raise RelaxError("For the BMRB, the spin '%s' must be named." % spin_id) 255 if not hasattr(spin, 'isotope') or spin.isotope == None: 256 raise RelaxError("For the BMRB, the spin isotope type of '%s' must be specified." % spin_id) 257 if not hasattr(spin, 'element') or spin.element == None: 258 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) 259 260 # The molecule/residue/spin info. 261 mol_name_list.append(mol_name) 262 res_num_list.append(res_num) 263 res_name_list.append(res_name) 264 atom_name_list.append(spin.name) 265 266 # CSA values. 267 if hasattr(spin, 'csa'): 268 csa_list.append(spin.csa * 1e6) # In ppm. 269 else: 270 csa_list.append(None) 271 272 # Interatomic distances. 273 interatoms = return_interatom_list(spin_hash=spin._hash) 274 for i in range(len(interatoms)): 275 # No relaxation mechanism. 276 if not interatoms[i].dipole_pair: 277 continue 278 279 # Add the interatomic distance. 280 if hasattr(interatoms[i], 'r'): 281 r_list.append(interatoms[i].r) 282 else: 283 r_list.append(None) 284 285 # Stop adding. 286 break 287 288 # The nuclear isotope. 289 if hasattr(spin, 'isotope'): 290 isotope_list.append(int(spin.isotope.strip(string.ascii_letters))) 291 else: 292 isotope_list.append(None) 293 294 # The element. 295 if hasattr(spin, 'element'): 296 element_list.append(spin.element) 297 else: 298 element_list.append(None) 299 300 # Diffusion tensor. 301 if hasattr(spin, 'local_tm'): 302 local_tm_list.append(spin.local_tm) 303 else: 304 local_tm_list.append(None) 305 if hasattr(spin, 'local_tm_err'): 306 local_tm_err_list.append(spin.local_tm_err) 307 else: 308 local_tm_err_list.append(None) 309 310 # Model-free parameter values. 311 s2_list.append(spin.s2) 312 s2f_list.append(spin.s2f) 313 s2s_list.append(spin.s2s) 314 te_list.append(spin.te) 315 tf_list.append(spin.tf) 316 ts_list.append(spin.ts) 317 if spin.rex == None: 318 rex_list.append(None) 319 else: 320 rex_list.append(spin.rex / (2.0*pi*rex_frq**2)) 321 322 # Model-free parameter errors. 323 params = ['s2', 's2f', 's2s', 'te', 'tf', 'ts', 'rex'] 324 for param in params: 325 # The error list. 326 err_list = locals()[param+'_err_list'] 327 328 # Append the error. 329 if hasattr(spin, param+'_err'): 330 # The value. 331 val = getattr(spin, param+'_err') 332 333 # Scaling. 334 if param == 'rex' and val != None: 335 val = val / (2.0*pi*rex_frq**2) 336 337 # Append. 338 err_list.append(val) 339 340 # Or None. 341 else: 342 err_list.append(None) 343 344 345 # Opt stats. 346 chi2_list.append(spin.chi2) 347 348 # Model-free model. 349 model_list.append(to_bmrb_model(spin.model)) 350 351 # Check that spin data is present. 352 if not len(model_list): 353 raise RelaxError("No model-free data could be found for any of the currently selected spins.") 354 355 # Convert the molecule names into the entity IDs. 356 entity_ids = zeros(len(mol_name_list), int32) 357 mol_names = get_molecule_names() 358 for i in range(len(mol_name_list)): 359 for j in range(len(mol_names)): 360 if mol_name_list[i] == mol_names[j]: 361 entity_ids[i] = j+1 362 363 364 # Create Supergroup 2 : The citations. 365 ###################################### 366 367 # Generate the citations saveframe. 368 bmrb_write_citations(star) 369 370 371 # Create Supergroup 3 : The molecular assembly saveframes. 372 ########################################################## 373 374 # Generate the entity saveframe. 375 mol_res_spin.bmrb_write_entity(star) 376 377 378 # Create Supergroup 4: The experimental descriptions saveframes. 379 ################################################################# 380 381 # Generate the method saveframes. 382 bmrb_write_methods(star) 383 384 # Generate the software saveframe. 385 software_ids, software_labels = bmrb_write_software(star) 386 387 388 # Create Supergroup 5 : The NMR parameters saveframes. 389 ###################################################### 390 391 # Generate the CSA saveframe. 392 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) 393 394 395 # Create Supergroup 6 : The kinetic data saveframes. 396 #################################################### 397 398 # Generate the relaxation data saveframes. 399 relax_data.bmrb_write(star) 400 401 402 # Create Supergroup 7 : The thermodynamics saveframes. 403 ###################################################### 404 405 # Get the relax software id. 406 for i in range(len(software_ids)): 407 if software_labels[i] == 'relax': 408 software_id = software_ids[i] 409 410 # Generate the model-free data saveframe. 411 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) 412 413 414 # Create Supergroup 8 : The structure determination saveframes. 415 ############################################################### 416 417 # Generate the diffusion tensor saveframes. 418 diffusion_tensor.bmrb_write(star) 419 420 421 # Write the contents to the STAR formatted file. 422 star.write()

423 424

425 - def calculate(self, spin_id=None, scaling_matrix=None, verbosity=1, sim_index=None):

426 """Calculation of the model-free chi-squared value. 427 428 @keyword spin_id: The spin identification string. 429 @type spin_id: str 430 @keyword scaling_matrix: The per-model list of diagonal and square scaling matrices. 431 @type scaling_matrix: list of numpy rank-2, float64 array or list of None 432 @keyword verbosity: The amount of information to print. The higher the value, the greater the verbosity. 433 @type verbosity: int 434 @keyword sim_index: The optional MC simulation index. 435 @type sim_index: int 436 """ 437 438 # Test if sequence data is loaded. 439 if not exists_mol_res_spin_data(): 440 raise RelaxNoSequenceError 441 442 # Determine the model type. 443 model_type = determine_model_type() 444 445 # Test if diffusion tensor data exists. 446 if model_type != 'local_tm' and not diffusion_tensor.diff_data_exists(): 447 raise RelaxNoTensorError('diffusion') 448 449 # Test if the PDB file has been loaded. 450 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere' and not hasattr(cdp, 'structure'): 451 raise RelaxNoPdbError 452 453 # Loop over the spins. 454 for spin, id in spin_loop(spin_id, return_id=True): 455 # Skip deselected spins. 456 if not spin.select: 457 continue 458 459 # Test if the model-free model has been setup. 460 if not spin.model: 461 raise RelaxNoModelError 462 463 # Test if the nuclear isotope type has been set. 464 if not hasattr(spin, 'isotope'): 465 raise RelaxSpinTypeError 466 467 # Test if the model-free parameter values exist. 468 unset_param = are_mf_params_set(spin) 469 if unset_param != None: 470 raise RelaxNoValueError(unset_param) 471 472 # Test if the CSA value has been set. 473 if not hasattr(spin, 'csa') or spin.csa == None: 474 raise RelaxNoValueError("CSA") 475 476 # Test the interatomic data. 477 interatoms = return_interatom_list(spin_hash=spin._hash) 478 for interatom in interatoms: 479 # No relaxation mechanism. 480 if not interatom.dipole_pair: 481 continue 482 483 # Test if unit vectors exist. 484 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere' and not hasattr(interatom, 'vector'): 485 raise RelaxNoVectorsError 486 487 # Test if multiple unit vectors exist. 488 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): 489 raise RelaxMultiVectorError 490 491 # The interacting spin. 492 if spin_id != interatom.spin_id1: 493 spin_id2 = interatom.spin_id1 494 else: 495 spin_id2 = interatom.spin_id2 496 spin2 = return_spin(spin_id=spin_id2) 497 498 # Test if the nuclear isotope type has been set. 499 if not hasattr(spin2, 'isotope'): 500 raise RelaxSpinTypeError 501 502 # Test if the interatomic distance has been set. 503 if not hasattr(interatom, 'r') or interatom.r == None: 504 raise RelaxNoValueError("interatomic distance", spin_id=spin_id, spin_id2=spin_id2) 505 506 # Skip spins where there is no data or errors. 507 if not hasattr(spin, 'ri_data') or not hasattr(spin, 'ri_data_err'): 508 continue 509 510 # Make sure that the errors are strictly positive numbers. 511 for ri_id in cdp.ri_ids: 512 # Alias. 513 err = spin.ri_data_err[ri_id] 514 515 # Checks. 516 if err != None and err == 0.0: 517 raise RelaxError("Zero error for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (id, ri_id)) 518 elif err != None and err < 0.0: 519 raise RelaxError("Negative error of %s for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (err, id, ri_id)) 520 521 # Create the initial parameter vector. 522 param_vector = assemble_param_vector(spin=spin, sim_index=sim_index) 523 524 # The relaxation data optimisation structures. 525 data = relax_data_opt_structs(spin, sim_index=sim_index) 526 527 # The spin data. 528 ri_data = [array(data[0])] 529 ri_data_err = [array(data[1])] 530 num_frq = [data[2]] 531 num_ri = [data[3]] 532 ri_labels = [data[4]] 533 frq = [data[5]] 534 remap_table = [data[6]] 535 noe_r1_table = [data[7]] 536 gx = [periodic_table.gyromagnetic_ratio(spin.isotope)] 537 if sim_index == None: 538 csa = [spin.csa] 539 else: 540 csa = [spin.csa_sim[sim_index]] 541 542 # The interatomic data. 543 interatoms = return_interatom_list(spin_hash=spin._hash) 544 for i in range(len(interatoms)): 545 # No relaxation mechanism. 546 if not interatoms[i].dipole_pair: 547 continue 548 549 # The surrounding spins. 550 if spin_id != interatoms[i].spin_id1: 551 spin_id2 = interatoms[i].spin_id1 552 else: 553 spin_id2 = interatoms[i].spin_id2 554 spin2 = return_spin(spin_id=spin_id2) 555 556 # The data. 557 if sim_index == None: 558 r = [interatoms[i].r] 559 else: 560 r = [interatoms[i].r_sim[sim_index]] 561 562 # Vectors. 563 if model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere': 564 xh_unit_vectors = [interatoms[i].vector] 565 else: 566 xh_unit_vectors = [None] 567 568 # Gyromagnetic ratios. 569 gh = [periodic_table.gyromagnetic_ratio(spin2.isotope)] 570 571 # Count the number of model-free parameters for the residue index. 572 num_params = [len(spin.params)] 573 574 # Repackage the parameter values as a local model (ignore if the diffusion tensor is not fixed). 575 param_values = [assemble_param_vector(model_type='mf')] 576 577 # Package the diffusion tensor parameters. 578 if model_type == 'local_tm': 579 diff_params = [spin.local_tm] 580 diff_type = 'sphere' 581 else: 582 # Diff type. 583 diff_type = cdp.diff_tensor.type 584 585 # Spherical diffusion. 586 if diff_type == 'sphere': 587 diff_params = [cdp.diff_tensor.tm] 588 589 # Spheroidal diffusion. 590 elif diff_type == 'spheroid': 591 diff_params = [cdp.diff_tensor.tm, cdp.diff_tensor.Da, cdp.diff_tensor.theta, cdp.diff_tensor.phi] 592 593 # Ellipsoidal diffusion. 594 elif diff_type == 'ellipsoid': 595 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] 596 597 # Initialise the model-free function. 598 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) 599 600 # Chi-squared calculation. 601 try: 602 chi2 = mf.func(param_vector) 603 except OverflowError: 604 chi2 = 1e200 605 606 # Global chi-squared value. 607 if model_type == 'all' or model_type == 'diff': 608 cdp.chi2 = cdp.chi2 + chi2 609 else: 610 spin.chi2 = chi2

611 612

613 - def create_mc_data(self, data_id=None):

614 """Create the Monte Carlo Ri data. 615 616 @keyword data_id: The spin identification string, as yielded by the base_data_loop() generator method. 617 @type data_id: str 618 @return: The Monte Carlo simulation data. 619 @rtype: list of floats 620 """ 621 622 # Initialise the MC data structure. 623 mc_data = [] 624 625 # Get the spin container and global spin index. 626 spin = return_spin(spin_id=data_id) 627 global_index = find_index(data_id) 628 629 # Skip deselected spins. 630 if not spin.select: 631 return 632 633 # Test if the model is set, returning None to signal that the spin should be skipped. 634 if not hasattr(spin, 'model') or not spin.model: 635 return None 636 637 # Loop over the relaxation data. 638 for ri_id in cdp.ri_ids: 639 # Back calculate the value. 640 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]) 641 642 # Append the value. 643 mc_data.append(value) 644 645 # Return the data. 646 return mc_data

647 648

649 - def data_init(self, data, sim=False):

650 """Initialise the spin specific data structures. 651 652 @param data: The spin ID string from the _base_data_loop_spin() method. 653 @type data: str 654 @keyword sim: The Monte Carlo simulation flag, which if true will initialise the simulation data structure. 655 @type sim: bool 656 """ 657 658 # Get the spin container. 659 spin = return_spin(spin_id=data) 660 661 # Loop over the data structure names. 662 for name in self._PARAMS.loop(scope='spin'): 663 # Blacklisted data structures. 664 if name in ['ri_data', 'ri_data_bc', 'ri_data_err']: 665 continue 666 667 # Data structures which are initially empty arrays. 668 list_data = [ 'params' ] 669 if name in list_data: 670 init_data = [] 671 672 # Set everything else initially to None or False. 673 init_data = None 674 if self._PARAMS.type(name) == bool: 675 init_data = False 676 if name == 'select': 677 init_data = True 678 679 # If the name is not in the spin container, add it. 680 if not hasattr(spin, name): 681 setattr(spin, name, init_data)

682 683

684 - def deselect(self, sim_index=None, model_info=None):

685 """Deselect models or simulations. 686 687 @keyword sim_index: The optional Monte Carlo simulation index. If None, then models will be deselected, otherwise the given simulation will. 688 @type sim_index: None or int 689 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 690 @type model_info: int 691 """ 692 693 # Determine the model type. 694 model_type = determine_model_type() 695 696 # Local models. 697 if model_type == 'mf' or model_type == 'local_tm': 698 # Get the spin. 699 spin = return_spin_from_index(global_index=model_info) 700 701 # Spin deselection. 702 if sim_index == None: 703 spin.select = False 704 705 # Simulation deselection. 706 else: 707 spin.select_sim[sim_index] = False 708 709 # Global models. 710 else: 711 # Global model deselection. 712 if sim_index == None: 713 raise RelaxError("Cannot deselect the global model.") 714 715 # Simulation deselection. 716 else: 717 # Deselect. 718 cdp.select_sim[sim_index] = False

719 720

721 - def duplicate_data(self, pipe_from=None, pipe_to=None, model_info=None, global_stats=False, verbose=True):

722 """Duplicate the data specific to a single model-free model. 723 724 @keyword pipe_from: The data pipe to copy the data from. 725 @type pipe_from: str 726 @keyword pipe_to: The data pipe to copy the data to. 727 @type pipe_to: str 728 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 729 @type model_info: int 730 @keyword global_stats: The global statistics flag. 731 @type global_stats: bool 732 @keyword verbose: A flag which if True will cause info about each spin to be printed out as the sequence is generated. 733 @type verbose: bool 734 """ 735 736 # Arg tests. 737 if model_info == None: 738 raise RelaxError("The model_info argument cannot be None.") 739 740 # First create the pipe_to data pipe, if it doesn't exist, but don't switch to it. 741 if not pipes.has_pipe(pipe_to): 742 pipes.create(pipe_to, pipe_type='mf', switch=False) 743 744 # Get the data pipes. 745 dp_from = pipes.get_pipe(pipe_from) 746 dp_to = pipes.get_pipe(pipe_to) 747 748 # Duplicate all non-sequence specific data. 749 for data_name in dir(dp_from): 750 # Skip the container objects. 751 if data_name in ['diff_tensor', 'mol', 'interatomic', 'structure', 'exp_info']: 752 continue 753 754 # Skip special objects. 755 if search('^_', data_name) or data_name in dp_from.__class__.__dict__: 756 continue 757 758 # Get the original object. 759 data_from = getattr(dp_from, data_name) 760 761 # The data already exists. 762 if hasattr(dp_to, data_name): 763 # Get the object in the target pipe. 764 data_to = getattr(dp_to, data_name) 765 766 # The data must match! 767 if data_from != data_to: 768 raise RelaxError("The object " + repr(data_name) + " is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 769 770 # Skip the data. 771 continue 772 773 # Duplicate the data. 774 setattr(dp_to, data_name, deepcopy(data_from)) 775 776 # Diffusion tensor comparison. 777 if hasattr(dp_from, 'diff_tensor'): 778 # Duplicate the tensor if it doesn't exist. 779 if not hasattr(dp_to, 'diff_tensor'): 780 setattr(dp_to, 'diff_tensor', deepcopy(dp_from.diff_tensor)) 781 782 # Otherwise compare the objects inside the container. 783 else: 784 # Loop over the modifiable objects. 785 for data_name in dp_from.diff_tensor._mod_attr: 786 # Get the original object. 787 data_from = None 788 if hasattr(dp_from.diff_tensor, data_name): 789 data_from = getattr(dp_from.diff_tensor, data_name) 790 791 # Get the target object. 792 if data_from and not hasattr(dp_to.diff_tensor, data_name): 793 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.") 794 elif data_from: 795 data_to = getattr(dp_to.diff_tensor, data_name) 796 else: 797 continue 798 799 # The data must match! 800 if data_from != data_to: 801 raise RelaxError("The object " + repr(data_name) + " is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 802 803 # Structure comparison. 804 if hasattr(dp_from, 'structure'): 805 # Duplicate the tensor if it doesn't exist. 806 if not hasattr(dp_to, 'structure'): 807 setattr(dp_to, 'structure', deepcopy(dp_from.structure)) 808 809 # Otherwise compare the objects inside the container. 810 else: 811 # Modifiable object checks. 812 compare_objects(dp_from.structure, dp_to.structure, pipe_from, pipe_to) 813 814 # Tests for the model and molecule containers. 815 if len(dp_from.structure.structural_data) != len(dp_from.structure.structural_data): 816 raise RelaxError("The number of structural models is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 817 818 # Loop over the models. 819 for i in range(len(dp_from.structure.structural_data)): 820 # Alias. 821 model_from = dp_from.structure.structural_data[i] 822 model_to = dp_to.structure.structural_data[i] 823 824 # Model numbers. 825 if model_from.num != model_to.num: 826 raise RelaxError("The structure models are not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 827 828 # Molecule number. 829 if len(model_from.mol) != len(model_to.mol): 830 raise RelaxError("The number of molecules is not consistent between the pipes " + repr(pipe_from) + " and " + repr(pipe_to) + ".") 831 832 # Loop over the models. 833 for mol_index in range(len(model_from.mol)): 834 # Modifiable object checks. 835 compare_objects(model_from.mol[mol_index], model_to.mol[mol_index], pipe_from, pipe_to) 836 837 # No sequence data, so skip the rest. 838 if dp_from.mol.is_empty(): 839 return 840 841 # Duplicate the sequence data if it doesn't exist. 842 if dp_to.mol.is_empty(): 843 sequence.copy(pipe_from=pipe_from, pipe_to=pipe_to, preserve_select=True, verbose=verbose) 844 845 # Duplicate the interatomic data if it doesn't exist. 846 if dp_to.interatomic.is_empty(): 847 interatomic.copy(pipe_from=pipe_from, pipe_to=pipe_to, verbose=verbose) 848 849 # Determine the model type of the original data pipe. 850 pipes.switch(pipe_from) 851 model_type = determine_model_type() 852 853 # Sequence specific data. 854 spin, spin_id = return_spin_from_index(global_index=model_info, pipe=pipe_from, return_spin_id=True) 855 if model_type == 'mf' or (model_type == 'local_tm' and not global_stats): 856 # Duplicate the spin specific data. 857 for name in dir(spin): 858 # Skip special objects. 859 if search('^__', name): 860 continue 861 862 # Get the object. 863 obj = getattr(spin, name) 864 865 # Skip methods. 866 if isinstance(obj, MethodType): 867 continue 868 869 # Duplicate the object. 870 new_obj = deepcopy(getattr(spin, name)) 871 setattr(dp_to.mol[spin._mol_index].res[spin._res_index].spin[spin._spin_index], name, new_obj) 872 873 # Duplicate the relaxation active spins which have not been copied yet. 874 interatoms = interatomic.return_interatom_list(spin_hash=spin._hash) 875 for interatom in interatoms: 876 # No relaxation mechanism. 877 if not interatom.dipole_pair: 878 continue 879 880 # The interacting spin. 881 if spin_id != interatom.spin_id1: 882 spin_id2 = interatom.spin_id1 883 else: 884 spin_id2 = interatom.spin_id2 885 spin2 = return_spin(spin_id=spin_id2) 886 887 # Duplicate the spin specific data. 888 mol_index, res_index, spin_index = return_spin_indices(spin_id=spin_id2) 889 dp_to.mol[mol_index].res[res_index].spin[spin_index] = deepcopy(spin2) 890 891 # Other data types. 892 else: 893 # Duplicate all the spin specific data. 894 dp_to.mol = deepcopy(dp_from.mol)

895 896

897 - def eliminate(self, name, value, args, sim=None, model_info=None):

898 """Model-free model elimination, parameter by parameter. 899 900 @param name: The parameter name. 901 @type name: str 902 @param value: The parameter value. 903 @type value: float 904 @param args: The c1 and c2 elimination constant overrides. 905 @type args: None or tuple of float 906 @keyword sim: The Monte Carlo simulation index. 907 @type sim: int 908 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 909 @type model_info: int 910 @return: True if the model is to be eliminated, False otherwise. 911 @rtype: bool 912 """ 913 914 # Default values. 915 c1 = 50.0 * 1e-9 916 c2 = 1.5 917 918 # Depack the arguments. 919 if args != None: 920 c1, c2 = args 921 922 # Determine the model type. 923 model_type = determine_model_type() 924 925 # Can't handle this one yet! 926 if model_type != 'mf' and model_type != 'local_tm': 927 raise RelaxError("Elimination of the global model is not yet supported.") 928 929 # Get the spin and it's id string. 930 spin, spin_id = return_spin_from_index(global_index=model_info, return_spin_id=True) 931 932 # Get the tm value. 933 if model_type == 'local_tm': 934 tm = spin.local_tm 935 else: 936 tm = cdp.diff_tensor.tm 937 938 # No tm value set, so skip the tests (no elimination). 939 if tm == None: 940 return False 941 942 # Local tm. 943 if name == 'local_tm' and value >= c1: 944 if sim == None: 945 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)) 946 else: 947 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)) 948 return True 949 950 # Internal correlation times. 951 if match('t[efs]', name) and value >= c2 * tm: 952 if sim == None: 953 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)) 954 else: 955 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)) 956 return True 957 958 # Accept model. 959 return False

960 961

962 - def get_param_names(self, model_info=None):

963 """Return a vector of parameter names. 964 965 @keyword model_info: The model information from model_loop(). This index 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 @return: The vector of parameter names. 968 @rtype: list of str 969 """ 970 971 # Determine the model type. 972 model_type = determine_model_type() 973 974 # Get the spin ids. 975 if model_type == 'mf' or model_type == 'local_tm': 976 # Get the spin and it's id string. 977 spin, spin_id = return_spin_from_index(global_index=model_info, return_spin_id=True) 978 else: 979 spin_id = None 980 981 # Assemble and return the parameter names. 982 return assemble_param_names(model_type, spin_id=spin_id)

983 984

985 - def get_param_values(self, model_info=None, sim_index=None):

986 """Return a vector of parameter values. 987 988 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 989 @type model_info: int 990 @keyword sim_index: The Monte Carlo simulation index. 991 @type sim_index: int 992 @return: The vector of parameter values. 993 @rtype: list of str 994 """ 995 996 # Test if the model-free models have been set up. 997 for spin in spin_loop(): 998 # Skip deselected spins. 999 if not spin.select: 1000 continue 1001 1002 # Not setup. 1003 if not spin.model: 1004 raise RelaxNoModelError 1005 1006 # Determine the model type. 1007 model_type = determine_model_type() 1008 1009 # Set the spin container (to None if the model is global). 1010 if model_type == 'mf' or model_type == 'local_tm': 1011 spin = return_spin_from_index(global_index=model_info) 1012 else: 1013 spin = None 1014 1015 # Assemble the parameter values and return them. 1016 return assemble_param_vector(spin=spin, sim_index=sim_index, model_type=model_type)

1017 1018

1019 - def grid_search(self, lower=None, upper=None, inc=None, scaling_matrix=None, constraints=True, verbosity=1, sim_index=None):

1020 """The model-free grid search function. 1021 1022 @keyword lower: The per-model lower bounds of the grid search which must be equal to the number of parameters in the model. 1023 @type lower: list of lists of numbers 1024 @keyword upper: The per-model upper bounds of the grid search which must be equal to the number of parameters in the model. 1025 @type upper: list of lists of numbers 1026 @keyword inc: The per-model increments for each dimension of the space for the grid search. The number of elements in the array must equal to the number of parameters in the model. 1027 @type inc: list of lists of int 1028 @keyword scaling_matrix: The per-model list of diagonal and square scaling matrices. 1029 @type scaling_matrix: list of numpy rank-2, float64 array or list of None 1030 @keyword constraints: If True, constraints are applied during the grid search (eliminating parts of the grid). If False, no constraints are used. 1031 @type constraints: bool 1032 @keyword verbosity: A flag specifying the amount of information to print. The higher the value, the greater the verbosity. 1033 @type verbosity: int 1034 @keyword sim_index: The index of the simulation to apply the grid search to. If None, the normal model is optimised. 1035 @type sim_index: int 1036 """ 1037 1038 # Minimisation. 1039 self.minimise(min_algor='grid', lower=lower, upper=upper, inc=inc, scaling_matrix=scaling_matrix, constraints=constraints, verbosity=verbosity, sim_index=sim_index)

1040 1041

1042 - def map_bounds(self, param, spin_id=None):

1043 """Create bounds for the OpenDX mapping function. 1044 1045 @param param: The name of the parameter to return the lower and upper bounds of. 1046 @type param: str 1047 @param spin_id: The spin identification string. 1048 @type spin_id: str 1049 @return: The upper and lower bounds of the parameter. 1050 @rtype: list of float 1051 """ 1052 1053 # Diffusion tensor bounds. 1054 if self._PARAMS.scope(param) == 'global': 1055 return diffusion_tensor.map_bounds(param) 1056 1057 # Get the spin. 1058 spin = return_spin(spin_id=spin_id) 1059 1060 # {s2, s2f, s2s}. 1061 if search('^s2', param): 1062 return [0.0, 1.0] 1063 1064 # {local tm, te, tf, ts}. 1065 elif search('^t', param) or param == 'local_tm': 1066 return [0.0, 1e-8] 1067 1068 # Rex. 1069 elif param == 'rex': 1070 return [0.0, 30.0 / (2.0 * pi * cdp.spectrometer_frq[cdp.ri_ids[0]])**2] 1071 1072 # Interatomic distances. 1073 elif param == 'r': 1074 return [1.0 * 1e-10, 1.1 * 1e-10] 1075 1076 # CSA. 1077 elif param == 'csa': 1078 return [-100 * 1e-6, -300 * 1e-6]

1079 1080

1081 - def minimise(self, min_algor=None, min_options=None, func_tol=None, grad_tol=None, max_iterations=None, constraints=False, scaling_matrix=None, verbosity=0, sim_index=None, lower=None, upper=None, inc=None):

1082 """Model-free minimisation function. 1083 1084 Three categories of models exist for which the approach to minimisation is different. These 1085 are: 1086 1087 Single spin optimisations: The 'mf' and 'local_tm' model types which are the 1088 model-free parameters for single spins, optionally with a local tm parameter. These 1089 models have no global parameters. 1090 1091 Diffusion tensor optimisations: The 'diff' diffusion tensor model type. No spin 1092 specific parameters exist. 1093 1094 Optimisation of everything: The 'all' model type consisting of all model-free and all 1095 diffusion tensor parameters. 1096 1097 1098 @keyword min_algor: The minimisation algorithm to use. 1099 @type min_algor: str 1100 @keyword min_options: An array of options to be used by the minimisation algorithm. 1101 @type min_options: array of str 1102 @keyword func_tol: The function tolerance which, when reached, terminates optimisation. Setting this to None turns of the check. 1103 @type func_tol: None or float 1104 @keyword grad_tol: The gradient tolerance which, when reached, terminates optimisation. Setting this to None turns of the check. 1105 @type grad_tol: None or float 1106 @keyword max_iterations: The maximum number of iterations for the algorithm. 1107 @type max_iterations: int 1108 @keyword constraints: If True, constraints are used during optimisation. 1109 @type constraints: bool 1110 @keyword scaling_matrix: The per-model list of diagonal and square scaling matrices. 1111 @type scaling_matrix: list of numpy rank-2, float64 array or list of None 1112 @keyword verbosity: The amount of information to print. The higher the value, the greater the verbosity. 1113 @type verbosity: int 1114 @keyword sim_index: The index of the simulation to optimise. This should be None if normal optimisation is desired. 1115 @type sim_index: None or int 1116 @keyword lower: The per-model lower bounds of the grid search which must be equal to the number of parameters in the model. This optional argument is only used when doing a grid search. 1117 @type lower: list of lists of numbers 1118 @keyword upper: The per-model upper bounds of the grid search which must be equal to the number of parameters in the model. This optional argument is only used when doing a grid search. 1119 @type upper: list of lists of numbers 1120 @keyword inc: The per-model increments for each dimension of the space for the grid search. The number of elements in the array must equal to the number of parameters in the model. This argument is only used when doing a grid search. 1121 @type inc: list of lists of int 1122 """ 1123 1124 # Test if sequence data is loaded. 1125 if not exists_mol_res_spin_data(): 1126 raise RelaxNoSequenceError 1127 1128 # Test if the model-free model has been setup, and that the nuclear isotope types have been set. 1129 for spin in spin_loop(): 1130 # Skip deselected spins. 1131 if not spin.select: 1132 continue 1133 1134 # Not setup. 1135 if not spin.model: 1136 raise RelaxNoModelError 1137 1138 # Test if the nuclear isotope type has been set. 1139 if not hasattr(spin, 'isotope'): 1140 raise RelaxSpinTypeError 1141 1142 # Containers for the model-free data and optimisation parameters. 1143 data_store = Data_container() 1144 opt_params = Data_container() 1145 1146 # Add the imported parameters. 1147 data_store.h_bar = h_bar 1148 data_store.mu0 = mu0 1149 opt_params.min_algor = min_algor 1150 opt_params.min_options = min_options 1151 opt_params.func_tol = func_tol 1152 opt_params.grad_tol = grad_tol 1153 opt_params.max_iterations = max_iterations 1154 1155 # Add the keyword args. 1156 opt_params.verbosity = verbosity 1157 1158 # Determine the model type. 1159 data_store.model_type = determine_model_type() 1160 if not data_store.model_type: 1161 return 1162 1163 # Model type for the back-calculate function. 1164 if min_algor == 'back_calc' and data_store.model_type != 'local_tm': 1165 data_store.model_type = 'mf' 1166 1167 # Test if diffusion tensor data exists. 1168 if data_store.model_type != 'local_tm' and not diffusion_tensor.diff_data_exists(): 1169 raise RelaxNoTensorError('diffusion') 1170 1171 # Tests for the PDB file and unit vectors. 1172 if data_store.model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere': 1173 # Test if the structure file has been loaded. 1174 if not hasattr(cdp, 'structure'): 1175 raise RelaxNoPdbError 1176 1177 # Test if unit vectors exist. 1178 for spin, spin_id in spin_loop(return_id=True): 1179 # Skip deselected spins. 1180 if not spin.select: 1181 continue 1182 1183 # Get the interatomic data container. 1184 interatoms = return_interatom_list(spin_hash=spin._hash) 1185 1186 # Unit vectors. 1187 for i in range(len(interatoms)): 1188 # No relaxation mechanism. 1189 if not interatoms[i].dipole_pair: 1190 continue 1191 1192 # Check for the vectors. 1193 if not hasattr(interatoms[i], 'vector'): 1194 raise RelaxNoVectorsError 1195 1196 # Test if the model-free parameter values are set for minimising diffusion tensor parameters by themselves. 1197 if data_store.model_type == 'diff': 1198 # Loop over the sequence. 1199 for spin in spin_loop(): 1200 unset_param = are_mf_params_set(spin) 1201 if unset_param != None: 1202 raise RelaxNoValueError(unset_param) 1203 1204 # Print out. 1205 if verbosity >= 1: 1206 if data_store.model_type == 'mf': 1207 print("Only the model-free parameters for single spins will be used.") 1208 elif data_store.model_type == 'local_mf': 1209 print("Only a local tm value together with the model-free parameters for single spins will be used.") 1210 elif data_store.model_type == 'diff': 1211 print("Only diffusion tensor parameters will be used.") 1212 elif data_store.model_type == 'all': 1213 print("The diffusion tensor parameters together with the model-free parameters for all spins will be used.") 1214 1215 # Test if the CSA and interatomic distances have been set. 1216 for spin, spin_id in spin_loop(return_id=True): 1217 # Skip deselected spins. 1218 if not spin.select: 1219 continue 1220 1221 # CSA value. 1222 if not hasattr(spin, 'csa') or spin.csa == None: 1223 raise RelaxNoValueError("CSA") 1224 1225 # Get the interatomic data container. 1226 interatoms = return_interatom_list(spin_hash=spin._hash) 1227 1228 # Interatomic distances. 1229 count = 0 1230 for i in range(len(interatoms)): 1231 # No relaxation mechanism. 1232 if not interatoms[i].dipole_pair: 1233 continue 1234 1235 # Check for the distances. 1236 if not hasattr(interatoms[i], 'r') or interatoms[i].r == None: 1237 raise RelaxNoValueError("interatomic distance", spin_id=interatoms[i].spin_id1, spin_id2=interatoms[i].spin_id2) 1238 1239 # Count the number of interactions. 1240 count += 1 1241 1242 # Too many interactions. 1243 if count > 1: 1244 raise RelaxError("The spin '%s' has %s dipolar relaxation interactions defined, but only a maximum of one is currently supported." % (spin_id, count)) 1245 1246 # Number of spins, minimisation instances, and data sets for each model type. 1247 if data_store.model_type == 'mf' or data_store.model_type == 'local_tm': 1248 num_data_sets = 1 1249 data_store.num_spins = 1 1250 elif data_store.model_type == 'diff' or data_store.model_type == 'all': 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_data_sets = 0 1257 data_store.num_spins = 1 1258 1259 # Get the Processor box singleton (it contains the Processor instance) and alias the Processor. 1260 processor_box = Processor_box() 1261 processor = processor_box.processor 1262 1263 # Loop over the models. 1264 for index in self.model_loop(): 1265 # Get the spin container if required. 1266 if data_store.model_type == 'diff' or data_store.model_type == 'all': 1267 spin_index = None 1268 spin, data_store.spin_id = None, None 1269 elif min_algor == 'back_calc': 1270 spin_index = opt_params.min_options[0] 1271 spin, data_store.spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 1272 else: 1273 spin_index = index 1274 spin, data_store.spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True) 1275 1276 # Individual spin stuff. 1277 if spin and (data_store.model_type == 'mf' or data_store.model_type == 'local_tm') and not min_algor == 'back_calc': 1278 # Skip deselected spins. 1279 if not spin.select: 1280 continue 1281 1282 # Skip spins missing relaxation data or errors. 1283 if not hasattr(spin, 'ri_data') or not hasattr(spin, 'ri_data_err'): 1284 continue 1285 1286 # Skip spins missing the dipolar interaction. 1287 if spin and (data_store.model_type == 'mf' or data_store.model_type == 'local_tm'): 1288 interatoms = return_interatom_list(spin_hash=spin._hash) 1289 if not len(interatoms): 1290 continue 1291 1292 # Parameter vector and diagonal scaling. 1293 if min_algor == 'back_calc': 1294 # Create the initial parameter vector. 1295 opt_params.param_vector = assemble_param_vector(spin=spin, model_type=data_store.model_type) 1296 1297 # Diagonal scaling. 1298 data_store.scaling_matrix = None 1299 1300 else: 1301 # Create the initial parameter vector. 1302 opt_params.param_vector = assemble_param_vector(spin=spin, sim_index=sim_index) 1303 1304 # The number of parameters. 1305 num_params = len(opt_params.param_vector) 1306 1307 # Diagonal scaling. 1308 data_store.scaling_matrix = scaling_matrix[index] 1309 if data_store.scaling_matrix is not None: 1310 opt_params.param_vector = dot(inv(data_store.scaling_matrix), opt_params.param_vector) 1311 1312 # Store the grid search options. 1313 opt_params.lower, opt_params.upper, opt_params.inc = None, None, None 1314 if lower != None: 1315 opt_params.lower = lower[index] 1316 if upper != None: 1317 opt_params.upper = upper[index] 1318 if inc != None: 1319 opt_params.inc = inc[index] 1320 1321 # Scaling of values for the set function. 1322 if match('^[Ss]et', min_algor): 1323 opt_params.min_options = dot(inv(data_store.scaling_matrix), opt_params.min_options) 1324 1325 # Linear constraints. 1326 if constraints: 1327 opt_params.A, opt_params.b = linear_constraints(num_params, model_type=data_store.model_type, spin=spin, scaling_matrix=data_store.scaling_matrix) 1328 else: 1329 opt_params.A, opt_params.b = None, None 1330 1331 # Get the data for minimisation. 1332 minimise_data_setup(data_store, min_algor, num_data_sets, opt_params.min_options, spin=spin, sim_index=sim_index) 1333 1334 # Setup the minimisation algorithm when constraints are present. 1335 if constraints and not match('^[Gg]rid', min_algor): 1336 algor = opt_params.min_options[0] 1337 else: 1338 algor = min_algor 1339 1340 # Initialise the function to minimise (for back-calculation and LM minimisation). 1341 if min_algor == 'back_calc' or match('[Ll][Mm]$', algor) or match('[Ll]evenburg-[Mm]arquardt$', algor): 1342 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) 1343 1344 # Levenberg-Marquardt minimisation. 1345 if match('[Ll][Mm]$', algor) or match('[Ll]evenburg-[Mm]arquardt$', algor): 1346 # Total number of ri. 1347 number_ri = 0 1348 for k in range(len(ri_data_err)): 1349 number_ri = number_ri + len(ri_data_err[k]) 1350 1351 # Reconstruct the error data structure. 1352 lm_error = zeros(number_ri, float64) 1353 index = 0 1354 for k in range(len(ri_data_err)): 1355 lm_error[index:index+len(ri_data_err[k])] = ri_data_err[k] 1356 index = index + len(ri_data_err[k]) 1357 1358 opt_params.min_options = opt_params.min_options + (mf.lm_dri, lm_error) 1359 1360 # Back-calculation. 1361 if min_algor == 'back_calc': 1362 return mf.calc_ri() 1363 1364 # Parallelised grid search for the diffusion parameter space. 1365 if match('^[Gg]rid', min_algor) and data_store.model_type == 'diff': 1366 # Print out. 1367 print("Parallelised diffusion tensor grid search.") 1368 1369 # Loop over each grid subdivision. 1370 for subdivision in grid_split(divisions=processor.processor_size(), lower=opt_params.lower, upper=opt_params.upper, inc=opt_params.inc, verbosity=verbosity): 1371 # Set the points. 1372 opt_params.subdivision = subdivision 1373 1374 # Grid search initialisation. 1375 command = MF_grid_command() 1376 1377 # Pass in the data and optimisation parameters. 1378 command.store_data(deepcopy(data_store), deepcopy(opt_params)) 1379 1380 # Set up the model-free memo and add it to the processor queue. 1381 memo = MF_memo(model_free=self, model_type=data_store.model_type, spin=spin, sim_index=sim_index, scaling_matrix=data_store.scaling_matrix) 1382 processor.add_to_queue(command, memo) 1383 1384 # Execute the queued elements. 1385 processor.run_queue() 1386 1387 # Exit this method. 1388 return 1389 1390 # Normal grid search (command initialisation). 1391 if search('^[Gg]rid', min_algor): 1392 command = MF_grid_command() 1393 1394 # Minimisation of all other model types (command initialisation). 1395 else: 1396 command = MF_minimise_command() 1397 1398 # Pass in the data and optimisation parameters. 1399 command.store_data(deepcopy(data_store), deepcopy(opt_params)) 1400 1401 # Set up the model-free memo and add it to the processor queue. 1402 memo = MF_memo(model_free=self, model_type=data_store.model_type, spin=spin, sim_index=sim_index, scaling_matrix=data_store.scaling_matrix) 1403 processor.add_to_queue(command, memo) 1404 1405 # Execute the queued elements. 1406 processor.run_queue()

1407 1408

1409 - def model_desc(self, model_info=None):

1410 """Return a description of the model. 1411 1412 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1413 @type model_info: int 1414 @return: The model description. 1415 @rtype: str 1416 """ 1417 1418 # Determine the model type. 1419 model_type = determine_model_type() 1420 1421 # Global models. 1422 if model_type == 'all': 1423 return "Global model - all diffusion tensor parameters and spin specific model-free parameters." 1424 elif model_type == 'diff': 1425 return "Diffusion tensor model." 1426 1427 # Spin specific model. 1428 else: 1429 # Get the spin container. 1430 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 1431 1432 # Return the description. 1433 return "Model-free model of spin '%s'." % spin_id

1434 1435

1436 - def model_loop(self):

1437 """Generator method for looping over the models (global or local). 1438 1439 If the model type is 'all' or 'diff', then this yields the single value of zero. Otherwise 1440 the global spin index is yielded. 1441 1442 1443 @return: The model index. This index is zero for the global models or equal to the global spin 1444 index (which covers the molecule, residue, and spin indices). 1445 @rtype: int 1446 """ 1447 1448 # Determine the model type. 1449 model_type = determine_model_type() 1450 1451 # Global model. 1452 if model_type == 'all' or model_type == 'diff': 1453 yield 0 1454 1455 # Spin specific models. 1456 else: 1457 # Loop over the spins. 1458 global_index = -1 1459 for spin in spin_loop(): 1460 # Increment the global spin index. 1461 global_index = global_index + 1 1462 1463 # Yield the spin index. 1464 yield global_index

1465 1466

1467 - def model_statistics(self, model_info=None, spin_id=None, global_stats=None):

1468 """Return the k, n, and chi2 model statistics. 1469 1470 k - number of parameters. 1471 n - number of data points. 1472 chi2 - the chi-squared value. 1473 1474 1475 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1476 @type model_info: int 1477 @keyword spin_id: The spin identification string. Either this or the instance keyword argument must be supplied. 1478 @type spin_id: None or str 1479 @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. 1480 @type global_stats: None or bool 1481 @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). 1482 @rtype: tuple of (int, int, float) 1483 """ 1484 1485 # Bad argument combination. 1486 if model_info == None and spin_id == None: 1487 raise RelaxError("Either the model_info or spin_id argument must be supplied.") 1488 elif model_info != None and spin_id != None: 1489 raise RelaxError("The model_info arg " + repr(model_info) + " and spin_id arg " + repr(spin_id) + " clash. Only one should be supplied.") 1490 1491 # Determine the model type. 1492 model_type = determine_model_type() 1493 1494 # Determine if local or global statistics will be returned. 1495 if global_stats == None: 1496 if model_type in ['mf', 'local_tm']: 1497 global_stats = False 1498 else: 1499 global_stats = True 1500 1501 # Statistics for a single residue. 1502 if not global_stats: 1503 # Get the SpinContainer. 1504 if spin_id: 1505 spin = return_spin(spin_id) 1506 else: 1507 spin = return_spin_from_index(model_info) 1508 1509 # Skip deselected residues. 1510 if not spin.select: 1511 return None, None, None 1512 1513 # Missing data sets. 1514 if not hasattr(spin, 'ri_data'): 1515 return None, None, None 1516 1517 # Count the number of parameters. 1518 param_vector = assemble_param_vector(spin=spin) 1519 k = len(param_vector) 1520 1521 # Count the number of data points. 1522 n = len(spin.ri_data) 1523 1524 # The chi2 value. 1525 chi2 = spin.chi2 1526 1527 # Global stats. 1528 elif global_stats: 1529 # Count the number of parameters. 1530 param_vector = assemble_param_vector() 1531 k = len(param_vector) 1532 1533 # Count the number of data points. 1534 n = 0 1535 chi2 = 0 1536 for spin in spin_loop(): 1537 # Skip deselected residues. 1538 if not spin.select: 1539 continue 1540 1541 # Skip residues with no relaxation data. 1542 if not hasattr(spin, 'ri_data') or not len(spin.ri_data): 1543 continue 1544 1545 n = n + len(spin.ri_data) 1546 1547 # Local tm models. 1548 if model_type == 'local_tm': 1549 chi2 = chi2 + spin.chi2 1550 1551 # The chi2 value. 1552 if model_type != 'local_tm': 1553 if not hasattr(cdp, 'chi2'): 1554 raise RelaxError("Global statistics are not available, most likely because the global model has not been optimised.") 1555 chi2 = cdp.chi2 1556 1557 # Return the data. 1558 return k, n, chi2

1559 1560

1561 - def model_type(self):

1562 """Return the type of the model, either being 'local' or 'global'. 1563 1564 @return: The model type, one of 'local' or 'global'. 1565 @rtype: str 1566 """ 1567 1568 # Determine the model type. 1569 model_type = determine_model_type() 1570 1571 # Global models. 1572 if model_type in ['all', 'diff']: 1573 return 'global' 1574 1575 # Local models. 1576 else: 1577 return 'local'

1578 1579

1580 - def num_instances(self):

1581 """Function for returning the number of instances. 1582 1583 @return: The number of instances used for optimisation. Either the number of spins if 1584 the local optimisations are setup ('mf' and 'local_tm'), or 1 for the global 1585 models. 1586 @rtype: int 1587 """ 1588 1589 # Test if sequence data exists. 1590 if not exists_mol_res_spin_data(): 1591 return 0 1592 1593 # Determine the model type. 1594 model_type = determine_model_type() 1595 1596 # Sequence specific data. 1597 if model_type == 'mf' or model_type == 'local_tm': 1598 return count_spins() 1599 1600 # Other data. 1601 elif model_type == 'diff' or model_type == 'all': 1602 return 1 1603 1604 # Should not be here. 1605 else: 1606 raise RelaxFault

1607 1608

1609 - def overfit_deselect(self, data_check=True, verbose=True):

1610 """Deselect spins which have insufficient data to support minimisation. 1611 1612 @keyword data_check: A flag to signal if the presence of base data is to be checked for. 1613 @type data_check: bool 1614 @keyword verbose: A flag which if True will allow printouts. 1615 @type verbose: bool 1616 """ 1617 1618 # Print out. 1619 if verbose: 1620 print("\nOver-fit spin deselection:") 1621 1622 # Test if sequence data exists. 1623 if not exists_mol_res_spin_data(): 1624 raise RelaxNoSequenceError 1625 1626 # Is structural data required? 1627 need_vect = False 1628 if hasattr(cdp, 'diff_tensor') and (cdp.diff_tensor.type == 'spheroid' or cdp.diff_tensor.type == 'ellipsoid'): 1629 need_vect = True 1630 1631 # Loop over the sequence. 1632 deselect_flag = False 1633 spin_count = 0 1634 for spin, spin_id in spin_loop(return_id=True): 1635 # Skip deselected spins. 1636 if not spin.select: 1637 continue 1638 1639 # The interatomic data. 1640 interatoms = interatomic.return_interatom_list(spin_hash=spin._hash) 1641 1642 # Loop over the interatomic data. 1643 dipole_relax = False 1644 for i in range(len(interatoms)): 1645 # No dipolar relaxation mechanism. 1646 if not interatoms[i].dipole_pair: 1647 continue 1648 1649 # The surrounding spins. 1650 if spin_id != interatoms[i].spin_id1: 1651 spin_id2 = interatoms[i].spin_id1 1652 else: 1653 spin_id2 = interatoms[i].spin_id2 1654 spin2 = return_spin(spin_id2) 1655 1656 # Dipolar relaxation flag. 1657 dipole_relax = True 1658 1659 # No relaxation mechanism. 1660 if not dipole_relax or not hasattr(spin, 'csa') or spin.csa == None: 1661 warn(RelaxDeselectWarning(spin_id, 'an absence of relaxation mechanisms')) 1662 spin.select = False 1663 deselect_flag = True 1664 continue 1665 1666 # Data checks. 1667 if data_check: 1668 # The number of relaxation data points (and for infinite data). 1669 data_points = 0 1670 inf_data = False 1671 if hasattr(cdp, 'ri_ids') and hasattr(spin, 'ri_data'): 1672 for id in cdp.ri_ids: 1673 if id in spin.ri_data and spin.ri_data[id] != None: 1674 data_points += 1 1675 1676 # Infinite data! 1677 if isInf(spin.ri_data[id]): 1678 inf_data = True 1679 1680 # Infinite data. 1681 if inf_data: 1682 warn(RelaxDeselectWarning(spin_id, 'infinite relaxation data')) 1683 spin.select = False 1684 deselect_flag = True 1685 continue 1686 1687 # Relaxation data must exist! 1688 if not hasattr(spin, 'ri_data'): 1689 warn(RelaxDeselectWarning(spin_id, 'missing relaxation data')) 1690 spin.select = False 1691 deselect_flag = True 1692 continue 1693 1694 # Require 3 or more relaxation data points. 1695 elif data_points < 3: 1696 warn(RelaxDeselectWarning(spin_id, 'insufficient relaxation data, 3 or more data points are required')) 1697 spin.select = False 1698 deselect_flag = True 1699 continue 1700 1701 # Require at least as many data points as params to prevent over-fitting. 1702 elif hasattr(spin, 'params') and spin.params and len(spin.params) > data_points: 1703 warn(RelaxDeselectWarning(spin_id, 'over-fitting - more parameters than data points')) 1704 spin.select = False 1705 deselect_flag = True 1706 continue 1707 1708 # Test for structural data if required. 1709 for i in range(len(interatoms)): 1710 # No dipolar relaxation mechanism. 1711 if not interatoms[i].dipole_pair: 1712 continue 1713 1714 # Check the unit vectors. 1715 if need_vect: 1716 if not hasattr(interatoms[i], 'vector') or interatoms[i].vector is None: 1717 warn(RelaxDeselectWarning(spin_id, 'missing structural data')) 1718 spin.select = False 1719 deselect_flag = True 1720 continue 1721 1722 # Increment the spin number. 1723 spin_count += 1 1724 1725 # No spins selected, so fail hard to prevent the user from going any further. 1726 if spin_count == 0: 1727 warn(RelaxWarning("No spins are selected therefore the optimisation or calculation cannot proceed.")) 1728 1729 # Final printout. 1730 if verbose and not deselect_flag: 1731 print("No spins have been deselected.")

1732 1733

1734 - def print_model_title(self, prefix=None, model_info=None):

1735 """Print out the model title. 1736 1737 @keyword prefix: The starting text of the title. This should be printed out first, followed by the model information text. 1738 @type prefix: str 1739 @keyword model_info: The model information from model_loop(). 1740 @type model_info: unknown 1741 """ 1742 1743 # Determine the model type. 1744 model_type = determine_model_type() 1745 1746 # Local models. 1747 if model_type == 'mf' or model_type == 'local_tm': 1748 spin, spin_id = return_spin_from_index(global_index=model_info, return_spin_id=True) 1749 text = "%sSpin '%s'" % (prefix, spin_id) 1750 1751 # Global models. 1752 else: 1753 text = prefix + "Global model" 1754 1755 # The printout. 1756 subsection(file=sys.stdout, text=text, prespace=2)

1757 1758

1759 - def set_error(self, index, error, model_info=None):

1760 """Set the parameter errors. 1761 1762 @param index: The index of the parameter to set the errors for. 1763 @type index: int 1764 @param error: The error value. 1765 @type error: float 1766 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1767 @type model_info: int 1768 """ 1769 1770 # Parameter increment counter. 1771 inc = 0 1772 1773 # Determine the model type. 1774 model_type = determine_model_type() 1775 1776 # Get the parameter object names. 1777 param_names = self.data_names(set='params', scope='spin') 1778 1779 1780 # Diffusion tensor parameter errors. 1781 #################################### 1782 1783 if model_type == 'diff' or model_type == 'all': 1784 # Spherical diffusion. 1785 if cdp.diff_tensor.type == 'sphere': 1786 # Return the parameter array. 1787 if index == 0: 1788 cdp.diff_tensor.set(param='tm', value=error, category='err') 1789 1790 # Increment. 1791 inc = inc + 1 1792 1793 # Spheroidal diffusion. 1794 elif cdp.diff_tensor.type == 'spheroid': 1795 # Return the parameter array. 1796 if index == 0: 1797 cdp.diff_tensor.set(param='tm', value=error, category='err') 1798 elif index == 1: 1799 cdp.diff_tensor.set(param='Da', value=error, category='err') 1800 elif index == 2: 1801 cdp.diff_tensor.set(param='theta', value=error, category='err') 1802 elif index == 3: 1803 cdp.diff_tensor.set(param='phi', value=error, category='err') 1804 1805 # Increment. 1806 inc = inc + 4 1807 1808 # Ellipsoidal diffusion. 1809 elif cdp.diff_tensor.type == 'ellipsoid': 1810 # Return the parameter array. 1811 if index == 0: 1812 cdp.diff_tensor.set(param='tm', value=error, category='err') 1813 elif index == 1: 1814 cdp.diff_tensor.set(param='Da', value=error, category='err') 1815 elif index == 2: 1816 cdp.diff_tensor.set(param='Dr', value=error, category='err') 1817 elif index == 3: 1818 cdp.diff_tensor.set(param='alpha', value=error, category='err') 1819 elif index == 4: 1820 cdp.diff_tensor.set(param='beta', value=error, category='err') 1821 elif index == 5: 1822 cdp.diff_tensor.set(param='gamma', value=error, category='err') 1823 1824 # Increment. 1825 inc = inc + 6 1826 1827 1828 # Model-free parameter errors for the model type 'all'. 1829 ####################################################### 1830 1831 if model_type == 'all': 1832 # Loop over the spins. 1833 for spin in spin_loop(): 1834 # Skip deselected spins. 1835 if not spin.select: 1836 continue 1837 1838 # Loop over the residue specific parameters. 1839 for param in param_names: 1840 # Return the parameter array. 1841 if index == inc: 1842 setattr(spin, param + "_err", error) 1843 1844 # Increment. 1845 inc = inc + 1 1846 1847 1848 # Model-free parameters for the model types 'mf' and 'local_tm'. 1849 ################################################################ 1850 1851 if model_type == 'mf' or model_type == 'local_tm': 1852 # Get the spin container. 1853 spin = return_spin_from_index(model_info) 1854 1855 # Skip deselected residues. 1856 if not spin.select: 1857 return 1858 1859 # Loop over the residue specific parameters. 1860 for param in param_names: 1861 # Return the parameter array. 1862 if index == inc: 1863 setattr(spin, param + "_err", error) 1864 1865 # Increment. 1866 inc = inc + 1

1867 1868

1869 - def set_param_values(self, param=None, value=None, index=None, spin_id=None, error=False, force=True):

1870 """Set the model-free parameter values. 1871 1872 @keyword param: The parameter name list. 1873 @type param: list of str 1874 @keyword value: The parameter value list. 1875 @type value: list 1876 @keyword index: The index for parameters which are of the list-type. This is unused. 1877 @type index: None or int 1878 @keyword spin_id: The spin identification string, only used for spin specific parameters. 1879 @type spin_id: None or str 1880 @keyword error: A flag which if True will allow the parameter errors to be set instead of the values. 1881 @type error: bool 1882 @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. 1883 @type force: bool 1884 """ 1885 1886 # Checks. 1887 is_str_list(param, 'parameter name') 1888 1889 # Separate out the diffusion tensor parameters from the model-free parameters. 1890 diff_params = [] 1891 diff_vals = [] 1892 mf_params = [] 1893 mf_vals = [] 1894 for i in range(len(param)): 1895 # Diffusion tensor parameter. 1896 if self._PARAMS.scope(param[i]) == 'global': 1897 if error: 1898 diff_params.append(param[i] + '_err') 1899 else: 1900 diff_params.append(param[i]) 1901 diff_vals.append(value[i]) 1902 1903 # Model-free parameter. 1904 else: 1905 mf_params.append(param[i]) 1906 mf_vals.append(value[i]) 1907 1908 # Set the diffusion tensor parameters. 1909 if diff_params: 1910 diffusion_tensor.set(value=diff_vals, param=diff_params) 1911 1912 # Set the model-free parameters. 1913 for i in range(len(mf_params)): 1914 # Check if it is a model-free parameter. 1915 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'): 1916 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'))) 1917 1918 # The error object name. 1919 if error: 1920 mf_params[i] += '_err' 1921 1922 # Set the parameter. 1923 for spin in spin_loop(spin_id): 1924 setattr(spin, mf_params[i], mf_vals[i])

1925 1926

1927 - def set_selected_sim(self, select_sim, model_info=None):

1928 """Set all simulation selection flags. 1929 1930 @param select_sim: The selection flags. 1931 @type select_sim: bool 1932 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 1933 @type model_info: int 1934 """ 1935 1936 # Determine the model type. 1937 model_type = determine_model_type() 1938 1939 # Global model. 1940 if model_type == 'all' or model_type == 'diff': 1941 cdp.select_sim = select_sim 1942 1943 # Spin specific model. 1944 else: 1945 # Get the spin container. 1946 spin = return_spin_from_index(model_info) 1947 1948 # Skip if deselected. 1949 if not spin.select: 1950 return 1951 1952 # Set the simulation flags. 1953 spin.select_sim = deepcopy(select_sim)

1954 1955

1956 - def set_update(self, param, spin):

1957 """Function to update the other model-free parameters. 1958 1959 @param param: The name of the parameter which has been changed. 1960 @type param: str 1961 @param spin: The SpinContainer object. 1962 @type spin: SpinContainer 1963 """ 1964 1965 # S2f parameter. 1966 if param == 's2f': 1967 # Update S2 if S2s exists. 1968 if hasattr(spin, 's2s') and spin.s2s != None: 1969 spin.s2 = spin.s2f * spin.s2s 1970 1971 1972 # S2s parameter. 1973 if param == 's2s': 1974 # Update S2 if S2f exists. 1975 if hasattr(spin, 's2f') and spin.s2f != None: 1976 spin.s2 = spin.s2f * spin.s2s

1977 1978

1979 - def sim_init_values(self):

1980 """Initialise the Monte Carlo parameter values.""" 1981 1982 # Determine the model type. 1983 model_type = determine_model_type() 1984 1985 # Get the parameter object names. 1986 param_names = self.data_names(set='params', scope='spin') 1987 1988 # Get the minimisation statistic object names. 1989 min_names = self.data_names(set='min', scope='spin') 1990 1991 # List of diffusion tensor parameters. 1992 if model_type == 'diff' or model_type == 'all': 1993 # Spherical diffusion. 1994 if cdp.diff_tensor.type == 'sphere': 1995 diff_params = ['tm'] 1996 1997 # Spheroidal diffusion. 1998 elif cdp.diff_tensor.type == 'spheroid': 1999 diff_params = ['tm', 'Da', 'theta', 'phi'] 2000 2001 # Ellipsoidal diffusion. 2002 elif cdp.diff_tensor.type == 'ellipsoid': 2003 diff_params = ['tm', 'Da', 'Dr', 'alpha', 'beta', 'gamma'] 2004 2005 2006 # Test if Monte Carlo parameter values have already been set. 2007 ############################################################# 2008 2009 # Diffusion tensor parameters and non spin specific minimisation statistics. 2010 if model_type == 'diff' or model_type == 'all': 2011 # Loop over the parameters. 2012 for object_name in diff_params: 2013 # Name for the simulation object. 2014 sim_object_name = object_name + '_sim' 2015 2016 # Test if the simulation object already exists. 2017 if hasattr(cdp.diff_tensor, sim_object_name): 2018 raise RelaxError("Monte Carlo parameter values have already been set.") 2019 2020 # Loop over the minimisation stats objects. 2021 for object_name in min_names: 2022 # Name for the simulation object. 2023 sim_object_name = object_name + '_sim' 2024 2025 # Test if the simulation object already exists. 2026 if hasattr(cdp, sim_object_name): 2027 raise RelaxError("Monte Carlo parameter values have already been set.") 2028 2029 # Spin specific parameters. 2030 if model_type != 'diff': 2031 for spin in spin_loop(): 2032 # Skip deselected spins. 2033 if not spin.select: 2034 continue 2035 2036 # Loop over all the parameter names. 2037 for object_name in param_names: 2038 # Name for the simulation object. 2039 sim_object_name = object_name + '_sim' 2040 2041 # Test if the simulation object already exists. 2042 if hasattr(spin, sim_object_name): 2043 raise RelaxError("Monte Carlo parameter values have already been set.") 2044 2045 2046 # Set the Monte Carlo parameter values. 2047 ####################################### 2048 2049 # Loop over the global minimisation stats objects. 2050 for object_name in min_names: 2051 # Skip non-existent objects. 2052 if not hasattr(cdp, object_name): 2053 continue 2054 2055 # Name for the simulation object. 2056 sim_object_name = object_name + '_sim' 2057 2058 # Create the simulation object. 2059 setattr(cdp, sim_object_name, []) 2060 2061 # Get the simulation object. 2062 sim_object = getattr(cdp, sim_object_name) 2063 2064 # Loop over the simulations. 2065 for j in range(cdp.sim_number): 2066 # Get the object. 2067 object = getattr(cdp, object_name) 2068 2069 # Copy and append the data. 2070 sim_object.append(deepcopy(object)) 2071 2072 # Diffusion tensor parameters and non spin specific minimisation statistics. 2073 if model_type == 'diff' or model_type == 'all': 2074 # Set up the number of simulations. 2075 cdp.diff_tensor.set_sim_num(cdp.sim_number) 2076 2077 # Loop over the parameters, setting the initial simulation values to those of the parameter value. 2078 for object_name in diff_params: 2079 for j in range(cdp.sim_number): 2080 cdp.diff_tensor.set(param=object_name, value=deepcopy(getattr(cdp.diff_tensor, object_name)), category='sim', sim_index=j) 2081 2082 # Spin specific parameters. 2083 if model_type != 'diff': 2084 for spin in spin_loop(): 2085 # Skip deselected spins. 2086 if not spin.select: 2087 continue 2088 2089 # Loop over all the data names. 2090 for object_name in param_names: 2091 # Name for the simulation object. 2092 sim_object_name = object_name + '_sim' 2093 2094 # Create the simulation object. 2095 setattr(spin, sim_object_name, []) 2096 2097 # Get the simulation object. 2098 sim_object = getattr(spin, sim_object_name) 2099 2100 # Loop over the simulations. 2101 for j in range(cdp.sim_number): 2102 # Copy and append the data. 2103 sim_object.append(deepcopy(getattr(spin, object_name))) 2104 2105 # Loop over all the minimisation object names. 2106 for object_name in min_names: 2107 # Name for the simulation object. 2108 sim_object_name = object_name + '_sim' 2109 2110 # Create the simulation object. 2111 setattr(spin, sim_object_name, []) 2112 2113 # Get the simulation object. 2114 sim_object = getattr(spin, sim_object_name) 2115 2116 # Loop over the simulations. 2117 for j in range(cdp.sim_number): 2118 # Copy and append the data. 2119 sim_object.append(deepcopy(getattr(spin, object_name)))

2120 2121

2122 - def sim_return_chi2(self, index=None, model_info=None):

2123 """Return the simulation chi-squared values. 2124 2125 @keyword index: The optional simulation index. 2126 @type index: int 2127 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2128 @type model_info: int 2129 @return: The list of simulation chi-squared values. If the index is supplied, only a single value will be returned. 2130 @rtype: list of float or float 2131 """ 2132 2133 # Determine the model type. 2134 model_type = determine_model_type() 2135 2136 # Single instance. 2137 if model_type == 'all' or model_type == 'diff': 2138 return cdp.chi2_sim 2139 2140 # Multiple instances. 2141 else: 2142 # Get the spin container. 2143 spin = return_spin_from_index(model_info) 2144 2145 # Return the list. 2146 return spin.chi2_sim

2147 2148

2149 - def sim_return_param(self, index, model_info=None):

2150 """Return the array of simulation parameter values. 2151 2152 @param index: The index of the parameter to return the array of values for. 2153 @type index: int 2154 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2155 @type model_info: int 2156 @return: The array of simulation parameter values. 2157 @rtype: list of float 2158 """ 2159 2160 # Parameter increment counter. 2161 inc = 0 2162 2163 # Determine the model type. 2164 model_type = determine_model_type() 2165 2166 # Get the parameter object names. 2167 param_names = self.data_names(set='params', scope='spin') 2168 2169 2170 # Diffusion tensor parameters. 2171 ############################## 2172 2173 if model_type == 'diff' or model_type == 'all': 2174 # Spherical diffusion. 2175 if cdp.diff_tensor.type == 'sphere': 2176 # Return the parameter array. 2177 if index == 0: 2178 return cdp.diff_tensor.tm_sim 2179 2180 # Increment. 2181 inc = inc + 1 2182 2183 # Spheroidal diffusion. 2184 elif cdp.diff_tensor.type == 'spheroid': 2185 # Return the parameter array. 2186 if index == 0: 2187 return cdp.diff_tensor.tm_sim 2188 elif index == 1: 2189 return cdp.diff_tensor.Da_sim 2190 elif index == 2: 2191 return cdp.diff_tensor.theta_sim 2192 elif index == 3: 2193 return cdp.diff_tensor.phi_sim 2194 2195 # Increment. 2196 inc = inc + 4 2197 2198 # Ellipsoidal diffusion. 2199 elif cdp.diff_tensor.type == 'ellipsoid': 2200 # Return the parameter array. 2201 if index == 0: 2202 return cdp.diff_tensor.tm_sim 2203 elif index == 1: 2204 return cdp.diff_tensor.Da_sim 2205 elif index == 2: 2206 return cdp.diff_tensor.Dr_sim 2207 elif index == 3: 2208 return cdp.diff_tensor.alpha_sim 2209 elif index == 4: 2210 return cdp.diff_tensor.beta_sim 2211 elif index == 5: 2212 return cdp.diff_tensor.gamma_sim 2213 2214 # Increment. 2215 inc = inc + 6 2216 2217 2218 # Model-free parameters for the model type 'all'. 2219 ################################################# 2220 2221 if model_type == 'all': 2222 # Loop over the spins. 2223 for spin in spin_loop(): 2224 # Skip deselected spins. 2225 if not spin.select: 2226 continue 2227 2228 # Loop over the spin specific parameters. 2229 for param in param_names: 2230 # Return the parameter array. 2231 if index == inc: 2232 return getattr(spin, param + "_sim") 2233 2234 # Increment. 2235 inc = inc + 1 2236 2237 2238 # Model-free parameters for the model types 'mf' and 'local_tm'. 2239 ################################################################ 2240 2241 if model_type == 'mf' or model_type == 'local_tm': 2242 # Get the spin container. 2243 spin = return_spin_from_index(model_info) 2244 2245 # Skip deselected spins. 2246 if not spin.select: 2247 return 2248 2249 # Loop over the spin specific parameters. 2250 for param in param_names: 2251 # Return the parameter array. 2252 if index == inc: 2253 return getattr(spin, param + "_sim") 2254 2255 # Increment. 2256 inc = inc + 1

2257 2258

2259 - def sim_return_selected(self, model_info=None):

2260 """Return the array of selected simulation flags for the spin. 2261 2262 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2263 @type model_info: int 2264 @return: The array of selected simulation flags. 2265 @rtype: list of int 2266 """ 2267 2268 # Determine the model type. 2269 model_type = determine_model_type() 2270 2271 # Single instance. 2272 if model_type == 'all' or model_type == 'diff': 2273 return cdp.select_sim 2274 2275 # Multiple instances. 2276 else: 2277 # Get the spin container. 2278 spin = return_spin_from_index(model_info) 2279 2280 # Skip if deselected. 2281 if not spin.select: 2282 return 2283 2284 # Return the list. 2285 return spin.select_sim

2286 2287

2288 - def skip_function(self, model_info=None):

2289 """Skip certain data. 2290 2291 @keyword model_info: The model information from model_loop(). This index is zero for the global models or equal to the global spin index (which covers the molecule, residue, and spin indices). 2292 @type model_info: int 2293 @return: True if the data should be skipped, False otherwise. 2294 @rtype: bool 2295 """ 2296 2297 # Determine the model type. 2298 model_type = determine_model_type() 2299 2300 # Sequence specific data. 2301 if (model_type == 'mf' or model_type == 'local_tm') and not return_spin_from_index(model_info).select: 2302 return True 2303 2304 # Don't skip. 2305 return False

Source Code for Module specific_analyses.model_free.api