206
Anal. Chem. 1981, 53,206-213
explain the sluggish response near the electrode, since concentrations averaged over 25 pm change more slowly than those within 3 pm of the electrode as shown in Figure 6. A reasonable conclusion is that the wave properties of the light prevent a resolution of 3 pm because the light diffracts near the electrode surface, allowing light from a wider region than 3 pm to enter the slit. This effect should be less severe for shorter wavelengths, but apparently the difference was not observable for the two wavelengths (515 and 633 nm) used here. One would predict that a shorter electrode would yield better results, since the diffracted light would travel a shorter distance and the slit would sample a smaller increment of the diffusion layer. It was observed experimentally (21)that a shorter electrode gave better results, but the short (0.04 cm) electrode was exceedingly difficult to align, and the small scale made mechanical operations very difficult. However, it can be concluded that the wave properties of the light ultimately degrade resolution when attempting to monitor events on a scale of only a few tens of wavelengths. The same effect also limits the response time of the absorbance, since a wider effective sampling window leads to sluggish response. In addition, there is a trade-off between path length and spatial resolution, since longer electrodes worsen the diffraction problem. Even without spatial resolution, however, a longer electrode could be used to enhance sensitivity and lower detection limits. In conclusion, a spectroelectrochemical experiment using light parallel to the electrode represents a geometry with some attractive features. High sensitivity results from long optical path length, and regions within -20 pm of the electrode surface can be monitored. Concentration vs. distance profiles can be constructed with diffraction of the light being the
ultimate limit on spatial resolution and response time. Applications of the technique to problems involving reactions of electrogeneratedspecies and more complex mass transport phenomena are in progress.
LITERATURE CITED Winograd, N.; Kuwana, T. In “ElectroanalyticalChemistry”; Bard, A. J., Ed.; Marcel Dekker: New York, 1974; Vol. 4. Kuwana, T. Ber. Bunsenges. Phys. Chem. 1973, 77, 858. Heineman, W. R. Anal. Chem. 1978, 50, 390A. Winograd, N.; Kuwana, T. J. Elecfroanal. Chem. 1976, 23, 333. McCreery, R. L.; Prulksma, R.; Fagan. R. Anal. Chem. 1079, 57, 749. Muller, R. H. Adv. Electrochem. Ehsctrochem. Eng. 1073, 9 , 281. Srinivasan, V. S. A&. Electrochem. Elecfrochem. Eng. 1973, 9 , 389. IkeshoJi,T.; Seklne, T. Denkl Kagaku 1977, 45, 575. Pruiksma, R.; McCreery, R. L. Anal. Chem. 1979, 51, 2253. Tyson, J. F.; West, T. S. Talanta 1980, 27, 335. Referemce 9, p 2254. Hecht, E. “Schaum Outline Serles, Optics”; McGraw-HIII: New York, 1975; p 197. Adams, R. N. “Electrochemlstryat Solid Electrodes”;Marcel Dekker: New York, 1989; p 218. Bishop, E.; Hartshorn, L. Q. Analyst (London) 1971, 96, 26. Feldberg. S. W. In “Electroanalytlcal Chemistry”; Bard, A. J., Ed.; Marcel Dekker: New York, 1970; Vol. 3. Dovlchl, N. J.; Harrls, J. M. Anal. Chem. 1978, 57, 728. Reference 12, pp 159-181. Mackev. L. Ph.D. Thesis. The Ohlo State Universltv, Columbus. OH. 1975, p 135. McIntyre, J. D. E. A&. Elecfrochem. Electrochem. Eng. 1073, 9 , 61. Srlnivasan, V. S. A&. Elecfrochem, Electrochem. Eng. 1973, 8, 331. Pruiksma, R. Ph.D. Dissertation, The Ohlo State University, Columbus, OH, 1980.
RJXEIVFDfor review September 10,1980. Accepted November 18,1980. This work was supported by The National Science Foundation (CHE 7828068) and The National Institute of Mental Health (MH-28412).
Linearization Methods for First-Order Kinetic Analysis Lowell M. Schwartz Department of Chemistry, Unlversity of Massachusetts- Boston, Boston, Massachusetts 02 125
Three linearization methods are presented for determlnlng the parameters of a first-order klnetlc decay superposed on a flnlte sum of power terms. Although primary consideration is glven to the estimation of the rate constant and Its statistical uncertalnty, the preexponential factor and the coefflcients In the polynomial are discussed as well. The three methods depend upon the freedom of the analyst to preprogram the tlme schedule of recording data polnts, and recording at equally spaced tlme Intervals Is a convenient and suitable program. Special cases of these three general methods have been previously suggested by Guggenhelm, ref 1, by Gutfreund and Sturtevant, ref 7, by Kezdy, Jaz, and Bruylants, ref 2, by Mangelsdorf, ref 3, by Swinbourne, ref 4, and by Glick, Brubacher, and Leggett, ref 8.
In this paper the phrase “kinetic analysis” refers to the procedure whereby measurements are made of the decay variable z vs. time t , and these data are used to estimate the magnitudes and the statistical uncertainties of unknown parameters. The simplest exponential decay scheme is char-
acterized by the well-known nonlinear expression z = y exp(-lzt) (1) where the unknown parameters are y, the preexponential factor, and k , the rate constant, or l/k, the first-order lifetime. We will assume here that the analyst is quite certain of the functional form of the time dependence of z. He knows that eq 1 or one of the generalizations to be discussed is the proper model. While the problem of deducing the kinetic order or the functional form of the model is sometimes of legitimate concern in kinetic investigations, this problem will not be addressed here. Linearization of eq 1 is not a prerequisite for calculating the parameters since several alternative procedures are available, but this is no doubt the most popular technique. Logarithmic transformation yields In z = In y - kt (2) which is a straight line on In z vs. t coordinates with the rate constant becoming the negative slope and the preexponential factor being calculated by exponentiatingthe intercept on the In z axis. In many instances the primary decay variable z is not measured directly but rather some other physical property
0003-2700/81/0353-0206$01.00/00 1981 American Chemlcal Society
ANALYTICAL CHEMISTRY, VOL. 53, NO. 2, FEBRUARY 1981
related to the primary variable is measured. A simple case is when the primary variable is a chemical concentration but measurement is made of the absorption of incident radiation a t a particular wavelength and this property is assumed proportional to the concentration via Beer’s law. In this case the relationship between measurement #J and primary variable z is #J = KZ, where K is the constant of proportionality. If the kinetic analysis seeks to evaluate the preexponential factor y in eq 1,K must be known a priori since measurement of #J rather than z leads to the new preexponentialfactor KY so that y alone cannot be found from #J vs. t data. A knowledge of K is not required to determine k however. If #J is proportional to z the logarithmictransformation still succeeds in linearizing the problem, but this is not generally true for other #J vs. z relationships. Even if #J is linearly related to z, i.e., of the form #J = KZ + #Ja,where 4- as well as K are constants, the #J vs. t relationship is 4 = KY exp(-kt) 4(3) and cannot be linearized by logarithmictransformation unless a nonzero #Jais known a priori and its value subtracted from each #J datum. Since #Jm is in reality the value of 4 at infinite time, it may be impossible or inconvenient to measure directly. Two linearizing procedures exist to circumvent the problem. One is due to Guggenheim (1) and the other was reported independently in three accounts Kezdy et al. (2),Mangelsdorf (31, and Swinbourne (4). This second procedure will be denoted by the abbreviation KMS for the first authors’ names and in order of publication date. Both procedures are special cases of the more general treatment to be offered here. Rather than a simple linear relationship between #J and z, consider the possibility that a time-dependent relationship exists of the form
+
where
is a polynomial of degree M with constant coefficients cp In these equations K is assumed to be known a priori so that the constant KY will henceforth be written simply as y. M will be assumed to be known also but the coefficients cj will be treated as unknown parameters whose values may be estimated from the data along with y and k. Equation 3 can be seen to be the particular case of a zero-degree polynomial Po(t) = co = #Jm. The literature of physical and biological science abounds with applications of this form. In most cases the analysis seeks to determine lz alone with #Jm and y existing as nuisance parameters but some applications seek these parameters as well. In constant-current coulometric analysis involving an uncomplicated reduction of a single analyte, the instantaneous cell current i(t) is expected to vary with time according to the Lingane equation (5) i(t) = io exp(-kt) (5) The analyst wishes to know the totalcharge qm= io/k required to reduce the analyte quantitatively. Therefore, both the preexponential factor io and the rate constant must be determined from i vs. t measurements. Furthermore, if statistical uncertainty is to be estimated for qmin this way, the variances of io and 12 as well as their covariance must be calculated. Alternatively the analyst may choose to record the instantaneous output q ( t ) of an integrator circuit for which the model equation becomes q ( t ) = q m [ l- exp(-kt)]
(6) which has the form of eq 3 with qmidentified with both -y
207
and Thus, in principle, calculating either of these parameters serves to determine qm,but evaluation of k is of secondary importance. Frequently, the primary reduction current of eq 5 is superposed on a constant background current ib, the continuous faradaic current (6). Then the total instantaneous current i(t) + ib is of the form of eq 3 and the instantaneous integrator output has the form q ( t ) = qm[l- exp(-kt)] i$ (7) which can be identified with eq 4 with the polynomial of degree M = 1, y = -qm, co = qm,and c1 = ib. Other examples of kinetic analyses which have the form of eq 4 with M = 1 arise in studies of enzyme catalysis. Gutfreund and Sturtevant (7) describe a modified Guggenheim method and Glick, Brubacher, and Leggett (8) describe a modified KMS method of data treatment to handle these cases. Both modifications are special applications of the general procedures which follow. While most applications of kinetic analysis are interested in the exponential term in eq 4, the methods to be described here are useful in cases where 4 vs. t is regarded as a signal P&) of primary interest which is superposed on a nuisance transient expressed by the exponential term. Once k is evaluated by any one of the following methods, the problem of finding the remaining parameters is linear and thus easier to attack. A later section will deal with this problem.
+
THREE GENERAL TIME-LAG LINEARIZATION PROCEDURES
Both the Guggenheim and KMS methods rely on the freedom of the analyst to select the discrete times t, at which he records the discrete measurements4,. These methods both require that data points be recorded in pairs separated by a fixed time-lag 7. While this requirement does not imply that all t, are equally spaced in time, it is frequently convenient to do so by means of an automatic recording device. Then the fixed time-lag requirement can be satisfied by treating measurements in pairs whose members are separated in time by a certain number of spacings. The generalized methods have similar requirementsand options. They require that data be recorded in accordance with a particular time schedule which depends on the degree M of the polynomial. But the requirement can always be satisfied by selecting an array of points from the whole set of points recorded at a fixed time spacing. The linearizations come about because of the property of the exponential function g(8) exp(-kO) that any linear combination of these in the form Zlv& + 17) factors into a product of two functions such that one is a function o f t and the other of T . This property is expressed by Zlv& + / T ) = g ( t ) Z p l exp(-Zk.r) (8) and where it is understood that 7 and vl are constank3 and that 1 are integers. Method 1. Consider a kinetic analysis where the model equation involves a polynomial of known degree, say, P2(t) so that
4 = yg(t) + co + c1t + C 2 t 2 (9) and data points 4, vs. t, are recorded. Corresponding to the four parameters y, co, cl, and c2 which appear explicitly, we will seek four Coefficients vl, 1 = 1-4, of a linear combination A2 of four data points & , 1 = 1-4, separated by a fixed time-lag 7 with the requirement that all terms stemming from Pzvanish in Az. The linear combination of points is 4
A2
Cvdt = yxlvg(tn) + 2,vJ’2(t,) 1=1
(10)
and when the measurement times t 2 = tl + 7,t3 = tl + 27 and
208
ANALYTICAL CHEMISTRY, VOL. 53, NO. 2, FEBRUARY 1981
+ 37 are substituted, this becomes A2 = Cvg(t1 + ( I - 1 ) ~ +)u1[co + cltl + ~ 2 t ? ] t 1=1 t 4 = tl
4
+ c1(tl + 7)+ c2(t12 + 2t17 + 791 4+ Cl(t1 + 27) + C2(tl2 + 4tlT + 4T2)] + u ~ [ c O+ cl(t1 + 37) + cp(&t2+ 6t17 + 9 ~ ~ ) (11) ]
U~[CO
U3[Co
In order to make all the polynomial terms vanish, the strategy is to make each of the coefficients of powers of T vanish separately and this leads to (VI
+ u2 + u3 + vq)(co + cltl + czt12) = 0 for TO (v2 + 2v3 + 3v4)(c1 + 2c2tl) = 0 for 7’ (v2 + 4u3 + 9u4)(c2) = 0 for
which is the Guggenheim (I) procedure. Therefore, it seems appropriate to refer to the method of eq 14 as “generalized Guggenheim”. The modified Guggenheim method of Gutfreund and Sturtevant (7)is eq 14 with M = 1. Method 2. An alternative procedure can be derived by starting with eq 13 which expresses the relationship between a linear combination AMof a group of M 2 data points and tl, the time of measurement of the first of these points. Consider two such groups, both characterized by time-lag T , one combined into and starting at time tl(l) and another combined into AM(^) and starting a t time t1(2).The time interval (to be called a time “displacement”) between these two starting times will be denoted by #, Le., t1(2)- t;’) = 7’. With this notation, the ratio is
+
A M ( ~ ) / A ~=( ~exp(-W) )
T~
But since tl and c j are arbitrary, the requirement reduces to a set of simultaneous linear equations in uI u2 + u3 + u4 = -v1
+ 2u3 + 3u4 = 0 u2 + 4u3 + 9u4 = 0
vz
which can be solved for u2, u3, u4 only relative to some arbitrary value of ul. The choice of ul = 1 leads, with the aid of Cramer’s rule, to the solution vector ( u l , u2, us, v4) = (1,-3, 3, -1). By following this same procedure for model equations with other polynomial degrees, there results for Po@)the solution vector ( u l , us) = (1,-1) and for Pl(t) the vector (vl, u2, u3) = (1,-2, 1). In general the terms stemming from a polynomial of degree M a r e made to vanish by taking a linear combination of L = M + 2 data points with coefficient values
which are binomial coefficients with alternating signs. When this is done, the nonvanishing terms in AM are a linear combination of time-lagged exponentials M+2
AM = 7 C vg(tl + ( l - UT) 1=1
which according to the property embodied in eq 8 factors conveniently into M+2
AM = y exp(-ktl) C u@-l 1=1
M+2
AM
= In [y C v@-’] I= 1
- kt1
(13)
(14)
Recalling that tl is the time of measurement of the first of a group of L = M 2 data points, the kinetic analysis proceeds by seeking the straight line In AM vs. tl of eq 14 and so requires the recording of at least 2L 4,, vs. t,, data points in order to form at least two such groups. If statistical uncertainties of the kinetic parameters are to be estimated from the data, many more than two groups are needed. As with the simple linearization of eq 2, the slope of the In AMvs. tl line determines k and the intercept on the In AM axis determines y since ul, T , and k can be used to evaluate the summation which multiplies it. When the polynomial is zero degree (constant term only), the model equation is 4 = y exp(-kt) + co which is the same as eq 3. The linearization eq 14 then utilizes (ul, v2) = (1,-1) and so is In ($1 - Cpz) = In [y(exp(-k7) - l)] - k t ,
+
E’
(15)
which can be seen by writing eq 13 for each linear combination and then dividing the two. Both AM(^) and Ay(l) are easily calculated from recorded data points and known coefficients from eq 12 so that the ratio yields E’ and, hence, k. Again Statisticalconsiderationssuggest that more than a single pair of (M + 2)-member groups should be recorded with the starting times within each pair displaced by 7’. A series of such pairs might be treated by plotting AM(^) vs. which according to eq 15 yields a straight line passing through the origin with slope E’. Equation 15 might be regarded as an extension of a simple procedure applied to eq 1which is simply an exponential term with no superposed polynomial: Measure a series of zi1) at times t i l )and ~ “ ( ~at1 times t i 2 )= t i 1 )t T’. Then calculate the corresponding set of ratios zi2)/zi1) = E,,’ or plot z i 2 ) vs. z,,(I) to extract the slope. However, we are not aware of a reported application to a case involving at least a zero-degree polynomial and would be happy to be informed of such cases by readers. Meanwhile, we will refer to eq 15 as “method 2”. Method 3. As with method 1, we will explain method 3 by illustrating the development of equations specific to a model with degree two polynomial and then generalize the result to degree M. Here we seek a linear combination K 2 of three data points with the first recorded a t time tl and the others time lagged by T such that all polynomial terms vanish except those that multiply T ~ By . analogy with eq 11 3
+ clt1 + c2tl2] + Uz[co + C l ( t 1 + + C 2 ( t l 2 + 2t17 + + V3[Co + Cl(t1 + 2 7 ) + C2(tl2 + 4t17 + 4T2)]
K2 = ~ C ~ g (+t (1I - 1 ) ~+)V1[Co 1=1
where E = exp(-k7). Linearization of AMwith respect to the independent variable tl is achieved by taking logarithms In
E
7)
T2)]
leads to a pair of simultaneous equations which is identical with the set for the linear combination Al and has solution vector (ul, u2, us) = (1, -2, 1). With these coefficients K2 becomes 3
KZ(l) TCVg(t1 + (I - 1)T) + c2 1=1
where C2= 2c2?. A second linear combination KJ2)is identical in all respects with Kz(’) except that the three points begin T’ later than tl, i.e., at tl + T / . The coefficients U I are the same and so 3
K2(2)= yCvg(t1 +T’ 1=1
+ (I - 1 ) ~ )+ C2
Subtracting C2from both and KJ2),taking the ratio of these results, and utilizing the factorization property of eq 8 yield
ANALYTICAL CHEMISTRY, VOL. 53,
This rearranges to
K2(2)= C2(l -E’)+ E’K,(l) which is linear in K2@)vs. K 2 ( l )and thus a plot on these coordinates is a straight line from which k can be calculated from the slope E’and the coefficient c2 from the intercept on the Kk2)axis. This result is generalized to any degree polynomial model equation by recognizing that the ( M + 1)member linear combination KM is identical with the ( M 1)-member combination AM-^ so that the linearization can be written
+
AM-^(') = c M ( 1 - E? + E’AM-“’)
(16)
and where reference to eq 12 with L = M - 1 provides the binomial coefficients for the linear combinations except for the special case hl= 4 when M = 0. The general expression for C M is M+l
CM =
cflM
1=1
~ [ ( l l- ) M
M L. 1
(17)
and several examples are Co = co, C1 = -c17, C2 = 2c2r2,and c3 = -6c3r3. A zero-degree polynomial leads to
4(’) = cO[l-- exp(-k~’)] + #l) exp(-k7’) which is the KMS linearization, and so we will refer to the procedure of eq 16 as “generalized KMS”. The modification suggested by Glick et al. (8)is an application of eq 16 with M = 1 and with 7’ = 7.
COMPARISONS This section summarizes some of the factors which enter a decision of selecting a procedure best suited to a particular application. While this discussion will be limited to the three linearizations described above, it should be kept in mind that other general techniques exist: nonlinear regression (9) which may be used with both equally spaced and unequally spaced data, and Cornell’s method of partial totals (IO)which may be used only with equally spaced data. The generalized Guggenheim and generalized KMS methods both produce two-parameter lines. Both slopes yield the rate constant k but the Guggenheim intercept yields the preexponential factor while the KMS intercept yields cM,the highest order coefficient in P&). Method 2 produces a single-parameter line whose slope yields k . An analyst seeking to determine a parameter in eq 4 which is not calculated directly from a straight-line parameter is obliged to perform a secondary calculation which will be discussed later. It is useful now to introduce a common notation for the lines to be drawn or calculated according to the various methods. This is the familiar yi vs. x , for data points plotted on an ordinate axis y and an abscissa axis x . With degree M polynomial, the generalized Guggenheim implies y = AM and x = tl. Method 2 implies y = AM@)and x = AM(1),and the generalized KMS implies y = and x = AM-^(^). It will be helpful if other terminology and nomenclature conventions are made clear also. To distinguish the two types of data points being discussed, we will refer to 4, vs. t, data as “recorded data” and use N to denote the total number of these gathered in a single kinetic run. Recorded data are then processed into y I vs. xi “plotted data” and I denotes the total number of these in a run. We will further assume that the experimenter takes the convenient recourse of recording data a t equally spaced time intervals At. Fixed numbers of At intervals are used to setup the equally spaced intervals within linear combinations and the equally spaced intervals 7’ between initial times tl of x1 and yl. At is a “spacing”, r is a “time-lag”, and 7’ is a “displacement”. The integer p will
NO. 2,
FEBRUARY 1981
209
denote the number of spacings per time lag, i.e., 7 = pAt, and the integer p ’ will denote the number of spacings per displacement, Le., 7’ = p’At. For a given method and a given number of points, N , many options may exist for setting the time lags and displacements. Each such option involves selecting particular p and p’integers, and we will refer to these selections as “programs”. The large number of possible programs for each of the three methods makes a full consideration a rather formidable undertaking for this single article. Therefore, for the purpose of reducing the total number of alternatives we will impose the restriction that every recorded point be used once but only once in some linear combination. Except for purposes of this condensation there is no a priori reason why in certain cases a program might omit some recorded points or might use some points in two or more linear combinations. An illustration will serve to show the possible programs existing even with the restriction imposed and this can be seen in Figure 1which diagrams all the programs for a polynomial degree M = 2 analysis when N = 24 points are recorded. The ( M N = (2,24) configurations lead to different sets of programs all adhering to the same restriction. Since N is frequently under the control of the experiment, this number may be adjusted to make more or fewer programs possible. In some circumstances the analyst wishing to minimize calculational effort might seek only a crude estimate of the rate constant k by establishing the straight-line graphically, i.e., by plotting the points y1vs. x i and fitting a line by eye. Without use of statistical principles there are few criteria to help select the method or the program here. The choice perhaps could be based on that which produced the greatest number of plotted points from the same number of recorded points. Figure 1 shows that methods 1,2, and 3 produce 6-, 3-, and &point lines, respectively, for that configuration so that the generalized Guggenheim method is best under this criterion. This same conclusion will hold generally regardless of the ( M N configuration, but the choice of program is more difficult to make. Program 1in Figure 1maximizes the range of plotted points along the x = tl axis but at the expense of lesser precision in the A, values which are calculated by subtracting closely spaced recorded points. The other extreme is program 4 which reverses the effect of these two considerations, and programs 2 and 3 are compromises. In the spirit of minimizing calculational effort in determining a crude k alone, method 2 is attractive because of ita one-parameter feature. It is not necessary to draw a line but simply to calculate a series of E’values corresponding to eq 15 and to select some central value from these, say, the median or the arithmetic mean. Then k is calculated from the central E ’. The following discussion applies to circumstances where both kinetic parameter values and statistical uncertainties of those values are required from the data. Toward this end it is generally necessary to have some a priori information about the measurement error or “noise” in the 4, data. For example, the error in 4, may have a uniform variance independent of the magnitude of 4, or the standard deviation of the scatter may be proportional to the magnitude of 4, or the variance may depend on the magnitude in some other functional way. This variance function var 4, is difficult to determine from the kinetic data especially if the kinetic parameters are unknown as well. If runs can be made in cases where the form of the model equation is known and all the parameters can be fixed as well, then a study of the scatter of recorded points about the known 4 vs. t curve will, in principle, lead to the variance function. More typically, however, the analyst relies on his experience with the experimental system for the nature of the statistical scatter of his data.
ANALYTICAL CHEMISTRY, VOL. 53, NO. 2, FEBRUARY 1981
210
RECORDED POINT I N D E X n PROGRAM NUMBER
2
4
6
8
I-r-T-Vy-I--T-T-
12
10
14 I
I
16
I
I
18
I
1
20
I
I
I
22 I
I
24 = N 1
J
J=6At
T=At, 7'=4At ]T=At, T'= I 2 A t 7=2At, J'=At
LLL
J=At, T'=3At ],=At, f'=12At J=2At, 7'=At
v)
7=2At,
7'=6At
J=4At,
T'=At
7=4At, 7'=12At J L
I
2
_
_
4
I
_
L
I
I
8
6
I
I
10
I
I
12
I
I
1
14
I
~
16
~
I
18
I
I
20
I
I
I
1
22
1
24=N
RECORDED POINT I N D E X n
Figure 1. A diagram showing possible programs for treating N = 24 recorded data points in a kinetlc model involving a polynomial of degree
M = 2. Each horizontal line connecting three or four tic marks represents a linear combination of those recorded data points with the tics.
For introduction of these statistical considerationsinto the model, eq 4 is rewritten as &) = @(t)+ e = y exp(-kt) PM(t)+ e (18) where 3(t)denotes the expectation of 4 as a function of time and e is a random variable with expectation value of zero representing the measurement error or noise. Each recorded value is a sample taken from I$ at time t, and this discrete function has two components, 9, which is mathematical (nonrandom) and en which is random. The variance function var 4, = var e, may be a function of 3,. The variable t, which is the time of data recording will be assumed to be a nonrandom variable since its value is usually measured with high precision or preset by a timing device. The method of least squares is the most popular procedure for finding straight-line parameters. Consider the use of least squares with a generalized Guggenheim line obtained from recorded data characterized by the variance function var 4,. According to eq 14, the line is plotted as yi = In Ami vs. xi = ti. Ordinary unweighted least squares requires that the xi be nonrandom or at least of significantly higher precision than yi and that yi variances be uniform. The requirement on xi is usually met in kinetic analysis but the yi requirement is not. An arbitrary var 4, will not propagate a uniform variance to yi. This is particularly true for uniform var 4, for which the logarithmic transformation makes var yi nonuniform as has been noted many times elsewhere. Nevertheless, a generalized Guggenheim line may be calculated by weighted least squares with weighting factors applied to each point. This technique will be discussed in the next section. Neither weighted nor unweighted least-squares methods are appropriate for method 2 or generalized KMS lines because in these cases both xi and
+
+,
+,, vertically aligned
yi quantities are linear combinations of 4, and so both are random variables with similar variances. WHEN MEASUREMENT ERRORS ARE SMALL AND UNIFORM
Frequently experimentalconditions are such that the errors are small relative to 9, and independent of a, Le., conditions of high signal-to-noise and uniform variance. In these circumstancesweighted least squares applied to Guggenheim lines and other statistical methods applied to method 2 and KMS lines become tractable and amenable to direct comparisons. Here the comparisons will he made on the basis of the relative variance estimates var k of the predicted rate constants. This is the criterion of statistical "efficiency" which is particularly useful for analysts whose prime consideratiop is to minimize the statistical uncertainty of the estimated k . The condition of small e, relative to 3, allows the use of propagation-of-variance which is an approximation technique valid in this limit. The uniform variances var 4, will be denoted in the usual way by the constant u2. If, following convention, the weighting factors for a least-squares Guggenheim line are taken as the reciprocals of var yi, propagation-of-variance leads to e,
w;l
=
= var(1n AM)i a 2 ~ / [ A ~ ( @ ) l ?
(19)
where L a2L = u 2 c y p
(20)
for L-term linear combinations. The absence of covariance terms in eq 20 stems from the further assumption that the measurement errors of the various recorded points are sta-
ANALYTICAL CHEMISTRY, VOL. 53, NO. 2, FEBRUARY 1981
tistically independent. This circumstance, which will be assumed throughout this paper, will generally be true in experiments where the spacing At is not too small. The functional notation in the denominator of eq 19 means that AM should be calculated as a linear combination of the expectations a,, rather than the actual measurements $J,, and this is because the measurements have a random component which introduces unwarranted statistical scatter into the w,. Use of 0,in place of 4,,, however, leads to the problem that 0, are not known until the straight line parameters and k are calculated but these depend on wi. Two procedures are possible and both are iterative: (1) minimize the weighted s u m of squares with respect to y and k with these parameters written explicitly into the weighting factors, a nonlinear regression problem, or (2) proceed in the usual way via normal equations but iterate with respect to the weighting factors by starting with 4,, in place of a,, and then successively refining w ias and k are calculated. Experience shows that the latter iteration tends to converge with less trouble th,an the former and so the formulas to be derived here for var k will be based on the latter. The variance es_timatefor the weighted least-squares slope is the same as for k since these parameters have the same value but opposite sign. Thus (11)
and this variance can be calculated from the weighting factors once the iteration has converged. However, for the purpose of comparison of these variance estimates with those of the other methods, we substitute for wi and x i and derive an ?xpression for a nondimensionalreduced variance V, y2 (var k)/(a2kz) in terms of k~ (or kAt) and the program integers N, L,p . An expanded version of this section is available as supplementary material on microform. This material describes this derivation, the resulting variance expression,and similar developments for method 2 and generalized KMS equations. For method 2 and generalized KMS lines, it is noted that both x ; and yi are the same linear combination of 4,, and so with uniform var &, var x i = var yi = u22pf which is a constant for any given linearization method and integer L (or M). There are many procedures in the statistics literature for estimating the parameters of a straight line where both variables are subject to uniform errors. An entrance to this literature might begin with surveys by Madansky (12) or by Kendall and Stuart (13). One such method is Wald‘s method of averages (14) which Schwartz and Gelb (15)apply to a KMS line with M = 0 and derive a formula for var k (their eq 7) under the same conditions assumed here. This derivation will now be extended to any generalized KMS line. Wald‘s method is to divide the total of I plotted points into two groups of equal or nearly equal size. In these experiments where N can usually be selected at will, it should be no trouble to choose N such that I is even and each group contains Ig = I / 2 points. If the fitted line is written y = a + bx, then Wald’s method is to set the slope b through the centroid of each group. A variaGce estimate of this slope leads to a variance estimate for k. However, as opposed to generalized Guggenheim programs all of which condense to a single variance equation, the method 2 and KMS programs do not, but require four different expression. Referring to the method 2 and generalized KMS programs illustrated in Figure 1,we will characterize some of these as “alternating” and others as “sequential” based on the order of the initial times tl of successive linear combinations x, or yi. For example, in program 5 the successive initial times correspond to the plotted coordinate in the order xyxyxy which is alternating, whereas
211
program 6 initial times are in the order xxxyyy which is Sequential. Similarly, programs 7,8,10, and 12 are alternating and programs 9 and 11are sequential. Thus an alternating program is defined as one in which x and y strictly alternate and a sequential program involves some repetitions, at least x x and yy. The var k expressions differ for alternating and sequential programs and they also differ for p = 1 and p 1 2.
The Wald method sets the slope of the line through the centroids of two groups of plotted points, and while this is appropriate for a generalized KMS line, it is not for a method 2 line which is inherently forced to pass through the origin of x , y coordinates. To treat method 2 data in a manner analogous to Wald, we will adopt an ad hoc procedure of setting the slope by passing the line through the centroid of the group consisting of all plotted points. The slope then is b = ji/f where the means j i and X include all I points. With this modification of Wald’s procedure, the reduced variance estimates V, for k calculated for method 2 can be derived. In addition to the modified Wald procedure we have derived expressions for a second ad hoc but statistically reasonable way of treating method 2 data. Equation 15 shows that the ratio y i / x i represents a value E;! The mean E’of all these Ei’Flds an wtimate k. Although both xiand yi have the same variance, the variance of the ratio is not uniform with respect to the index i. Therefore, each Ei’must be weighted in the mean so that E’= (2iwgi/xi)/(Z;wi) = exp(-k.r’) where the weighting factors are taken equal to the reciprocals of var &/xi). Under the conditions assumed here, the weighting factors turn out to be
The variance estimate for such a weighted mean is (2;wJ-l and proceeding as described above for the other methods yields reduced variances. Because E’ is the weighted mean of the slopes drawn from the origin through plotted x;,y; points, this procedure will be referred to 89 the “weighted mean slope” procedure. To study the relative statistical efficiencies predicted by these four procedures, we have evaluated the reduced variances for a variety of programs covering ranges of parameters likely to be used in practice. In all cases, we fixed to = 0 and NkAt = 3 which means that 4, data would be recorded for three firsborder lifetimes. Then we varied (kat)-‘, the number of spacings per lifetime, from 8 to about 100, varied M from 0 to 2, and varied p from 1up to a maximum of 36 where this was ppsible for a given (MN configuration. Considering the four k-calculating methods and the use of alternating vs. sequential programs for all these methods except the Guggenheim, over 300 evaluations were done and these revealed the following trends and conclusions. (1) For generalized Guggenheim lines, V, decreases monotornically with increasing p for fiied ( M N configurations. The same is true for all other methods with p 1 2. However, for many, but not all cases, V , for p = 1is less than V, for p = 2, 3, etc. Nevertheless, there always exists some high p which yields a smaller V, than p = 1 for a given (M,N). (2) In all cases when an alternating program is compared with a sequential program both for the same (M,N,p) and method, the sequential program yields a substantially smaller V,. (3) When several k-calculating methods yield programs having the same set of integers (M,N,p) the generalized KMS-Wald line with sequential program always yields the smallest V, when this is a feasible program for a particular p . The generalized Guggenheim method always yields the next higher V,. (4) When one
212
ANALYTICAL CHEMISTRY, VOL. 53, NO. 2, FEBRUARY 1981
surveys all programs considered here for a given M , the smallest V, is realized by a generalized KMS-Wald line with large p , large N , and a sequential program. Conclusion 4 sumarizes the important outcome of this Gmited study of the relative statistical efficiencies of seven k-calculating methods. For a fixed total run time, the least statistical uncertainty can be obtained by maximizing N (or by minimizing the spacing At), by maximizing p , and by using the KMS-Wald method with sequential program, which is illustrated in Figure 1 by program 13 for a (M,N) = (2,24) configuration. These conclusions are valid only for kinetic experiments with small and uniform measurement errors and for conditions such that these errors are statistically independent. This latter restriction places a lower limit on At and thus an upper limit on N. In principle, the measurement error or noise e in eq 18 has a t least some component which is a continuous function of time, and this component can be characterized by an autocovariance function which expresses the covariance between points separated in time by At. As At approaches zero, the value of this function approaches the variance of the continuous component of the noise, and this value cannot be zero. Therefore, the condition of statistical independence which implies zero covariancesbetween d,, must eventually be violated for adjacent 4, when At is sufficiently small or N is sufficiently large. While this section considered only the single criterion of statistical efficiency in selecting an optimal program, it is worthwhile to mention another criterion and that is of digital computing requirements. By this we mean the extent of storage capacity and/or the amount of real time required to compute an estimate k from dn data. This criterion is especially important in experiments where conditions dictate that large computer storage cannot be devoted to accumulating 4, values or that a k estimate must be made as early as possible in the kinetic run. By examination of all the programs in Figure 1, it can be seen that program 13, which is statistically most efficient, is one of the least efficient in terms or computing requirements. Considering that at least two ( x L , y,) points must be calculated before a first evaluation can be made of the slope and, hence, of k in program 13 at least 22 of the 24 dnpoints must be stored and 92% of the run time must pass before these data can be collected. In contrast, of some statistically inefficient programs, the alternating program 8 requires only 12 4, points for a two-point line and the alternating program 5 requires only 8, 4, points and one-third of the run time since in method 2, k can be estimated from a single (xi, yJ point. The trade-off between statistical efficiency and computing storage also enters into the decision of how large a number N should be used in any given program. ILLUSTRATIVE EXAMPLE
To help readers become familar with these linearization methods, we offer a particular example as an illustration. Kinetic data are simulated by adding random normal noise with zero mean and 0.01 standard deviation to each of 24 values of 4, = 160[exp(-0.0672tn)- 11 + 12.7t, + 0.080t: taken with time spacings of At = 2. The recorded data are shown in Table I and since these conform to a ( M , N ) = (2,24) configuration, the programs diagrammed in Figure 1 are appropriate for kinetic analyses. To conserve journal space, we detail here only 3 of the 16 calculations made and these three are shown in Table 11. The results of the remaining calculations including statistical uncertainty estimates are given as part of a microform supplement to this article. ESTIMATING KINETIC PARAMETERS OTHER THAN K
While all thTe linearization methods yield estimates of the rate constant k from the line slope, the generalized Guggen-
Table I. Simulated Recorded Dataa Used for the Illustrative Example tn 2 4 6
@n
5.5957 14.366 25.992 40.169 56.679 75.363 95.928 118.29
8
10 12 14 16 a
11
tn 18 20 22 24 26 28 30 32
On
tn
On
142.25 167.72 194.60 222.78 252.14 282.68 314.31 346.96
34 36 38 40 42 44 46 48
380.57 415.12 450.58 486.88 524.05 561.99 600.76 640.28
Simulated from the model
= 160[exp(-0.0672tn) -
+ 12.7tn + 0.080tn2 + e n , where en is a random nor-
mal variate with u = 0.01 and zero mean.
Table 11. Three Kinetic Analyses of Recorded Data Shown in Table I Generalized Guggenheim Method: Program 4 of Figure 1 (p = 6 ) plotted points weighting factors x i = tli y i = In wi,lo6 b 2 4 6
3.164 3.028 2.898 2.7 67 2.626 2.494
8
10 12
0.280 0.214 0.164 0.127 0.095 0.073
Result: of weighted 15ast-squar:s analysis: -slope = k = 0.06686, E = exp(-kT) = 0.4483 intercept = ln[$(l- 3E + 3EZ - E 3 ) ] = 3.2976 preexponential factor +? = 161 Method 2: Program 7 of Figure 1 (p = 2 , p ' = 1) initial times plotted points tli(l)
tl,(2)
2 18 34
4 20 36
xi =
a
1.729 0.5533 0.2338
y i = A2i(2) a
1.658 0.4783 0.1654
Results of modified Wald linehanalysis: $ope through the origin = exp(-k T ' ) = 0.9145 k = 0.0447
Generalized KMS Method: Program 13 of Figure 1 (p = 4 , p ' = 12) initial times plotted points tli(l)
tli(2)
2 4 6 8
26 28 30 32
x. = Ali(l)
c
34.48 31.35 28.74 26.36
yi = Ali(2)
15.06 14.44 13.91 13.48
Results ,Of Wald line analysis:
