Elementary Vectors and Conformal Sums in Polyhedral Geometry and their Relevance for Metabolic Pathway Analysis

Connectivity of Parameter Regions of Multistationarity for Multisite Phosphorylation Networks

The parameter region of multistationarity of a reaction network contains all the parameters for which the associated dynamical system exhibits multiple steady states. Describing this region is challenging and remains an active area of research. In this paper, we concentrate on two biologically relevant families of reaction networks that model multisite phosphorylation and dephosphorylation of a substrate at n sites. For small values of n, it had previously been shown that the parameter region of multistationarity is connected. Here, we extend these results and provide a proof that applies to all values of n. Our techniques are based on the study of the critical polynomial associated with these reaction networks together with polyhedral geometric conditions of the signed support of this polynomial.

A hierarchy of metabolite exchanges in metabolic models of microbial species and communities

The metabolic network of an organism can be analyzed as a constraint-based model. This analysis can be biased, optimizing an objective such as growth rate, or unbiased, aiming to describe the full feasible space of metabolic fluxes through pathway analysis or random flux sampling. In particular, pathway analysis can decompose the flux space into fundamental and formally defined metabolic pathways. Unbiased methods scale poorly with network size due to combinatorial explosion, but a promising approach to improve scalability is to focus on metabolic subnetworks, e.g., cells' metabolite exchanges with each other and the environment, rather than the full metabolic networks. Here, we applied pathway enumeration and flux sampling to metabolite exchanges in microbial species and a microbial community, using models ranging from central carbon metabolism to genome-scale and focusing on pathway definitions that allow direct targeting of subnetworks such as metabolite exchanges (elementary conversion modes, elementary flux patterns, and minimal pathways). Enumerating growth-supporting metabolite exchanges, we found that metabolite exchanges from different pathway definitions were related through a hierarchy, and we show that this hierarchical relationship between pathways holds for metabolic networks and subnetworks more generally. Metabolite exchange frequencies, defined as the fraction of pathways in which each metabolite was exchanged, were similar across pathway definitions, with a few specific exchanges explaining large differences in pathway counts. This indicates that biological interpretation of predicted metabolite exchanges is robust to the choice of pathway definition, and it suggests strategies for more scalable pathway analysis. Our results also signal wider biological implications, facilitating detailed and interpretable analysis of metabolite exchanges and other subnetworks in fields such as metabolic engineering and synthetic biology.

A mathematically rigorous algorithm to define, compute and assess relevance of the probable dissociation constants in characterizing a biochemical network

Metabolism results from enzymatic- and non-enzymatic interactions of several molecules, is easily parameterized with the dissociation constant and occurs via biochemical networks. The dissociation constant is an empirically determined parameter and cannot be used directly to investigate in silico models of biochemical networks. Here, we develop and present an algorithm to define, compute and assess the relevance of the probable dissociation constant for every reaction of a biochemical network. The reactants and reactions of this network are modelled by a stoichiometry number matrix. The algorithm computes the null space and then serially generates subspaces by combinatorially summing the spanning vectors that are non-trivial and unique. This is done until the terms of each row either monotonically diverge or form an alternating sequence whose terms can be partitioned into subsets with almost the same number of oppositely signed terms. For a selected null space-generated subspace the algorithm utilizes several statistical and mathematical descriptors to select and bin terms from each row into distinct outcome-specific subsets. The terms of each subset are summed, mapped to the real-valued open interval [Formula: see text] and used to populate a reaction-specific outcome vector. The p1-norm for this vector is then the probable dissociation constant for this reaction. These steps are continued until every reaction of a modelled network is unambiguously annotated. The assertions presented are complemented by computational studies of a biochemical network for aerobic glycolysis. The fundamental premise of this work is that every row of a null space-generated subspace is a valid reaction and can therefore, be modelled as a reaction-specific sequence vector with a dimension that corresponds to the cardinality of the subspace after excluding all trivial- and redundant-vectors. A major finding of this study is that the row-wise sum or the sum of the terms contained in each reaction-specific sequence vector is mapped unambiguously to a positive real number. This means that the probable dissociation constants, for all reactions, can be directly computed from the stoichiometry number matrix and are suitable indicators of outcome for every reaction of the modelled biochemical network. Additionally, we find that the unambiguous annotation for a biochemical network will require a minimum number of iterations and will determine computational complexity.

ReDirection: an R-package to compute the probable dissociation constant for every reaction of a user-defined biochemical network

Biochemical networks integrate enzyme-mediated substrate conversions with non-enzymatic complex formation and disassembly to accomplish complex biochemical and physiological functions. The choice of parameters and constraints used in most of these studies is numerically motivated and network-specific. Although sound in theory, the outcomes that result depart significantly from the intracellular milieu and are less likely to retain relevance in a clinical setting. There is a need for a computational tool which is biochemically relevant, mathematically rigorous, and unbiased, and can ascribe functionality to and generate potentially testable hypotheses for a user-defined biochemical network. Here, we present "ReDirection," an R-package which computes the probable dissociation constant for every reaction of a biochemical network directly from a null space-generated subspace of the stoichiometry number matrix of the modeled network. "ReDirection" delineates this subspace by excluding all trivial and redundant or duplicate occurrences of non-trivial vectors, combinatorially summing the vectors that remain and verifying that the upper or lower bounds of the sequence of terms formed by each row of this subspace belong to the open real-valued intervals - ∞ , - 1 or 1 , ∞ or whether the number of terms that are differently signed are almost equal. "ReDirection" iterates these steps until these bounds are consistent and unambiguous for all reactions of the modeled biochemical network. Thereafter, "ReDirection" filters the terms from each row of this subspace, bins them to outcome-specific subsets, sums and maps this to an outcome-specific reaction vector, and computes the p1-norm, which is the probable dissociation constant for a reaction. "ReDirection" works on first principles, does not discriminate between enzymatic and non-enzymatic reactions, offers a biochemically relevant and mathematically rigorous environment to explore user-defined biochemical networks under baseline and perturbed conditions, and can be used to address empirically intractable biochemical problems. The utility and relevance of "ReDirection" are highlighted by numerical studies on stoichiometric number models of biochemical networks of galactose metabolism and heme and cholesterol biosynthesis. "ReDirection" is freely available and accessible from the comprehensive R archive network (CRAN) with the URL (https://cran.r-project.org/package=ReDirection).

Modeling the metabolic dynamics at the genome-scale by optimized yield analysis

The hybrid cybernetic model (HCM) approach is a dynamic modeling framework that integrates enzyme synthesis and activity regulation. It has been widely applied in bioreaction engineering, particularly in the simulation of microbial growth in different mixtures of carbon sources. In a HCM, the metabolic network is decomposed into elementary flux modes (EFMs), whereby the network can be reduced into a few pathways by yield analysis. However, applying the HCM approach on conventional genome-scale metabolic models (GEMs) is still a challenge due to the high computational demands. Here, we present a HCM strategy that introduced an optimized yield analysis algorithm (opt-yield-FBA) to simulate metabolic dynamics at the genome-scale without the need for EFMs calculation. The opt-yield-FBA is a flux-balance analysis (FBA) based method that can calculate optimal yield solutions and yield space for GEM. With the opt-yield-FBA algorithm, the HCM strategy can be applied to get the yield spaces and avoid the computational burden of EFMs, and it can therefore be applied for developing dynamic models for genome-scale metabolic networks. Here, we illustrate the strategy by applying the concept to simulate the dynamics of microbial communities.

What makes a reaction network "chemical"?

BACKGROUND

Reaction networks (RNs) comprise a set X of species and a set [Formula: see text] of reactions [Formula: see text], each converting a multiset of educts [Formula: see text] into a multiset [Formula: see text] of products. RNs are equivalent to directed hypergraphs. However, not all RNs necessarily admit a chemical interpretation. Instead, they might contradict fundamental principles of physics such as the conservation of energy and mass or the reversibility of chemical reactions. The consequences of these necessary conditions for the stoichiometric matrix [Formula: see text] have been discussed extensively in the chemical literature. Here, we provide sufficient conditions for [Formula: see text] that guarantee the interpretation of RNs in terms of balanced sum formulas and structural formulas, respectively.

RESULTS

Chemically plausible RNs allow neither a perpetuum mobile, i.e., a "futile cycle" of reactions with non-vanishing energy production, nor the creation or annihilation of mass. Such RNs are said to be thermodynamically sound and conservative. For finite RNs, both conditions can be expressed equivalently as properties of the stoichiometric matrix [Formula: see text]. The first condition is vacuous for reversible networks, but it excludes irreversible futile cycles and-in a stricter sense-futile cycles that even contain an irreversible reaction. The second condition is equivalent to the existence of a strictly positive reaction invariant. It is also sufficient for the existence of a realization in terms of sum formulas, obeying conservation of "atoms". In particular, these realizations can be chosen such that any two species have distinct sum formulas, unless [Formula: see text] implies that they are "obligatory isomers". In terms of structural formulas, every compound is a labeled multigraph, in essence a Lewis formula, and reactions comprise only a rearrangement of bonds such that the total bond order is preserved. In particular, for every conservative RN, there exists a Lewis realization, in which any two compounds are realized by pairwisely distinct multigraphs. Finally, we show that, in general, there are infinitely many realizations for a given conservative RN.

CONCLUSIONS

"Chemical" RNs are directed hypergraphs with a stoichiometric matrix [Formula: see text] whose left kernel contains a strictly positive vector and whose right kernel does not contain a futile cycle involving an irreversible reaction. This simple characterization also provides a concise specification of random models for chemical RNs that additionally constrain [Formula: see text] by rank, sparsity, or distribution of the non-zero entries. Furthermore, it suggests several interesting avenues for future research, in particular, concerning alternative representations of reaction networks and infinite chemical universes.

EFMlrs: a Python package for elementary flux mode enumeration via lexicographic reverse search

Background

Elementary flux mode (EFM) analysis is a well-established, yet computationally challenging approach to characterize metabolic networks. Standard algorithms require huge amounts of memory and lack scalability which limits their application to single servers and consequently limits a comprehensive analysis to medium-scale networks. Recently, Avis et al. developed mplrs—a parallel version of the lexicographic reverse search (lrs) algorithm, which, in principle, enables an EFM analysis on high-performance computing environments (Avis and Jordan. mplrs: a scalable parallel vertex/facet enumeration code. arXiv:1511.06487, 2017). Here we test its applicability for EFM enumeration.

Results

We developed EFMlrs, a Python package that gives users access to the enumeration capabilities of mplrs. EFMlrs uses COBRApy to process metabolic models from sbml files, performs loss-free compressions of the stoichiometric matrix, and generates suitable inputs for mplrs as well as efmtool, providing support not only for our proposed new method for EFM enumeration but also for already established tools. By leveraging COBRApy, EFMlrs also allows the application of additional reaction boundaries and seamlessly integrates into existing workflows.

Conclusion

We show that due to mplrs’s properties, the algorithm is perfectly suited for high-performance computing (HPC) and thus offers new possibilities for the unbiased analysis of substantially larger metabolic models via EFM analyses. EFMlrs is an open-source program that comes together with a designated workflow and can be easily installed via pip.

Supplementary Information

The online version contains supplementary material available at 10.1186/s12859-021-04417-9.

Understanding FBA Solutions under Multiple Nutrient Limitations

Genome-scale stoichiometric modeling methods, in particular Flux Balance Analysis (FBA) and variations thereof, are widely used to investigate cell metabolism and to optimize biotechnological processes. Given (1) a metabolic network, which can be reconstructed from an organism’s genome sequence, and (2) constraints on reaction rates, which may be based on measured nutrient uptake rates, FBA predicts which reactions maximize an objective flux, usually the production of cell components. Although FBA solutions may accurately predict the metabolic behavior of a cell, the actual flux predictions are often hard to interpret. This is especially the case for conditions with many constraints, such as for organisms growing in rich nutrient environments: it remains unclear why a certain solution was optimal. Here, we rationalize FBA solutions by explaining for which properties the optimal combination of metabolic strategies is selected. We provide a graphical formalism in which the selection of solutions can be visualized; we illustrate how this perspective provides a glimpse of the logic that underlies genome-scale modeling by applying our formalism to models of various sizes.

Unlocking Elementary Conversion Modes: ecmtool Unveils All Capabilities of Metabolic Networks

The metabolic capabilities of cells determine their biotechnological potential, fitness in ecosystems, pathogenic threat levels, and function in multicellular organisms. Their comprehensive experimental characterization is generally not feasible, particularly for unculturable organisms. In principle, the full range of metabolic capabilities can be computed from an organism's annotated genome using metabolic network reconstruction. However, current computational methods cannot deal with genome-scale metabolic networks. Part of the problem is that these methods aim to enumerate all metabolic pathways, while computation of all (elementally balanced) conversions between nutrients and products would suffice. Indeed, the elementary conversion modes (ECMs, defined by Urbanczik and Wagner) capture the full metabolic capabilities of a network, but the use of ECMs has not been accessible until now. We explain and extend the theory of ECMs, implement their enumeration in ecmtool, and illustrate their applicability. This work contributes to the elucidation of the full metabolic footprint of any cell.

•

Elementary conversion modes (ECMs) specify all metabolic capabilities of any organism

•

Ecmtool computes all ECMs from a reconstructed metabolic network

•

ECM enumeration enables metabolic characterization of larger networks than ever

•

Focusing on ECMs between relevant metabolites even enables genome-scale enumeration

Understanding the metabolic capabilities of cells is of profound importance. Microbial metabolism shapes global cycles of elements and cleans polluted soils. Human and pathogen metabolism affects our health. Recent advances allow for automatic reconstruction of reaction networks for any organism, which is already used in synthetic biology, (food) microbiology, and agriculture to compute optimal yields from resources to products. However, computational tools are limited to optimal states or subnetworks, leaving many capabilities of organisms hidden. Our program, ecmtool, creates a blueprint of any organism's metabolic functionalities, drastically improving insights obtained from genome sequences. Ecmtool may become essential in exploratory research, especially for studying cells that are not culturable in laboratory conditions. Ideally, elementary conversion mode enumeration will someday be a standard step after metabolic network reconstruction, achieving the metabolic characterization of all known organisms.

MCS2: minimal coordinated supports for fast enumeration of minimal cut sets in metabolic networks

Motivation

Constraint-based modeling of metabolic networks helps researchers gain insight into the metabolic processes of many organisms, both prokaryotic and eukaryotic. Minimal cut sets (MCSs) are minimal sets of reactions whose inhibition blocks a target reaction in a metabolic network. Most approaches for finding the MCSs in constrained-based models require, either as an intermediate step or as a byproduct of the calculation, the computation of the set of elementary flux modes (EFMs), a convex basis for the valid flux vectors in the network. Recently, Ballerstein et al. proposed a method for computing the MCSs of a network without first computing its EFMs, by creating a dual network whose EFMs are a superset of the MCSs of the original network. However, their dual network is always larger than the original network and depends on the target reaction. Here we propose the construction of a different dual network, which is typically smaller than the original network and is independent of the target reaction, for the same purpose. We prove the correctness of our approach, minimal coordinated support (MCS²), and describe how it can be modified to compute the few smallest MCSs for a given target reaction.

Results

We compare MCS² to the method of Ballerstein et al. and two other existing methods. We show that MCS² succeeds in calculating the full set of MCSs in many models where other approaches cannot finish within a reasonable amount of time. Thus, in addition to its theoretical novelty, our approach provides a practical advantage over existing methods.

Availability and implementation

MCS² is freely available at https://github.com/RezaMash/MCS under the GNU 3.0 license.

Supplementary information

Supplementary data are available at Bioinformatics online.

Flux tope analysis: studying the coordination of reaction directions in metabolic networks

Motivation

Elementary flux mode (EFM) analysis allows an unbiased description of metabolic networks in terms of minimal pathways (involving a minimal set of reactions). To date, the enumeration of EFMs is impracticable in genome-scale metabolic models. In a complementary approach, we introduce the concept of a flux tope (FT), involving a maximal set of reactions (with fixed directions), which allows one to study the coordination of reaction directions in metabolic networks and opens a new way for EFM enumeration.

Results

A FT is a (nontrivial) subset of the flux cone specified by fixing the directions of all reversible reactions. In a consistent metabolic network (without unused reactions), every FT contains a 'maximal pathway', carrying flux in all reactions. This decomposition of the flux cone into FTs allows the enumeration of EFMs (of individual FTs) without increasing the problem dimension by reaction splitting. To develop a mathematical framework for FT analysis, we build on the concepts of sign vectors and hyperplane arrangements. Thereby, we observe that FT analysis can be applied also to flux optimization problems involving additional (inhomogeneous) linear constraints. For the enumeration of FTs, we adapt the reverse search algorithm and provide an efficient implementation. We demonstrate that (biomass-optimal) FTs can be enumerated in genome-scale metabolic models of B.cuenoti and E.coli, and we use FTs to enumerate EFMs in models of M.genitalium and B.cuenoti.

Availability and implementation

The source code is freely available at https://github.com/mpgerstl/FTA.

Supplementary information

Supplementary data are available at Bioinformatics online.

OptMDFpathway: Identification of metabolic pathways with maximal thermodynamic driving force and its application for analyzing the endogenous CO2 fixation potential of Escherichia coli

Constraint-based modeling techniques have become a standard tool for the in silico analysis of metabolic networks. To further improve their accuracy, recent methodological developments focused on integration of thermodynamic information in metabolic models to assess the feasibility of flux distributions by thermodynamic driving forces. Here we present OptMDFpathway, a method that extends the recently proposed framework of Max-min Driving Force (MDF) for thermodynamic pathway analysis. Given a metabolic network model, OptMDFpathway identifies both the optimal MDF for a desired phenotypic behavior as well as the respective pathway itself that supports the optimal driving force. OptMDFpathway is formulated as a mixed-integer linear program and is applicable to genome-scale metabolic networks. As an important theoretical result, we also show that there exists always at least one elementary mode in the network that reaches the maximal MDF. We employed our new approach to systematically identify all substrate-product combinations in Escherichia coli where product synthesis allows for concomitant net CO₂ assimilation via thermodynamically feasible pathways. Although biomass synthesis cannot be coupled to net CO₂ fixation in E. coli we found that as many as 145 of the 949 cytosolic carbon metabolites contained in the genome-scale model iJO1366 enable net CO₂ incorporation along thermodynamically feasible pathways with glycerol as substrate and 34 with glucose. The most promising products in terms of carbon assimilation yield and thermodynamic driving forces are orotate, aspartate and the C4-metabolites of the tricarboxylic acid cycle. We also identified thermodynamic bottlenecks frequently limiting the maximal driving force of the CO₂-fixing pathways. Our results indicate that heterotrophic organisms like E. coli hold a possibly underestimated potential for CO₂ assimilation which may complement existing biotechnological approaches for capturing CO₂. Furthermore, we envision that the developed OptMDFpathway approach can be used for many other applications within the framework of constrained-based modeling and for rational design of metabolic networks.

When analyzing metabolic networks, one often searches for metabolic pathways with certain (desired) properties, for example, conversion routes that maximize the yield of a product from a given substrate. While those problems can be solved with established methods of constraint-based modeling, no algorithm is currently available for genome-scale models to identify the pathway that has the highest possible thermodynamic driving force among all solutions with predefined stoichiometric properties. This gap is closed with our new approach OptMDFpathway which is based on the recently introduced concept of Max-min Driving Force (MDF). OptMDFpathway offers various applications, especially in the context of metabolic design of cell factories. To demonstrate the power and usefulness of OptMDFpathway, we employed it to analyze the endogenous CO₂ fixation potential of Escherichia coli. While E. coli cannot assimilate CO₂ into biomass, net CO₂ fixation can take place along linear pathways from substrate to product and we show that thermodynamically feasible pathways with net CO₂ assimilation exist for 145 (34) products when choosing glycerol (glucose) as substrate. Our results indicate that heterotrophic organisms like E. coli hold a possibly underestimated potential for CO₂ assimilation which may complement existing biotechnological approaches for capturing CO₂.

A mathematical framework for yield (vs. rate) optimization in constraint-based modeling and applications in metabolic engineering

BACKGROUND

The optimization of metabolic rates (as linear objective functions) represents the methodical core of flux-balance analysis techniques which have become a standard tool for the study of genome-scale metabolic models. Besides (growth and synthesis) rates, metabolic yields are key parameters for the characterization of biochemical transformation processes, especially in the context of biotechnological applications. However, yields are ratios of rates, and hence the optimization of yields (as nonlinear objective functions) under arbitrary linear constraints is not possible with current flux-balance analysis techniques. Despite the fundamental importance of yields in constraint-based modeling, a comprehensive mathematical framework for yield optimization is still missing.

RESULTS

We present a mathematical theory that allows one to systematically compute and analyze yield-optimal solutions of metabolic models under arbitrary linear constraints. In particular, we formulate yield optimization as a linear-fractional program. For practical computations, we transform the linear-fractional yield optimization problem to a (higher-dimensional) linear problem. Its solutions determine the solutions of the original problem and can be used to predict yield-optimal flux distributions in genome-scale metabolic models. For the theoretical analysis, we consider the linear-fractional problem directly. Most importantly, we show that the yield-optimal solution set (like the rate-optimal solution set) is determined by (yield-optimal) elementary flux vectors of the underlying metabolic model. However, yield- and rate-optimal solutions may differ from each other, and hence optimal (biomass or product) yields are not necessarily obtained at solutions with optimal (growth or synthesis) rates. Moreover, we discuss phase planes/production envelopes and yield spaces, in particular, we prove that yield spaces are convex and provide algorithms for their computation. We illustrate our findings by a small example and demonstrate their relevance for metabolic engineering with realistic models of E. coli.

CONCLUSIONS

We develop a comprehensive mathematical framework for yield optimization in metabolic models. Our theory is particularly useful for the study and rational modification of cell factories designed under given yield and/or rate requirements.

From elementary flux modes to elementary flux vectors: Metabolic pathway analysis with arbitrary linear flux constraints

Elementary flux modes (EFMs) emerged as a formal concept to describe metabolic pathways and have become an established tool for constraint-based modeling and metabolic network analysis. EFMs are characteristic (support-minimal) vectors of the flux cone that contains all feasible steady-state flux vectors of a given metabolic network. EFMs account for (homogeneous) linear constraints arising from reaction irreversibilities and the assumption of steady state; however, other (inhomogeneous) linear constraints, such as minimal and maximal reaction rates frequently used by other constraint-based techniques (such as flux balance analysis [FBA]), cannot be directly integrated. These additional constraints further restrict the space of feasible flux vectors and turn the flux cone into a general flux polyhedron in which the concept of EFMs is not directly applicable anymore. For this reason, there has been a conceptual gap between EFM-based (pathway) analysis methods and linear optimization (FBA) techniques, as they operate on different geometric objects. One approach to overcome these limitations was proposed ten years ago and is based on the concept of elementary flux vectors (EFVs). Only recently has the community started to recognize the potential of EFVs for metabolic network analysis. In fact, EFVs exactly represent the conceptual development required to generalize the idea of EFMs from flux cones to flux polyhedra. This work aims to present a concise theoretical and practical introduction to EFVs that is accessible to a broad audience. We highlight the close relationship between EFMs and EFVs and demonstrate that almost all applications of EFMs (in flux cones) are possible for EFVs (in flux polyhedra) as well. In fact, certain properties can only be studied with EFVs. Thus, we conclude that EFVs provide a powerful and unifying framework for constraint-based modeling of metabolic networks.