1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25 """Target functions for relaxation dispersion."""
26
27
28 from copy import deepcopy
29 from math import pi
30 from numpy import complex64, dot, float64, int16, sqrt, zeros
31
32
33 from lib.dispersion.b14 import r2eff_B14
34 from lib.dispersion.cr72 import r2eff_CR72
35 from lib.dispersion.dpl94 import r1rho_DPL94
36 from lib.dispersion.it99 import r2eff_IT99
37 from lib.dispersion.lm63 import r2eff_LM63
38 from lib.dispersion.lm63_3site import r2eff_LM63_3site
39 from lib.dispersion.m61 import r1rho_M61
40 from lib.dispersion.m61b import r1rho_M61b
41 from lib.dispersion.mp05 import r1rho_MP05
42 from lib.dispersion.mmq_cr72 import r2eff_mmq_cr72
43 from lib.dispersion.ns_cpmg_2site_3d import r2eff_ns_cpmg_2site_3D
44 from lib.dispersion.ns_cpmg_2site_expanded import r2eff_ns_cpmg_2site_expanded
45 from lib.dispersion.ns_cpmg_2site_star import r2eff_ns_cpmg_2site_star
46 from lib.dispersion.ns_mmq_3site import r2eff_ns_mmq_3site_mq, r2eff_ns_mmq_3site_sq_dq_zq
47 from lib.dispersion.ns_mmq_2site import r2eff_ns_mmq_2site_mq, r2eff_ns_mmq_2site_sq_dq_zq
48 from lib.dispersion.ns_r1rho_2site import ns_r1rho_2site
49 from lib.dispersion.ns_r1rho_3site import ns_r1rho_3site
50 from lib.dispersion.ns_matrices import r180x_3d
51 from lib.dispersion.tp02 import r1rho_TP02
52 from lib.dispersion.tap03 import r1rho_TAP03
53 from lib.dispersion.tsmfk01 import r2eff_TSMFK01
54 from lib.errors import RelaxError
55 from lib.float import isNaN
56 from target_functions.chi2 import chi2
57 from specific_analyses.relax_disp.variables import EXP_TYPE_CPMG_DQ, EXP_TYPE_CPMG_MQ, EXP_TYPE_CPMG_PROTON_MQ, EXP_TYPE_CPMG_PROTON_SQ, EXP_TYPE_CPMG_SQ, EXP_TYPE_CPMG_ZQ, EXP_TYPE_R1RHO, MODEL_B14, MODEL_B14_FULL, MODEL_CR72, MODEL_CR72_FULL, MODEL_DPL94, MODEL_IT99, MODEL_LIST_CPMG, MODEL_LIST_FULL, MODEL_LIST_MMQ, MODEL_LIST_MQ_CPMG, MODEL_LIST_R1RHO, MODEL_LM63, MODEL_LM63_3SITE, MODEL_M61, MODEL_M61B, MODEL_MP05, MODEL_MMQ_CR72, MODEL_NOREX, 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, MODEL_TAP03, MODEL_TP02, MODEL_TSMFK01
58
59
61 - 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):
62 """Relaxation dispersion target functions for optimisation.
63
64 Models
65 ======
66
67 The following analytic models are currently supported:
68
69 - 'No Rex': The model for no chemical exchange relaxation.
70 - 'LM63': The Luz and Meiboom (1963) 2-site fast exchange model.
71 - 'LM63 3-site': The Luz and Meiboom (1963) 3-site fast exchange model.
72 - 'CR72': The reduced Carver and Richards (1972) 2-site model for all time scales with R20A = R20B.
73 - 'CR72 full': The full Carver and Richards (1972) 2-site model for all time scales.
74 - 'IT99': The Ishima and Torchia (1999) 2-site model for all time scales with skewed populations (pA >> pB).
75 - 'TSMFK01': The Tollinger et al. (2001) 2-site very-slow exchange model, range of microsecond to second time scale.
76 - 'B14': The Baldwin (2014) 2-site exact solution model for all time scales with R20A = R20B..
77 - 'B14 full': The Baldwin (2014) 2-site exact solution model for all time scales.
78 - 'M61': The Meiboom (1961) 2-site fast exchange model for R1rho-type experiments.
79 - 'DPL94': The Davis, Perlman and London (1994) 2-site fast exchange model for R1rho-type experiments.
80 - 'M61 skew': The Meiboom (1961) on-resonance 2-site model with skewed populations (pA >> pB) for R1rho-type experiments.
81 - 'TP02': The Trott and Palmer (2002) 2-site exchange model for R1rho-type experiments.
82 - 'TAP03': The Trott, Abergel and Palmer (2003) 2-site exchange model for R1rho-type experiments.
83 - 'MP05': The Miloushev and Palmer (2005) off-resonance 2-site exchange model for R1rho-type experiments.
84
85 The following numerical models are currently supported:
86
87 - '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.
88 - 'NS CPMG 2-site 3D full': The full numerical solution for the 2-site Bloch-McConnell equations for CPMG data using 3D magnetisation vectors.
89 - '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.
90 - 'NS CPMG 2-site star full': The full numerical solution for the 2-site Bloch-McConnell equations for CPMG data using complex conjugate matrices.
91 - 'NS CPMG 2-site expanded': The numerical solution for the 2-site Bloch-McConnell equations for CPMG data expanded using Maple by Nikolai Skrynnikov.
92 - '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.
93 - '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.
94 - '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.
95 - 'NS R1rho 2-site': The numerical solution for the 2-site Bloch-McConnell equations for R1rho data with R20A = R20B.
96 - '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.
97 - 'NS R1rho 3-site': The numerical solution for the 3-site Bloch-McConnell equations for R1rho data with R20A = R20B = R20C.
98
99
100 Indices
101 =======
102
103 The data structures used in this target function class consist of many different index types. These are:
104
105 - Ei: The index for each experiment type.
106 - Si: The index for each spin of the spin cluster.
107 - Mi: The index for each magnetic field strength.
108 - Oi: The index for each spin-lock offset. In the case of CPMG-type data, this index is always zero.
109 - Di: The index for each dispersion point (either the spin-lock field strength or the nu_CPMG frequency).
110
111
112 @keyword model: The relaxation dispersion model to fit.
113 @type model: str
114 @keyword num_param: The number of parameters in the model.
115 @type num_param: int
116 @keyword num_spins: The number of spins in the cluster.
117 @type num_spins: int
118 @keyword num_frq: The number of spectrometer field strengths.
119 @type num_frq: int
120 @keyword exp_types: The list of experiment types. The dimensions are {Ei}.
121 @type exp_types: list of str
122 @keyword values: The R2eff/R1rho values. The dimensions are {Ei, Si, Mi, Oi, Di}.
123 @type values: rank-4 list of numpy rank-1 float arrays
124 @keyword errors: The R2eff/R1rho errors. The dimensions are {Ei, Si, Mi, Oi, Di}.
125 @type errors: rank-4 list of numpy rank-1 float arrays
126 @keyword missing: The data structure indicating missing R2eff/R1rho data. The dimensions are {Ei, Si, Mi, Oi, Di}.
127 @type missing: rank-4 list of numpy rank-1 int arrays
128 @keyword frqs: The spin Larmor frequencies (in MHz*2pi to speed up the ppm to rad/s conversion). The dimensions are {Ei, Si, Mi}.
129 @type frqs: rank-3 list of floats
130 @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}.
131 @type frqs_H: rank-3 list of floats
132 @keyword cpmg_frqs: The CPMG frequencies in Hertz. This will be ignored for R1rho experiments. The dimensions are {Ei, Mi}.
133 @type cpmg_frqs: rank-2 list of floats
134 @keyword spin_lock_nu1: The spin-lock field strengths in Hertz. This will be ignored for CPMG experiments. The dimensions are {Ei, Mi}.
135 @type spin_lock_nu1: rank-2 list of floats
136 @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}.
137 @type chemical_shifts: rank-3 list of floats
138 @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}.
139 @type offset: rank-4 list of floats
140 @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}.
141 @type tilt_angles: rank-5 list of floats
142 @keyword r1: The R1 relaxation rates. This is only used for off-resonance R1rho models. The dimensions are {Ei, Si, Mi}.
143 @type r1: rank-3 list of floats
144 @keyword relax_times: The experiment specific fixed time period for relaxation (in seconds). The dimensions are {Ei, Mi}.
145 @type relax_times: rank-2 list of floats
146 @keyword scaling_matrix: The square and diagonal scaling matrix.
147 @type scaling_matrix: numpy rank-2 float array
148 @keyword recalc_tau: A flag which if True will cause tau_CPMG to be recalculated to remove user input truncation.
149 @type recalc_tau: bool
150 """
151
152
153 if model not in MODEL_LIST_FULL:
154 raise RelaxError("The model '%s' is unknown." % model)
155 if values == None:
156 raise RelaxError("No values have been supplied to the target function.")
157 if errors == None:
158 raise RelaxError("No errors have been supplied to the target function.")
159 if missing == None:
160 raise RelaxError("No missing data information has been supplied to the target function.")
161 if model in [MODEL_DPL94, MODEL_TP02, MODEL_TAP03, MODEL_MP05]:
162 if chemical_shifts == None:
163 raise RelaxError("Chemical shifts must be supplied for the '%s' R1rho off-resonance dispersion model." % model)
164 if r1 == None:
165 raise RelaxError("R1 relaxation rates must be supplied for the '%s' R1rho off-resonance dispersion model." % model)
166
167
168 self.model = model
169 self.num_params = num_params
170 self.num_spins = num_spins
171 self.num_frq = num_frq
172 self.exp_types = exp_types
173 self.values = values
174 self.errors = errors
175 self.missing = missing
176 self.frqs = frqs
177 self.frqs_H = frqs_H
178 self.cpmg_frqs = cpmg_frqs
179 self.spin_lock_nu1 = spin_lock_nu1
180 self.chemical_shifts = chemical_shifts
181 self.offset = offset
182 self.tilt_angles = tilt_angles
183 self.r1 = r1
184 self.relax_times = relax_times
185 self.scaling_matrix = scaling_matrix
186
187
188 self.back_calc = deepcopy(values)
189
190
191 self.experiment_type_setup()
192
193
194 self.num_disp_points = []
195 self.num_offsets = []
196 for ei in range(len(values)):
197 self.num_disp_points.append([])
198 self.num_offsets.append([])
199 for si in range(len(values[ei])):
200 self.num_disp_points[ei].append([])
201 self.num_offsets[ei].append([])
202 for mi in range(len(values[ei][0])):
203 self.num_disp_points[ei][si].append([])
204 self.num_offsets[ei][si].append([])
205
206
207 if len(offset[ei][si][mi]):
208 self.num_offsets[ei][si][mi] = len(self.offset[ei][si][mi])
209 else:
210 self.num_offsets[ei][si][mi] = 0
211
212
213 for oi in range(self.num_offsets[ei][si][mi]):
214 self.num_disp_points[ei][si][mi].append([])
215 if cpmg_frqs != None and len(cpmg_frqs[ei][mi][oi]):
216 self.num_disp_points[ei][si][mi][oi] = len(self.cpmg_frqs[ei][mi][oi])
217 elif spin_lock_nu1 != None and len(spin_lock_nu1[ei][mi][oi]):
218 self.num_disp_points[ei][si][mi][oi] = len(self.spin_lock_nu1[ei][mi][oi])
219 else:
220 self.num_disp_points[ei][si][mi][oi] = 0
221
222
223 self.scaling_flag = False
224 if self.scaling_matrix != None:
225 self.scaling_flag = True
226
227
228 self.end_index = []
229
230
231 self.end_index.append(self.num_exp * self.num_spins * self.num_frq)
232 if model in [MODEL_B14_FULL, MODEL_CR72_FULL, MODEL_NS_CPMG_2SITE_3D_FULL, MODEL_NS_CPMG_2SITE_STAR_FULL]:
233 self.end_index.append(2 * self.num_exp * self.num_spins * self.num_frq)
234
235
236 self.end_index.append(self.end_index[-1] + self.num_spins)
237 if model in [MODEL_IT99, MODEL_LM63_3SITE, MODEL_MMQ_CR72, MODEL_NS_MMQ_2SITE]:
238 self.end_index.append(self.end_index[-1] + self.num_spins)
239 elif model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]:
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 elif model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]:
243 self.end_index.append(self.end_index[-1] + self.num_spins)
244 self.end_index.append(self.end_index[-1] + self.num_spins)
245 self.end_index.append(self.end_index[-1] + self.num_spins)
246
247
248 if model in [MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL]:
249
250 self.Rr = zeros((2, 2), complex64)
251
252
253 self.Rex = zeros((2, 2), complex64)
254
255
256 self.RCS = zeros((2, 2), complex64)
257
258
259 self.R = zeros((2, 2), complex64)
260
261
262 if model in [MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL]:
263 self.r180x = r180x_3d()
264
265
266 if model in [MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL, MODEL_NS_MMQ_2SITE]:
267 self.M0 = zeros(2, float64)
268 if model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]:
269 self.M0 = zeros(3, float64)
270 if model in [MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL]:
271 self.M0 = zeros(7, float64)
272 self.M0[0] = 0.5
273 if model in [MODEL_NS_R1RHO_2SITE]:
274 self.M0 = zeros(6, float64)
275 if model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]:
276 self.M0 = zeros(9, float64)
277
278
279 if model in [MODEL_B14, MODEL_B14_FULL, 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]:
280
281 self.power = []
282 for ei in range(self.num_exp):
283 self.power.append([])
284 for mi in range(self.num_frq):
285 self.power[ei].append(zeros(self.num_disp_points[ei][si][mi][0], int16))
286 for di in range(self.num_disp_points[ei][si][mi][0]):
287
288 if isNaN(self.relax_times[ei][mi]):
289 self.power[ei][mi][di] = 0.0
290
291
292 else:
293 self.power[ei][mi][di] = int(round(self.cpmg_frqs[ei][mi][0][di] * self.relax_times[ei][mi]))
294
295
296 self.tau_cpmg = []
297 for ei in range(len(values)):
298 self.tau_cpmg.append([])
299 for mi in range(len(values[ei][0])):
300 self.tau_cpmg[ei].append(zeros(self.num_disp_points[ei][si][mi][0], float64))
301 for di in range(self.num_disp_points[ei][si][mi][0]):
302
303 if recalc_tau:
304 self.tau_cpmg[ei][mi][di] = 0.25 * self.relax_times[ei][mi] / self.power[ei][mi][di]
305 else:
306 self.tau_cpmg[ei][mi][di] = 0.25 / self.cpmg_frqs[ei][mi][0][di]
307
308
309 if spin_lock_nu1 != None:
310 self.spin_lock_omega1 = []
311 self.spin_lock_omega1_squared = []
312 for ei in range(len(values)):
313 self.spin_lock_omega1.append([])
314 self.spin_lock_omega1_squared.append([])
315 for mi in range(len(values[ei][0])):
316 self.spin_lock_omega1[ei].append([])
317 self.spin_lock_omega1_squared[ei].append([])
318 for oi in range(len(values[ei][0][mi])):
319 self.spin_lock_omega1[ei][mi].append([])
320 self.spin_lock_omega1_squared[ei][mi].append([])
321 self.spin_lock_omega1[ei][mi][oi] = 2.0 * pi * self.spin_lock_nu1[ei][mi][oi]
322 self.spin_lock_omega1_squared[ei][mi][oi] = self.spin_lock_omega1[ei][mi][oi] ** 2
323
324
325 if model in [MODEL_B14, MODEL_B14_FULL, 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]:
326 self.inv_relax_times = 1.0 / relax_times
327
328
329 if model == MODEL_NS_MMQ_2SITE:
330 self.m1 = zeros((2, 2), complex64)
331 self.m2 = zeros((2, 2), complex64)
332 elif model in [MODEL_NS_MMQ_3SITE, MODEL_NS_MMQ_3SITE_LINEAR]:
333 self.m1 = zeros((3, 3), complex64)
334 self.m2 = zeros((3, 3), complex64)
335 elif model == MODEL_NS_R1RHO_2SITE:
336 self.matrix = zeros((6, 6), float64)
337 elif model in [MODEL_NS_R1RHO_3SITE, MODEL_NS_R1RHO_3SITE_LINEAR]:
338 self.matrix = zeros((9, 9), float64)
339
340
341 if model == MODEL_NOREX:
342 self.func = self.func_NOREX
343 if model == MODEL_LM63:
344 self.func = self.func_LM63
345 if model == MODEL_LM63_3SITE:
346 self.func = self.func_LM63_3site
347 if model == MODEL_CR72_FULL:
348 self.func = self.func_CR72_full
349 if model == MODEL_CR72:
350 self.func = self.func_CR72
351 if model == MODEL_IT99:
352 self.func = self.func_IT99
353 if model == MODEL_TSMFK01:
354 self.func = self.func_TSMFK01
355 if model == MODEL_B14:
356 self.func = self.func_B14
357 if model == MODEL_B14_FULL:
358 self.func = self.func_B14_full
359 if model == MODEL_NS_CPMG_2SITE_3D_FULL:
360 self.func = self.func_ns_cpmg_2site_3D_full
361 if model == MODEL_NS_CPMG_2SITE_3D:
362 self.func = self.func_ns_cpmg_2site_3D
363 if model == MODEL_NS_CPMG_2SITE_EXPANDED:
364 self.func = self.func_ns_cpmg_2site_expanded
365 if model == MODEL_NS_CPMG_2SITE_STAR_FULL:
366 self.func = self.func_ns_cpmg_2site_star_full
367 if model == MODEL_NS_CPMG_2SITE_STAR:
368 self.func = self.func_ns_cpmg_2site_star
369 if model == MODEL_M61:
370 self.func = self.func_M61
371 if model == MODEL_M61B:
372 self.func = self.func_M61b
373 if model == MODEL_DPL94:
374 self.func = self.func_DPL94
375 if model == MODEL_TP02:
376 self.func = self.func_TP02
377 if model == MODEL_TAP03:
378 self.func = self.func_TAP03
379 if model == MODEL_MP05:
380 self.func = self.func_MP05
381 if model == MODEL_NS_R1RHO_2SITE:
382 self.func = self.func_ns_r1rho_2site
383 if model == MODEL_NS_R1RHO_3SITE:
384 self.func = self.func_ns_r1rho_3site
385 if model == MODEL_NS_R1RHO_3SITE_LINEAR:
386 self.func = self.func_ns_r1rho_3site_linear
387 if model == MODEL_MMQ_CR72:
388 self.func = self.func_mmq_CR72
389 if model == MODEL_NS_MMQ_2SITE:
390 self.func = self.func_ns_mmq_2site
391 if model == MODEL_NS_MMQ_3SITE:
392 self.func = self.func_ns_mmq_3site
393 if model == MODEL_NS_MMQ_3SITE_LINEAR:
394 self.func = self.func_ns_mmq_3site_linear
395
396
397 - def calc_B14_chi2(self, R20A=None, R20B=None, dw=None, pA=None, kex=None):
398 """Calculate the chi-squared value of the Baldwin (2014) 2-site exact solution model for all time scales.
399
400
401 @keyword R20A: The R2 value for state A in the absence of exchange.
402 @type R20A: list of float
403 @keyword R20B: The R2 value for state B in the absence of exchange.
404 @type R20B: list of float
405 @keyword dw: The chemical shift differences in ppm for each spin.
406 @type dw: list of float
407 @keyword pA: The population of state A.
408 @type pA: float
409 @keyword kex: The rate of exchange.
410 @type kex: float
411 @return: The chi-squared value.
412 @rtype: float
413 """
414
415
416 pB = 1.0 - pA
417 k_BA = pA * kex
418 k_AB = pB * kex
419
420
421 chi2_sum = 0.0
422
423
424 for ei in range(self.num_exp):
425
426 for si in range(self.num_spins):
427
428 for mi in range(self.num_frq):
429
430 r20_index = mi + si*self.num_frq
431
432
433 dw_frq = dw[si] * self.frqs[ei][si][mi]
434
435
436 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ:
437 aliased_dw = dw_frq
438 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ:
439 aliased_dw = dw_frq
440
441
442 r2eff_B14(r20a=R20A[r20_index], r20b=R20B[r20_index], pA=pA, pB=pB, dw=aliased_dw, kex=kex, k_AB=k_AB, k_BA=k_BA, ncyc=self.power[ei][mi], 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])
443
444
445 for di in range(self.num_disp_points[ei][si][mi][0]):
446 if self.missing[ei][si][mi][0][di]:
447 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di]
448
449
450 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0])
451
452
453 return chi2_sum
454
455
456 - def calc_CR72_chi2(self, R20A=None, R20B=None, dw=None, pA=None, kex=None):
457 """Calculate the chi-squared value of the Carver and Richards (1972) 2-site exchange model on all time scales.
458
459 @keyword R20A: The R2 value for state A in the absence of exchange.
460 @type R20A: list of float
461 @keyword R20B: The R2 value for state B in the absence of exchange.
462 @type R20B: list of float
463 @keyword dw: The chemical shift differences in ppm for each spin.
464 @type dw: list of float
465 @keyword pA: The population of state A.
466 @type pA: float
467 @keyword kex: The rate of exchange.
468 @type kex: float
469 @return: The chi-squared value.
470 @rtype: float
471 """
472
473
474 chi2_sum = 0.0
475
476
477 for si in range(self.num_spins):
478
479 for mi in range(self.num_frq):
480
481 r20_index = mi + si*self.num_frq
482
483
484 dw_frq = dw[si] * self.frqs[0][si][mi]
485
486
487 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])
488
489
490 for di in range(self.num_disp_points[0][si][mi][0]):
491 if self.missing[0][si][mi][0][di]:
492 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
493
494
495 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
496
497
498 return chi2_sum
499
500
502 """Calculate the chi-squared value of the 'NS CPMG 2-site' models.
503
504 @keyword R20A: The R2 value for state A in the absence of exchange.
505 @type R20A: list of float
506 @keyword R20B: The R2 value for state B in the absence of exchange.
507 @type R20B: list of float
508 @keyword dw: The chemical shift differences in ppm for each spin.
509 @type dw: list of float
510 @keyword pA: The population of state A.
511 @type pA: float
512 @keyword kex: The rate of exchange.
513 @type kex: float
514 @return: The chi-squared value.
515 @rtype: float
516 """
517
518
519 pB = 1.0 - pA
520 k_BA = pA * kex
521 k_AB = pB * kex
522
523
524 self.M0[1] = pA
525 self.M0[4] = pB
526
527
528 chi2_sum = 0.0
529
530
531 for si in range(self.num_spins):
532
533 for mi in range(self.num_frq):
534
535 r20_index = mi + si*self.num_frq
536
537
538 dw_frq = dw[si] * self.frqs[0][si][mi]
539
540
541 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])
542
543
544 for di in range(self.num_disp_points[0][si][mi][0]):
545 if self.missing[0][si][mi][0][di]:
546 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
547
548
549 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
550
551
552 return chi2_sum
553
554
556 """Calculate the chi-squared value of the 'NS CPMG 2-site star' models.
557
558 @keyword R20A: The R2 value for state A in the absence of exchange.
559 @type R20A: list of float
560 @keyword R20B: The R2 value for state B in the absence of exchange.
561 @type R20B: list of float
562 @keyword dw: The chemical shift differences in ppm for each spin.
563 @type dw: list of float
564 @keyword pA: The population of state A.
565 @type pA: float
566 @keyword kex: The rate of exchange.
567 @type kex: float
568 @return: The chi-squared value.
569 @rtype: float
570 """
571
572
573 pB = 1.0 - pA
574 k_BA = pA * kex
575 k_AB = pB * kex
576
577
578 self.Rex[0, 0] = -k_AB
579 self.Rex[0, 1] = k_BA
580 self.Rex[1, 0] = k_AB
581 self.Rex[1, 1] = -k_BA
582
583
584 self.M0[0] = pA
585 self.M0[1] = pB
586
587
588 chi2_sum = 0.0
589
590
591 for si in range(self.num_spins):
592
593 for mi in range(self.num_frq):
594
595 r20_index = mi + si*self.num_frq
596
597
598 dw_frq = dw[si] * self.frqs[0][si][mi]
599
600
601 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])
602
603
604 for di in range(self.num_disp_points[0][si][mi][0]):
605 if self.missing[0][si][mi][0][di]:
606 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
607
608
609 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
610
611
612 return chi2_sum
613
614
615 - 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):
616 """Calculate the chi-squared value for the 'NS MMQ 3-site' models.
617
618 @keyword R20A: The R2 value for state A in the absence of exchange.
619 @type R20A: list of float
620 @keyword R20B: The R2 value for state B in the absence of exchange.
621 @type R20B: list of float
622 @keyword R20C: The R2 value for state C in the absence of exchange.
623 @type R20C: list of float
624 @keyword dw_AB: The chemical exchange difference between states A and B in rad/s.
625 @type dw_AB: float
626 @keyword dw_BC: The chemical exchange difference between states B and C in rad/s.
627 @type dw_BC: float
628 @keyword dwH_AB: The proton chemical exchange difference between states A and B in rad/s.
629 @type dwH_AB: float
630 @keyword dwH_BC: The proton chemical exchange difference between states B and C in rad/s.
631 @type dwH_BC: float
632 @keyword pA: The population of state A.
633 @type pA: float
634 @keyword kex_AB: The rate of exchange between states A and B.
635 @type kex_AB: float
636 @keyword kex_BC: The rate of exchange between states B and C.
637 @type kex_BC: float
638 @keyword kex_AC: The rate of exchange between states A and C.
639 @type kex_AC: float
640 @return: The chi-squared value.
641 @rtype: float
642 """
643
644
645 pC = 1.0 - pA - pB
646 pA_pB = pA + pB
647 pA_pC = pA + pC
648 pB_pC = pB + pC
649 k_BA = pA * kex_AB / pA_pB
650 k_AB = pB * kex_AB / pA_pB
651 k_CB = pB * kex_BC / pB_pC
652 k_BC = pC * kex_BC / pB_pC
653 k_CA = pA * kex_AC / pA_pC
654 k_AC = pC * kex_AC / pA_pC
655 dw_AC = dw_AB + dw_BC
656 dwH_AC = dwH_AB + dwH_BC
657
658
659 self.M0[0] = pA
660 self.M0[1] = pB
661 self.M0[2] = pC
662
663
664 chi2_sum = 0.0
665
666
667 for ei in range(self.num_exp):
668
669 for si in range(self.num_spins):
670
671 for mi in range(self.num_frq):
672
673 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp
674
675
676 dw_AB_frq = dw_AB[si] * self.frqs[ei][si][mi]
677 dw_AC_frq = dw_AC[si] * self.frqs[ei][si][mi]
678 dwH_AB_frq = dwH_AB[si] * self.frqs_H[ei][si][mi]
679 dwH_AC_frq = dwH_AC[si] * self.frqs_H[ei][si][mi]
680
681
682 aliased_dwH_AB = 0.0
683 aliased_dwH_AC = 0.0
684 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ:
685 aliased_dw_AB = dw_AB_frq
686 aliased_dw_AC = dw_AC_frq
687 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ:
688 aliased_dw_AB = dwH_AB_frq
689 aliased_dw_AC = dwH_AC_frq
690 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ:
691 aliased_dw_AB = dw_AB_frq + dwH_AB_frq
692 aliased_dw_AC = dw_AC_frq + dwH_AC_frq
693 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ:
694 aliased_dw_AB = dw_AB_frq - dwH_AB_frq
695 aliased_dw_AC = dw_AC_frq - dwH_AC_frq
696 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ:
697 aliased_dw_AB = dw_AB_frq
698 aliased_dw_AC = dw_AC_frq
699 aliased_dwH_AB = dwH_AB_frq
700 aliased_dwH_AC = dwH_AC_frq
701 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ:
702 aliased_dw_AB = dwH_AB_frq
703 aliased_dw_AC = dwH_AC_frq
704 aliased_dwH_AB = dw_AB_frq
705 aliased_dwH_AC = dw_AC_frq
706
707
708 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])
709
710
711 for di in range(self.num_disp_points[ei][si][mi][0]):
712 if self.missing[ei][si][mi][0][di]:
713 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di]
714
715
716 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0])
717
718
719 return chi2_sum
720
721
722 - 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):
723 """Calculate the chi-squared value for the 'NS MMQ 3-site' models.
724
725 @keyword r1rho_prime: The R1rho value for all states in the absence of exchange.
726 @type r1rho_prime: list of float
727 @keyword dw_AB: The chemical exchange difference between states A and B in rad/s.
728 @type dw_AB: float
729 @keyword dw_BC: The chemical exchange difference between states B and C in rad/s.
730 @type dw_BC: float
731 @keyword pA: The population of state A.
732 @type pA: float
733 @keyword kex_AB: The rate of exchange between states A and B.
734 @type kex_AB: float
735 @keyword kex_BC: The rate of exchange between states B and C.
736 @type kex_BC: float
737 @keyword kex_AC: The rate of exchange between states A and C.
738 @type kex_AC: float
739 @return: The chi-squared value.
740 @rtype: float
741 """
742
743
744 pC = 1.0 - pA - pB
745 pA_pB = pA + pB
746 pA_pC = pA + pC
747 pB_pC = pB + pC
748 k_BA = pA * kex_AB / pA_pB
749 k_AB = pB * kex_AB / pA_pB
750 k_CB = pB * kex_BC / pB_pC
751 k_BC = pC * kex_BC / pB_pC
752 k_CA = pA * kex_AC / pA_pC
753 k_AC = pC * kex_AC / pA_pC
754 dw_AC = dw_AB + dw_BC
755
756
757 chi2_sum = 0.0
758
759
760 for si in range(self.num_spins):
761
762 for mi in range(self.num_frq):
763
764 r20_index = mi + si*self.num_frq
765
766
767 dw_AB_frq = dw_AB[si] * self.frqs[0][si][mi]
768 dw_AC_frq = dw_AC[si] * self.frqs[0][si][mi]
769
770
771 for oi in range(self.num_offsets[0][si][mi]):
772
773 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])
774
775
776 for di in range(self.num_disp_points[0][si][mi][oi]):
777 if self.missing[0][si][mi][oi][di]:
778 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
779
780
781 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
782
783
784 return chi2_sum
785
786
788 """Check the experiment types and simplify data structures.
789
790 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.
791 """
792
793
794 self.num_exp = len(self.exp_types)
795
796
797 if self.model in MODEL_LIST_MMQ:
798
799 self.r2eff_ns_mmq = []
800
801
802 for ei in range(self.num_exp):
803
804 if self.exp_types[ei] in [EXP_TYPE_CPMG_SQ, EXP_TYPE_CPMG_PROTON_SQ, EXP_TYPE_CPMG_DQ, EXP_TYPE_CPMG_ZQ]:
805 if self.model == MODEL_NS_MMQ_2SITE:
806 self.r2eff_ns_mmq.append(r2eff_ns_mmq_2site_sq_dq_zq)
807 else:
808 self.r2eff_ns_mmq.append(r2eff_ns_mmq_3site_sq_dq_zq)
809
810
811 elif self.exp_types[ei] in [EXP_TYPE_CPMG_MQ, EXP_TYPE_CPMG_PROTON_MQ]:
812 if self.model == MODEL_NS_MMQ_2SITE:
813 self.r2eff_ns_mmq.append(r2eff_ns_mmq_2site_mq)
814 else:
815 self.r2eff_ns_mmq.append(r2eff_ns_mmq_3site_mq)
816
817
818 else:
819
820 if self.model != MODEL_NOREX and self.model in MODEL_LIST_CPMG and self.exp_types[0] != EXP_TYPE_CPMG_SQ:
821 raise RelaxError("The '%s' CPMG model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0]))
822 if self.model != MODEL_NOREX and self.model in MODEL_LIST_R1RHO and self.exp_types[0] != EXP_TYPE_R1RHO:
823 raise RelaxError("The '%s' R1rho model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0]))
824 if self.model != MODEL_NOREX and self.model in MODEL_LIST_MQ_CPMG and self.exp_types[0] != EXP_TYPE_CPMG_MQ:
825 raise RelaxError("The '%s' CPMG model is not compatible with the '%s' experiment type." % (self.model, self.exp_types[0]))
826
827
829 """Target function for the Baldwin (2014) 2-site exact solution model for all time scales, whereby the simplification R20A = R20B is assumed.
830
831 This assumes that pA > pB, and hence this must be implemented as a constraint.
832
833
834 @param params: The vector of parameter values.
835 @type params: numpy rank-1 float array
836 @return: The chi-squared value.
837 @rtype: float
838 """
839
840
841 if self.scaling_flag:
842 params = dot(params, self.scaling_matrix)
843
844
845 R20 = params[:self.end_index[0]]
846 dw = params[self.end_index[0]:self.end_index[1]]
847 pA = params[self.end_index[1]]
848 kex = params[self.end_index[1]+1]
849
850
851 return self.calc_B14_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)
852
853
855 """Target function for the Baldwin (2014) 2-site exact solution model for all time scales.
856
857 This assumes that pA > pB, and hence this must be implemented as a constraint.
858
859
860 @param params: The vector of parameter values.
861 @type params: numpy rank-1 float array
862 @return: The chi-squared value.
863 @rtype: float
864 """
865
866
867 if self.scaling_flag:
868 params = dot(params, self.scaling_matrix)
869
870
871 R20 = params[:self.end_index[1]].reshape(self.num_spins*2, self.num_frq)
872 R20A = R20[::2].flatten()
873 R20B = R20[1::2].flatten()
874 dw = params[self.end_index[1]:self.end_index[2]]
875 pA = params[self.end_index[2]]
876 kex = params[self.end_index[2]+1]
877
878
879 return self.calc_B14_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)
880
881
883 """Target function for the reduced Carver and Richards (1972) 2-site exchange model on all time scales.
884
885 This assumes that pA > pB, and hence this must be implemented as a constraint. For this model, the simplification R20A = R20B is assumed.
886
887
888 @param params: The vector of parameter values.
889 @type params: numpy rank-1 float array
890 @return: The chi-squared value.
891 @rtype: float
892 """
893
894
895 if self.scaling_flag:
896 params = dot(params, self.scaling_matrix)
897
898
899 R20 = params[:self.end_index[0]]
900 dw = params[self.end_index[0]:self.end_index[1]]
901 pA = params[self.end_index[1]]
902 kex = params[self.end_index[1]+1]
903
904
905 return self.calc_CR72_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)
906
907
909 """Target function for the full Carver and Richards (1972) 2-site exchange model on all time scales.
910
911 This assumes that pA > pB, and hence this must be implemented as a constraint.
912
913
914 @param params: The vector of parameter values.
915 @type params: numpy rank-1 float array
916 @return: The chi-squared value.
917 @rtype: float
918 """
919
920
921 if self.scaling_flag:
922 params = dot(params, self.scaling_matrix)
923
924
925 R20 = params[:self.end_index[1]].reshape(self.num_spins*2, self.num_frq)
926 R20A = R20[::2].flatten()
927 R20B = R20[1::2].flatten()
928 dw = params[self.end_index[1]:self.end_index[2]]
929 pA = params[self.end_index[2]]
930 kex = params[self.end_index[2]+1]
931
932
933 return self.calc_CR72_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)
934
935
937 """Target function for the Davis, Perlman and London (1994) fast 2-site off-resonance exchange model for R1rho-type experiments.
938
939 @param params: The vector of parameter values.
940 @type params: numpy rank-1 float array
941 @return: The chi-squared value.
942 @rtype: float
943 """
944
945
946 if self.scaling_flag:
947 params = dot(params, self.scaling_matrix)
948
949
950 R20 = params[:self.end_index[0]]
951 phi_ex = params[self.end_index[0]:self.end_index[1]]
952 kex = params[self.end_index[1]]
953
954
955 chi2_sum = 0.0
956
957
958 for si in range(self.num_spins):
959
960 for mi in range(self.num_frq):
961
962 r20_index = mi + si*self.num_frq
963
964
965 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2
966
967
968 for oi in range(self.num_offsets[0][si][mi]):
969
970 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])
971
972
973 for di in range(self.num_disp_points[0][si][mi][oi]):
974 if self.missing[0][si][mi][oi][di]:
975 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
976
977
978 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
979
980
981 return chi2_sum
982
983
985 """Target function for the Ishima and Torchia (1999) 2-site model for all timescales with pA >> pB.
986
987 @param params: The vector of parameter values.
988 @type params: numpy rank-1 float array
989 @return: The chi-squared value.
990 @rtype: float
991 """
992
993
994 if self.scaling_flag:
995 params = dot(params, self.scaling_matrix)
996
997
998 R20 = params[:self.end_index[0]]
999 dw = params[self.end_index[0]:self.end_index[1]]
1000 pA = params[self.end_index[1]]
1001 tex = params[self.end_index[1]+1]
1002
1003
1004 pB = 1.0 - pA
1005
1006
1007 chi2_sum = 0.0
1008
1009
1010 for si in range(self.num_spins):
1011
1012 for mi in range(self.num_frq):
1013
1014 r20_index = mi + si*self.num_frq
1015
1016
1017 dw_frq = dw[si] * self.frqs[0][si][mi]
1018
1019
1020 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])
1021
1022
1023 for di in range(self.num_disp_points[0][si][mi][0]):
1024 if self.missing[0][si][mi][0][di]:
1025 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1026
1027
1028 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1029
1030
1031 return chi2_sum
1032
1033
1035 """Target function for the Luz and Meiboom (1963) fast 3-site exchange model.
1036
1037 @param params: The vector of parameter values.
1038 @type params: numpy rank-1 float array
1039 @return: The chi-squared value.
1040 @rtype: float
1041 """
1042
1043
1044 if self.scaling_flag:
1045 params = dot(params, self.scaling_matrix)
1046
1047
1048 R20 = params[:self.end_index[0]]
1049 phi_ex_B = params[self.end_index[0]:self.end_index[1]]
1050 phi_ex_C = params[self.end_index[1]:self.end_index[2]]
1051 kB = params[self.end_index[2]]
1052 kC = params[self.end_index[2]+1]
1053
1054
1055 rex_B = phi_ex_B / kB
1056 rex_C = phi_ex_C / kC
1057 quart_kB = kB / 4.0
1058 quart_kC = kC / 4.0
1059
1060
1061 chi2_sum = 0.0
1062
1063
1064 for si in range(self.num_spins):
1065
1066 for mi in range(self.num_frq):
1067
1068 r20_index = mi + si*self.num_frq
1069
1070
1071 frq2 = self.frqs[0][si][mi]**2
1072 rex_B_scaled = rex_B[si] * frq2
1073 rex_C_scaled = rex_C[si] * frq2
1074
1075
1076 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])
1077
1078
1079 for di in range(self.num_disp_points[0][si][mi][0]):
1080 if self.missing[0][si][mi][0][di]:
1081 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1082
1083
1084 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1085
1086
1087 return chi2_sum
1088
1089
1091 """Target function for the Luz and Meiboom (1963) fast 2-site exchange model.
1092
1093 @param params: The vector of parameter values.
1094 @type params: numpy rank-1 float array
1095 @return: The chi-squared value.
1096 @rtype: float
1097 """
1098
1099
1100 if self.scaling_flag:
1101 params = dot(params, self.scaling_matrix)
1102
1103
1104 R20 = params[:self.end_index[0]]
1105 phi_ex = params[self.end_index[0]:self.end_index[1]]
1106 kex = params[self.end_index[1]]
1107
1108
1109 chi2_sum = 0.0
1110
1111
1112 for si in range(self.num_spins):
1113
1114 for mi in range(self.num_frq):
1115
1116 r20_index = mi + si*self.num_frq
1117
1118
1119 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2
1120
1121
1122 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])
1123
1124
1125 for di in range(self.num_disp_points[0][si][mi][0]):
1126 if self.missing[0][si][mi][0][di]:
1127 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1128
1129
1130 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1131
1132
1133 return chi2_sum
1134
1135
1137 """Target function for the Meiboom (1961) fast 2-site exchange model for R1rho-type experiments.
1138
1139 @param params: The vector of parameter values.
1140 @type params: numpy rank-1 float array
1141 @return: The chi-squared value.
1142 @rtype: float
1143 """
1144
1145
1146 if self.scaling_flag:
1147 params = dot(params, self.scaling_matrix)
1148
1149
1150 R20 = params[:self.end_index[0]]
1151 phi_ex = params[self.end_index[0]:self.end_index[1]]
1152 kex = params[self.end_index[1]]
1153
1154
1155 chi2_sum = 0.0
1156
1157
1158 for si in range(self.num_spins):
1159
1160 for mi in range(self.num_frq):
1161
1162 r20_index = mi + si*self.num_frq
1163
1164
1165 phi_ex_scaled = phi_ex[si] * self.frqs[0][si][mi]**2
1166
1167
1168 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])
1169
1170
1171 for di in range(self.num_disp_points[0][si][mi][0]):
1172 if self.missing[0][si][mi][0][di]:
1173 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1174
1175
1176 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1177
1178
1179 return chi2_sum
1180
1181
1183 """Target function for the Meiboom (1961) R1rho on-resonance 2-site model for skewed populations (pA >> pB).
1184
1185 @param params: The vector of parameter values.
1186 @type params: numpy rank-1 float array
1187 @return: The chi-squared value.
1188 @rtype: float
1189 """
1190
1191
1192 if self.scaling_flag:
1193 params = dot(params, self.scaling_matrix)
1194
1195
1196 R20 = params[:self.end_index[0]]
1197 dw = params[self.end_index[0]:self.end_index[1]]
1198 pA = params[self.end_index[1]]
1199 kex = params[self.end_index[1]+1]
1200
1201
1202 chi2_sum = 0.0
1203
1204
1205 for si in range(self.num_spins):
1206
1207 for mi in range(self.num_frq):
1208
1209 r20_index = mi + si*self.num_frq
1210
1211
1212 dw_frq = dw[si] * self.frqs[0][si][mi]
1213
1214
1215 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])
1216
1217
1218 for di in range(self.num_disp_points[0][si][mi][0]):
1219 if self.missing[0][si][mi][0][di]:
1220 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1221
1222
1223 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1224
1225
1226 return chi2_sum
1227
1228
1230 """Target function for the Miloushev and Palmer (2005) R1rho off-resonance 2-site model.
1231
1232 @param params: The vector of parameter values.
1233 @type params: numpy rank-1 float array
1234 @return: The chi-squared value.
1235 @rtype: float
1236 """
1237
1238
1239 if self.scaling_flag:
1240 params = dot(params, self.scaling_matrix)
1241
1242
1243 R20 = params[:self.end_index[0]]
1244 dw = params[self.end_index[0]:self.end_index[1]]
1245 pA = params[self.end_index[1]]
1246 kex = params[self.end_index[1]+1]
1247
1248
1249 pB = 1.0 - pA
1250
1251
1252 chi2_sum = 0.0
1253
1254
1255 for si in range(self.num_spins):
1256
1257 for mi in range(self.num_frq):
1258
1259 r20_index = mi + si*self.num_frq
1260
1261
1262 dw_frq = dw[si] * self.frqs[0][si][mi]
1263
1264
1265 for oi in range(self.num_offsets[0][si][mi]):
1266
1267 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])
1268
1269
1270 for di in range(self.num_disp_points[0][si][mi][oi]):
1271 if self.missing[0][si][mi][oi][di]:
1272 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
1273
1274
1275 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
1276
1277
1278 return chi2_sum
1279
1280
1282 """Target function for the CR72 model extended for MQ CPMG data.
1283
1284 @param params: The vector of parameter values.
1285 @type params: numpy rank-1 float array
1286 @return: The chi-squared value.
1287 @rtype: float
1288 """
1289
1290
1291 if self.scaling_flag:
1292 params = dot(params, self.scaling_matrix)
1293
1294
1295 R20 = params[:self.end_index[0]]
1296 dw = params[self.end_index[0]:self.end_index[1]]
1297 dwH = params[self.end_index[1]:self.end_index[2]]
1298 pA = params[self.end_index[2]]
1299 kex = params[self.end_index[2]+1]
1300
1301
1302 pB = 1.0 - pA
1303 k_BA = pA * kex
1304 k_AB = pB * kex
1305
1306
1307 chi2_sum = 0.0
1308
1309
1310 for ei in range(self.num_exp):
1311
1312 for si in range(self.num_spins):
1313
1314 for mi in range(self.num_frq):
1315
1316 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp
1317
1318
1319 dw_frq = dw[si] * self.frqs[ei][si][mi]
1320 dwH_frq = dwH[si] * self.frqs_H[ei][si][mi]
1321
1322
1323 aliased_dwH = 0.0
1324 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ:
1325 aliased_dw = dw_frq
1326 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ:
1327 aliased_dw = dwH_frq
1328 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ:
1329 aliased_dw = dw_frq + dwH_frq
1330 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ:
1331 aliased_dw = dw_frq - dwH_frq
1332 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ:
1333 aliased_dw = dw_frq
1334 aliased_dwH = dwH_frq
1335 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ:
1336 aliased_dw = dwH_frq
1337 aliased_dwH = dw_frq
1338
1339
1340 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])
1341
1342
1343 for di in range(self.num_disp_points[ei][si][mi][0]):
1344 if self.missing[ei][si][mi][0][di]:
1345 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di]
1346
1347
1348 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0])
1349
1350
1351 return chi2_sum
1352
1353
1355 """Target function for no exchange.
1356
1357 @param params: The vector of parameter values.
1358 @type params: numpy rank-1 float array
1359 @return: The chi-squared value.
1360 @rtype: float
1361 """
1362
1363
1364 if self.scaling_flag:
1365 params = dot(params, self.scaling_matrix)
1366
1367
1368 R20 = params
1369
1370
1371 chi2_sum = 0.0
1372
1373
1374 for ei in range(self.num_exp):
1375
1376 for si in range(self.num_spins):
1377
1378 for mi in range(self.num_frq):
1379
1380 r20_index = mi + si*self.num_frq
1381
1382
1383 for oi in range(self.num_offsets[ei][si][mi]):
1384
1385 for di in range(self.num_disp_points[ei][si][mi][oi]):
1386 self.back_calc[ei][si][mi][oi][di] = R20[r20_index]
1387
1388
1389 for di in range(self.num_disp_points[ei][si][mi][oi]):
1390 if self.missing[ei][si][mi][oi][di]:
1391 self.back_calc[ei][si][mi][oi][di] = self.values[ei][si][mi][oi][di]
1392
1393
1394 chi2_sum += chi2(self.values[ei][si][mi][oi], self.back_calc[ei][si][mi][oi], self.errors[ei][si][mi][oi])
1395
1396
1397 return chi2_sum
1398
1399
1401 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations.
1402
1403 @param params: The vector of parameter values.
1404 @type params: numpy rank-1 float array
1405 @return: The chi-squared value.
1406 @rtype: float
1407 """
1408
1409
1410 if self.scaling_flag:
1411 params = dot(params, self.scaling_matrix)
1412
1413
1414 R20 = params[:self.end_index[0]]
1415 dw = params[self.end_index[0]:self.end_index[1]]
1416 pA = params[self.end_index[1]]
1417 kex = params[self.end_index[1]+1]
1418
1419
1420 return self.calc_ns_cpmg_2site_3D_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)
1421
1422
1424 """Target function for the full numerical solution for the 2-site Bloch-McConnell equations.
1425
1426 @param params: The vector of parameter values.
1427 @type params: numpy rank-1 float array
1428 @return: The chi-squared value.
1429 @rtype: float
1430 """
1431
1432
1433 if self.scaling_flag:
1434 params = dot(params, self.scaling_matrix)
1435
1436
1437 R20 = params[:self.end_index[1]].reshape(self.num_spins*2, self.num_frq)
1438 R20A = R20[::2].flatten()
1439 R20B = R20[1::2].flatten()
1440 dw = params[self.end_index[1]:self.end_index[2]]
1441 pA = params[self.end_index[2]]
1442 kex = params[self.end_index[2]+1]
1443
1444
1445 return self.calc_ns_cpmg_2site_3D_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)
1446
1447
1449 """Target function for the numerical solution for the 2-site Bloch-McConnell equations using the expanded notation.
1450
1451 @param params: The vector of parameter values.
1452 @type params: numpy rank-1 float array
1453 @return: The chi-squared value.
1454 @rtype: float
1455 """
1456
1457
1458 if self.scaling_flag:
1459 params = dot(params, self.scaling_matrix)
1460
1461
1462 R20 = params[:self.end_index[0]]
1463 dw = params[self.end_index[0]:self.end_index[1]]
1464 pA = params[self.end_index[1]]
1465 kex = params[self.end_index[1]+1]
1466
1467
1468 pB = 1.0 - pA
1469 k_BA = pA * kex
1470 k_AB = pB * kex
1471
1472
1473 chi2_sum = 0.0
1474
1475
1476 for si in range(self.num_spins):
1477
1478 for mi in range(self.num_frq):
1479
1480 r20_index = mi + si*self.num_frq
1481
1482
1483 dw_frq = dw[si] * self.frqs[0][si][mi]
1484
1485
1486 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])
1487
1488
1489 for di in range(self.num_disp_points[0][si][mi][0]):
1490 if self.missing[0][si][mi][0][di]:
1491 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1492
1493
1494 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1495
1496
1497 return chi2_sum
1498
1499
1501 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations using complex conjugate matrices.
1502
1503 This is the model whereby the simplification R20A = R20B is assumed.
1504
1505
1506 @param params: The vector of parameter values.
1507 @type params: numpy rank-1 float array
1508 @return: The chi-squared value.
1509 @rtype: float
1510 """
1511
1512
1513 if self.scaling_flag:
1514 params = dot(params, self.scaling_matrix)
1515
1516
1517 R20 = params[:self.end_index[0]]
1518 dw = params[self.end_index[0]:self.end_index[1]]
1519 pA = params[self.end_index[1]]
1520 kex = params[self.end_index[1]+1]
1521
1522
1523 return self.calc_ns_cpmg_2site_star_chi2(R20A=R20, R20B=R20, dw=dw, pA=pA, kex=kex)
1524
1525
1527 """Target function for the full numerical solution for the 2-site Bloch-McConnell equations using complex conjugate matrices.
1528
1529 @param params: The vector of parameter values.
1530 @type params: numpy rank-1 float array
1531 @return: The chi-squared value.
1532 @rtype: float
1533 """
1534
1535
1536 if self.scaling_flag:
1537 params = dot(params, self.scaling_matrix)
1538
1539
1540 R20 = params[:self.end_index[1]].reshape(self.num_spins*2, self.num_frq)
1541 R20A = R20[::2].flatten()
1542 R20B = R20[1::2].flatten()
1543 dw = params[self.end_index[1]:self.end_index[2]]
1544 pA = params[self.end_index[2]]
1545 kex = params[self.end_index[2]+1]
1546
1547
1548 return self.calc_ns_cpmg_2site_star_chi2(R20A=R20A, R20B=R20B, dw=dw, pA=pA, kex=kex)
1549
1550
1552 """Target function for the combined SQ, ZQ, DQ and MQ CPMG numeric solution.
1553
1554 @param params: The vector of parameter values.
1555 @type params: numpy rank-1 float array
1556 @return: The chi-squared value.
1557 @rtype: float
1558 """
1559
1560
1561 if self.scaling_flag:
1562 params = dot(params, self.scaling_matrix)
1563
1564
1565 R20 = params[:self.end_index[0]]
1566 dw = params[self.end_index[0]:self.end_index[1]]
1567 dwH = params[self.end_index[1]:self.end_index[2]]
1568 pA = params[self.end_index[2]]
1569 kex = params[self.end_index[2]+1]
1570
1571
1572 pB = 1.0 - pA
1573 k_BA = pA * kex
1574 k_AB = pB * kex
1575
1576
1577 self.M0[0] = pA
1578 self.M0[1] = pB
1579
1580
1581 chi2_sum = 0.0
1582
1583
1584 for ei in range(self.num_exp):
1585
1586 for si in range(self.num_spins):
1587
1588 for mi in range(self.num_frq):
1589
1590 r20_index = mi + ei*self.num_frq + si*self.num_frq*self.num_exp
1591
1592
1593 dw_frq = dw[si] * self.frqs[ei][si][mi]
1594 dwH_frq = dwH[si] * self.frqs_H[ei][si][mi]
1595
1596
1597 aliased_dwH = 0.0
1598 if self.exp_types[ei] == EXP_TYPE_CPMG_SQ:
1599 aliased_dw = dw_frq
1600 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_SQ:
1601 aliased_dw = dwH_frq
1602 elif self.exp_types[ei] == EXP_TYPE_CPMG_DQ:
1603 aliased_dw = dw_frq + dwH_frq
1604 elif self.exp_types[ei] == EXP_TYPE_CPMG_ZQ:
1605 aliased_dw = dw_frq - dwH_frq
1606 elif self.exp_types[ei] == EXP_TYPE_CPMG_MQ:
1607 aliased_dw = dw_frq
1608 aliased_dwH = dwH_frq
1609 elif self.exp_types[ei] == EXP_TYPE_CPMG_PROTON_MQ:
1610 aliased_dw = dwH_frq
1611 aliased_dwH = dw_frq
1612
1613
1614 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])
1615
1616
1617 for di in range(self.num_disp_points[ei][si][mi][0]):
1618 if self.missing[ei][si][mi][0][di]:
1619 self.back_calc[ei][si][mi][0][di] = self.values[ei][si][mi][0][di]
1620
1621
1622 chi2_sum += chi2(self.values[ei][si][mi][0], self.back_calc[ei][si][mi][0], self.errors[ei][si][mi][0])
1623
1624
1625 return chi2_sum
1626
1627
1629 """Target function for the combined SQ, ZQ, DQ and MQ 3-site MMQ CPMG numeric solution.
1630
1631 @param params: The vector of parameter values.
1632 @type params: numpy rank-1 float array
1633 @return: The chi-squared value.
1634 @rtype: float
1635 """
1636
1637
1638 if self.scaling_flag:
1639 params = dot(params, self.scaling_matrix)
1640
1641
1642 R20 = params[:self.end_index[0]]
1643 dw_AB = params[self.end_index[0]:self.end_index[1]]
1644 dw_BC = params[self.end_index[1]:self.end_index[2]]
1645 dwH_AB = params[self.end_index[2]:self.end_index[3]]
1646 dwH_BC = params[self.end_index[3]:self.end_index[4]]
1647 pA = params[self.end_index[4]]
1648 kex_AB = params[self.end_index[4]+1]
1649 pB = params[self.end_index[4]+2]
1650 kex_BC = params[self.end_index[4]+3]
1651 kex_AC = params[self.end_index[4]+4]
1652
1653
1654 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)
1655
1656
1658 """Target function for the combined SQ, ZQ, DQ and MQ 3-site linearised MMQ CPMG numeric solution.
1659
1660 @param params: The vector of parameter values.
1661 @type params: numpy rank-1 float array
1662 @return: The chi-squared value.
1663 @rtype: float
1664 """
1665
1666
1667 if self.scaling_flag:
1668 params = dot(params, self.scaling_matrix)
1669
1670
1671 R20 = params[:self.end_index[0]]
1672 dw_AB = params[self.end_index[0]:self.end_index[1]]
1673 dw_BC = params[self.end_index[1]:self.end_index[2]]
1674 dwH_AB = params[self.end_index[2]:self.end_index[3]]
1675 dwH_BC = params[self.end_index[3]:self.end_index[4]]
1676 pA = params[self.end_index[4]]
1677 kex_AB = params[self.end_index[4]+1]
1678 pB = params[self.end_index[4]+2]
1679 kex_BC = params[self.end_index[4]+3]
1680
1681
1682 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)
1683
1684
1686 """Target function for the reduced numerical solution for the 2-site Bloch-McConnell equations for R1rho data.
1687
1688 @param params: The vector of parameter values.
1689 @type params: numpy rank-1 float array
1690 @return: The chi-squared value.
1691 @rtype: float
1692 """
1693
1694
1695 if self.scaling_flag:
1696 params = dot(params, self.scaling_matrix)
1697
1698
1699 r1rho_prime = params[:self.end_index[0]]
1700 dw = params[self.end_index[0]:self.end_index[1]]
1701 pA = params[self.end_index[1]]
1702 kex = params[self.end_index[1]+1]
1703
1704
1705 pB = 1.0 - pA
1706 k_BA = pA * kex
1707 k_AB = pB * kex
1708
1709
1710 chi2_sum = 0.0
1711
1712
1713 for si in range(self.num_spins):
1714
1715 for mi in range(self.num_frq):
1716
1717 r20_index = mi + si*self.num_frq
1718
1719
1720 dw_frq = dw[si] * self.frqs[0][si][mi]
1721
1722
1723 for oi in range(self.num_offsets[0][si][mi]):
1724
1725 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])
1726
1727
1728 for di in range(self.num_disp_points[0][si][mi][oi]):
1729 if self.missing[0][si][mi][oi][di]:
1730 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
1731
1732
1733 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
1734
1735
1736 return chi2_sum
1737
1738
1740 """Target function for the R1rho 3-site numeric solution.
1741
1742 @param params: The vector of parameter values.
1743 @type params: numpy rank-1 float array
1744 @return: The chi-squared value.
1745 @rtype: float
1746 """
1747
1748
1749 if self.scaling_flag:
1750 params = dot(params, self.scaling_matrix)
1751
1752
1753 r1rho_prime = params[:self.end_index[0]]
1754 dw_AB = params[self.end_index[0]:self.end_index[1]]
1755 dw_BC = params[self.end_index[1]:self.end_index[2]]
1756 pA = params[self.end_index[2]]
1757 kex_AB = params[self.end_index[2]+1]
1758 pB = params[self.end_index[2]+2]
1759 kex_BC = params[self.end_index[2]+3]
1760 kex_AC = params[self.end_index[2]+4]
1761
1762
1763 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)
1764
1765
1767 """Target function for the R1rho 3-site numeric solution linearised with kAC = kCA = 0.
1768
1769 @param params: The vector of parameter values.
1770 @type params: numpy rank-1 float array
1771 @return: The chi-squared value.
1772 @rtype: float
1773 """
1774
1775
1776 if self.scaling_flag:
1777 params = dot(params, self.scaling_matrix)
1778
1779
1780 r1rho_prime = params[:self.end_index[0]]
1781 dw_AB = params[self.end_index[0]:self.end_index[1]]
1782 dw_BC = params[self.end_index[1]:self.end_index[2]]
1783 pA = params[self.end_index[2]]
1784 kex_AB = params[self.end_index[2]+1]
1785 pB = params[self.end_index[2]+2]
1786 kex_BC = params[self.end_index[2]+3]
1787
1788
1789 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)
1790
1791
1793 """Target function for the Trott, Abergel and Palmer (2003) R1rho off-resonance 2-site model.
1794
1795 @param params: The vector of parameter values.
1796 @type params: numpy rank-1 float array
1797 @return: The chi-squared value.
1798 @rtype: float
1799 """
1800
1801
1802 if self.scaling_flag:
1803 params = dot(params, self.scaling_matrix)
1804
1805
1806 R20 = params[:self.end_index[0]]
1807 dw = params[self.end_index[0]:self.end_index[1]]
1808 pA = params[self.end_index[1]]
1809 kex = params[self.end_index[1]+1]
1810
1811
1812 pB = 1.0 - pA
1813
1814
1815 chi2_sum = 0.0
1816
1817
1818 for si in range(self.num_spins):
1819
1820 for mi in range(self.num_frq):
1821
1822 r20_index = mi + si*self.num_frq
1823
1824
1825 dw_frq = dw[si] * self.frqs[0][si][mi]
1826
1827
1828 for oi in range(self.num_offsets[0][si][mi]):
1829
1830 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])
1831
1832
1833 for di in range(self.num_disp_points[0][si][mi][oi]):
1834 if self.missing[0][si][mi][oi][di]:
1835 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
1836
1837
1838 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
1839
1840
1841 return chi2_sum
1842
1843
1845 """Target function for the Trott and Palmer (2002) R1rho off-resonance 2-site model.
1846
1847 @param params: The vector of parameter values.
1848 @type params: numpy rank-1 float array
1849 @return: The chi-squared value.
1850 @rtype: float
1851 """
1852
1853
1854 if self.scaling_flag:
1855 params = dot(params, self.scaling_matrix)
1856
1857
1858 R20 = params[:self.end_index[0]]
1859 dw = params[self.end_index[0]:self.end_index[1]]
1860 pA = params[self.end_index[1]]
1861 kex = params[self.end_index[1]+1]
1862
1863
1864 pB = 1.0 - pA
1865
1866
1867 chi2_sum = 0.0
1868
1869
1870 for si in range(self.num_spins):
1871
1872 for mi in range(self.num_frq):
1873
1874 r20_index = mi + si*self.num_frq
1875
1876
1877 dw_frq = dw[si] * self.frqs[0][si][mi]
1878
1879
1880 for oi in range(self.num_offsets[0][si][mi]):
1881
1882 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])
1883
1884
1885 for di in range(self.num_disp_points[0][si][mi][oi]):
1886 if self.missing[0][si][mi][oi][di]:
1887 self.back_calc[0][si][mi][oi][di] = self.values[0][si][mi][oi][di]
1888
1889
1890 chi2_sum += chi2(self.values[0][si][mi][oi], self.back_calc[0][si][mi][oi], self.errors[0][si][mi][oi])
1891
1892
1893 return chi2_sum
1894
1895
1897 """Target function for the the Tollinger et al. (2001) 2-site very-slow exchange model, range of microsecond to second time scale.
1898
1899 @param params: The vector of parameter values.
1900 @type params: numpy rank-1 float array
1901 @return: The chi-squared value.
1902 @rtype: float
1903 """
1904
1905
1906 if self.scaling_flag:
1907 params = dot(params, self.scaling_matrix)
1908
1909
1910 R20A = params[:self.end_index[0]]
1911 dw = params[self.end_index[0]:self.end_index[1]]
1912 k_AB = params[self.end_index[1]]
1913
1914
1915 chi2_sum = 0.0
1916
1917
1918 for si in range(self.num_spins):
1919
1920 for mi in range(self.num_frq):
1921
1922 r20a_index = mi + si*self.num_frq
1923
1924
1925 dw_frq = dw[si] * self.frqs[0][si][mi]
1926
1927
1928 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])
1929
1930
1931 for di in range(self.num_disp_points[0][si][mi][0]):
1932 if self.missing[0][si][mi][0][di]:
1933 self.back_calc[0][si][mi][0][di] = self.values[0][si][mi][0][di]
1934
1935
1936 chi2_sum += chi2(self.values[0][si][mi][0], self.back_calc[0][si][mi][0], self.errors[0][si][mi][0])
1937
1938
1939 return chi2_sum
1940