Towards kinetic modeling of genome-scale metabolic networks without sacrificing stoichiometric, thermodynamic and physiology constraints

Mathematical modeling is an essential tool for a comprehensive understanding of cell metabolism and its interactions with the environmental and process conditions. Recent developments in the construction and analysis of stoichiometric models made it possible to define limits on steady-state metabolic behavior using flux balance analysis. However, detailed information about enzyme kinetics and enzyme regulation is needed to formulate kinetic models that can accurately capture the dynamic metabolic responses. The use of mechanistic enzyme kinetics is a difficult task due to uncertainty in the kinetic properties of enzymes. Therefore, the majority of recent works consider only the mass action kinetics for the reactions in the metabolic networks. In this work, we applied the ORACLE framework and constructed a large-scale, mechanistic kinetic model of optimally grown E. coli. We investigated the complex interplay between stoichiometry, thermodynamics, and kinetics in determining the flexibility and capabilities of metabolism. Our results indicate that enzyme saturation is a necessary consideration in modeling metabolic networks and it extends the feasible ranges of the metabolic fluxes and metabolite concentrations. Our results further suggest that the enzymes in metabolic networks have evolved to function at different saturation states to ensure greater flexibility and robustness of the cellular metabolism

Data, Parameters & Nonlinearities: Development and Applications of Large-scale Dynamic Models of Metabolism

Dynamic nonlinear models of metabolism offer a significant advantage as compared to constraint-based stoichiometric descriptions. However, progress in the development of large-scale nonlinear models has been hindered by both structural and quantitative uncertainties. In particular, the knowledge about kinetic rate laws and their parameters is till today still very limited when compared to the number of stoichiometric reactions known to be present in a large-scale metabolic model. In addition, strategies to systematically identify and implement large-scale dynamic models for metabolism are still lacking. In this contribution, we propose a novel methodology for development of dynamic nonlinear models for metabolism. Using the ORACLE (Optimization and Risk Analysis of Complex Living Entities) framework, we integrate thermodynamics and available omics and kinetic data into a large-scale stoichiometric model. The resulting set of log-linear kinetic models is used to compute kinetic parameters of the involved enzymatic reactions such as the maximal velocities and Michaelis constants. These kinetic parameters are in turn used to compute populations of stable, nonlinear, dynamic models sharing the same stable steady-state as the log-linear ones. The computed models offer unprecedented possibilities for system analysis, e.g. to study the responses of metabolism upon large perturbations; to investigate time course evolutions in and around the steady state; and to identify multiple steady-states and their basins of attraction. We illustrate the features of the generated models in the case of optimally grown E. coli, where our analysis of the estimated maximal reaction rates highlights the significance of network thermodynamics in constraining the variability of these quantities

Identification of Feasible Metabolic Fluxes and Metabolite Concentrations using Large-scale Kinetic Models

The constraints imposed during modeling must satisfy biologically representative phenotypes of the studied organism. Simultaneously, identification and analysis of these constraints enhances our understanding of the evolution/operational paradigms of the organism. It was postulated by Varma and Palsson[1] that it is possible to define limits on metabolic behavior using flux balance analysis, but in order to accurately capture the metabolic responses, detailed information about enzyme kinetics and their regulation is needed. Since development of mechanistic kinetic models is a difficult task due to uncertainty in kinetic properties of enzymes, a substantial number of recent works consider only the mass action (MA) term in their model formulation. As kinetics is one of crucial factors in governing the metabolic capabilities of a cell, i.e. realizable metabolic flux and concentration states, considering only the mass action term does not necessarily provide a realistic description of the feasible space of fluxes and concentrations. In this work, using the ORACLE[2] framework, we constructed a large-scale mechanistic kinetic model of optimally grown E. coli that considers the enzyme saturations as observed in biological systems. Using this model, we performed an analysis of the complex interplay between stoichiometry, thermodynamics, and kinetics in determining flexibility and capabilities of metabolic networks. Our analysis indicates that enzyme saturation is an important and necessary consideration in modeling metabolic networks. Extended ranges of feasibility, both in the space of metabolic fluxes and metabolite concentrations, of kinetic models involving the enzyme saturation suggests that the enzymes in metabolic networks have evolved to function at different saturation states so as to ensure higher flexibility and robustness of the cell

Differential Metabolism of Medium-Chain Fatty Acids in Differentiated Human-Induced Pluripotent Stem Cell-Derived Astrocytes.

Medium-chain triglyceride (MCT) ketogenic diets increase ketone bodies, which are believed to act as alternative energy substrates in the injured brain. Octanoic (C8:0) and decanoic (C10:0) acids, which produce ketone bodies through β-oxidation, are used as part of MCT ketogenic diets. Although the ketogenic role of MCT is well-established, it remains unclear how the network metabolism underlying β-oxidation of these medium-chain fatty acids (MCFA) differ. We aim to elucidate basal β-oxidation of these commonly used MCFA at the cellular level. Human-induced pluripotent stem cell-derived (iPSC) astrocytes were incubated with [U-13C]-C8:0 or [U-13C]-C10:0, and the fractional enrichments (FE) of the derivatives were used for metabolic flux analysis. Data indicate higher extracellular concentrations and faster secretion rates of β-hydroxybutyrate (βHB) and acetoacetate (AcAc) with C8:0 than C10:0, and an important contribution from unlabeled substrates. Flux analysis indicates opposite direction of metabolic flux between the MCFA intermediates C6:0 and C8:0, with an important contribution of unlabeled sources to the elongation in the C10:0 condition, suggesting different β-oxidation pathways. Finally, larger intracellular glutathione concentrations and secretions of 3-OH-C10:0 and C6:0 were measured in C10:0-treated astrocytes. These findings reveal MCFA-specific ketogenic properties. Our results provide insights into designing different MCT-based ketogenic diets to target specific health benefits

Computational Modeling and Analysis of Insulin Induced Eukaryotic Translation Initiation

Insulin, the primary hormone regulating the level of glucose in the bloodstream, modulates a variety of cellular and enzymatic processes in normal and diseased cells. Insulin signals are processed by a complex network of biochemical interactions which ultimately induce gene expression programs or other processes such as translation initiation. Surprisingly, despite the wealth of literature on insulin signaling, the relative importance of the components linking insulin with translation initiation remains unclear. We addressed this question by developing and interrogating a family of mathematical models of insulin induced translation initiation. The insulin network was modeled using mass-action kinetics within an ordinary differential equation (ODE) framework. A family of model parameters was estimated, starting from an initial best fit parameter set, using 24 experimental data sets taken from literature. The residual between model simulations and each of the experimental constraints were simultaneously minimized using multiobjective optimization. Interrogation of the model population, using sensitivity and robustness analysis, identified an insulin-dependent switch that controlled translation initiation. Our analysis suggested that without insulin, a balance between the pro-initiation activity of the GTP-binding protein Rheb and anti-initiation activity of PTEN controlled basal initiation. On the other hand, in the presence of insulin a combination of PI3K and Rheb activity controlled inducible initiation, where PI3K was only critical in the presence of insulin. Other well known regulatory mechanisms governing insulin action, for example IRS-1 negative feedback, modulated the relative importance of PI3K and Rheb but did not fundamentally change the signal flow

Modeling The Unfolded Protein Response And Its Connection To Cellular Stress

The Endoplasmic Reticulum (ER) plays a pivotal role in cellular functioning. Proteins undergo modifications in the ER. The cell monitors protein folding and has a quality control. The incorrectly folded proteins are tagged for degradation or sent back to the refolding cycle. Upon accumulation of incorrectly folded proteins in the cell, due to any alterations in the cellular homeostasis, cellular stress is created. This leads to a cascade of events, termed the Unfolded Protein Response (UPR). Cellular stresses range from hypoxia, nutrient deprivation to external stimulus such as heat or radiation. It is vital to understand the complexities involved in this system, to predict the role of UPR in cancer, diabetes and other diseases associated with misfolded proteins. In the current study, we have created an initial model of the UPR. This model is based on mass action kinetics. We solve this model deterministically, to look at the concentration profiles of the proteins and protein complexes, upon UPR induction. This kind of study, along with the mathematical and statistical tools developed in the Varner lab provide tools to identify fragile elements of the UPR network. These fragile components of the system become targets for intervention by drugs and inform the development of experiments to discover new pharmaceutical treatments

Systems Analysis Of Core Architectures Regulating Cellular Responses Under Stress In Eukaryotes

Eukaryotes have developed evolutionarily conserved mechanisms to respond to diverse ranges of internal and external perturbations e.g., changes in oxygen/nutrient levels, temperature oscillations, protein folding load, viral/bacterial attacks and graft implants. Alterations and malfunctions within these core regulatory architectures provide a peek into shifts towards cancerous and disease states. These core architectures include signal integrators and actuators, positive and negative feedback mechanisms which mediate the effectiveness and ultimate outcome of the response. Contemporary modeling approaches in the era of genomics revolution and high-throughput technology presents a unique systems level insight into many areas of biology, ecology, developmental biology and immunology. In my research, I have used a combination of bottom-up analysis (signaling networks based analysis) and a top-down analysis (using microarray/high throughput experimental data) to investigate stress responses in eukaryotes. Using the bottom-up analysis scheme, we assembled a series of molecular modules describing different aspects of the cellular response to stress. Some of these modules included the Unfolded Protein Response (UPR), Hypoxic Response (HR) and Tumor Angiogenesis, Epithelial to Mesenchymal Transition (EMT), Translation Initiation and Renal Allograft Failure. These modules were modeled using mass action kinetics with an Ordinary Differential Equation (ODE) based framework to investigate the internal regula- tory cores. For example, in UPR we identified the differential negative feedback of activating transcription factor 4 (ATF4) as the key in the adaptation phase via regulation of binding immunoglobulin protein (BiP). Similarly in HR, we identified the role of activator protein 1 (AP1) in mediating the autocrine response via vascular endothelial growth factor (VEGF) and interleukin 8 (IL8) signaling modules. Model generation was done using UNIVERSAL, an in house software freely available at google code. We addressed issues pertaining to uncertainties within these models by developing POETs, a multi-objective optimizing algorithm which allowed us to train our models with experimental data from the literature. POETs presented us with an advantage by generating an ensemble of models consistent with experimental data. The diversity within these ensembles were used to study different operational paradigms within these modules. For example, in EMT we identified the differential modes of crosstalk between mitogen-activated protein kinase (MAPK) and SMADs in mediating the cellular transformation. These different modes of operation suggest insights into different diseases and irregularities in cellular adaptation. We subsequently analyzed these models using parameter independent structural analysis tools like extreme pathways and parameter dependent tools such as fragility and robustness to identify targets relevant to therapeutic interventions. These configurations represent experimentally testable hypothesis and potentially new strategies to manipulate the cellular responses. At an intermediate (length-scale) level, we developed multiscale modeling strategies by inductively extrapolating the consequences of cell signaling to tissue/organ function. We employed this (signaling assisted multiscale modeling (SAMM)) strategy to investigate tumor growth and angiogenesis. Using a top down strategy, we used microarray datasets to investigate cellular signaling, identify malfunctions and create predictive models to infer patient outcome in case of hypoxia induced tumor growth and angiogenesis

Exploring Valine Metabolism in Astrocytic and Liver Cells: Lesson from Clinical Observation in TBI Patients for Nutritional Intervention

The utilization of alternative energy substrates to glucose could be beneficial in traumatic brain injury (TBI). Recent clinical data obtained in TBI patients reported valine, β-hydroxyisobutyrate (ibHB) and 2-ketoisovaleric acid (2-KIV) as three of the main predictors of TBI outcome. In particular, higher levels of ibHB, 2-KIV, and valine in cerebral microdialysis (CMD) were associated with better clinical outcome. In this study, we investigate the correlations between circulating and CMD levels of these metabolites. We hypothesized that the liver can metabolize valine and provide a significant amount of intermediate metabolites, which can be further metabolized in the brain. We aimed to assess the metabolism of valine in human-induced pluripotent stem cell (iPSC)-derived astrocytes and HepG2 cells using 13C-labeled substrate to investigate potential avenues for increasing the levels of downstream metabolites of valine via valine supplementation. We observed that 94 ± 12% and 84 ± 16% of ibHB, and 94 ± 12% and 87 ± 15% of 2-KIV, in the medium of HepG2 cells and in iPSC-derived astrocytes, respectively, came directly from valine. Overall, these findings suggest that both ibHB and 2-KIV are produced from valine to a large extent in both cell types, which could be of interest in the design of optimal nutritional interventions aiming at stimulating valine metabolism

Investigation of network flexibility and identification of reaction essentiality of human host cells infected by S. flexneri

We have developed a consistently reduced core model for HeLa based on the genome scale human reconstruction. This reduced model includes the central carbon pathways, the electron transport chains, all the necessary transport and exchange reactions and the essential compartments for HeLa. In addition, the pyruvate/glutamate metabolism and the purine/pyrimidine catabolism for HeLa were incorporated based on available experimental observations. The core HeLa model consists of 273 reactions and 202 metabolites. The tFBA analysis of the resulting model, where we integrated experimental information about the metabolites concentrations and fluxes, allowed us to: (i) identify the thermodynamically feasible direction profiles of the system; (ii) classify the reactions according to how close (or far) they operate from their thermodynamic equilibrium. This analysis also shed light on the flexibility of the network in terms of the feasible ranges of the metabolite levels and Gibbs energies of reactions. Furthermore, we used the ORACLE (Optimization and Risk Analysis of Complex Living Entities) framework to develop sets of mechanistic kinetics models. ORACLE allows us, despite incomplete information about kinetic properties of the network, to integrate thermodynamics, and available omics and kinetic data. Using these models, we identified and ranked the metabolites according to their impact on the HeLa physiology when drained from the system. This analysis allowed us to uncover the ‘weak’ and ‘strong’ links in host-pathogen interactions. The identification of such links of the metabolic interface between host and pathogen networks could potentially be validated in vitro and can help us understand the mechanisms of infection and provide insights in the fight against intracellular pathogens

Large-scale dynamic models of metabolic networks

The development and application of methods quantifying the metabolic fluxes such as Flux Balance Analysis (FBA) has been one of the driving forces behind the successful growth of metabolic engineering. However, applications of FBA have been limited by the fact that FBA does not take into account kinetic properties of the network and therefore it cannot be used to identify rate-limiting steps and comprehend time course evolutions of the system. Dynamic mathematical descriptions of the metabolism offer a large advantage compared to constraint-based stoichiometric models, but unfortunately their development comes with inevitable difficulties due to: (i) structural uncertainties, such as incomplete knowledge about stoichiometry or about kinetic laws of the enzymes, and (ii) quantitative/parametric uncertainties such as lack of knowledge concerning kinetic parameters. Recent developments and vast resources of curated genome scale metabolic networks address to a great extent the issues around stoichiometric uncertainty. However, the knowledge about kinetic rate laws and in particular their parameters is to these days still limited. In this contribution, starting from large-scale stoichiometric models we use the ORACLE (Optimization and Risk Analysis of Complex Living Entities) framework that integrates available information in a set of stable log-linear kinetic models sharing the same steady state. These models are used to compute kinetic parameters of the enzymatic mechanisms in the metabolic network, e.g. the maximal velocities, Vmax, and Michaelis constants, Km. Using these kinetic parameters, we systematically develop populations of stable dynamic models having the same steady-state as the log-linear ones. The estimated parameters are comparable to the experimental information as seen in BRENDA and other databases. These non-linear estimations about the stable state can be therefore used to analyze properties of the system upon large perturbations and investigate time course evolutions in and around this steady state. We demonstrate the capabilities of the proposed approach by building a dynamic E. coli core model that includes ca. 200 metabolites and more than 400 reactions. We discuss the strengths and limitations of this approach and possible avenues for development