Kinetic modelling of metabolic pathways

Reposted from biostars

http://bioinformatics.knowledgeblog.org/2011/06/21/kinetic-modelling-of-metabolic-pathways/

Authors

Natalie J Stanford^1,2 , Kieran Smallbone^2,3
1. Integrative Systems Biology Doctoral Training Centre; 2. Manchester Centre for Integrative Systems Biology; and 3. School of Mathematics, University of Manchester, M13 9PL, UK.

Abstract

We show how to build a kinetic model of a metabolic pathway. We provide the example of building a model of glycerol synthesis using Copasi. However, the techniques required remain the same for any pathway, using any software.

Introduction

A mathematical description of a kinetic metabolic model may be given in differential equation form as $latex x^\prime = N v (x, y, p)$ $latex x(0) = x_0$ and this may be used as a guide as to the data required to create and parameterise a kinetic model. First, N is the stoichiometric matrix, which may be derived from the topology of the model. x denotes metabolite concentrations. y denotes boundary metabolites, whose concentrations are not allowed to vary, but do affect the reaction rates. Initial concentrations for both x and y must be defined, though note that only concentrations x will change over time. Finally, v denotes reaction rates; these are dependent on kinetic mechanisms, concentrations x and y and parameters p.

Stoichiometry

The first stage in kinetic modelling is defining the pathway of interest and its boundaries, which are captured mathematically in the stoichiometric matrix. As an example, we shall use a model of glycerol synthesis containing only two reactions [1], though the same procedure can be used for more complex systems. The model’s two reactions are glycerol 3-phosphate dehydrogenase (GPD) – which has reactants DHAP and NADH and products glycerol 3-phosphate (G3P) and NAD – and glycerol 3-phosphatase (GPP), which has reactant G3P and product glycerol (Gly) and phosphate (Phi). We shall show here how to implement the model using Copasi [2], but the techniques we describe will be applicable to any software. To enter the stoichiometric knowledge in Copasi, navigate to
Model > Biochemical > Reactions
Double click on the equation box, which should take you to the reaction interface (see figure 2). In the chemical equation box enter
DHAP + NADH = G3P + NAD
Click commit and navigate back to Model > Biochemical > Reactions

Figure 1: Navigating to the reaction interface

Figure 2: The reaction interface.

Note: You may notice that Copasi automatically defines a rate law. This is fine, and will be addressed later in the tutorial. Double click the new reaction box, and enter the following in the chemical equation box G3P -> Gly + Phi Note: You may have noticed that in reaction 1 we use = to signify the conversion between substrates and products, whilst reaction 2 uses ->. The former is used to represent a reversible chemical reaction. The latter signifies an irreversible reaction. It is important to note that Copasi will not let you assign pre-defined rate laws for reversible reactions that are irreversible and vice versa. However, if you choose to define your own rate law, these constraints can be violated. This will be discussed in more detail when we address entering rate laws. Navigating now to
Model > Biochemical > Species

Figure 3: The species interface

we see that the six metabolic species have already been added to the model. However, we do not want to simulate changes in concentration of the boundary metabolites in our system. To address this, change Simulation Type in DHAP from reactions to fixed, and repeat this for the other boundary metabolites, Gly, NAD, NADH and Phi.

Figure 4: Changing the metabolite boundaries

The stoichiometry matrix N is now found at
Model > Mathematical > Matrices
As expected, it contains only G3P, which is produced by the first reaction and consumed by the second.

Figure 5: The stoichiometric matrix.

Metabolites

In order to build a kinetic model, concentration values must be provided for all metabolites. Typically these will come from metabolomics measurements, or databases such as HMDB [3]. In Table 1, the concentrations for all metabolites in our model are set out. These may be added by accessing the Species menu in Copasi, and entering the concentration as the initial concentration (see figure 6).

Figure 6: Setting the initial concentration.

ADP, ATP and F16BP are included in the model as modifiers of GPD (Reaction 1). This is represented in Copasi by changing the reaction’s definition to
DHAP + NADH = G3P + NAD; ADP ATP F16BP

metabolite	concentration (mM)
ADP	2.17
ATP	2.37
DHAP	0.59
F16BP	6.01
G3P	0
Gly	15.1
NAD	1.45
NADH	1.87
Phi	1

Table 1: Metabolite concentrations used in the model. It is important to ensure that consistent units are used throughout the model. Navigating to Model, we change Time, Quantity Unit and Volume Unit to min, mmol and l, respectively. In this way, all concentrations are measured in units of mmol/l = mM, and all fluxes in units mM/min.

Figure 7: Model unit interface.

Rate equations

Kinetic rate laws v may be derived from knowledge of the mechanism underlying the enzymatic process. In many cases these are not known, and we must instead resort to approximations such as linlog [4, 5], power law [6, 7] or convenience kinetics [8]. For our model, the rate of GPD is known to be
$latex v_1 = \frac{\frac{Vf1}{K1nadh \cdot K1dhap} \left(NADH \cdot DHAP – \frac{NAD \cdot G3P}{Keq1} \right)} { \left(1 + \frac{F16BP}{K1f16bp} + \frac{ATP}{K1atp} + \frac{ADP}{K1adp} \right) \left( 1 + \frac{NADH}{K1nadh} + \frac{NAD}{K1nad} \right) \left(1+ \frac{DHAP}{K1dhap} + \frac{G3P}{K1g3p} \right)}$
whilst the rate of GPP is given by
$latex v_2 = \frac{\frac{V2 \cdot G3P}{K2g3p}} {\left( 1+ \frac{G3P}{K2g3p} \right) \left( 1+ \frac{Phi}{K2phi} \right)}$
Copasi has a wide variety of inbuilt kinetic laws that may be accessed via the Functions menu. However, the forms above must be entered manually; this is achieved by returning to the Reactions menu and selecting GPD (Reaction 1). Click New Rate Law and enter the following
Vf1/(K1nadh*K1dhap)*(NADH*DHAP-NAD*G3P/Keq1)/((1+F16BP/K1f16bp+ATP/K1atp+ADP/K1adp)*(1+NADH/K1nadh+NAD/K1nad)*(1+DHAP/K1dhap+G3P/K1g3p))
Once the equation has been entered you then need to associate each identifier as a Substrate (NADH, DHAP), Product (NAD, G3P), Modifier (F16BP, ATP, ADP) or Parameter(all other identifiers). This can be seen in figure 8

Figure 8: Entering a new rate law.

Repeat for GPP (Reaction 2)
V2*G3P/K2g3p/((1+G3P/K2g3p)*(1+Phi/K2phi))
once again identifying Substrate (G3P), Product (Gly, Phi) or Parameter (all other identifiers). Returning to the Reactionsmenu, these formulae need to be selected again under Rate Law for each reaction. Where there are multiple substrates, products or modifiers within a reaction it is important to check that the correct identifier is associated with the correct metabolite. These are circled in figure 9

Figure 9: Selecting the rate law.

Kinetic parameters

The kinetic formulae above lead to a range of kinetic parameters. Ideally, these would be measured using enzymatic assay, When this is not possible, databases such as Brenda [9] or Sabio-RK [10] could be used. In Table 2, the parameter values in our model are set out. These can be modified in the same interface highlighted in figure 9.

parameter	value	units
V1	147	mM/min
Keq	10000	1
K1adp	2	mM
K1atp	0.73	mM
K1dhap	0.54	mM
K1f16bp	4.8	mM
K1g3p	1.2	mM
K1nad	0.93	mM
K1nadh	0.023	mM
V2	53	mM/min
K2g3p	3.5	mM
K2phi	1	mM

Table 2: Parameter values used in the model.

Analysis

Having built a model of the system, we would like to use mathematical tools to further our understanding of its behaviour. Copasi provides a number of tools to carry out these tasks under the Tasks menu. In the original paper [1], the authors investigate the response of flux through the pathway to changes in boundary metabolite concentrations. This analysis may be performed running the Sensitivities task, with Function: Concentration Fluxes of Reactions and Variable: Initial Concentrations. Its results are presented in Table 3.

metabolite	response
ADP	-0.141
ATP	-0.420
DHAP	0.475
F16BP	-0.162
G3P	0
Gly	0
NAD	-0.0159
NADH	0.0260
Phi	-0.0741

Table 3: Scaled response coeffcients of flux to metabolite concentrations.

Systems biology standards

Describing mathematical models as above is unwieldy and error-prone, and naturally leads to difficulties in reproduction of results. Thus researchers have developed SBML (the Systems Biology Markup Language [11]). A computer-readable format for representing models of biological processes. Exchanging our models in SBML format (rather than e.g. Copasi’s own xml format) not only allows us to describe our model in unambiguous terms. It also gives us access to over 200 software tools specialising in many elements from data integration to visualisation [12]. SBML can be combined with MIRIAM [13] to annotate the entities of those models; for example, by marking-up the molecule G3P as CHEBI:15978 allows its unambiguous identification and automatically links to many additional sources of information. This may be achieved in Copasi under the Species submenu Annotation. BioModels.net [14] is a modelling repository containing hundreds of models marked up with the above standards. If you are having trouble reproducing the glycerol model described above, it is available from the repository at BIOMD0000000076.

References

1. G. Cronwright and J. Rohwer, “Metabolic control analysis of glycerol synthesis in Saccharomyces cerevisiae,” Applied and Environmental Microbiology, vol. 68, pp. 4448–4456, 2002.
2. S. Hoops, S. Sahle, R. Gauges, C. Lee, J. Pahle, N. Simus, M. Singhal, L. Xu, P. Mendes, and U. Kummer, “COPASI – a COmplex PAthway SImulator,” Bioinformatics, vol. 22, pp. 3067–74, Dec. 2006.
3. D. S. Wishart, D. Tzur, C. Knox, R. Eisner, A. C. Guo, N. Young, D. Cheng, K. Jewell, D. Arndt, S. Sawhney, C. Fung, L. Nikolai, M. Lewis, M.A. Coutouly, I. Forsythe, P. Tang, S. Shrivastava, K. Jeroncic, P. Stothard, G. Amegbey, D. Block, D. D. Hau, J. Wagner, J. Miniaci, M. Clements, M. Gebremedhin, N. Guo, Y. Zhang, G. E. Duggan, G. D. Macinnis, A. M. Weljie, R. Dowlatabadi, F. Bamforth, D. Clive, R. Greiner, L. Li, T. Marrie, B. D. Sykes, H. J. Vogel, and L. Querengesser, “HMDB: the Human Metabolome Database,” Nucleic Acids Res, vol. 35, pp. D521–6, Jan. 2007.
4. D. Visser and J. J. Heijnen, “Dynamic simulation and metabolic re-design of a branched pathway using linlog kinetics,” Metabolic engineering, vol. 5, pp. 164–76, July 2003.
5. K. Smallbone, E. Simeonidis, D. S. Broomhead, and D. B. Kell, “Something from nothing: bridging the gap between constraint-based and kinetic modelling.,” The FEBS journal, vol. 274, pp. 5576–85, Nov. 2007.
6. M. A. Savageau, Biochemical Systems Analysis: A Study of Function and Design in Molecular Biology. Addison-Wesley Pub (Sd), 1976.
7. J. J. Heijnen, “Approximative kinetic formats used in metabolic network modeling,” Biotechnology and bioengineering, vol. 91, pp. 534–45, Sept. 2005.
8. W. Liebermeister and E. Klipp, “Bringing metabolic networks to life: convenience rate law and thermodynamic constraints,” Theor Biol Med Model, vol. 3, p. 41, 2006.
9. M. Scheer, A. Grote, A. Chang, I. Schomburg, C. Munaretto, M. Rother, C. Sohngen, M. Stelzer, J. Thiele, and D. Schomburg, “BRENDA, the enzyme information system in 2011,” Nucleic Acids Research, vol. 39, pp. D670–D676, Nov. 2010.
10. U. Leser, F. Naumann, B. Eckman, U. Wittig, M. Golebiewski, R. Kania, O. Krebs, S. Mir, A. Weidemann, S. Anstein, J. Saric, and I. Rojas, Data Integration in the Life Sciences, vol. 4075 of Lecture Notes in Computer Science. Berlin, Heidelberg: Springer Berlin Heidelberg, 2006.
11. M. Hucka, A. Finney, H. M. Sauro, H. Bolouri, J. C. Doyle, H. Kitano, A. P. Arkin, B. J. Bornstein, D. Bray, A. Cornish-Bowden, A. A. Cuellar, S. Dronov, E. D. Gilles, M. Ginkel, V. Gor, I. I. Goryanin, W. J. Hedley, T. C. Hodgman, J.-H. Hofmeyr, P. J. Hunter, N. S. Juty, J. L. Kasberger, A. Kremling, U. Kummer, N. Le Nov`re, L. M.e Loew, D. Lucio, P. Mendes, E. Minch, E. D. Mjolsness, Y. Nakayama, M. R. Nelson, P. F. Nielsen, T. Sakurada, J. C. Schaff, B. E. Shapiro, T. S. Shimizu, H. D. Spence, J. Stelling, K. Takahashi, M. Tomita, J. Wagner, and J. Wang, “The systems biology markup language (SBML): a medium for representation and exchange of biochemical network models,” Bioinformatics, vol. 19, pp. 524–531, Mar. 2003.
12. SBML software guide
13. N. Le Novère, A. Finney, M. Hucka, U. S. Bhalla, F. Campagne, J. Collado-Vides, E. J. Crampin, M. Halstead, E. Klipp, P. Mendes, P. Nielsen, H. Sauro, B. Shapiro, J. L. Snoep, H. D. Spence, and B. L. Wanner, “Minimum information requested in the annotation of biochemical models (MIRIAM),” Nature Biotechnology, vol. 23, pp. 1509–15, Dec. 2005.
14. “BioModels database – A database of annotated published models”

Tags: copasi, metabolic, modelling, simulation

• Citation

  Please cite this article as: Natalie Stanford,Kieran Smallbone (2011) Kinetic modelling of metabolic pathways. The Bioinformatics Knowledgeblog. http://bioinformatics.knowledgeblog.org/2011/06/21/kinetic-modelling-of-metabolic-pathways/

READ MORE...