CC BY
48
ЭЛЕКТРОНИКА
UDC 517.958:541.14
KINFITSIM (VERSION 2.1) - A POWERFUL TOOL FOR KINETIC SIMULATION OF ANY REACTION MECHANISM AND FITTING OF ANY NUMBER OF PAIRS OF THEORETICAL AND EXPERIMENTAL DATA SETS
SVIR I.B., KLYMENKO O.V, OLEINICKA.I., PLATZ M.S.
A new version of the software package KinFitSim (version 2.1) for kinetic simulation and experimental data fitting is presented. In addition to the possibilities offered by the previous version of KinFitSim, which was described before (Computers and Chemistry, V.26, 2002, pp. 379386), new features have been added to the program that include a robust automatic regime for kinetic simulation which simplifies the operation with KinFitSim for the user. KinFitSim 2.1 allows simultaneous fitting of several experimental curves corresponding to the same reaction mechanism that permits better quality of the results due to the better use of input information. The program can be used for both research and educational purposes.
1. Introduction
Homogeneous reactions are of constant interest among chemists and biochemists due to the diversity of situations where they may take place. These range from simple laboratory experiments to extremely complex reactions occurring in living organisms and ecological systems. Despite their different nature all homogeneous chemical reactions can be represented by (systems of) stoichiometric equations. This general approach allows for relatively simple construction of mathematical models of homogeneous chemical processes which is implemented in a number of existing software packages offering kinetic simulation capabilities. The most known and frequently used software are briefly outlined below.
Encora [1] is a program for enzymatic kinetic parameter fitting using progress curve analysis. It was developed by R.J.W. Slats et al. at Delft University of Technology. The program can handle a limited number of reaction mechanisms which constrains its applicability. Differential equations are solved using the fourth order Runge-Kutta method which means that the program would fail to solve a stiff system of ordinary differential equations. The Nelder-Mead method is used for fitting.
EZ-Fit [2] is designed to simplify nonlinear regression analysis of enzyme kinetic data. After entering the data the program may suggests a fit or allow the user to select an equation from the menu of commonly used inhibition
models. Thus the program is not able to simulate an arbitrary kinetic scheme.
KINETICS [3] uses the Gear method for the numerical integration ofmultistep rate equations. Systems containing up to 100 reversible reactions and 50 chemical compounds can be examined. Each individual reaction can contain as many as 5 reactants and 5 products, allowing most systems of interest to be examined. However, the program does not permit fitting experimental data.
KINSIM and FITSIM [4] are perhaps the most popular programs for kinetic simulation and data fitting. This is due to their universal nature: KINSIM can model any system of reactions. However, these programs have severe limitations because they were written for the MS-DOS operating system. The drawback is the absence of a modern user-friendly interface. The mechanism must be written using a separate text editor and compiled before use.
The Chemical Kinetic Simulator (CKS) [5] is based on a stochastic simulation technique, has advanced user interface and a comprehensive manual for users. However, it lacks tools for kinetic parameter determination which makes it useful only for purely theoretical observations or teaching.
Facsimile [6] models complex steady state and time dependent processes. Features are available to allow users to solve both differential and algebraic equations and fit parameters to experimental data. Kinetic simulations are performed using the Gear method with automatic time stepping. The disadvantages of Facsimile are complex model writing procedure and the need to use command language.
As follows from its name, KinFitSim [7] was designed to carry out the kinetic simulation of any given reaction mechanisms and fit experimental data in order to determine the values of kinetic parameters governing the behaviour of a system under study. Both tasks can be performed with ease due to intrinsic straightforwardness of the operation with this software package. With KinFitSim, there are no limitations on the number of chemical steps in the mechanism or the number of species taking part in it. The program has been successfully employed in scientific research in the area of photochemistry [8].
The new version of the program has a number of important features such as the automatic simulation regime and simultaneous fitting of several datasets which make it an essential tool for the investigation of homogeneous chemical and biochemical processes. These features are discussed in the following sections.
2. Kinetic simulation
The core of the KinFitSim package is represented by the kinetic simulator (KS). The Kinetic Simulator parses the entered reaction mechanism and generates its mathematical model which is a system of ordinary differential equations (ODEs). KS comprises various algorithms and methods aimed at efficient simulation of reaction mechanisms. The latter can be performed in either the automatic regime or by manual setting of simulation parameters depending on the user’s choice.
22
РИ, 2004, № 4
The simulation ofa reaction mechanism has never been as easy as with the automatic simulation regime: the only requirement to the user is to specify the reaction scheme in the conventional format and set the values of initial concentrations, rate constants and the desired length of the simulation interval. Then the automatic simulation procedure analyses the mathematical model and solves it numerically ensuring that the results are physically meaningful and precise. This is achieved by applying appropriate methods for the solution of systems of ordinary differential equations and choosing the necessary accuracy of numerical integration. To this end, a specially designed algorithm of discrimination between stiff and non-stiff systems of ODEs [9] is implemented in KinFitSim. This permits the program to automatically select the appropriate numerical integration method, i.e. the Gear method for stiff systems and the Adams-Moulton method for non-stiff ones. In some cases it may happen that according to a priori information an ODE system appears to be non-stiffbut then becomes stiff in the course of the simulation (i.e. for some t>0). This may lead to erroneous solutions because of wrong choice of the solution method. However, the algorithm implemented in the automatic simulation subroutine of KinFitSim is capable of overcoming such situations and choosing an appropriate method in any case.
The accuracy of a numerical method (defined by the limit of local integration error on integration time step) is also chosen automatically based on the previously proven theorem and criterion [9] which ensure that the resulting concentration-time curves do not suffer from non-physical oscillations. This means that KS begins calculations with reasonably high accuracy and then augments it if necessary. The importance of this issue is demonstrated by the following example.
Thus the automatic simulation regime makes KinFitSim an extremely convenient and robust instrument for modelling homogeneous reaction mechanisms which can be exploited by the user with any level of mathematical skill.
At the same time an interested user can benefit from the manual simulation regime where he/she can learn more about the numerical methods used for the solution of ODEs by choosing and tuning them at their discretion. The program can be conveniently switched to the manual simulation regime where the user can choose between the Runge-Kutta and Gear methods (if the latter is chosen and the ODE system is non-stiff the Adams-Moulton procedure is applied), set the desired accuracy of integration and other simulation parameters [7]. By manually setting the simulation method parameters the user will be able to observe the effects of different simulation parameters on the results of computation. These settings are especially important when a stiff system of OD Es is to be solved. In this case the program may fail when running in a manual simulation regime due to numerical solutions of stiff ODE systems being very sensitive to applied numerical methods and their accuracy. Therefore in this case the program suggests the user to apply the automatic regime which is a recommended option for any ODE system regardless of its stiffness.
3. Fitting experimental data
The fitting simulator (FS) is designed for automatic fitting of theoretically predicted time dependent distributions to those obtained experimentally. This is achieved by iteratively running KS with various values of kinetic parameters in order to try and minimise the difference between experimental and theoretical curves.
In the new version of the KinFitSim package, FS offers a new possibility of fitting multiple pairs of curves. This is useful in a situation when a reaction under study is characterised by several different quantities, which can be independently measured. These quantities can be, for example, solution absorbance at different wavelengths.
The importance of such a feature is shown by the following example. Consider a reaction scheme:
A —^ B, (1)
B + C k2 > D , (2)
D + E k3 > A + F (3)
with initial concentrations [A]0=1M, [C]q=0.5M,
[E]0=0.7M, [B]0=[D]0=[F]0=0 and rate constants k1=10s-1, k2=k3=10cm3mol-1s-1. By simulating reactions
(1)-(3) with KS and adding 2% Gaussian noise to the concentrations of A, B, D and F we obtain spiky curves which replicate experimental ones measured with random fluctuations (see Figure 1, which shows a snapshot of the program window designed for loading and editing experimental data).
Next we assume that the values ofthe rate constants are unknown to us and try to fit the ‘experimental’ curves from Figure 1 using FS given an initial guess for the rate constants of k1=1s-1, k2=1cm3mol-1s-1,
k3=20cm3mol-1s-1. A usual way of doing this would be to fit the experimental curve for [A] to the simulated one and repeat this for the rest of the curves. The results of this procedure are shown in Table 1, where the Nelder-Mead minimisation algorithm [10] was used with an accuracy of 10-4.
It is clear that the results cannot be regarded as satisfactory because the resulting values of the rate constants are far from their true values (the deviation approaches 400% for kil). This result was obtained because even though each experimentally measured curve contains all the information about the reaction mechanism it has different sensitivity to different kinetic parameters. For instance, the concentration distribution of species A is much more sensitive to k1 than it is to k2 and k3, which is observed in Table, i.e. the value of k1 resulting from the fitting is much closer to the true one than those of k2 and k3. Similar conclusions can be drawn from the other three fitting results which, however, are not as straightforward.
Fitting results for the reaction mechanism (1)-(3)
Fitted curve ki / s-1 k2 / cm3mol'1s'1 k3 / cm3mol'1s'1
[A] 8.6099 1.0708 1.3051
[B] 5.2051 0.1428 2.3245
[D] 49.986 5.6705 9.7073
[F] 33.130 6.8927 9.7240
РИ, 2004, № 4
23
Let us now fit all the ‘experimental’ concentration distributions simultaneously instead of fitting them one by one. In this case the fitting results are far better than those given in Table 1: ki=9.990s_1, k2= 10.138cm3mol-1s-1, k3=9.922cm3mol-1s-1. The deviation from the true values of the rate constants does not exceed 1.5% and can be regarded as excellent given the significant amplitude of noise on the ‘experimental’ curves. Such a drastic difference in the fitting results with the previous approach is due to the fact that by fitting all the data simultaneously we even out the lack of sensitivity of certain curves to certain kinetic parameters (rate constants). Therefore the fitting procedure converges to the true solution in all kinetic parameters.
Consequently, the possibility of fitting multiple data sets is a very important feature which allows obtaining significantly more accurate fitting results provided that several different experimental quantities characterising the same kinetic process are available.
Apart from this new feature the FS has become much more reliable and robust compared to previous versions of the program. This was achieved by adding the Nelder-Mead minimisation method to the FS and optimising the work of methods which were implemented in FS previously (these are the direct search algorithm, Gradient search, Gauss-Newton and Marquardt methods). The Nelder-Mead method usually requires less computation compared to Gradient search, Gauss-Newton and Marquardt methods and therefore is typically more efficient than other methods.
Although the recommended way of fitting data is the automatic fitting procedure the user can perform this task manually by adjusting the values of kinetic parameters, running KS and visually evaluating the quality of fit. This may be a good idea to find a reasonable initial fit between experimental and theoretical curves manually before running the automatic procedure,
especially when the user is in possession of any additional information about possible values of kinetic parameters. This will reduce the time necessary to achieve a good fit by the FS.
4. Program interface
KinFitSim has a well designed interface which is simple to understand for an untrained user (see Figure 2). The new version of the KinFitSim software includes a comprehensive user manual. Context help is available in most program windows which facilitates obtaining necessary information immediately when it is needed.
The kinetic mechanism can be saved together with the simulation data and user’s comments for later use. This facilitates repeated simulation of the mechanism and storing research archives.
The output of the program is also fully controllable by the user: it may be not only the concentrations of reacting species but also any algebraic functions of concentrations that permit direct comparison of simulated curves to experimental data. This feature significantly broadens the range of applicability of the KinFitSim software. The simulation or fitting results can be conveniently presented in the form of graphical images and reports containing comprehensive information about the initial data, solution process and its results.
5. Conclusions
The new version of KinFitSim (2.1) has many advantages in comparison with previous versions of KinFitSim and other programs that are currently available. The principal advantages are the automatic simulation regime in KS, which allows simplification of the kinetic simulation for users of any mathematical skill, and a possibility to fit several pairs of theoretical and experimental curves at the same time in FS.
Fig. 1. ‘Spiked’ concentration distributions of species A, B, D and F in the program window for loading experimental data
24
РИ, 2004, № 4
The overall quality of the program has also been improved which results in faster simulation and fitting and various improvements in the interface. All this makes KinFitSim software a valuable tool for research and chemical education.
Acknowledgements
The authors wish to thank the users of all previous versions of KinFitSim for their valuable comments and suggestions, which have allowed us to develop this new improved version.
References: 1. Straathof A.J.J. Development of a computer program for analysis of enzyme kinetics by progress curve fitting // Journal of Molecular Catalysis B: Enzymatic. 2001. V. 11. P. 991-998. 2. http: //www.jlc.net/~fperrell/webps04.htm.
3. http://www.ari.net/ars/kin.html. 4. Barshop B.A., Wrenn R.F., Frieden C. Analysis of numerical methods for computer simulation of kinetic processes: development of KINSIM — a flexible, portable system // Analytical Biochemistry. 1983. V. 130. P. 134-145. 5. http://www.almaden.ibm.com/st/ msim. 6. http://www.esm-software.com/facsimile. 7. Svir I.B, Klymenko O. V, Platz M.S. ‘KINFITSIM’ - a software to fit kinetic data to a user selected mechanism // Computers and Chemistry. 2002. V. 26. P. 379-386. 8. Tippmann E., Platz M, Svir I., Klymenko O. Evidence for Specific Solvation of Two Halocarbene Amides // Journal of the American Chemical Society 2004. V.126, P. 5750-5762. 9. Klymenko O. V. Methods for the determination of stiffness and required accuracy of numerical integration of mathematical models describing homogeneous chemical processes //
Radioelektronika i Informatika. 2004. V. 3. P. 42-47. 10. Himmelblau D.M. Applied nonlinear programming. Moscow: Mir, 1975. 534 p.
Поступила в редколлегию 24.11.2004
Рецензент: д-р техн. наук, проф. Хаханов В.И.
Svir Irina Borisovna, Chief Scientist, Dr. Tech. Sci., Head of the Mathematical and Computer Modelling Laboratory, KhNURE. Scientific interests: mathematical modelling and numerical methods in research. Address: KhNURE, Kharkov, 61166, 14 Lenin Avenue. Email:
svir@kture. kharkov .ua.
Klymenko Oleksiy Victorovich, PhD, Senior Scientist of the Mathematical and Computer Modelling Laboratory, KhNURE. Scientific interests: mathematical modelling and numerical methods in research, mathematical physics, programming. Address: KhNURE, Kharkov, 61166, 14 Lenin Avenue. Email: klymenko@kture.kharkov.ua.
Oleinick Alexander Igorevich, Scientist of the Mathematical and Computer Modelling Laboratory, KhNURE. Scientific interests: mathematical modelling and numerical methods in research, mathematical physics, programming. Address: KhNURE, Kharkov, 61166, 14 Lenin Avenue. Email: oleinick@kture.kharkov.ua.
Platz Matthew S., Distinguished Professor of Chemistry, Department of Chemistry, Newman and Wolfrom Laboratory, The Ohio State University. Scientific interests: organic chemistry and photochemistry. Address: 100 West 18th Avenue, Columbus, OH 43210-1183, USA.
Fig. 2. A snapshot of KinFitSim windows including “Help”, “Report” about fitting results and graphical plots after fitting
РИ, 2004, № 4
25
