J . Phys. Chem. 1989, 93, 8294-8299
8294
0 86 eV
0 96 eV
I
0.52 eV
Figure 13. Summary of energetics in Rb. sphaeroides RCs obtained from analysis of the decay of the 3Pstate. As discussed in the text, these are likely to be close to the equilibrium energies and may not reflect the energetics on a very short time scale which may be important for the initial charge separation step (cf. ref 104 and 105).
However, energetics alone will not suffice; the shapes and displacements of potential surfaces as well as an accurate calculation of the electronic coupling are also required. The subtlety of the problem is reflected in the question of the origin of unidirectional
electron transfer along the L-side of the RC, illustrated in Figure 1. As pointed out in the analysis of Michel-Beyerle et a1.,20it appears that slight differences in distances on the M- and L-sides and differential electrostatic stabilization due to the protein environment combine to produce unidirectional electron transfer, although this remains an open question. This subject can now be approached by recombinant DNA methods82-86and by measurements of electric field effects on reaction rates and pathways. Note Added in Proof. Zinth and co-workers have recently published time-resolved kinetic data which purport to show a two-step formation of P+BPheo- in Rb. sphaeroides RCs at room temperature [Holzapfel, W.; Finkele, U.; Kaiser, W.; Oesterhelt, D.; Scheer, H.; Stilz, H. U.; Zinth, W. Chem. Phys. Lett. 1989, 160, 1-71. This result differs from earlier quantitative analyses of the kinetics [cf. ref 12 and 131 and the electric field induced fluorescence anisotropy [cf. ref 281 at low tmeperature. Thus, the mechanism of the initial charge separation step is far from settled. Acknowledgment. Parts of this work were supported by grants from the NSF Biophysics Program and a Presidential Young Investigator Award to S.G.B. We thank Drs. Bixon, Holten, Jortner, Parson, and Small for providing manuscripts prior to publication and Drs. Bixon, Fischer, Friesner, Holten, Jortner, Marcus, Norris, Parson, Plato, Ulstrup, and Woodbury for helpful discussions.
ARTICLES Experimental Deslgn Optimization and Confidence Interval Estimation in Muttiexponentlal NMR Relaxation Measurements John D. Decatur and Thomas C. Farrar* Department of Chemistry, University of Wisconsin, Madison, Wisconsin 53706 (Received: April 3, 1989)
Numerical simulations were used to estimate the precision of parameters responsible for NMR spin-lattice relaxation in a coupled two-spin system where two relaxation mechanisms, dipole-dipoleand chemical shift anisotropy,are present. Parameter confidence intervals in nonlinear least-squares fitting, necessary for multiexponential relaxation, can be obtained without restrictive assumptions. The simulations allow one to investigate the effect of experimental variables on parameter precision without performing many lengthy experiments. It was found that the results of several pulse experiments, including the traditional inversion recovery, must be fit simultaneously to achieve reasonable precision. Precision depends upon many factors including the nucleus observed, how the recovery curves are sampled, the sample temperature stability, and the number of simultaneously estimated parameters.
Introduction
Measurement of spin-lattice relaxation times is a powerful method for obtaining information about the structure and dynamics of molecules in solution. If the relaxation is due to a single relaxation mechanism, the inversion recovery curve for either a coupled or decoupled spin system can usually be characterized, to a high degree of approximation, by a single-exponential time constant, T I . ] When two or more relaxation processes are present, recovery curves may be m~ltiexponential.~~~ At the high fields (1) Campbell, 1. D.; Freeman, R. J . Mugn. Reson. 1973, 11, 143-162. (2) Shimizu, H. J . Chem. Phys. 1964, 40(11), 3357-3364. (3) Mackor, E. L.; MacLean, C . Prog. Nucl. Magn. Reson. Spectrosc. 1967, 3, 129-157.
0022-3654/89/2093-8294$01.50/0
commonly used, the contribution to relaxation from the chemical shift anisotropy (CSA) may be as large as or even exceed the dipolar contribution. For PH03*- at 4.7 T (200 MHz 'H), interference terms between dipolar and CSA relaxation cause very noticeable effects in both the observed spectrum and in the inversion recovery curves. Previous work in this laboratory has demonstrated that these effects can be exploited to yield a simultaneous determination of all molecular parameters responsible for nuclear spin-rela~ation,~ For PHOj2-, these include the chemical shift anisotropy of the proton and the phosphorus, the H-P bond distance, and the molecular correlation time. (4) Farrar, T. C.; Locker, I. C . J . Chem. Phys. 1987,87(6), 3281-3287. See also: 2.Phys. Chem. (Munich) 1987, I S I , 25-33.
0 1989 American Chemical Society
Multiexponential N M R Relaxation Measurements
The Journal of Physical Chemistry, Vol. 93, No. 26, 1989 8295 design of N M R relaxation experiments for simple spin systems described by single-exponential recovery. The numerical Monte Carlo method is ideally suited for the more complicated case of two spin 1/2 nuclei where both dipolar and CSA interactions contribute to relaxation. This paper addresses the effect of pulse errors, radio frequency field imperfections, inclusion of additional simultaneously estimated parameters, temperature, and other experimental conditions upon the size of these confidence intervals. The set of N M R pulse experiments that produces the smallest uncertainty in the fitted parameters is also described. Theory The theory of longitudina1,relaxation in a coupled two spin system relaxed by dipolar and CSA mechanisms in the weak coupling limit has been presented el~ewhere.~.~ The relaxation of the four lines is described by four coupled differential equations:
d -(I) = R(1) dt where I is a vector of line intensities (which are simply differences in the populations of the energy levels that are given by Nonlinear least-squares methods are necessary for fitting multiexponential experimental recovery curves. In any leastsquares analysis, statistical fluctuations in measured data result in uncertainty in the fitted parameters. Error analysis in nonlinear fitting is not, however, as straightforward as in the linear case. Many of the nonlinear fitting routines available in the literature output confidence intervals but, in general, they rely on one or more assumptions: (1) the measurement error is normally distributed and a quantitative relationship exists between the covariance matrix and confidence intervals, (2) the uncertainties in the fitted parameters do not extend outside a region in which the model could be replaced by a suitable linearized one, and (3) the variance space is symmetric and can be estimated by extrapolating from the curvature at the minimum.s-6 The true confidence intervals of the fitted parameters may be calculated by the ratio of variances method.7 These methods, however, have no predictive utility because they give no informaton that could lead to improved confidence limits. A conceptually simple alternative is a Monte Carlo simulation of experimental data.5 The procedure is as follows: The experimental data are fit to obtain the set of ”best fit” parameters, Ai,,. With these parameters, “perfect data” are generated and random numbers are added to it from an appropriate distribution to mimic the measurement errors in the apparatus. With such random draws, synthetic data sets are constructed with exactly the same number of measured points and the same number of independent variables as the actual data set. For each synthetic data set, the least-squares minimization is performed, giving simulated measured parameters Ail, ..., A , If enough data sets are simulated, the desired probability distribution in M dimensions can be mapped out, here M is the number of independent parameters. This method, although costly in computer time, allows one to optimize both the computer-analysis method and the experimental conditions according to any desired criterion without actually performing many lengthy experiments. Several authors have discussed least-squares analysis of onedimensional relaxation in NMR.8-” Most recently, Weiss and Ferretti12 have reported analytical expressions for the optimal (5) Press, W. H.; Flannery, B. P.; Teukolsky, S. A. Vetterling, W. T. Numerical Recipes, The Art of Scientific Computing; Cambridge University
Press: New York, 1986. (6) Johnson, M. Biophys. J. 1983,44, 101-106. (7) Lampton, M.; Margon, B., Bowyer, S. Astrophys. J . 1976, 208, 177-1 90. (8) Hanssum, H.; Ruterjans, H. J . Magn. Reson. 1980,39, 65-78. (9) Celmins, A. J . Magn. Reson. 1982,50, 373-381. (10) LeiDert. T. K.:Marauardt. D. W. J. M a m . Reson. 1976.24.181-199. ( 1 l j W e k G.H.f Gupia, R. K.; Ferretti, .f A.; Becker, E: D. J . Magn. Reson. 1980,37, 369-319. ( 1 2 ) Weiss, G. H.; Ferretti, J. A. Frog. Nucl. Magn. Reson. Spectrosc. 1988,20, 317-335.
1’ = PI1 - P22
(2)
12 = P33 - P44
(3)
- P33 = P22 - P44
(4)
13 14
=
PlI
where the pii are the populations of the energy levels. For PH03”, which has a positive spin coupling constant (JpH = +565 Hz), Zland Z2 refer to the phosphorus low- and highfrequency transitions and I3and 1, refer to the proton low- and high-frequency transitions, respectively. The relaxation matrix, R, may be calculated from the theory developed by Redfield.” Elements of the relaxation matrix for the line intensity basis are listed are Table I. Elements of R depend upon the phosphorus CSA, Aup = (ull- ql)p, the proton CSA, AuH = (uli- c ~ )the~ , molecular correlation time, T ~ ,the dipolar parameter, p = YHYph/2rHp3, the two Larmor frequencies wH and wp,and two random-field spectral densities A and B (to account for residual intermolecular dipolar effects from solvent deuterium and cations). The solution to eq 1 is
+
+
+
I = cleXlrvl c2eX2‘v2 c3ex3‘v4
I(m)
(6)
In the least-squares analysis program, the relaxation matrix R is calculated based on estimated values of the Aui, p , and T ~ . This is diagonalized to give eigenvalues XI, Xi, ..., X4 and eigenvectors vl, v2, ..., v4. One eigenvalue is always zero because the total population of all four energy levels does not change. With the initial conditions and the equilibrium line intensities, the constants ~ 1 - ~are 4 obtained, and one may then calculate the predicted line intensities as a function of time. The criterion for the best fit is the set of parameters that minimizes the x2 function (7)
where the yi are measured line intensities, thefi’(zi) are calculated intensities depending on the fitted parameters, zi, and ui are estimates of the standard deviation of yi. The ai may be replaced by Wi,a general weighting function. This allows weighting more strongly sections of data that are sensitive to particular parameters. Experimental Section Computations were performed on a microVAX and a VAX 8650. The nonlinear least-squares minimization program uses the Gauss method of minimization. The uniform deviate (between 0 and 1, constant probability) random numbers were generated by an algorithm using two linear congruential generators designed (13) Redfield, A. G. Advances in Magnetic Resonance; Waugh, J. S., Ed.; Academic: New York, 1965; Vol. 1.
Decatur and Farrar
8296 The Journal of Physical Chemistry, Vol. 93, No. 26, 1989
TABLE 11: Parameter Values Determined from Experimental Data (Second Row of Column Heading) and Their Percent Uncertainties' As Determined by Simulation of 31PData expt AUH, ppm A
B C D E F G
H I
fHp,
A
i Cs,
X 10"
Abp, ppm
A,
rad/s
B, rad/s
(31Pcq)c (31Pinv)
('Hi",)
includedb
-3.6
1.514
8.5
-123.8
0.2
0.3
100
-97.7
-96
1 132 1-3 1-3 1-3 1-3 1-3 1,2,4d 1,2,4'
104 5.99 4.05 4.26 9.4 4.37 6.53 5.67 9.8
8.16 0.50 0.37 0.39 0.43
39 3.06 2.32 2.34 10.7
76 18.3 14.5 15.4 15 9.47
fix
fix
fix
fix
fix
0.12 0.21 0.27
2.44 0.97 8.75
fix
fix
0.15 0.16 0.15 0.15
0.41 0.5 0.39 0.31
fix fix fix fix
fix
160 49.4 37.5 37.7 74 6.22
fix
fix
65 4.78 3.62 3.66 8.6 0.78 1.27 1.52 6.35
20.0 54.2
4.11 7.78
fix
fix
fix
0.19
0.67
14.4
16 1.43 5.1
'Percent uncertainty defined as Iu/xI X 100%. bIdentity of recovery curves fit simultaneously. Numbers refer to NMR pulse experiments described in the Experimental Section. CEquilibriummagnetizations based on 100 for both nuclei. dFor experiment 4 (INEPT), the initial conditions require four additional parameters, one for each line. In this case they were also fixed to their ideal values. CLiked except the four initial conditions were allowed to vary (variances not shown). to eliminate any sequential correlations in the random numbers.5 Uniform deviates were then used to construct a normal (Gaussian) distribution characterized by an adjustable variance. Normally distributed random numbers were then added to the perfect synthetic intensity data simulating the measurement errors in our NMR spectrometer. Replicate measurements of the equilibrium line intensities can then be used to establish the magnitude of the variance of the noise to be added. In all simulations, except where noted, the standard deviation of the noise was 1%. No error was added to the discrete delay times used to sample the recovery. A set of 25 delay times was used to sample each recovery. When more than one recovery curve was fit simultaneously, the same set of delays was used for each recovery curve. For each simulation, the least-squares minimization process was repeated 250 times in order to adequately map out the probability distribution. The PHOj2- sample and N M R spectrometer have been described b e f ~ r e .PHOj2~ was dissolved in a 1:3 mixture of D 2 0 and ethylene glycol-d6 to make a 0.2 M solution. The basicity was adjusted to approximately pH 12. Excellent temperature control, required for the 36 h necessary to perform these experiments, was provided by a home-built digital controller based on an IBM PC-XT.I4 The temperature was stable to *0.040 O C for up to 48 h or more. Several different N M R pulse experiments were simulated corresponding to different boundary conditions for eq 6 . They are as follows: (1) broad-band (nonselective) inversion of jlP followed by observation of jlP; (2) broad-band inversion of 'H followed by observation of 3'P; (3) broad-band inversion of both jlP and IH followed by observation of 31P;(4) selective inversion of (INEPT pulse sequence) of a single 'H line followed by observation of jlP; (5) broad-band inversion of 31Pfollowed by observation of 'H; (6) broad-band inversion of 'H followed by observation of 'H; (7) broad-band inversion of both 31Pand 'H followed by observation of 'H. All population inversions were accomplished by composite p ~ 1 s e s . l ~Peak heights were used for intensities.
Results and Discussion Typical simulated recovery curves resulting from the different NMR pulse experiments are shown in Figure I . The parameters minimizing the x2 function for actual experimental data are shown in Table 11. A typical result of a simulation using the above experimental best fit parameters is shown for rHPvs Aup in Figure 2. Only two dimensions of the multiparameter space can be shown at a time. The shape of the distribution indicates that these two parameters are highly correlated. Although most available fitting routines output a correlation matrix that indicates the degree of correlation between parameters, Monte Carlo simulation provides a graphical illustration. Correlation between parameters affects the choice of the confidence region shape. When little or no (14) Farrar, T. C.; Sidky, E. Y.; Decatur, J. D. J . Ma@. Reson., in press. ( I 5 ) Levitt, M . H. Prog. Nucl. Magn. Reson. Spectrosc. 1986, 18,61-122.
P. OBSERLE
INVERT
a
7
/--
,/
so!
INVERT
C
P
P,
OBSERVE
H
I:CL
i. OBSERLE
IYVEQT
d
11
I60
i
l o o m
1
--
I
0
3 0 r
I eo I
0
e
__L__i~_
2.c
3.0
-AI"
c
Figure 1. Simulated PH03" recovery curves for different boundary conditions and observing different nuclei. Parameters values were those given in Table 11. Key: (a) experiment 1, (b) experiment 2, (c) experiment 5, (d) experiment 6.
t.
,
y
,
,
i , ,
so0
, ,
,
,
, ,I I 55c
,
,1
F H 30NO D1ST4r\(lE
Figure 2. Result of Monte Carlo simulation for T H and ~ Aup when recovery curves from experiment 2 only were fitted. The projection of the elliptical confidence region onto each axis determines the confidence interval for each parameter. Correlation increases the parameter confidence interval (-) relative to those when they are uncorrelated (--).
correlation is present, circular confidence regions define the probability distribution. When correlation is present, ellipses are drawn that include 68%, 90%, ... of the points corresponding to the desired percent confidence level. In more than two dimensions, ellipsoids are constructed. The projection of the ellipse onto each parameter axis determines the confidence interval. The effect of correlation is to increase the confidence interval over that determined from a single dimension. In such cases, an accurate,
Multiexponential NMR Relaxation Measurements
The Journal of Physical Chemistry, Vol. 93, No. 26, 1989 8297
TABLE 111: Systematic Error Introduced by Fixing Initial Magnetizations to Incorrect Valueso conditn Aa,, ppm fHp, A T ~ s, X 10" Aap, ppm A, rad/s B, rad/s b -3.60 1.514 8.50 -123.8 0.20 0.30 C -3.44 1.512 8.96 -1 19.3 0.13 0.30 -113.2 0.085 0.0489 1.489 8.94 d -3.38
("p,) 100 99.98 99.43
('lpinv)
-100 -99.5 fix
('Hinv) -100
fix fix
Values are averages of 250 simulations. Pulse experiments 1-3 were fit simultaneously in each case. bCorrect values. CSyntheticdata generated with (Hinv) -0.92(Hq). During fitting (Hinv)was fixed to ideal value ( ( H i n v )= -l.O(Hq)). dSynthetic data generated with (Pin,) = -0.9(P,) and (Hinv)= -l.O(Hq). During fitting (Pin")fixed to ideal value. (Hinv)also fixed to ideal value. independent measurement of one of the parameters would greatly increase the accuracy and precision of the other parameters: Correlation results when two or more parameters are partially redundant and the available experimental data cannot separate the parameters fully. The power of this simulation method lies in the ability to investigate what factors will reduce confidence limits or correlation between parameters. For simplicity, we used the one-dimensional standard deviation of the parameters as the criterion for which set of pulse experiments or experimental conditions produce the smallest uncertainty in the fitted parameters. The number of boundary conditions or the number of N M R pulse experiments, which are fit simultaneously, is very important. Table I1 shows results for different boundary conditions. Fitting the recovery curve from a single boundary condition, regardless of which one, produces unacceptably poor confidence intervals, indicating a shallow minimum in the xz space; Le., there are too many combinations of three exponentials capable of describing the recovery curve from a single boundary condition. Simultaneously fitting two or three boundary conditions under the same experimental conditions greatly improves the precision. For example, fitting any two of the experiments simultaneously reduces the standard deviation on Aap from over 65% (row A) to less than 4.78% (row B). Including a third experiment provides a marginal reduction to 3.6% (row C). This is probably only due to the 50% increase in the number of data points, analogous to signal averaging. In all cases, 25 points were used to define each recovery curve. The signal to noise (S/N) ratio of the experimenta data is crucial. The above numbers were based on data with 1% standard deviation. Reducing the standard deviation of the noise added onto the data to 0.5% produced a corresponding decrease of one-half in the confidence intervals of all parameters. Even further reduction in noise provides a greater than linear increase in the precision of the parameters as shown for proton data in Table 111. It is well-known that incomplete inversion pulses or incorrect values for the equilirium magnetization lead to erroneous results in ordinary inversion recovery experiments.I6 The so-called three-parameter fit, which allows fitting of both the inversion and equilibrium values, has been shown to give more accurate TI values than a single parameter fits1' Allowing the boundary conditions (inversion and equilibrium magnetizations) in our experiments to be fit would be ideal. In the broad-band inversion experiments, both lines in each spectrum can be assumed to be equally inverted because the transmitter is placed in the middle of the doublet. Three additional parameters are then required: two for the inverted magnetization of each nucleus and one for the equilibrium value. An equilibrium value is required for only the lines of one of the spins since the equilibrium value for the other spin is fixed by the ratio of their Larmor frequencies. Estimating additional parameters simultaneously, however, can induce correlations and reduce the precision of the estimates of all parameters. Consequently, it is important to examine the benefits of including these parameters, which contain no explicit chemical information. When boundary conditions were allowed to vary in the simulations, several effects were noticed (Table 11). First, allowing the equilibrium and inversion values for the phosphorus to vary, while keeping the proton inversion value constant (row D), resulted in small variances in these parameters. Furthermore, the confidence intervals on the other parameters
were not significantly increased. This indicates the model is sensitive to the equilibrium and phosphorus inversion parameters: they are uncorrelated with other parameters and are easily fit. Allowing the inverted magnetization of the unobserved nucleus, the proton, to be fit, however, resulted in a large variance in its value (row E). Moreover, it induced much larger confidence intervals into all other fitted parameters. The correlation matrix reveals strong correlation between the proton inversion value and most other parameters. Clearly, it would be desirable to *fixn the proton inversion value. Composite inversion pulses produce nearly ideal inversion, but it is difficult to measure precisely its value while the spectrometer is set to observe phosphorus. How large are the induced errors if the proton inversion value is fixed to an incorrect value? Systematic error is not detected by least-squares analysis, and confidence intervals are unaffected. Systematic error, or bias, is introduced, however, into all parameters. Table 111shows results where "perfect data" were generated with the proton inversion value incorrect by 8% and then subsequently fit with the proton value fixed to its "ideal" value. The average of the 250 best fit values of the other parameters then contain errors up to 5%. Results for an error in the fixed phosphorus inversion value are shown in Table 111. The resulting errors are slightly greater. Fortunately, a 10%error in the inversion value is much larger than is observed experimentally when inversion pulses are composite. Composite pulses typically achieve better than 98% inversion. Systematic error resulting from imperfect boundary conditions can be expected to be small if inversion pulses contain only a few percent error. Increased precision is possible when an independent estimate of a parameter is available. If there is good reason to assume a parameter will not change over different experimental conditions, such as temperature, the average of several determinations will yield a more precise value. In either case, the parameter can be removed from the simultaneous estimation. When rHPand T, (which can be related to temperature by an Arrhenius relationship) are known and held constant, a substantial increase in the precision is realized (row F, Table 11). In the present case, the precision of Aap increases by more than a factor of 10. Under conditions for which the random field terms, A and B, can be held constant, an increase in precision is also noticed (row G, Table 11). It is advantageous to use a solvent system, if possible, that is free of magnetic moments, thus eliminating dipolar intermolecular relaxation contributions. Other interactions also described by random-field contributions, however, must also be absent in order to remove A and B.I8 Because the INEPT experiment provides unique boundary conditions and increased signal intensity, it is tempting to include it as one of the set of initial conditions. The INEPT sequence involves several pulses for each nucleus. This results in a boundary condition quite unrelated to the boundary conditions of the other experiments. For example, for all of the other boundary conditions phosphorus inversion pulses were fit with a single parameter regardless of whether a proton inversion pulse was also applied or which nucleus was subsequently observed. In addition, the initial intensities of the four lines, when prepared by INEPT, all have different values. To fit the boundary conditions requires four additional parameters. One would desire a constraint that relates the intensities of the four lines, but in the general case of unequal pulse errors on each nucleus, we know of no such constraint.
(16) Crouch, R.; Hurlert, S.;Ragouzeos, A. J . Magn. Reson. 1982, 49, 371-382. (17) Sass, M.; Ziessow, D.J . Magn.Reson. 1977, 25, 263-276.
12, 79.
(18) Vold, R. L.; Vold, R. R. Prog. Nucl. Magn. Reson. Spectrosc. 1978,
8298
Decatur and Farrar
The Journal of Physical Chemistry, Vol. 93, No. 26, 1989
TABLE I V Percent Uncertainties in Parameters Determined from Simulations of 'H Data" with Different Amounts of Added Noise B, rad/s ("P,) Aup, ppm A , rad/s noise std dev expt included Ab,, ppm rHp, A T ~ s, X 10" ("Pin,) 1 .O%
0.258% 0.0662% I .O%
35.4 6.95 1.77 9.02
5-7 5-7 5-7 1-3'
3.07 0.74 0.19 0.43
20.4 4.68 1.18 8.92
197 47.0 12.3 66.6
45.0 6.37 1.61 1.52
115 32.4 8.47 15.1
1.28 0.03 0.008 0.153
'For quantitative comparison between 'H and 31P,all recovery curves were sampled with nonlinear T spacings up to line (see text).
5tnUll for
data similarly sampled, for comparison.
('Hi"") 0.3 1 0.08 0.02 13.1
8.78 2.24 0.57 0.263
the slower recovering
TABLE V Effect of Linear and Quadratic r-Time Spacing on the Percent Uncertainty of Parameters"
linear spacing If
1tmtll
3tnui1 5tnu11 7tnu11 9t,,11
A'H,
ppm
18.1 9.52 10.1 11.2 15.4
THp.
A
0.887 0.360 0.343 0.407 0.530
T'c, S
X
quadratic spacing 10''
21.9 9.84 10.5 11.5 15.3
ppm 20. I 7.69 7.77 8.91 12.8
AUp,
ppm 17.8 9.28 9.02 9.49 9.66
A'H~
rHp,
A
1 .oo
0.453 0.430 0.420 0.417
X 10'' 21.8 9.68 8.92 9.23 9.40
Tc, S
AUp,
ppm
20.3 7.99 7.52 7.73 7.84
"All simulations included experiments 1-3.
The results of including the INEPT experiment are shown in Table 11. Comparing simultaneously fitting experiments 1-3 and 1,2, and INEPT, each with all parameters variable, shows a small increase in the precision in the estimates of most the parameters when INEPT is included (row I). When all initial intensities are fixed to their ideal values, inclusion of INEPT provides the best precision of any experimental combination (row H). The multiple pulses required, however, produce initial intensities after the INEPT sequence that are probably less ideal than for the other experiments. The fact that the total time involved for the composite pulse trains and evolution times, is relatively long, and may be a nonnegligible fraction of the relaxation time, is an additional problem. The nucleus whose magnetization is sampled after the inversion pulses affects the parameter precision. When observing different nuclei, differences in S/N must be considered. For a fixed spectral width (Hz) and given field strength: (S/N)H/(S/N)p c: ( 7 H / r p ) 3 = 15.08. For comparison, consider simultaneously fitting experiments 1-3 versus experiments 5-7. These experiments correspond to identical sets of boundary conditions, or inversion pulses, but subsequently observing a different nucleus, phosphorus or proton. The proton results, shown in Table IV, are for simulations in which 1%, 0.258%, and 0.0662% (the theoretical improvement over phosphorus data) noise has been added to synthetic data. Shown also is the result from phosphorus data with 1 % noise added. Clearly, precision depends upon the relative S/N of the proton and phoshorus data. If the added noise for the two nuclei is equal, phosphorus data provides much better precision. If the noise added to proton is reduced to its theoretical amount (0.0662%), proton data is superior. For 0.258% added noise for proton and 1% added noise for phosphorus, the precision obtainable from each nucleus is similar. This implied superiority of proton data has not, however, been achieved in practice. Although the spectral S / N is better for proton, the reproducibility of replicate equilibrium magnetizations has been equal for the two nuclei. This suggests that other factors such as the stability of the spectrometer and/or the Bo homogeneity are responsibe for the small (less than 1%) scatter i n replicates rather than simply thermal noise in the coil. This may be true, however, only for large, relatively concentrated samples that provide large S/N ratios. Our sample provided phosphorus spectral S/N greater than 500 after 64 scans. One could intuitively explain the greater precision obtainable from phosphorus data for equal S/N by considering the relative size of Aup and AuH. Aup is more than 14 times as great as AuH in absolute terms (Hz) and 34.5 times larger in parts per million, leading to greater interference effects in the phosphorus spectrum. Indeed, visually comparing the recovery curves of experiments 1-3 and 5-7 one sees greater differential effects between the two lines in the phosphorus sectrum than in the proton spectrum. The equivalent sensitivity of the phosphorus data is fortuitous from an experimental viewpoint. The coil design of most probes
is such that the proton coil is longer than the phosphorus (Xnucleus) coil and thus covers a greater volume of the sample. Since the radio frequency homogeneity falls off rapidly outside the coil, that part of the sample outside the phosphorus coil would not have been completely inverted following the 180' phosphorus pulse. That same part of the sample, however, would be detected by the proton observe pulse. The actual boundary condition for that part of the sample would then be far from ideal. The opposite experiment, Le., proton inversion followed by phosphorus observe actually benefits from this coil arrangement since proton magnetization is completely inverted spacially beyond what is detected by the phosphorus coil. This is crucial because for this experiment it is difficult to measure the inverted proton magnetization under the same conditions as it was created. Reverse probes, which are becoming widely available for reverse-detection 2D experiments, have the proton coil inside a larger X coil. These would be useful for experiments in which the proton is observed. The manner in which the recovery curves are sampled can affect the attainable precision. In a recent paper on the optimal design of relaxation exeriments, Weiss and Ferretti12 investigated various spacings of the sampling times ( T times) on the precision of TI for single-exponential recoveries. Following their ideas, the spacings investigated were restricted to those following the form ti = 21 + (tf - t J ( i - 1 n- 1 where n is the number of times, and tl and ti are the final and initial times, respectively. Two values of r were examined: r = 1 for linear spacing and r = 2 for quadratic spacing. The results from fitting experiments 1-3 simultaneously with n = 25 for each experiment are shown in Table V. For these and all subsequent simulations, all boundary conditions were simultaneously fit with the other parameters. Several effects were noticed. For linear spacing, the value of tf strongly influences the precision of all parameters. The values of tf examined ranged from 1 to 9 times the point where the magnetization equals zero (tn,,ll) for the more slowly recovering line in the recovery curve for experiment 1. Maximum precision was found when tf was near 3t,,ll. For quadratic spacing, precision was much less dependent upon t p Maximum precision was found when tf was near 5tnull. The best precision obtainable by both methods was very similar. Use of quadratic spacing is to be recommended due to its relative insensitivity to tf and thus to the experimenter's guess for fnull. These findings are very similar to the results of Weiss and Ferretti. The temperature dependence of molecular parameters is of particular interest. The isotropic chemical shifts of many nuclei, e.g., I3C, 31P,and I9F, are known to have a large temperature dependence. For PHO,", which has C3, symmetry, there are two components of the chemical shift tensor u,, and uI. These components may be temperature dependent, which may lead to a temperature-dependent CSA. The precision of fitted parameters may, however, also be temperature dependent. In order to assess
-)
Multiexponential N M R Relaxation Measurements
t-
z Q
KO
" 9
40
L
3 y 20
3.00
10.00
I1.oo
-loa CORRELATION TIME
Figure 3. Percent uncertainty of parameters at different correlation times. Key: Abp (X), A ~ H rHp (A),2, (+I, MNV) (0).
(oh
the temperature dependence of the CSA, we need to know what happens to the confidence intervals at different temperatures. The value of the molecular correlation time changes substantially with temperature. The precision is expected to depend on the magnitude of all the parameters and thus may change with temperature. The results of simulations1at different correlation times are illustrated in Figure 3. In order to equally sample all recovery curves, 5t,,,i for each correlation time was chosen as tp Optimum precision for all parameters is found with correlation times such that 2.0 X lo4 I T , I 3.0 X The center of this range corresponds to uHrc= 1. At progressively shorter correlation times (UT, < 1) precision decreases drastically and then increases again. This is probably related to the relative magnitudes of the preexponential constants, cI, cq, ..., c4. Often their magnitudes are such that the recovery curves are essentially described by a biexponential. In the region of poor precision, however, the three exponential components contribute substantially. Thus, care must be taken to adjust the correlation time, through temperature or solvent viscosity, to a value that produces reasonable precision. The effect of temperature instability on parameter precision may be investigated. Once an activation energy for reorientation has been determined, the uncertainty in T , may be related to temperature instability. For our sample, E, = 7.51 kcal/mol. For AT = 0.05 K at 273 K, ATc/Tc= 0.0025. Instead of adding noise only to the perfect intensity data, random noise can be added to the correlation time used to generate each synthetic data set. Thus, each of the 250 simulated recovery curves contain not only intensities modulated by noise but the correlation time responsible for each recovery curve modulated by noise as well. The standard deviation of noise added to the correlation time can be chosen to equal that produced by the temperature uncertainty. The results indicate that uncertainty caused by temperature variation becomes significant when it approaches the magnitude of the AT, already present resulting from noise on the intensities. This additional uncertainty in I , also induces uncertainty in all other parameters. To achieve 1% uncertainty in AT^, AT,/T, induced from temperature variation should be small, 0.25% or less. This requires temperature stability of 0.05 K or better. This neglects the additional uncertainty caused by spinspin relaxation effects. The line widths and thus peak heights depend upon temperature and the effects of temperature variation on precision will enter by modulating the peak heights as well as 7,. Since the magnitude of the CSA interaction depends on the magnetic field strength, the precision obtainable can be expected to be field dependent. Simulations were carried out at 2.35, 4.7, and 8.46 T because these fields are available. The precision was
The Journal of Physical Chemistry, Vol. 93, No. 26, 1989 8299 similar at all three fields and showed a dependence on correlation time. Finally, a general weighting function Wi, used in the x2 minimization (eq 7), allows different sections of the recovery curves to be weighted more strongly. Inspection of Figure 1 indicates that the two lines in the recovery are indistinguishable at small T times. Weighting shorter T times less strongly than longer T times might naively be expected to enhance the differential information present when the two recovery curves diverge. The resulting precision, however, was much reduced. Short times are important, as in the single-exponential case, because there the rate of change of the magnetization is greatest. Short times are also important for the estimation of the initial inverted magnetization. Conclusions
Numerical simulation is a powerful tool to estimate confidence intervals and optmize experimental design in NMR relaxation studies when the longitudinal spin relaxation is multiexponential. It is clear from the work reported here that in order to obtain accurate information about the molecular parameters upon which the relaxation depends, it is essential to simultaneously fit data from at least two independent sets of boundary conditions. Inclusion of the INEPT sequence would greatly improve the attainable precision if the initial intensities for all four lines did not need to be simultaneously estimated. In an INEPT experiment, however, additional parameters must simultaneously be estimated, and these estimates reduce the precision on all parameters resulting in only marginal improvement when including the INEPT experiment. The nucleus observed is also important. In PHO:-, longitudinal relaxation recovery curves for phosphorus, which has a large CSA, show more pronounced differential relaxation than the proton, which has a rather small CSA. The phosphorus recovery curves are consequently more sensitive to changes in the molecular parameters and the phosphorus data result in greater precision of fitted parameters. Sampling the recovery curves quadratically provides optimum precision and lower sensitivity to the experimenters' guess for the rate of recovery. Sources of systematic error, such as incorrect pulse widths, may also be estimated. If composite inversion pulses are used, errors can be expected to be less than 2%. For a typical, modern NMR spectrometer, such as the one used in our laboratory, the drift in the measured intensity of a given spectral line over the time required to collect all of the data is less than OS%, and the temperature stability over this same time period (up to 36 h) is better than *0.035 "C. Under these conditions if all parameters including the inversions (boundary conditions) and the A and B terms are fit, it is possible to obtain simultaneouly, at a given magnetic field value and a given temperature, values for the two CSA's and the molecular correlation time to a precision of about 4%. The bond distance is obtained to a precision of about 0.2%. If the boundary conditions are known and stable, as they have been in the experiments that we have conducted, they may be fixed. In this event, the precision on the two CSA's and the correlation time increases to about 2% and the precision on the bond distance improves to about 0.1%. Acknowledgment. We gratefully acknowledge the support of the National Science Foundation for the support of this research (NSF Grant No. CHE-8802373). Registry No. PO,HZ-, 15477-76-6.
