specific_analyses.model

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

74 - def __init__(self):

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

93 94

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

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

134 135

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

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

174 175

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

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

412 413

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

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

600 601

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

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

636 637

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

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

671 672

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

674 """Deselect models or simulations. 675 676 @keyword sim_index: The optional Monte Carlo simulation index. If None, then models will be deselected, otherwise the given simulation will. 677 @type sim_index: None or int 678 @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). 679 @type model_info: int 680 """ 681 682 # Determine the model type. 683 model_type = determine_model_type() 684 685 # Local models. 686 if model_type == 'mf' or model_type == 'local_tm': 687 # Get the spin. 688 spin = return_spin_from_index(model_info) 689 690 # Spin deselection. 691 if sim_index == None: 692 spin.select = False 693 694 # Simulation deselection. 695 else: 696 spin.select_sim[sim_index] = False 697 698 # Global models. 699 else: 700 # Global model deselection. 701 if sim_index == None: 702 raise RelaxError("Cannot deselect the global model.") 703 704 # Simulation deselection. 705 else: 706 # Deselect. 707 cdp.select_sim[sim_index] = False

708 709

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

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

884 885

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

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

949 950

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

952 """Return a vector of parameter names. 953 954 @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). 955 @type model_info: int 956 @return: The vector of parameter names. 957 @rtype: list of str 958 """ 959 960 # Determine the model type. 961 model_type = determine_model_type() 962 963 # Get the spin ids. 964 if model_type == 'mf' or model_type == 'local_tm': 965 # Get the spin and it's id string. 966 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 967 else: 968 spin_id = None 969 970 # Assemble and return the parameter names. 971 return assemble_param_names(model_type, spin_id=spin_id)

972 973

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

975 """Return a vector of parameter values. 976 977 @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). 978 @type model_info: int 979 @keyword sim_index: The Monte Carlo simulation index. 980 @type sim_index: int 981 @return: The vector of parameter values. 982 @rtype: list of str 983 """ 984 985 # Test if the model-free models have been set up. 986 for spin in spin_loop(): 987 # Skip deselected spins. 988 if not spin.select: 989 continue 990 991 # Not setup. 992 if not spin.model: 993 raise RelaxNoModelError 994 995 # Determine the model type. 996 model_type = determine_model_type() 997 998 # Set the spin container (to None if the model is global). 999 if model_type == 'mf' or model_type == 'local_tm': 1000 spin = return_spin_from_index(model_info) 1001 else: 1002 spin = None 1003 1004 # Assemble the parameter values and return them. 1005 return assemble_param_vector(spin=spin, sim_index=sim_index, model_type=model_type)

1006 1007

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

1009 """The model-free grid search function. 1010 1011 @keyword lower: The per-model lower bounds of the grid search which must be equal to the number of parameters in the model. 1012 @type lower: list of lists of numbers 1013 @keyword upper: The per-model upper bounds of the grid search which must be equal to the number of parameters in the model. 1014 @type upper: list of lists of numbers 1015 @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. 1016 @type inc: list of lists of int 1017 @keyword scaling_matrix: The per-model list of diagonal and square scaling matrices. 1018 @type scaling_matrix: list of numpy rank-2, float64 array or list of None 1019 @keyword constraints: If True, constraints are applied during the grid search (eliminating parts of the grid). If False, no constraints are used. 1020 @type constraints: bool 1021 @keyword verbosity: A flag specifying the amount of information to print. The higher the value, the greater the verbosity. 1022 @type verbosity: int 1023 @keyword sim_index: The index of the simulation to apply the grid search to. If None, the normal model is optimised. 1024 @type sim_index: int 1025 """ 1026 1027 # Minimisation. 1028 self.minimise(min_algor='grid', lower=lower, upper=upper, inc=inc, scaling_matrix=scaling_matrix, constraints=constraints, verbosity=verbosity, sim_index=sim_index)

1029 1030

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

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

1068 1069

1070 - 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):

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

1396 1397

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

1399 """Return a description of the model. 1400 1401 @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). 1402 @type model_info: int 1403 @return: The model description. 1404 @rtype: str 1405 """ 1406 1407 # Determine the model type. 1408 model_type = determine_model_type() 1409 1410 # Global models. 1411 if model_type == 'all': 1412 return "Global model - all diffusion tensor parameters and spin specific model-free parameters." 1413 elif model_type == 'diff': 1414 return "Diffusion tensor model." 1415 1416 # Spin specific model. 1417 else: 1418 # Get the spin container. 1419 spin, spin_id = return_spin_from_index(model_info, return_spin_id=True) 1420 1421 # Return the description. 1422 return "Model-free model of spin '%s'." % spin_id

1423 1424

1425 - def model_loop(self):

1426 """Generator method for looping over the models (global or local). 1427 1428 If the model type is 'all' or 'diff', then this yields the single value of zero. Otherwise 1429 the global spin index is yielded. 1430 1431 1432 @return: The model index. This index is zero for the global models or equal to the global spin 1433 index (which covers the molecule, residue, and spin indices). 1434 @rtype: int 1435 """ 1436 1437 # Determine the model type. 1438 model_type = determine_model_type() 1439 1440 # Global model. 1441 if model_type == 'all' or model_type == 'diff': 1442 yield 0 1443 1444 # Spin specific models. 1445 else: 1446 # Loop over the spins. 1447 global_index = -1 1448 for spin in spin_loop(): 1449 # Increment the global spin index. 1450 global_index = global_index + 1 1451 1452 # Yield the spin index. 1453 yield global_index

1454 1455

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

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

1548 1549

1550 - def model_type(self):

1551 """Return the type of the model, either being 'local' or 'global'. 1552 1553 @return: The model type, one of 'local' or 'global'. 1554 @rtype: str 1555 """ 1556 1557 # Determine the model type. 1558 model_type = determine_model_type() 1559 1560 # Global models. 1561 if model_type in ['all', 'diff']: 1562 return 'global' 1563 1564 # Local models. 1565 else: 1566 return 'local'

1567 1568

1569 - def num_instances(self):

1570 """Function for returning the number of instances. 1571 1572 @return: The number of instances used for optimisation. Either the number of spins if 1573 the local optimisations are setup ('mf' and 'local_tm'), or 1 for the global 1574 models. 1575 @rtype: int 1576 """ 1577 1578 # Test if sequence data exists. 1579 if not exists_mol_res_spin_data(): 1580 return 0 1581 1582 # Determine the model type. 1583 model_type = determine_model_type() 1584 1585 # Sequence specific data. 1586 if model_type == 'mf' or model_type == 'local_tm': 1587 return count_spins() 1588 1589 # Other data. 1590 elif model_type == 'diff' or model_type == 'all': 1591 return 1 1592 1593 # Should not be here. 1594 else: 1595 raise RelaxFault

1596 1597

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

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

1721 1722

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

1724 """Print out the model title. 1725 1726 @keyword prefix: The starting text of the title. This should be printed out first, followed by the model information text. 1727 @type prefix: str 1728 @keyword model_info: The model information from model_loop(). 1729 @type model_info: unknown 1730 """ 1731 1732 # Determine the model type. 1733 model_type = determine_model_type() 1734 1735 # Local models. 1736 if model_type == 'mf' or model_type == 'local_tm': 1737 spin, spin_id = return_spin_from_index(global_index=model_info, return_spin_id=True) 1738 text = "%sSpin '%s'" % (prefix, spin_id) 1739 1740 # Global models. 1741 else: 1742 text = prefix + "Global model" 1743 1744 # The printout. 1745 subsection(file=sys.stdout, text=text, prespace=2)

1746 1747

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

1749 """Set the parameter errors. 1750 1751 @param index: The index of the parameter to set the errors for. 1752 @type index: int 1753 @param error: The error value. 1754 @type error: float 1755 @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). 1756 @type model_info: int 1757 """ 1758 1759 # Parameter increment counter. 1760 inc = 0 1761 1762 # Determine the model type. 1763 model_type = determine_model_type() 1764 1765 # Get the parameter object names. 1766 param_names = self.data_names(set='params', scope='spin') 1767 1768 1769 # Diffusion tensor parameter errors. 1770 #################################### 1771 1772 if model_type == 'diff' or model_type == 'all': 1773 # Spherical diffusion. 1774 if cdp.diff_tensor.type == 'sphere': 1775 # Return the parameter array. 1776 if index == 0: 1777 cdp.diff_tensor.set(param='tm', value=error, category='err') 1778 1779 # Increment. 1780 inc = inc + 1 1781 1782 # Spheroidal diffusion. 1783 elif cdp.diff_tensor.type == 'spheroid': 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='theta', value=error, category='err') 1791 elif index == 3: 1792 cdp.diff_tensor.set(param='phi', value=error, category='err') 1793 1794 # Increment. 1795 inc = inc + 4 1796 1797 # Ellipsoidal diffusion. 1798 elif cdp.diff_tensor.type == 'ellipsoid': 1799 # Return the parameter array. 1800 if index == 0: 1801 cdp.diff_tensor.set(param='tm', value=error, category='err') 1802 elif index == 1: 1803 cdp.diff_tensor.set(param='Da', value=error, category='err') 1804 elif index == 2: 1805 cdp.diff_tensor.set(param='Dr', value=error, category='err') 1806 elif index == 3: 1807 cdp.diff_tensor.set(param='alpha', value=error, category='err') 1808 elif index == 4: 1809 cdp.diff_tensor.set(param='beta', value=error, category='err') 1810 elif index == 5: 1811 cdp.diff_tensor.set(param='gamma', value=error, category='err') 1812 1813 # Increment. 1814 inc = inc + 6 1815 1816 1817 # Model-free parameter errors for the model type 'all'. 1818 ####################################################### 1819 1820 if model_type == 'all': 1821 # Loop over the spins. 1822 for spin in spin_loop(): 1823 # Skip deselected spins. 1824 if not spin.select: 1825 continue 1826 1827 # Loop over the residue specific parameters. 1828 for param in param_names: 1829 # Return the parameter array. 1830 if index == inc: 1831 setattr(spin, param + "_err", error) 1832 1833 # Increment. 1834 inc = inc + 1 1835 1836 1837 # Model-free parameters for the model types 'mf' and 'local_tm'. 1838 ################################################################ 1839 1840 if model_type == 'mf' or model_type == 'local_tm': 1841 # Get the spin container. 1842 spin = return_spin_from_index(model_info) 1843 1844 # Skip deselected residues. 1845 if not spin.select: 1846 return 1847 1848 # Loop over the residue specific parameters. 1849 for param in param_names: 1850 # Return the parameter array. 1851 if index == inc: 1852 setattr(spin, param + "_err", error) 1853 1854 # Increment. 1855 inc = inc + 1

1856 1857

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

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

1914 1915

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

1917 """Set all simulation selection flags. 1918 1919 @param select_sim: The selection flags. 1920 @type select_sim: bool 1921 @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). 1922 @type model_info: int 1923 """ 1924 1925 # Determine the model type. 1926 model_type = determine_model_type() 1927 1928 # Global model. 1929 if model_type == 'all' or model_type == 'diff': 1930 cdp.select_sim = select_sim 1931 1932 # Spin specific model. 1933 else: 1934 # Get the spin container. 1935 spin = return_spin_from_index(model_info) 1936 1937 # Skip if deselected. 1938 if not spin.select: 1939 return 1940 1941 # Set the simulation flags. 1942 spin.select_sim = deepcopy(select_sim)

1943 1944

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

1946 """Function to update the other model-free parameters. 1947 1948 @param param: The name of the parameter which has been changed. 1949 @type param: str 1950 @param spin: The SpinContainer object. 1951 @type spin: SpinContainer 1952 """ 1953 1954 # S2f parameter. 1955 if param == 's2f': 1956 # Update S2 if S2s exists. 1957 if hasattr(spin, 's2s') and spin.s2s != None: 1958 spin.s2 = spin.s2f * spin.s2s 1959 1960 1961 # S2s parameter. 1962 if param == 's2s': 1963 # Update S2 if S2f exists. 1964 if hasattr(spin, 's2f') and spin.s2f != None: 1965 spin.s2 = spin.s2f * spin.s2s

1966 1967

1968 - def sim_init_values(self):

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

2109 2110

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

2112 """Return the simulation chi-squared values. 2113 2114 @keyword index: The optional simulation index. 2115 @type index: int 2116 @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). 2117 @type model_info: int 2118 @return: The list of simulation chi-squared values. If the index is supplied, only a single value will be returned. 2119 @rtype: list of float or float 2120 """ 2121 2122 # Determine the model type. 2123 model_type = determine_model_type() 2124 2125 # Single instance. 2126 if model_type == 'all' or model_type == 'diff': 2127 return cdp.chi2_sim 2128 2129 # Multiple instances. 2130 else: 2131 # Get the spin container. 2132 spin = return_spin_from_index(model_info) 2133 2134 # Return the list. 2135 return spin.chi2_sim

2136 2137

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

2139 """Return the array of simulation parameter values. 2140 2141 @param index: The index of the parameter to return the array of values for. 2142 @type index: int 2143 @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). 2144 @type model_info: int 2145 @return: The array of simulation parameter values. 2146 @rtype: list of float 2147 """ 2148 2149 # Parameter increment counter. 2150 inc = 0 2151 2152 # Determine the model type. 2153 model_type = determine_model_type() 2154 2155 # Get the parameter object names. 2156 param_names = self.data_names(set='params', scope='spin') 2157 2158 2159 # Diffusion tensor parameters. 2160 ############################## 2161 2162 if model_type == 'diff' or model_type == 'all': 2163 # Spherical diffusion. 2164 if cdp.diff_tensor.type == 'sphere': 2165 # Return the parameter array. 2166 if index == 0: 2167 return cdp.diff_tensor.tm_sim 2168 2169 # Increment. 2170 inc = inc + 1 2171 2172 # Spheroidal diffusion. 2173 elif cdp.diff_tensor.type == 'spheroid': 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.theta_sim 2181 elif index == 3: 2182 return cdp.diff_tensor.phi_sim 2183 2184 # Increment. 2185 inc = inc + 4 2186 2187 # Ellipsoidal diffusion. 2188 elif cdp.diff_tensor.type == 'ellipsoid': 2189 # Return the parameter array. 2190 if index == 0: 2191 return cdp.diff_tensor.tm_sim 2192 elif index == 1: 2193 return cdp.diff_tensor.Da_sim 2194 elif index == 2: 2195 return cdp.diff_tensor.Dr_sim 2196 elif index == 3: 2197 return cdp.diff_tensor.alpha_sim 2198 elif index == 4: 2199 return cdp.diff_tensor.beta_sim 2200 elif index == 5: 2201 return cdp.diff_tensor.gamma_sim 2202 2203 # Increment. 2204 inc = inc + 6 2205 2206 2207 # Model-free parameters for the model type 'all'. 2208 ################################################# 2209 2210 if model_type == 'all': 2211 # Loop over the spins. 2212 for spin in spin_loop(): 2213 # Skip deselected spins. 2214 if not spin.select: 2215 continue 2216 2217 # Loop over the spin specific parameters. 2218 for param in param_names: 2219 # Return the parameter array. 2220 if index == inc: 2221 return getattr(spin, param + "_sim") 2222 2223 # Increment. 2224 inc = inc + 1 2225 2226 2227 # Model-free parameters for the model types 'mf' and 'local_tm'. 2228 ################################################################ 2229 2230 if model_type == 'mf' or model_type == 'local_tm': 2231 # Get the spin container. 2232 spin = return_spin_from_index(model_info) 2233 2234 # Skip deselected spins. 2235 if not spin.select: 2236 return 2237 2238 # Loop over the spin specific parameters. 2239 for param in param_names: 2240 # Return the parameter array. 2241 if index == inc: 2242 return getattr(spin, param + "_sim") 2243 2244 # Increment. 2245 inc = inc + 1

2246 2247

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

2249 """Return the array of selected simulation flags for the spin. 2250 2251 @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). 2252 @type model_info: int 2253 @return: The array of selected simulation flags. 2254 @rtype: list of int 2255 """ 2256 2257 # Determine the model type. 2258 model_type = determine_model_type() 2259 2260 # Single instance. 2261 if model_type == 'all' or model_type == 'diff': 2262 return cdp.select_sim 2263 2264 # Multiple instances. 2265 else: 2266 # Get the spin container. 2267 spin = return_spin_from_index(model_info) 2268 2269 # Skip if deselected. 2270 if not spin.select: 2271 return 2272 2273 # Return the list. 2274 return spin.select_sim

2275 2276

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

2278 """Skip certain data. 2279 2280 @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). 2281 @type model_info: int 2282 @return: True if the data should be skipped, False otherwise. 2283 @rtype: bool 2284 """ 2285 2286 # Determine the model type. 2287 model_type = determine_model_type() 2288 2289 # Sequence specific data. 2290 if (model_type == 'mf' or model_type == 'local_tm') and not return_spin_from_index(model_info).select: 2291 return True 2292 2293 # Don't skip. 2294 return False

Source Code for Module specific_analyses.model_free.api