Positive feedback loops generate bistability

Transcription of a gene is the process by which RNA polymerase (RNAp) produce mRNA that corresponds to that gene’s coding sequence. The mRNA is then translated into a protein. The transcription rate is controlled by the quality of RNAp-binding sites on promoter, a regulatory region of DNA that precedes the gene. Whereas RNAp acts on virtually all of the genes, changes in the expression of specific genes are due to transcription factors, a kind of protein. Transcription factors affect the transcription rate by binding specific sites in the promoters of the regulated genes. When bound, they change the probability per unit time that RNAp binds the promoter and produces an mRNA molecule. The transcription factors thus affect the rate at which RNAp initiates transcription of the gene. The transcription factors can act as activators that increase the transcription rate of a gene, or as repressors that reduce the transcription rate. If the gene product P regulates its own transcription process, it would be called autoregulation or autoregulatory protein. If protein P acts as an activator, then it is a positive autoregulation. And if protein P acts as a repressor, then it is a negative autoregulation. Figure 1

Schematic positive autoregulation.

Typically, transcription and translation processes, binding of transcription factor to promoter and translating of mRNA to protein, occur much faster than protein accumulation does. The former timescale is in order of second or minute. The later timescale is in order of many minutes to hours. Thus, the production rate of a protein could be assumed in a quasi-steady-state description of promoter activity, transcription factor binding-unbinding events.

The dynamics of protein P thus are described by its production rate kb+ksA(P) size 12{k rSub { size 8{b} } +k rSub { size 8{s} } A ( P ) } {} and degradation/dilution rate kdP size 12{k rSub { size 8{d} } P} {}

dPdt=kb+ksA(P)−kdP size 12{ { { ital "dP"} over { ital "dt"} } =k rSub { size 8{b} } +k rSub { size 8{s} } A ( P ) - k rSub { size 8{d} } P} {} (eq1)

where P denotes the protein concentration. k_b, k_s, and k_d are basal, maximal product and degradation rates, respectively. The function A(P) size 12{A ( P ) } {} represents autoregulation and typically is a monotonic, S-shaped function. A useful function that describes many real gene input functions is called Hill function. The Hill input function for positive autoregulation is a curve that rises from zero and approaches a maximal saturated level

A(P)=PnKn+Pn size 12{A ( P ) = { {P rSup { size 8{n} } } over {K rSup { size 8{n} } +P rSup { size 8{n} } } } } {} (eq2)

where parameter K is activation coefficient with unit of concentration. It defines the concentration of active X needed to significantly activate expression. The half-maximal expression is reached, A(P)=1/2 size 12{A ( P ) = {1} slash {2} } {} when P=X. The Hill coefficient n represents the cooperativity of protein P on its own promoter. Typically, Hill coefficient n represents the number of protein molecule in complex, it is monomer, dimer, trimer, … when n=1, n=2, n=3, .., respectively. The Hill coefficient n governs the steepness of the input function. If n=1, the Hill input function turns out to be the Michaelis-Menten function which is a hyperbolic curve. With n>1, Hill function is a S-shaped function. And the larger is n, the steeper the function. Typically, the input function are moderately steep, with n=1-4. Figure 2

the input function A(P). With n=1, the Hill input function turns out to be the Michaelis-Menten function which is a hyperbolic curve. With n>1, Hill function A(P) is a S-shaped function.

The activation coefficient K can be changed by mutations that alter the DNA sequence of the binding site of P in the promoter of gene p. Even a change of a single DNA letter in the binding site can strengthen or weaken the chemical bonds between transcription factor and promoter and change K. The parameter K can also be varied if the position of the binding site is changed. Similarly, the maximal activity k_s can be tuned by mutations in the RNA binding site or many other factors.

The steady state concentration of protein P can be obtained by analytically or numerically solving the equation dP/dt=0 size 12{ { ital "dP"} slash { ital "dt"} =0} {}, or equivalently using rate balance plot to find out the steady state at intersection of production rate with degradation rate. The former method gives out exact solution of steady state of P, while the later gives an intuition of how the steady state of P depending on parameter changes. Here, the rate balance plot is used to illustrate how the concentration of P changes according to the parameters change, especially the maximal rate, k_s, and activation coefficient K.

Without cooperativity, n=1, the feedback function A(P) size 12{A ( P ) } {} simply is a hyperbolic curve. The Figure 3A shows the plot of degradation rate, blue line, and production rate, red lines, with maximal rate k_s equally changing. The response curve of P corresponding to intersection point of two rates with maximal rate changing is shown in Figure 3B. It shows that the steady state of protein P concentration is increasing from basal activity, P=kb/kd size 12{P= {k rSub { size 8{b} } } slash {k rSub { size 8{d} } } } {} at ks=0 size 12{k rSub { size 8{s} } =0} {}, with increment of k_s. Initially, the steady state of protein P slowly increases at small maximal rate, k_s, then almost linearly increases at large maximal rate, k_s. The small maximal rate, k_s, produces small concentration of protein P, at P<K, the feedback function A(P) size 12{A ( P ) } {} is weak or suppressed A(P)<1/2 size 12{A ( P ) < {1} slash {2} } {}. Thus, the steady state of P, which is approximated by P≈kb+ksA(P)/kd size 12{P approx left (k rSub { size 8{b} } +k rSub { size 8{s} } A ( P ) right )/k rSub { size 8{d} } } {}, with A(P)<1/2 size 12{A ( P ) < {1} slash {2} } {}, slowly increases. However, at the large enough maximal rate, k_s, to produce protein P concentration larger than activation coefficient, P>K, hence the feedback function A(P) size 12{A ( P ) } {} reaches saturation, A(P)~1 size 12{A ( P ) "~" 1} {}, thus the steady state of P now is approximated by P≈kb+ks/kd size 12{P approx left (k rSub { size 8{b} } +k rSub { size 8{s} } right )/k rSub { size 8{d} } } {} linearly increasing with maximal rate, k_s.

Similarly, the Figure 3 C and D show rate plot and response curve, respectively. The production P decreases quickly from maximal activity then reaches saturation at basal activity as long as activation coefficient K increases (Figure 3D). The reasoning is that, the larger activation coefficient K, the larger production is required to activate the feedback function, A(P)>1/2 size 12{A ( P ) > {1} slash {2} } {}, the later activation of feedback function is (Figure 3C). If the activation coefficient K is larger than maximal activity of production P, the product P would be at low level since the feedback function is never activated.

(C)(A)(D)(B)The rate balance plots and bifurcation diagrams of autoregulation with Hill coefficient n=1. (A, B) The rate balance plots and bifurcation diagrams of autoregulation, respectively, with respect to maximal rate k_s. (C, D) The rate balance plots and bifurcation diagrams of autoregulation, respectively, with respect to activation coefficient K. Typically used parameters: kb=0.1 size 12{k rSub { size 8{b} } =0 "." 1} {}, kd=1 size 12{k rSub { size 8{d} } =1} {}, n=3 size 12{n=3} {}, ks=1 size 12{k rSub { size 8{s} } =1} {}and K=0.5 size 12{K=0 "." 5} {}.

With cooperativity of protein P at its promoter, n>1, the feedback function A(P)=Pn/Kn+Pn;n>1 size 12{A ( P ) = {P rSup { size 8{n} } } slash { left (K rSup { size 8{n} } +P rSup { size 8{n} } right )} ;n>1} {} now is an S-shaped function. Thus production rate, kb+ksA(P) size 12{k rSub { size 8{b} } +k rSub { size 8{s} } A ( P ) } {}, is also an S-shaped curve. The production rate with respect to maximal rate k_s is shown in Figure 3A, red lines. The blue line is degradation rate, kdP size 12{k rSub { size 8{d} } P} {}. Now, there might be more than one intersection, two or three intersections, between production and degradation rates. Meaning that, the autoregulation could generate bistability, product P have two stable steady states (filled circles) and one unstable steady state (unfilled circle) at a certain range of parameters, Figure 4A and C. The product P can access those stable steady states by tuning parameters. For instance, in Figure 4B, as maximal rate k_s increases, the product P is low until k_s exceeds some critical value, ks2 size 12{k rSub { size 8{s} } rSup { size 8{2} } } {}, at which point the product P increases abruptly to a high value. Then if k_s decreases, the product P does not drop down until k_s crosses another critical value, ks1 size 12{k rSub { size 8{s} } rSup { size 8{1} } } {}. Notice that, for k_s value between ks1 size 12{k rSub { size 8{s} } rSup { size 8{1} } } {}and ks2 size 12{k rSub { size 8{s} } rSup { size 8{2} } } {}, the product P is bistable- that is, it has two stable steady state values (on the upper and lower braches-the solid lines) separated by an unstable steady state (on the intermediate branch – the dashed line). This sort of reversible bistability is also referred as hysteresis.

(B)(C)(A)(D)The rate balance plots and bifurcation diagrams of autoregulation with Hill coefficient n>1. (A, B) The production (red) and degradation (blue) rate plots and bifurcation diagrams of autoregulation, respectively, with respect to maximal rate k_s. (C, D) The production (red) and degradation (blue) rate plots and bifurcation diagrams of autoregulation, respectively, with respect to activation coefficient K. In rate balance plot, filled and un-filled circle represent stable and unstable steady state, respectively. In bifurcation diagram, the solid and dashed lines represent stable and unstable steady state, respectively. Typically used parameters: kb=0.1 size 12{k rSub { size 8{b} } =0 "." 1} {}, kd=1 size 12{k rSub { size 8{d} } =1} {}, n=3 size 12{n=3} {}, ks=1 size 12{k rSub { size 8{s} } =1} {}and K=0.5 size 12{K=0 "." 5} {}.

The signal-response curves in Figure 4B would be called ‘one-parameter bifurcation diagram’. The parameter is maximal rate, k_s. The steady state response, on the Y axis, is concentration of P with respect to maximal rate, k_s. At the critical points, ks1 size 12{k rSub { size 8{s} } rSup { size 8{1} } } {} and ks2 size 12{k rSub { size 8{s} } rSup { size 8{2} } } {}, the concentration of P changes abruptly from low to high concentration and vice versa. Such points are called bifurcation points, in this case, ‘saddle –node bifurcation points - SN’.

By tuning other parameter, such as activation coefficient K, the autoregulation also could generate bistability as shown by one-parameter bifurcation diagram in Figure 4D. The bistable region is in between K1 size 12{K rSup { size 8{1} } } {} and K2 size 12{K rSup { size 8{2} } } {}, which are saddle –node bifurcation points.

In fact, the two parameters, maximal rate and activation coefficient, would be tuned therefore the bistable region might changed accordingly. The Figure 5 shows bistable region of autoregulation with changing of both parameters. The bistability disappears through emergent of two saddle node points at the cusp when both of them are small. As the activation coefficient K increases, the maximal rate also need to increases to maintain the bistability. The reasoning is that, increasing of activation coefficient K causes decrement of the concentration of P, thus the maximal rate should increase in order to maintain the concentration of P.

The two-parameter bifurcation diagram of autoregulation showing the region of bistability, the region inside the cusp.

The reversible bistability or hysteresis of autoregulation has been well reported in many experiments such as the lac-operon in bacteria, the activation of M-phase-promoting factor in frog egg extracts, and the autocatalytic form.

Beside reversible bistability, the autoregulation could also generate irreversible bistability, meaning that the concentration of product P changes abruptly and irreversibly as parameter crosses a certain threshold. The irreversible bistability can be achieved by introducing an additive term to the production rate, kb+ks1A(P)+ks2S size 12{k rSub { size 8{b} } +k rSub { size 8{s1} } A ( P ) +k rSub { size 8{s2} } S} {}, the additive S term could be interpreted as an input signal or an additional transcription factor binding gene promoter. Thus, the dynamics of product P now is

dPdt=kb+ks1S+ks2A(P)−kdP size 12{ { { ital "dP"} over { ital "dt"} } =k rSub { size 8{b} } +k rSub { size 8{s1} } S+k rSub { size 8{s2} } A ( P ) - k rSub { size 8{d} } P} {} (eq3)

In fact, one gene possibly is regulated by multiple transcription factors (refs). By introducing additive transcription factor S, the autoregulation possibly generate irreversible bistability. For instance, an example is shown in the Figure 6B, the concentration of product P is low until the input signal S crosses a certain threshold S₁, at which point the product concentration increases abruptly to high level. Once jump on, the concentration of P would never drop down even the input signal is washed out, S=0. However, introducing an additional input signal S into production rate does not guarantee irreversible bistability. Figure 6C, D show regions of mono stability (un-shaded), reversible bistability (filled) and irreversible bistability (patterned) with respect to k_s and K, respectively.

The irreversible bistability is utilized in the processes characterized by a point-of-no-return such as frog oocyte maturation in response to progesterone, the frog oocyte can not returned to original state once matured; apoptosis is another decision that must be a one-way switch.

(D)(C)(B)(A)The rate balance plots (A) and bifurcation diagrams (B) of autoregulation with additional transcription factor S. The irreversible bistability could be generated by autoregulation. Typically used parameters: kb=0.1 size 12{k rSub { size 8{b} } =0 "." 1} {}, ks1=0.1 size 12{k rSub { size 8{s1} } =0 "." 1} {}, ks=1 size 12{k rSub { size 8{s} } =1} {}, kd=1 size 12{k rSub { size 8{d} } =1} {}, n=3 size 12{n=3} {}, K=0.5 size 12{K=0 "." 5} {}.

Proteins undergo a huge number of post translational modifications. However, only a subset of these modifications is reversible, covalent modifications such as acetylation, fatty acid acylation, glycosylation and phosphorylation. These modifications affect the activity, life span, or cellular location of the modified proteins. Of the post translational modifications mentioned above, phosphorylation is an important and ubiquitous one. The phosphorylation reaction in the cell is a reversible one where kinases catalyze the addition of the phosphoryl group and phosphatases catalyze the removal of the phosphoryl group, Figure 7. Phosphorylation and dephosphorylation reactions under the control of kinases and phosphatases, can occur in less than a second or over a span of hours which makes this system ideal as a regulatory process. Phosphorylation and dephosphorylation reactions can be part of a cascade of reactions which can amplify a signal which has an extracellular origin such as hormones and growth factors. The kinases and phosphatases are enzymes, a kind of protein.

Reversible phosphorylation of protein. Kinase E1 catalyzes the addition of phosphoryl group (yellow circle) and phosphatase E2 catalyzes the removal of the phosphoryl group. In turn, kinase E1 is regulated by phosphorylated protein P.

The auto-phosphorylation protein, the protein is phosphorylated/dephosphorylated by an enzyme and that enzyme is also phosphorylated/dephosphorylated by the protein, have been reported in many other sources (fission, budding yeast cell cycle, .. refs). The phosphorylation and dephosphorylation reactions of phosphorylation are typically govern by the Michaelis-Menten equation. We assume the total protein concentration, PT size 12{P rSub { size 8{T} } } {}, is constantly produced and auto degraded. And the protein exists in two forms, free form P and phosphorylated form P size 12{P rSup { size 8{*} } } {}, thus P+P=PT size 12{P+P rSup { size 8{*} } =P rSub { size 8{T} } } {}. And assume that, the free form P size 12{P} {}is active form which is in our interest. Thus, dynamics of active form or free form P size 12{P} {} can be described as the following,

dPTdt=ksS−kdPT size 12{ { { ital "dP" rSub { size 8{T} } } over { ital "dt"} } =k rSub { size 8{s} } S - k rSub { size 8{d} } P rSub { size 8{T} } } {} (eq4)

dPdt=ksS+VapE(PT−P)Kap+(PT−P)−VipPKip+P−kdP size 12{ { { ital "dP"} over { ital "dt"} } =k rSub { size 8{s} } S+V rSub { size 8{ ital "ap"} } { {E ( P rSub { size 8{T} } - P ) } over {K rSub { size 8{ ital "ap"} } + ( P rSub { size 8{T} } - P ) } } - V rSub { size 8{ ital "ip"} } { {P} over {K rSub { size 8{ ital "ip"} } +P} } - k rSub { size 8{d} } P} {} (eq5)

dEdt=VaeP(1−E)Kae+(1−E)−VieEKie+E size 12{ { { ital "dE"} over { ital "dt"} } =V rSub { size 8{ ital "ae"} } { {P ( 1 - E ) } over {K rSub { size 8{ ital "ae"} } + ( 1 - E ) } } - V rSub { size 8{ ital "ie"} } { {E} over {K rSub { size 8{ ital "ie"} } +E} } } {} (eq6)

The first and second terms on the right hand side of the first equation are production and degradation rate of protein P with product rate k_s and degradation rate k_d, respectively. The first and second terms on the right hand side of the second and the third equation are dephosphorylation and phosphorylation reaction with maximal rate V_ap, V_ip, V_ae, and V_ie, kinetic constant K_ap, K_ip, K_ae, and K_ie, respectively.

The Figure 8A shows the rate plots with assumption that the phosphorylation process occurs much faster than protein synthesis, thus the total concentration PT size 12{P rSub { size 8{T} } } {} is in equilibrium state, PT=ksS/kd size 12{P rSub { size 8{T} } = {k rSub { size 8{s} } S} slash {k rSub { size 8{d} } } } {}; and the enzyme E activation is much faster than that of P, thus the enzyme E is also in equilibrium state, dE/dt=0 size 12{ { ital "dE"} slash { ital "dt"} =0} {}. The auto-phosphorylation is also a positive feedback loop, so that it is able to generate bistability as shown in Figure 8A (rate balance plot) and B (bifurcation diagram plot).

The rate balance plots (A), one parameter bifurcation diagrams (B) and two parameters bifurcation diagram (C) of autophosphorylation with additional transcription factor S. The reversible bistability could be generated by autophosphorylation. Red line and blue line represent production and degradation rates, respectively. used parameters: ks=1ks=1 size 12{k rSub { size 8{s} } =1} {}, kd=1kd=1 size 12{k rSub { size 8{d} } =1} {}, Vap=1Vap=1 size 12{V rSub { size 8{ ital "ap"} } =1} {}, Vip=0.5Vip=0.5 size 12{V rSub { size 8{ ital "ip"} } =0 "." 5} {}, Vae=1Vae=1 size 12{V rSub { size 8{ ital "ae"} } =1} {}, Vie=0.5Vie=0.5 size 12{V rSub { size 8{ ital "ie"} } =0 "." 5} {}, Kap=Kip=Kae=Kie=0.01