specific_analyses.relax_disp.parameters

1 ############################################################################### 2 # # 3 # Copyright (C) 2004-2014 Edward d'Auvergne # 4 # Copyright (C) 2009 Sebastien Morin # 5 # # 6 # This file is part of the program relax (http://www.nmr-relax.com). # 7 # # 8 # This program is free software: you can redistribute it and/or modify # 9 # it under the terms of the GNU General Public License as published by # 10 # the Free Software Foundation, either version 3 of the License, or # 11 # (at your option) any later version. # 12 # # 13 # This program is distributed in the hope that it will be useful, # 14 # but WITHOUT ANY WARRANTY; without even the implied warranty of # 15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # 16 # GNU General Public License for more details. # 17 # # 18 # You should have received a copy of the GNU General Public License # 19 # along with this program. If not, see <http://www.gnu.org/licenses/>. # 20 # # 21 ############################################################################### 22 23 # Module docstring. 24 """Functions relating to the parameters of the relaxation dispersion models.""" 25 26 # Python module imports. 27 from copy import deepcopy 28 from numpy import array, float64, identity, zeros 29 import sys 30 31 # relax module imports. 32 from lib.errors import RelaxError, RelaxNoSequenceError 33 from lib.text.sectioning import subsection 34 from pipe_control import pipes 35 from pipe_control.mol_res_spin import exists_mol_res_spin_data, return_spin 36 from specific_analyses.relax_disp.disp_data import count_spins, generate_r20_key, has_exponential_exp_type, loop_cluster, loop_exp_frq, return_value_from_frq_index 37 from specific_analyses.relax_disp.variables import MODEL_LIST_MMQ, MODEL_M61B, MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR, MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR 38 39

40 -def assemble_param_vector(spins=None, key=None, sim_index=None):

41 """Assemble the dispersion relaxation dispersion curve fitting parameter vector. 42 43 @keyword spins: The list of spin data containers for the block. 44 @type spins: list of SpinContainer instances 45 @keyword key: The key for the R2eff and I0 parameters. 46 @type key: str or None 47 @keyword sim_index: The optional MC simulation index. 48 @type sim_index: int 49 @return: An array of the parameter values of the dispersion relaxation model. 50 @rtype: numpy float array 51 """ 52 53 # Initialise. 54 param_vector = [] 55 56 # Loop over the parameters of the cluster. 57 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 58 # Get the value. 59 value = get_value(key=key, spins=spins, sim_index=sim_index, param_name=param_name, spin_index=spin_index, r20_key=r20_key) 60 61 # Add to the vector. 62 param_vector.append(value) 63 64 # Convert all None values to 0.0. 65 for i in range(len(param_vector)): 66 if param_vector[i] == None: 67 param_vector[i] = 0.0 68 69 # Return a numpy array. 70 return array(param_vector, float64)

71 72

73 -def assemble_scaling_matrix(spins=None, key=None, scaling=True):

74 """Create and return the scaling matrix. 75 76 @keyword spins: The list of spin data containers for the block. 77 @type spins: list of SpinContainer instances 78 @keyword key: The key for the R2eff and I0 parameters. 79 @type key: str or None 80 @keyword scaling: A flag which if False will cause the identity matrix to be returned. 81 @type scaling: bool 82 @return: The diagonal and square scaling matrix. 83 @rtype: numpy diagonal matrix 84 """ 85 86 # Initialise. 87 scaling_matrix = identity(param_num(spins=spins), float64) 88 i = 0 89 param_index = 0 90 91 # No diagonal scaling. 92 if not scaling: 93 return scaling_matrix 94 95 # Loop over the parameters of the cluster. 96 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 97 # Transversal relaxation rate scaling. 98 if param_name in ['r2', 'r2a', 'r2b']: 99 scaling_matrix[param_index, param_index] = 10 100 101 # The pA.pB.dw**2, phi_ex_B, phi_ex_C and pA.dw**2 parameters. 102 elif param_name in ['phi_ex', 'phi_ex_B', 'phi_ex_C', 'padw2']: 103 scaling_matrix[param_index, param_index] = 1 104 105 # Chemical shift difference between states A and B scaling. 106 elif param_name in ['dw', 'dw_AB', 'dw_AC', 'dw_BC', 'dwH', 'dwH_AB', 'dwH_AC', 'dwH_BC']: 107 scaling_matrix[param_index, param_index] = 1 108 109 # The population of state X. 110 elif param_name in ['pA', 'pB', 'pC']: 111 scaling_matrix[param_index, param_index] = 1 112 113 # Exchange rate scaling. 114 elif param_name in ['kex', 'kex_AB', 'kex_AC', 'kex_BC', 'k_AB', 'kB', 'kC']: 115 scaling_matrix[param_index, param_index] = 10000 116 117 # Time of exchange scaling. 118 elif param_name == 'tex': 119 scaling_matrix[param_index, param_index] = 1e-4 120 121 # Return the scaling matrix. 122 return scaling_matrix

123 124

125 -def copy(pipe_from=None, pipe_to=None):

126 """Copy dispersion parameters from one data pipe to another, averaging values for clusters. 127 128 @param pipe_from: The data pipe to copy the value from. This defaults to the current data pipe. 129 @type pipe_from: str 130 @param pipe_to: The data pipe to copy the value to. This defaults to the current data pipe. 131 @type pipe_to: str 132 """ 133 134 # The current data pipe. 135 pipe_orig = pipes.cdp_name() 136 if pipe_from == None: 137 pipe_from = pipe_orig 138 if pipe_to == None: 139 pipe_to = pipe_orig 140 141 # Test that the pipes exist. 142 pipes.test(pipe_from) 143 pipes.test(pipe_to) 144 145 # Test that the pipes are not the same. 146 if pipe_from == pipe_to: 147 raise RelaxError("The source and destination pipes cannot be the same.") 148 149 # Test if the sequence data for pipe_from is loaded. 150 if not exists_mol_res_spin_data(pipe_from): 151 raise RelaxNoSequenceError(pipe_from) 152 153 # Test if the sequence data for pipe_to is loaded. 154 if not exists_mol_res_spin_data(pipe_to): 155 raise RelaxNoSequenceError(pipe_to) 156 157 # Switch to the destination data pipe. 158 pipes.switch(pipe_to) 159 160 # Loop over the clusters. 161 for spin_ids in loop_cluster(): 162 # Initialise some variables. 163 model = None 164 pA = 0.0 165 pB = 0.0 166 pC = 0.0 167 kex = 0.0 168 kex_AB = 0.0 169 kex_AC = 0.0 170 kex_BC = 0.0 171 k_AB = 0.0 172 kB = 0.0 173 kC = 0.0 174 tex = 0.0 175 count = 0 176 spins_from = [] 177 spins_to = [] 178 selected_cluster = False 179 180 # Loop over the spins, summing the parameters to be averaged. 181 for id in spin_ids: 182 # Get the spins, then store them. 183 spin_from = return_spin(id, pipe=pipe_from) 184 spin_to = return_spin(id, pipe=pipe_to) 185 spins_from.append(spin_from) 186 spins_to.append(spin_to) 187 188 # Skip deselected spins. 189 if not spin_from.select or not spin_to.select: 190 continue 191 192 # The first printout. 193 if not selected_cluster: 194 subsection(file=sys.stdout, text="Copying parameters for the spin block %s"%spin_ids, prespace=2) 195 196 # Change the cluster selection flag. 197 selected_cluster = True 198 199 # The model. 200 if not model: 201 model = spin_from.model 202 203 # Check that the models match for all spins of the cluster. 204 if spin_from.model != model: 205 raise RelaxError("The model '%s' of spin '%s' from the source data pipe does not match the '%s' model of previous spins of the cluster." % (spin_from.model, id, model)) 206 if spin_to.model != model: 207 raise RelaxError("The model '%s' of spin '%s' from the destination data pipe does not match the '%s' model of previous spins of the cluster." % (spin_from.model, id, model)) 208 209 # Sum the source parameters. 210 if 'pA' in spin_from.params: 211 pA += spin_from.pA 212 if 'pB' in spin_from.params: 213 pB += spin_from.pB 214 if 'pC' in spin_from.params: 215 pC += spin_from.pC 216 if 'kex' in spin_from.params: 217 kex += spin_from.kex 218 if 'kex_AB' in spin_from.params: 219 kex_AB += spin_from.kex_AB 220 if 'kex_AC' in spin_from.params: 221 kex_AC += spin_from.kex_AC 222 if 'kex_BC' in spin_from.params: 223 kex_BC += spin_from.kex_BC 224 if 'k_AB' in spin_from.params: 225 k_AB += spin_from.k_AB 226 if 'kB' in spin_from.params: 227 kB += spin_from.kB 228 if 'kC' in spin_from.params: 229 kC += spin_from.kC 230 if 'tex' in spin_from.params: 231 tex += spin_from.tex 232 233 # Increment the spin count. 234 count += 1 235 236 # The cluster is not selected, so move to the next. 237 if not selected_cluster: 238 continue 239 240 # Average parameters. 241 if pA != 0.0: 242 pA = pA / count 243 print("Averaged pA value: %.15f" % pA) 244 if pB != 0.0: 245 pB = pB / count 246 print("Averaged pA value: %.15f" % pA) 247 if pC != 0.0: 248 pC = pC / count 249 print("Averaged pC value: %.15f" % pC) 250 if kex != 0.0: 251 kex = kex / count 252 print("Averaged kex value: %.15f" % kex) 253 if kex_AB != 0.0: 254 kex_AB = kex_AB / count 255 print("Averaged k_AB value: %.15f" % kex_AB) 256 if kex_AC != 0.0: 257 kex_AC = kex_AC / count 258 print("Averaged k_AC value: %.15f" % kex_AC) 259 if kex_BC != 0.0: 260 kex_BC = kex_BC / count 261 print("Averaged k_BC value: %.15f" % kex_BC) 262 if k_AB != 0.0: 263 k_AB = k_AB / count 264 print("Averaged k_AB value: %.15f" % k_AB) 265 if kB != 0.0: 266 kB = kB / count 267 print("Averaged kB value: %.15f" % kB) 268 if kC != 0.0: 269 kC = kC / count 270 print("Averaged kC value: %.15f" % kC) 271 if tex != 0.0: 272 tex = tex / count 273 print("Averaged tex value: %.15f" % tex) 274 275 # Loop over the spins, this time copying the parameters. 276 for i in range(len(spin_ids)): 277 # Alias the containers. 278 spin_from = spins_from[i] 279 spin_to = spins_to[i] 280 281 # Skip deselected spins. 282 if not spin_from.select or not spin_to.select: 283 continue 284 285 # The R20 parameters. 286 if 'r2' in spin_from.params: 287 spin_to.r2 = deepcopy(spin_from.r2) 288 if 'r2a' in spin_from.params: 289 spin_to.r2a = deepcopy(spin_from.r2a) 290 if 'r2b' in spin_from.params: 291 spin_to.r2b = deepcopy(spin_from.r2b) 292 293 # The averaged parameters. 294 if 'pB' in spin_from.params and 'pC' not in spin_from.params: 295 spin_to.pA = pA 296 spin_to.pB = pB 297 spin_to.pC = 1.0 - pA - pB 298 elif 'pA' in spin_from.params: 299 spin_to.pA = pA 300 spin_to.pB = 1.0 - pA 301 if 'kex' in spin_from.params: 302 spin_to.kex = kex 303 if 'kex_AB' in spin_from.params: 304 spin_to.kex_AB = kex_AB 305 if 'kex_AC' in spin_from.params: 306 spin_to.kex_AC = kex_AC 307 if 'kex_BC' in spin_from.params: 308 spin_to.kex_BC = kex_BC 309 if 'k_AB' in spin_from.params: 310 spin_to.k_AB = k_AB 311 if 'kB' in spin_from.params: 312 spin_to.kB = kB 313 if 'kC' in spin_from.params: 314 spin_to.kC = kC 315 if 'tex' in spin_from.params: 316 spin_to.tex = tex 317 318 # All other spin specific parameters. 319 for param in spin_from.params: 320 if param in ['r2', 'pA', 'pB', 'pC', 'kex', 'kex_AB', 'kex_AC', 'kex_BC', 'k_AB', 'kB', 'kC', 'tex']: 321 continue 322 323 # Copy the value. 324 setattr(spin_to, param, deepcopy(getattr(spin_from, param))) 325 326 # Switch back to the original data pipe. 327 pipes.switch(pipe_orig)

328 329

330 -def disassemble_param_vector(param_vector=None, key=None, spins=None, sim_index=None):

331 """Disassemble the parameter vector. 332 333 @keyword param_vector: The parameter vector. 334 @type param_vector: numpy array 335 @keyword key: The key for the R2eff and I0 parameters. 336 @type key: str or None 337 @keyword spins: The list of spin data containers for the block. 338 @type spins: list of SpinContainer instances 339 @keyword sim_index: The optional MC simulation index. 340 @type sim_index: int 341 """ 342 343 # Initialise parameters if needed. 344 for spin in spins: 345 # Skip deselected spins. 346 if not spin.select: 347 continue 348 349 # The R2 parameter. 350 if 'r2' in spin.params: 351 if sim_index != None: 352 spin.r2_sim[sim_index] = {} 353 else: 354 spin.r2 = {} 355 356 # The R2A parameter. 357 if 'r2a' in spin.params: 358 if sim_index != None: 359 spin.r2a_sim[sim_index] = {} 360 else: 361 spin.r2a = {} 362 363 # The R2B parameter. 364 if 'r2b' in spin.params: 365 if sim_index != None: 366 spin.r2b_sim[sim_index] = {} 367 else: 368 spin.r2b = {} 369 370 # Loop over the parameters of the cluster, setting the values. 371 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 372 set_value(value=param_vector[param_index], key=key, spins=spins, sim_index=sim_index, param_name=param_name, spin_index=spin_index, r20_key=r20_key)

373 374

375 -def get_param_names(spins=None):

376 """Generate a list of dispersion parameter names for the given spins. 377 378 @keyword spins: The list of spin data containers for the block. 379 @type spins: list of SpinContainer instances 380 """ 381 382 # Initialise the structure. 383 names = [] 384 385 # Loop over the parameters. 386 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 387 # Set the initial text. 388 param_text = param_name 389 390 # The parameters with additional details. 391 if param_name in ['r2', 'r2a', 'r2b']: 392 param_text += " (%s)" % r20_key 393 394 # Append the text. 395 names.append(param_text) 396 397 # Return the structure. 398 return names

399 400

401 -def get_value(key=None, spins=None, sim_index=None, param_name=None, spin_index=None, r20_key=None):

402 """Return the value for the given parameter. 403 404 @keyword key: The key for the R2eff and I0 parameters. 405 @type key: str or None 406 @keyword spins: The list of spin data containers for the block. 407 @type spins: list of SpinContainer instances 408 @keyword sim_index: The optional MC simulation index. 409 @type sim_index: int 410 @keyword param_name: The parameter name. 411 @type param_name: str 412 @keyword spin_index: The spin index (for the cluster). 413 @type spin_index: int 414 @keyword r20_key: The unique R20 parameter key. 415 @type r20_key: str 416 @return: The parameter value. 417 @rtype: float 418 """ 419 420 # Default value of 0.0. 421 value = 0.0 422 423 # Spin specific parameters. 424 if spin_index != None: 425 # Get the simulation value. 426 if sim_index != None: 427 # Get the simulation object. 428 sim_obj = getattr(spins[spin_index], param_name+'_sim') 429 430 # R20 parameter. 431 if r20_key != None: 432 if r20_key in sim_obj[sim_index].keys(): 433 value = sim_obj[sim_index][r20_key] 434 435 # All other parameters. 436 else: 437 value = sim_obj[sim_index] 438 439 # Get the normal value. 440 else: 441 # Get the object. 442 obj = getattr(spins[spin_index], param_name) 443 444 # R20 parameter. 445 if r20_key != None: 446 if r20_key in obj.keys(): 447 value = obj[r20_key] 448 449 # All other parameters. 450 else: 451 value = getattr(spins[spin_index], param_name) 452 453 # Cluster specific parameters - use the parameter value from the first spin. 454 else: 455 # Set the simulation value. 456 if sim_index != None: 457 value = getattr(spins[0], param_name+'_sim')[sim_index] 458 459 # Set the normal value. 460 else: 461 value = getattr(spins[0], param_name) 462 463 # The R2eff model parameters. 464 if key != None: 465 if not key in value.keys(): 466 value = 0.0 467 else: 468 value = value[key] 469 470 # Return the value. 471 return value

472 473

474 -def linear_constraints(spins=None, scaling_matrix=None):

475 """Set up the relaxation dispersion curve fitting linear constraint matrices A and b. 476 477 Standard notation 478 ================= 479 480 The different constraints used within different models are:: 481 482 0 <= R2 <= 200 483 0 <= R2A <= 200 484 0 <= R2B <= 200 485 pC <= pB <= pA <= 1 486 pA >= 0.85 (the skewed condition, pA >> pB) 487 pB >= 0 488 pC >= 0 489 phi_ex >= 0 490 phi_ex_B >= 0 491 phi_ex_C >= 0 492 padw2 >= 0 493 dw >= 0 494 0 <= kex <= 2e6 495 0 <= k_AB <= 2e6 496 0 <= kB <= 2e6 497 0 <= kC <= 2e6 498 tex >= 0 499 k_AB >= 0 500 501 502 Matrix notation 503 =============== 504 505 In the notation A.x >= b, where A is a matrix of coefficients, x is an array of parameter values, and b is a vector of scalars, these inequality constraints are:: 506 507 | 1 0 0 | | R2 | | 0 | 508 | | | | | | 509 |-1 0 0 | | R2 | | -200 | 510 | | | | | | 511 | 1 0 0 | | R2A | | 0 | 512 | | | | | | 513 |-1 0 0 | | R2A | | -200 | 514 | | | | | | 515 | 1 0 0 | | R2B | | 0 | 516 | | | | | | 517 |-1 0 0 | | R2B | | -200 | 518 | | | | | | 519 | 1 0 0 | | pA | | 0.5 | 520 | | | | | | 521 |-1 0 0 | | pA | | -1 | 522 | | | | | | 523 | 1 0 0 | | pA | | 0.85 | 524 | | | | | | 525 | 1 0 0 | | pB | | 0 | 526 | | | | | | 527 | 1 0 0 | | phi_ex | | 0 | 528 | | | | | | 529 | 1 0 0 | . | phi_ex_B | >= | 0 | 530 | | | | | | 531 | 1 0 0 | | phi_ex_C | | 0 | 532 | | | | | | 533 | 1 0 0 | | padw2 | | 0 | 534 | | | | | | 535 | 1 0 0 | | dw | | 0 | 536 | | | | | | 537 | 1 0 0 | | kex | | 0 | 538 | | | | | | 539 |-1 0 0 | | kex | | -2e6 | 540 | | | | | | 541 | 1 0 0 | | kB | | 0 | 542 | | | | | | 543 |-1 0 0 | | kB | | -2e6 | 544 | | | | | | 545 | 1 0 0 | | kC | | 0 | 546 | | | | | | 547 |-1 0 0 | | kC | | -2e6 | 548 | | | | | | 549 | 1 0 0 | | tex | | 0 | 550 | | | | | | 551 | 1 0 0 | | k_AB | | 0 | 552 553 554 @keyword spins: The list of spin data containers for the block. 555 @type spins: list of SpinContainer instances 556 @keyword scaling_matrix: The diagonal, square scaling matrix. 557 @type scaling_matrix: numpy diagonal matrix 558 """ 559 560 # Initialisation (0..j..m). 561 A = [] 562 b = [] 563 n = param_num(spins=spins) 564 zero_array = zeros(n, float64) 565 j = 0 566 567 # Loop over the parameters of the cluster. 568 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 569 # Effective transversal relaxation rate. 570 if param_name == 'r2eff': 571 A.append(zero_array * 0.0) 572 A.append(zero_array * 0.0) 573 A[j][param_index] = 1.0 574 A[j+1][param_index] = -1.0 575 b.append(0.0) 576 b.append(-200.0 / scaling_matrix[param_index, param_index]) 577 j += 2 578 579 # Initial intensity. 580 elif param_name == 'i0': 581 A.append(zero_array * 0.0) 582 A[j][param_index] = 1.0 583 b.append(0.0) 584 j += 1 585 586 # The transversal relaxation rates (0 <= r2 <= 200). 587 elif param_name in ['r2', 'r2a', 'r2b']: 588 A.append(zero_array * 0.0) 589 A.append(zero_array * 0.0) 590 A[j][param_index] = 1.0 591 A[j+1][param_index] = -1.0 592 b.append(0.0) 593 b.append(-200.0 / scaling_matrix[param_index, param_index]) 594 j += 2 595 596 # The pA.pB.dw**2, phi_ex_B, phi_ex_C and pA.dw**2 parameters (phi_ex* >= 0 and padw2 >= 0). 597 elif param_name in ['phi_ex', 'phi_ex_B', 'phi_ex_C', 'padw2']: 598 A.append(zero_array * 0.0) 599 A[j][param_index] = 1.0 600 b.append(0.0) 601 j += 1 602 603 # Chemical exchange difference (dw >= 0). 604 elif param_name in ['dw', 'dw_AB', 'dw_AC', 'dw_BC', 'dwH', 'dwH_AB', 'dwH_AC', 'dwH_BC']: 605 if not spins[0].model in MODEL_LIST_MMQ + [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]: 606 A.append(zero_array * 0.0) 607 A[j][param_index] = 1.0 608 b.append(0.0) 609 j += 1 610 611 # The population of state A and B. 612 elif param_name == 'pA': 613 # First the pA <= 1 constraint (which also rearranged is pB >= 0). 614 A.append(zero_array * 0.0) 615 A[j][param_index] = -1.0 616 b.append(-1.0 / scaling_matrix[param_index, param_index]) 617 j += 1 618 619 # The skewed condition (pA >> pB). 620 if spins[0].model == MODEL_M61B: 621 A.append(zero_array * 0.0) 622 A[j][param_index] = 1.0 623 b.append(0.85 / scaling_matrix[param_index, param_index]) 624 j += 1 625 626 # Otherwise use the pA >= pB constraint. 627 else: 628 A.append(zero_array * 0.0) 629 A[j][param_index] = 1.0 630 b.append(0.5 / scaling_matrix[param_index, param_index]) 631 j += 1 632 633 # The population of state C (the pB parameter is only present when pC exists). 634 elif param_name == 'pB': 635 # Find pA. 636 for param_name2, param_index2, spin_index2, r20_key2 in loop_parameters(spins=spins): 637 if param_name2 == 'pA': 638 # First pC >= 0 (which rearranged is -pA - pB >= -1). 639 A.append(zero_array * 0.0) 640 A[j][param_index2] = -1.0 641 A[j][param_index] = -1.0 642 b.append(-1.0 / scaling_matrix[param_index, param_index]) 643 j += 1 644 645 # Then the pA >= pC constraint. 646 A.append(zero_array * 0.0) 647 A[j][param_index2] = 1.0 648 A[j][param_index] = -1.0 649 b.append(0.0) 650 j += 1 651 break 652 653 # Exchange rates and times (0 <= k <= 2e6). 654 elif param_name in ['kex', 'kex_AB', 'kex_AC', 'kex_BC', 'k_AB', 'kB', 'kC']: 655 A.append(zero_array * 0.0) 656 A.append(zero_array * 0.0) 657 A[j][param_index] = 1.0 658 A[j+1][param_index] = -1.0 659 b.append(0.0) 660 b.append(-2e6 / scaling_matrix[param_index, param_index]) 661 j += 2 662 663 # Exchange times (tex >= 0). 664 elif param_name in ['tex']: 665 A.append(zero_array * 0.0) 666 A[j][param_index] = 1.0 667 b.append(0.0) 668 j += 1 669 670 # Convert to numpy data structures. 671 A = array(A, float64) 672 b = array(b, float64) 673 674 # Return the matrices. 675 return A, b

676 677

678 -def loop_parameters(spins=None):

679 """Generator function for looping of the parameters of the cluster. 680 681 @keyword spins: The list of spin data containers for the block. 682 @type spins: list of SpinContainer instances 683 @return: The parameter name, the parameter index (for the parameter vector), the spin index (for the cluster), and the R20 parameter key (for R20, R20A, and R20B parameters stored as dictionaries). 684 @rtype: str, int, int, str 685 """ 686 687 # The parameter index. 688 param_index = -1 689 690 # The R2eff model. 691 if cdp.model_type == 'R2eff': 692 # Loop over the spins. 693 for spin_index in range(len(spins)): 694 # Skip deselected spins. 695 if not spins[spin_index].select: 696 continue 697 698 # Yield the two parameters. 699 params = ['r2eff', 'i0'] 700 for i in range(2): 701 # First increment the indices. 702 param_index += 1 703 704 # Yield the data. 705 yield params[i], param_index, spin_index, None 706 707 # All other models. 708 else: 709 # First the R2 parameters (one per spin per field strength). 710 for spin_index in range(len(spins)): 711 # Skip deselected spins. 712 if not spins[spin_index].select: 713 continue 714 715 # The R2 parameter. 716 if 'r2' in spins[0].params: 717 for exp_type, frq in loop_exp_frq(): 718 param_index += 1 719 yield 'r2', param_index, spin_index, generate_r20_key(exp_type=exp_type, frq=frq) 720 721 # The R2A parameter. 722 if 'r2a' in spins[0].params: 723 for exp_type, frq in loop_exp_frq(): 724 param_index += 1 725 yield 'r2a', param_index, spin_index, generate_r20_key(exp_type=exp_type, frq=frq) 726 727 # The R2B parameter. 728 if 'r2b' in spins[0].params: 729 for exp_type, frq in loop_exp_frq(): 730 param_index += 1 731 yield 'r2b', param_index, spin_index, generate_r20_key(exp_type=exp_type, frq=frq) 732 733 # Then the chemical shift difference parameters 'phi_ex', 'phi_ex_B', 'phi_ex_C', 'padw2', 'dw', 'dw_AB', 'dw_BC', 'dw_AB' (one per spin). 734 for spin_index in range(len(spins)): 735 # Skip deselected spins. 736 if not spins[spin_index].select: 737 continue 738 739 # Yield the data. 740 if 'phi_ex' in spins[spin_index].params: 741 param_index += 1 742 yield 'phi_ex', param_index, spin_index, None 743 if 'phi_ex_B' in spins[spin_index].params: 744 param_index += 1 745 yield 'phi_ex_B', param_index, spin_index, None 746 if 'phi_ex_C' in spins[spin_index].params: 747 param_index += 1 748 yield 'phi_ex_C', param_index, spin_index, None 749 if 'padw2' in spins[spin_index].params: 750 param_index += 1 751 yield 'padw2', param_index, spin_index, None 752 if 'dw' in spins[spin_index].params: 753 param_index += 1 754 yield 'dw', param_index, spin_index, None 755 if 'dw_AB' in spins[spin_index].params: 756 param_index += 1 757 yield 'dw_AB', param_index, spin_index, None 758 if 'dw_BC' in spins[spin_index].params: 759 param_index += 1 760 yield 'dw_BC', param_index, spin_index, None 761 if 'dw_AC' in spins[spin_index].params: 762 param_index += 1 763 yield 'dw_AC', param_index, spin_index, None 764 765 # Then a separate block for the proton chemical shift difference parameters for the MQ models (one per spin). 766 for spin_index in range(len(spins)): 767 # Skip deselected spins. 768 if not spins[spin_index].select: 769 continue 770 771 if 'dwH' in spins[spin_index].params: 772 param_index += 1 773 yield 'dwH', param_index, spin_index, None 774 if 'dwH_AB' in spins[spin_index].params: 775 param_index += 1 776 yield 'dwH_AB', param_index, spin_index, None 777 if 'dwH_BC' in spins[spin_index].params: 778 param_index += 1 779 yield 'dwH_BC', param_index, spin_index, None 780 if 'dwH_AC' in spins[spin_index].params: 781 param_index += 1 782 yield 'dwH_AC', param_index, spin_index, None 783 784 # All other parameters (one per spin cluster). 785 for param in spins[0].params: 786 if not param in ['r2', 'r2a', 'r2b', 'phi_ex', 'phi_ex_B', 'phi_ex_C', 'padw2', 'dw', 'dw_AB', 'dw_BC', 'dw_AB', 'dwH', 'dwH_AB', 'dwH_BC', 'dwH_AB']: 787 param_index += 1 788 yield param, param_index, None, None

789 790

791 -def param_conversion(key=None, spins=None, sim_index=None):

792 """Convert Disassemble the parameter vector. 793 794 @keyword key: The key for the R2eff and I0 parameters. 795 @type key: str or None 796 @keyword spins: The list of spin data containers for the block. 797 @type spins: list of SpinContainer instances 798 @keyword sim_index: The optional MC simulation index. 799 @type sim_index: int 800 """ 801 802 # Loop over the parameters of the cluster. 803 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 804 # Get the value. 805 value = get_value(key=key, spins=spins, sim_index=sim_index, param_name=param_name, spin_index=spin_index, r20_key=r20_key) 806 807 # The pA to pB to pC conversion. 808 if param_name == 'pA': 809 # 3-site exchange. 810 if spins[0].model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]: 811 # Get the pB value. 812 pB = get_value(key=key, spins=spins, sim_index=sim_index, param_name='pB', spin_index=spin_index, r20_key=r20_key) 813 814 # Set the pC value. 815 pC = 1.0 - value - pB 816 set_value(value=pC, key=key, spins=spins, sim_index=sim_index, param_name='pC', spin_index=spin_index) 817 818 # 2-site exchange. 819 else: 820 pB = 1.0 - value 821 set_value(value=pB, key=key, spins=spins, sim_index=sim_index, param_name='pB', spin_index=spin_index) 822 823 # The kex to tex conversion. 824 if param_name == 'kex': 825 tex = 1.0 / value 826 set_value(value=tex, key=key, spins=spins, sim_index=sim_index, param_name='tex', spin_index=spin_index) 827 828 # The kex to k_AB and k_BA conversion. 829 if param_name == 'kex' and 'pA' in spins[0].params: 830 # Get pA value. 831 pA = get_value(key=key, spins=spins, sim_index=sim_index, param_name='pA', spin_index=spin_index) 832 833 # Calculate k_AB value and set it. 834 k_AB = value * (1.0 - pA) 835 set_value(value=k_AB, key=key, spins=spins, sim_index=sim_index, param_name='k_AB', spin_index=spin_index) 836 837 # Calculate k_BA value and set it. 838 k_BA = value * pA 839 set_value(value=k_BA, key=key, spins=spins, sim_index=sim_index, param_name='k_BA', spin_index=spin_index) 840 841 # The tex to kex conversion. 842 if param_name == 'tex': 843 kex = 1.0 / value 844 set_value(value=kex, key=key, spins=spins, sim_index=sim_index, param_name='kex', spin_index=spin_index)

845 846

847 -def param_index_to_param_info(index=None, spins=None):

848 """Convert the given parameter array index to parameter identifying information. 849 850 The parameter index will be converted to the parameter name, the relevant spin index in the cluster, and relevant exponential curve key. 851 852 @keyword index: The index of the parameter array. 853 @type index: int 854 @keyword spins: The list of spin data containers for the block. 855 @type spins: list of SpinContainer instances 856 @return: The parameter name, the spin index (for the cluster), and the frequency index (for parameters with different values per spectrometer field strength). 857 @rtype: str, int, int 858 """ 859 860 # Loop over the parameters, yielding when a match is found. 861 for param_name, param_index, spin_index, r20_key in loop_parameters(spins=spins): 862 if param_index == index: 863 return param_name, spin_index, r20_key

864 865

866 -def param_num(spins=None):

867 """Determine the number of parameters in the model. 868 869 @keyword spins: The list of spin data containers for the block. 870 @type spins: list of SpinContainer instances 871 @return: The number of model parameters. 872 @rtype: int 873 """ 874 875 # Initialise the number. 876 num = 0 877 878 # The R2eff model. 879 if cdp.model_type == 'R2eff': 880 # Count the selected spins. 881 spin_num = count_spins(spins) 882 883 # Exponential curves (with clustering). 884 if has_exponential_exp_type(): 885 return 2 * spin_num 886 887 # Fixed time period experiments (with clustering). 888 return 1 * spin_num 889 890 # Check the spin cluster. 891 for spin in spins: 892 # Skip deselected spins. 893 if not spin.select: 894 continue 895 896 if len(spin.params) != len(spins[0].params): 897 raise RelaxError("The number of parameters for each spin in the cluster are not the same.") 898 899 # Count the number of R20 parameters. 900 r20_params = ['r2', 'r2a', 'r2b'] 901 for spin in spins: 902 # Skip deselected spins. 903 if not spin.select: 904 continue 905 906 for i in range(len(spin.params)): 907 if spin.params[i] in r20_params: 908 for exp_type, frq in loop_exp_frq(): 909 num += 1 910 911 # Count the number of spin specific parameters for all spins. 912 spin_params = ['phi_ex', 'phi_ex_B', 'phi_ex_C', 'padw2', 'dw', 'dwH'] 913 for spin in spins: 914 # Skip deselected spins. 915 if not spin.select: 916 continue 917 918 for i in range(len(spin.params)): 919 if spin.params[i] in spin_params: 920 num += 1 921 922 # Count all other parameters, but only for a single spin. 923 all_params = r20_params + spin_params 924 for spin in spins: 925 # Skip deselected spins. 926 if not spin.select: 927 continue 928 929 for i in range(len(spin.params)): 930 if not spin.params[i] in all_params: 931 num += 1 932 break 933 934 # Return the number. 935 return num

936 937

938 -def set_value(value=None, key=None, spins=None, sim_index=None, param_name=None, spin_index=None, r20_key=None):

939 """Return the value for the given parameter. 940 941 @keyword value: The parameter value to set. 942 @type value: float 943 @keyword key: The key for the R2eff and I0 parameters. 944 @type key: str or None 945 @keyword spins: The list of spin data containers for the block. 946 @type spins: list of SpinContainer instances 947 @keyword sim_index: The optional MC simulation index. 948 @type sim_index: int 949 @keyword param_name: The parameter name. 950 @type param_name: str 951 @keyword spin_index: The spin index (for the cluster). 952 @type spin_index: int 953 @keyword r20_key: The unique R20 parameter key. 954 @type r20_key: str 955 """ 956 957 # Spin specific parameters. 958 if spin_index != None: 959 # Set the simulation value. 960 if sim_index != None: 961 # Get the simulation object. 962 obj = getattr(spins[spin_index], param_name+'_sim') 963 964 # R20 parameter. 965 if r20_key != None: 966 obj[sim_index][r20_key] = value 967 968 # All other parameters. 969 else: 970 if key != None: 971 obj[sim_index][key] = value 972 else: 973 obj[sim_index] = value 974 975 # Set the normal value. 976 else: 977 # Get the object. 978 obj = getattr(spins[spin_index], param_name) 979 980 # R20 parameter. 981 if r20_key != None: 982 obj[r20_key] = value 983 984 # Set the normal value. 985 else: 986 if key != None: 987 obj[key] = value 988 else: 989 setattr(spins[spin_index], param_name, value) 990 991 # Cluster specific parameters. 992 else: 993 # Set the same parameter value for all spins in the cluster. 994 for spin in spins: 995 # Skip deselected spins. 996 if not spin.select: 997 continue 998 999 # Set the simulation value. 1000 if sim_index != None: 1001 sim_obj = getattr(spin, param_name+'_sim') 1002 sim_obj[sim_index] = value 1003 1004 # Set the normal value. 1005 else: 1006 setattr(spin, param_name, value)

1007

Source Code for Module specific_analyses.relax_disp.parameters