Taking Variable Correlation into Consideration during Parameter Estimation

Taking Variable Correlation into Consideration during Parameter Estimation

T.J. Santos¹ and J.C. Pinto2,^*

¹Departamento de Estatística - Universidade Tiradentes - Aracaju, SE - Brazil

²Programa de Engenharia Química / COPPE - Universidade Federal do Rio de Janeiro

Cidade Universitária - CP 68502 - Rio de Janeiro 21945-970 RJ - Brazil

(Received: April, 24, 1997; Accepted: August 18, 1997)

Abstract - Variable correlations are usually neglected during parameter estimation. Very frequently these are gross assumptions and may potentially lead to inadequate interpretation of final estimation results. For this reason, variable correlation and model parameters are sometimes estimated simultaneously in certain parameter estimation procedures. It is shown, however, that usually taking variable correlation into consideration during parameter estimation may be inadequate and unnecessary, unless independent experimental analysis of measurement procedures is carried out.

Keywords: Parameter estimation, maximum likelihood, correlation.

Introduction

During the development of mathematical models, one often has to evaluate certain variables which cannot be measured by independent experimentation and which are called the model parameters. The model parameters are usually evaluated by fitting model responses to experimental responses. Such a procedure is called parameter estimation. Parameter estimation procedures generally lead to the minimization of a scalar function which measures the inadequacy of model responses. This scalar function is normally presented as a weighted sum of model deviations from experimental results in the form:

(1)

where is the response of the j^th dependent model variable at the i^th experiment, which depends on the independent experimental conditions xi and on the model parameters b .

The literature regarding parameter estimation has been largely concentrated on the development of numerical techniques for fast and reliable minimization of nonlinear quadratic scalar objective functions such as that presented by Equation (1). Bard (1970, 1974) and Edgar and Himmelblau (1988) presented a very complete compilation of numerical methods for minimization of nonlinear quadratic scalar functions. Special techniques have also been designed to exploit particular functional forms of problems which are commonly found in chemical engineering. Britt and Luecke (1973), Anderson et al (1978) and Schittkowski (1994) presented numerical procedures for parameter estimation of nonlinear algebraic equations and analysed the performance of the proposed procedures in thermodynamic problems. Froment (1975) presented a review of parameter estimation procedures applied to heterogeneous kinetic problems. A significant amount of work has been devoted to the development of numerical procedures for parameter estimation of distributed systems described by ordinary and partial differential equations (for example, see Leal et al. 1978, Hosten 1979, Frenklach and Miller 1985, Rutzler 1987, Özgülsen et al. 1992, Lohmann et al. 1992, Mangold et al. 1994). More recently special attention has also been given to parameter estimation for systems of algebraic-differential equations (for example, see Biegler et al. 1986, and Bilardello et al. 1993).

The problem of estimating parameters is often presented simultaneously with the problem of designing and analysing experiments. In this case, experiments are designed and parameters are estimated in order to allow model discrimination, to improve parameter accuracy and to reduce parameter correlation (for example, see Blau et al. 1972, Hosten 1974, Hosten and Emig 1975, Pritchard and Bacon 1978, Ferraris et al. 1984, Rippin 1988, Espie and Macchietto 1989, Pinto et al. 1990, 1991, Doví 1993). For experimental design and nalysis to make sense, a somewhat rigorous statistical interpretation of the objective function to be minimized and of the final parameter estimates is required. One possible way to introduce a rigorous statistical interpretation into the parameter estimation problem is the definition of the maximum likelihood function.

Let us assume that both the dependent and independent variables are subject to experimental errors. In this case, the experimental results may be seen as random variables which can be described by a certain joint probability density function

(2)

which shows the probability of obtaining and in experiment i, given the real unknown experimental values and and a measure of the experimental error and . The maximum likelihood estimation technique consists of maximizing Equation (2), given the model constraints

, (3)

i = 1, ..., NE, j = 1, ..., NY

where it is assumed that the model is perfect and that the true unknown values of the dependent and independent variables are related through the model equations. As pointed out by Bard (1974) if the model equations are perfect and the experiments are well done, it seems plausible to accept that the experimental results actually observed are the most likely ones, so that it seems reasonable to vary the model parameters in order to maximize the probability of observing the experimental results obtained. It is important to emphasize that Equations (2-3) depend on NP+NE*NX parameters, where NP is the dimension of the vector of model parameters and NE*NX is the number of independent variables whose true values at each experiment cannot be evaluated independently. For this reason, it is usual to assume that the independent variables are known with great precision, which means that the additional constraints are imposed on the problem:

i = 1...NE (4)

When the experimental deviations are distributed normally and the model is assumed to be perfect, maximizing the joint probability density function is equal to minimizing the function

(5)

where Z is a vector which contains the dependent and independent variables and V_Z is the covariance matrix of the Z deviations shown below

(6)

(7)

which has the form of Equation (1) when the independent variables are assumed to be free of error.

When the experiments are carried out independently, it is generally written

(8)

It is also usual to assume that dependent and independent variables are not correlated, so that

(9)

Equations (5-9) show that the experimental errors are very important to describe the final parameter estimates. Mezaki et al. (1973) and Pritchard and Bacon (1975) analysed how parameter estimates change when the actual response variables are replaced by predictor variables in the problem formulation. Although using the least-squares technique may be justifiable when actual responses are used for parameter estimation, the technique is generally improper when transformations of the actual responses are used to formulate the objective function. Nevertheless, weighting procedures may be implemented successfully in these cases if explicit transformations are used. Graham and Stevenson (1972) analysed the influence of experimental errors on the convergence of a sequential experimental design and observed that large experimental errors may cause a lack of uniform convergence. Ricker (1984) compared the performance of traditional least-squares and maximum likelihood algorithms and concluded that the additional computer efforts of maximum likelihood methods would not necessarily lead to significant gains in parameter accuracy. Kim et al. (1991) observed, however, that maximum likelihood methods which take into account the errors of all experimental variables can provide both better parameter estimates and more accurate model predictions. More recently, Straja (1993) considered that experimental errors in standard first-order kinetic experiments should follow a log-normal distribution and analysed the differences observed when both log-normal and normal distributions were used to formulate the maximum-likelihood objective function.

Very often the covariance matrix of the experimental errors is totally or partially unknown. This is particularly true when the experimental errors are correlated. In this case, the unknown components of the covariance matrix are usually considered as additional parameters to be estimated with the model parameters. Pritchard and Bacon (1977), Pritchard et al. (1977) and Stewart et al. (1992) present some functional forms for adequate description of the unknown covariance components and analyse the performance of some estimation procedures. Although a large number of numerical procedures has been presented for parameter estimation of models described by differential equations, most of them have just neglected the fact that the experimental errors are naturally correlated in these systems, as an error in the boundary conditions causes deviations in all other experimental values.

According to the formulation presented, it has been implictly assumed that the deviations observed between model predictions and experiments can be described by the experimental error distribution. Actually this is not true, even when the model is perfect. The reason for this is due to the parameter estimation procedure, which is based on uncertain data. As experimental data are uncertain, parameter estimates and model predictions are also uncertain. Moreover, as parameter estimates are based on all experimental data included in the objective function, model predictions for a certain experimental condition depend on all the other experimental values, which means that the covariance matrix of model predictions is dense even if the covariance matrix of experimental errors is diagonal. This is an aspect that has been neglected in previous studies. Our main objective here is analysing the importance of taking variable correlation into consideration during parameter estimation when we have and when we do not have an independent experimental evaluation of measurement correlations.

Theoretical Framework

Let a mathematical model be described by

(10)

where the dependent variables are assumed to be explicit functions of the independent variables and of the model parameters. This assumption does not imply any significant restriction on the problem formulation, as it is assumed that numerical techniques may be used to solve the model when it is not explicit. When Y*, X* and ? are subject to small errors, the following equation may be written:

(11)

where B_x and Bb are the sensitivity matrixes for perturbations on the dependent variables and on the model parameters, and are respectively given by

(12)

(13)

While the deviations D Y* may be assumed to be equal to the experimental errors e _y when the model is perfect, the terms D X* and D b are related to the parameter estimation technique. As shown by Bard, 1974,

where

(15)

(16)

(17)

Combining Equations (11) and (14),

where

(19)

which shows how model uncertainties are correlated with experimental errors during parameter estimation. It may be seen that even when the model is perfect, prediction error distributions are much more complex than the experimental error distributions.

Equations (14) and (18) may be averaged as

where

(22)

so that

(24)

which is the covariance matrix of prediction errors. So, this is the matrix which should be used in Equation (5) to weight the differences observed between model and experimental results. Observe that, due to the matrix operations, V_Z* is full even when VZe is diagonal. To make this point clear, let us assume that the independent variables are free of error and that the experiments are independent of each other. In this case,

(25)

(26)

(27)

Therefore, predictions Yi* and Yj* are correlated and the objective function (25) does not represent the distribution of deviations observed between model and experiments. It is then proposed that the matrix defined by Equation (24) be used as the weighting function of the maximum likelihood estimation procedure.

A problem in the formulation of the parameter estimation procedure is the fact that V_Z* depends on the parameter estimates, which are not known beforehand. If the minimization problem is posed in terms of finding the roots of the set of nonlinear equations

(28)

then V_Z* may be updated during the iterative procedure. However, using minimization algorithms to find the minimum of F is a much more common and reliable procedure to estimate b and X* (see Bard 1974 and Edgar and Himmelblau 1988). In this case, V_Z* should not be changed during minimization, as the numerical scheme might become unstable due to the frequent modifications of F, which is the function to be minimized. Therefore, the parameter estimation procedure must be iterative, in the sense that after finding the minimum of F_k, V_Z*k is computed and a new maximum likelihood objective function F_k+1 is defined and minimized. The iterative procedure must be repeated if certain tolerance criteria, which measure the rate of convergence, are not satisfied. It seems that the most natural tolerance is the Euclidean norm of parameter variation

(29)

although additional convergence criteria may be imposed to both F and V_Z*.

We do not intend to develop here the detailed conditions that guarantee the convergence of the proposed numerical scheme. However, it may be interesting to show some results regarding convergence of the proposed algorithm in simple problems. First of all, convergence on F may always be attained if a simple linear transformatin of F_k+1 is used in the form

(30)

where

(31)

and c _k is the vector of parameter estimates after k iterations. Thus, F_k(c _k) = F_k+1(c _k). As F _k+1 is positive and decreases during the minimization procedure, F necessarily converges as k increases. It is important to emphasize that the linear transformation in Equation (30) does not change the parameter estimates and the statistical interpretation of the final results.

Regarding the parameter estimates, it is interesting to analyse how the algorithm behaves when the model is linear. Let us assume that the independent variables are free of error, so that F_k has the general form

(32)

where

(33)

If

(34)

then it may be shown that the minimum of F_k+1 is at

(35)

Equations (26) and (27) may then be written as

(36)

(37)

When k is equal to zero, Equations (36) and (37) reduce to

(38)

(39)

In order to compute W₁, the Matrix Inversion Lemma may be used (see Ogata 1987), leading to

(40)

Using Equation (33)

(41)

so that

(42)

According to Equations (35), (36) and (37), it may be seen that convergence is attained after just one iteration. Therefore, results are the same when either experimental or model prediction errors are used to define the maximum likelihood function. These results indicate that taking the variable correlation which is induced by the model structure and the parameter estimation procedure into consideration is not important at all. In a more general nonlinear problem, the results are similar to the ones already presented, as a linear approximation may always be used successfully in the neighborhood of the minimum. The same may be said when the independent variables are subject to experimental errors. Therefore, induced model and parameter estimation correlations may be completely neglected.

As F_k+1 follows a chi-square distribution when the deviations are normal, W_k may be corrected as suggested by Anderson et al. (1978) in the form

(43)

so that F _k+1 assumes its most probable value. This technique is very useful when experimental errors are not known and does not change any of the results presented above. However, as it will be shown below, it should be used with caution when experimental errors are known.

Numerical Results

Examples 1, 2 and 3 presented below illustrate how variable correlation may change the final parameter estimation results. In the first and second examples simulated problems are presented to show how the parameter estimates change when variables are assumed to be correlated. The third example regards the characterization of the activity of a Ziegler-Natta catalyst used in propylene polymerization reactions. The main objective of the third example is to show how variable correlation may change the experimental sequence during a sequential experimental plan.

Example 1 - Parameter Estimation for a Dynamic Kinetic Model

Reagent A is converted to the product B in a batch stirred tank reactor. The reaction is promoted by a catalyst and may be described by a first-order reaction rate expression. Temperature is kept constant and the catalyst activity presents an exponential decay. Then, the mass balance may be written as:

(44)

which may be integrated as

(45)

where x is conversion, t is time, k_O is the initial catalyst activity and k_d is the catalyst decay rate constant.

An experimental conversion dynamic profile is available, as shown in Table 1, and parameter estimates must be obtained. (Experimental results shown in Table 1 were generated at random, assuming that k_O and k_d were equal to 1.35 and 0.35, respectively.) It is assumed that the standard deviation of the independent variable is equal to 0.01 h. The triplicates available when t is equal to 1 hour are used to evaluate the experimental standard conversion deviation, which is then found to be equal to 0.029. Experimental errors are assumed to follow the normal distribution.

Different procedures were used to evaluate model parameters. Table 2 shows results obtained when the independent variables are assumed to be free of error, while Table 3 shows results obtained when independent variables are assumed to be subject to errors. In all cases, parameters were obtained with ESTIMA, a software for parameter estimation and experimental design which uses the Gauss-Newton procedure to find the minimum of F (see Noronha et al. 1993). The initial guesses of k_O and k_d were always set to 2.3 and 1.5, while the relative tolerances were set to 1.x10^-4.

Table 1:
Experimental Conversion Dynamic Profile

t (h)	0.001	0.001	1.0	1.0	1.0	2.0	3.0
x	0.01	0.02	0.652	0.710	0.682	0.875	0.905

Table 2:
Parameter Estimates When Independent Variables Are Free of Error - Example 1

Internal Iterations	Parameter kO	Parameter kd	Standard Deviation of kO	Standard Deviation of kd	Parameter Correlation
Model predictions in different experiments are independent - Independent variables are free of error - Experimental errors are known
4	1.37395	0.366795	0.885694	1.0087	0.927534
6	1.35799	0.344143	0.884152	1.0224	0.928354
7	1.36298	0.352161	0.879434	1.01378	0.927774
8	1.3614	0.349552	0.883499	1.01605	0.928054
9	1.36182	0.350245	0.882276	1.01543	0.927973
Model predictions in different experiments are independent - Independent variables are free of error - Experimental errors are unknown
4	1.37395	0.366795	0.106725	0.121548	0.927534
6	1.35799	0.344143	0.881573	1.01942	0.928354
7	1.36298	0.352161	0.87902	1.0133	0.927774
8	1.3614	0.349552	0.883452	1.016	0.928054
9	1.36182	0.350245	0.882273	1.01543	0.927973
Model predictions in different experiments are correlated - Independent variables are free of error - Experimental errors are known
4	1.37395	0.366795	0.885694	1.0087	0.927534
5	1.37391	0.366734	0.885427	1.00866	0.927524
Model predictions in different experiments are correlated - Independent variables are free of error - Experimental errors are unknown
4	1.37395	0.366795	0.885694	1.0087	0.927534
5	1.37391	0.366734	0.885427	1.00866	0.927524

Table 3:
Parameter Estimates When Independent Variables are Subject to Experimental Error - Example 1

Internal Iterations	Parameter kO	Parameter kd	Standard Deviation of kO	Standard Deviation of kd	Parameter Correlation
Model predictions in different experiments are independent - Independent variables are not free of error - Experimental errors are known
4	1.37335	0.365524	0.398876	0.486362	0.935603
6	1.36678	0.355408	0.426711	0.50237	0.929629
8	1.3653	0.355127	0.405547	0.490821	0.933469
Model predictions in different experiments are independent - Independent variables are not free of error - Experimental errors are unknown
4	1.37335	0.365524	0.0464095	0.0565885	0.935603
6	1.36678	0.355408	0.0201973	0.0237785	0.929629
7	1.36712	0.357642	0.406009	0.491571	0.933595
Model predictions in different experiments are correlated - Independent variables are not free of error - Experimental errors are known
4	1.37335	0.365524	0.435074	0.500487	0.928394
5	1.37335	0.365524	0.435074	0.500487	0.928394
Model predictions in different experiments are correlated - Independent variables are not free of error - Experimental errors are unknown
4	1.37335	0.365524	0.435074	0.500487	0.928394
5	1.37335	0.365524	0.435074	0.500487	0.928394

The first two sets of numerical results presented in Tables 2 and 3 were obtained when the iterative procedure was used to update the covariance matrix of prediction errors. In these cases, the covariance matrix of prediction errors was assumed to be diagonal (in spite of some huge prediction variable correlations), although, as already discussed, this matrix is known to be full. Numerical results show fast convergence to final parameter estimates, which are not identical but are very close to the results obtained when the full covariance matrix is used, as shown in the following two sets of numerical results. Particularly, the second set of numerical results in Tables 2 and 3 show estimates obtained when Equation (43) is used to update experimental errors. It is interesting to observe how the iterative procedure of the second set of results converges smoothly to the results obtained in the first set. However, as parameter estimation would probably not be continued after the first iteration, it is important to observe how parameter uncertainties were underpredicted after applying the correction procedure. In this case, the final value of F after the first iteration is approximately equal to 0.08. If the data obtained when t is equal to zero are discarded (as conversions were set to zero on purpose, to avoid handling negative conversions), then, according to the c ² distribution, the 99% confidence interval of F, with three degrees of freedom, is equal to 0.0717 < F < 12.84. Therefore, the value of F obtained is not unlikely and correction should be carried out with caution to avoid underestimating actual parameter uncertainties. Similar results are obtained when the F-test is used to compare the experimental errors measured with the replicates and with the value of F. The third and fourth sets of numerical results of Tables 2 and 3 also confirm that induced model and parameter estimation correlations do not modify the final parameter estimation results.

In the cases discussed above, experimental measurements were assumed to be independent. As discussed by Pritchard and Bacon (1977), Pritchard et al. (1977) and Stewart et al. (1992) results may be improved if variable correlation is

also estimated during the parameter estimation procedure. Figures 1 and 2 show numerical results obtained when the experimental values obtained at similar experimental conditions are assumed to be correlated and when all experiments are assumed to be correlated respectively. Therefore, the covariance matrixes of experimental errors are assumed to have the forms

Obviously, variable correlations are not necessarily equal, although this simplifying assumption is made here in order to allow an easier interpretation of final numerical results. Figures 1 and 2 show numerical results obtained when the variable correlation factor (a ) varies within the interval [-1,+1]. It may be seen from Figures 1 and 2 that parameter estimates are not very sensitive to a , although large changes are observed in the stardard deviations of the final parameter estimates. More important, if a is taken as an additional parameter, F does not present a minimum. Even when a is kept constant, the numerical scheme may fail in a large range of reasonable values of a .

In order to interpret these results properly, it is necessary to remember that Equation (5) has a mimimum when V_Z is positive definite, which means that the eigenvalues of V_Z must be positive. In the first case, the eigenvalues are positive when a is in the range [-0.5,+1], while in the second case, the eigenvalues are positive when a is in the range [-0.15,+1]. When a is outside these intervals, F does not have a minimum and its stationary point (the point where the gradient of F is equal to zero) is a saddle. As F increases monotonically with a in both cases, the simultaneous estimation of a necessarily leads to the failure of the numerical scheme.

Example 2 - Parameter Estimation for a Linear Model

An interesting point regarding the simultaneous estimation of model parameters and variable correlations is the fact that the parameter estimates depend nonlinearly on variable correlation factors. It is well known that nonlinear functions may present complex behavior, such as multiple solutions, and that singularities arise at certain special points (called bifurcation points) where the number and / or stability of solutions change (see Iooss and Joseph 1980 for details). As may be seen below, these singular responses are likely to occur when a is varied in estimation problems, even when the model is linear.

Let us assume that the experimental data shown in Table 4 are available and that a best straight line is to be estimated. The maximum likelihood method must be used. The independent variables are assumed to be free of error and the following covariance matrix describes the experimental errors of the dependent variables, where a must be estimated:

x -1 -0.5 0 0.5 1 y -0.9 -0.4 -0.1 0.6 0.9

x	-1	-0.5	0	0.5	1
y	-0.9	-0.4	-0.1	0.6	0.9

Table 4: Experimental Data for the Linear Model

Figure 1:
Simultaneous Estimation of Model Parameters and Variable Correlation - Similar Experiments Are Correlated

Figure 2:
Simultaneous Estimation of Model Parameters and Variable Correlation - All Experiments Are Correlated

According to the discussion presented before, a must be a real number in the interval [-1,+1]; otherwise, a minimum does not exist. In this case, parameter estimates are the solutions of the following set of algebraic equations:

(46)

(47a)

(47b)

Figures 3 and 4 show how parameter estimates depend on a . The existence of discontinuities may be clearly seen from these figures. It is important to emphasize that the discontinuities are observed in the range of acceptable values of a . The discontinuities constitute the borders of the basins of attractions of three different minima, which are located at a =0.085, a =0.625 and a =0.970. The discontinuities are observed when a is equal to -1, 0.18, 0.93 and 1. At these points, the system of linear equations (47) is not well defined, as its discriminant is equal to zero and infinite solutions are possible. Figure 5 shows that the models described by each of the minima obtained may be considered significantly worse than the model obtained when it is assumed that there is no correlation between the variables. Therefore, a smaller value of F does not necessarily imply that the model has been improved, as correlation terms may artificially add negative terms to the objective function.

Figure 3:
Simultaneous Estimation of Model Parameters and Variable Correlation - Existence of Multiple Solutions in the Usual Straight Line Estimation

Figure 4: Simultaneous Estimation of Model Parameters and Variable Correlation

- Zoom of Figure 3

Figure 5:
Simultaneous Estimation of Model Parameters and Variable Correlation - Model Responses at the Different Minima

Example 3 - Parameter Estimation for Actual Ziegler-Natta Olefin Polymerization

In order to investigate the kinetics of a certain Ziegler-Natta catalyst, a series of kinetic experiments were carried out. The main objective of carrying out this series of experiments was to characterize the catalyst activity by performing semi-batch slurry propylene (PP) polymerizations in a stirred tank reactor. The polymerization takes place at isothermal conditions and it is verified experimentally that reaction rate is approximately constant during the whole batch. At such experimental conditions, the following model may be written

(48)

where Q is the polymer yield, A and D E are the kinetic parameters to be evaluated, Cat is the catalyst mass, [Cm] is the monomer concentration, R is the universal gas constant, T is the reaction temperature, T_O is a reference temperature and t is the batch time. The monomer concentration in the reaction medium [Cm] is not measured, but calculated as a function of reactor temperature and pressure (P). It is assumed that both gas and liquid phases are in equilibrium and that both phases may be described by the Soave-Redlich-Kwong (SRK) equation of state. A full description of the experimental apparatus and procedure is beyond the scope of this paper and has been presented in detail by Santos, 1995.

This type of olefin polymerization is very sensitive to the presence of impurities (inhibitors), so that one may wonder whether experimental errors of consecutive experiments may be correlated, as feed lines and reactor internals cannot be completely cleaned after reaction is over due to experimental constraints. So, it is proposed that the error correlation of consecutive experiments and the other two process parameters be estimated simultaneously. In order to improve the quality of the final parameter estimates and reduce the costs of carrying out the experiments, a sequential experimental design procedure is used to define the most adequate experimental conditions for parameter estimation.

Sequential experimental design procedures for parameter estimation have been discussed in detail in previous papers (see Pinto et al. 1990, 1991). Different design criteria have been proposed and most of them are based on properties of the confidence region of the parameter estimates, which may be described by the posterior covariance matrix of parameter estimates

(49)

which represents the expected covariance matrix of the parameter estimates if the additional experiment (k+1) is carried out. Equation (49) may be derived if, among other assumptions (see Bard 1974), experimental errors of different experiments are assumed to be independent, which is not true in the case under study. Therefore, a new equation must be derived from Equation (20). If the Gauss-Newton approximation is used (see Bard 1974) and if the independent variables are free of error, then

(50)

where contains the covariances of experimental errors of experiments i and j. Equation (50) may be presented in an iterative form which resembles Equation (49). However, as individual sensitivity matrixes Bb i and covariance matrixes appear explicitly in the iterative form, Equation (50) is certainly more convenient. It is important to notice that the lines and columns i,j = 1...k are known, as it is assumed that k experiments have already been performed. Lines and columns k+1 depend on the experimental conditions to be chosen.

Based on the discussion presented before, it was assumed that the covariance matrix of experimental errors had the tridiagonal form

(51)

where s _y² was evaluated from independent replicates as 0.1 mol² of PP. Then, an initial set of experimental results was used to allow the estimation of A, D E and a . These experiments are presented in Table 5 and were carried out by introducing small perturbations at usual operation conditions. After estimating the parameters assuming either that a was equal to zero or that a was not equal to zero, a sequential experimental design procedure was used with both relative volume and relative b -trace design criteria (see Pinto et al. 1991). Equation (50) was used to calculate the posterior covariance matrix of parameter estimates. The experimental grid of possible experiments was built by allowing temperature to vary within the range 40 ^oC £ T £ 75 ^oC, pressure to vary within the range 7 kgf/cm²£ P £ 9 kgf/cm² and batch time to vary within the range 1 h £ t £ 3h. It was assumed that the amount of catalyst used in each experiment would always be equal to 2 mmol. Results obtained are shown in Tables 6 and 7 and Figures 6 and 7. It is important to notice that T_O was made equal to the average of the experimental temperature data.

Table 5:
Initial Set of Experimental Conditions - Example 3

Experiment	Cat (mmol)	P (kgf/cm2)	T (K)	t (h)	Q (PPmol)
1	4.62	8	63	2	8.30
2	5.12	7.5	62	2.5	9.82
3	4.01	7	60	2	5.67

Table 6:
Sequential Experimental Design and Parameter Estimates - Example 3

Experiment	Parameter Estimates					Selected Experiment
Results Available	a	A (l/mmol h)	D E (cal/mol K)	T (oC)	P (kgf/cm2)	t (h)	Cat (mmol)
Relative b -Trace Design With Varying a
3	-0.601	340.8 ± 7.14	11400 ± 10300	75	9	3	2
4	-0.601	420.6 ± 2.85	14380 ± 1280	40	8	3	2
5	-0.202	308.3 ± 12.4	13950 ± 1210	50	7	1	2
6	-0.102	277.7 ± 13.4	13870 ± 1150	40	7	3	2
Relative Volume Design With Varying a
3	-0.601	340.8 ± 7.14	11400 ± 10300	75	9	3	2
4	-0.601	420.6 ± 2.85	14380 ± 1280	75	7	3	2
5	-0.551	478.4 ± 3.6	14590 ± 800	40	9	3	2
6	-0.102	361.9 ± 12.9	14170 ± 1170	40	7	1	2
Relative b -Trace Design With a =0
3	0	413.5 ± 16.3	9602 ± 10300	40	9	3	2
4	0	247.6 ± 15.7	12800 ± 2080	75	9	3	2
5	0	309 ± 14.1	13790 ± 1130	40	9	3	2
6	0	253.3 ± 12.2	13540 ± 920	75	8	3	2
Relative Volume Design With a =0
3	0	413.5± 16.3	9602 ± 10300	75	9	3	2
4	0	418.7 ± 15.7	14160 ± 1370	75	9	3	2
5	0	476.0 ± 15.7	14450 ± 1150	75	9	3	2

1 Border between regions of convergence and divergence

² Minimum of the maximum likelihood objective function

T (oC) P (kgf/cm2) Cat (mmol) t (h) Q (PPmol) 75 9 2.02 3 10.09 40 9 2.14 3 2.87 40 8 2.56 3 2.87 75 9 1.78 3 9.21 50 7 1.90 1 0.90 40 9 2.58 3 3.47

T (oC)	P (kgf/cm2)	Cat (mmol)	t (h)	Q (PPmol)
75	9	2.02	3	10.09
40	9	2.14	3	2.87
40	8	2.56	3	2.87
75	9	1.78	3	9.21
50	7	1.90	1	0.90
40	9	2.58	3	3.47

Table 7: Experimental Results - Example 3

Figure 6:
Simultaneous Estimation of Model Parameters and Variable Correlation - Example 3

Figure 7:
Experimental Design Sensitivity to Variable Correlation

Figure 6 shows that, when a is estimated simultaneously with the other process parameters in the first stages of the sequential plan, the objective function does not have a minimum. This leads to a very optimistic evaluation of parameter uncertainties, although parameter estimates do not change significantly, as they are not very sensitive to a in the range -1 £ a £ 0. Table 6 also shows that the sequential experimental plan depends on whether a is assumed to be equal or different from zero. It is interesting to observe that, when a is estimated, parameter uncertainties increase when a minimum appears in the objective function. Another interesting point is that the sequential plans are very similar, if they are compared as a whole. The selected experiments are concentrated around the extreme temperature conditions, at conditions which maximize polymer production (larger pressures and batch times). This is caused by the larger relative parameter uncertainties observed in the activation energy. As pointed out by Pinto et al. (1990) it may be observed in the last lines of Table 6 that the volume design tends to concentrate experiments in a very narrow range of experimental conditions, which may bias the parameter estimation procedure.

Figure 7 shows that the experimental design may be relatively insensitive to changes in a . Figure 7 shows how the selected experiment changes when a is allowed to vary if the relative b -trace design criterion is used to design the seventh experiment. The "turbulence" observed in this figure is generated when one or more of the characteristic values of Vb becomes negative. Therefore, in the case analysed, the most important effect caused by a in the sequential plan is changing parameter estimates and parameter uncertainties during parameter estimation.

If all experiments available are used to estimate the parameters, the results obtained are A = 245 ± 11.7 l/mmol h, D E = 13750 ± 745 cal/mol K when a is allowed to vary (a =0.07), and A = 245 ± 11.2 l/mmol h, D E = 13740 ± 750 cal/mol K when a is equal to zero. It can probably be concluded that variable correlation between successive experiments is equal to zero.

Conclusions

The importance of taking variable correlation into consideration during parameter estimation was analysed. It was shown that variable correlation introduced by the model structure during the parameter estimation procedure may be neglected. It was also shown that the existence of experimental variable correlation may cause changes of the final parameter estimates and parameter uncertainties and of sequential experimental designs, although both parameter estimates and uncertainties and experimental designs are often insensitive to variable correlations. It was shown that, when experimental variable correlation is not known, the simultaneous estimation of model parameters and variable correlation may lead to meaningless parameter estimates and parameter uncertainties, as the objective function often does not have a minimum or has multiple minima when variable correlation is allowed to vary. In the first case, the practical minimum is the border which separates the regions where V_y is positive definite and where it is not. At such a point, prediction of parameter uncertainties is usually too optimistic and parameter correlation is equal to ± 1. For all these reasons, it seems that variable correlation should only be considered when it can be evaluated experimentally and independently from the parameter estimation procedure and when it is observed that parameter estimates are sensitive to changes in the variable correlation.

Nomenclature

a Straight line angular coefficient

A Kinetic constant

b Sensitivity, straight line linear coefficient

B Sensitivity matrix

Cat Mass of catalyst

[Cm] Monomer concentration

f Explicit model constraints

F Objective function

g Implicit model constraints

G Matrix of cross second derivatives

H Hessian matrix of second derivatives

I Identity matrix

k Kinetic constant

NE Number of experiments

NF Number of degrees of freedom

NX Number of independent variables

NY Number of dependent variables

P Probability density function

Q Polymer yield

R Universal gas constant

t Time

T Temperature

V Covariance matrix

x Independent variables, conversion

W Matrix of weighting factors

y Dependent variable

z Augmented vector of variables

0 Zero matrix

Greek

a Variable correlation

b Parameter

c Augmented vector of model parameters

D Deviation

D E Activation energy

s Standard deviation

X Defined by Equation (34)

Special Subscripts

o Reference state

Special Superscripts

e Experiment

m Model

* Real unknown value

T Transpose

Acknowledgement

We thank CNPq - Conselho Nacional de Desenvolvimento Científico e Tecnológico - for providing scholarships for both authors and for supporting part of this research. We also thank Polibrasil S.A. for supporting part of this research.

References

Anderson, T.F.; Abrams, D.S. and Grens II, E.A., Evaluation of Parameters for Nonlinear Thermodynamic Models. AIChE Jl., 24, 20-29, (1978).

Bard Y., Comparison of Gradient Methods for the Solution of Nonlinear Parameter Estimation Problems. SIAM Jl. Numer. Anal., 7, 157-186, (1970).

Bard, Y., Nonlinear Parameter Estimation, Academic Press, San Diego, (1974).

Biegler, L.T.; Damiano, J.J. and Blau, G.E., Nonlinear Parameter Estimation: A Case Study Comparison. AIChE Jl., 32, 29-45, (1986).

Bilardello, P.; Joulia, X.; Le Lann, J.M.; Delmas, H. and Koehret, B., General Strategy for Parameter Estimation in Differential-Algebraic Systems. Comp. Chem. Engng., 17, 517-525, (1993).

Blau, G.E.; Klimpel, R.R. and Steiner, E.C., Nonlinear Parameter Estimation and Model Distinguishability of Physicochemical Models at Chemical Equilibrium. Can. Jl. Chem.

Engng., 50, 399-409, (1972).

Britt, H.I. and Luecke, R.H., The Estimation of Parameters in Nonlinear, Implicit Models. Technometrics, 15, 233-247, (1973).

Doví, V.G.; Reverberi, A.P. and Maga, L., Optimal Design of Sequential Experiments for Error-in-Variables Models. Comp. Chem. Engng., 17, 111-115, (1993).

Edgar, T.F. and Himmelblau, D.M., Optimization of Chemical Processes. McGraw-Hill, New York, (1988).

Espie, D. and Macchietto, S., The Optimal Design of Dynamic Experiments. AIChE Jl., 35, 223-229, (1989).

Ferraris, G.B.; Forzatti, P.; Emig, G. and Hofmann, H., Sequential Experimental Design for Model Discrimination in the Case of Multiple Responses. Chem. Engng. Sci., 39, 81-85, (1984).

Frenklach, M. and Miller, D.L., Statistically Rigorous Parameter Estimation in Dynamic Modeling Using Approximate Empirical Models. AIChE Jl., 31, 498-500, (1985).

Froment, G.F., Model Discrimination and Parameter Estimation in Heterogeneous Catalysis. AIChE Jl., 21, 1041-1057, (1975).

Graham, R.J. and Stevenson, F.D., Effects of Experimental Error on Parameter Estimation and Convergence of a Sequential Experimental Design. Ind. Engng. Chem. Proc. Des. Dev., 11, 164-169, (1972).

Hosten, L.H., A Sequential Experimental Design Procedure for Precise Parameter Estimation Based Upon the Shape of the Joint Confidence Region. Chem. Engng. Sci., 29, 2247-2252, (1974).

Hosten, L.H., A Comparative Study of Short Cut Procedures for Parameter Estimation in Differential Equations. Comp. Chem. Engng., 3, 117-126, (1979).

Hosten, L.H. and Emig, G., Sequential Experimental Design Procedures for Precise Parameter Estimation in Ordinary Differential Equations. Chem. Engng. Sci., 30, 1357-1364, (1975).

Iooss, G. and Joseph, D.D., Elementary Stability and Bifurcation Theory. Springer, New York, (1980).

Kim, I.-W.; Edgar, T.F. and Bell, N.H., Parameter Estimation for a Laboratory Water-Gas-Shift Reactor Using a Nonlinear Error-in-Variables Method. Comp. Chem. Engng., 15, 361-367, (1991).

Leal, D.J.; Georgantzis, G. and Roberts, P.D., Parameter Estimation in Uncertain Models of Nonlinear Dynamic Systems. Electron. Lett., 14, 718-720, (1978).

Lohmann, T.; Bock, H.G. and Schlöder, J.P., Numerical Methods for Parameter Estimation and Optimal Experiment Design in Chemical Reaction Systems. Ind. Engng. Chem. Res., 31, 54-57, (1992).

Mangold, M.; Lauschke, G.; Schaffner, J.; Zeitz, M. and Gilles, E.-D., State and Parameter Estimation for Adsorption Columns by Nonlinear Distributed Parameter State Observers. Jl. Proc. Cont., 4, 163-172, (1994).

Mezaki, R.; Draper, N.R. and Johnson, R.A., On the Violation of Assumptions in Nonlinear Least Squares by Interchange of Response and Predictor Variables. Ind. Engng. Chem. Fundam., 12, 251-254, (1973).

Noronha, F.B.; Pinto, J.C.; Monteiro, J.L.; Lobão, M.W. and Santos, T.J., ESTIMA - Um Pacote Computacional para Estimação de Parâmetros e Projeto de Experimentos. Internal Report, PEQ / COPPE, Universidade Federal do Rio de Janeiro (in Portuguese), (1993).

Ogata, K., Discrete-Time Control Systems. Prentice-Hall, Englenwood Cliffs, NJ, (1987).

Özgülsen, F.; Adomaitis, R.A. and Çinar, A., A Numerical Method for Determining Optimal Parameter Values in Forced Periodic Operation. Chem. Engng. Sci., 47, 605-613, (1992).

Pinto, J.C.; Lobão, M.W. and Monteiro, J.L., Sequential Experimental Design for Parameter Estimation: A Different Approach. Chem. Engng. Sci., 45, 883-892, (1990).

Pinto, J.C.; Lobão, M.W. and Monteiro, J.L., Sequential Experimental Design for Parameter Estimation: Analysis of Relative Deviations. Chem. Engng. Sci., 46, 3129-3138, (1991).

Pritchard, D.J. and Bacon, D.W., Statistical Assessment of Chemical Kinetic Models. Chem. Engng. Sci., 30, 567-575, (1975).

Pritchard, D.J. and Bacon, D.W., Accounting for Heteroscedasticity in Experimental Design. Technometrics, 19, 109-115, (1977)

Pritchard, D.J. and Bacon, D.W., Prospects for Reducing Correlations Among Parameter Estimates in Kinetic Models. Chem. Engng. Sci., 33, 1539-1543, (1978).

Pritchard, D.J.; Downie, J. and Bacon, D.W., Further Consideration of Heteroscedasticity in Fitting Kinetic Models," Technometrics, 19, 227-236, (1977).

Ricker, N.L., Comparison of Methods for Nonlinear Parameter Estimation. Ind. Engng. Chem. Proc. Des. Dev., 23, 283-286, (1984).

Rippin, D.W.T., Statistical Methods for Experimental Planning in Chemical Engineering. Chem. Engng. Sci., 12, 109-116, (1988).

Rutzler, W., Nonlinear and Adaptive Parameter Estimation Methods for Tubular Reactors. Ind. Engng. Chem. Res., 26, 325-333, (1987).

Santos, T.J., Estimação de Parâmetros e Projeto de Experimentos em Sistemas Dinâmicos. Internal Report, PEQ / COPPE, Universidade Federal do Rio de Janeiro, (1995).

Schittkowski, K., Parameter Estimation in Systems of Nonlinear Equations. Numer. Math., 68, 129-142, (1994).

Stewart, W.E.; Caracotsios, M. and Sorensen, J.P., Parameter Estimation from Multiresponse Data. AIChE Jl., 38, 641-650, (1992).

Straja, S., Parameters Estimation for a First Order Irreversible Reaction. A Stochastic Approach. Chem. Engng. Comm., 124, 165-175, (1993).

* To whom correspondence should be addressed.

To whom correspondence should be addressed.

Publication Dates

• Publication in this collection
  09 Oct 1998
• Date of issue
  Mar 1998

History

• Accepted
  18 Aug 1997
• Received
  21 Apr 1997