Anal. Chem. 1982, 5 4 , 1525-1528
1525
Monte Carlo Simulation for the Study of Error Propagation in the Double Known Addition Method with Ion-§elective Electrodes C. E. Efstathiou and T. P. HadJiloannou" Laboratory of Analytical Chemistry, University of Athens, Athens, Greece
The error propagation In tlhe double known addltlon method used with Ion-selectlve electrodes Is studled with a Monte Carlo slmulatlon. Uncertalmtles associated wlth five experlmental quantltles, three potential readlngs and two concentratlon Increments, expressed In terms of standard devlatlon and relatlve standard devlatlon are propagated affectlng the reproduclblllty of the final result. Results from 500 slmulated measurements are used ito estlmate the relatlve standard devlatlon of the flnal resullt In each case under a varlety of condltlons. Graphlcal prwentatlon of results Indicates the expected uncertainty of analytlcal results obtalned by the double known addltlon method. Experlmental data obtained wlth a number of determlnstlons wlth a fluoride Ion selective electrode, uslng the double known addltlon method, reveal relatlve standard deVlatlON proflles In agreement with those predlcted with the stated technlque.
The standard known addition method is a commonly used technique with many of the existing instrumental methods of analysis. The effect of the increment size on the attained precision, for a variety of transfer functions of the sensors in use, has been thoroughly iatudied by Ratzlaff (I). The standard known addition method is widely used in potentiometry with ion-selective electrodes (2). This method is commonly in use whenewer speed and experimental convenience are required or when a reliable calibration curve cannot be obtained because1 the chemical background of the samples is either unknown, variable, or impossible to be duplicated artificially. An accurate value of the prelogarithmic factor, usually called slope, S, should be known in order to apply the standard known addition method and obtain reasonably accurate results. Generally, it is not a good practice to use the theoretical slope, provided by the Nernst equation. On the other hand, the recommended cdculation of slope from a recently obtained calibration curve is not a cure either, because the slope can be affected by the chemical background of the actual sample. An erroneous value of slope will bias all results with a systematic error. An interesting approach which eliminates the necessity of using the slope is the double known addition method (DKAM). This method ha(3been criticized and has not been widely accepted because the mathematics involved are rather complicated and time-consuming since approximation techniques should be used to calculate the final result; also, the accuracy of the method is poor being 5 to 10 times worse than that of other ion-selective ellectrode methods commonly used (3, 4). In light of the wide use of programmable calculators, the aforementioned difficulty of complicated mathematics is eliminated. A hand-held programmable calculator will give an acurate, within f0.1% or even better, result, in a few seconds using for instance the Newton-Raphson iterative
method (5). A simplified reiterative procedure for the calculations involved in the DKAM has also been proposed (6). Recently, methods employing nomograms (7) and numerical-graphical techniques (8)have been developed to facilitate calculations in the DKAM. The purpose of this work is to reexamine the accuracy of the DKAM and demonstrate a statistical procedure to obtain more quantitative information on the accuracy of the method in terms of relative standard deviation of the analytical results under a variety of conditions. Experimental data obtained with a number of determinations of fluoride in synthetic water samples with a fluoride ion selective electrode, using the DKAM, reveal relative standard deviation profiles in close agreement with those predicted with the Monte Carlo simulation. THEORY
In the double known addition method two successive known additions of solutions containing the ion being measured are made to the sample. Three separate readings of the potential are made: Eo,El, and E2. Eo is taken on the initial sample with the unknown concentration C,, El is taken after the first known addition which will result in an increase of the concentration by ACl, and E2 is taken after the second known addition which will result in a further increase of the concentration by AC,. These potentials can be expressed as follows: Eo E ' + S log C, (1) El = E ' + S log (C, AC1) (2) E2 = E'+ S log (C, + AC, ACZ) (3) Equations 1-3 are valid if the following assumptions are made: (i) The electrode slope, S, remains constant over the concentration range C, to (C, + AC1 + AC,). (ii) There is no potential drift and the term E' remains stable during all measurements. (iii) Small volumes of sufficiently concentrated solutions are added, so that the volume change on making additions can be neglected. (iv) There is a strong ionic strength background and the activity coefficients remain essentially constant. (v) No chemical interferences are encountered. The last two assumptions allow us to deal with analytical concentrations instead of activities. Combining eq 1-3, we obtain
+
+
Equation 4 cannot be solved explicitly for C,, since it cannot be rearranged to any expression of the form C, = F(Eo,El, E29 ACi, AC2). Error Propagation. All five experimental quantities, Eo, El, E2, ACl, and AC,, are associated with an intrinsic exper-
0003-2700/82/0354-1525$01.25/00 1982 American Chemical SocEety
1526
ANALYTICAL CHEMISTRY, VOL. 54, NO. 9, AUGUST 1982 0.21
U,
>
204.8
205.1
205.4
0.16
205.7
0.3t
0,3t
n
c 0 N C E N T R AT I O N S ,arbitrary units POTENTIALS,
mV
Figure 1. Computergeneratedprobability histograms of 500 simulated data required to estimate a single RSDg based on the hypothetical transfer function: E(mV) = 100 59.15 log C(C in arbitrary units). Intrinsic uncertalnties for ail potential measurements (in terms of standard deviation) are taken equal to 0.1 mV and intrinsic uncertainties for all concentration increments (in terms of relative standard devlatlon) are taken equal to 0.7%. Unknown concentration, C,, is taken equal to 10.0 and concentration increments, AC, and AC2, equal to 10.0 and 40.0. Dashed llnes represent the calculated normal (Gaussian) curves.
+
imental uncertainty or error. This uncertainty may be expressed in terms of standard deviation, uE, for all potentials and a relative standard deviation (percent), RSDAc,for all concentration increments. Nonrandom (determinate) errors on the five experimental quantities are not considered. The reasoning of using the same relative standard deviation for both concentration increments is that the main source of uncertainty in AC1 and AC2,from a practical point of view, may be confined to the uncertainty introduced by the volumetric instrument (e.g., syringe or micropipet) used to transfer volumes of equal size from standard solutions of different concentrations. The titer of these solutions can be known much more precisely and cannot be considered as another source of uncertainty. For simplicity, UE is considered constant for all potential measurements but it is reasonable to say that this uncertainty may be larger for Eo (more dilute solution) and smaller for E2(more concentrated solution). The intrinsic uncertainty of AC1 will be propagated to El and combined with its intrinsic uncertainty. Similarly, intrinsic uncertainties of both AC1 and AC2will be propagated to E2and combined with its intrinsic uncertainty. Taking the aforementioned error propagation routes in consideration along with the fact that no direct solution of eq 4 exists, then obviously the uncertainty propagated to the unknown concentration C,, expressed in terms of relative standard deviation (percent), RSDc,, cannot be calculated directly. This is a typical case where a numerical method, such as a Monte Carlo simulation, may be applied to estimate RSDc, values, given the other uncertainties, UE and RSDAc It is expected that RSDc, values will also depend on the size of the increments ACI and ACz. Monte Carlo simulations have been recommended by Schwartz (9) for solving problems of this nature frequently encountered in analytical chemistry. SIMULATED DATA Monte Carlo simulations are always carried out with the aid of a computer. According to this procedure, in the present
case, a number of simulated measurements takes place with an equal number of sets of all five quantities. All arithmetic values of each quantity should be normally distributed around given mean values with given standard deviation, cE,for Eo, relative standard deviation (percent), RSDAc,for both AC, and AC2, and standard deviations for El and Ez dependent on both Q and RSDAc values. If a sufficiently large number of simulated measurements takes place, then the relative standard deviation of all calculated values of C, will be a good estimate of RSDc,. Equation 1 is considered as the hypothetical transfer function of the ion-selective electrode in use and the following numerical values are used: E’ = 100 mV, S = 59.15 mV/ concentration decade, and C, = 10 arbitrary concentration units (any set of numerical values may be used for E’and C, without affecting the results). In calculation of the RSDc, for a given set of aE, RSDAc, AC,, and AC2the following steps are taking place: (i) Normally distributed random values (ACJN and (AC,), are generated (see Appendix A) with mean values AC1 and AC2and standard deviations RSDApAC,/100 and RSD,pAC2/100, respectively. (ii) By use of the real value of C, and the generated values (ACJN and (AC2INthepotentials Eo,El, and Ezare calculated, using eq 1-3. Eois the theoretical value but potentials El and Ez are not because they are calculated using (ACJN and (AC,), (iii) Normally distributed random values EO)^, and (E2)N are generated with mean values of the potentials Eo,E,, and E2,respectively, and standard deviation UE. (iv) Equation 4 is solved, using the Newton-Raphson iterative procedure (see Appendix B). The values (EO)N,(E~)N, and (E2)N are used in place of Eo, El, and E,. Through this step, a simulated “experimental”value ( C x ) N is computed. (v) Steps i-iv me repeated 500 times. RSDc, is calculated from all ( C , ) N values found. In Figure 1,typical probability histograms are shown of all simulated data required to estimate a single RSD, value. It can be seen from Figure 1 that the overall precision of potential values deteriorates from Eoto E,. It is also clear that (Cr)Nvalues are normally distributed around the “real” one.
ANALYTICAL CHEMISTRY, VOL. 54, NO. 9, AUGUST 1982
1527
The normality of all distributions has been verified with the x2 test. REAL DATA
A fluoride ion selective electrode is used to obtain some real data with the DKAM to compare them with those predicted by the Monte Carlo simulation. This electrode combines high potential stability and fast. response and it is less influenced by interferences and probably is the most widely used ionselective electrode. A synthetic water sample (SWS) containing 2.000 ppm F (as NaF), 50.0 ppm Ca (as CaCl,), and 10.0 ppm Mg (as MgSO,) was analyzed for its fluoride content with the DKAM after mixing with an equal volume of total ionic strength adjustment buffer (TISAB). Two 50-pL aliquots of fluoride working standard solutions with concentrations ranging from 200 to 10 000 ppm F were added in 50.00 mL of 1+1 SWS-TISAB mixture. For comparison of the FISDczsof the real data with those obtained by the Monte Carlo simulation, a good estimate of the actual RSDAcand uE i ~ necessary. i RSDAcis easily calculated by repetitive weighings (to the nearest 0.01 mg) of 50-~L water aliquots transferred by the 50-bL syringe used. Long-term potential drift prevents an analogous direct estimation of UE from potentials obtained by repetitive immersings of the fluoride-selwtive electrode in the same fluoride solution, but it may be estimated as follows: Sequential readings of potentials EA and EB are obtained respectively in 1+1SWS-TISAB mixtureri containing 1ppm and 2 ppm F. The standard deviation of itlhe difference AE (AI3 = EA - EB), urn, is equal to (UE? + u,F:)l/’ where Q, and uE, are the standard deviations of potential readings EA and EB. Assuming that uE is constant through the entire working range of fluoride concentration (tsE = UE, = Q,), then Q is equal to urn/Z1I2. urn may be calculated from all AE values which are not subject to long-term drift. EXPERIMENTAL SECTION Instrumentation. An Orion combination fluoride-selective electrode (Model 96-09) is used. The outer chamber of this electrode was filled with Orion 90-00-01equitransferent solution. Emf values were measured with an Orion Ion-Analyzer (Model 801 digital pH/pIon meter). ,411 solutions are measured at ambient temperature in 50-mL polypropylene beakers with constant magnetic stirring. Reagents. All solutions were prepared with deionized-distilled water from reagent grade materials. Sodium Fluoride. Stock; solution, 10000 ppm in fluoride. Dissolve 22.101 g of NaF (Merck, p.a.1 in water and dilute to 1 L. Necessary working fluoride solutions ranging from 200 to 5000 ppm F are prepared by appropriate dilutions of the stock fluoride solution. All fluoride solutions are kept in polyethylene bottles. Total Ionic Strength Adjustment Buffer (TISAB). Prepare this solution as described in the literature (IO). Synthetic Water Sample, .2.000 p p m F-50.0ppm Ca-10.0 p p m Mg. In a 2-L volumetric flask add with the following order: 25.00 mL of a 4000 ppm calcium solution (as GaCl2),1.8 L of water, 25.00 mL of a 800 ppm magnesium solution (as MgS04), 10.00 mL of the 400 ppm fluoride working solution, and dilute to the mark. Store in a polyethylene bottle. Measurements. Pipet 50 CH) mL of 1+1SWS-TISAB mixture into the polypropylene beaher. Immerse the fluoride-selective electrode and start the magnetic stirrer. Read the potential (Eo) after 90 h 5 s. Immediately after the potential reading transfer a 50-pL aliquot of the appropriate working fluoride solution with a 50-pL syringe (Hamilton Go., Reno, NV) and read the potential (E,) after 90 f 5 s. Repeat the same with the next appropriate working solution (E2). It is essential that all measurements are performed in a similar fashion with constant delays between potential measurements to compensate any long-term potential drift. All potential measurements are estimated to the nearest 0.05 mV by observing the relative appearance of the last digit (0.1 mV)
0
5
15
20
Flgure 2. Effect of the ratios ACIICxand ACpIAC1on the relative standard deviation of the analytical results for monovalent ions $31 = 59.15 mVlconcentration decade). Standard deviation (intrinsic)of all potential measurements, 0.1 mV, and relative standard deviation of concentratlon increments, 0.5 % .
of the Ion-Analyzer. RESULTS AND DISCUSSION RSD, values estimated by Monte Carlo simulations were reproducible within hts% when 500 simulations were performed. We found it more convenient to present results in the form of plots of RSDc, vs. AC2/AC1 for various ACl/C, values. A representative plot of this kind is shown in Figure 2.
For the plots of Figure 2 a reasonable UE value of 0.1 mV has been associated with the intrinsic uncertainty of all three potential measurements whereas an RSDAcvalue of 0.5% has been associated with AC1 and ACz values, which is a representative reproducibility with most micropipets and microsyringes used when small volumes of solutions are to be transferred. The arithmetic figures of the RSD, axis are valid for ion-selective electrodes responding to monovalent ions (IS1 = 59.15 mV/concentration decade). These figures should be doubled in case of divalent ions (IS1 = 29.58 mV/concentration decade) or generally multiplied by the factor 59.15/S, where S is the slope of the electrode actually used. As a typical example it can be seen from Figure 2 that a relative standard deviation of about 5% should be expected on the calculated C, values of a monovalent ion with the DKAM; if the standard deviation of the measured potentials is 0.1 mV, the relative standard deviation of concentration increments is 0.5%, when the first increment increases the unknown concentration of the ion by 50% (ACl/Cx = 0.5) and the second concentration increment is four times larger than the first one (AC2/AC1 = 4). Obviously, increasing ACl/C, and ACz/AC1 will provide more accurate results down to a limit value of 1.5% for RSD,. This limit is set by the uncertainty associated with AC1 and ACz values. However, practical reasons such as the requirement of constant ionic strength and the limited linear range of electrode response will prevent the unlimited increase of AC1/C, and AC2/AC1. It can be seen from Figure 2 that the previously, for simplicity, recommended value of AC2/AC1 = 1 (3) is not the proper one for minimizing the error. In Figure 3 the effect of u~ on RSDc, is shown. I t can be
1528
ANALYTICAL CHEMISTRY, VOL. 54, NO. 9, AUGUST 1982
Table I. Comparison of RSDc, Values Obtained Experimentally with Those Obtained by the Monte Carlo Simulation expt no. 1 2 3 4 5 6 7
0.200 0.200 0.200 0.500 0.500 0.500 0.500 1.00 1.00 1.00 1.00 5.00 5.00 10.0
8
Flgure 3. Effect of various standard deviations of all potential measurements on the relative standard deviation of the analytical results for monovalent ions (IS(= 59.15 mV/concentratlon decade). AC ,/CX ratio Is equal to 1 for all curves. An RSDA, equal to 0.5% is assumed for all curves except the lower one where RSDA, is taken equal to zero.
seen that RSDczvalues are practically linearly related to uE values (but the minimum value of RSDc, depends on the actual RSD,c value). This fact extends the usability of the plots in Figure 2. In Table I experimentally obtained values of RSD, are compared with those predicted by Monte Carlo simulations. CONCLUSIONS The double known addition method with ion-selective electrodes has not obtained the recognition that it certainly deserves. It should be emphasized that errors associated with this method are random in nature compared with the systematic errors which may occur with the (single) standard known addition method. That means that multiple measurements can improve the accuracy of the DKAM. The contemporary trend to incorporate calculating power in modern potentiometric analyzers could increase the accuracy of this method. More accurate measurements of potentials may be obtained and the tedious calculations involved in this method could be left to the properly programmed built-in microcomputer. APPENDIX A Generator of Normally Distributed Random Numbers. Random numbers possessing normal (Gaussian)distribution, with given mean, p, and standard deviation, u,, can be generated by using eq 5
x =p
+
--,g
C r i - (n/2) i=l
(n/ 12)l I 2
where rj is a uniformly distributed random number whose value lies between 0 and 1 provided by the high-level programming language in use (BASIC). n is frequently selected equal to 12 (11). The underlying basis for eq 5 is the central limit theorem which states that if rl, r2,..., r,, are identically distributed numbers with finite mean p and variance c2,then Ci=lnritends to be normally distributed with mean n p and variance nu2, as n (12). APPENDIX B Iterative Solution of Equation 4. Equation 4 is solved numerically using the Newton-Raphson iterative procedure (5). Starting with an initial guess value for C,, C,,o, a new one, C,,l, is calculated which is a better approximationof the actual
-
RSDc,, % concns, ppm F Monte -AC, AC, C, exptl a Carlo b
9 10 11 12 13 14
1.00 2.00 4.00 0.500 1.00 2.50 5.00 1.00 2.00 5.00 10.00 5.00 10.0 10.0
1.031 1.002 0.994 0.985 0.988 0.994 0.996 1.O2Oc 1.007 1.012c 0.997 0.988 1.006 0.985c
8.73 4.61 5.22 6.03 3.99 3.39 2.11 3.25 1.94 1.99 1.41 2.67 1.44 1.90
7.15 5.11 4.39 7.33 4.47 2.85 2.11 3.92 2.51 1.67 1.38 2.40 1.65 2.44
a Mean C, and RSDc, values from 15 determinations with the DKAM on mixtures 1t 1 SWS-TISAB containing 1.000 ppm F. RSDc, values calculated from results of 500 simulated measurements. U E and RSDACvalues used were 0.055 mV and 0.35%,respectively. These results tested with Student's t test at 95%confidence level denote significant deviation from the true value.
C, value. Generally the new value C,,++,is calculated from the previous one, C,,j Cz,i+l = Cx,i - q(Cz,J/q'(Cz,J (6) i = 0, 1, 2, where q (C,,J is the arithmetic value of the left term of eq 4 for C, = C,,t and q'(C,J is the arithmetic value of the first derivative of this term with respect to C,, for C, = CXi. The first derivative is given by q'(Cz) = (AC1 + A C J ( E 1 - Eo) - ACI(E, - Eo) (7) C, C, + AC1 + AC2 C, AC1
...
'[
-
~
_
_
I
_
_
]
_ _ _ l _ _ l
+
The iterations are continued until the following condition is met: 10Ol(C,,i+l - C,,J/C,,il 5 (8) where t is the required degree of precision (percent) of the numerical solution for C., In the present study c was set equal to 0.05. LITERATURE CITED (1) Ratzlaff, K. L. Anal. Chem. 1979, 51, 232-235. (2) Durst, R. A. In "Ion-Selectlve Electrodes"; Durst, R. A., Ed.; Natlonal Bureau of Standards: Washington, DC, 1969;Special Publlcatlon 31 4, Chapter 11, pp 381-385. (3) Orion Research Inc.,Newsletter 1970, II ( 7 4 ,34-35. (4) Orion Research Inc., Newsletter 1970, I1 (9-lo),43. (5) Dence, J. 8. "Mathematical Technlques In Chemistry"; Wlley: New York, 1975;p 54. (6) Slmpson, R. J. I n "Ion-Selectlve Electrode Methodology"; Covington, A. K., Ed.; CRC Press: Boca Raton, FL, 1979;Vol. I, Chapter 3,p 51. (7) Zhao, Z.;Zhou, X.; Sun, T. Hua Hsueh Tung Pa0 1981,(3),21-22, 58. (8) Kaplan, B. Ya.; Varavko, T. N. Zh. Anal. Khim. 1981, 36 (3),
591-592. (9) Schwartz, L. M. Anal. Chem. 1975, 4 7 , 963-964. (10) "Manual of Methods for Chemlcal Analysis of Water and Wastes"; U.S. Envlronmental Protectlon Agency: Washington, DC, 1974;p 66. (11) Gottfried. B. S. "Programming wlth Basic", Schaum's Outline Series; McGraw-HIII: New York, 1975b p 155. (12) Green, J. R.; Margerlson, D. Statistical Treatment of Experimental Data"; Elsevler: Amsterdam, 1977;p 55.
RECEIVED for review October 29,1981. Accepted April 5,1982. This work was supported in part by the Greek National Institue of Research.
