New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Understanding bistability in complex enzymedriven reaction networks

Communicated by Avner Friedman, Ohio State University, Columbus, OH, April 6, 2006 (received for review February 1, 2006)
Abstract
Much attention has been paid recently to bistability and switchlike 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 switchlike behavior results from a wellstudied mechanism for the action of human dihydrofolate reductase, an important anticancer target.
Even at the most basic levels of cell operation, there is evidence, both theoretical and experimental, that important metabolic pathways have the capacity to exhibit bistability, i.e., to exhibit two very different stable steady states, with switching between them driven by response to chemical signaling (1–9). Indeed, evolution would seem to favor cells that can switch quickly between radically different operational states as circumstances change, either in the cell’s environment or within the cell itself.
Origins of Bistability
Bistability is often associated with the presence of a pathway in which the product of one reaction either inhibits or promotes occurrence of another reaction, perhaps far removed along the pathway from the first reaction. It should not be supposed, however, that connections between bistability and the presence of feedback are straightforward. In an important and very general context, we argue here that there are very stringent (and somewhat odd) structural requirements a reaction network must meet if it is to engender bistability; an apparent presence of feedback is far from enough. At the same time, we also argue that bistability can result from very simple chemistry in which feedback of information from one reaction to another is not readily apparent. Indeed, taken with very traditional (e.g., massaction) kinetics, common enzymatic mechanisms for a single overall reaction, involving just one or two substrates, already carry the capacity for bistability in simple physicochemical circumstances.
There are important implications of the fact that certain chemically simple enzymecatalyzed reaction networks can give rise to bistability. If two experiments under identical operating conditions can give rise to two very different steadystate rates of conversion of reactants to products, then it is easy to presume that the fundamental chemistry in the two experiments is different, that in one there might be new reactions or enzymedegradation not present in the other. In the absence of some understanding that the capacity for bistability can be rooted in very simple chemistry, it is not hard to imagine how, in at least some instances, changes in cell behavior might be attributed to variations in enzyme transcription activity (or even to enzyme malformation) when, in fact, those changes have their origin in the reaction network itself.
If this understanding is to be both broad and deep, one must ultimately confront the fact that there is an overwhelming array of distinct enzymecatalyzed biochemical networks that might be key operatives in various aspects of cell function. Each has its own structure, and, for each, the governing system of (nonlinear) equations will be daunting. Some networks will be sources of bistability, while others will not, no matter what values the various kinetic parameters take. Our aim here is to provide general tools for drawing precise connections between reaction network structure and the capacity for bistability. We intend these tools to be applicable to quite complex networks of enzymecatalyzed networks and to be sufficiently subtle as to distinguish between very similar networks having rather different capacities for bistability.
We will be more precise in A Motivational Example, but it is worth stating here that when we say that a reaction network has the capacity for bistability, we mean that, in a certain very simple physicochemical context, there exist combinations of parameter values (e.g., rate constants and substrate supply rates) such that the governing equations admit at least two distinct stable steady states.
A Motivational Example
For the purpose of subsequent discussion, consider a “toy” spatially homogenous cell in which only one overall reaction, S1 + S2 → P, occurs. The two substrates S1 and S2 are supplied to the cell at fixed rates, while the two substrates and the product, P, diffuse from the cell (or are degraded) at rates proportional to their concentrations within the cell. The overall conversion of substrates to product is catalyzed by a single enzyme, E, which remains within the cell.
Now imagine two copies of the cell, both supplied with substrates at precisely the same (fixed) rates, and suppose that, at steady state, the two copies are observed to exhibit radically different production rates of P. One might posit reasonable explanations of the disparity related to differences in the supply or condition of the enzyme within the two cells.
It should be understood, however, that there is another explanation, having nothing at all to do with the supply of the enzyme or its condition. Indeed, imagine that in the two cells there is a permanent fixed charge of healthy enzyme that is identical in both copies. Suppose too that the mechanism for enzyme catalysis is a classical one, shown in Eq. 1 , corresponding to unordered (random) substrate binding (10). [There is evidence that this is the mechanism driving at least some cyclindependent kinasecatalyzed phosphorylation reactions at the heart of the cell cycle (11).]
Taken with the usual mass action kinetics, the (nonlinear) differential equations governing the species concentrations within the toy cell become those shown in Eq. 2 . The overdot represents timedifferentiation, k _{1}–k _{9} are rate constants, ξ_{S1}, ξ_{S2}, ξ_{P} are mass transport (or degradation) coefficients, and F _{S1}, F _{S2} are supply rates of the substrates. (Expressions in Eq. 2 , such as F _{S1} − ξ_{S1} c _{S1}, can be identified with diffusive terms, such as ξ_{S} _{1} (c _{S1} ^{0} − c _{S1}), where c _{S1} ^{0} is a fixed concentration of S1 in a medium external to the cell.) Although the biochemical mechanism is simple, the induced differential equations are nevertheless complex.
We ask now whether these equations have the capacity to explain the existence of the rather different steady states described, one characterized by high productivity of P and the other characterized by low productivity. More precisely, we ask whether there are combinations of parameter values (i.e., rate constants, mass transfer coefficients and substrate supply rates) such that Eq. 2 is consistent with the existence of two steady states, each compatible with the same fixed availability of enzyme (i.e., compatible with an equation, in suitably chosen units, such as c _{E} + c _{ES1} + c _{ES2} + c _{ES1S2} = 1.) In fact, the answer is yes, ^{‖}
^{‖}There are many combinations of parameters that give rise to bistability. Fig. 1 was drawn for the following choice (see Supporting Appendix): k _{1} = 93.43, k _{2} = 2539, k _{3} = 481.6, k _{4} = 1,183, k _{5} = 1,556, k _{6} = 121,192, k _{7} = 0.02213, k _{8} = 1,689, k _{9} = 85,842, ξ_{S1} = 1, ξ_{Σ2} = 1, ξ_{P} = 1, F _{S1} = 2,500, F _{S2} = 1,500. All initial conditions were chosen such that c _{E} + c _{ES1} + c _{ES2} + c _{ES1S2} = 1. Thereafter, the reactions themselves conserve the total amount of enzyme. Ho and Li (12) have done some rudimentary bifurcation studies.
and, for one such selection of parameters, we show some substrate–product composition trajectories in Fig. 1. There are two stable steady states, one characterized by a high productivity of P and the other by a substantially lower one. The steady state actually visited depends on the initial conditions within the cell. Switching between steady states would result, for example, from a signal in the form of a temporary disturbance in a substrate supply rate. (In terms of the extracellular medium picture alluded to earlier, such a disturbance might correspond to a temporary perturbation in, for example, the extracellular concentration of S1.)Consideration of this very simple example is meant to make an important point: The capacity for bistability is already present in certain biochemical reactions of the most elementary kind. The presence of apparent “feedback loops” in the overall biochemistry is not a necessary component of switching phenomena, given that sources of bistability can lurk behind the fine mechanistic details of even a single overall reaction. Although the toy cell picture was invoked merely to indicate the capacity for bistability in a simple situation, it should be noted that the governing equations are, in structural terms, reflective of those commonly used to model more sophisticated aspects of cell behavior (see, for example, refs. 4 and 13).
With these ideas in mind, we aim to provide a rigorous conceptual basis for understanding the relationship between the detailed structure of massaction biochemical reaction networks and their capacity for bistability. That relationship is quite subtle, as Table 1 indicates. In each entry we show the mechanism for enzyme catalysis at the underlying massaction level, the overall reaction(s), and the capacity for bistability in the same elementary context discussed earlier. That is, substrates and inhibitors are supplied at fixed rates; total concentrations of enzyme(s) of the various kinds are fixed; and substrates, inhibitors, and products are removed (or are degraded) at rates proportional to their current concentrations. As indicated in Origins of Bistability, when we say that a particular entry “has the capacity for bistability,” we mean that there are combinations of parameters (rate constants, mass transfer coefficients, and supply rates) such that the induced massaction differential equations admit more than one stable rest point. (Parameter values are indicated in Supporting Appendix, which is published as supporting information on the PNAS web site.)
The first six entries in Table 1 show extremely basic mechanisms for enzyme catalysis of a single overall reaction. Even at this elementary level, the capacity for bistability is already present in some of the mechanisms, notably entries 4 and 6. Moreover, the first six entries tell us that, with respect to the capacity for bistability, the mechanistic details of enzyme catalysis are of real consequence.
In Table 1, entries 7–9, we show just how delicate the relationship between reaction network structure and the capacity for bistability can be. In entries 7–9, there are multiple overall reactions, with each overall reaction promoted by its own enzyme. For each overall reaction, the catalytic mechanism depicted is of the kind that cannot, by itself, give rise to bistability. (Those mechanisms appear earlier as entries 1 and 5.) Yet, note that entry 8 does have the capacity for bistability. Note also that entry 8 sits between two remarkably similar examples, entries 7 and 9, that do not have the capacity for bistability.
Table 1 is intended to demonstrate, through just a few examples, the subtle relationship between network structure and the capacity for bistability. The full range of distinct enzymecatalyzed reaction network structures that might present themselves for consideration in various biological contexts is, of course, highly daunting. And each new reaction network on its own would, more often than not, be very difficult to study. The equations governing entry 9, for example, involve 18 variables (species concentrations) and many parameters (e.g., rate constants and mass transfer coefficients)!
If there is to be some general and deep understanding of the behavior of intricate networks of enzymepromoted reactions (an understanding appropriate to the complexities of cell biology), it is clear that supportive theory must, simultaneously, be sufficiently powerful to accommodate great complexity in the mathematics and sufficiently subtle as to draw distinctions between networks as similar as those shown in Table 1. In fact, there already exists a body of work (14–17) that serves to draw rigorous qualitative connections between reaction network structure and the variety of dynamics that might be exhibited. The next section contains a recent theorem along these same lines that sheds considerable light on precisely where in reaction network structure sources of bistability might lie. That same theorem also indicates just how enzymedriven reaction networks can carry the seeds of bistability in their finest mechanistic details.
What Is Required of a Reaction Network If It Is to Engender Bistability?
All reaction networks considered hereafter are taken to be governed by massaction kinetics. That is, reaction networks, like those shown in Table 1, are taken to be descriptions of elementary chemical events, not coarsegrained overall reactions. There are several reasons for this: First, coarsegrained models in enzyme kinetics (Michaelis–Menten being the simplest) derive, as approximations, from more fundamental massactiongoverned elementary reactions, and one can never be certain that specious dynamical phenomena aren’t introduced as artifacts of the approximation. Second, the approximations themselves are sometimes wrenched from the context in which they were derived and applied, without care, in rather different contexts. Third, the fine mechanistic details have real consequences (18), as should be clear from the differences in behavior across the first six entries in Table 1. And, finally, massaction kinetics bears a transparent, coherent, and precise relationship to the network’s architecture, a relationship that ultimately makes possible very powerful statements about the connection between reaction network structure and varieties of behavior that the induced differential equations might exhibit.
We proceed by way of a concrete problem that has nothing to do with enzyme catalysis but that is better suited to the exposition. Subsequently, however, we shall return to enzymedriven reaction networks, in particular those in Table 1. Consider the reaction network shown in Eq. 3 , and suppose also that the network is, for the purposes of the example, operative within a homogeneous cell. All species are fed to the cell at fixed rates, F _{A}, F _{B}, …, F _{G}, and all species are removed from the cell (or are degraded) at rates proportional to their current concentrations within the cell.
The proportionality constants are ξ_{A}, …, ξ_{G}. The massaction differential equations are shown in Eq. 4 , where k _{1}, …, k _{8} are reaction rate constants. Note that the system (Eq. 4 ) contains 22 parameters (rate constants, species supply rates, and mass transfer coefficients).
We now ask whether there is some combination of the 22 parameters such that Eq. 4 admits more than one positive steady state (that is, at which all species concentrations are positive). We are asking whether the network shown in Eq. 3 has the capacity for multiple steady states in the context described.
This question is difficult, but recent theory provides an almost immediate answer. First, we need to explain how a reaction network (in particular, Eq. 3 ) gives rise to what we call its speciesreaction(SR) graph. Then we shall state a theorem that indicates how the SR graph gives immediate and very subtle information about the network’s capacity for multiple steady states. We need a small amount of vocabulary introduced by way of our example: The species of the network in Eq. 3 are, of course, A, B, C, …, G. The objects, such as A + B and F, appearing at the heads and tails of the reaction arrows are called the complexes of the network; thus, the complexes in the network (Eq. 3 ) are A + B, F, C + G, A, C + D, B, C + E, and D. The reactions of the network are evident.
We depict in Fig. 2 the SR graph for the network (Eq. 3 ). Its construction is simple and is, in fact, reminiscent of reaction diagrams drawn in biochemistry: Note that there is a symbol for each of the species, and, within boxes, a symbol for each of the reactions. (Reversible reaction pairs are drawn within the same box.) If a species appears within a reaction, then an arc is drawn from the species symbol to the reaction symbol, and the arc is labeled with the name of the complex in which the species appears. For example, species A appears within the reaction(s) A + B ⇄ F. Thus, an arc is drawn from A to reactions A + B ⇄ F and labeled with the complex A + B. The SR graph is completed once the species nodes are connected to the reaction nodes in the manner described. (If a species appears in both complexes of a reaction, as in A + B ⇄ 2A, then two arcs are drawn, each labeled by a different complex.)
Before we indicate how the diagram gives information about bistability, we need a little more terminology. By a complex pair(cpair) in the SR graph, we mean a pair of arcs that are adjacent to the same reaction symbol and that carry the same complex label. Thus, for example, the two arcs labeled C + E constitute a cpair because they are adjacent to the same reaction symbol and carry the same complex label, C + E. In Fig. 2 there are four cpairs, which we have colored blue (C + E), green (C + G), red (C + D) and purple (A + B). Because of their importance in what follows, we have colored cpairs to make them more evident. For readers without access to color, when we refer to a color in the text we also indicate the corresponding complex label [e.g., blue (C + E)]. White is the default arc color and will not be used to indicate cpairs.
Note that there are cycles, labeled cycle 1 and cycle 2, in Fig. 2. (By a cycle we mean a closed path in which no arc or vertex is traversed twice. In fact, there is a third unlabeled outer cycle, traversing species A–C–D–B–A.) We need to describe three kinds of cycles: oddcycles, evencycles, and 1cycles. These classifications are not mutually exclusive. A cycle can, for example, be both an oddcycle and a 1cycle.
By an oddcycle in an SR graph we mean a cycle containing an odd number of cpairs. Similarly, an evencycle contains an even number of cpairs. Note that cycle 1 contains only one cpair, the purple (A + B) one, so cycle 1 is an oddcycle. Similarly, cycle 2 contains only the red (C + D) cpair, so cycle 2 also is an oddcycle. Finally, the large outer cycle contains only the purple (A + B) cpair, so it too is an oddcycle.
To indicate what a 1cycle is, we first observe that every arc has at one end a species node. Moreover, that species has a stoichiometric coefficient (usually an implied “1”) in the complex that labels the arc. In this way, we can associate uniquely with each arc a stoichiometric coefficient. A 1cycle in an SR graph is a cycle such that the stoichiometric coefficient associated with each of its arcs is a 1. Note that if, for a particular reaction network, every (nonzero) stoichiometric coefficient is a 1 then, in that network’s SR graph, every cycle will be a 1cycle, which will be the case more often than not, and it is the case with the network shown in Eq. 3 .
Finally, we say that two cycles split a cpair if both arcs of the cpair are among the arcs of the two cycles, and one of the arcs is contained in just one of the cycles. In Fig. 2, cycles 1 and 2 split the red (C + D) cpair: Both red arcs are among the arcs of the two cycles, but cycle 1 contains only one of the red (C + D) arcs. The large outer cycle and cycle 1 also split the red cpair, as do the large outer cycle and cycle 2.
The following theorem (19–21) extends considerably an earlier theorem of Schlosser and Feinberg (17). (The theorem below is a special case of theorem 1.1 in ref. 21; the connection between the two is given by lemma 7.1 in ref. 21.)
Theorem.
Consider a reaction network for which the SR graph has both of the following properties. (i) Each cycle in the graph is a 1cycle, an oddcycle, or both. (ii) No cpair is split by two evencycles. For such a reaction network, the corresponding massaction differential equations **
**We mean here the differential equations formulated for the homogeneous cell, as in the passage from the network in Eq. 3 to Eq. 4 , having the largest number of free parameters, i.e., with all species permitted to pass through the cell boundary. When, as in the toy cell example, certain species are entrapped or, more precisely, do not give rise to terms of the form −ξ_{A} c _{A} in the model differential equations, the theorem nevertheless remains true with a technical modification: In this case, the hypothesis serves to deny the capacity of a network to engender two nondegenerate stoichiometrically compatible positive steady states (22). In subsequent discussion of Table 1, then, preclusion of multiple steady states should be understood to carry, implicitly, this technical qualification.
cannot admit more than one positive steady state, no matter what (positive) values the rate constants, effluent coefficients, or species supply rates take.The theorem indicates that, for a massaction network to have the capacity for more than one steady state, its SR graph must satisfy quite stringent conditions. In particular, when every stoichiometric coefficient is a 1 (the common situation), two evencycles must split a cpair. In fact, the theorem tells us immediately that the complex system shown in Eq. 4 cannot admit more than one positive steady state, no matter what positive values are assigned to the 22 parameters. In this case, each cycle in Fig. 2 is both an oddcycle and a 1cycle. Although there is a split cpair [the red (C + D) one], it is not split by two evencycles.
In light of the stringent necessary conditions for bistability that the theorem provides, it is significant (and perhaps essential to biology itself) that certain common mechanisms for enzyme catalysis actually meet those conditions. In fact, the theorem is already sufficiently delicate as to discriminate correctly between those networks in Table 1 that do have the capacity for multiple steady states and those that do not. (Note that when conditions i and ii are not satisfied, the theorem gives no information. When the table indicates that a network does have the capacity for bistability, it means that parameter values, indicated in Supporting Appendix, were found to elicit bistability. In such instances, the network must violate either condition i or ii.) Space limitations preclude the drawing of all nine SR graphs, but we can draw some and describe others.
The SR graphs for entries 1–3 in Table 1 each contain just one cycle. Because all stoichiometric coefficients are 1, each such cycle is a 1cycle. Moreover, because each SR graph contains only one cycle, there is no possibility of a split cpair. Thus, the conditions of the theorem are satisfied, and bistability is precluded, no matter what values the parameters take.
In Fig. 3 we show the SR graph for entry 4. Here, the situation is very different. Every cycle is again a 1cycle, but now there are several cycles. In particular, the red (E + S) cpair is split by two evencycles: The red cpair is split by the large outer cycle containing two cpairs, red (E + S) and yellow (EI + S), and a smaller cycle containing two cpairs, orange (ES + I) and blue (E + I), and just one red (E + S) arc. In this case, the hypothesis of the theorem is not satisfied, and the theorem stands silent. In fact, entry 4 does have the capacity to admit more then one steady state (see Supporting Appendix).
The SR graph for Table 1, entry 5, contains just one cycle, and it is a 1cycle; thus, bistability is precluded. The SR graph for entry 6 is similar to one drawn in the example below, where we discuss the unordered binding of two substrates to dihydrofolate reductase. All cycles are again 1cycles, but there is a cpair split by two evencycles. Bistability is not precluded. Indeed, as we indicated earlier, there are parameter values for which bistability is exhibited.
We turn now to consideration of the strikingly similar entries 7–9 in Table 1, of which only the “middle case” entry 8 has the capacity for multiple steady states. Again, space limitations preclude inclusion of SR graphs for all three entries, so we present only the simplest of these (for entry 7) in Fig. 4. Then we can indicate just how the disparate entry 8 is subtly different from entries 7 and 9. There are several cycles in Fig. 4, but it will be important for us to focus on the large outer one traversing (in order) species S1, E1, ES1S2, S2, E2, with a return to S1 via the two yellow (E2 + 2S1) arcs. Note that this cycle is not a 1cycle, for the topmost yellow arc has an associated stoichiometric coefficient of 2. On the other hand, that same cycle is an oddcycle, because it contains three cpairs, two red (S1 + E1, S2 + E2) and one yellow (E2 + 2S1). (Note that the white arcs containing the complex label E1S1S2 do not constitute a cpair, because they are not adjacent to the same reaction node.) In fact, both conditions of the theorem are satisfied, and bistability is quickly precluded.
The situation for entry 8 is different. The SR graph is similar but larger, with an additional “layer” containing a third “red” (S3 + E3) cpair of arcs (analogous to the red S1 + E1 and S2 + E2 arcs shown in Fig. 4). Again, there will be a large outer cycle, with an arc having a stoichiometric coefficient of 2. Thus, that cycle is not a 1cycle, but for entry 8 it will not be an oddcycle either. (There are four cpairs.) In this case, hypothesis i of the theorem is not satisfied, and bistability is not precluded. As indicated in Table 1, network 8 does have the capacity for bistability (see Supporting Appendix).
For entry 9 the situation returns to what it was for entry 7: The problematic cycle containing a 2 once again become an oddcycle (with five cpairs), and bistability is again precluded.
Example: A Mechanism Underlying the Operation of a Classical AntiCancer Target, Dihydrofolate Reductase (DHFR)
DHFR is a crucial enzyme along the pathway for thymine synthesis. Because thymine is, in turn, a necessary ingredient in DNA synthesis, its absence thwarts cell division. For this reason, strong inhibitors of DHFR, such as methotrexate, provide the backbone of some classical chemotherapy regimens. DHFR promotes the overall reaction shown below as Eq. 5 .
But there is much important mechanistic detail hidden behind this overall reaction. Appleman et al. (23) have proposed a mechanism at the level of elementary reactions for the operation of human DHFR and have measured or inferred massaction rate constants for the various steps. The comprehensive kinetics (at 20°C, pH 7.65) is shown in Fig. 5. First and secondorder rate constants are, respectively, in units of sec^{−1} and μM^{−1}·sec^{−1}. Despite the presence of 13 reversible reactions, the mechanism is quite ordinary: The first group of reactions (group A) merely represent binding of substrates to the enzyme, followed by the chemical transformation of substrates to products. The second group (group B) merely represents unbinding of the products from the enzyme, while the third group (group C) represents rebinding of products to enzyme–substrate complexes. Apart from binding–unbinding steps, there is only one step representing chemical transformation, ENHH2F ENH4F, where E is DHFR and NH is NADPH.
Although the rebinding steps in Group C provide a degree of feedback in the form of product inhibition, it should be noted that the steps in Group A are already reminiscent of the far simpler unordered twosubstrate mechanism 6 in Table 1, a mechanism that has the capacity for bistability and for which composition trajectories were drawn in Fig. 1. Indeed, Fig. 6 depicts a fragment of the SR graph for the putative DHFR mechanism, drawn for the subnetwork in blocks A and B. Cycle 1 (passing through species NH, EH2F, and E) contains two cpairs, colored purple (NH + EH2F) and green (NH + E). Cycle 2 (passing through E, ENH, and H2F) also contains two cpairs, colored yellow (H2F + ENH) and red (H2F + E). Note that these evencycles split the red cpair so that the conditions of the theorem are not met. In this case, the capacity for multiple steady states is not precluded.
In fact, simulations based on the rate constants of Appleman et al. (23) indicate a simple situation in which sharp switchlike behavior would be elicited: Consider a 3.5ml stirred cell having inlet and outlet ports, each covered by a membrane that denies passage of DHFR but that permits smaller substrate and product molecules to pass freely. A solution containing H2F and NH in concentrations of 100 μM and 400 μM, respectively, is fed continuously through the inlet port at a constant flow rate of F ml/min, and mixture is drawn from the cell through the outlet port, also at F ml/min. (The concentrations of substrates and products in the effluent mixture are presumed identical to those in the cell.) The cell contains entrapped DHFR at a total concentration of 15 nM. (This is the total concentration of all DHFRcontaining entities, with or without bound substrates and products.)
Computations indicate that, for values of F of <0.530 ml/min and >0.577 ml/min, there is a unique steadystate composition within the cell (and within the effluent stream). For values of F in between, however, there are three steady states, two of them stable. The computed locus of steadystate concentrations of tetrahydrofolate vs. F is shown in Fig. 7. Note that, for reactant supply rates of <0.53 ml/min, the steadystate conversion of dihydrofolate, supplied at 100 μM, is nearly complete, for the resulting concentration of tetrahydrofolate is ≈96 μM. In a sequence of simulated experiments in which the supply rate is elevated in small increments (see Fig. 7), the resulting steadystate effluent concentrations of tetrahydrofolate remain fairly constant (in the vicinity of 96 μM) until the flow rate exceeds 0.577 ml/min. At that point, the steadystate conversion of dihydrofolate to tetrahydrofolate drops markedly, with an effluent tetrahydrofolate concentration of ≈65 μM. A subsequent reduction of the flow rate in small increments results in an elevation of the steadystate tetrahydrofolate concentration until the flow rate drops just below 0.53 ml/min, at which point it increases discontinuously from 80 to 97 μM. Note that these pronounced transitions take place within an extremely narrow flow rate range. In this sense, the simulations indicate a sharp and dramatic switchlike response to small changes in the reactant supply rate.
Discussion
The theorem presented here gives precise and delicate conditions a mass action reaction network must satisfy if it is to engender more than one stationary steady state in a simple physicochemical setting. Despite the apparent severity of those conditions, it is noteworthy that they are met by classical and quite ordinary mechanisms for enzyme catalysis for a single overall reaction, including the mechanism posited for the action of human dihydrofolate reductase by Appleman et al. (23). Indeed, in many instances, the basis for switchlike behavior might lie not in broad apparent feedback across a pathway but, rather, in the fine details of enzyme catalysis at the level of individual overall reactions. If so, the cell’s arsenal of biochemical switches, evolved for response to changing circumstances, might be broader than is generally supposed.
Acknowledgments
We thank T. Mitchison and J. Chalmers for discussions. This work was supported by National Science Foundation Grant BES 0425459 and Agreement 0112050.
Footnotes
 ^{¶}To whom correspondence should be addressed. Email: feinberg.14{at}osu.edu

Author contributions: G.C. and M.F. designed research; G.C., Y.T., and M.F. performed research; G.C., Y.T., and M.F. contributed analytic tools; and G.C. and M.F. wrote the paper.

^{‖}There are many combinations of parameters that give rise to bistability. Fig. 1 was drawn for the following choice (see Supporting Appendix): k _{1} = 93.43, k _{2} = 2539, k _{3} = 481.6, k _{4} = 1,183, k _{5} = 1,556, k _{6} = 121,192, k _{7} = 0.02213, k _{8} = 1,689, k _{9} = 85,842, ξ_{S1} = 1, ξ_{Σ2} = 1, ξ_{P} = 1, F _{S1} = 2,500, F _{S2} = 1,500. All initial conditions were chosen such that c _{E} + c _{ES1} + c _{ES2} + c _{ES1S2} = 1. Thereafter, the reactions themselves conserve the total amount of enzyme. Ho and Li (12) have done some rudimentary bifurcation studies.

↵**We mean here the differential equations formulated for the homogeneous cell, as in the passage from the network in Eq. 3 to Eq. 4 , having the largest number of free parameters, i.e., with all species permitted to pass through the cell boundary. When, as in the toy cell example, certain species are entrapped or, more precisely, do not give rise to terms of the form −ξ_{A} c _{A} in the model differential equations, the theorem nevertheless remains true with a technical modification: In this case, the hypothesis serves to deny the capacity of a network to engender two nondegenerate stoichiometrically compatible positive steady states (22). In subsequent discussion of Table 1, then, preclusion of multiple steady states should be understood to carry, implicitly, this technical qualification.

Conflict of interest statement: No conflicts declared.
 Abbreviations:
 SR,
 speciesreaction;
 cpair,
 complex pair;
 DHFR,
 dihydrofolate reductase.
Abbreviations:
 © 2006 by The National Academy of Sciences of the USA
References

↵
 Lisman J. E.

↵
 Murray A. W. ,
 Kirschner M. W.

↵
 Gonze D. ,
 Goldbeter A.

↵
 Chen K. C. ,
 CsikaszNagy A. ,
 Gyorffy B. ,
 Val J. ,
 Novak B. ,
 Tyson J. J.

↵
 Cross F. R. ,
 Archambault V. ,
 Miller M. ,
 Klovstad M.

↵
 Sha W. ,
 Moore J. ,
 Chen K. ,
 Lassaletta A. D. ,
 Yi C.S. ,
 Tyson J. J. ,
 Sible J. C.
 ↵
 ↵

↵
 Markevich N. I. ,
 Hoek J. B. ,
 Kolodenko B. N.

↵
 Segel I. H.

↵
 Clare P. M. ,
 Poorman R. A. ,
 Kelley L. C. ,
 Watenpaugh K. D. ,
 Bannow B. A. ,
 Leach K. L.
 ↵

↵
 Lee E. ,
 Salic A. ,
 Kruger R. ,
 Heinrich R. ,
 Kirschner M. W.
 ↵

↵
 Feinberg M.
 ↵
 ↵
 ↵

↵
 Craciun G.
 ↵
 ↵

↵
 Craciun G. ,
 Feinberg M.

↵
 Appleman J. R. ,
 Beard W. A. ,
 Delcamp T. J. ,
 Prendergast N. J. ,
 Freisheim J. H. ,
 Blakley R. L.
Citation Manager Formats
Sign up for Article Alerts
Jump to section
You May Also be Interested in
More Articles of This Classification
Physical Sciences
Applied Mathematics
Biological Sciences
Cell Biology
Related Content
 No related articles found.
Cited by...
 A structural classification of candidate oscillators and multistationary systems
 Cellcycle transitions: a common role for stoichiometric inhibitors
 Parameterfree methods distinguish Wnt pathway models and guide design of experiments
 Models of signalling networks  what cell biologists can gain from them and give to them
 Pattern formation of Rho GTPases in single cell wound healing
 Signal processing in cellular clocks
 Structural Sources of Robustness in Biochemical Reaction Networks
 Building biological memory by linking positive feedback loops
 Underexplored Niches in Research on Plant Pathogenic Bacteria
 Subnetwork analysis reveals dynamic features of complex (bio)chemical networks
 Understandable Complexity
 Understandable Complexity