Abortive Initiation as a Bottleneck for Transcription in the Early Drosophila Embryo

Abortive Initiation as a Bottleneck for Transcription in the Early Drosophila Embryo

Alexander S. Serov Corresponding author. E-mail: alexander.serov@pasteur.fr Alexander J. Levine Department of Chemistry And Biochemistry, University of California, Los Angeles, CA 90095, USA Madhav Mani Engineering Sciences and Applied Mathematics, Northwestern University, 2145 Sheridan Road, Evanston, IL 60208, USA

Gene transcription is a critical step in gene expression. The currently accepted physical model of transcription of MacDonald, Gibbs & Pipkin (1969) predicts the existence of a physical limit on the maximal rate of transcription, which does not depend on the transcribed gene. This limit appears as a result of polymerase “traffic jams” forming in the bulk of the 1D DNA chain at high polymerase concentrations. Recent experiments have, for the first time, allowed one to access live gene expression dynamics in the Drosophila fly embryo in vivo under the conditions of heavy polymerase load and test the predictions of the model (Garcia et al., 2013). Our analysis of the data shows that the maximal rate of transcription is indeed the same for the Hunchback, Snail and Knirps gap genes, and modified gene constructs in nuclear cycles 13, 14, but the experimentally observed value of the maximal transcription rate corresponds to only 40  of the one predicted by this model. We argue that such a decrease must be due to a slower polymerase elongation rate in the vicinity of the promoter region. This effectively shifts the bottleneck of transcription from the bulk to the promoter region of the gene. We suggest a quantitative explanation of the difference by taking into account abortive transcription initiation. Our calculations based on the independently measured abortive initiation constant in vitro confirm this hypothesis and find quantitative agreement with MS2 fluorescence live imaging data in the early fruit fly embryo. If our explanation is correct, then the transcription rate cannot be increased by replacing “slow codons” in the bulk with synonymous codons, and experimental efforts must be focused on the promoter region instead. This study extends our understanding of transcriptional regulation, re-examines physical constraints on the kinetics of transcription and re-evaluates the validity of the standard totally asymmetric simple exclusion process model of transcription. Several ideas on experimental validation of our hypothesis are also suggested.


Keywords: transcription, Drosophila fly embryonic development, pattern formation, k-TASEP, extended particles

1 Introduction

Pattern formation in the developing organism relies on precise spatiotemporal control of gene expression Gilbert (2013). Gene expression occurs as the result of a complex dynamical process, the first step of which is transcription Alberts et al. (2002). During transcription, a multi-protein polymerase complex moves sequentially down a strand of DNA Alberts et al. (2002) and produces an mRNA molecule, which, after various downstream processes, produces a particular protein. While it is clear in some situations that the observed rate of protein production is selected by the developmental program, that rate is also subject to physical limits resulting from its molecular mechanisms Bialek (2012).

Recently, techniques have been developed that allow us to track mRNA levels over time and across nuclear cycles of a fruit fly early development at the individual cell level Garcia et al. (2013), Bothma et al. (2014, 2015). Consequently, we have have access to at least parts of the fly’s developmental program (as monitored via the transcription rates of a handful of specific genes) from the first cell division of the fertilized egg through the formation of a complex and spatially patterned embryo.

These data provide a novel window on the unfolding of the fly’s developmental pattern. It is tempting to interpret the spatiotemporal changes in the rate of mRNA transcript production solely in terms of evolutionary optimization of embryonic development. The development process, however, relies on physical processes, which impose their own constraints. In this manuscript, we examine the mRNA production data for three genes, known to be a part of the head-to-tail patterning systems active during early fly development Gilbert (2013), with an eye towards determining whether the physical dynamics of transcription impose limitations on the developmental dynamics. To put this another way: we ask whether the observed developmental program runs in a dynamical regime where it saturates the physical constraints of the underlying biochemical machinery.

There is a distinct reason to imagine that the gene expression machinery does, in fact, saturate some dynamical upper bounds. On the one hand, there may exist selection pressure on the speed of embryo growth. On the other hand, there are speed limits to transcription inherent in the mechanism itself. In particular, the gene is a linear array of base pairs that must be read sequentially; the polymerase complex is doing the reading by stochastically making steps along the gene, and while doing so, it takes up a finite-length footprint on it (Fig. 1). A number of such polymerase complexes may simultaneously run along the gene Garcia et al. (2013), that in a first approximation can be assumed to interact only through steric repulsion (excluded volume). During maximal protein production, the number of such complexes is large enough (up to 100, Garcia et al. (2013)) that their combined footprint takes up a significant fraction of the gene’s length (up to 80 Garcia et al. (2013)). At such high coverages, one might well expect that the steric interactions between stochastically moving polymerase complexes set a maximal rate for gene expression. These basic dynamics are reminiscent of traffic jams observed on congested roadways, and has been well studied in a variety of contexts Schütz (2001), Schmittmann and Zia (1995), as we describe below.

Figure 1: Schematic representation of the abortive k-TASEP model with extended particles of length . The abortive k-TASEP model features an additional slow first site (shown in red) compared to the standard k-TASEP model. From the biological point of view, this model represents sequential elongation of polymerase complexes with footprint  along a gene of length 

Here we use earlier developed models to determine from theoretical considerations alone what that maximal mRNA production rate should be and ask whether this rate is indeed observed in the data. It is reasonable to imagine that Drosophila embryonic development is precisely the right place to look for nature to saturate this limiting value of protein production, since optimal production likely leads to the most efficient method for rapidly making a fly. Our findings are, at first, perplexing. We indeed find that for three different genes that are transcribed during development, the rates of transcriptional initiation production reach the same maximal value (within 15  of each other), suggesting the presence of a physically-imposed production speed limit that developmental biology saturates, at least in these extreme circumstances. The observed speed limit, however, is only around 40  that predicted by basic theoretical considerations, to be described below. A likely explanation of this result is that we have misidentified the rate-limiting step in the mRNA production. One cannot discount the existence of other bottlenecks in the system, and we note here that abortive initiation of transcription has been independently observed in polymerase dynamics in the proximity of the promoter in in vitro setups Revyakin et al. (2006), Margeat et al. (2006), Kapanidis et al. (2006). In effect, it was reported that the polymerase complex gets stuck for some time at the start of a gene. When we take into account this aspect of transcription and perform calculations based on an independently measured abortive initiation time constant Revyakin et al. (2006), we find that the traffic jams located at the start of the gene further limit the polymerase throughput and bring the predicted maximal rate inline with the observed one.

The remainder of this manuscript is organized as follows. First we briefly review the basic physics of traffic jams in one-dimensional models (see also Schmittmann and Zia, 1995, Zia et al., 2011, Chou et al., 2011). We then introduce a modified model based on our knowledge of abortive transcription initiation, and discuss our calculations and simulations. The modified model is a special case of what is known in the literature as an inhomogeneous TASEP model with slow (or defective) sites Zia et al. (2011), Dong et al. (2009), Chou and Lakatos (2004). In the next section we present the data and show that these cannot be satisfactory explained by the original model, but are quantitatively consistent with a transcriptional bottleneck associated with the initiation of the transcription process observed in the modified model. In the conclusions we discuss the implication of the demonstrated effect for developmental biology and for the statistical physics of non-equilibrium systems, and propose several directions for future research.

2 The statistical mechanics of traffic jams

Recognizing that the collective dynamics of many polymerases walking on DNA can be thought of as that of processive stochastic walkers with simple steric (excluded volume) interactions, we briefly review this long-standing model of non-equilibrium statistical physics. The physics of this one-dimensional model provides perhaps the simplest example of non-equilibrium steady states in a strongly interacting many-body system. This model, originally proposed by MacDonald and collaborators MacDonald et al. (1968), MacDonald and Gibbs (1969) and commonly referred to as a Totally Asymmetric Simple Exclusion Process (TASEP), has been studied as part of a general interest in the physics of driven diffusive systems Schmittmann and Zia (1995) and as a model for biologically relevant dynamics Zia et al. (2011), Chowdhury et al. (2005), Frey et al. (2004), with the applications to ribosome motion along mRNA. We refer the interested reader to Refs. Schmittmann and Zia (1995), Zia et al. (2011), Chou et al. (2011) for extensive reviews of the literature on TASEP.

In brief, it is essential to recall that the basic problem is defined by a particle injection rate  at the start (e.g., the left end of the one dimensional track on which the particles move) and an extraction rate at the right end (Fig. 1). The injection rate  may be interpreted as the attempt frequency to add a particle to the left, which is successful only when that space is unoccupied. The extraction rate  may be similarly interpreted. Once injected, particles move stochastically to the right. We will consider only the case of Poisson walkers taking single steps on a lattice so that, without interparticle interactions, a single stepping rate  determines both the mean velocity of the walker and all its velocity fluctuations van Kampen (1992). We also restrict our analysis to the case in which particles cannot be added or removed in the bulk (i.e. no Langmuir kinetics), although novel dynamics have been studied in this case Parmeggiani et al. (2003). This restriction is demanded by the experimentally observed processivity of the polymerase complex on the DNA. The model contains two more control parameters: the length of the track  (i.e., the total number of base pairs in the gene) and the length of the walking particles  (i.e., the number of lattice sites taken up by one walker). In most TASEP studies the walkers are assumed to be “point particles” (taking up one site) on a lattice constant in length, but, as we will see below, the finite length of the polymerase complex plays a significant role in understanding the collective dynamics of transcription in the developing embryo and cannot be ignored.

The main phenomenology of TASEP dynamics of point-like particles can be characterized by a phase diagram defined by the injection and extraction rates and respectively Schütz (2001), Zia et al. (2011). It is remarkable that this one-dimensional system with short-ranged interactions can develop long-range order such that the bulk properties of the system can be controlled by these boundary conditions at the endpoints. This is one of the many features distinguishing non-equilibrium steady states from the better understood thermodynamic phases of equilibrium states. A mean-field analysis of the steady states leads to the prediction of three distinct phases Schütz (2001). For  sufficiently small, one finds a stable shock front where the particle density jumps from a bulk value defined by the bulk jump rate  to the one controlled by the boundary extraction/injection rates . The position of this stable shock within the system depends on the relative values of  and . For sufficiently high values of both and , one observes a different phase, the maximal current regime. For  () there is also a high (low) density phase in which the current through the system is below the maximal value. These mean-field results are supported by exact solutions to the steady-state density profile obtained by Liggett Liggett (1975) and by numerical simulations.

In the case of extended particles (the k-TASEP model), the conclusions of the mean-field approach stay valid, but additional factors appear in now approximate formulas for the mean elongation rate  and particle density . For instance, in the low-density (LD) phase for low injection attempt rate  and high exit attempt rate  one obtains MacDonald and Gibbs (1969), Shaw et al. (2003):


When the injection attempt rate  increases to its critical value while  stays high (), the system experiences a phase transition to the maximal current phase that is characterized by the maximal particle density  and the maximal particle elongation rate :


For more details on the phase diagram of the k-TASEP model, the reader is referred to Refs Shaw et al. (2003), Schütz (2001). Assuming the footprint of a polymerase complex on the DNA is base pairs (bp) long (as determined by electron crystallography, (Rice1993, Selby et al., 1997, Poglitsch et al., 1999)), and its bulk jump rate on the DNA is around  (Garcia et al., 2013), in the maximal current regime, the theory predicts the maximal polymerase injection rate . New polymerase complexes cannot be recruited on the given gene faster than this rate, which is defined by how fast the injected particle clears the injection site (the first  sites of the lattice). On a gene  base pairs long (like Hunchback gene, (Garcia et al., 2013)), the maximal number of polymerases is .

3 Abortive k-TASEP Model

In this paper we suggest that abortive initiation of transcription can be the rate-limiting step of transcription and should be added to the classical k-TASEP model, when used to explain experimental data in early embryo development. One can expect that in such model the bottleneck of the transcription process will be located in the promoter region at the 5’ end of the gene, rather than in its bulk as in the original model. Abortive initiation is an effect wherein polymerase complexes loaded on a promoter region of a gene will several times engage in transcription cycles, but will abort and release short RNA products returning to the initial position (Kapanidis et al., 2006, Revyakin et al., 2006, Margeat et al., 2006). However once a product of 9 to 11 nucleotides is synthesized, the polymerase will break its interactions with the promoter and will proceed with mRNA synthesis. It is currently accepted that the transition from the abortive initiation state to processive RNA synthesis occurs through a “DNA scrunching” mechanism, when the mechanical stress stored in the “scrunched” DNA during abortive transcription allows the polymerase to break interactions with the promoter. However other mechanisms (namely, “transient excursions”, “inchworming”) have been proposed as well (Kapanidis et al., 2006).

To assess the effect of abortive initiation on the transcription rate, we modify the standard k-TASEP model by adding a slow site at the beginning of the lattice. This modification will cause the polymerase to stochastically pause at the first site imitating the effect observed in the experiment. Experimental studies have revealed that the time spent by a polymerase in the abortive initiation stage in the promoter region is an exponentially-distributed random variable with the average pause time  (Revyakin et al., 2006), which we will use in the simulations below. In what follows, we will call this modification of the k-TASEP model an abortive k-TASEP model.

The proposed abortive k-TASEP model is a special case of an inhomogeneous k-TASEP model with only one slow site located at the first site of the lattice. The inhomogeneous TASEP model has attracted significant interest in recent years, and some important analytical and numerical results have been obtained for point particles  (see Cook et al., 2013, Poker et al., 2015, Dhiman and Gupta, 2016, Xiao et al., 2016, and references therein), as well as for extended particles Shaw et al. (2003, 2004b, 2004a), Dong et al. (2007), Klumpp and Hwa (2008), Dong et al. (2009), Zia et al. (2011), Brackley et al. (2011). The standard way to treat a problem with an individual slow site like ours would consist in using mean-field approximations for the sublattices to the left and to the right of the slow site and stitching the two solutions together by the flux conservation equation on the slow site Kolomeisky (1998), Shaw et al. (2004a). However, when the only slow site is located at the beginning of the lattice, the problem can be easily reduced to the standard k-TASEP model (a similar problem with point particles was recently considered by Xiao et al. (2016)). Indeed, the slow site can be considered as a boundary condition for the remaining  sites of the lattice that can be treated as an ordinary -site k-TASEP model. Ignoring entry-rule-related boundary effects, one can then, in the first approximation, define the injection attempt rate  of the rest of the lattice through the pause time  at the first site: . For long genes , one can safely assume  and using Eq. (1), (2) get the following expressions for the normalized polymerase number in the steady state and normalized polymerase elongation rate:


These dependencies are shown in Figs. (a)a,LABEL:sub@fig:SimulationInjectionRate as solid red lines for the biologically relevant range of values of . Note that the normalized polymerase number  and the normalized injection rate  decrease from  and  for  to  and  for . For the experimentally reported value of  Revyakin et al. (2006), the corresponding ratios are  and .

Figure 2: Numerical simulations of the maximal polymerase number in the steady state  LABEL:sub@fig:SimulationMaxNumber, the polymerase injection rate  LABEL:sub@fig:SimulationInjectionRate, and the evolution of polymerase number  LABEL:sub@fig:SimulationMeanEvolution as functions of the mean time  spent in the abortive initiation regime. Each circle in LABEL:sub@fig:SimulationMaxNumber,LABEL:sub@fig:SimulationInjectionRate and each curve in LABEL:sub@fig:SimulationMeanEvolution represent the average of 20 simulations; the dashed black line marks the maximal loading rate  predicted by the original k-TASEP model (Eq. (2)). Independent experiments report that the mean duration of the abortive initiation phase is around  (Revyakin et al., 2006), which gives and , or and in the absolute values. Black arrows in LABEL:sub@fig:SimulationMeanEvolution indicate the established steady state of the simulation where  is evaluated.  corresponds to the slope of the curves in LABEL:sub@fig:SimulationMeanEvolution, with  labeled in the plot. Eq. (3) and plots LABEL:sub@fig:SimulationMaxNumber and LABEL:sub@fig:SimulationInjectionRate can be used to assess the mean time  polymerases spend in the abortive initiation regime based on experimental data for each gene, nc and construct. For instance, for the Hunchback bac construct (Fig. 3), one gets in nc 13 and in nc 14 from the mean , and in nc 13 and in nc 14 from the mean injection rate  (dashed green lines).

Note that while describes the steady-state of the system, the expression for  describes the current related to a shock wave of density propagating in the system before the steady state is established. The shock-wave description Schütz (2001) is possible thanks to the fact that the interface separating two distinct phases in the TASEP model is microscopic, i.e. the transition occurs on a scale of several lattice sites Janowsky and Lebowitz (1992, 1994). We verify that the shock-wave approach and the approximate formulas (3) ignoring boundary-layer effects correctly describe our abortive k-TASEP model by performing Monte Carlo simulations with the numerical values of parameters , and  provided above. We use the “complete entry, incremental exit” rule implying that polymerases are loaded on the gene only if the first  sites are free, while they leave the gene one site at a time with the bulk rate  (Fig. 1(Lakatos and Chou, 2003). To focus on the effect of the first slow site, we set the polymerase injection rate  to the value  much greater than the inverse time a polymerase spends on the first site of the lattice  for all  in the range 1–10 . The slow site is modeled as a stochastic process, in which at each time step  a polymerase can either remain at the first  sites or continue transcription with a certain finite probability. The escape probability  (the probability to leave the slow site) is calculated based on the chosen values of  and . The simulations were performed with the rejection-free kinetic Monte Carlo algorithm (KMC, an event-based algorithm, also known as residence-time or Bortz-Kalos-Lebowitz algorithm, see (Bortz et al., 1975)). The theoretical prediction (3) is found to perfectly describe simulation results (shown as blue circles in Figs (a)a,LABEL:sub@fig:SimulationInjectionRate), which validates our analytical approach. Fig. (c)c completes the simulation results illustrating the mean evolution of the system from the non-steady state to the steady state.

4 Comparison to Experimental Data

In the present section we compare the predictions of the abortive k-TASEP model to experimental results inferred from MS2 fluorescence data after FISH calibration, as described in (Garcia et al., 2013). Each data entry (fluorescence trace) contains the number of RNA polymerase II complexes  currently loaded (transcribing) on a selected gene as a function of time  recorded in one nucleus of one Drosophila fly embryo. The whole observation period is divided into nuclear cycles (nc), at the end of each of which the number of cells in the embryo doubles, and all polymerase complexes are unloaded from the gene. We define the time point  to correspond to the start of nc 10 (Garcia et al., 2013). In this analysis we will consider only the nuclear cycles 13 (starting around ) and 14 (starting around ), since the genes we analyze are mainly expressed in these nuclear cycles.

The following three genes were considered in these study: Hunchback (Hb), Snail (Sn) and Knirps (Kn). All considered genes have one primary and one shadow enhancer that influence the spatio-temporal dynamics of expression of these genes. The data presented below were obtained for each of these genes in three modifications: wild-type (labeled bac), a gene with no primary enhancer (no primary) and a gene with no shadow enhancer (no shadow). Table 1 summarizes the number of embryos and fluorescence traces analyzed for each gene-enhancer construct combination. The time resolution of the measurements is 37 .

Hunchback Snail Knirps
Bac 1965 (11) 547 (5) 326 (5)
No primary 979 (9) 607 (4) 84 (2)
No shadow 1279 (8) 159 (1) 60 (1)
Table 1: Number of fluorescence traces analyzed for each gene-enhancer construct combination. The number of embryos is given in the parentheses

Figs (c)c,LABEL:sub@fig:SlopeHistogramsNC14 show the distribution of polymerase loading rates  as measured in our experiments for bac Hb in nc 13 and 14. The maximal theoretical polymerase injection rate  predicted by the original k-TASEP model is marked by dashed black lines. One can see that  appears to envelope all observed experimental values, suggesting that Hb transcription in nc 13 and 14 cannot occur faster than . The steady-state polymerase number  for the Hunchback bac gene is shown in Figs (a)a,LABEL:sub@fig:MaxNumberHistogramsNC14, with the k-TASEP prediction  shown by a dotted black line. Once again we see that  provides a correct overall limit for all observed experimental values, suggesting that the original k-TASEP model indeed defines the upper boundary on the number of polymerase complexes loaded on the gene in nc 13 and 14. Moreover, Figs 4 and 5 demonstrate that the actual mean transcription rate  does not depend on a particular gene, construct or location of the nucleus on the anterior-to-posterior (AP) axis of the embryo, suggesting gene-independent nature of the transcription rate constraint.

Figure 3: LABEL:sub@fig:MaxNumberHistogramsNC13,LABEL:sub@fig:MaxNumberHistogramsNC14: Probability density function (PDF) of the steady-state number of polymerases  loaded on the HunchBack bac gene in nuclear cycles 13 LABEL:sub@fig:MaxNumberHistogramsNC13 and 14 LABEL:sub@fig:MaxNumberHistogramsNC14 (all AP positions combined). All data entries are grouped in bins of width 5, and the y-axis is proportional to the number of fluorescence traces falling in the current bin. The theoretically predicted upper limit  is shown by a dashed black line. The mean values of the distributions in the figure are: (nc 13) and (nc 14). LABEL:sub@fig:SlopeHistogramsNC13,LABEL:sub@fig:SlopeHistogramsNC14: PDF of initial polymerase injection rates  for the HunchBack bac enhancer construct in nuclear cycles 13 LABEL:sub@fig:SlopeHistogramsNC13 and 14 LABEL:sub@fig:SlopeHistogramsNC14 (all AP positions combined). All slopes are grouped in bins of width 1, and the y-axis is proportional to the number of fluorescence traces falling in the current bin. The dashed vertical line shows the maximal polymerase injection rate  as predicted by Eq. (2). The mean values of the shown distributions are: (nc 13) and (nc 14). In all figures the total number of traces shown is 682 for nc 13, and 1283 for nc 14. One can see that both  and  correctly predict the maximal observed values, although the average values are lower. The few outliers above  and  are likely due to noisy experimental data

Figure 4: Mean evolution of active polymerase number on the Hunchback (blue), Snail (yellow) and Knirps (red) genes. Each curve corresponds to a specific construct (bac, no primary or no shadow) of a particular gene (Hunchback, Snail or Knirps) in one of the nuclear cycles (13 or 14) and is averaged over the anterior-posterior (AP) position in the embryo. The curves are combined in the same plot to demonstrate that the initial injection rate  does not depend on any of these parameters (some noise is present). The thick black curve shows the average over all the parameters. The dashed straight line shows the maximal theoretical injection rate  as predicted by the standard k-TASEP model (Eq. (2)). The solid straight line shows the initial injection rate  as predicted by the abortive k-TASEP model for . Eqs (3) can be used to adjust the  parameter and fit the initial slope more accurately (see main text). All fluorescence traces used throughout this article required synchronization of the beginning of nuclear cycles between different embryos and nuclei. The desynchronization is due to recording started at different moments in different embryos and to the limited temporal resolution of the data. All traces shorter than 10  were also removed, since we could not reliably synchronize them with the other traces. The steady state of the system seems to be evolving with time, showing a gradual decrease in the active polymerase number. Under these conditions, we have taken the peak number of polymerases as an estimate of  in all presented data

Figure 5: Initial injection rate  as a function of the AP position for different constructs (bac, no primary, no shadow) of the Hunchback gene. The error bars show 95  confidence intervals as estimated by bootstrapping (Efron, 1979). Each curve was slightly shifted sideways to allow all error bars to be visible. The dashed line shows the injection rate predicted by the abortive TASEP model () for . This prediction fits well the experimental data lying slightly below. Narrow AP windows and low number of data sets available to us for the Snail and Knirps genes (Table 1) did not allow us to draw any reliable conclusions on the dependence of  on AP for those genes

It is however clear that the mean values of the distributions in Fig. 3 lie well below the predictions of the unmodified k-TASEP model MacDonald and Gibbs (1969), Shaw et al. (2003). The histograms in Fig. 3 clearly demonstrate a downshift from the predicted numbers, indicating that on average the 1D DNA chain is used significantly lower than its maximum polymerase flow capacity for a homogeneous lattice (with, for instance,  and  for Hb bac, Fig. 3). In this paper we interpret this decrease as that on average no traffic jams are observed in the bulk of the DNA chain for the selected genes and constructs, and that for them abortive transcription initiation is the real bottleneck of transcription (other interpretations are also possible, and we return to this discussion later). Our speculations are supported by the abortive k-TASEP model that for independently measured value of  (Revyakin et al., 2006) predicts a similar decrease:  and . The slight difference between the predicted and the observed values of  and  may be (besides other reasons such as biological variability) attributed to a slightly different value of  for Drosophila embryos in vivo compared to previous in vitro experiments in Ref. Revyakin et al. (2006). If one assumes the validity of the abortive k-TASEP model, the distributions in Fig. 3 can be explained by the following values of : for nc 13 and for nc 14 from the polymerase number , and for nc 13 and for nc 14 from the mean injection rate . These values are shown by dashed green lines in Figs (a)a,LABEL:sub@fig:SimulationInjectionRate.

Figure 6: Relative standard deviation of the steady state number of polymerases () and normalized initial injection slope () calculated numerically as a function of mean abortive transcription initiation duration . Experimentally calculated values for the Hunchback bac construct are shown by dashed horizontal lines:  (nc 13),  (nc 14),  (nc 13),  (nc 14). The experimental lines show a constant value, since  dependence was not available from the experiment. Note that the noise level of the experimental data is much higher than the noise predicted by the numerical model, which suggests the existence of additional noise sources in experimental data (see main text)

To finalize the comparison between the simulations and the experimental data, it is interesting to compare the widths of the calculated distributions for  and  to the widths of experimental distributions in Fig. 3. Fig. 6 shows standard deviations (STD) of and calculated for different values of mean abortive initiation time  and plotted along the experimentally measured values. It can be seen that experimental distributions are about 3 times wider than the calculated ones, which suggests that the noise sources included into the numerical abortive k-TASEP model (namely, stochastic elongation and stochastic abortive initiation duration) are responsible for only about one third of the variability of the experimental values. Other noise sources (such as stochastic promoter activation, heterogeneous or dynamic elongation rates, instrumental noise or slope fitting) must account for the rest of the width of the distributions (for further discussion of noise sources in protein synthesis see (Tkacik et al., 2008)).

5 Discussion

Our analysis of transcriptional dynamics in the Drosophila fly embryo exhibits a common value of steady-state polymerase number  and polymerase injection rate  in the beginning of nuclear cycles 13 and 14 for several analyzed genes and constructs. Although such universal character can be explained by the standard k-TASEP model, the observed values of  and  are significantly lower than the theoretically predicted ones. In the present article we have demonstrated that this difference in slopes can be explained by the abortive initiation process that effectively decreases the polymerase elongation rate in the promoter region of the gene. We have modified the k-TASEP model adding a slow first site, and provided simple formulas and performed Monte Carlo simulations to estimate  and  in it. Our results for independently measured abortive transcription initiation duration  give values much closer to the experimental ones, and suggest that the promoter region may be the real bottleneck of the transcription process rather than polymerase “traffic jams” downstream of the gene.

It is however important to mention that our hypothesis is not the only possible explanation for the observed phenomenon. It seems to be impossible to explain the observed decrease in  and  without slow sites near the start of the gene, but the distribution of slow sites may take different forms. We have only analyzed the simplest option of only one first slow site. Other options may include several distributed slow sites, clusters of slow sites Chou and Lakatos (2004), or inhomogeneous jump rates (quenched disorder) near the start of the gene Shaw et al. (2004b), Zia et al. (2011), Shaw et al. (2003). Interestingly, the effect of quenched disorder on the particle current in the TASEP model has been shown to mainly reduce to the effect of several slowest sites Shaw et al. (2004b, 2003). Although the physical description of systems with different distributions of slow sites would be similar, the biological processes behind the slow sites may be completely different. For instance, the quenched disorder model describes the fact that different nucleotides for the building mRNA chain may have different concentrations in the cell, so that the polymerase complex has to wait longer for some nucleotides than for the others thus modulating the elongation speed. It is interesting to mention here that, for example, for some of the proteins of Escherichia coli, the slow codons were reported to be preferentially located within the first 25 codons Chen and Inouye (1990). Regardless of the exact distribution of slow sites, an important practical conclusion of this study is that if the bottleneck of transcription is indeed located in the promoter region, the overall transcription rate cannot be increased by modifications in the bulk of the gene, such as replacement of slow codons with synonymous ones Shaw et al. (2004b), Zia et al. (2011), Shaw et al. (2003). Experimental research of higher transcription rates must then focus on acceleration of the initiation dynamic or on substitution of slow codons in the promoter region of the gene.

From the physical point of view, it seems to be impossible to distinguish different configurations of slow sites having only access to the steady-state dynamics of the experimental system. However, one can probably search for the answer in the non-steady state dynamics, or the slow evolution of the steady-states seen in Fig. 4. For instance, if slow sites are distributed along the whole gene, the propagating density wave must slow down each time a new, slower site is discovered. Surprisingly, the curves in Fig. 4 don’t seem to exhibit this gradual or stepwise slow-down, suggesting again that the slowest sites are located near the promoter region. This feature, as well as a local decrease of transcription speed around 1–2  into nc 13 (Fig. 4, middle pannel), or the gradual decrease of polymerase density after the maximum has been reached (Fig. 4), cannot be explained by the abortive k-TASEP model, and require further investigation.

From the point of view of molecular biology, one could conceive additional experiments aimed at identification of the role of the abortive transcription process in the decrease of the transcription rate. On the one hand, if one could alter the molecular mechanisms behind the abortive initiation and slow it down, the corresponding decrease in the transcription rate and steady-state polymerase number could be observed experimentally and compared to the predictions of Eq. (3) for the new values of . If on the other hand, the abortive initiation could be accelerated, at some point one can expect that the transcription bottleneck may shift to other involved mechanisms, so that the transcription rate may be determined by polymerase “traffic jams” in the bulk of the gene or, for instance, promoter activation kinetics. In either case, the curves in Fig. 4 must reflect the changes.

We have also compared the width of the distributions of  and  in our simulations and in the experimental data. Much wider experimental distributions indicate that the noise sources taken into account in our abortive k-TASEP model cannot satisfactory explain the full variability of experimental data, and that additional noise sources should be investigated in future studies. As an interesting development, one could try to extract additional information about the transcriptional bottlenecks and noise sources of the system by analyzing the noise, identifying noise signatures seen in the experimental data and shapes of the experimental distributions (i.e., in Fig. 3) and comparing them to those predicted by the abortive k-TASEP model or other noise sources.

Finally, the question of whether the system of polymerases on the DNA actually reaches a steady non-equilibrium state under physiological conditions of nc 13 and 14 also requires additional investigation. For instance, Fig. 4 is inconclusive about whether the transcription rate stays constant for some time after reaching the peak value and before starting to decrease. If it is not the case, the validity of the steady-state k-TASEP model as a model of transcription under physiological conditions should be critically re-evaluated.


  • Alberts et al. (2002) Alberts, B., Johnson, A., Lewis, J., Raff, M., Roberts, K., and Walter, P., 2002. Molecular Biology of the Cell. Garland Science, New York, 4 edition.
  • Bialek (2012) Bialek, W., 2012. Biophysics: Searching for Principles. Princeton University Press, Princeton, Oxford.
  • Bortz et al. (1975) Bortz, A., Kalos, M., and Lebowitz, J. L., 1975. A new algorithm for Monte Carlo simulation of Ising spin systems. J. Comput. Phys., 17(1):10–18. http://dx.doi.org/10.1016/0021-9991(75)90060-1.
  • Bothma et al. (2014) Bothma, J. P., Garcia, H. G., Esposito, E., Schlissel, G., Gregor, T., and Levine, M., 2014. Dynamic regulation of eve stripe 2 expression reveals transcriptional bursts in living Drosophila embryos. Proc. Natl. Acad. Sci. U. S. A., 111(29):10598–603. http://dx.doi.org/10.1073/pnas.1410022111.
  • Bothma et al. (2015) Bothma, J. P., Garcia, H. G., Ng, S., Perry, M. W., Gregor, T., and Levine, M., 2015. Enhancer additivity and non-additivity are determined by enhancer strength in the Drosophila embryo. Elife, 4(AUGUST2015):1–14. http://dx.doi.org/10.7554/eLife.07956.001.
  • Brackley et al. (2011) Brackley, C. A., Romano, M. C., and Thiel, M., 2011. The Dynamics of Supply and Demand in mRNA Translation. PLoS Comput. Biol., 7(10):e1002203. http://dx.doi.org/10.1371/journal.pcbi.1002203.
  • Chen and Inouye (1990) Chen, G. F. T. and Inouye, M., 1990. Suppression of the negative effect of minor arginine codons on gene expression; preferential usage of minor codons within the first 25 codons of the Escherichia coli genes. Nucleic Acids Res., 18(6):1465–1473. http://dx.doi.org/10.1093/nar/18.6.1465.
  • Chou and Lakatos (2004) Chou, T. and Lakatos, G., 2004. Clustered bottlenecks in mRNA translation and protein synthesis. Phys. Rev. Lett., 93(19):1–4. http://dx.doi.org/10.1103/PhysRevLett.93.198101.
  • Chou et al. (2011) Chou, T., Mallick, K., and Zia, R. K. P., 2011. Non-equilibrium statistical mechanics: from a paradigmatic model to biological transport. Reports Prog. Phys., 74(11):116601. http://dx.doi.org/10.1088/0034-4885/74/11/116601.
  • Chowdhury et al. (2005) Chowdhury, D., Schadschneider, A., and Nishinari, K., 2005. Physics of transport and traffic phenomena in biology: From molecular motors and cells to organisms. Phys. Life Rev., 2(4):318–352. http://dx.doi.org/10.1016/j.plrev.2005.09.001.
  • Cook et al. (2013) Cook, L. J., Dong, J. J., and Lafleur, A., 2013. Interplay between finite resources and a local defect in an asymmetric simple exclusion process. Phys. Rev. E - Stat. Nonlinear, Soft Matter Phys., 88(4):1–8. http://dx.doi.org/10.1103/PhysRevE.88.042127.
  • Dhiman and Gupta (2016) Dhiman, I. and Gupta, A. K., 2016. Origin and dynamics of a bottleneck-induced shock in a two-channel exclusion process. Phys. Lett. Sect. A Gen. At. Solid State Phys., 380(24):2038–2044. http://dx.doi.org/10.1016/j.physleta.2016.04.031.
  • Dong et al. (2007) Dong, J. J., Schmittmann, B., and Zia, R. K. P., 2007. Inhomogeneous exclusion processes with extended objects: The effect of defect locations. Phys. Rev. E - Stat. Nonlinear, Soft Matter Phys., 76(5):31–33. http://dx.doi.org/10.1103/PhysRevE.76.051113.
  • Dong et al. (2009) Dong, J. J., Zia, R. K. P., and Schmittmann, B., 2009. Understanding the edge effect in TASEP with mean-field theoretic approaches. J. Phys. A Math. Theor., 42(1):015002. http://dx.doi.org/10.1088/1751-8113/42/1/015002.
  • Efron (1979) Efron, B., 1979. Bootstrap Methods: Another Look at the Jackknife. Ann. Stat., 7(1):1–26. http://dx.doi.org/10.1214/aos/1176344552.
  • Frey et al. (2004) Frey, E., Parmeggiani, A., and Franosch, T., 2004. Collective phenomena in intracellular processes. Genome informatics, 15(1):46–55. http://dx.doi.org/10.11234/gi1990.15.46.
  • Garcia et al. (2013) Garcia, H. G., Tikhonov, M., Lin, A., and Gregor, T., 2013. Quantitative imaging of transcription in living Drosophila embryos links polymerase activity to patterning. Curr. Biol., 23(21):2140–5. http://dx.doi.org/10.1016/j.cub.2013.08.054.
  • Gilbert (2013) Gilbert, S. F., 2013. Developmental Biology. Sinauer Associates, Inc., Sunderland, MA, 10 edition.
  • Janowsky and Lebowitz (1994) Janowsky, S. A. and Lebowitz, J. L., 1994. Exact results for the asymmetric simple exclusion process with a blockage. J. Stat. Phys., 77(1-2):35–51. http://dx.doi.org/10.1007/BF02186831.
  • Janowsky and Lebowitz (1992) Janowsky, S. A. and Lebowitz, J. L., 1992. Finite-size effects and shock fluctuations in the asymmetric simple-exclusion process. Phys. Rev. A, 45(2):618–625. http://dx.doi.org/10.1103/PhysRevA.45.618.
  • Kapanidis et al. (2006) Kapanidis, A. N., Margeat, E., Ho, S. O., Kortkhonjia, E., Weiss, S., and Ebright, R. H., 2006. Initial transcription by RNA polymerase proceeds through a DNA-scrunching mechanism. Science (80-. )., 314(5802):1144–1147. http://dx.doi.org/10.1126/science.1131399.
  • Klumpp and Hwa (2008) Klumpp, S. and Hwa, T., 2008. Stochasticity and traffic jams in the transcription of ribosomal RNA: Intriguing role of termination and antitermination. Proc. Natl. Acad. Sci. U. S. A., 105(47):18159–18164. http://dx.doi.org/10.1073/pnas.0806084105.
  • Kolomeisky (1998) Kolomeisky, A. B., 1998. Asymmetric simple exclusion model with local inhomogeneity. J. Phys. A. Math. Gen., 31(4):1153–1164. http://dx.doi.org/10.1088/0305-4470/31/4/006.
  • Lakatos and Chou (2003) Lakatos, G. and Chou, T., 2003. Totally asymmetric exclusion processes with particles of arbitrary size. J. Phys. A. Math. Gen., 36(8):2027–2041. http://dx.doi.org/10.1088/0305-4470/36/8/302.
  • Liggett (1975) Liggett, T. M., 1975. Ergodic theorems for the asymmetric simple exclusion process. Trans. Am. Math. Soc., 213:237–61. http://dx.doi.org/10.1090/S0002-9947-1975-0410986-7.
  • MacDonald and Gibbs (1969) MacDonald, C. T. and Gibbs, J. H., 1969. Concerning the kinetics of polypeptide synthesis on polyribosomes. Biopolymers, 7(5):707–725. http://dx.doi.org/10.1002/bip.1969.360070508.
  • MacDonald et al. (1968) MacDonald, C. T., Gibbs, J. H., and Pipkin, A. C., 1968. Kinetics of biopolymerization on nucleic acid templates. Biopolymers, 6(1):1–25. http://dx.doi.org/10.1002/bip.1968.360060102.
  • Margeat et al. (2006) Margeat, E., Kapanidis, A. N., Tinnefeld, P., Wang, Y., Mukhopadhyay, J., Ebright, R. H., and Weiss, S., 2006. Direct observation of abortive initiation and promoter escape within single immobilized transcription complexes. Biophys. J., 90(4):1419–1431. http://dx.doi.org/10.1529/biophysj.105.069252.
  • Parmeggiani et al. (2003) Parmeggiani, A., Franosch, T., and Frey, E., 2003. Phase coexistence in driven one-dimensional transport. Phys. Rev. Lett., 90(February):086601. http://dx.doi.org/10.1103/PhysRevLett.90.086601.
  • Poglitsch et al. (1999) Poglitsch, C. L., Meredith, G. D., Gnatt, a. L., Jensen, G. J., Chang, W. H., Fu, J., and Kornberg, R. D., 1999. Electron crystal structure of an RNA polymerase II transcription elongation complex. Cell, 98(6):791–798. http://dx.doi.org/S0092-8674(00)81513-5[pii].
  • Poker et al. (2015) Poker, G., Margaliot, M., and Tuller, T., 2015. Sensitivity of mRNA Translation. Sci. Rep., 5:12795. http://dx.doi.org/10.1038/srep12795.
  • Revyakin et al. (2006) Revyakin, A., Liu, C., Ebright, R. H., and Strick, T. R., 2006. Abortive Initiation and Productive Initiation by RNA Polymerase Involve DNA Scrunching. Science (80-. )., 314(5802):1139–43. http://dx.doi.org/10.1126/science.1131398.
  • Schmittmann and Zia (1995) Schmittmann, B. and Zia, R. Statistical Mechanics of Driven Diffusive Systems. In Domb, C. and Lebowitz, J. L., editors, Phase Transitions Crit. Phenomena. Vol. 17, pages 3–214. Academic Press, London, 1995. http://dx.doi.org/10.1016/S1062-7901(06)80014-5.
  • Schütz (2001) Schütz, G. M. Exactly Solvable Models for Many-Body Systems Far from Equilibrium. In Domb, C. and Lebowitz, J. L., editors, Phase Transitions Crit. Phenom., chapter 1, pages 3–255. Academic Press, San Diego, 2001. http://dx.doi.org/10.1016/S1062-7901(01)80015-X.
  • Selby et al. (1997) Selby, C. P., Drapkin, R., Reinberg, D., and Sancar, A., 1997. RNA polymerase II stalled at a thymine dimer: footprint and effect on excision repair. Nucleic Acids Res., 25(4):787–793. http://dx.doi.org/10.1093/nar/25.4.787.
  • Shaw et al. (2003) Shaw, L. B., Zia, R. K. P., and Lee, K. H., 2003. Totally asymmetric exclusion process with extended objects: A model for protein synthesis. Phys. Rev. E, 68(2):021910. http://dx.doi.org/10.1103/PhysRevE.68.021910.
  • Shaw et al. (2004a) Shaw, L. B., Kolomeisky, A. B., and Lee, K. H., 2004a. Local inhomogeneity in asymmetric simple exclusion processes with extended objects. J. Phys. A, 37(6):2105–2113. http://dx.doi.org/10.1088/0305-4470/37/6/010.
  • Shaw et al. (2004b) Shaw, L. B., Sethna, J. P., and Lee, K. H., 2004b. Mean-field approaches to the totally asymmetric exclusion process with quenched disorder and large particles. Phys. Rev. E - Stat. Nonlinear, Soft Matter Phys., 70(2 1):1–7. http://dx.doi.org/10.1103/PhysRevE.70.021901.
  • Tkacik et al. (2008) Tkacik, G., Gregor, T., and Bialek, W., 2008. The role of input noise in transcriptional regulation. PLoS One, 3(7):e2774. http://dx.doi.org/10.1371/journal.pone.0002774.
  • van Kampen (1992) van Kampen, N., 1992. Stochastic Processes in Physics and Chemistry. North-Holland Personal Library, Amsterdam.
  • Xiao et al. (2016) Xiao, S., Chen, X., and Liu, Y., 2016. Totally asymmetric simple exclusion process with a single defect site on boundaries. Int. J. Mod. Phys. B, 30(14):1650083. http://dx.doi.org/10.1142/S0217979216500831.
  • Zia et al. (2011) Zia, R. K. P., Dong, J. J., and Schmittmann, B., 2011. Modeling Translation in Protein Synthesis with TASEP: A Tutorial and Recent Developments. J. Stat. Phys., 144(2):405–428. http://dx.doi.org/10.1007/s10955-011-0183-1.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description