Author response: Dynamics of co-substrate pools can constrain and regulate metabolic fluxes

Abstract

Article Figures and data Abstract Editor's evaluation eLife digest Introduction Results Discussion Materials and methods Appendix 1 Appendix 2 Appendix 3 Appendix 4 Appendix 5 Appendix 6 Appendix 7 Appendix 8 Appendix 9 Appendix 10 Appendix 11 Data availability References Decision letter Author response Article and author information Metrics Abstract Cycling of co-substrates, whereby a metabolite is converted among alternate forms via different reactions, is ubiquitous in metabolism. Several cycled co-substrates are well known as energy and electron carriers (e.g. ATP and NAD(P)H), but there are also other metabolites that act as cycled co-substrates in different parts of central metabolism. Here, we develop a mathematical framework to analyse the effect of co-substrate cycling on metabolic flux. In the cases of a single reaction and linear pathways, we find that co-substrate cycling imposes an additional flux limit on a reaction, distinct to the limit imposed by the kinetics of the primary enzyme catalysing that reaction. Using analytical methods, we show that this additional limit is a function of the total pool size and turnover rate of the cycled co-substrate. Expanding from this insight and using simulations, we show that regulation of these two parameters can allow regulation of flux dynamics in branched and coupled pathways. To support these theoretical insights, we analysed existing flux measurements and enzyme levels from the central carbon metabolism and identified several reactions that could be limited by the dynamics of co-substrate cycling. We discuss how the limitations imposed by co-substrate cycling provide experimentally testable hypotheses on specific metabolic phenotypes. We conclude that measuring and controlling co-substrate dynamics is crucial for understanding and engineering metabolic fluxes in cells. Editor's evaluation This manuscript presents an important mathematical analysis of metabolic "co-substrates" and how their cycling can affect metabolic fluxes. Through mathematical analysis of simple network motifs, it shows the impact of co-substrate cycling on constraining metabolic fluxes. The combination of mathematical modeling and comparisons with existing data from previous studies offers convincing support for the potential biological relevance of co-substrate cycling. The work will be of interest to researchers who study microbial metabolism and metabolic engineering. https://doi.org/10.7554/eLife.84379.sa0 Decision letter Reviews on Sciety eLife's review process eLife digest Metabolism powers individual cells and ultimately the body. It comprises a sequence of chemical reactions that cells use to break down substances and generate energy. These reactions are catalyzed by enzymes, which are proteins that speed up the rate of the reaction. Many reactions also involve co-substrates, which are themselves transformed by individual reactions but are eventually converted back into their original form in a series of steps. This process is known as co-substrate cycling. Scientists have long been interested in understanding what controls the rate at which metabolic reactions and metabolic pathways convert a substance into a final product. This is a difficult subject to study because of the complexity of the metabolic pathways, with their branched, linear or coupled structures. In the past, researchers have looked at the influence of enzymes on the rate of a metabolic pathway, but less has been known about the effect of co-substrate cycling. To find out more, West, Delattre et al. developed a series of mathematical models to describe different types of metabolic pathways in terms of the number of metabolites that enter and leave it, including the influence of co-substrates. They found that co-substrate cycling, when involved in a metabolic reaction, limits the speed with which the reaction happens. This is distinct from the limit that enzymes impose on the speed of the reaction. It depends on the total amount of co-substrates in the cell: changing the number of co-substrates in the cell influences the speed at which the metabolic reaction takes place. This study has increased our understanding of how metabolic pathways work, and what controls the speed at which reactions take place. It opens up a new potential method for explaining how cells control metabolic reaction rates and how metabolic substrates can be directed across different pathways. This research is likely to inspire future research into the influence of co-substrates in different cell types and conditions. Introduction Dynamics of cell metabolism directly influences individual and population-level cellular responses. Examples include metabolic oscillations underpinning the cell cycle (Papagiannakis et al., 2017; Murray et al., 2007) and metabolic shifts from respiration to fermentation, as observed in cancer phenotypes (Warburg, 1956; Diaz-Ruiz et al., 2009; Carmona-Fontaine et al., 2013) and cell-to-cell cross-feeding (Ponomarova et al., 2017; Campbell et al., 2015; Großkopf et al., 2016). Predicting or conceptualising these physiological responses using dynamical models is difficult due to the large size and high connectivity of cellular metabolism. Despite this complexity, however, it is possible that cellular metabolism features certain ‘design principles’ that determine the overall dynamics. There is ongoing interest in finding such simplifying principles. A key concept for understanding the dynamics of any metabolic system is that of ‘reaction flux’, which is a measure of the rate of biochemical conversion in a given reaction. To identify possible limitations on reaction fluxes, early studies focused on linear pathways involving ATP production and studied their dynamics under the optimality assumption of maximisation of overall pathway flux under limited enzyme levels available to the pathway (Heinrich et al., 1991). The resulting theory predicted a trade-off between pathway flux vs. yield (i.e. rate of ATP generation vs. amount of ATP generated per metabolite consumed by the pathway) in linear pathways (Heinrich and Hoffmann, 1991). This theory is subsequently used to explain the emergence of different metabolic phenotypes (Pfeiffer et al., 2001). In related studies, models pertaining to flux optimisation and enzyme levels being a key limitation are used to explain the structure of different metabolic pathways (Flamholz et al., 2013), and the metabolic shifting from respiration to fermentative pathways under increasing glycolysis rates (Großkopf et al., 2016; Basan et al., 2015; Majewski and Domach, 1990). There are, however, increasing number of studies suggesting that enzyme levels alone might not be sufficient to explain observed flux levels. For example, it was shown that the maximal value of the apparent activities (ka⁢p⁢pm⁢a⁢x) of an enzyme, derived using measured enzyme levels and fluxes under different conditions, was a good estimate for the specific activity of that enzyme in vitro (kc⁢a⁢t) (Davidi et al., 2016). However, individual estimates from each condition (i.e. individual ka⁢p⁢p values) were commonly lower than the specific activity – suggesting that the flux is limited by something other than enzyme levels under those conditions. Other studies have shown that metabolic flux changes, caused by perturbations in media conditions, are not explained solely by changes in expression levels of enzymes (Chubukov et al., 2013; Gerosa et al., 2015). Another conceptual framework emphasized the importance of cyclic reaction motifs, particularly those involving so-called co-substrate pairs, such as ATP / ADP or NAD(P)H / NAD(P)+, as a key to understanding metabolic system dynamics (Reich and Sel‘kov, 1981). This framework is linked to the idea of considering the supply and demand structures around specific metabolites as regulatory blocks within metabolism (Hofmeyr and Cornish-Bowden, 2000). For example, the total pool of ATP and its derivates (the ‘energy charge’) is suggested as a key determinant of physiological cell states (Atkinson, 1968). Inspired by these ideas, theoretical studies have shown that metabolic systems featuring metabolite cycling together with allosteric regulation can introduce switch-like and bistable dynamics (Okamoto and Hayashi, 1983; Hervagault and Cimino, 1989), and that metabolite cycling motifs introduce total co-substrate level as an additional control element in metabolic control analysis (Hofmeyr et al., 1986; Sauro, 1994). Specific analyses of ATP cycling in the glycolysis pathway, sometimes referred to as a ‘turbo-design’, and metabolite cycling with autocatalysis, as seen for example in the glyoxylate cycle, have shown that these features constrain pathway fluxes (Koebmann et al., 2002; Teusink et al., 1998; van Heerden et al., 2014; Hatakeyama and Furusawa, 2017; Barenholz et al., 2017; Kurata, 2019). Taken together, these studies indicate that metabolite cycling, in general, and co-substrate cycling specifically, could provide a key ‘design feature’ in cell metabolism, imposing certain constraints or dynamical properties to it. Towards better understanding the role of co-substrate cycling in cell metabolism dynamics, we undertook here an analytical and simulation-based mathematical study together with analyses of measured fluxes. We created models of enzymatic reaction systems featuring co-substrate cycling, abstracted from real metabolic systems such as glycolysis, nitrogen-assimilation, and central carbon metabolism. We found that co-substrate cycling introduces a fundamental constraint on reaction flux. In the case of single reaction and short linear pathways, we were able to derive a mathematical expression of the constraint, showing that it relates to the pool size and turnover rate of the co-substrate. Analysing measured fluxes, we find that several of the co-substrate featuring reactions in central carbon metabolism carry lower fluxes than expected from the kinetics of their primary enzymes, suggesting that these reactions might be limited by co-substrate cycling. In addition to its possible constraining role, we show that co-substrate cycling can also act as a regulatory element, where control of co-substrate pool size can allow control of flux dynamics across connected or branching pathways. Together, these findings show that co-substrate cycling can act both as a constraint and a regulatory element in cellular metabolism. The resulting theory provides testable hypotheses on how to manipulate metabolic fluxes and cell physiology through the control of co-substrate pool sizes and turnover dynamics, and can be expanded to explain dynamic measurements of metabolite concentrations in different perturbation experiments. Results Co-substrate cycling represents a ubiquitous motif in metabolism with co-substrate pools acting as ‘conserved moieties’ Certain pairs of metabolites can be interconverted via different reactions in the cell, thereby resulting in their ‘cycling’. This cycling creates interconnections within metabolism, spanning either multiple reactions in a single, linear pathway, or multiple pathways that are independent or are branching from common metabolites. For example, in glycolysis, ATP is consumed in reactions mediated by the enzymes glucose hexokinase and phosphofructokinase, and is produced by the downstream reactions mediated by phosphoglycerate and pyruvate kinase (Appendix 1—figure 1A). In the nitrogen assimilation pathway, the NAD+ / NADH pair is cycled by the enzymes glutamine oxoglutarate aminotransferase and glutamate dehydrogenase (Appendix Dynamics of co-substrate pools can constrain and regulate metabolic fluxes - Appendix 1—figure 1B). Many other cycling motifs can be identified, involving either metabolites from the central carbon metabolism or metabolites that are usually referred to as co-substrates. Examples for the latter include NADPH, FADH2, GTP, and Acetyl-CoA and their corresponding alternate forms, while examples for the former include the tetrahydrofolate (THF) / 5,10-Methylene-THF and glutamate / α-ketoglutarate (akg) pairs involved in one-carbon transfer and in amino acid biosynthesis pathways, respectively (Appendix 1—figure 1C and D). For some of these metabolites, their cycling can connect many reactions in the metabolic network. Taking ATP (NADH) as an example, there are 265 (118) and 833 (601) reactions linked to the cycling of this metabolite in the genome-scale metabolic models of Escherichia coli and human respectively models iJO1366 (Orth et al., 2011) and Recon3d (Brunk et al., 2018). We notice here that many of the co-substrate involving cycling reactions can be abstracted as a simplified motif as shown in (Figure 1A). This abstract representation highlights the fact that the total pool-size involving all the different forms of a cycled metabolite can become a conserved quantity. This would be the case even when we consider biosynthesis or environmental uptake of co-substrates, as the total concentration of a cycled metabolite across its different forms at steady state would then be given by a constant defined by the ratio of the influx and outflux rates (see Appendices 2 and 3). In other words, the cycled metabolite would become a ‘conserved moiety’ for the rest of the metabolic system and can have a constant ‘pool size’. Supporting this, temporal measurement of specific co-substrate pool sizes shows that ATP and GTP pools are constant under stable metabolic conditions, but can rapidly change in response to external perturbations, possibly through inter-conversions among pools rather than through biosynthesis (Walther et al., 2010). Figure 1 Download asset Open asset Motif, time-series and threshold in a single co-substrate involving reaction. (A) Cartoon representation of a single irreversible reaction with co-substrate cycling (see Appendices for other reaction schemes). The co-substrate is considered to have two forms A0 and A1. (B) Concentrations of the metabolites M0 (red) and M1 (green), and the A0/A1 ratio (blue) are shown as a function of time. At t=500, the parameters are switched from the white dot in panel (C) (where a steady state exists) to the black dot (where we see continual build-up of M0 and decline of A0 without steady state). (C & D) Heatmap of the steady state concentration of M0 as a function of the total co-substrate pool size (Atot) and inflow flux (kin). White area shows the region where there is no steady state. On both panels, the dashed line indicates the limitation from the primary enzyme, kin<Vmax,E0, and the solid line indicates the limitation from co-substrate cycling, kin<AtotVmax,Ea/(KM,Ea+Atot). In panel (C), there is a range of Atot values for which the first limitation is more severe than the second. In contrast, in panel (D), the second limitation is always more severe than the first. In (B & C) the parameters used for the primary enzyme (for the reaction converting M0 into M1) are picked from within a physiological range (see Supplementary file 1) and are set to: Etot=0.01 mM, kcat=100/s, KM,E0=KM,Ea=50μM, while kout is set to 0.1/s. The Etot and kcat for the co-substrate cycling enzyme are 1.2 times those for the primary enzyme. In panel (D) the parameters are the same except for the Etot and kcat of the co-substrate cycling enzyme, which are set to 0.7 times those for the primary enzyme. Co-substrate cycling introduces a limitation on reaction flux To explore the effect of co-substrate cycling on pathway fluxes, we first consider a didactic case of a single reaction. This reaction converts an arbitrary metabolite M0 to M1 and involves co-substrate cycling (Figure 1A). For co-substrate cycling, we consider additional ‘background’ enzymatic reactions that are independent of M0 and can also convert the co-substrate (denoted Ea on Figure 1A). We use either irreversible or reversible enzyme dynamics to build an ordinary differential equation (ODE) kinetic model for this reaction system and solve for its steady states analytically (see Methods and Appendix 3). In the case of using irreversible enzyme kinetics, we obtain that the steady state concentration of the two metabolites, M0 and M1 (denoted as m0 and m1) are given by: (1) m0=α kinKM,E0(Vmax,E0−kin)(Vmax,EaAtot−kin(KM,Ea+Atot))m1=kinkout where kin and kout denote the rate of in-flux of M0, and out-flux of M1, either in-and-out of the cell or from other pathways, and Atot denotes the total pool size of the cycled metabolite (with the different forms of the cycled metabolite indicated as A0 and A1 in Figure 1A). The parameters Vmax,E0 and Vmax,Ea are the maximal rates (i.e. Vmax=kcat⁢Etot) for the enzymes catalysing the conversion of A0 and M0 into A1 and M1 (enzyme E0), and the turnover of A1 into A0 (enzyme Ea), respectively, while the parameters KM,E0 and KM,Ea are the individual or combined Michaelis-Menten coefficients for these enzymes’ substrates (i.e. for A0 and M0 and A1, respectively). The term α is (in this case where all reactions are irreversible) equal to Vmax,Ea-kin, and in general is a positive expression comprising kin, and the Michaelis-Menten coefficients and the Vmax parameters of the background enzymes in the model (see Appendix 3, Equations 7; 9; 11). The steady states for the model with all enzymatic conversions being reversible, and for a model with degradation and synthesis of A0 and A1, are given in Appendix 3. The steady state solutions of these alternative models are structurally akin to (1), and do not alter the qualitative conclusions we make in what follows. A key property of (1) is that it contains terms in the denominator that involve a subtraction. The presence of these terms introduces a limit on the parameter values for the system to attain a positive steady state. Specifically, we obtain the following conditions for positive steady states to exist: (2) kin<Vmax,E0andkin<AtotVmax,EaKM,Ea+Atot. Additionally, the ‘shape’ of (1) indicates a ‘threshold effect’ on the steady state value of m0, where it would rise towards infinity as kin increases towards the lower one among the limits given in (2) (see Figure 1B). Why does (1) show this specific form, leading to these limits? We find that this is a direct consequence of the steady state condition, where metabolite production and consumption rates need to be the same at steady state. In the case of co-substrate cycling, the production rate of M0 is given by kin, while its consumption rate is a function of the Vmax,E0 and the concentration of A0. In turn, the concentration of A0 is determined by its re-generation rate (which is a function of KM,Ea and Vmax,Ea) and the pool size (Atot). This explains the inequalities given in (2) and shows that a cycled co-substrate creates the same type of limitation (mathematically speaking) on the flux of a reaction it is involved in, as that imposed by the enzyme catalysing that reaction (E0 in this example) (see Figure 1C & D). We also show that considering the system shown in Figure 1A as an enzymatic reaction without co-substrate cycling leads to only the constraint kin<Vmax,E0, while considering it as a non-enzymatic reaction with co-substrate cycling only, leads to only the constraint kin<AtotVmax,Ea/(KM,Ea+Atot) becoming the sole limitation on the system (see Appendix 3). In other words, the two limitations act independently. To conclude this section, we re-iterate its main result. The flux of a reaction involving co-substrate cycling is limited either by the kinetics of the primary enzyme mediating that reaction, or by the turnover rate of the co-substrate. The latter is determined by the co-substrate pool size and the kinetics of the enzyme(s) mediating its turnover. Co-substrate cycling causes a flux limit on linear metabolic pathways We next considered a generalised, linear pathway model with n+1 metabolites and arbitrary locations of reactions for co-substrate cycling, for example as seen in upper glycolysis (Appendix 1—figure 1A). In this model, we only consider intra-pathway metabolite cycling, i.e. the co-substrate is consumed and re-generated solely by the reactions of the pathway. Here, we show results for this model with 5 metabolites as an illustration (Figure 2A), while the general case is presented in Appendix 4. Figure 2 Download asset Open asset Motif, time-series and thresholds for the linear pathway model with n=4. (A) Cartoon representation of a chain of reversible reactions with co-substrate cycling occurring solely inter-pathway. The co-substrate is considered to have two forms A0 and A1. (B) Heatmap of the steady state concentration of M0 as a function of the total metabolite pool size (Atot) and inflow rate constant (kin). White area shows the region where there is no steady state. The dashed and solid lines indicate the limitations arising from primary enzyme (E1 in this case) and co-substrate cycling, respectively, as in Figure 1. (C) Concentrations of M0-4, and A0/A1 ratio as a function of time (with colours as indicated in the inset). At t=1000 s, the parameters are switched from the white dot in panel (B) (where a steady state exists) to the black dot (where we see build-up of all substrates that are produced before the first co-substrate cycling reaction, and continued decline of A0). The parameters used are picked from within a physiological range (see Supplementary file 1) and are set to: Etot=0.01 mM, kcat=100/s, KM=50μm, for all reactions, and kout=0.1/s. We find the same kind of threshold dynamics as in the single reaction case. When kin is above a threshold value, the metabolite M0 accumulates towards infinity and the system does not have a steady state (Figure 2B). A numerical analysis, as well as our analytical solution, reveals that the accumulation of metabolites applies to all metabolites upstream of the first reaction with co-substrate cycling (Figure 2C and Appendix 4). Additionally, metabolites downstream of the cycling reaction accumulate to a steady state level that does not depend on kin (Figure 2C and Appendix 3—figure 1). In other words, pathway output cannot be increased further by increasing kin beyond the threshold. Finally, as kin increases, the cycled metabolite pool shifts towards one form and the ratio of the two forms approaches zero (Figure 2C). An analytical expression for the threshold for kin, like shown in (2), could not be derived for linear pathways with n>3, but our analytical study indicates that (i) the threshold is always linked to Atot and enzyme kinetic parameters, and (ii) the concentration of all metabolites upstream (downstream) to the reaction coupled to metabolite cycling will accumulate towards infinity (a fixed value) as kin approaches the threshold (see Appendix 4). In Figure 2, we illustrate these dynamics with simulations for a system with n=4. We also considered several variants of this generalised linear pathway model, corresponding to biologically relevant cases as shown in Appendix 1—figure 1. These included (i) intra-pathway cycling of two different metabolites, as seen with ATP and NADH in combined upper glycolysis and fermentation pathways (Appendix 5), (ii) different stoichiometries for consumption and re-generation reactions of the cycled metabolite, as seen in upper glycolysis (Appendix 6), and (iii) cycling of one metabolite interlinked with that of another, as seen in nitrogen assimilation (Appendix 7). The results in the Appendices confirm that all these cases display similar threshold dynamics, where the threshold point is a function of the co-substrate pool size and the enzyme kinetics. Cycled metabolite related limit could be relevant for specific reactions from central metabolism Based on flux values that are either experimentally measured or predicted by flux balance analysis (FBA), many reactions from the central carbon metabolism of the model organism Escherichia coli are shown to have lower flux than expected from the kinetics of their immediate enzymes (i.e. Vmax) (Davidi et al., 2016). This finding is based on calculating Vmax from in vitro measured kcat values of specific enzymes and their in vivo levels based on proteomics studies in E. coli (see Materials and methods). The flux and enzyme concentration data were from other studies which measured them during the exponential phase in E. coli growing on minimal media supplemented with various carbon sources (Schmidt et al., 2016; Gerosa et al., 2015). If we consider measured fluxes for each reaction as a proxy for kin (notice that these two would be equal at steady state), we can conclude from the fact that there were no observed substrate accumulation in these reactions, as an indication for the analysed reactions carrying fluxes below the first limit identified above in (2). There could be several explanations for this observation of measured fluxes being lower than the limit set by measured enzyme kinetics and level. One simple explanation could be that there is a discrepancy between in vitro measured enzyme kinetics and in vivo realised ones. Alternatively, this discrepancy can be low, but the lower flux could be arising because there are additional limiting factors other than the enzymes mediating the main reaction. Among such additional limiting factors, substrate limitation and thermodynamic effects are shown to partially explain observed lower fluxes in some reactions (Davidi et al., 2016; see also below results). Here, we highlight that the presented theory shows that an additional possible limitation could be the co-substrate pool size and turnover dynamics. To explore this possibility, we re-analysed the flux values compiled previously (Davidi et al., 2016; Gerosa et al., 2015) and focused solely on reactions that are linked to ATP, NADH, or NADPH pools (see Materials and methods and Supplementary file 1). The resulting dataset contained fluxes, substrate concentrations, and enzyme levels for 45 different reactions determined under 7 different conditions along with turnover numbers and kinetic constants of the corresponding enzymes. In total, we gathered 49 combinations of enzyme-flux-kcat values with full experimental data and 259 combinations with only FBA-predicted flux values. We compared the flux values that would be expected from the primary enzyme limit identified above, under all conditions analysed (Figure 3A), and in addition checked whether the saturation effect of the primary substrate could explain the difference (Appendix 8—figure 1). We found that in both cases, about 80% of these reactions carry flux lower than what is expected from enzyme kinetics (Appendix 8—figure 2), suggesting that the limits imposed by co-factor dynamics might be constraining the flux further. The low number of the cases where the flux exceeds the limit might be due to uncertainties in measurement of flux, enzyme or substrate level. Figure 3 Download asset Open asset Measured and FBA-predicted flux values are typically lower than the calculated primary enzyme threshold. (A) Measured and FBA-predicted flux values (from Davidi et al., 2016; Gerosa et al., 2015) plotted against the calculated primary enzyme kinetic threshold (first part of eq. (1)). Notice that there are 7 points for each reaction, corresponding to the different experimental conditions under which measurements or FBA modelling was done (see Supplementary file 1 for data, along with reaction names and metabolites involved). (B–D) Measured flux values under different experimental conditions (from Gerosa et al., 2015) for select reactions plotted against the corresponding co-substrate pool size. Panels B to D show reactions for phosphoglycerate kinase (PGK), malate dehydrogenase (MDH), and glucose-6-phosphate dehydrogenase (G6PDH). Each point on these panels is a separate flux measurement under a different environmental condition, where the co-substrate pool size is also measured. Error bars represent standard deviations of flux and metabolite measurements as they appear in the dataset from Gerosa et al. Point colours represent co-substrate type and are as shown in the legend to panel A. Lines show the best linear fit with the corresponding normalised RMSE shown in the panel title. It is also possible that observed lower fluxes are due to thermodynamic limitations. This is very difficult to analyse without more data, as calculating reaction thermodynamics requires knowledge of concentrations for all substrates and products, as well as enzyme Michaelis-Menten constants in both forward and backward directions. This information is currently not available except for few of the reactions among the ones we analysed. Nevertheless, to give as much insight as possible on the thermodynamic effect, we analysed the physiological Gibbs free energy (the Δr⁢G′ is calculated assuming that all reactants are at 1 mM and pH = 7) against the normalized flux – v/(E0⋅kc⁢a⁢t) (Appendix 8—figure 3). This shows that although in few cases, such as malate dehydrogenase (MDH), the normalised flux seems to be greatly reduced by the thermodynamic barrier, the general picture is that there is little correlation between reaction flux and thermodynamics. We have also checked the relation between fluxes and co-substrate pool sizes. Co-substrate pool sizes do change between different conditions, and we note that such changes cannot be due to flux changes in co-substrate utilising reactions. But, on the other hand, changes in pool size can affect flux in those reactions, where co-substrate dynamics is limiting (as predicted by the theory). For both measured and FBA-predicted fluxes, we find that several reactions show significant correlation between flux and co-substrate pool size (see Figure 3B–D, see also Appendix 8—table 1 and Appendix 8—figure 4).