Stochasticity in protein levels drives colinearity of gene order in metabolic operons of Escherichia coli

Abstract

In bacterial genomes, gene order is not random. This is most evident when looking at operons, these often encoding enzymes involved in the same metabolic pathway or proteins from the same complex. Is gene order within operons nonrandom, however, and if so why? We examine this issue using metabolic operons as a case study. Using the metabolic network of Escherichia coli, we define the temporal order of reactions. We find a pronounced trend for genes to appear in operons in the same order as they are needed in metabolism (colinearity). This is paradoxical as, at steady state, enzymes abundance should be independent of order within the operon. We consider three extensions of the steady-state model that could potentially account for colinearity: (1) increased productivity associated with higher expression levels of the most 5' genes, (2) a faster metabolic processing immediately after up-regulation, and (3) metabolic stalling owing to stochastic protein loss. We establish the validity of these hypotheses by employing deterministic and stochastic models of enzyme kinetics. The stochastic stalling hypothesis correctly and uniquely predicts that colinearity is more pronounced both for lowly expressed operons and for genes that are not physically adjacent. The alternative models fail to find any support. These results support the view that stochasticity is a pervasive problem to a cell and that gene order evolution can be driven by the selective consequences of fluctuations in protein levels.

Figure 1. The concept of colinearity between positional sequence of operonic genes and enzymatic steps in a hypothetical metabolic pathway.

The first gene arrangement is perfectly colinear, whereas in the second case, two of the six intraoperonic gene pairs have colinear enzymatic orders; therefore, the degree of colinearity is one-third. m₀, m₁, …, m₄ denote metabolites.

Figure 2. Model of operon expression.

Reaction scheme for modelling gene expression of a polycistronic operon. We constructed a model of a four-gene operon and a linear metabolic pathway, containing 5 metabolites and 4 enzymes encoded by the operonic genes. Operon expression was modelled following the “read-through” operon model of Swain . Transcription is modelled as reversible binding of RNAP to promoter (D) with rates: f ₀ (association) and b ₀ (dissociation). Isomerization from closed to open complex and initiation of transcription are approximated as a first-order process (with rate k ₀). Only the leader region of the mRNA, M, is tracked in the model, which is made by transcribing polymerase, T, at rate v ₀. mRNA molecules are degraded with rate mf ₀, and diluted with rate D. Ribosomes compete with degradosomes for leader mRNA and bind reversibly (rates mf ₁ for association and mb ₁ for dissociation). Translation is started from the mC ² state with rate k ₁, which then frees M for further ribosome or degradasome binding. Enzymes are translated in the mT state with rate v ₁, and decay and dilute with rate α (α = D+k_degr). In case of a “read-through” operon, only the first cistron has a ribosome binding site; thus, a translating ribosome, mT ₂, releases enzyme E ₁ before translating the next protein (in the state mT ₃). The translation rate parameter was fine tuned to achieve realistic time delays between the appearances of consecutive gene products (approximately 60 s, on average). See Table S1 for parameter values and constants.

Figure 3. The relative advantage of colinearity decreases with elapsed time after operon induction according to deterministic simulations of the model.

Bars represent the relative excess of end product synthesized by a colinear operon (ABCD) compared to an operon with anti-colinear arrangement (DCBA). The operon is induced at t = 0. Cell generation time is set to 60 min. See Table S4 for simulation results with different parameter values.

Figure 4. Stochastic simulation results.

(A) Temporal fluctuations of enzyme molecule numbers in a lowly and a highly expressed operon (according to stochastic simulations of the model; see Table S1 for model parameters). (B) Calculated average relative advantage of colinearity for a highly and a lowly expressed operon after 50 cell generations (180,000 s) according to stochastic simulations of the model. Mean values of 1,000 repeated simulations and 95% confidence intervals are shown. Colinearity in a lowly expressed operon has a significantly higher advantage than in a highly expressed operon (p<2.2×10⁻¹⁶, Brunner-Munzel test; a rank-based heteroscedastic method to compare two groups [43]). Error bar indicates 95% confidence interval.

Figure 5. Distribution of the number of colinear gene pairs in randomised samples of highly (A) and lowly (B) expressed operons.

Red line indicates the number of colinear gene pairs observed in the E. coli genome ([A] 67/143, [B] 127/178 gene pairs). mRNA levels measured under aerobic glucose minimal condition were used to define highly and lowly expressed operons (see Methods). Intraoperonic gene orders were randomised 100,000 times.