Kinetic Modelling of Metabolic Pathways

Kinetic modelling is a powerful computational tool employed in the study of biochemical systems, offering a dynamic perspective on the intricate processes governing cellular functions. At its core, kinetic modelling involves the mathematical representation of biological reactions, capturing the rates at which molecules interact and transform within a system. Kinetic ordinary differential equation (ODE) models have been traditionally used for dynamic optimization of culture conditions in a bioreactor. Regardless of the type of the models, the process of model building is typically iterative, combining wet-lab experiments and in-silico analysis and optimization. Kinetic models provide a quantitative framework to understand the temporal evolution of complex biological networks. This approach allows researchers to simulate and analyse the behaviour of biochemical pathways, providing insights into the dynamics of cellular processes.

A major hurdle lies in the precise determination of kinetic parameters. Obtaining comprehensive and accurate data on reaction rates, enzyme kinetics, and molecular interactions can be labour-intensive and subject to variations in experimental conditions. Additionally, the sheer complexity of biological systems poses a challenge in developing models that accurately represent interconnected pathways and regulatory networks. Bridging the gap between theoretical predictions and experimental validations remains a persistent challenge.[1]

For our model, we used the enzyme database BRENDA[2]. The BRENDA (BRaunschweig ENzyme DAtabase) Enzyme Database stands as a comprehensive and authoritative resource in the field of enzymology. Developed and maintained by the Department of Bioinformatics at the Technische Universität Braunschweig, Germany, BRENDA provides a wealth of information on enzyme functional data, including kinetics, molecular and cellular properties, and environmental parameters. Other available databases include Swiss-Prot/UniProt[3], KEGG (Kyoto Encyclopedia of Genes and Genomes)[4], and many others.

Kinetic Modelling in general

What We Did and Why?

During our iGEM cycle, we developed a comprehensive kinetic model for the fatty acid synthesis pathway in yeast. Our model encompasses the intricate interactions within the glycolysis pathway and encodes data for all kinetic parameters associated with each enzymatic reaction involved in the pathway.

The primary objective of constructing this kinetic model was to compute the theoretical yield of palmitate, the initial target product for our project. Through time course simulations of the complete model, we aimed to gain insights into the metabolite flow within the reaction network. This approach enables us to identify key metabolites and parameters influencing palmitate yield, providing a deeper understanding of the underlying metabolic dynamics.

Additionally, our model serves as a foundation for estimating the production cost of palmitate. By incorporating data on metabolite inputs, media, bioreactor costs, and other external factors, we aimed to provide a holistic view of the economic feasibility of large-scale palmitate production. Optimization strategies were explored to enhance palmitate yield, followed by a revised cost estimate, contributing to the economic viability of our proposed biofuel production process.

Notably, our work fills a critical gap as there was no existing kinetic model for fatty acid synthesis in yeast. We believe that, with refinement and experimental validation, our model can serve as a novel tool for understanding and optimising fatty acid synthesis. While our current data is based on Saccharomyces cerevisiae and Homo sapiens, we envision adapting our model with experimental data from Yarrowia lipolytica, a more efficient lipogenic yeast. This extension holds significant promise for the biofuel and energy industries, providing a platform for informed decision-making and large-scale biofuel production.

Our theoretical approach, driven by the constraints of our iGEM cycle, lays the groundwork for future experimental validation and refinement. The proposed model stands as a valuable resource for researchers and industry professionals interested in advancing sustainable biofuel production using yeast, particularly Yarrowia lipolytica. The integration of experimental data from this promising lipogenic yeast will further enhance the accuracy and applicability of our model, contributing to the broader goal of developing environmentally friendly and economically viable biofuel production processes.

Through our project, we have created a preliminary version of a potential, small scale dynamic metabolic model, which we aim to validate and modify using experimental data in the future.

Procedure and Results

Data collection

The kinetic parameters for the glycolysis and fatty acid synthesis pathways in Saccharomyces cerevisiae and some from Homo sapiens were systematically gathered from the BRENDA database, ensuring a comprehensive dataset for subsequent model construction.

Enzyme Kinetics

Most of the rate laws of the reactions in our model were approximated to linear, single substrate Michaelis-Menten kinetics. A few reactions use the Mass-Action law instead. More information on the approximations for these rate laws is mentioned in the Approximations and Assumptions section below in this page itself.

Our approaches

We initially attempted to use the software COPASI[5] to build the model. COPASI carries out simulation and analysis of biochemical network models. Models in COPASI are based on reactions that convert a set of species into another set of species. COPASI automatically converts the reaction network to a set of differential equations or to a system of stochastic reaction events — the user does not have to write down the equations explicitly, the software does that. COPASI can import and export models in the SBML format (levels 1 to 3). Simulation can be performed either with stochastic kinetics or with differential equations, and the software easily allows switching between them. Despite the easy to use GUI interface of COPASI, we faced some fundamental issues while building the model, and we were unable to troubleshoot and resolve them due to the inability to see the code behind the algorithms.

Thus, we decided to write an SBML[6] file, which encapsulates our model. This file can be loaded on GUI softwares like COPASI, and other programming languages and softwares such as Python and MATLAB. SBML (Systems Biology Markup Language) is a standardised format facilitating the exchange of computational models in systems biology. It employs an XML-based syntax to represent biological processes, allowing for model sharing and comparison across different software tools. SBML promotes collaboration and model integration, serving as a crucial tool for encoding mathematical models of diverse biological systems and pathways, enhancing transparency and reproducibility in computational biology.

Data compilation

We added all the reactions for the glycolysis pathway and the fatty acid synthesis pathway for palmitate in our SBML file. Firstly, we defined all the units used for the parameters and species. Our model has three compartments, the cytoplasm, mitochondria and extracellular volume. All the species along with their compartments and initial concentrations were defined before the reactions.The parameters Kcat and Km were added for each reaction, taken from BRENDA. Kcat represents the turnover number of an enzyme, indicating the number of substrate molecules converted to product per unit time when an enzyme is fully saturated with substrate. Km is the Michaelis constant, representing the substrate concentration at which the reaction velocity is half of Vmax, which is the maximum reaction velocity.

Time course simulation and Results

Once all the reactions and their rate laws were complete in the SBML file, we validated the SBML file using an online SBML validator[7]. Next, the model was loaded on COPASI, and the concentrations of the metabolites over a time of 15 minutes was observed using its “Time Course” function.

We obtained a graph for the temporal evolution of palmitate according to the reactions and parameters defined in our model over a time of 15 minutes. The final yield we obtained for palmitate was 0.7624 mM, and the graph was seen to have an initial increasing phase, followed by saturation.

Optimisation, Cost analysis and Results

Our primary goal in optimising our computational model was to enhance the yield of palmitate, a key component of our biofuel production, while maintaining consistent initial conditions. This optimization initiative is pivotal in addressing the current challenge of elevated costs associated with biofuels in comparison to conventional fuels. Recognizing this economic hurdle, we devised an algorithm capable of maximising a user-defined product within the model.

The algorithm, detailed in the Software section, will dynamically vary a designated parameter—in our case, the catalytic constant Kcat for the fatty acid synthesis reaction—over a symmetric range to achieve the highest palmitate yield. The result is a refined concentration-versus-time graph, illustrating the enhanced palmitate production attained through the optimisation process.

In parallel, we conducted a meticulous cost analysis, incorporating external factors such as reactor costs, electricity, and labour expenses. Initially performed using our baseline model, this assessment will be then extended to the optimised model. The Entrepreneurship section provides comprehensive insights into the implications of our cost analysis for the project's viability, shedding light on the economic aspects and future prospects of our biofuel initiative.

Temporal activity of the palmitic acid concentration

Approximations and Assumptions

Enzyme Kinetics

Most of the reactions in our model are approximated to single substrate, linear Michaelis-Menten kinetics, following the rate law:

This equation illustrates how the rate of product formation ([P]) is influenced by the substrate concentration ([S]). At low substrate concentrations (when [S] << Km), the rate is approximately proportional to [S]. As [S] increases, the rate approaches Vmax, demonstrating saturation of the enzyme's active sites. The Km value reflects the affinity of the enzyme for its substrate, with lower Km values indicating higher affinity.

Some of our reactions follow the Mass-Action law. For an irreversible reaction following the mass action law with only one substrate (let's call it [A]), the general form of the differential equation is:

This equation simplifies to describe how the rate of an irreversible single-substrate reaction is directly proportional to the concentration of the substrate raised to the power of the reaction order.

Parameter data

The data for the kinetic parameters included in all our reactions are obtained from BRENDA, for Saccharomyces cerevisiae, and the data for some reactions are obtained from Homo sapiens, due to lack of availability of data for any yeast. Our project involves Yarrowia lipolytica, the lipogenic yeast, closely related to Saccharomyces cerevisiae. Due to lack of availability of data for this yeast, we resorted to using data from related organisms. As a further implementation for our project, we intend to use experimental data from Yarrowia lipolytica in order to modify our model into a novel fatty acid synthesis model for the lipogenic yeast, which would be invaluable for research in the biofuel sector.