How catalytic mechanisms reveal themselves in multiple steady-state data: I. Basic principles

Multistationarity in the Space of Total Concentrations for Systems that Admit a Monomial Parametrization

We apply tools from real algebraic geometry to the problem of multistationarity of chemical reaction networks. A particular focus is on the case of reaction networks whose steady states admit a monomial parametrization. For such systems, we show that in the space of total concentrations multistationarity is scale invariant: If there is multistationarity for some value of the total concentrations, then there is multistationarity on the entire ray containing this value (possibly for different rate constants)-and vice versa. Moreover, for these networks it is possible to decide about multistationarity independent of the rate constants by formulating semi-algebraic conditions that involve only concentration variables. These conditions can easily be extended to include total concentrations. Hence, quantifier elimination may give new insights into multistationarity regions in the space of total concentrations. To demonstrate this, we show that for the distributive phosphorylation of a protein at two binding sites multistationarity is only possible if the total concentration of the substrate is larger than either the total concentration of the kinase or the total concentration of the phosphatase. This result is enabled by the chamber decomposition of the space of total concentrations from polyhedral geometry. Together with the corresponding sufficiency result of Bihan et al., this yields a characterization of multistationarity up to lower-dimensional regions.

Stochastic multistationarity in a model of the hematopoietic stem cell differentiation network

A central issue in the analysis of multi-stable systems is that of controlling the relative size of the basins of attraction of alternative states through suitable choices of system parameters. We are interested here mainly in the stochastic version of this problem, that of shaping the stationary probability distribution of a Markov chain so that various alternative modes become more likely than others. Although many of our results are more general, we were motivated by an important biological question, that of cell differentiation. In the mathematical modeling of cell differentiation, it is common to think of internal states of cells (quanfitied by activation levels of certain genes) as determining the different cell types. Specifically, we study here the "PU.1/GATA-1 circuit" which is involved in the control of the development of mature blood cells from hematopoietic stem cells (HSCs). All mature, specialized blood cells have been shown to be derived from multipotent HSCs. Our first contribution is to introduce a rigorous chemical reaction network model of the PU.1/GATA-1 circuit, which incorporates current biological knowledge. We then find that the resulting ODE model of these biomolecular reactions is incapable of exhibiting multistability, contradicting the fact that differentiation networks have, by definition, alternative stable steady states. When considering instead the stochastic version of this chemical network, we analytically construct the stationary distribution, and are able to show that this distribution is indeed capable of admitting a multiplicity of modes. Finally, we study how a judicious choice of system parameters serves to bias the probabilities towards different stationary states. We remark that certain changes in system parameters can be physically implemented by a biological feedback mechanism; tuning this feedback gives extra degrees of freedom that allow one to assign higher likelihood to some cell types over others.

Dimensionality reduction of bistable biological systems

Time hierarchies, arising as a result of interactions between system's components, represent a ubiquitous property of dynamical biological systems. In addition, biological systems have been attributed switch-like properties modulating the response to various stimuli across different organisms and environmental conditions. Therefore, establishing the interplay between these features of system dynamics renders itself a challenging question of practical interest in biology. Existing methods are suitable for systems with one stable steady state employed as a well-defined reference. In such systems, the characterization of the time hierarchies has already been used for determining the components that contribute to the dynamics of biological systems. However, the application of these methods to bistable nonlinear systems is impeded due to their inherent dependence on the reference state, which in this case is no longer unique. Here, we extend the applicability of the reference-state analysis by proposing, analyzing, and applying a novel method, which allows investigation of the time hierarchies in systems exhibiting bistability. The proposed method is in turn used in identifying the components, other than reactions, which determine the systemic dynamical properties. We demonstrate that in biological systems of varying levels of complexity and spanning different biological levels, the method can be effectively employed for model simplification while ensuring preservation of qualitative dynamical properties (i.e., bistability). Finally, by establishing a connection between techniques from nonlinear dynamics and multivariate statistics, the proposed approach provides the basis for extending reference-based analysis to bistable systems.

Concordant chemical reaction networks

We describe a large class of chemical reaction networks, those endowed with a subtle structural property called concordance. We show that the class of concordant networks coincides precisely with the class of networks which, when taken with any weakly monotonic kinetics, invariably give rise to kinetic systems that are injective - a quality that, among other things, precludes the possibility of switch-like transitions between distinct positive steady states. We also provide persistence characteristics of concordant networks, instability implications of discordance, and consequences of stronger variants of concordance. Some of our results are in the spirit of recent ones by Banaji and Craciun, but here we do not require that every species suffer a degradation reaction. This is especially important in studying biochemical networks, for which it is rare to have all species degrade.

Parameter-free model discrimination criterion based on steady-state coplanarity

We introduce a procedure for deciding when a mass-action model is incompatible with observed steady-state data that does not require any parameter estimation. Thus, we avoid the difficulties of nonlinear optimization typically associated with methods based on parameter fitting. Instead, we borrow ideas from algebraic geometry to construct a transformation of the model variables such that any set of steady states of the model under that transformation lies on a common plane, irrespective of the values of the model parameters. Model rejection can then be performed by assessing the degree to which the transformed data deviate from coplanarity. We demonstrate our method by applying it to models of multisite phosphorylation and cell death signaling. Our framework offers a parameter-free perspective on the statistical model selection problem, which can complement conventional statistical methods in certain classes of problems where inference has to be based on steady-state data and the model structures allow for suitable algebraic relationships among the steady-state solutions.

Multistationarity in mass action networks with applications to ERK activation

Ordinary Differential Equations (ODEs) are an important tool in many areas of Quantitative Biology. For many ODE systems multistationarity (i.e. the existence of at least two positive steady states) is a desired feature. In general establishing multistationarity is a difficult task as realistic biological models are large in terms of states and (unknown) parameters and in most cases poorly parameterized (because of noisy measurement data of few components, a very small number of data points and only a limited number of repetitions). For mass action networks establishing multistationarity hence is equivalent to establishing the existence of at least two positive solutions of a large polynomial system with unknown coefficients. For mass action networks with certain structural properties, expressed in terms of the stoichiometric matrix and the reaction rate-exponent matrix, we present necessary and sufficient conditions for multistationarity that take the form of linear inequality systems. Solutions of these inequality systems define pairs of steady states and parameter values. We also present a sufficient condition to identify networks where the aforementioned conditions hold. To show the applicability of our results we analyse an ODE system that is defined by the mass action network describing the extracellular signal-regulated kinase (ERK) cascade (i.e. ERK-activation).

Emergence of switch-like behavior in a large family of simple biochemical networks

Bistability plays a central role in the gene regulatory networks (GRNs) controlling many essential biological functions, including cellular differentiation and cell cycle control. However, establishing the network topologies that can exhibit bistability remains a challenge, in part due to the exceedingly large variety of GRNs that exist for even a small number of components. We begin to address this problem by employing chemical reaction network theory in a comprehensive in silico survey to determine the capacity for bistability of more than 40,000 simple networks that can be formed by two transcription factor-coding genes and their associated proteins (assuming only the most elementary biochemical processes). We find that there exist reaction rate constants leading to bistability in ∼90% of these GRN models, including several circuits that do not contain any of the TF cooperativity commonly associated with bistable systems, and the majority of which could only be identified as bistable through an original subnetwork-based analysis. A topological sorting of the two-gene family of networks based on the presence or absence of biochemical reactions reveals eleven minimal bistable networks (i.e., bistable networks that do not contain within them a smaller bistable subnetwork). The large number of previously unknown bistable network topologies suggests that the capacity for switch-like behavior in GRNs arises with relative ease and is not easily lost through network evolution. To highlight the relevance of the systematic application of CRNT to bistable network identification in real biological systems, we integrated publicly available protein-protein interaction, protein-DNA interaction, and gene expression data from Saccharomyces cerevisiae, and identified several GRNs predicted to behave in a bistable fashion.

Switch-like behavior is found across a wide range of biological systems, and as a result there is significant interest in identifying the various ways in which biochemical reactions can be combined to yield a switch-like response. In this work we use a set of mathematical tools from chemical reaction network theory that provide information about the steady-states of a reaction network irrespective of the values of network rate constants, to conduct a large computational study of a family of model networks consisting of only two protein-coding genes. We find that a large majority of these networks (∼90%) have (for some set of parameters) the mathematical property known as bistability and can behave in a switch-like manner. Interestingly, the capacity for switch-like behavior is often maintained as networks increase in size through the introduction of new reactions. We then demonstrate using published yeast data how theoretical parameter-free surveys such as this one can be used to discover possible switch-like circuits in real biological systems. Our results highlight the potential usefulness of parameter-free modeling for the characterization of complex networks and to the study of network evolution, and are suggestive of a role for it in the development of novel synthetic biological switches.

The capacity for multistability in small gene regulatory networks

BACKGROUND

Recent years have seen a dramatic increase in the use of mathematical modeling to gain insight into gene regulatory network behavior across many different organisms. In particular, there has been considerable interest in using mathematical tools to understand how multistable regulatory networks may contribute to developmental processes such as cell fate determination. Indeed, such a network may subserve the formation of unicellular leaf hairs (trichomes) in the model plant Arabidopsis thaliana.

RESULTS

In order to investigate the capacity of small gene regulatory networks to generate multiple equilibria, we present a chemical reaction network (CRN)-based modeling formalism and describe a number of methods for CRN analysis in a parameter-free context. These methods are compared and applied to a full set of one-component subnetworks, as well as a large random sample from 40,680 similarly constructed two-component subnetworks. We find that positive feedback and cooperativity mediated by transcription factor (TF) dimerization is a requirement for one-component subnetwork bistability. For subnetworks with two components, the presence of these processes increases the probability that a randomly sampled subnetwork will exhibit multiple equilibria, although we find several examples of bistable two-component subnetworks that do not involve cooperative TF-promoter binding. In the specific case of epidermal differentiation in Arabidopsis, dimerization of the GL3-GL1 complex and cooperative sequential binding of GL3-GL1 to the CPC promoter are each independently sufficient for bistability.

CONCLUSION

Computational methods utilizing CRN-specific theorems to rule out bistability in small gene regulatory networks are far superior to techniques generally applicable to deterministic ODE systems. Using these methods to conduct an unbiased survey of parameter-free deterministic models of small networks, and the Arabidopsis epidermal cell differentiation subnetwork in particular, we illustrate how future experimental research may be guided by network structure analysis.

The effects of reversibility and noise on stochastic phosphorylation cycles and cascades

The phosphorylation-dephosphorylation cycle is a common motif in cellular signaling networks. Previous work has revealed that, when driven by a noisy input signal, these cycles may exhibit bistable behavior. Here, a recently introduced theorem on network bistability is applied to prove that the existence of bistability is dependent on the stochastic nature of the system. Furthermore, the thermodynamics of simple cycles and cascades is investigated in the stochastic setting. Because these cycles are driven by the ATP hydrolysis potential, they may operate far from equilibrium. It is shown that sufficient high ATP hydrolysis potential is necessary for the existence of a bistable steady state. For the single-cycle system, the ensemble average behavior follows the ultrasensitive response expected from analysis of the corresponding deterministic system, but with significant fluctuations. For the two-cycle cascade, the average behavior begins to deviate from the expected response of the deterministic system. Examination of a two-cycle cascade reveals that the bistable steady state may be either propagated or abolished along a cascade, depending on the parameters chosen. Likewise, the variance in the response can be maximized or minimized by tuning the number of enzymes in the second cycle.

Subnetwork analysis reveals dynamic features of complex (bio)chemical networks

In analyzing and mathematical modeling of complex (bio)chemical reaction networks, formal methods that connect network structure and dynamic behavior are needed because often, quantitative knowledge of the networks is very limited. This applies to many important processes in cell biology. Chemical reaction network theory allows for the classification of the potential network behavior-for instance, with respect to the existence of multiple steady states-but is computationally limited to small systems. Here, we show that by analyzing subnetworks termed elementary flux modes, the applicability of the theory can be extended to more complex networks. For an example network inspired by cell cycle control in budding yeast, the approach allows for model discrimination, identification of key mechanisms for multistationarity, and robustness analysis. The presented methods will be helpful in modeling and analyzing other complex reaction networks.

Understanding bistability in complex enzyme-driven reaction networks

Much attention has been paid recently to bistability and switch-like behavior that might be resident in important biochemical reaction networks. There is, in fact, a great deal of subtlety in the relationship between the structure of a reaction network and its capacity to engender bistability. In common physicochemical settings, large classes of extremely complex networks, taken with mass action kinetics, cannot give rise to bistability no matter what values the rate constants take. On the other hand, bistable behavior can be induced in those same settings by certain very simple and classical mass action mechanisms for enzyme catalysis of a single overall reaction. We present a theorem that distinguishes between those mass action networks that might support bistable behavior and those that cannot. Moreover, we indicate how switch-like behavior results from a well-studied mechanism for the action of human dihydrofolate reductase, an important anti-cancer target.