Influence of the Experimental Errors and Their Propagation on the

ARTICLE pubs.acs.org/IECR

Influence of the Experimental Errors and Their Propagation on the Accuracy of Identified Kinetics Parameters: Oxygen and Temperature Effects on Ascorbic Acid Oxidation during Storage Caroline Penicaud,† Philippe Bohuon,‡ Stephane Peyron,† Nathalie Gontard,† and Valerie Guillard*,† †

UMR 1208 IATE Agropolymers Engineering and Emerging Technologies, Universite Montpellier 2, CIRAD, INRA, Montpellier SupAgro, CC 023 Place Eugene Bataillon F-34095 Montpellier Cedex 5, France ‡ UMR Qualisud, Food Process Engineering Research Unit, Montpellier SupAgro, CIRAD, Universite Montpellier 2, 1101 avenue Agropolis, CS 24501, F-34093 Montpellier Cedex 5, France ABSTRACT: The oxidation of ascorbic acid was chosen as a case study to investigate the effect of experimental errors and their propagation on the identification of the parameters of a kinetics law involving both ascorbic acid and oxygen content, on a typical food storage temperature range (833 °C) and under conditions representing anaerobic storage (0% O2) to storage under air (21% O2). A nonlinear method was used to identify kinetic parameters, both by incremental identifications of intermediate parameters and by direct identification of relevant parameters. Confidence intervals on all these identified parameters were evaluated using Monte Carlo simulations. The propagation error during the incremental stages of parameters identification induced a high uncertainty on the values obtained, thus a direct identification was preferred. The accuracy of the identified parameters was also shown to be strongly dependent on the experimental procedure carried out to obtain experimental data, but the time required to perform an ideal procedure is unrealizable. A methodology based on Monte Carlo simulations was proposed to optimize a priori the experimental procedure and data treatment to improve the characterization of such a kinetics law.

1. INTRODUCTION The determination of kinetics parameters is crucial to investigate chemical reactions. Such determinations are based on a comparison between the experimental data and the mathematical models (kinetic laws) that represent the reaction. The resulting values of the reaction rate constant are dependent on the accuracy of the experimental data, the pertinence of the structure of the mathematical model chosen, and the identification procedure employed. The experimental values necessarily go with errors, resulting in uncertainties on the identified reaction rate constants. The identification procedure can also introduce uncertainty on identified values of constant rates. Although these statements are generally taken into account in chemistry,14 they are rarely considered in food science and the influence of the experimental and numerical errors on the accuracy of the identified parameters is scarcely evaluated. Moreover, the choice of the mathematical model to be used can be controversial.5,6 In the case of ascorbic acid oxidation, the reaction is generally recognized to follow a first-order kinetics law, with respect to ascorbic acid,7 which assumes that the oxygen content in the medium does not impact the kinetics rate. This is surprising if we consider that some of the authors using this kinetics law simultaneously observed that the rate of oxidation was dependent on oxygen’s availability.810 On the other hand, a few authors considered the reaction rate as a first-order kinetics, with respect only to oxygen, considering that the ascorbic acid content in the medium does not affect the reaction rate constant.1113 Despite the wide number of investigations carried out on ascorbic acid oxidation, only three studies have pointed out a joined effect of ascorbic acid and oxygen contents in food, which is the case that seems the most likely. Singh et al.14 studied the kinetics r 2011 American Chemical Society

of ascorbic acid oxidation in an infant formula during storage at 7.2 °C. Below the oxygen content limit of 8.71 mg L1 in the product, the reaction rate was assumed to follow a second-order kinetics (order one with respect to ascorbic acid and order one with respect to oxygen). Eison-Perchonok and Downes15 also described the kinetics of ascorbic acid oxidation in a buffered model system by a second-order kinetics (order one with respect to ascorbic acid and order one with respect to oxygen) on the temperature range of 3355 °C. Patkai et al.16 analyzed the decrease of vitamin C content of model solutions simulating real citrus juices, on the temperature range of 80100 °C. The rate of vitamin C consumption was considered to be first order, with respect to vitamin C, and the kinetics constant rate was expressed as a polynomial equation involving the initial oxygen content, the temperature, and the duration of heat treatment. In these previous studies, the typical law used for modeling the reaction rate constant of ascorbic acid degradation is r = k[AA]α[O2]βd , where k is the reaction rate constant and α and β are the orders of the reaction, with respect to ascorbic acid and oxygen, respectively. Three input parameters are required for modeling ascorbic acid oxidation, according to this law: k, α, and β. They are deduced from experimental measurements and are consequently not true values but rather values going with uncertainties that must be taken into account in a prediction step. In practice, that means that model input variables and parameters can be represented by probabilistic distributions due to observed variability. Received: May 20, 2011 Accepted: December 14, 2011 Revised: December 8, 2011 Published: December 14, 2011 1131

dx.doi.org/10.1021/ie201087h | Ind. Eng. Chem. Res. 2012, 51, 1131–1142

Industrial & Engineering Chemistry Research As previously illustrated in a complex system (i.e., cheese ripening), the integration and the mathematically process of uncertainty permitted to provide reliable upper and lower probability bounds for the prediction.17 Unfortunately, this statistical treatment is scarcely used in food science and engineering and was never applied to the prediction of oxidation kinetic law. This is an important methodological issue since oxidation reactions are major causes of food degradation, especially during storage. Therefore, the objective of this paper was to investigate the effect of experimental errors and their propagation on the uncertainty of identified parameters in the case of oxidation kinetic law. Oxidation of ascorbic acid was chosen as a case study and the most complete equation was used for modeling the reaction rate: r = k[AA]α[O2]βd . A typical food storage temperature range (833 °C) and atmospheric conditions from anaerobic (0% O2) to aerobic storage (21% O2) were investigated. The partial orders of reaction were experimentally determined over the entire temperature range, which was never done previously in the case of ascorbic acid degradation. Confidence intervals on identified kinetics parameters were evaluated by means of Monte Carlo simulations to take into account both experimental error and propagation error due to the identification procedure. All these results lay the foundations to propose a methodology to improve both experimental and modeling procedures in order to enhance accuracy on identified parameters of such oxidation reaction with the lowest experimental cost.

2. MATERIALS AND METHODS 2.1. Materials. 2.1.1. Reagents. Ultrapure water (conductivity = 18.2 MΩ cm) was used in all experiments. Sodium azide (purity of g99.0%), metaphosphoric acid analytical (33.5%36.5%), ferric chloride hexahydrate (98%102%), and ascorbic acid (purity of g99.0%) were supplied by Fluka. Hydrochloric acid (5 N) was purchased from Prolabo. Glacial acetic acid (99%100%), absolute ethanol (analytic grade), orthophenanthroline (99%), and ammonium acetate (purity g98%) were obtained from SigmaAldrich. 2.1.2. Controlled Environmental System. In order to control the oxygen content, experiments were carried out into a glovebox (ATLabo, Montpellier, France). The oxygen content was manually adjusted and controlled with a luminescence dissolvedoxygen probe (LDO probe, Hach Lange, Noisy Le Grand, France). A fixed oxygen concentration in the glovebox atmosphere was kept with a precision of (1% (v/v). Temperature control was carried out at (2.0 °C, using water circulation and a fan to homogenize the temperature in the glovebox. The glovebox was darkened in order to protect samples from exposure to light. 2.2. Methods. 2.2.1. Ascorbic Acid Experimental Assay. Ascorbic acid content was determined as described by Penicaud et al.18 Briefly, samples were diluted six times (v/v) in an extraction solution and stored at 4 °C until assay. The extraction solution was an aqueous metaphosphoric acid 3% (w/v)acetic acid 8% (v/v) solution, as recommended in the official method.19 The assay, carried out on microplates, was based on the redox reaction that takes place between ascorbic acid and Fe(III), yielding dehydroascorbic acid and Fe(II). Fe(II) reacts with 1,10-phenanthroline, originating the reddish orange Fe(phen)32+ complex (ferroin). This complex was spectrophotometrically monitored at 505 nm, using a

ARTICLE

microplate reader (Multiskan Spectrum, Thermo Electron Corporation, Saint-Herblain, France). 2.2.2. Experimental Acquisition of Ascorbic Acid Oxidation Kinetics. Fifteen (15) plastic jars initially containing ultrapure water (100.0 g) and sodium azide (0.2 g) in order to avoid microbial growth were subjected to, and maintained under, magnetic stirring (400 rpm) during the experiment to the required oxygen concentration, which was measured with an LDO probe (Hach Lange, Noisy Le Grand, France). The temperature was also measured by the same probe. The pH of the solutions was 5.1. For each reaction time, the ascorbic acid (0.1 g precisely weighted) was introduced in a jar; this was done in duplicate or triplicate (i.e., for each reaction time, two or three jars were used). That means that seven different reaction times were analyzed to obtain one kinetics curve. One experiment (at 20.0 °C under 21% oxygen) was done with half the initial content of ascorbic acid, to confirm that the kinetics constant rate was independent of the initial ascorbic acid concentration. The time at which the ascorbic acid was added to the water represented time zero for the reaction. After the ascorbic acid had been allowed to react for the desired time, 5 mL were extracted and placed into 25 mL of extraction solution and two or three repetitions assays were made on each extraction solution sample. The reaction jar was weighed before and after the reaction, which allowed assay values correction due to evaporation. Finally, a point on a kinetics curve represents the average of at least four measurements (two samples with two assays for each sample) to, at most, nine measurements (three samples with three assays for each sample). Experiments were carried out at 20.0 °C ( 2.0 °C under 2.5%, 5%, 7.5%, 10%, 15%, and 21% oxygen levels ((1%) in the glovebox, at 8.0, 13.7, and 27.4 °C (( 2.0 °C) under 10%, 15% and 21% oxygen levels ((1%), and at 33.0 °C ( 2.0 °C under 5%, 10%, and 21% oxygen levels ((1%). The total duration of one experiment (i.e., seizure of one kinetics curve) was one week, except for measurements at 8.0 °C, which required two weeks. 2.2.3. Ascorbic Acid Oxidation Kinetics Model. The rate law for the dependence of concentration (expressed in mol L1) on time is assumed to be d½AA ¼ k½AAα ½O2 βd ð1Þ dt where k represents the reaction rate constant (molβ Lβ s1), and α and β are the partial orders of reaction with respect to ascorbic acid and oxygen, respectively. In the case where oxygen is present in nonlimiting amounts (i.e., the oxygen content is a constant), eq 1 becomes

d½AA ¼ kapp ½AAα dt

ð2Þ

with kapp ¼ k½O2 βd

ð3Þ

The linear form of eq 3, with a logarithmic transformation, gives ln kapp ¼ β ln½O2 d þ ln k

ð4Þ

As it is generally described in literature in the case of aerobic oxidation, it is assumed that ascorbic acid follows first-order 1132

dx.doi.org/10.1021/ie201087h |Ind. Eng. Chem. Res. 2012, 51, 1131–1142

Industrial & Engineering Chemistry Research

ARTICLE

Figure 1. Two strategies for the identification of kinetics parameters: incremental identification (strategy a) and simultaneous identification (strategy b). Each modeling stage is presented, together with the corresponding equation.

degradation kinetics,7 so that α = 1 and eq 2 can be transformed as

d½AA ¼ kapp ½AA dt

ð5Þ

And, by integration with respect to time, ½AAt ¼ ½AAt ¼ 0 expð kapp tÞ

ð6Þ

Equation 6 can be linearized as ln ½AAt ¼ ln½AAt ¼ 0 kapp t

ð7Þ

The rate constant k varied with the system’s absolute temperature, T (K), according to the Arrhenius law, as follows: Ea 1 1 ð8Þ k ¼ kref exp R T Tref where kref is the rate constant (molβ Lβ s1) at the reference temperature Tref (K), Ea is the apparent activation energy ( J mol1), and R is the gas constant (R = 8.314 J mol1 K1). The linear form of eq 8, with a logarithmic transformation, gives Ea 1 1 ln k ¼ ln kref ð9Þ R T Tref 2.2.4. Identification of Kinetics Parameters. The kinetics parameters were identified by comparison between the experimental data measured as described in section 2.2.2 and theoretical data calculated by solving adequate model equation from those presented in section 2.2.3. The strategies which may be used to identify the kinetics parameters are summarized in Figure 1. The first strategy was to identify incrementally the different parameters (Figure 1a), the second one was to identify simultaneously all the parameters (Figure 1b).

2.2.4.1. Incremental Identifications (Strategy a in Figure 1). Stage 1. First, the apparent kinetics constant rate kapp was determined from eq 6 for each experiment carried out under specific conditions of temperature and oxygen level. The kapp value was iteratively adjusted to the goodness-of-fit merit function with the GaussNewton minimization procedure,20 using the Matlab software (The Mathworks Inc., Natick, MA, USA). This merit function was the root-mean-square error (RMSE) between the ascorbic acid contents (mol L1) actually measured (Yexp) and ascorbic acid contents (mol L1) calculated from eq 6 (Ypred): sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 n RMSE ¼ ð10Þ ðYpred Yexp Þ2 n 1 i¼1

∑

where n is the number of data. Stage 2. Second, the objective was to characterize parameters of eq 3. For doing this, eq 3 was reparameterized to introduce a reference kinetics constant rate. The basis for this arises from the application of eq 3 at two dissolved oxygen levels, [O2]d,1 and [O2]d,2: k1app ¼ k½O2 βd, 1

ð11Þ

k2app ¼ k½O2 βd, 2

ð12Þ

By dividing eq 11 by eq 12, and assuming that [O2]d,2 is the reference chosen, eq 13 is obtained: !β ½O2 d ref ð13Þ kapp ¼ kapp ½O2 ref d where kapp is the apparent kinetics constant rate (s1) identified in the previous step, [O2]d the measured dissolved oxygen content into the medium (mol L1), and kref app the rate constant 1133



ARTICLE

(s1) at the reference dissolved oxygen concentration [O2]ref d (mol L1). This dissolved-oxygen content of the reference 1 (i.e., approximately the [O2]ref d was chosen to be 6 mg L average of all oxygen levels investigated). From eq 13, the kref app and β values were simultaneously adjusted as described above, by iterations to goodness-of-fit merit function (eq 10) with Ypred the predicted apparent constant rate (kapp,pred) and Yexp the previously identified apparent constant rate (kapp,id). Stage 3. Finally, the third and last step was to assess the influence of temperature on kinetics parameters. The Arrhenius law (eq 8) was combined with the expression of k deduced from eq 3 in the case of the reference dissolved-oxygen content: k¼

kref app β ð½O2 ref d Þ

ð14Þ

which finally resulted in eq 15, assuming that β remains constant with temperature: Ea 1 1 ref β ref kapp ¼ ð½O2 d Þ kref exp ð15Þ R T Tref where β is the partial order of reaction with respect to oxygen, R the gas constant (R = 8.314 J mol1 K1), and T the absolute temperature (K). The reference temperature Tref was chosen to be in the middle of the studied temperature range, at 293.15 K (20.0 °C). From eq 15, the reference kinetics constant rate kref and activation energy Ea values were simultaneously adjusted as described above by iterations to the goodness-of-fit merit function (eq 10) with Ypred the predicted apparent reference constant rate (kref app,pred) and Yexp the previously identified apparent constant rate (kref app,id). 2.2.4.2. Simultaneous Identification (Strategy b in Figure 1). An alternative to the incremental identifications was to simultaneously identify the parameters values of the partial order of reaction β, the reference rate constant kref and the activation energy Ea (Figure 1b). The complete rate law for concentrations dependence on time was obtained by combination of eqs 1 and 8: d½AA Ea 1 1 ¼ kref exp ð16Þ ½AAα ½O2 βd dt R T Tref Assuming that ascorbic acid follows first-order degradation kinetics (i.e., α = 1), that oxygen concentration was constant during each experiment (experimentally measured) and that order of reaction with respect to oxygen, β, was constant whatever the temperature (this will be discussed in the Results section), eq 16 could be integrated with respect to time: Ea 1 1 ½AAt ¼ ½AAt ¼ 0 exp kref exp ½O2 βd t R T Tref ð17Þ The values of the three parameters—β, kref, and Ea—were iteratively adjusted to the goodness-of-fit merit function with the LevenbergMarcquardt minimization procedure,20 using the Matlab software (The Mathworks, Inc., Natick, MA, USA). The merit function was also eq 10 with Yexp representing the ascorbic acid contents (mol L1) actually measured and Ypred the ascorbic acid contents (mol L1) calculated from eq 17. 2.2.5. Evaluation of the Uncertainty on the Identified Kinetics Parameters. The errors on identified parameters (kapp, kref app, β, kref, and Ea) were determined via Monte Carlo simulations.21 The basic purpose is the sampling in silico of large numbers of

data sets to simulate the duplication of experimental data. Kinetics parameters can be identified from these data sets and thus large numbers of values can be obtained. Statistical treatment of these identified values of kinetics parameters allows the evaluation of the uncertainty on the parameters. 2.2.5.1. Experimental Data Generated In Silico. In silico ~ exp) was generated by superposition of a experimental data (X pseudo-random noise on the experimental data really measured (Xexp). The noise reflected the uncertainty (Ux) on the experimental data really measured: ~ exp ¼ Xexp þ Ux δ X

ð18Þ

~ exp being the experimental data generated in silico (ascorbic with X acid concentrations, dissolved oxygen content, temperature, or parameters kapp,id, kref app,id and β), Xexp is the experimental value really measured (ascorbic acid concentrations, dissolved oxygen content and temperature) or the identified value (parameters kapp,id, kref app,id, and β), Ux is the uncertainty on Xexp (see below), δ is a random number whose elements are normally distributed with mean 0 and variance 1. These normally distributed random errors were considered to be independent (in which case their covariance is zero). As shown in Figure 1, in the first stage of strategy a (identification of kapp), random scatter drawn was added on ascorbic acid concentrations [AA]t=0 and [AA]t. In the second stage (identification of kref app and β), random scatter drawn was added on kapp,id and the dissolved oxygen content [O2]d. In the third stage (identification of kref and Ea) random scatter drawn was added on kref app,id, β, and the temperature T. In the case of strategy b, random scatter drawn was added on ascorbic acid concentrations [AA]t=0, [AA]t, the dissolved-oxygen content [O2]d and the temperature T. That means that, for each identification, whatever the strategy or the stage, each set of experimental data generated in silico included joined effects of all adequate errors. These in silico experimental data were generated using the Matlab software (The Mathworks, Inc., Natick, MA, USA). 2.2.5.2. Errors Ux Considered in Simulations. Experimental Parameters. Initial ascorbic acid content [AA]t=0 (mol L1) could be calculated from weightings of introduced ascorbic acid and water: mAA Fwater 103 ð19Þ ½AAt ¼ 0 ¼ mwater þ mAA MAA where mAA and mwater are the weights (kg) of ascorbic acid and water, respectively, Fwater is the water volumic weight (kg m3) and MAA is the molecular weight of ascorbic acid (kg mol1). Thus, assuming that (i) uncertainties in the values of Fwater and MAA are null, and (ii) the uncertainties are additive, the uncertainty u[AA]t=0 on initial ascorbic acid content could be deduced: !2 u½AAt ¼ 0 u2 þ u2mAA umAA 2 ¼ þ mwater ð20Þ ½AAt ¼ 0 mAA ðmwater þ mAA Þ2 where uα is the relative uncertainty on the experimental quantity α. Uncertainties on weightings were known from calibration certificates of the balances used, the relative uncertainty in this case was thus equal to 0.04%. Assuming that experimental measurements follow a normal law, what is commonly admitted, the distribution of the values is centered on the mean (x) and the 1134



ARTICLE

variability within the population can be expressed by the relative uncertainty (ux), equivalent to the standard deviation.22,23 A level of confidence of 95% was chosen to quantify the experimental variability; thus, the relative uncertainty was multiplied by a covering factor of 1.96 to obtain the relative error.22 In the case of initial ascorbic acid content, relative error at 95% confidence level was rounded up 0.1% and thus enlarged uncertainty U[AA]t=0 was defined as U½AAt ¼ 0 ¼ 0:1% ½AAt ¼ 0

!2

¼

u½AAassay ½AAassay

!2

þ

4um ∑m

2 þ

uVpipette Vpipette

!2

ð22Þ

The relative uncertainty on the assay was known to be 3%,18 and uncertainties on the balances and pipet were known from calibration certificates, so the relative uncertainty on ascorbic acid content was calculated to be 3%. This was enlarged by a factor of 1.96 to obtain an error at 95% confidence, so the relative error on ascorbic acid content was set at 6%: U½AAt ¼ 6% ½AAt

MSE ¼

ð23Þ

The dissolved-oxygen content and temperature were measured as a function of time during the experiments; therefore, errors on both parameters were taken as extreme variations observed: U½O2 d ¼ 0:4 mg L1

ð24Þ

UT ¼ 2:0°C

ð25Þ

Identified Parameters. The determinations of the parameters kapp,id, kref app,id, and β were joined to a statistical evaluation of their standard deviation via the Monte Carlo method; therefore, errors on these identified values to consider in the further stages of the identification process were taken as their standard deviation multiplied by a factor 1.96 to obtain confidence intervals at 95% confidence levels. The peculiar case of β will be discussed in the Results section. 2.2.5.3. Identification of Parameters. For a given operating condition, a set of experimental data was generated in silico, as described in Section 2.2.5.1. Each set of in silico experimental data was used to identify the kinetics parameters by comparison with the theoretical data calculated from the adequate equation (eq 6, 13, 15 or 17, as shown in Figure 1). The equation used to calculate the theoretical data was expressed as Y = f(X) and errors could potentially be allocated both to X and Y; thus, both errors on X and Y must be taken into account. This was not the case in stage 1 of the incremental approach (strategy a) since the variable X was the time and no error was considered on it. Consequently, the merit function was eq 10 with Yexp being the ascorbic acid contents generated in silico, as described in Section 2.2.5.1, and Ypred represents the ascorbic acid contents calculated from eq 6. In contrast, in Stages 2 and 3 of the incremental approach (strategy a), errors could be

1 n 2 L n 1 i¼1 i

∑

ð26Þ

with ~ exp , i Þ2 L2i ¼ ðYpred, i Y~ exp , i Þ2 þ ðXexp , i X

ð21Þ

Ascorbic acid content as a function of time, [AA]t, was obtained by assay corrected by a factor which included dilutions carried out to perform assay and water evaporation during measurement. Thus, [AA]t resulted from one assay, four weighings (on the same balance) and one pipeting, its relative uncertainty could be expressed as u½AAt ½AAt

allocated both to X and Y; thus, merit functions used in the minimization procedure were related both to X and Y as follows:

where n is the number of data, Ypred,i are the theoretical data ~ exp are the associated to Xexp,i from eq 6, 13, or 15, and Y~ exp and X experimental data generated in silico from the Monte Carlo method. In practical terms, in Stage 2, Xexp,i represented the actually measured [O2]d, Ypred,i the kapp corresponding to ~ exp were the experimental [O2]d the actually measured [O2]d, X generated in silico, as described in section 2.2.5.1, and Y~ exp represented the identified kapp values generated in silico, as described in section 2.2.5.1. In Stage 3, the variables were T ~ exp) and kref ~ (Xexp,i and X app (Ypred,i and Y exp). In the simultaneous identification (strategy b), it was the same case as Stage 1 of the incremental approach (strategy a), since the variable X was the time and no error was considered on it. Consequently, the merit function was eq 10 with Yexp being the ascorbic acid contents generated in silico (Y~ exp), as described in section 2.2.5.1, and Ypred represents the ascorbic acid contents calculated from eq 17. For each operating condition and whatever the strategy used for the identification, the generation of in silico data followed by the identification of kinetics parameters was repeated m = 2000 times, resulting in m values of the kinetics parameters. 2.2.5.4. Statistical Treatment of Identified Parameters. It was checked that the m identified parameters values were normally distributed using a Lilliefors test carried out with the Matlab software (The Mathworks, Inc., Natick, MA, USA) at a 95% confidence level. If the test confirmed the normality, the results were expressed as the mean value of the 2000 identified parameters associated with the corresponding standard deviation. All these standard deviations had to be enlarged by a factor 1.96 to obtain a 95% confidence interval for parameters. If testing did not confirm the normality, the confidence interval at 95% level was obtained by applying its general definition (by repetition in 95% of the cases, the obtained value will be in the calculated interval; in 5% of the cases, however, it does not), i.e., by eliminating 5% of the m values (the 2.5% lowest and the 2.5% highest) to keep only 95% of the values, which thus define the limits of the confidence interval at the 95% confidence level.

3. RESULTS AND DISCUSSION The measurement of ascorbic acid degradation kinetic was done in stirred water in order to prevent any limitation due to the diffusion step. Dissolved oxygen concentration was measured and proved to be constant during each experiment ((0.4 mg L1 on the entire experiment, due to temperature variations). At first, the objective was to evaluate the reaction orders α and β experimentally, as well as the influence of the temperature on them. This could be achieved only by using the incremental approach (strategy a in Figure 1). Indeed, this approach permits to assess the influence of the temperature on each parameter separately, which is not the case with the simultaneous identification (strategy b in Figure 1). 1135



ARTICLE

Figure 2. Linear representations of ascorbic acid oxidation kinetics in water at 20 °C under different oxygen levels: (a) evolution of ascorbic acid content, as a function of time for different surrounding oxygen atmospheres, and (b) apparent kinetics constant rates as a function of dissolved oxygen in the media. On both graphs, data points represent experimental data, error bars represent confidence intervals at the 95% confidence level, and straight lines represent linear regressions from data.

3.1. Identification of Parameters Using an Incremental Approach. 3.1.1. Stage 1. The kinetics order, with respect to

ascorbic acid, was supposed to be equal to 1, as referenced in the literature.7 To validate this hypothesis, plots of ln([AA]t/ [AA]t=0) versus time t (eq 7) were drawn for each experiment (under specific conditions of temperature and oxygen level), as shown, for example, in Figure 2a for experiments carried out at 20.0 °C under different oxygen contents in the atmosphere (from 2.5% to 21%). The straight lines obtained with corresponding determination coefficients (R2) >0.95 in all cases confirmed that the assumption of a reaction order, with respect to ascorbic acid, equal to 1 can be considered as being correct. This statement was systematically assessed for each level of oxygen and temperature tested (other data not shown). In a first stage (Figure 1a), apparent constant kinetics rates (kapp) of ascorbic acid oxidation in stirred aqueous systems under different oxygen contents in the atmosphere were identified from experimental data by using kinetics law in its exponential form (eq 6). This form was preferred to the linearized one, because the transformation from the exponential form to the linearized form changed the data but also the error structure related to the data. This may lead to a bias in estimations, because then some assumptions that underlie the regression (such as homoscedaticity, independence, and normality of the residuals) are violated.24 Values obtained from the nonlinear regression are reported in Table 1. Corresponding confidence intervals at the 95% confidence level were estimated from the Monte Carlo method. Errors on experimental measurements, i.e., [AA]t and [AA]t=0, were taken into account to generate random experimental data from which the kinetics rate kapp was identified. In this way, 2000 kapp values were produced and proven to be normally distributed (confirmed via a Lilliefors test), the average of the normal distribution was kapp (the initially identified value), and the confidence interval at the 95% confidence level from this normal distribution was the standard deviation on kapp multiplied by a factor of 1.96. The identified apparent rate constants were increased with higher oxygen contents at each temperature, except for low temperatures (8.0 ( 2.0 °C and 13.7 ( 2.0 °C), where change in kapp values was not so significant, probably because kinetics are so slowed at these temperatures that variations in oxygen level do

Table 1. Apparent Rate Constant (kapp) for the Fit of [ AA]t = [ AA]t=0 exp(kapp t) (eq 6) at Several Temperatures T and Oxygen Contents [O2]d, Together with the Root Mean Square Error between Experimental and Predicted [AA]ta temperature,

[O2]d

kapp

root mean square error,

T (°C)

(mg L1)

( 107 s 1)

RMSE ( 104 mol L1)

8.0 ( 2.0

5.9 ( 0.4 8.2 ( 0.4

0.90 ( 0.44 2.26 ( 0.46

1.52 1.85

11.6 ( 0.4

1.64 ( 0.44

1.88

13.7 ( 2.0

20.0 ( 2.0

27.4 ( 2.0

33.0 ( 2.0

5.7 ( 0.4

5.36 ( 1.14

1.80

7.8 ( 0.4

3.68 ( 1.24

0.96

9.7 ( 0.4

8.18 ( 1.40

2.33

1.1 ( 0.4 2.1 ( 0.4

1.90 ( 1.68 2.16 ( 1.14

1.28 0.88

3.5 ( 0.4

4.44 ( 1.60

1.14

4.4 ( 0.4

6.41 ( 1.32

1.60

6.3 ( 0.4

8.82 ( 1.10

1.86

8.5 ( 0.4

13.9 ( 1.0

3.26

4.0 ( 0.4

5.71 ( 0.94

1.70

5.8 ( 0.4 7.6 ( 0.4

9.07 ( 1.04 14.9 ( 0.98

2.39 2.39

2.0 ( 0.4

9.99 ( 1.24

1.86

3.7 ( 0.4

17.8 ( 1.10

4.98

7.5 ( 0.4

36.9 ( 2.4

6.50

a

Data are presented as the mean values with their errors at the 95% confidence level (experimental variations for T and [O2]d, Monte Carlo simulations for kapp).

not limit the apparent kinetics rate. Singh et al.14 also stated that, at 7.2 °C, the apparent ascorbic acid oxidation kinetics rate constants in infant formula under dark conditions for different initial oxygen concentrations remained unchanged. Nevertheless, according to results of the present study listed in Table 1, 1136



ARTICLE

Table 2. First-Order Kinetics Constant Rates (kapp) of Ascorbic Acid Oxidation Found in the Literature with Their Errors at 95% Confidence Level (When Determined)a for Different Food Beverages under Air and over a Temperature Range of 237 °C food beverage

a

kapp (s1)

temperature, T (°C)

7

reference 8

orange-carrot pasteurized juice

2

5.4 10

fresh blood orange juice high-pressured orange juice

3 5

4.6 108 1.4 107

27 28

thermally pasteurized orange juice

5

2.2 107

28

high-pressured orange juice

5

6.7 108 ( 3.3 109

29


5

1.3 107 ( 2.4 108

29

( 3.9 10

26


10

2.2 107

28


10

2.7 107

28

high-pressured orange juice thermally pasteurized orange juice

10 10

9.9 108 ( 7.1 109 2.0 107 ( 2.5 108

29 28

orange-carrot pasteurized juice

10

1.1 106 ( 2.2 107

26

fresh blood orange juice

10.6

7.8 108

27


15

1.2 107

27


15

4.6 107

28


15

5.3 107

28

high-pressured orange juice thermally pasteurized orange juice

15 15

2.8 107 ( 2.7 108 4.0 107 ( 2.2 108

28 28


20

1.3 107

27

orange juice

28

4.6 108 ( 3.5 108

25

lemon juice

28

1.1 107 ( 1.9 108

25

grapefruit juice

28

7.8 108 ( 2.2 108

25

tangerine juice

28

7.6 108 ( 1.1 108

25


30

7.9 107 ( 1.3 107

29


30

8.3 107 ( 1.5 107

29

orange juice

37

3.1 107 ( 4.3 108

25

lemon juice

37

3.0 107 ( 1.3 108

25

grapefruit juice

37

2.7 107 ( 6.5 108

25

tangerine juice

37

2.8 107 ( 3.5 108

25

The confidence interval at 95% level was calculated from the standard deviation presented by the authors multiplied by a factor of 1.96.

it was hard to establish clear correlations between the kapp results obtained at a same dissolved oxygen level but for different temperatures, since oxygen solubility directly depends on temperature and thus the same oxygen level in the glovebox induced different oxygen concentrations into the reacting media. Values obtained under about 21% O2 (corresponding to the maximum dissolved oxygen level at each temperature), i.e., under air, were varying between 1.64 107 ( 0.44 107 s1 at 8.0 ( 2.0 °C and 36.9 107 ( 2.4 107 s1 at 33.0 ( 2.0 °C and could be compared to the ones found in the literature (Table 2), since kinetics constant rates are generally determined under air. It can be observed that values are in the same order of magnitude than the ones reported from other authors for degradation of ascorbic acid, even if slightly higher. For instance, at 13.7 ( 2.0 °C, a value of 8.18 107 ( 1.40 107 s1 was determined whereas an average value (mean of five values) of 3.6 107 ( 1.6 107 s1 were reported in the literature at 15 °C (Table 2). The same statement can be done at each temperature. This slight difference may be due to the influence of oxygen, since care was taken in this

study to supply oxygen continuously in the medium by permanent stirring in order to avoid perturbation due to oxygen diffusion in the medium. Note that this precaution was not systematically taken in the previous studies listed in Table 2. Subsequently, fluctuations in oxygen content in the reacting media and especially rate-limiting oxygen diffusion could have slightly reduced the kinetics rates. After examination of Table 2, it should be noticed that, in most of the studies available in the literature, the uncertainty on the kapp value was not systematically given and even determined. However, this value is crucial for using kinetic model in decision-making system. 3.1.2. Stage 2. In a second stage, the apparent rate constant (kapp) of ascorbic acid oxidation were plotted against dissolved oxygen concentration expressed on logarithmic scales (see eq 4) for each temperature. Figure 2b presents an example of the result obtained for a temperature of 20.0 ( 2.0 °C. A straight line could be observed, whose slope was β, representing the reaction order with respect to oxygen. Intersections between the straight line and the y-axis enabled to calculate the kinetics constant rate k. 1137



ARTICLE

Table 3. Rate Constant kref app and Reaction Order with Respect to Dissolved Oxygen β at Several Temperatures T ref β for the Fit of kapp = kref app([O2]d/[O2]d ) (eq 13), Together with the Root Mean Square Error (RMSE) between Experimental and Predicted kapp Valuesa kref app T (°C)

( 10

7

root mean square error,

molβ Lβ s1)

8.0 ( 2.0 13.7 ( 2.0

1.3 ( 0.6 4.4 ( 2.4

β

RMSE ( 108 s1)

0.6 ( 0.8 1.0 ( 0.8

4.88 14.9

20.0 ( 2.0

8.6 ( 1.4

1.2 ( 0.4

4.65

27.4 ( 2.0

10 ( 2

1.6 ( 0.8

3.89

33.0 ( 2.0

29 ( 4

1.0 ( 0.4

2.24

a

Data are presented as the mean values with their errors at the 95% confidence level (experimental variations for T and Monte Carlo simulations for kref app and β). 2

The straightness of the line (R = 0.99) confirmed the adequacy of the proposed model at 20.0 ( 2.0 °C (Figure 2b). Similarly good correlations (R2 > 0.90) were also observed for the other temperatures investigated (data not shown). As explained above for the estimation of kapp, the parameter determination from a linearized representation is not accurate for statistical reasons. To improve kref app estimation, nonlinear regressions were thus performed at each temperature from adjustment between the experimental data and the nonlinear relation (eq 3), reparameterized as described by eq 13 (Figure 1a). Resulting identified kref app and β parameters are shown in Table 3 with their confidence intervals at the 95% confidence level. These errors were obtained from Monte Carlo simulations taking into account uncertainties previously identified on rate constant kapp and on measured [O2]d. The 2000 values of parameters kref app and β were normally distributed (confirmed by Lilliefors tests), the averages of the normal distributions were kref app and β (the initially identified values) and the confidence intervals at the 95% confidence level from these normal distributions corresponded to the standard deviations on kref app and β multiplied by a factor of 1.96. Table 3 shows that the reference kinetics constant rates kref app significantly increased with temperature, as expected, from 1.3 107 ( 0.6 107 s1 at 8.0 °C to 29 107 ( 4 107 s1 at 33.0 °C. In contrast, the partial order of reaction with respect to oxygen (β) was not significantly different from one temperature to another. Obviously, the Arrhenius law cannot be validated for β, probably because of too high uncertainties on β concomitant to the high-temperature uncertainty ((2.0 °C). The both effects lead to overlapped results for β vs T. To formalize the temperature effect on the kinetics parameter, it is essential to know the relationship between the parameter and the temperature, the equation that models this relationship, and the number of parameters needed to identify in this equation: three items that could not be obtained for β vs T. As a consequence, any mathematical law relating β to the temperature could be identified and β was thus averaged to 1.1 ( 0.9 in the following. It was the first time that this value of β was experimentally determined. In the literature, when partial order with respect to oxygen in ascorbic acid oxidation reaction was assessed, it was supposed to be equal to 1.1115 These studies were carried out under different conditions than in the present work, but our experimental value of 1.1 is not significantly different from 1 and then is consistent with assumptions previously made by other authors. However, the wide uncertainty on

Figure 3. Arrhenius plot describing the temperature dependence of rate constants for ascorbic acid oxidation from 8 °C to 33 °C. Dots represent the experimental temperatures on the x-axis and the identified rate constants on the y-axis, error bars represent confidence intervals at the 95% confidence level, and straight lines represent linear regressions from data.

Table 4. Activation Energies with Their Errors at 95% Confidence Level (When Determined)a Reported in the Literature for Ascorbic Acid Loss in Food Beverages on a Storage Temperature Range activation energy, food beverage fresh blood orange juice high-pressured orange juice

Ea (kJ mol1) 42 61


44


69


53

orange juice

105 ( 35

lemon juice

53 ( 8

grapefruit juice

77 ( 10

tangerine juice

79 ( 6

reference 27 28 29 25

a

The confidence interval at the 95% confidence level was calculated from the standard deviation presented by the authors, multiplied by a factor of 1.96.

the β value obtained in the present case must be evidenced and will be investigated in the following. The range of possible values found for β as a function of the temperature was thus considered as the range of variation of this parameter. 3.1.3. Stage 3. In the third and last stage of parameter identification (Figure 1a), the effect of temperature was tentatively modeled using the Arrhenius law. Because β was considered to be constant, regardless of the temperature, the activation energy thus referred only to kref and could be estimated from eq 8. The plot of the linearized form of the reparameterized version of the Arrhenius law (eq 9) is shown in Figure 3. Uncertainties on both temperature and kinetics constant rates are added to the graph. The adequacy of the Arrhenius law to describe the influence of temperature on the kinetics constant rate was confirmed (R2 = 0.93). For the same statistical reasons than those cited above for incremental estimations of kapp, kref app and β, parameter determination from the linearized form 1138



ARTICLE

Figure 4. Histograms of the Monte Carlo simulations to estimate the confidence interval at the 95% confidence level (IC 95%) on the activation energy for ascorbic acid oxidation (results of 2000 simulations). Two strategies were assessed (Figure 1): incremental identification on experimental data (strategy a) versus simultaneous identification on simulated data (strategy b).

of the Arrhenius law is not desirable. Thus, the regression was performed by adjustment between experimental data and the nonlinear reparameterized relation described by eq 15. The activation energy (Ea) obtained was 93 kJ mol1 and the rate constant at the reference temperature (20.0 °C) kref was 1 102 mol1.1 L1.1 s1 with a RMSE value of 4.8 107 mol1.1 L1.1 s1. Literature data about activation energies for ascorbic acid oxidation are presented in Table 4. The value obtained in this study was not significantly different from that obtained by Burdurlu et al.25 in orange juice (105 ( 36 kJ mol1) but higher than all other values. These reported values were determined under air with simplified kinetics laws involving only ascorbic acid without any influence of oxygen, which could explain the observed differences. In the study of Eison-Perchonok and Downes,15 activation energies were determined as a function of the oxygen level in the atmosphere; the values were 40 ( 2, 57 ( 2, and 69 ( 2 kJ mol1 at 10%, 15%, and 21% oxygen, respectively. These values are lower than those obtained in this study, but the temperature range (3355 °C) can result in a change in the reaction mechanisms and, consequently, the activation energy. Standard deviations of kref and Ea were obtained from Monte Carlo simulations, which took into account the uncertainties previously identified on kinetics constant rates kref app, partial order of reaction with respect to oxygen β, and measured temperatures T. In this case, the 2000 values of kref and Ea resulting of the Monte Carlo simulations were not normally distributed and the confidence intervals were thus estimated to 95% of obtained values, as shown for Ea on the frequency histogram in Figure 4a. The confidence intervals obtained for Ea and kref were dramatically large, asymmetrical, and contained the zero value (248 to 995 kJ mol1 for Ea and 1 105 to 357 mol1.1 L1.1 s1 for kref), which completely withdrawn these two values. Such a large confidence interval rigorously indicates that identification was not possible under these conditions. Considering β to be a constant, regardless of the temperature, with values ranging from 0.1 to 2.0, is a hypothesis that can have important influence on this result, because of this wide range of possible values. Indeed, it cannot be guaranteed that the best values for Ea and kref were obtained for a given value of β. On the other hand, as mentioned previously, any law could be established to relate β to the

temperature; therefore, it was impossible to perform an iterative approach on β to consider variations of β as a function of the temperature T. The uncertainty on the temperature is also high ((2.0 °C). As the uncertainties on kref app, β, and T were combined together to evaluate the uncertainty in Ea and kref, the joined uncertainties on β and T could explain the identification problems. This proves that the error propagation taken into account via the Monte Carlo procedure was very important, especially when incremental identifications of parameters are obliged. Although this procedure is more rigorous to determine the confidence interval of Arrhenius parameters, it is scarcely used by the scientific community. Many authors establish their kinetics parameters from linear regression using the linear form of the Arrhenius law and logarithmic plot of their data. The confidence interval (when it is determined, which is not systematically the case) is calculated from the model error without taking into account the probability distribution of the data in the interval of experimental errors. If the parameters and their confidence intervals were calculated in this study using this simplified approach by linear regression with the Microsoft Excel program, the values gave 78 ( 24 kJ mol1 for Ea and 1 102 ( 4 103 mol1.1 L1.1 s1 for kref. The mean value of Ea was shifted, compared with that obtained by nonliner regression (93 kJ mol1), proving the sensitivity of this parameter to the method used to estimate it. Moreover, the confidence intervals obtained using linear regression with Microsoft Excel were considerably restricted, compared to those obtained by taking into account errors and their propagation during the numerical treatment of the data. A significant underestimation of confidence intervals on identified parameters could be obtained by using a simplified calculation of the confidence interval. 3.2. Improvement of the Accuracy of the Identified Parameters. The results obtained here proved that a specific caution should be taken during the experimentation and the parameters identification to avoid wide uncertainties on measured values, which could trigger the value of confidence interval on the identified parameters. It must be underlined that the choice of a 95% confidence level on all parameters (thus a covering factor of 1.96 on the uncertainties) has an impact on the accuracy of the 1139


Industrial & Engineering Chemistry Research identified parameters. A higher level of confidence would imply larger errors on the data. This does not seem useful since a 95% level of confidence was enough to represent the experimental uncertainty (assessed for example by carrying out several assays on ascorbic acid concentration: each additional assay lay between the boundaries of the 95% confidence interval previously determined). Moreover, the 95% level of confidence proved the very important impact of the errors on the accuracy of the identified parameters. In contrast, a lower confidence interval would not depict the error on experimental data correctly. Consequently, this level of 95% appeared to be the best compromise between the reality to represent and the objective of the paper, and was thus kept in the following. To minimize the effects of error propagation, one should ideally: (a) improve accuracy on the experimental parameters (ascorbic acid assay, dissolved oxygen, and temperature controls); (b) improve the method of determination of the kinetics parameters degradation. For instance, in this study, at low temperatures (814 °C), even with kinetics carried out over two weeks, only ∼20% of the ascorbic acid was lost, especially when low levels of oxygen were imposed (

Influence of the Experimental Errors and Their Propagation on the

Recommend Documents