target_functions.relax

59 - def __init__(self, model=None, num_params=None, num_spins=None, num_frq=None, exp_types=None, values=None, errors=None, missing=None, frqs=None, frqs_H=None, cpmg_frqs=None, spin_lock_nu1=None, chemical_shifts=None, offset=None, tilt_angles=None, r1=None, relax_times=None, scaling_matrix=None, recalc_tau=True):

60 """Relaxation dispersion target functions for optimisation. 61 62 Models 63 ====== 64 65 The following analytic models are currently supported: 66 67 - 'No Rex': The model for no chemical exchange relaxation. 68 - 'LM63': The Luz and Meiboom (1963) 2-site fast exchange model. 69 - 'LM63 3-site': The Luz and Meiboom (1963) 3-site fast exchange model. 70 - 'CR72': The reduced Carver and Richards (1972) 2-site model for all time scales with R20A = R20B. 71 - 'CR72 full': The full Carver and Richards (1972) 2-site model for all time scales. 72 - 'IT99': The Ishima and Torchia (1999) 2-site model for all time scales with skewed populations (pA >> pB). 73 - 'TSMFK01': The Tollinger et al. (2001) 2-site very-slow exchange model, range of microsecond to second time scale. 74 - 'M61': The Meiboom (1961) 2-site fast exchange model for R1rho-type experiments. 75 - 'DPL94': The Davis, Perlman and London (1994) 2-site fast exchange model for R1rho-type experiments. 76 - 'M61 skew': The Meiboom (1961) on-resonance 2-site model with skewed populations (pA >> pB) for R1rho-type experiments. 77 - 'TP02': The Trott and Palmer (2002) 2-site exchange model for R1rho-type experiments. 78 - 'TAP03': The Trott, Abergel and Palmer (2003) 2-site exchange model for R1rho-type experiments. 79 - 'MP05': The Miloushev and Palmer (2005) off-resonance 2-site exchange model for R1rho-type experiments. 80 81 The following numerical models are currently supported: 82 83 - 'NS CPMG 2-site 3D': The reduced numerical solution for the 2-site Bloch-McConnell equations for CPMG data using 3D magnetisation vectors with R20A = R20B. 84 - 'NS CPMG 2-site 3D full': The full numerical solution for the 2-site Bloch-McConnell equations for CPMG data using 3D magnetisation vectors. 85 - 'NS CPMG 2-site star': The reduced numerical solution for the 2-site Bloch-McConnell equations for CPMG data using complex conjugate matrices with R20A = R20B. 86 - 'NS CPMG 2-site star full': The full numerical solution for the 2-site Bloch-McConnell equations for CPMG data using complex conjugate matrices. 87 - 'NS CPMG 2-site expanded': The numerical solution for the 2-site Bloch-McConnell equations for CPMG data expanded using Maple by Nikolai Skrynnikov. 88 - 'NS MMQ 2-site': The numerical solution for the 2-site Bloch-McConnell equations for combined proton-heteronuclear SQ, ZD, DQ, and MQ CPMG data with R20A = R20B. 89 - 'NS MMQ 3-site linear': The numerical solution for the 3-site Bloch-McConnell equations linearised with kAC = kCA = 0 for combined proton-heteronuclear SQ, ZD, DQ, and MQ CPMG data with R20A = R20B = R20C. 90 - 'NS MMQ 3-site': The numerical solution for the 3-site Bloch-McConnell equations for combined proton-heteronuclear SQ, ZD, DQ, and MQ CPMG data with R20A = R20B = R20C. 91 - 'NS R1rho 2-site': The numerical solution for the 2-site Bloch-McConnell equations for R1rho data with R20A = R20B. 92 - 'NS R1rho 3-site linear': The numerical solution for the 3-site Bloch-McConnell equations linearised with kAC = kCA = 0 for R1rho data with R20A = R20B = R20C. 93 - 'NS R1rho 3-site': The numerical solution for the 3-site Bloch-McConnell equations for R1rho data with R20A = R20B = R20C. 94 95 96 Indices 97 ======= 98 99 The data structures used in this target function class consist of many different index types. These are: 100 101 - Ei: The index for each experiment type. 102 - Si: The index for each spin of the spin cluster. 103 - Mi: The index for each magnetic field strength. 104 - Oi: The index for each spin-lock offset. In the case of CPMG-type data, this index is always zero. 105 - Di: The index for each dispersion point (either the spin-lock field strength or the nu_CPMG frequency). 106 107 108 @keyword model: The relaxation dispersion model to fit. 109 @type model: str 110 @keyword num_param: The number of parameters in the model. 111 @type num_param: int 112 @keyword num_spins: The number of spins in the cluster. 113 @type num_spins: int 114 @keyword num_frq: The number of spectrometer field strengths. 115 @type num_frq: int 116 @keyword exp_types: The list of experiment types. The dimensions are {Ei}. 117 @type exp_types: list of str 118 @keyword values: The R2eff/R1rho values. The dimensions are {Ei, Si, Mi, Oi, Di}. 119 @type values: rank-4 list of numpy rank-1 float arrays 120 @keyword errors: The R2eff/R1rho errors. The dimensions are {Ei, Si, Mi, Oi, Di}. 121 @type errors: rank-4 list of numpy rank-1 float arrays 122 @keyword missing: The data structure indicating missing R2eff/R1rho data. The dimensions are {Ei, Si, Mi, Oi, Di}. 123 @type missing: rank-4 list of numpy rank-1 int arrays 124 @keyword frqs: The spin Larmor frequencies (in MHz*2pi to speed up the ppm to rad/s conversion). The dimensions are {Ei, Si, Mi}. 125 @type frqs: rank-3 list of floats 126 @keyword frqs_H: The proton spin Larmor frequencies for the MMQ-type models (in MHz*2pi to speed up the ppm to rad/s conversion). The dimensions are {Ei, Si, Mi}. 127 @type frqs_H: rank-3 list of floats 128 @keyword cpmg_frqs: The CPMG frequencies in Hertz. This will be ignored for R1rho experiments. The dimensions are {Ei, Mi}. 129 @type cpmg_frqs: rank-2 list of floats 130 @keyword spin_lock_nu1: The spin-lock field strengths in Hertz. This will be ignored for CPMG experiments. The dimensions are {Ei, Mi}. 131 @type spin_lock_nu1: rank-2 list of floats 132 @keyword chemical_shifts: The chemical shifts in rad/s. This is only used for off-resonance R1rho models. The ppm values are not used to save computation time, therefore they must be converted to rad/s by the calling code. The dimensions are {Ei, Si, Mi}. 133 @type chemical_shifts: rank-3 list of floats 134 @keyword offset: The structure of spin-lock or hard pulse offsets in rad/s. This is only currently used for off-resonance R1rho models. The dimensions are {Ei, Si, Mi, Oi}. 135 @type offset: rank-4 list of floats 136 @keyword tilt_angles: The spin-lock rotating frame tilt angle. This is only used for off-resonance R1rho models. The dimensions are {Ei, Si, Mi, Oi, Di}. 137 @type tilt_angles: rank-5 list of floats 138 @keyword r1: The R1 relaxation rates. This is only used for off-resonance R1rho models. The dimensions are {Ei, Si, Mi}. 139 @type r1: rank-3 list of floats 140 @keyword relax_times: The experiment specific fixed time period for relaxation (in seconds). The dimensions are {Ei, Mi}. 141 @type relax_times: rank-2 list of floats 142 @keyword scaling_matrix: The square and diagonal scaling matrix. 143 @type scaling_matrix: numpy rank-2 float array 144 @keyword recalc_tau: A flag which if True will cause tau_CPMG to be recalculated to remove user input truncation. 145 @type recalc_tau: bool 146 """ 147 148 # Check the args. 149 if model not in MODEL_LIST_FULL: 150 raise RelaxError("The model '%s' is unknown." % model) 151 if values == None: 152 raise RelaxError("No values have been supplied to the target function.") 153 if errors == None: 154 raise RelaxError("No errors have been supplied to the target function.") 155 if missing == None: 156 raise RelaxError("No missing data information has been supplied to the target function.") 157 if model in [MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05]: 158 if chemical_shifts == None: 159 raise RelaxError("Chemical shifts must be supplied for the '%s' R1rho off-resonance dispersion model." % model) 160 if r1 == None: 161 raise RelaxError("R1 relaxation rates must be supplied for the '%s' R1rho off-resonance dispersion model." % model) 162 163 # Store the arguments. 164 self.model = model 165 self.num_params = num_params 166 self.num_spins = num_spins 167 self.num_frq = num_frq 168 self.exp_types = exp_types 169 self.values = values 170 self.errors = errors 171 self.missing = missing 172 self.frqs = frqs 173 self.frqs_H = frqs_H 174 self.cpmg_frqs = cpmg_frqs 175 self.spin_lock_nu1 = spin_lock_nu1 176 self.chemical_shifts = chemical_shifts 177 self.offset = offset 178 self.tilt_angles = tilt_angles 179 self.r1 = r1 180 self.relax_times = relax_times 181 self.scaling_matrix = scaling_matrix 182 183 # Create the structure for holding the back-calculated R2eff values (matching the dimensions of the values structure). 184 self.back_calc = deepcopy(values) 185 186 # Check the experiment types, simplifying the data structures as needed. 187 self.experiment_type_setup() 188 189 # Determine the number of dispersion points (for each experiment type, magnetic field strength, and offset) and offsets. 190 self.num_disp_points = [] 191 self.num_offsets = [] 192 for ei in range(len(values)): 193 self.num_disp_points.append([]) 194 self.num_offsets.append([]) 195 for si in range(len(values[ei])): 196 self.num_disp_points[ei].append([]) 197 self.num_offsets[ei].append([]) 198 for mi in range(len(values[ei][0])): 199 self.num_disp_points[ei][si].append([]) 200 self.num_offsets[ei][si].append([]) 201 202 # The offset data. 203 if len(offset[ei][si][mi]): 204 self.num_offsets[ei][si][mi] = len(self.offset[ei][si][mi]) 205 else: 206 self.num_offsets[ei][si][mi] = 0 207 208 # The dispersion points. 209 for oi in range(self.num_offsets[ei][si][mi]): 210 self.num_disp_points[ei][si][mi].append([]) 211 if cpmg_frqs != None and len(cpmg_frqs[ei][mi][oi]): 212 self.num_disp_points[ei][si][mi][oi] = len(self.cpmg_frqs[ei][mi][oi]) 213 elif spin_lock_nu1 != None and len(spin_lock_nu1[ei][mi][oi]): 214 self.num_disp_points[ei][si][mi][oi] = len(self.spin_lock_nu1[ei][mi][oi]) 215 else: 216 self.num_disp_points[ei][si][mi][oi] = 0 217 218 # Scaling initialisation. 219 self.scaling_flag = False 220 if self.scaling_matrix != None: 221 self.scaling_flag = True 222 223 # Initialise the post spin parameter indices. 224 self.end_index = [] 225 226 # The spin and frequency dependent R2 parameters. 227 self.end_index.append(self.num_exp * self.num_spins * self.num_frq) 228 if model in [MODEL_CR72_FULL, MODEL_NS_CPMG_2SITE_3D_FULL, MODEL_NS_CPMG_2SITE_STAR_FULL]: 229 self.end_index.append(2 * self.num_exp * self.num_spins * self.num_frq) 230 231 # The spin and dependent parameters (phi_ex, dw, padw2). 232 self.end_index.append(self.end_index[-1] + self.num_spins) 233 if model in [MODEL_IT99, MODEL_LM63_3SITE, MODEL_MMQ_CR72, MODEL_NS_MMQ_2SITE]: 234 self.end_index.append(self.end_index[-1] + self.num_spins) 235 elif model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]: 236 self.end_index.append(self.end_index[-1] + self.num_spins) 237 self.end_index.append(self.end_index[-1] + self.num_spins) 238 elif model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]: 239 self.end_index.append(self.end_index[-1] + self.num_spins) 240 self.end_index.append(self.end_index[-1] + self.num_spins) 241 self.end_index.append(self.end_index[-1] + self.num_spins) 242 243 # Set up the matrices for the numerical solutions. 244 if model in [MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL]: 245 # The matrix that contains only the R2 relaxation terms ("Redfield relaxation", i.e. non-exchange broadening). 246 self.Rr = zeros((2, 2), complex64) 247 248 # The matrix that contains the exchange terms between the two states A and B. 249 self.Rex = zeros((2, 2), complex64) 250 251 # The matrix that contains the chemical shift evolution. It works here only with X magnetization, and the complex notation allows to evolve in the transverse plane (x, y). 252 self.RCS = zeros((2, 2), complex64) 253 254 # The matrix that contains all the contributions to the evolution, i.e. relaxation, exchange and chemical shift evolution. 255 self.R = zeros((2, 2), complex64) 256 257 # Pi-pulse propagators. 258 if model in [MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL]: 259 self.r180x = r180x_3d() 260 261 # This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 262 if model in [MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL, MODEL_NS_MMQ_2SITE]: 263 self.M0 = zeros(2, float64) 264 if model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]: 265 self.M0 = zeros(3, float64) 266 if model in [MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL]: 267 self.M0 = zeros(7, float64) 268 self.M0[0] = 0.5 269 if model in [MODEL_NS_R1RHO_2SITE]: 270 self.M0 = zeros(6, float64) 271 if model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]: 272 self.M0 = zeros(9, float64) 273 274 # Special CPMG-type data structures. 275 if model in [MODEL_MMQ_CR72, MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL, MODEL_NS_CPMG_2SITE_EXPANDED, MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL, MODEL_NS_MMQ_2SITE, MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR, MODEL_TSMFK01]: 276 # The number of CPMG blocks. 277 self.power = [] 278 for ei in range(self.num_exp): 279 self.power.append([]) 280 for mi in range(self.num_frq): 281 self.power[ei].append(zeros(self.num_disp_points[ei][si][mi][0], int16)) 282 for di in range(self.num_disp_points[ei][si][mi][0]): 283 # Missing data for an entire field strength. 284 if isNaN(self.relax_times[ei][mi]): 285 self.power[ei][mi][di] = 0.0 286 287 # Normal value. 288 else: 289 self.power[ei][mi][di] = int(round(self.cpmg_frqs[ei][mi][0][di] * self.relax_times[ei][mi])) 290 291 # The tau_cpmg times. 292 self.tau_cpmg = [] 293 for ei in range(len(values)): 294 self.tau_cpmg.append([]) 295 for mi in range(len(values[ei][0])): 296 self.tau_cpmg[ei].append(zeros(self.num_disp_points[ei][si][mi][0], float64)) 297 for di in range(self.num_disp_points[ei][si][mi][0]): 298 # Recalculate the tau_cpmg times to avoid any user induced truncation in the input files. 299 if recalc_tau: 300 self.tau_cpmg[ei][mi][di] = 0.25 * self.relax_times[ei][mi] / self.power[ei][mi][di] 301 else: 302 self.tau_cpmg[ei][mi][di] = 0.25 / self.cpmg_frqs[ei][mi][0][di] 303 304 # Convert the spin-lock data to rad.s^-1. 305 if spin_lock_nu1 != None: 306 self.spin_lock_omega1 = [] 307 self.spin_lock_omega1_squared = [] 308 for ei in range(len(values)): 309 self.spin_lock_omega1.append([]) 310 self.spin_lock_omega1_squared.append([]) 311 for mi in range(len(values[ei][0])): 312 self.spin_lock_omega1[ei].append([]) 313 self.spin_lock_omega1_squared[ei].append([]) 314 for oi in range(len(values[ei][0][mi])): 315 self.spin_lock_omega1[ei][mi].append([]) 316 self.spin_lock_omega1_squared[ei][mi].append([]) 317 self.spin_lock_omega1[ei][mi][oi] = 2.0 * pi * self.spin_lock_nu1[ei][mi][oi] 318 self.spin_lock_omega1_squared[ei][mi][oi] = self.spin_lock_omega1[ei][mi][oi] ** 2 319 320 # The inverted relaxation delay. 321 if model in [MODEL_MMQ_CR72, MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL, MODEL_NS_CPMG_2SITE_EXPANDED, MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL, MODEL_NS_MMQ_2SITE, MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR, MODEL_NS_R1RHO_2SITE, MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]: 322 self.inv_relax_times = 1.0 / relax_times 323 324 # Special storage matrices for the multi-quantum CPMG N-site numerical models. 325 if model == MODEL_NS_MMQ_2SITE: 326 self.m1 = zeros((2, 2), complex64) 327 self.m2 = zeros((2, 2), complex64) 328 elif model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]: 329 self.m1 = zeros((3, 3), complex64) 330 self.m2 = zeros((3, 3), complex64) 331 elif model == MODEL_NS_R1RHO_2SITE: 332 self.matrix = zeros((6, 6), float64) 333 elif model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]: 334 self.matrix = zeros((9, 9), float64) 335 336 # Set up the model. 337 if model == MODEL_NOREX: 338 self.func = self.func_NOREX 339 if model == MODEL_LM63: 340 self.func = self.func_LM63 341 if model == MODEL_LM63_3SITE: 342 self.func = self.func_LM63_3site 343 if model == MODEL_CR72_FULL: 344 self.func = self.func_CR72_full 345 if model == MODEL_CR72: 346 self.func = self.func_CR72 347 if model == MODEL_IT99: 348 self.func = self.func_IT99 349 if model == MODEL_TSMFK01: 350 self.func = self.func_TSMFK01 351 if model == MODEL_NS_CPMG_2SITE_3D_FULL: 352 self.func = self.func_ns_cpmg_2site_3D_full 353 if model == MODEL_NS_CPMG_2SITE_3D: 354 self.func = self.func_ns_cpmg_2site_3D 355 if model == MODEL_NS_CPMG_2SITE_EXPANDED: 356 self.func = self.func_ns_cpmg_2site_expanded 357 if model == MODEL_NS_CPMG_2SITE_STAR_FULL: 358 self.func = self.func_ns_cpmg_2site_star_full 359 if model == MODEL_NS_CPMG_2SITE_STAR: 360 self.func = self.func_ns_cpmg_2site_star 361 if model == MODEL_M61: 362 self.func = self.func_M61 363 if model == MODEL_M61B: 364 self.func = self.func_M61b 365 if model == MODEL_DPL94: 366 self.func = self.func_DPL94 367 if model == MODEL_TP02: 368 self.func = self.func_TP02 369 if model == MODEL_TAP03: 370 self.func = self.func_TAP03 371 if model == MODEL_MP05: 372 self.func = self.func_MP05 373 if model == MODEL_NS_R1RHO_2SITE: 374 self.func = self.func_ns_r1rho_2site 375 if model == MODEL_NS_R1RHO_3SITE: 376 self.func = self.func_ns_r1rho_3site 377 if model == MODEL_NS_R1RHO_3SITE_LINEAR: 378 self.func = self.func_ns_r1rho_3site_linear 379 if model == MODEL_MMQ_CR72: 380 self.func = self.func_mmq_CR72 381 if model == MODEL_NS_MMQ_2SITE: 382 self.func = self.func_ns_mmq_2site 383 if model == MODEL_NS_MMQ_3SITE: 384 self.func = self.func_ns_mmq_3site 385 if model == MODEL_NS_MMQ_3SITE_LINEAR: 386 self.func = self.func_ns_mmq_3site_linear

387 388

389 - def calc_CR72_chi2(self, R20A=None, R20B=None, dw=None, pA=None, kex=None):

390 """Calculate the chi-squared value of the 'NS CPMG 2-site star' models. 391 392 @keyword R20A: The R2 value for state A in the absence of exchange. 393 @type R20A: list of float 394 @keyword R20B: The R2 value for state B in the absence of exchange. 395 @type R20B: list of float 396 @keyword dw: The chemical shift differences in ppm for each spin. 397 @type dw: list of float 398 @keyword pA: The population of state A. 399 @type pA: float 400 @keyword kex: The rate of exchange. 401 @type kex: float 402 @return: The chi-squared value. 403 @rtype: float 404 """ 405 406 # Initialise. 407 chi2_sum = 0.0 408 409 # Loop over the spins. 410 for si in range(self.num_spins): 411 # Loop over the spectrometer frequencies. 412 for mi in range(self.num_frq): 413 # The R20 index. 414 r20_index = mi + si*self.num_frq 415 416 # Convert dw from ppm to rad/s. 417 dw_frq = dw[si] * self.frqs[0][si][mi] 418 419 # Back calculate the R2eff values. 420 r2eff_CR72(r20a=R20A[r20_index], r20b=R20B[r20_index], pA=pA, dw=dw_frq, kex=kex, cpmg_frqs=self.cpmg_frqs[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 421 422 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 423 for di in range(self.num_disp_points[0][si][mi][0]): 424 if self.missing[0][si][mi][0][di]: 425 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 426 427 # Calculate and return the chi-squared value. 428 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 429 430 # Return the total chi-squared value. 431 return chi2_sum

432 433

434 - def calc_ns_cpmg_2site_3D_chi2(self, R20A=None, R20B=None, dw=None, pA=None, kex=None):

435 """Calculate the chi-squared value of the 'NS CPMG 2-site' models. 436 437 @keyword R20A: The R2 value for state A in the absence of exchange. 438 @type R20A: list of float 439 @keyword R20B: The R2 value for state B in the absence of exchange. 440 @type R20B: list of float 441 @keyword dw: The chemical shift differences in ppm for each spin. 442 @type dw: list of float 443 @keyword pA: The population of state A. 444 @type pA: float 445 @keyword kex: The rate of exchange. 446 @type kex: float 447 @return: The chi-squared value. 448 @rtype: float 449 """ 450 451 # Once off parameter conversions. 452 pB = 1.0 - pA 453 k_BA = pA * kex 454 k_AB = pB * kex 455 456 # This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 457 self.M0[1] = pA 458 self.M0[4] = pB 459 460 # Chi-squared initialisation. 461 chi2_sum = 0.0 462 463 # Loop over the spins. 464 for si in range(self.num_spins): 465 # Loop over the spectrometer frequencies. 466 for mi in range(self.num_frq): 467 # The R20 index. 468 r20_index = mi + si*self.num_frq 469 470 # Convert dw from ppm to rad/s. 471 dw_frq = dw[si] * self.frqs[0][si][mi] 472 473 # Back calculate the R2eff values. 474 r2eff_ns_cpmg_2site_3D(r180x=self.r180x, M0=self.M0, r20a=R20A[r20_index], r20b=R20B[r20_index], pA=pA, pB=pB, dw=dw_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_times[0][mi], tcp=self.tau_cpmg[0][mi], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0], power=self.power[0][mi]) 475 476 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 477 for di in range(self.num_disp_points[0][si][mi][0]): 478 if self.missing[0][si][mi][0][di]: 479 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 480 481 # Calculate and return the chi-squared value. 482 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 483 484 # Return the total chi-squared value. 485 return chi2_sum

486 487

488 - def calc_ns_cpmg_2site_star_chi2(self, R20A=None, R20B=None, dw=None, pA=None, kex=None):

489 """Calculate the chi-squared value of the 'NS CPMG 2-site star' models. 490 491 @keyword R20A: The R2 value for state A in the absence of exchange. 492 @type R20A: list of float 493 @keyword R20B: The R2 value for state B in the absence of exchange. 494 @type R20B: list of float 495 @keyword dw: The chemical shift differences in ppm for each spin. 496 @type dw: list of float 497 @keyword pA: The population of state A. 498 @type pA: float 499 @keyword kex: The rate of exchange. 500 @type kex: float 501 @return: The chi-squared value. 502 @rtype: float 503 """ 504 505 # Once off parameter conversions. 506 pB = 1.0 - pA 507 k_BA = pA * kex 508 k_AB = pB * kex 509 510 # Set up the matrix that contains the exchange terms between the two states A and B. 511 self.Rex[0, 0] = -k_AB 512 self.Rex[0, 1] = k_BA 513 self.Rex[1, 0] = k_AB 514 self.Rex[1, 1] = -k_BA 515 516 # This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 517 self.M0[0] = pA 518 self.M0[1] = pB 519 520 # Chi-squared initialisation. 521 chi2_sum = 0.0 522 523 # Loop over the spins. 524 for si in range(self.num_spins): 525 # Loop over the spectrometer frequencies. 526 for mi in range(self.num_frq): 527 # The R20 index. 528 r20_index = mi + si*self.num_frq 529 530 # Convert dw from ppm to rad/s. 531 dw_frq = dw[si] * self.frqs[0][si][mi] 532 533 # Back calculate the R2eff values. 534 r2eff_ns_cpmg_2site_star(Rr=self.Rr, Rex=self.Rex, RCS=self.RCS, R=self.R, M0=self.M0, r20a=R20A[r20_index], r20b=R20B[r20_index], dw=dw_frq, inv_tcpmg=self.inv_relax_times[0][mi], tcp=self.tau_cpmg[0][mi], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0], power=self.power[0][mi]) 535 536 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 537 for di in range(self.num_disp_points[0][si][mi][0]): 538 if self.missing[0][si][mi][0][di]: 539 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 540 541 # Calculate and return the chi-squared value. 542 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 543 544 # Return the total chi-squared value. 545 return chi2_sum

546 547

548 - def calc_ns_mmq_3site_chi2(self, R20A=None, R20B=None, R20C=None, dw_AB=None, dw_BC=None, dwH_AB=None, dwH_BC=None, pA=None, pB=None, kex_AB=None, kex_BC=None, kex_AC=None):

549 """Calculate the chi-squared value for the 'NS MMQ 3-site' models. 550 551 @keyword R20A: The R2 value for state A in the absence of exchange. 552 @type R20A: list of float 553 @keyword R20B: The R2 value for state B in the absence of exchange. 554 @type R20B: list of float 555 @keyword R20C: The R2 value for state C in the absence of exchange. 556 @type R20C: list of float 557 @keyword dw_AB: The chemical exchange difference between states A and B in rad/s. 558 @type dw_AB: float 559 @keyword dw_BC: The chemical exchange difference between states B and C in rad/s. 560 @type dw_BC: float 561 @keyword dwH_AB: The proton chemical exchange difference between states A and B in rad/s. 562 @type dwH_AB: float 563 @keyword dwH_BC: The proton chemical exchange difference between states B and C in rad/s. 564 @type dwH_BC: float 565 @keyword pA: The population of state A. 566 @type pA: float 567 @keyword kex_AB: The rate of exchange between states A and B. 568 @type kex_AB: float 569 @keyword kex_BC: The rate of exchange between states B and C. 570 @type kex_BC: float 571 @keyword kex_AC: The rate of exchange between states A and C. 572 @type kex_AC: float 573 @return: The chi-squared value. 574 @rtype: float 575 """ 576 577 # Once off parameter conversions. 578 pC = 1.0 - pA - pB 579 pA_pB = pA + pB 580 pA_pC = pA + pC 581 pB_pC = pB + pC 582 k_BA = pA * kex_AB / pA_pB 583 k_AB = pB * kex_AB / pA_pB 584 k_CB = pB * kex_BC / pB_pC 585 k_BC = pC * kex_BC / pB_pC 586 k_CA = pA * kex_AC / pA_pC 587 k_AC = pC * kex_AC / pA_pC 588 dw_AC = dw_AB + dw_BC 589 dwH_AC = dwH_AB + dwH_BC 590 591 # This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 592 self.M0[0] = pA 593 self.M0[1] = pB 594 self.M0[2] = pC 595 596 # Initialise. 597 chi2_sum = 0.0 598 599 # Loop over the experiment types. 600 for ei in range(self.num_exp): 601 # Loop over the spins. 602 for si in range(self.num_spins): 603 # Loop over the spectrometer frequencies. 604 for mi in range(self.num_frq): 605 # The R20 index. 606 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp 607 608 # Convert dw from ppm to rad/s. 609 dw_AB_frq = dw_AB[si] * self.frqs[ei][si][mi] 610 dw_AC_frq = dw_AC[si] * self.frqs[ei][si][mi] 611 dwH_AB_frq = dwH_AB[si] * self.frqs_H[ei][si][mi] 612 dwH_AC_frq = dwH_AC[si] * self.frqs_H[ei][si][mi] 613 614 # Alias the dw frequency combinations. 615 aliased_dwH_AB = 0.0 616 aliased_dwH_AC = 0.0 617 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ: 618 aliased_dw_AB = dw_AB_frq 619 aliased_dw_AC = dw_AC_frq 620 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ: 621 aliased_dw_AB = dwH_AB_frq 622 aliased_dw_AC = dwH_AC_frq 623 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ: 624 aliased_dw_AB = dw_AB_frq + dwH_AB_frq 625 aliased_dw_AC = dw_AC_frq + dwH_AC_frq 626 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ: 627 aliased_dw_AB = dw_AB_frq - dwH_AB_frq 628 aliased_dw_AC = dw_AC_frq - dwH_AC_frq 629 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ: 630 aliased_dw_AB = dw_AB_frq 631 aliased_dw_AC = dw_AC_frq 632 aliased_dwH_AB = dwH_AB_frq 633 aliased_dwH_AC = dwH_AC_frq 634 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ: 635 aliased_dw_AB = dwH_AB_frq 636 aliased_dw_AC = dwH_AC_frq 637 aliased_dwH_AB = dw_AB_frq 638 aliased_dwH_AC = dw_AC_frq 639 640 # Back calculate the R2eff values for each experiment type. 641 self.r2eff_ns_mmq[ei](M0=self.M0, m1=self.m1, m2=self.m2, R20A=R20A[r20_index], R20B=R20B[r20_index], R20C=R20C[r20_index], pA=pA, pB=pB, pC=pC, dw_AB=aliased_dw_AB, dw_AC=aliased_dw_AC, dwH_AB=aliased_dwH_AB, dwH_AC=aliased_dwH_AC, k_AB=k_AB, k_BA=k_BA, k_BC=k_BC, k_CB=k_CB, k_AC=k_AC, k_CA=k_CA, inv_tcpmg=self.inv_relax_times[ei][mi], tcp=self.tau_cpmg[ei][mi], back_calc=self.back_calc[ei][si][mi][0], num_points=self.num_disp_points[ei][si][mi][0], power=self.power[ei][mi]) 642 643 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 644 for di in range(self.num_disp_points[ei][si][mi][0]): 645 if self.missing[ei][si][mi][0][di]: 646 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di] 647 648 # Calculate and return the chi-squared value. 649 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0]) 650 651 # Return the total chi-squared value. 652 return chi2_sum

653 654

655 - def calc_ns_r1rho_3site_chi2(self, r1rho_prime=None, dw_AB=None, dw_BC=None, pA=None, pB=None, kex_AB=None, kex_BC=None, kex_AC=None):

656 """Calculate the chi-squared value for the 'NS MMQ 3-site' models. 657 658 @keyword r1rho_prime: The R1rho value for all states in the absence of exchange. 659 @type r1rho_prime: list of float 660 @keyword dw_AB: The chemical exchange difference between states A and B in rad/s. 661 @type dw_AB: float 662 @keyword dw_BC: The chemical exchange difference between states B and C in rad/s. 663 @type dw_BC: float 664 @keyword pA: The population of state A. 665 @type pA: float 666 @keyword kex_AB: The rate of exchange between states A and B. 667 @type kex_AB: float 668 @keyword kex_BC: The rate of exchange between states B and C. 669 @type kex_BC: float 670 @keyword kex_AC: The rate of exchange between states A and C. 671 @type kex_AC: float 672 @return: The chi-squared value. 673 @rtype: float 674 """ 675 676 # Once off parameter conversions. 677 pC = 1.0 - pA - pB 678 pA_pB = pA + pB 679 pA_pC = pA + pC 680 pB_pC = pB + pC 681 k_BA = pA * kex_AB / pA_pB 682 k_AB = pB * kex_AB / pA_pB 683 k_CB = pB * kex_BC / pB_pC 684 k_BC = pC * kex_BC / pB_pC 685 k_CA = pA * kex_AC / pA_pC 686 k_AC = pC * kex_AC / pA_pC 687 dw_AC = dw_AB + dw_BC 688 689 # Initialise. 690 chi2_sum = 0.0 691 692 # Loop over the spins. 693 for si in range(self.num_spins): 694 # Loop over the spectrometer frequencies. 695 for mi in range(self.num_frq): 696 # The R20 index. 697 r20_index = mi + si*self.num_frq 698 699 # Convert dw from ppm to rad/s. 700 dw_AB_frq = dw_AB[si] * self.frqs[0][si][mi] 701 dw_AC_frq = dw_AC[si] * self.frqs[0][si][mi] 702 703 # Loop over the offsets. 704 for oi in range(self.num_offsets[0][si][mi]): 705 # Back calculate the R2eff values for each experiment type. 706 ns_r1rho_3site(M0=self.M0, matrix=self.matrix, r1rho_prime=r1rho_prime[r20_index], omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], r1=self.r1[si, mi], pA=pA, pB=pB, pC=pC, dw_AB=dw_AB_frq, dw_AC=dw_AC_frq, k_AB=k_AB, k_BA=k_BA, k_BC=k_BC, k_CB=k_CB, k_AC=k_AC, k_CA=k_CA, spin_lock_fields=self.spin_lock_omega1[0][mi][oi], relax_time=self.relax_times[0][mi], inv_relax_time=self.inv_relax_times[0][mi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 707 708 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 709 for di in range(self.num_disp_points[0][si][mi][oi]): 710 if self.missing[0][si][mi][oi][di]: 711 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 712 713 # Calculate and return the chi-squared value. 714 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 715 716 # Return the total chi-squared value. 717 return chi2_sum

718 719

720 - def experiment_type_setup(self):

721 """Check the experiment types and simplify data structures. 722 723 For the single experiment type models, the first dimension of the values, errors, and missing data structures will be removed to simplify the target functions. 724 """ 725 726 # The number of experiments. 727 self.num_exp = len(self.exp_types) 728 729 # The MMQ combined data type models. 730 if self.model in MODEL_LIST_MMQ: 731 # Alias the r2eff functions. 732 self.r2eff_ns_mmq = [] 733 734 # Loop over the experiment types. 735 for ei in range(self.num_exp): 736 # SQ, DQ and ZQ data types. 737 if self.exp_types[ei] in [EXP_TYPE_CPMG_SQ, EXP_TYPE_CPMG_PROTON_SQ, EXP_TYPE_CPMG_DQ, EXP_TYPE_CPMG_ZQ]: 738 if self.model == MODEL_NS_MMQ_2SITE: 739 self.r2eff_ns_mmq.append(r2eff_ns_mmq_2site_sq_dq_zq) 740 else: 741 self.r2eff_ns_mmq.append(r2eff_ns_mmq_3site_sq_dq_zq) 742 743 # MQ data types. 744 elif self.exp_types[ei] in [EXP_TYPE_CPMG_MQ, EXP_TYPE_CPMG_PROTON_MQ]: 745 if self.model == MODEL_NS_MMQ_2SITE: 746 self.r2eff_ns_mmq.append(r2eff_ns_mmq_2site_mq) 747 else: 748 self.r2eff_ns_mmq.append(r2eff_ns_mmq_3site_mq) 749 750 # The single data type models. 751 else: 752 # Check that the data is correct. 753 if self.model != MODEL_NOREX and self.model in MODEL_LIST_CPMG and self.exp_types[0] != EXP_TYPE_CPMG_SQ: 754 raise RelaxError("The '%s' CPMG model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0])) 755 if self.model != MODEL_NOREX and self.model in MODEL_LIST_R1RHO and self.exp_types[0] != EXP_TYPE_R1RHO: 756 raise RelaxError("The '%s' R1rho model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0])) 757 if self.model != MODEL_NOREX and self.model in MODEL_LIST_MQ_CPMG and self.exp_types[0] != EXP_TYPE_CPMG_MQ: 758 raise RelaxError("The '%s' CPMG model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0]))

759 760

761 - def func_CR72(self, params):

762 """Target function for the reduced Carver and Richards (1972) 2-site exchange model on all time scales. 763 764 This assumes that pA > pB, and hence this must be implemented as a constraint. For this model, the simplification R20A = R20B is assumed. 765 766 767 @param params: The vector of parameter values. 768 @type params: numpy rank-1 float array 769 @return: The chi-squared value. 770 @rtype: float 771 """ 772 773 # Scaling. 774 if self.scaling_flag: 775 params = dot(params, self.scaling_matrix) 776 777 # Unpack the parameter values. 778 R20 = params[:self.end_index[0]] 779 dw = params[self.end_index[0]:self.end_index[1]] 780 pA = params[self.end_index[1]] 781 kex = params[self.end_index[1]+1] 782 783 # Calculate and return the chi-squared value. 784 return self.calc_CR72_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)

785 786

787 - def func_CR72_full(self, params):

788 """Target function for the full Carver and Richards (1972) 2-site exchange model on all time scales. 789 790 This assumes that pA > pB, and hence this must be implemented as a constraint. 791 792 793 @param params: The vector of parameter values. 794 @type params: numpy rank-1 float array 795 @return: The chi-squared value. 796 @rtype: float 797 """ 798 799 # Scaling. 800 if self.scaling_flag: 801 params = dot(params, self.scaling_matrix) 802 803 # Unpack the parameter values. 804 R20A = params[:self.end_index[0]] 805 R20B = params[self.end_index[0]:self.end_index[1]] 806 dw = params[self.end_index[1]:self.end_index[2]] 807 pA = params[self.end_index[2]] 808 kex = params[self.end_index[2]+1] 809 810 # Calculate and return the chi-squared value. 811 return self.calc_CR72_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)

812 813

814 - def func_DPL94(self, params):

815 """Target function for the Davis, Perlman and London (1994) fast 2-site off-resonance exchange model for R1rho-type experiments. 816 817 @param params: The vector of parameter values. 818 @type params: numpy rank-1 float array 819 @return: The chi-squared value. 820 @rtype: float 821 """ 822 823 # Scaling. 824 if self.scaling_flag: 825 params = dot(params, self.scaling_matrix) 826 827 # Unpack the parameter values. 828 R20 = params[:self.end_index[0]] 829 phi_ex = params[self.end_index[0]:self.end_index[1]] 830 kex = params[self.end_index[1]] 831 832 # Initialise. 833 chi2_sum = 0.0 834 835 # Loop over the spins. 836 for si in range(self.num_spins): 837 # Loop over the spectrometer frequencies. 838 for mi in range(self.num_frq): 839 # The R20 index. 840 r20_index = mi + si*self.num_frq 841 842 # Convert phi_ex from ppm^2 to (rad/s)^2. 843 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2 844 845 # Loop over the offsets. 846 for oi in range(self.num_offsets[0][si][mi]): 847 # Back calculate the R2eff values. 848 r1rho_DPL94(r1rho_prime=R20[r20_index], phi_ex=phi_ex_scaled, kex=kex, theta=self.tilt_angles[0][si][mi][oi], R1=self.r1[si, mi], spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 849 850 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 851 for di in range(self.num_disp_points[0][si][mi][oi]): 852 if self.missing[0][si][mi][oi][di]: 853 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 854 855 # Calculate and return the chi-squared value. 856 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 857 858 # Return the total chi-squared value. 859 return chi2_sum

860 861

862 - def func_IT99(self, params):

863 """Target function for the Ishima and Torchia (1999) 2-site model for all timescales with pA >> pB. 864 865 @param params: The vector of parameter values. 866 @type params: numpy rank-1 float array 867 @return: The chi-squared value. 868 @rtype: float 869 """ 870 871 # Scaling. 872 if self.scaling_flag: 873 params = dot(params, self.scaling_matrix) 874 875 # Unpack the parameter values. 876 R20 = params[:self.end_index[0]] 877 dw = params[self.end_index[0]:self.end_index[1]] 878 pA = params[self.end_index[1]] 879 tex = params[self.end_index[2]] 880 881 # Once off parameter conversions. 882 pB = 1.0 - pA 883 884 # Initialise. 885 chi2_sum = 0.0 886 887 # Loop over the spins. 888 for si in range(self.num_spins): 889 # Loop over the spectrometer frequencies. 890 for mi in range(self.num_frq): 891 # The R20 index. 892 r20_index = mi + si*self.num_frq 893 894 # Convert dw from ppm to rad/s. 895 dw_frq = dw[si] * self.frqs[0][si][mi] 896 897 # Back calculate the R2eff values. 898 r2eff_IT99(r20=R20[r20_index], pA=pA, pB=pB, dw=dw_frq, tex=tex, cpmg_frqs=self.cpmg_frqs[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 899 900 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 901 for di in range(self.num_disp_points[0][si][mi][0]): 902 if self.missing[0][si][mi][0][di]: 903 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 904 905 # Calculate and return the chi-squared value. 906 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 907 908 # Return the total chi-squared value. 909 return chi2_sum

910 911

912 - def func_LM63_3site(self, params):

913 """Target function for the Luz and Meiboom (1963) fast 3-site exchange model. 914 915 @param params: The vector of parameter values. 916 @type params: numpy rank-1 float array 917 @return: The chi-squared value. 918 @rtype: float 919 """ 920 921 # Scaling. 922 if self.scaling_flag: 923 params = dot(params, self.scaling_matrix) 924 925 # Unpack the parameter values. 926 R20 = params[:self.end_index[0]] 927 phi_ex_B = params[self.end_index[0]:self.end_index[1]] 928 phi_ex_C = params[self.end_index[1]:self.end_index[2]] 929 kB = params[self.end_index[2]] 930 kC = params[self.end_index[2]+1] 931 932 # Once off parameter conversions. 933 rex_B = phi_ex_B / kB 934 rex_C = phi_ex_C / kC 935 quart_kB = kB / 4.0 936 quart_kC = kC / 4.0 937 938 # Initialise. 939 chi2_sum = 0.0 940 941 # Loop over the spins. 942 for si in range(self.num_spins): 943 # Loop over the spectrometer frequencies. 944 for mi in range(self.num_frq): 945 # The R20 index. 946 r20_index = mi + si*self.num_frq 947 948 # Convert phi_ex (or rex) from ppm^2 to (rad/s)^2. 949 frq2 = self.frqs[0][si][mi]**2 950 rex_B_scaled = rex_B[si] * frq2 951 rex_C_scaled = rex_C[si] * frq2 952 953 # Back calculate the R2eff values. 954 r2eff_LM63_3site(r20=R20[r20_index], rex_B=rex_B_scaled, rex_C=rex_C_scaled, quart_kB=quart_kB, quart_kC=quart_kC, cpmg_frqs=self.cpmg_frqs[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 955 956 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 957 for di in range(self.num_disp_points[0][si][mi][0]): 958 if self.missing[0][si][mi][0][di]: 959 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 960 961 # Calculate and return the chi-squared value. 962 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 963 964 # Return the total chi-squared value. 965 return chi2_sum

966 967

968 - def func_LM63(self, params):

969 """Target function for the Luz and Meiboom (1963) fast 2-site exchange model. 970 971 @param params: The vector of parameter values. 972 @type params: numpy rank-1 float array 973 @return: The chi-squared value. 974 @rtype: float 975 """ 976 977 # Scaling. 978 if self.scaling_flag: 979 params = dot(params, self.scaling_matrix) 980 981 # Unpack the parameter values. 982 R20 = params[:self.end_index[0]] 983 phi_ex = params[self.end_index[0]:self.end_index[1]] 984 kex = params[self.end_index[1]] 985 986 # Initialise. 987 chi2_sum = 0.0 988 989 # Loop over the spins. 990 for si in range(self.num_spins): 991 # Loop over the spectrometer frequencies. 992 for mi in range(self.num_frq): 993 # The R20 index. 994 r20_index = mi + si*self.num_frq 995 996 # Convert phi_ex from ppm^2 to (rad/s)^2. 997 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2 998 999 # Back calculate the R2eff values. 1000 r2eff_LM63(r20=R20[r20_index], phi_ex=phi_ex_scaled, kex=kex, cpmg_frqs=self.cpmg_frqs[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 1001 1002 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1003 for di in range(self.num_disp_points[0][si][mi][0]): 1004 if self.missing[0][si][mi][0][di]: 1005 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 1006 1007 # Calculate and return the chi-squared value. 1008 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 1009 1010 # Return the total chi-squared value. 1011 return chi2_sum

1012 1013

1014 - def func_M61(self, params):

1015 """Target function for the Meiboom (1961) fast 2-site exchange model for R1rho-type experiments. 1016 1017 @param params: The vector of parameter values. 1018 @type params: numpy rank-1 float array 1019 @return: The chi-squared value. 1020 @rtype: float 1021 """ 1022 1023 # Scaling. 1024 if self.scaling_flag: 1025 params = dot(params, self.scaling_matrix) 1026 1027 # Unpack the parameter values. 1028 R20 = params[:self.end_index[0]] 1029 phi_ex = params[self.end_index[0]:self.end_index[1]] 1030 kex = params[self.end_index[1]] 1031 1032 # Initialise. 1033 chi2_sum = 0.0 1034 1035 # Loop over the spins. 1036 for si in range(self.num_spins): 1037 # Loop over the spectrometer frequencies. 1038 for mi in range(self.num_frq): 1039 # The R20 index. 1040 r20_index = mi + si*self.num_frq 1041 1042 # Convert phi_ex from ppm^2 to (rad/s)^2. 1043 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2 1044 1045 # Back calculate the R2eff values. 1046 r1rho_M61(r1rho_prime=R20[r20_index], phi_ex=phi_ex_scaled, kex=kex, spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 1047 1048 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1049 for di in range(self.num_disp_points[0][si][mi][0]): 1050 if self.missing[0][si][mi][0][di]: 1051 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 1052 1053 # Calculate and return the chi-squared value. 1054 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 1055 1056 # Return the total chi-squared value. 1057 return chi2_sum

1058 1059

1060 - def func_M61b(self, params):

1061 """Target function for the Meiboom (1961) R1rho on-resonance 2-site model for skewed populations (pA >> pB). 1062 1063 @param params: The vector of parameter values. 1064 @type params: numpy rank-1 float array 1065 @return: The chi-squared value. 1066 @rtype: float 1067 """ 1068 1069 # Scaling. 1070 if self.scaling_flag: 1071 params = dot(params, self.scaling_matrix) 1072 1073 # Unpack the parameter values. 1074 R20 = params[:self.end_index[0]] 1075 dw = params[self.end_index[0]:self.end_index[1]] 1076 pA = params[self.end_index[1]] 1077 kex = params[self.end_index[1]+1] 1078 1079 # Initialise. 1080 chi2_sum = 0.0 1081 1082 # Loop over the spins. 1083 for si in range(self.num_spins): 1084 # Loop over the spectrometer frequencies. 1085 for mi in range(self.num_frq): 1086 # The R20 index. 1087 r20_index = mi + si*self.num_frq 1088 1089 # Convert dw from ppm to rad/s. 1090 dw_frq = dw[si] * self.frqs[0][si][mi] 1091 1092 # Back calculate the R1rho values. 1093 r1rho_M61b(r1rho_prime=R20[r20_index], pA=pA, dw=dw_frq, kex=kex, spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][0], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 1094 1095 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1096 for di in range(self.num_disp_points[0][si][mi][0]): 1097 if self.missing[0][si][mi][0][di]: 1098 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 1099 1100 # Calculate and return the chi-squared value. 1101 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 1102 1103 # Return the total chi-squared value. 1104 return chi2_sum

1105 1106

1107 - def func_MP05(self, params):

1108 """Target function for the Miloushev and Palmer (2005) R1rho off-resonance 2-site model. 1109 1110 @param params: The vector of parameter values. 1111 @type params: numpy rank-1 float array 1112 @return: The chi-squared value. 1113 @rtype: float 1114 """ 1115 1116 # Scaling. 1117 if self.scaling_flag: 1118 params = dot(params, self.scaling_matrix) 1119 1120 # Unpack the parameter values. 1121 R20 = params[:self.end_index[0]] 1122 dw = params[self.end_index[0]:self.end_index[1]] 1123 pA = params[self.end_index[1]] 1124 kex = params[self.end_index[1]+1] 1125 1126 # Once off parameter conversions. 1127 pB = 1.0 - pA 1128 1129 # Initialise. 1130 chi2_sum = 0.0 1131 1132 # Loop over the spins. 1133 for si in range(self.num_spins): 1134 # Loop over the spectrometer frequencies. 1135 for mi in range(self.num_frq): 1136 # The R20 index. 1137 r20_index = mi + si*self.num_frq 1138 1139 # Convert dw from ppm to rad/s. 1140 dw_frq = dw[si] * self.frqs[0][si][mi] 1141 1142 # Loop over the offsets. 1143 for oi in range(self.num_offsets[0][si][mi]): 1144 # Back calculate the R1rho values. 1145 r1rho_MP05(r1rho_prime=R20[r20_index], omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], spin_lock_fields=self.spin_lock_omega1[0][mi][oi], spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 1146 1147 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1148 for di in range(self.num_disp_points[0][si][mi][oi]): 1149 if self.missing[0][si][mi][oi][di]: 1150 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 1151 1152 # Calculate and return the chi-squared value. 1153 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 1154 1155 # Return the total chi-squared value. 1156 return chi2_sum

1157 1158

1159 - def func_mmq_CR72(self, params):

1160 """Target function for the CR72 model extended for MQ CPMG data. 1161 1162 @param params: The vector of parameter values. 1163 @type params: numpy rank-1 float array 1164 @return: The chi-squared value. 1165 @rtype: float 1166 """ 1167 1168 # Scaling. 1169 if self.scaling_flag: 1170 params = dot(params, self.scaling_matrix) 1171 1172 # Unpack the parameter values. 1173 R20 = params[:self.end_index[0]] 1174 dw = params[self.end_index[0]:self.end_index[1]] 1175 dwH = params[self.end_index[1]:self.end_index[2]] 1176 pA = params[self.end_index[2]] 1177 kex = params[self.end_index[2]+1] 1178 1179 # Once off parameter conversions. 1180 pB = 1.0 - pA 1181 k_BA = pA * kex 1182 k_AB = pB * kex 1183 1184 # Initialise. 1185 chi2_sum = 0.0 1186 1187 # Loop over the experiment types. 1188 for ei in range(self.num_exp): 1189 # Loop over the spins. 1190 for si in range(self.num_spins): 1191 # Loop over the spectrometer frequencies. 1192 for mi in range(self.num_frq): 1193 # The R20 index. 1194 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp 1195 1196 # Convert dw from ppm to rad/s. 1197 dw_frq = dw[si] * self.frqs[ei][si][mi] 1198 dwH_frq = dwH[si] * self.frqs_H[ei][si][mi] 1199 1200 # Alias the dw frequency combinations. 1201 aliased_dwH = 0.0 1202 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ: 1203 aliased_dw = dw_frq 1204 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ: 1205 aliased_dw = dwH_frq 1206 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ: 1207 aliased_dw = dw_frq + dwH_frq 1208 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ: 1209 aliased_dw = dw_frq - dwH_frq 1210 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ: 1211 aliased_dw = dw_frq 1212 aliased_dwH = dwH_frq 1213 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ: 1214 aliased_dw = dwH_frq 1215 aliased_dwH = dw_frq 1216 1217 # Back calculate the R2eff values. 1218 r2eff_mmq_cr72(r20=R20[r20_index], pA=pA, pB=pB, dw=aliased_dw, dwH=aliased_dwH, kex=kex, k_AB=k_AB, k_BA=k_BA, cpmg_frqs=self.cpmg_frqs[ei][mi][0], inv_tcpmg=self.inv_relax_times[ei][mi], tcp=self.tau_cpmg[ei][mi], back_calc=self.back_calc[ei][si][mi][0], num_points=self.num_disp_points[ei][si][mi][0], power=self.power[ei][mi]) 1219 1220 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1221 for di in range(self.num_disp_points[ei][si][mi][0]): 1222 if self.missing[ei][si][mi][0][di]: 1223 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di] 1224 1225 # Calculate and return the chi-squared value. 1226 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0]) 1227 1228 # Return the total chi-squared value. 1229 return chi2_sum

1230 1231

1232 - def func_NOREX(self, params):

1233 """Target function for no exchange. 1234 1235 @param params: The vector of parameter values. 1236 @type params: numpy rank-1 float array 1237 @return: The chi-squared value. 1238 @rtype: float 1239 """ 1240 1241 # Scaling. 1242 if self.scaling_flag: 1243 params = dot(params, self.scaling_matrix) 1244 1245 # Unpack the parameter values. 1246 R20 = params 1247 1248 # Initialise. 1249 chi2_sum = 0.0 1250 1251 # Loop over the experiment types. 1252 for ei in range(self.num_exp): 1253 # Loop over the spins. 1254 for si in range(self.num_spins): 1255 # Loop over the spectrometer frequencies. 1256 for mi in range(self.num_frq): 1257 # The R20 index. 1258 r20_index = mi + si*self.num_frq 1259 1260 # Loop over the offsets. 1261 for oi in range(self.num_offsets[ei][si][mi]): 1262 # The R2eff values as R20 values. 1263 for di in range(self.num_disp_points[ei][si][mi][oi]): 1264 self.back_calc[ei][si][mi][oi][di] = R20[r20_index] 1265 1266 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1267 for di in range(self.num_disp_points[ei][si][mi][oi]): 1268 if self.missing[ei][si][mi][oi][di]: 1269 self.back_calc[ei][si][mi][oi][di] = self.values[ei][si][mi][oi][di] 1270 1271 # Calculate and return the chi-squared value. 1272 chi2_sum += chi2(self.values[ei][si][mi][oi], self.back_calc[ei][si][mi][oi], self.errors[ei][si][mi][oi]) 1273 1274 # Return the total chi-squared value. 1275 return chi2_sum

1276 1277

1278 - def func_ns_cpmg_2site_3D(self, params):

1279 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations. 1280 1281 @param params: The vector of parameter values. 1282 @type params: numpy rank-1 float array 1283 @return: The chi-squared value. 1284 @rtype: float 1285 """ 1286 1287 # Scaling. 1288 if self.scaling_flag: 1289 params = dot(params, self.scaling_matrix) 1290 1291 # Unpack the parameter values. 1292 R20 = params[:self.end_index[0]] 1293 dw = params[self.end_index[0]:self.end_index[1]] 1294 pA = params[self.end_index[1]] 1295 kex = params[self.end_index[1]+1] 1296 1297 # Calculate and return the chi-squared value. 1298 return self.calc_ns_cpmg_2site_3D_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)

1299 1300

1301 - def func_ns_cpmg_2site_3D_full(self, params):

1302 """Target function for the full numerical solution for the 2-site Bloch-McConnell equations. 1303 1304 @param params: The vector of parameter values. 1305 @type params: numpy rank-1 float array 1306 @return: The chi-squared value. 1307 @rtype: float 1308 """ 1309 1310 # Scaling. 1311 if self.scaling_flag: 1312 params = dot(params, self.scaling_matrix) 1313 1314 # Unpack the parameter values. 1315 R20A = params[:self.end_index[0]] 1316 R20B = params[self.end_index[0]:self.end_index[1]] 1317 dw = params[self.end_index[1]:self.end_index[2]] 1318 pA = params[self.end_index[2]] 1319 kex = params[self.end_index[2]+1] 1320 1321 # Calculate and return the chi-squared value. 1322 return self.calc_ns_cpmg_2site_3D_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)

1323 1324

1325 - def func_ns_cpmg_2site_expanded(self, params):

1326 """Target function for the numerical solution for the 2-site Bloch-McConnell equations using the expanded notation. 1327 1328 @param params: The vector of parameter values. 1329 @type params: numpy rank-1 float array 1330 @return: The chi-squared value. 1331 @rtype: float 1332 """ 1333 1334 # Scaling. 1335 if self.scaling_flag: 1336 params = dot(params, self.scaling_matrix) 1337 1338 # Unpack the parameter values. 1339 R20 = params[:self.end_index[0]] 1340 dw = params[self.end_index[0]:self.end_index[1]] 1341 pA = params[self.end_index[1]] 1342 kex = params[self.end_index[1]+1] 1343 1344 # Once off parameter conversions. 1345 pB = 1.0 - pA 1346 k_BA = pA * kex 1347 k_AB = pB * kex 1348 1349 # Chi-squared initialisation. 1350 chi2_sum = 0.0 1351 1352 # Loop over the spins. 1353 for si in range(self.num_spins): 1354 # Loop over the spectrometer frequencies. 1355 for mi in range(self.num_frq): 1356 # The R20 index. 1357 r20_index = mi + si*self.num_frq 1358 1359 # Convert dw from ppm to rad/s. 1360 dw_frq = dw[si] * self.frqs[0][si][mi] 1361 1362 # Back calculate the R2eff values. 1363 r2eff_ns_cpmg_2site_expanded(r20=R20[r20_index], pA=pA, dw=dw_frq, k_AB=k_AB, k_BA=k_BA, relax_time=self.relax_times[0][mi], inv_relax_time=self.inv_relax_times[0][mi], tcp=self.tau_cpmg[0][mi], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0], num_cpmg=self.power[0][mi]) 1364 1365 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1366 for di in range(self.num_disp_points[0][si][mi][0]): 1367 if self.missing[0][si][mi][0][di]: 1368 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 1369 1370 # Calculate and return the chi-squared value. 1371 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 1372 1373 # Return the total chi-squared value. 1374 return chi2_sum

1375 1376

1377 - def func_ns_cpmg_2site_star(self, params):

1378 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations using complex conjugate matrices. 1379 1380 This is the model whereby the simplification R20A = R20B is assumed. 1381 1382 1383 @param params: The vector of parameter values. 1384 @type params: numpy rank-1 float array 1385 @return: The chi-squared value. 1386 @rtype: float 1387 """ 1388 1389 # Scaling. 1390 if self.scaling_flag: 1391 params = dot(params, self.scaling_matrix) 1392 1393 # Unpack the parameter values. 1394 R20 = params[:self.end_index[0]] 1395 dw = params[self.end_index[0]:self.end_index[1]] 1396 pA = params[self.end_index[1]] 1397 kex = params[self.end_index[1]+1] 1398 1399 # Calculate and return the chi-squared value. 1400 return self.calc_ns_cpmg_2site_star_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)

1401 1402

1403 - def func_ns_cpmg_2site_star_full(self, params):

1404 """Target function for the full numerical solution for the 2-site Bloch-McConnell equations using complex conjugate matrices. 1405 1406 @param params: The vector of parameter values. 1407 @type params: numpy rank-1 float array 1408 @return: The chi-squared value. 1409 @rtype: float 1410 """ 1411 1412 # Scaling. 1413 if self.scaling_flag: 1414 params = dot(params, self.scaling_matrix) 1415 1416 # Unpack the parameter values. 1417 R20A = params[:self.end_index[0]] 1418 R20B = params[self.end_index[0]:self.end_index[1]] 1419 dw = params[self.end_index[1]:self.end_index[2]] 1420 pA = params[self.end_index[2]] 1421 kex = params[self.end_index[2]+1] 1422 1423 # Calculate and return the chi-squared value. 1424 return self.calc_ns_cpmg_2site_star_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)

1425 1426

1427 - def func_ns_mmq_2site(self, params):

1428 """Target function for the combined SQ, ZQ, DQ and MQ CPMG numeric solution. 1429 1430 @param params: The vector of parameter values. 1431 @type params: numpy rank-1 float array 1432 @return: The chi-squared value. 1433 @rtype: float 1434 """ 1435 1436 # Scaling. 1437 if self.scaling_flag: 1438 params = dot(params, self.scaling_matrix) 1439 1440 # Unpack the parameter values. 1441 R20 = params[:self.end_index[0]] 1442 dw = params[self.end_index[0]:self.end_index[1]] 1443 dwH = params[self.end_index[1]:self.end_index[2]] 1444 pA = params[self.end_index[2]] 1445 kex = params[self.end_index[2]+1] 1446 1447 # Once off parameter conversions. 1448 pB = 1.0 - pA 1449 k_BA = pA * kex 1450 k_AB = pB * kex 1451 1452 # This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 1453 self.M0[0] = pA 1454 self.M0[1] = pB 1455 1456 # Initialise. 1457 chi2_sum = 0.0 1458 1459 # Loop over the experiment types. 1460 for ei in range(self.num_exp): 1461 # Loop over the spins. 1462 for si in range(self.num_spins): 1463 # Loop over the spectrometer frequencies. 1464 for mi in range(self.num_frq): 1465 # The R20 index. 1466 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp 1467 1468 # Convert dw from ppm to rad/s. 1469 dw_frq = dw[si] * self.frqs[ei][si][mi] 1470 dwH_frq = dwH[si] * self.frqs_H[ei][si][mi] 1471 1472 # Alias the dw frequency combinations. 1473 aliased_dwH = 0.0 1474 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ: 1475 aliased_dw = dw_frq 1476 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ: 1477 aliased_dw = dwH_frq 1478 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ: 1479 aliased_dw = dw_frq + dwH_frq 1480 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ: 1481 aliased_dw = dw_frq - dwH_frq 1482 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ: 1483 aliased_dw = dw_frq 1484 aliased_dwH = dwH_frq 1485 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ: 1486 aliased_dw = dwH_frq 1487 aliased_dwH = dw_frq 1488 1489 # Back calculate the R2eff values for each experiment type. 1490 self.r2eff_ns_mmq[ei](M0=self.M0, m1=self.m1, m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB, dw=aliased_dw, dwH=aliased_dwH, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_times[ei][mi], tcp=self.tau_cpmg[ei][mi], back_calc=self.back_calc[ei][si][mi][0], num_points=self.num_disp_points[ei][si][mi][0], power=self.power[ei][mi]) 1491 1492 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1493 for di in range(self.num_disp_points[ei][si][mi][0]): 1494 if self.missing[ei][si][mi][0][di]: 1495 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di] 1496 1497 # Calculate and return the chi-squared value. 1498 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0]) 1499 1500 # Return the total chi-squared value. 1501 return chi2_sum

1502 1503

1504 - def func_ns_mmq_3site(self, params):

1505 """Target function for the combined SQ, ZQ, DQ and MQ 3-site MMQ CPMG numeric solution. 1506 1507 @param params: The vector of parameter values. 1508 @type params: numpy rank-1 float array 1509 @return: The chi-squared value. 1510 @rtype: float 1511 """ 1512 1513 # Scaling. 1514 if self.scaling_flag: 1515 params = dot(params, self.scaling_matrix) 1516 1517 # Unpack the parameter values. 1518 R20 = params[:self.end_index[0]] 1519 dw_AB = params[self.end_index[0]:self.end_index[1]] 1520 dw_BC = params[self.end_index[1]:self.end_index[2]] 1521 dwH_AB = params[self.end_index[2]:self.end_index[3]] 1522 dwH_BC = params[self.end_index[3]:self.end_index[4]] 1523 pA = params[self.end_index[4]] 1524 kex_AB = params[self.end_index[4]+1] 1525 pB = params[self.end_index[4]+2] 1526 kex_BC = params[self.end_index[4]+3] 1527 kex_AC = params[self.end_index[4]+4] 1528 1529 # Calculate and return the chi-squared value. 1530 return self.calc_ns_mmq_3site_chi2(R20A=R20, R20B=R20, R20C=R20, dw_AB=dw_AB, dw_BC=dw_BC, dwH_AB=dwH_AB, dwH_BC=dwH_BC, pA=pA, pB=pB, kex_AB=kex_AB, kex_BC=kex_BC, kex_AC=kex_AC)

1531 1532

1533 - def func_ns_mmq_3site_linear(self, params):

1534 """Target function for the combined SQ, ZQ, DQ and MQ 3-site linearised MMQ CPMG numeric solution. 1535 1536 @param params: The vector of parameter values. 1537 @type params: numpy rank-1 float array 1538 @return: The chi-squared value. 1539 @rtype: float 1540 """ 1541 1542 # Scaling. 1543 if self.scaling_flag: 1544 params = dot(params, self.scaling_matrix) 1545 1546 # Unpack the parameter values. 1547 R20 = params[:self.end_index[0]] 1548 dw_AB = params[self.end_index[0]:self.end_index[1]] 1549 dw_BC = params[self.end_index[1]:self.end_index[2]] 1550 dwH_AB = params[self.end_index[2]:self.end_index[3]] 1551 dwH_BC = params[self.end_index[3]:self.end_index[4]] 1552 pA = params[self.end_index[4]] 1553 kex_AB = params[self.end_index[4]+1] 1554 pB = params[self.end_index[4]+2] 1555 kex_BC = params[self.end_index[4]+3] 1556 1557 # Calculate and return the chi-squared value. 1558 return self.calc_ns_mmq_3site_chi2(R20A=R20, R20B=R20, R20C=R20, dw_AB=dw_AB, dw_BC=dw_BC, dwH_AB=dwH_AB, dwH_BC=dwH_BC, pA=pA, pB=pB, kex_AB=kex_AB, kex_BC=kex_BC, kex_AC=0.0)

1559 1560

1561 - def func_ns_r1rho_2site(self, params):

1562 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations for R1rho data. 1563 1564 @param params: The vector of parameter values. 1565 @type params: numpy rank-1 float array 1566 @return: The chi-squared value. 1567 @rtype: float 1568 """ 1569 1570 # Scaling. 1571 if self.scaling_flag: 1572 params = dot(params, self.scaling_matrix) 1573 1574 # Unpack the parameter values. 1575 r1rho_prime = params[:self.end_index[0]] 1576 dw = params[self.end_index[0]:self.end_index[1]] 1577 pA = params[self.end_index[1]] 1578 kex = params[self.end_index[1]+1] 1579 1580 # Once off parameter conversions. 1581 pB = 1.0 - pA 1582 k_BA = pA * kex 1583 k_AB = pB * kex 1584 1585 # Chi-squared initialisation. 1586 chi2_sum = 0.0 1587 1588 # Loop over the spins. 1589 for si in range(self.num_spins): 1590 # Loop over the spectrometer frequencies. 1591 for mi in range(self.num_frq): 1592 # The R20 index. 1593 r20_index = mi + si*self.num_frq 1594 1595 # Convert dw from ppm to rad/s. 1596 dw_frq = dw[si] * self.frqs[0][si][mi] 1597 1598 # Loop over the offsets. 1599 for oi in range(self.num_offsets[0][si][mi]): 1600 # Back calculate the R2eff values. 1601 ns_r1rho_2site(M0=self.M0, matrix=self.matrix, r1rho_prime=r1rho_prime[r20_index], omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], r1=self.r1[si, mi], pA=pA, pB=pB, dw=dw_frq, k_AB=k_AB, k_BA=k_BA, spin_lock_fields=self.spin_lock_omega1[0][mi][oi], relax_time=self.relax_times[0][mi], inv_relax_time=self.inv_relax_times[0][mi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 1602 1603 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1604 for di in range(self.num_disp_points[0][si][mi][oi]): 1605 if self.missing[0][si][mi][oi][di]: 1606 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 1607 1608 # Calculate and return the chi-squared value. 1609 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 1610 1611 # Return the total chi-squared value. 1612 return chi2_sum

1613 1614

1615 - def func_ns_r1rho_3site(self, params):

1616 """Target function for the R1rho 3-site numeric solution. 1617 1618 @param params: The vector of parameter values. 1619 @type params: numpy rank-1 float array 1620 @return: The chi-squared value. 1621 @rtype: float 1622 """ 1623 1624 # Scaling. 1625 if self.scaling_flag: 1626 params = dot(params, self.scaling_matrix) 1627 1628 # Unpack the parameter values. 1629 r1rho_prime = params[:self.end_index[0]] 1630 dw_AB = params[self.end_index[0]:self.end_index[1]] 1631 dw_BC = params[self.end_index[1]:self.end_index[2]] 1632 pA = params[self.end_index[2]] 1633 kex_AB = params[self.end_index[2]+1] 1634 pB = params[self.end_index[2]+2] 1635 kex_BC = params[self.end_index[2]+3] 1636 kex_AC = params[self.end_index[2]+4] 1637 1638 # Calculate and return the chi-squared value. 1639 return self.calc_ns_r1rho_3site_chi2(r1rho_prime=r1rho_prime, dw_AB=dw_AB, dw_BC=dw_BC, pA=pA, pB=pB, kex_AB=kex_AB, kex_BC=kex_BC, kex_AC=kex_AC)

1640 1641

1642 - def func_ns_r1rho_3site_linear(self, params):

1643 """Target function for the R1rho 3-site numeric solution linearised with kAC = kCA = 0. 1644 1645 @param params: The vector of parameter values. 1646 @type params: numpy rank-1 float array 1647 @return: The chi-squared value. 1648 @rtype: float 1649 """ 1650 1651 # Scaling. 1652 if self.scaling_flag: 1653 params = dot(params, self.scaling_matrix) 1654 1655 # Unpack the parameter values. 1656 r1rho_prime = params[:self.end_index[0]] 1657 dw_AB = params[self.end_index[0]:self.end_index[1]] 1658 dw_BC = params[self.end_index[1]:self.end_index[2]] 1659 pA = params[self.end_index[2]] 1660 kex_AB = params[self.end_index[2]+1] 1661 pB = params[self.end_index[2]+2] 1662 kex_BC = params[self.end_index[2]+3] 1663 1664 # Calculate and return the chi-squared value. 1665 return self.calc_ns_r1rho_3site_chi2(r1rho_prime=r1rho_prime, dw_AB=dw_AB, dw_BC=dw_BC, pA=pA, pB=pB, kex_AB=kex_AB, kex_BC=kex_BC, kex_AC=0.0)

1666 1667

1668 - def func_TAP03(self, params):

1669 """Target function for the Trott, Abergel and Palmer (2003) R1rho off-resonance 2-site model. 1670 1671 @param params: The vector of parameter values. 1672 @type params: numpy rank-1 float array 1673 @return: The chi-squared value. 1674 @rtype: float 1675 """ 1676 1677 # Scaling. 1678 if self.scaling_flag: 1679 params = dot(params, self.scaling_matrix) 1680 1681 # Unpack the parameter values. 1682 R20 = params[:self.end_index[0]] 1683 dw = params[self.end_index[0]:self.end_index[1]] 1684 pA = params[self.end_index[1]] 1685 kex = params[self.end_index[1]+1] 1686 1687 # Once off parameter conversions. 1688 pB = 1.0 - pA 1689 1690 # Initialise. 1691 chi2_sum = 0.0 1692 1693 # Loop over the spins. 1694 for si in range(self.num_spins): 1695 # Loop over the spectrometer frequencies. 1696 for mi in range(self.num_frq): 1697 # The R20 index. 1698 r20_index = mi + si*self.num_frq 1699 1700 # Convert dw from ppm to rad/s. 1701 dw_frq = dw[si] * self.frqs[0][si][mi] 1702 1703 # Loop over the offsets. 1704 for oi in range(self.num_offsets[0][si][mi]): 1705 # Back calculate the R1rho values. 1706 r1rho_TAP03(r1rho_prime=R20[r20_index], omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], spin_lock_fields=self.spin_lock_omega1[0][mi][oi], spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 1707 1708 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1709 for di in range(self.num_disp_points[0][si][mi][oi]): 1710 if self.missing[0][si][mi][oi][di]: 1711 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 1712 1713 # Calculate and return the chi-squared value. 1714 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 1715 1716 # Return the total chi-squared value. 1717 return chi2_sum

1718 1719

1720 - def func_TP02(self, params):

1721 """Target function for the Trott and Palmer (2002) R1rho off-resonance 2-site model. 1722 1723 @param params: The vector of parameter values. 1724 @type params: numpy rank-1 float array 1725 @return: The chi-squared value. 1726 @rtype: float 1727 """ 1728 1729 # Scaling. 1730 if self.scaling_flag: 1731 params = dot(params, self.scaling_matrix) 1732 1733 # Unpack the parameter values. 1734 R20 = params[:self.end_index[0]] 1735 dw = params[self.end_index[0]:self.end_index[1]] 1736 pA = params[self.end_index[1]] 1737 kex = params[self.end_index[1]+1] 1738 1739 # Once off parameter conversions. 1740 pB = 1.0 - pA 1741 1742 # Initialise. 1743 chi2_sum = 0.0 1744 1745 # Loop over the spins. 1746 for si in range(self.num_spins): 1747 # Loop over the spectrometer frequencies. 1748 for mi in range(self.num_frq): 1749 # The R20 index. 1750 r20_index = mi + si*self.num_frq 1751 1752 # Convert dw from ppm to rad/s. 1753 dw_frq = dw[si] * self.frqs[0][si][mi] 1754 1755 # Loop over the offsets. 1756 for oi in range(self.num_offsets[0][si][mi]): 1757 # Back calculate the R1rho values. 1758 r1rho_TP02(r1rho_prime=R20[r20_index], omega=self.chemical_shifts[0][si][mi], offset=self.offset[0][si][mi][oi], pA=pA, pB=pB, dw=dw_frq, kex=kex, R1=self.r1[si, mi], spin_lock_fields=self.spin_lock_omega1[0][mi][oi], spin_lock_fields2=self.spin_lock_omega1_squared[0][mi][oi], back_calc=self.back_calc[0][si][mi][oi], num_points=self.num_disp_points[0][si][mi][oi]) 1759 1760 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1761 for di in range(self.num_disp_points[0][si][mi][oi]): 1762 if self.missing[0][si][mi][oi][di]: 1763 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di] 1764 1765 # Calculate and return the chi-squared value. 1766 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi]) 1767 1768 # Return the total chi-squared value. 1769 return chi2_sum

1770 1771

1772 - def func_TSMFK01(self, params):

1773 """Target function for the the Tollinger et al. (2001) 2-site very-slow exchange model, range of microsecond to second time scale. 1774 1775 @param params: The vector of parameter values. 1776 @type params: numpy rank-1 float array 1777 @return: The chi-squared value. 1778 @rtype: float 1779 """ 1780 1781 # Scaling. 1782 if self.scaling_flag: 1783 params = dot(params, self.scaling_matrix) 1784 1785 # Unpack the parameter values. 1786 R20A = params[:self.end_index[0]] 1787 dw = params[self.end_index[0]:self.end_index[1]] 1788 k_AB = params[self.end_index[1]] 1789 1790 # Initialise. 1791 chi2_sum = 0.0 1792 1793 # Loop over the spins. 1794 for si in range(self.num_spins): 1795 # Loop over the spectrometer frequencies. 1796 for mi in range(self.num_frq): 1797 # The R20 index. 1798 r20a_index = mi + si*self.num_frq 1799 1800 # Convert dw from ppm to rad/s. 1801 dw_frq = dw[si] * self.frqs[0][si][mi] 1802 1803 # Back calculate the R2eff values. 1804 r2eff_TSMFK01(r20a=R20A[r20a_index], dw=dw_frq, k_AB=k_AB, tcp=self.tau_cpmg[0][mi], back_calc=self.back_calc[0][si][mi][0], num_points=self.num_disp_points[0][si][mi][0]) 1805 1806 # For all missing data points, set the back-calculated value to the measured values so that it has no effect on the chi-squared value. 1807 for di in range(self.num_disp_points[0][si][mi][0]): 1808 if self.missing[0][si][mi][0][di]: 1809 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di] 1810 1811 # Calculate and return the chi-squared value. 1812 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0]) 1813 1814 # Return the total chi-squared value. 1815 return chi2_sum

Source Code for Module target_functions.relax_disp