16th European Symposiumon Computer Aided Process Engineering and 9th International Symposiumon Process SystemsEngineering W. Marquardt, C. Pantelides (Editors) © 2006 Published by Elsevier B.V.
Model discrimination and parameter estimation through sensitivity analysis Mauricio Sales-Cruz, Rafiqul Gani
CAPEC, Department of Chemical Engineering, Technical University of Denmark, DK2800 Lyngby, Denmark Abstract In this paper, a systematic procedure to support the building of process models through the use of a computer-aided modeling system, called ICAS-MoT, is presented. Specifically, the step-by-step procedure allows the decomposition of the model solution into several sub-problems: first, a sensitivity analysis helps in the parameter discrimination and gives the order of the parameters to be estimated, then the parameters are estimated using either least squares method or maximum likelihood method, and finally the applicability of the model parameters is evaluated by determining the confidence limits of the estimated parameters. The application of these computer-aided tools is highlighted through a complex kinetic parameter estimation problem.
Keywords: computer-aided modeling, model discrimination, parameter estimation, sensitivity analysis.
1. Introduction Commonly, model-based optimization problems are related to (kinetic) parameter estimation, optimal start-up/shutdown operation, time scheduling, reactor design, optimizing process control, analysis of dynamic systems and so on. To solve these problems, several altemative models are often proposed to explain the same data, and objective criteria are needed to choose among models. The altemative models may, in some cases, be nested. Nested models are constructed such that a simpler model can be obtained from a more complex model by eliminating one or more parameters from the more complex model. Thus the model discrimination problem may be reduced to the appropriateness of adding or removing some of the related parameters. The question of increased model complexity with or without additional model parameters needs to be investigated with respect to computational cost, accuracy and parameter uncertainty. In general, the more parameters contained in a model, the less reliable are the parameter estimates. Criteria to select among models must weigh the trade-off between increased information and decreased reliability. Moreover, the model refinement is an iterative procedure where, experimental and expert guided process of adding, deleting, or modifying assumptions until a model that satisfactorily explains the data is obtained is in general a difficult and time consuming task. The objective of this paper is to highlight the model parameter estimation and model discrimination features within the computer-aided modeling system called ICAS-MoT (Sales-Cruz and Gani, 2003). Specifically, model refinement and model parameter discrimination are done through a sensitivity analysis, together with a statistical analysis, which are tools included in ICAS-MoT. As a case study, the propane aromatization on zeolite catalyst is considered. The large number of reactions and reactive species involved gives a better process description but complicates the process
625
626
M. Sales-Cruz and R. Gani
modeling as reported by Katare et al. (2004a; 2004b) and Lukyanov et al. (1995). To refine the model, a sensitivity analysis is performed by testing different values of the parameters in a systematic manner so that the less sensitive parameters can be discriminated/eliminated. Numerical optimization methods (i.e. least square or maximum likelihood method) are used to estimate the parameters and statistical information is obtained to determine the reliability of the model parameters.
2. Model construction strategy For building and discriminating nested models, the first step involves specifying a mathematical formulation of a detailed model that describes the process. Then three stages should be followed: Stage I: model parameter discrimination, the strategy uses this model to perform parameter sensitivity analysis before any data are collected, in order to determine whether or not the model parameters can be identified so that some of them can be discriminated, breaking down the parameter estimation problem into sub-problems; Stage II: model refinement, the strategy involves designing experiments for parameter estimation (using least squares or maximum likelihood method) and model adequacy checking to determine whether or not the proposed model can be refined; and Stage III: model validation, experiments are designed for improving the precision of the parameters within that model, eventually arriving at a final, statisticallyverified model formulation.
3. Computer-aided modeling system There are many programs for fitting nonlinear models to experimental data, and the use of these software is now widespread. In this case, the model discrimination procedure (i.e. the model screening, comparison and improvement) is done systematically through our computer-aided modeling system, called ICAS-MoT (Sales-Cruz and Gani, 2003) that allows performing the following modeling steps in a fast, reliable and efficient way: (1) model formulation, (2) model solution, (3) sensitivity analysis, (4) model parameter optimization, and (5) statistical analysis of results. An important feature of ICAS-MoT is that the model developer does not need to write any programming codes to enter the model equations. Models are entered (imported) as text-files or XML-files, which are then internally translated. In the model analysis step ICAS-MoT orders the equations into lower triangular form (if feasible), generates the incidence matrix, verifies the degrees of freedom, and checks for singularity. After this interactive model analysis, the appropriate solver for the model equations is selected together with a corresponding strategy of solution. Then the sensitivity analysis option allows testing different values of the parameter, which can help in the model parameter discrimination task. For model parameter regression (Englezos and Kalogerakis, 2001), the least square method (by minimizing the sum of the square errors between the experimental values and the predicted ones) or maximum likelihood method (by maximizing the probability distribution function of the estimation error) can be used. For the statistical analysis of results, ICAS-MoT performs analysis of variance, correlation factors calculations, F significance test and confidential intervals of the parameters estimated (Montgomery and Runger, 1999).
4. Case study: Kinetic model for propane aromatization 4.1. Process description Propane aromatization is a complex heterogeneous reaction that involves many
Model Discrimination and Parameter Estimation Through Sensitivity Analysis
627
individual reactions species and many more chemical reaction (one species can participate in many reactions). Katare et al. (2004a, 2004b) have developed a reduced kinetic model by lumping all the isomers of the same carbon number into a single component. The schematic mechanism of the reactions considered in the kinetic model is shown in Fig. 1.
C2H4 H2 ~
H4 4 __~C2H6
aliphatics/ olefins
H+l
aliphatics/ olefins
1
Diolefins
H+~"~ LA . . . . tics]
\ CH3+ - H+ \~ CzH5+ ,v....................... C3H7+ ReactionIntermediates Fig. 1. Reaction pathways for propane aromatization over HZSM-5 (Katare et al., 2004a).
4.2. Model formulation In this work, the detailed model for propane aromatization was taken from Venkatasubramanian (2005), who using information about the relative rates of the reactions, stability of the intermediates, quantum chemical calculations and empirical correlations from the literature lumped the reactions into 33 reaction families, so as to reduced the number of model parameters to 13. The reaction rate model where the kinetic model is embedded, consists of a differential-algebraic equation system with 29 algebraic equations (representing the mass balance in the active sites) and 31 ordinary differential equations (representing the concentration of the chemical species in the gas phase), the 13 kinetic model parameters found in the model are classified in this paper into four parameter types (as shown in Table 1). Table 1. Model parameters classified according type. Parameter Description Reference value Type Protolysis of carbonium ions 200.0 a Carbonium ion desorption 11.3x 104 a Carbenium ion desorption 853.0 a Beta-scission 25.6x 104 a Aromatization 6.7x 107 a Carbonium ion dehydrogenation 94.0 a Alkane adsorption 47.0 b Hydride transfer reactions 100.0 b Olefin adsorption 8923.0 b Increase in adsorption enthalpy for alkanes 6.0 c with carbon number k4 Increase in activation energy for alkanes 6.0 c desorption with carbon number k~ Increase in activation energy for carbonium 2.0 c ion dehydrogenation with carbon number k12 Entropy factor to determine the equilibrium 18.0 d between beta-scission and oli~omerization a: First order rate constants (in terms of mol/g/h), b: second order rate constants (in me/g/h), c: Energy terms (in KJ/mol), and d: the entropy term (normalized by the universal gas constant). Parameter k0 kl k3 k5 k6 kl0 k7 k8 k9 k2
M. Sales-Cruz and R. Gani
628 The generic form of the mathematical model is as follows: dx
--~=f(x,y,O,t);
g(x,y,0,t) = 0
(1)
Where vector x represents the concentration of the chemical species in the gas phase, vector y represents the concentration in the catalyst sites and vector 0 represents the kinetic model parameters.
4.3. Sensitivity analysis A sensitivity analysis is performed to evaluate the effect of each parameter, by perturbing them with respect to their reference value (given in Table 1) and computing the output response in terms of the percentage of change of the CsH8 composition (the main reactive component in the reaction). Alternatively, plots of the analytical derivatives would also provide the same information. Figures 3 and 4 show the results of this sensitivity analysis: (a) parameters k12 (Fig. 4b) and k6 (Fig. 3a) are the least sensitive (i.e. the change in output response is less than 1%); (b) parameters k5 (Fig. 3a), k8 (Fig. 3b) and k9 (Fig. 3b) are also not very sensitive as changes in output responses vary by up to 2.0%; (c) parameters k2 and kll (Fig. 4a) are lightly more sensitive as changes output responses vary by up to 5%, (d) parameters kl (Fig. 3a), k4 (Fig. 4a) and klo (Fig. 3a) are sensitive as changes in output responses vary by up to 15%, (e) parameters ko (Fig. 3a) and k7 (Fig. 3b) are very sensitive as changes in output responses vary by up to 20%, and finally, (f) the most sensitive parameter is k3 (Fig. 3a), as changes in output response varies by up to 100%. Based on the above analysis, the model is rather insensitive (as changes in output responses within 5 %) to seven parameters k12, k6, ks, k8, kg, k2 and k11, but highly sensitive to the other six parameters (kl, k4, klo, ko, k7 and k3).
First order rate constants
100
~k 0 ....
k1
......
k3
.....
k5
.......
k6
.......
kl0
Second orde rate constants 15-
k7 ....
ks
10-
5-
0
-120
|
|
-80
-40
,
0
.r.=.7-
40
% parameter change
(a)
80
~,
0
120
!
-60
-40
L
........... |
-20
0
......
-. . . . .
2'0
.-,
4'0
% parameter change
(b)
Fig. 3. Sensitivity analysis: (a) first order rate constants and (b) second order rate constants.
60
Model Discrimination and Parameter Estimation Through Sensitivity Analysis 1413 12 11
Entropy parameter
k~
k12
-- -k 4
~0 9 ~" C,.)
0.8
Energy parameters
....
629
0.6
kll
l 1
8 7
t
0.4
/ l
6 5 4 3 2
iI
°
//
0.2
o'°"
1 0 -60
~o
-~o
0
|
|
20
40
0.0 60
-10
% parameterchange (a)
-5
0
5
l0
% parameterchange (b)
Fig. 4. Sensitivity analysis: (a) energy parameters and (b) entropy parameter.
4.4. Strategy for parameter estimation and model validation Based on the above analysis, the following strategy for model parameter regression is proposed for stage II (i.e. model parameter estimation). Using available experimental data, the parameter identification problem is split into six sub-problems following the order of sensitivity of each parameter: (1) first all the parameters are fixed at their reference value, then the least sensitive parameter k12 (and also the only one parameter of type d) is regressed using the least squares minimization method; (2) then, taking this new value of the parameter k12, the following four less sensitive parameters (ks, k6, ks, k9) a r e regressed; (3) with these five parameters (k12, ks, k6, ks, k9) regressed then the parameters corresponding to the energy terms (k2 and k11) are regressed; (4) with these seven new values (k12, ks, k6, 1:8, k9, k2 and k11), the parameters kl, k4 and klo are regressed; (5) having these ten parameters regressed (k12, ks, k6, k8, k9, k2, k11, kl, k4 and klo), then the parameters ko and k7 a r e regressed; and (6) finally, with these twelve parameters regressed, the most sensitive parameter k3 has to be regressed. For the model refinement, once all these parameters have been identified, a simultaneous regression must be done but with a reduction in the parameter bounds to determine whether or not any of the estimated (step by step) model parameter changes with purpose of selecting the best model parameters. And lastly, for stage III (model validation), experiments must be designed in order to improve the precision of the parameters estimated, eventually arriving at a final statistically-verified model formulation. Following the aforementioned strategy of solution for stage II, the regressed parameters are presented in Table 2. Afterwards, when the simultaneous regression of all parameters was done, as it was expected, only the regressed values of parameters k6 and k12 (the least sensitive ones) did not change. That means that the model could be reduced by eliminating these two parameters (]i76 and k12). To verify the accuracy of the model parameters, the following statistics of the regression was calculated: correlation coefficient = 0.972549, standard error of the estimates = 13.612893 and F significance = 5.206348e-002. In addition the 95% confidence limits were calculated and are reported in Table 2. It can be noted that the most uncertain parameters also have the widest bounds, the need for which is a topic of further investigation.
630
M. Sales-Cruz and R. Gani
Table 2. Model parameter regression. Parameters
Botmds
Step by step Simultaneous value value k0 102 < k 0 < 107 217.8642 5.6472X105 kl 104 < kl < 101° 132579.0 69122.67 k2 6 < k2 < 14 6.82674 6.0 k3 103 < k3 < 109 1025.72 1.1277x107 k4 6 < k4 < 12 8.648333 12.0 k5 10 3 < k5 < 108 238174.3121 255194.49 k6 107 < k6 < 1013 6x107 6x107 k7 10-3 < k7 < 102 4.7773 46.107 k8 10-3 < k8 < 102 125.87432 99.7586 k9 10-1< k9 < 104 7863.5726 8922.80 kl0 102< kl0 < 108 93.3911 125.475 kll 2 < kll < 6 3.78924 6.0 k12 18 < k12 < 25 18.0 18.0 * Confidential interval of the simultaneous regression
+95% CI* 59680.8614 7305.03699 0.6340933 1191778.36 1.2681866 26969.5194 6340933 4.87268997 10.54271 942.981283 13.2604761 0.6340933 2.4389342
4.5. Other case studies
The ICAS-MoT system for model discrimination and parameter estimation has also been tested with all the problems reported by Katare et al. (2004a). In all cases, better optimized parameters were obtained (detailed results can be obtained from the authors).
5. Conclusions In this paper, a computer-aided procedure for parameter discrimination and parameter estimation of nested models has been presented, illustrating the procedure with a complex case study related to the propane aromatization kinetics. In general, a systematic model building (i.e. model parameter discrimination, model refinement and model validation) permits a substantial saving in experimental (time and cost) effort. Specifically the parameter model regression was done based on a sensitivity analysis, obtaining helpful information about the parameter identifiability. Based on the model parameter discrimination results, additional experimental data are being selected for final model validation. Once the model is ready the ICAS-MoT system also provides a customized simulator of the process for future design/optimization studies. References Englezos, P., Kalogerakis, N., 2001, Applied Parameter Estimation for Chemical Engineers, Marcel Dekker, New York. Froment, G.F., 1975, AIChE J., 21 (6), 1041-1057. Katare, S., Bhan, A., Caruthers, J.M., Delgass, W.N., Venkatasubramanian, V., 2004a, Comp. Chem. Eng., 28, 2569-2581. Katare, S., Caruthers, J.M., Delgass, W.N., Venkatasubramanian, V., 2004b, Ind. Eng. Chem. Res., 43, 3484-3512. Lukyanov, D.B., Gnep, N.S., Guisnet, M.R., 1995, Ind. Eng. Chem. Res., 34, 516-523. Montgomery, D.C., Runger, G.C., 1999, Applied Statistics and Probability for Eng., Wiley, NY. Sales-Cruz, M. and R. Gani, 2003, Computer-Aided Chemical Engineering, vol. 16: Dynamic Model Development, Eds. S.P. Asprey and S. Macchietto, Elsevier, Amsterdam.
Model Discrimination and Parameter Estimation Through Sensitivity Analysis
631
Venkatasubramanian, V., Laboratory for Intelligent Process Systems (LIPS), Purdue University, personal communication, January, 2005.