Rare events in biochemical networks

Many biochemical networks exhibit multistability. Examples are networks that drive cell differentiation during embryonic development, the cell cycle, and apoptosis. These states are often highly stable, meaning that the probability that these networks switch spontaneously from one state to another is very small. Yet, why these states are so stable is poorly understood, especially given the recent experimental observations that gene expression tends to be very noisy. To elucidate the stability of these networks, we need to understand the pathways of switching from one state to another. Characterizing the switching pathways via experiments is, however, very difficult, precisely because these switching events are rare. Also conventional simulation techniques are of little use, since few, if any, switching events are likely to happen during the simulation. We have developed a new class of simulation techniques, called Forward Flux Sampling (FFS) [1], which makes it possible to simulate rare events efficiently in both equilibrium and non-equilibrium systems. We have applied the technique to both minimal and more detailed models of a paradigm switch, the phage lambda switch. This system is unique, because, together with the Lac system and the bacterial chemotaxis system, it is one of the best characterized model systems.

The phage-lambda switch
After the bacteriophage has injected its DNA into the host, the bacterium E. coli, the DNA is integrated into the genome of E. coli. The DNA of the phage is then replicated with that of the host – the phage is in the lysogenic state – until the DNA of the host is damaged. Upon DNA damage a switch is flipped and the phage enters the lytic state. New capsids are being formed, the DNA of the phage is excised from that of the host and put into the capsids, and the E. coli cell bursts open. The phage can now abandon the sinking ship, and find a new host. The switch from the lysogenic state to the lytic state is controlled by a genetic switch that consists of two genes that mutually repress each other (see Fig.1). The mutual repression leads to a bistable system, with one state in which the lambda repressor is strongly expressed and one in which cro is strongly expressed. Interestingly, the switch is extremely stable: the spontaneous switching rate from lysogeny to lysis is so low as to be almost undetectable, being less than 10-8 per cell per generation. We have studied a number of characteristics of this switch to elucidate its remarkable stability.

biochemical networks ten wolde Phage

Fig.1 The genetic switch of the phage lambda
The switch consists of two genes, cIcoding for the lambda repressor protein, and cro coding for the Cro protein. Both the lambda repressor and Cro can bind to three operator sites called OR3, OR2 and OR1. Lambda repressor has a high affinity for OR2 and OR1, but can also bind at higher concentrations to OR3. When lambda is bound to OR2 and OR1, the RNAP is recruited to the cI promoter and the expression of cI is activated, while when it is bound to OR3 binding of the RNAP is blocked, and the exression of cI is repressed; this leads to positive feedback at low concentrations and negative feedback at higher concentrations. Importantly, when the lambda repressor is bound to OR2 and OR1, RNAP cannot bind to the cro promoter, which means that when lambda is abundant, cro is repressed; the switch is in the lysogenic state. Now, when the DNA of the host is damaged, a protein from the host, RecA, is expressed. This protein cleaves the lambda repressor, as a result of which it can no longer bind the operator sites. Expresson of cI is now no longer activated and cro can be expressed. Moreover, the newly formed Cro proteins will bind OR3, which leads to the  repression of cI, while still allowing for cro expression. The concentration of lambda repressor now steadily decreases, while the concentration of Cro rises: the switch is flipped into the lytic state.

Overlapping regulatory domains strongly enhance switch stability
One distinguishing feature of the phage lambda switch is that the two genes share operators: both lambda repressor and Cro can bind to the operator sites OR1, OR2 and OR3 (Fig.1). Why do the two genes not have their own operator sites? To address this question, we compared the stability of a minimal model of the phage lambda switch, called the exclusive switch, to a minimal model of a hypothetical switch, called the general switch, in which the two genes have their own operator sites [2]. The simulations revealed that the exclusive switch can be orders of magnitude more stable than the general switch (Fig.2). To understand this, let’s imagine that the switch is in the X state, in which gene x, say the cI gene, is expressed. If, due to a rare fluctuation, gene y (cro) is expressed, then in the general switch the newly expressed Y protein can immediately bind to the promoter of x, switching off the expression of x. In the exclusive switch, the odds are that a newly formed Y protein cannot bind the promoter leading to the repression of x, because an x protein is still blocking it. Indeed, the exclusive switch requires a second, rare fluctuation, namely the unbinding of X from the promoter, before the switch can be flipped. Also an analysis of the transition path ensemble revealed that not only the difference in the number of proteins between X and Y, but also the state of the promoter is an important ingredient of the switching mechanism [1].

biochemical networks ten wolde GenExcSketchStab

Fig.2 Sketch of the architecture of the general switch, in which gene x and y each have their own operator, and the exclusive switch, in which gene x and y share an operator; in the latter, when protein x binds, y cannot bind the operator, and vice versa. On the right the stability of a symmetric switch, where state x and y are equally stable, as a function of N, the average number of copies of x and y. It is seen that the exclusive switch can be orders of magnitude more stable than the general switch.

Slow promoter-state fluctuations lower switch stability
Since the promoter state is important for flipping the switch, we expected that its dynamics critically affects the switching rate. We therefore computed the switching rate as a function of the time scale of the promoter state fluctuations. The simulations revealed that the switch stability decreases as the promoter dynamics becomes slower. This can be explained by noticing that once a minority species is produced and subsequently binds the promoter of the majority species to repress its expression, the concentration of the majority species will decrease more when the minority species is bound to the promoter for a longer time; with slower promoter state fluctuations, a rare event by which the minority species is produced is thus more likely to lead to the flipping of the switch [3].

biochemical networks ten wolde Phage PhaseDiag

Fig.3 The role of DNA looping and non-specific DNA binding
The lambda repressor proteins can induce a loop in the DNA by forming an octamer that bridges the OR sites and the OL sites some 2000 base pairs away (b). The panel on the right shows the phase diagram as a function of the strength of the DNA loop (y-axis) and the magnitude of non-specific binding (x-axis). The green circles indicate the region of bistability, while the blue squares indicate that only the cI state is stable and the red diamonds that only the Cro state is stable. Note that non-specific binding favors the Cro state, while DNA looping widens the region of bistability.

DNA looping increases stability phage lambda switch
An intriguing characteristic of the phage lambda switch is that the lambda repressors can induce a loop in the DNA (Fig.3). To elucidate its role in the flipping of the switch, we performed extensive simulations of a detailed model of the system [4]. While this system had been very well characterized experimentally, one other aspect besides DNA looping had not been resolved. This is non-specific binding of lambda repressor and Cro to the DNA, which effectively lowers the concentrations of these species, thereby affecting the switch stability. In the simulations, we therefore not only varied the strength of the DNA loop, but also the magnitude of non-specific binding. The simulations revealed that non-specific binding shifts the balance of the switch towards the Cro state (Fig.3). More importantly, DNA looping widens the region of bistability, thus increasing the space over which the system can function as a genuine switch (Fig.3). In addition, large-scale FFS simulations showed that DNA looping also dramatically increases the stability of the lysogenic state, lowering the switching rate by many orders of magnitude. This is because DNA looping allows the formation of a very stable complex of lambda repressor on the promoter (Fig.3b), firmly locking the switch in the lysogenic state.


  1. R. J. Allen , P. B. Warren and P. R. ten Wolde, Sampling rare switching events in biochemical networks, Physical Review Letters 94, 018104: 1–4 (2005)
  2. P. B. Warren and P. R. ten Wolde, Enhancement of the stability of genetic switches by overlapping upstream regulatory domains Physical Review Letters 92, 128101: 1–4 (2004)
  3. Marco J. Morelli , Sorin Tanase-Nicola , Rosalind J. Allen and Pieter Rein ten Wolde, Reaction coordinates for the flipping of genetic switches Biophysical Journal 94, 3413–3423 (2008)
  4. Marco J. Morelli , Pieter Rein ten Wolde and Rosalind J. Allen, DNA looping provides stability and robustness to the bacteriophage λ switch PNAS 106, 20: 8101-8106 (2009)