# Deciphering mRNA Sequence Determinants of Protein Production Rate

Juraj Szavits-Nossan, Luca Ciandrini, M. Carmen Romano

► To cite this version:

Juraj Szavits-Nossan, Luca Ciandrini, M. Carmen Romano. Deciphering mRNA Sequence Determinants of Protein Production Rate. Physical Review Letters, 2018, 120, pp.128101. 10.1103/Phys-RevLett.120.128101. hal-01585251

HAL Id: hal-01585251

https://hal.science/hal-01585251

Submitted on 2 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Deciphering mRNA Sequence Determinants of Protein Production Rate

Juraj Szavits-Nossan*

SUPA, School of Physics and Astronomy, University of Edinburgh,
Peter Guthrie Tait Road, Edinburgh EH9 3FD, United Kingdom

Luca Ciandrini

DIMNP UMR 5235, Université de Montpellier and CNRS, F-34095, Montpellier, France and
Laboratoire Charles Coulomb UMR5221, Université de Montpellier and CNRS, F-34095, Montpellier, France

M. Carmen Romano

SUPA, Institute for Complex Systems and Mathematical Biology,
Department of Physics, Aberdeen AB24 3UE, United Kingdom and
Institute of Medical Sciences, University of Aberdeen,
Foresterhill, Aberdeen AB24 3FX, United Kingdom

One of the greatest challenges in biophysical models of translation is to identify coding sequences features that affect the rate of translation and therefore the overall protein production in the cell. We propose an analytic method to solve a translation model based on the inhomogeneous totally asymmetric simple exclusion process, which allows us to unveil simple design principles of nucleotide sequences determining protein production rates. Our solution shows an excellent agreement when compared to numerical genome-wide simulations of *S. cerevisiae* transcript sequences and predicts that the first 10 codons, which is the ribosome footprint length on the mRNA, together with the value of the initiation rate, are the main determinants of protein production rate under physiological conditions. Finally, we interpret the obtained analytic results based on the evolutionary role of codons' choice for regulating translation rates and ribosome densities.

PACS numbers: 87.16.aj, 87.10.Mn, 05.60.-k

Translation is one of the major steps in protein biosynthesis. During this process, the nucleotide sequence of a messenger RNA (mRNA) is translated into a functional protein. Each nucleotide triplet, called codon, codes for a specific amino acid, the proteins' building block. There is experimental evidence that the rate at which a certain mRNA is translated depends on its specific codon sequence [1-4], especially in the case of eukaryotes. Identifying sequence features that determine protein production rate, also commonly referred to as translation rate or efficiency, is a fundamental open question in molecular biology [1, 3].

Translation is performed by molecular motors called ribosomes, which move unidirectionally along the mRNA. The amino acids are delivered to the ribosome by molecules called transfer RNAs (tRNAs), which are specific to the codon and the amino acid they deliver (Fig. 1). The dwelling time of a ribosome on a specific codon depends primarily on the abundance of the corresponding tRNA [5, 6]. Understanding how codon sequences determine protein production rates can potentially unlock many synthetic applications [7, 8].

The standard biophysical model of translation is known as the totally asymmetric simple exclusion process (TASEP), which captures the concurrent motion of ribosomes on the mRNA [9, 10]. In this model ribosomes progress along the mRNA codon by codon, provided that the codon a ribosome moves onto is not occupied by another ribosome. The protein production rate

<figure><img src="image_2.png" /><figcaption>FIG. 1: Sketch of the mRNA translation process involving initiation (a), elongation consisting of tRNA delivery (b) and translocation (c), followed by termination (d). We emphasise that this is an oversimplified scheme of the process and that actual ribosomes cover $l = 10$ codons. At each elongation step, the ribosome receives an amino acid from a tRNA that matches the codon occupied by the A-site of the ribosome (the “reading” site). After the amino acid is added to the growing polypeptide, the ribosome translocates one codon forward, and the process is repeated.</figcaption></figure>

can then be identified as the ribosomal current of this
driven lattice gas. Due to the net current, TASEP is not in equilibrium and its steady state is only known for a few special cases [11, 12]. Unfortunately, most biologically relevant variants of TASEP can be studied only numerically [13–15], and efficient methods of exploring a large number of parameters are lacking. The steady state of the TASEP with non-uniform hopping rates is a long outstanding problem in nonequilibrium statistical physics [15, 16]. An analytic prediction of protein production rates is also needed in order to interpret recently developed ribosome profiling experiments that are capable of monitoring ribosome positions along the mRNA [17]. Such analytical prediction is missing in previous models of ribosome dynamics implicitly or explicitly based on the TASEP [1, 2, 19, 20, 22], and it is key to decipher sequence determinants of protein production rates.

In this Letter, we develop a versatile analytic method to solve TASEP-based models of translation. Our analytic approach, integrated with simulations and experimental data, allows us to efficiently identify the main features of mRNA codon sequence that determine the rate of protein production. The genome-wide comparison of our analytic predictions with numerical simulations shows an excellent agreement for the model organism S. cerevisiae (baker's yeast).

**Stochastic model of mRNA translation.** We focus on the model for translation introduced in [1, 23, 24]. The mRNA is represented by a one-dimensional lattice consisting of L discrete sites (codons), where site 1 designates the start codon. Ribosomes are represented by particles that occupy $\ell = 10$ lattice sites, which is the ribosome footprint length measured in ribosome profiling experiments [17]. We identify the position $1 \le i \le L$ of a ribosome with the position of its A-site, which is located $d = 5$ lattice sites from the trailing end of the ribosome (Fig. 1a). A ribosome that waits for a tRNA at position i is labelled by $1_i$ and a ribosome that has already received the correct tRNA and is ready to move is labelled by $2_i$. The set of labels of all translating ribosomes on the lattice is called a configuration C of the system. For example, $C = 112_{13}$ denotes a lattice with 2 ribosomes, one at site 1 waiting for a tRNA and another one at site 13 that has already received the correct tRNA.

Ribosomes initiate translation at rate $\alpha$ by binding to the mRNA so that their A-site is at the start codon, provided the sites $1, \dots, \ell - d + 1$ are empty (Fig. 1a). A ribosome at site i makes the transition $1_i \to 2_i$ at rate $k_i$ dependent on tRNA abundances (Fig. 1b) and moves one site forward at rate $\gamma$ (Fig. 1c). Due to steric interactions between particles, the particle at site i can move only if there is no ribosome at site $i + \ell$. Termination occurs when the particle at site L receives the last amino acid, releases the final protein and detaches from the mRNA, which we integrate into a single step occurring at rate $\beta$ (Fig. 1d).

Our main goal is to compute the rate of protein production as a function of the parameters of the model, which are $\alpha, \beta, \gamma$ and $k_i$ for each of the L codons used. We assume that translation takes place under steady-state conditions, so that the ribosomal current is constant along the mRNA and is equal to the rate of protein production, which in turn is equal to rate at which ribosomes load onto the mRNA and initiate translation. To this end we define the codon occupation number $\tau_i$ to be equal to 1 if the i-th codon is occupied by an A-site, and 0 if otherwise. By definition, the exact steady-state ribosomal current J then reads

$$J = \alpha \sum_C \left[ \prod_{i=1}^{\ell} (1 - \tau_i(C)) \right] P(C), \quad (1)$$

where $\tau_i(C)$ denotes the i-th codon occupation number for configuration C, $P(C)$ denotes the steady-state probability that the lattice is in configuration C, and the summation goes over all configurations C. Other quantities of interest that we compute are the local and total particle densities $\rho_i = \langle \tau_i \rangle$ and $\rho = (1/L) \sum_{i=1}^L \rho_i$, respectively.

**Series expansion method for computing $P(C)$**. In order to find $P(C)$ one has to solve the steady-state master equation $MP = 0$, where P is a column vector whose $N$ elements are the steady-state probabilities $P(C)$ of being in configuration C, and $N$ denotes the total number of configurations. The transition rate matrix $M_{C,C'}$ is given by $W_{C' \to C}$ for $C \neq C'$ and $-e(C)$ for $C = C'$, where $W_{C' \to C}$ is the transition rate from $C'$ to C and $e(C) = \sum_{C'' \neq C} W_{C \to C''}$ is the total exit rate from C. The exact solution of the master equation can be formally written as

$$P(C) = \frac{\det M^{(p,p)}}{\sum_{q=1}^{N} \det M^{(q,q)}}, \quad (2)$$

where $p$ is the position of configuration C in the column vector P and $\det M^{(p,p)}$ is a determinant of the matrix obtained by removing $p$-th row and $p$-th column from M (see Supplemental Material for details). Unfortunately, calculating this determinant is feasible only for unrealistically small system sizes.

To circumvent this problem, we exploit the fact that $\det M^{(p,p)}$ is a multivariate polynomial in the variables $\alpha, k_1, \dots, k_L, \beta$ and $\gamma$ [25]. In the biological literature it is often assumed that the initiation rate $\alpha$ is a major limiting step of the translation process, mainly determined by the presence of secondary structures [26, 27]. For this reason we use $\alpha$ as an expansion parameter and we assume that $\alpha \ll \gamma, \beta, k_1, \dots, k_L$. By collecting terms with the same power of $\alpha$, we can rewrite $P(C)$ as an univariate polynomial $f(C)$ in the variable $\alpha$ with unknown coefficients $f_n(C)$

$$P(C) = \frac{f(C)}{\sum_C f(C)}, \quad f(C) = \sum_{n=0}^{K(C)} f_n(C) \alpha^n, \quad (3)$$
where the coefficients $f_n(C)$ depend on the transition rates $k_1, \dots, k_L$, $\gamma$ and $\beta$ (in order to ease the notation, we leave out the explicit dependence on those parameters). Since we expect that $f(C)$ can be well approximated by the first few terms, the value of $K(C)$ in Eq. (3) is irrelevant in our study.

In order to find the unknown coefficients $f_n(C)$, we insert Eq. (3) into the master equation $MP = 0$, collect all the terms with the same n-th power of $\alpha$ and equate their sum to zero. For a given power $n$, the resulting equation is similar to the original master equation in which $P(C)$ is replaced by $f_n(C)$ unless the coefficient multiplying $P(C)$ is $\alpha$, in which case $P(C)$ is replaced by $f_{n-1}(C)$. Starting with $n=0$, we note the equations for the coefficients $f_0(C)$ have the same form as the original master equation, but with $\alpha=0$, leading to the trivial solution $f_0(C) = 0$ if $C \neq \emptyset$. Since the terms $f(C)$ are not normalized, we have the freedom to choose any value for $f_0(\emptyset)$, which we set to 1.

For $n \ge 1$, the equation for $f_n(C)$ reads

$$e_0(C)f_n(C) + f_{n-1}(C) \sum_{C'} I_{C,C'} = \sum_{C'} I_{C',C} f_{n-1}(C') + \sum_{C'} (1 - I_{C',C}) W_{C' \to C} f_n(C'), \quad (4)$$

where $I_{C,C'} = 1$ if $W_{C \to C'} = \alpha$ and is 0 otherwise, and $e_0(C) = \sum_{C'} (1 - I_{C,C'}) W_{C \to C'}$ is the total exit rate from $C$ excluding the rate $\alpha$.

A key observation in our analysis is that $f_n(C) = 0$ whenever the number of particles in $C$ is larger than $n$. This follows from the result $f_0(C) = 1(0)$ if $C = \emptyset$ ($C \neq \emptyset$) in conjunction with the hierarchical structure of Eq. (4), which connects configurations differing in the number of particles by no more than one (a proof for $n=1$ is presented in the Supplemental Material). This allows us to write Eq. (4) taking into account configurations with only one particle and discarding all the others, which yields

$$f_1(1_1) = \frac{1}{k_1} f_0(\emptyset), \quad f_1(1_L) = \frac{\beta}{k_L} f_1(2_L) \quad (5a)$$

$$f_1(1_i) = \frac{\gamma}{k_i} f_1(2_{i-1}), \quad 2 \le i \le L, \quad (5b)$$

$$f_1(2_i) = \frac{k_i}{\gamma} f_1(1_i), \quad 1 \le i \le L-1. \quad (5c)$$

The solution to Eqs. (5) is given by

$$f_1(1_i) = \frac{1}{k_i}, \quad f_1(2_i) = \begin{cases} \frac{1}{\gamma} & i = 1, \dots, L-1 \\ \frac{1}{\beta} & i = L \end{cases} \quad (6)$$

The equations for $f_2(C)$ involving configurations with one and two particles are presented in the Supplemental Material.

Once we determine the coefficients $f_n(C)$ up to a desired order $n$, we can compute the steady-state average

of any observable $O(C)$ by inserting $P(C)$ from Eq. (3) and expanding $\langle O \rangle$ around $\alpha = 0$,

$$\langle O(C) \rangle = \frac{\sum_{n=0}^{K} \sum_{C} O(C) f_n(C) \alpha^n}{\sum_{n=0}^{K} \sum_{C} f_n(C) \alpha^n} = \sum_{n=0}^{\infty} c_n \alpha^n. \quad (7)$$

For example, the first three coefficients $c_0, c_1$ and $c_2$ are given by

$$c_0 = \frac{a_0}{b_0}, \quad c_1 = \frac{a_1 - c_0 b_1}{b_0}, \quad c_2 = \frac{a_2 - c_1 b_1 - c_0 b_2}{b_0}, \quad (8)$$

where $a_n$ and $b_n$ are defined as

$$a_n = \sum_C O(C) f_n(C), \quad b_n = \sum_C f_n(C). \quad (9)$$

Notice that the expansion in Eq. (7) is slightly different for the current $J$ due to an extra $\alpha$ in Eq. (1) and is given by $J = \sum_{n=0}^{\infty} c_n \alpha^{n+1}$. Using the expressions for $f_n(C)$ and $f_n(C)$ computed earlier yields

$$J = \alpha \left[ 1 - \sum_{i=1}^{\ell} \left( \frac{1}{k_i} + \frac{1}{\gamma_i} \right) \alpha + O(\alpha^2) \right], \quad (10)$$

$$\rho = \frac{1}{L} \sum_{i=1}^{L} \left( \frac{1}{k_i} + \frac{1}{\gamma_i} \right) \alpha + O(\alpha^2), \quad (11)$$

$$\rho_i = \left( \frac{1}{k_i} + \frac{1}{\gamma_i} \right) \alpha + O(\alpha^2), \quad (12)$$

where $\gamma_i = \gamma + (\beta - \gamma)\delta_{i,L}$. These equations constitute our main result. Equation (10) shows that if the initiation rate $\alpha$ is small compared to $k_1, \dots, k_L$ and $\gamma$, than the protein production rate $J$ depends predominantly on the initiation rate, along with the translocation and elongation rates of the first 10 codons, corresponding to the ribosome footprint $\ell$. The importance of the first 10 codons is a direct consequence of the excluded volume interactions: any ribosome already present in that region will prevent a new ribosome from binding the mRNA. It is also important to emphasize that (10)–(12) are exact series expansions around $\alpha = 0$; the approximation is made only when the series is truncated.

**Independent Particle Approximation (IPA).** Interestingly, the excluded volume interactions between particles have no effect on the first two terms in the series expansion of $J$. This motivates us to ask how the expansion in Eq. (10) would look like if we assumed that all particles are independent, i.e. not experiencing any exclusion interaction. In our model, the IPA amounts to replacing $P(C)$ with

$$P^{\text{IPA}}(C) = \frac{1}{Z_L} \prod_{j=1}^{N(C)} w_1(\theta_{X(j)}), \quad (13)$$

where $N(C)$ is the number of particles in a configuration $C$, $\theta$ is one of the two particle states 1 and 2, $X(i)$ is
the position of the i-th particle on the lattice and $Z_L$ is the normalization constant. The weights $w_1(1_i)$ and $w_1(2_i)$ for $i = 1, \dots, L$ are obtained by solving the master equation for a single particle and are given by $\alpha/k_i$ and $\alpha/\gamma_i$, respectively. The corresponding expressions for the current $J$ and local density $\rho_i$ read

$$J^{\text{IPA}} = \frac{\alpha}{\prod_{i=1}^{\ell} (1 + p_i)}, \quad \rho_i^{\text{IPA}} = \frac{p_i}{1 + p_i}, \qquad (14)$$

where $p_i = \alpha(1/k_i + 1/\gamma_i)$. We will use these results later in order to determine the importance of ribosome collisions in real genetic sequences. We also note that the IPA in (13) provides a good approximation to $f_n(C)$ for $n = N(C)$ when the particles in $C$ are far apart from each other.

**Application to mRNA translation in yeast.** We now apply our results to the transcriptome of S. cerevisiae using realistic model parameters. The values of $\alpha$ in the range $0.005 - 4.2 \text{ s}^{-1}$ with the median value of $0.09 \text{ s}^{-1}$ have been previously estimated in Ref. [1] using genome-wide experimental values of the ribosomal density. We assume that the rates $k_i$ are mainly proportional to the gene copy number of tRNAs delivering the corresponding amino acid [28]; the rates are normalized so that the average codon translation rate is equal to the experimental value of 10 codons/s [3]; the estimates of all elongation rates along with the distribution of $\alpha$ is presented in the Supplemental Material. The translocation rate $\gamma$ is fixed to $\gamma = 35 \text{ codons/s}$ [29], and termination is assumed to be fast and comparable to translocation [3], $\beta \approx \gamma$, so that $\gamma_i = \gamma \forall i$ in Eqs. (10)-(12).

In total, we analyzed 5836 gene sequences; for each gene we calculated $J$, $\rho$ and $\rho_i$ for $1 \le i \le L$ up to and including the second order of the perturbative expansion at the corresponding physiological value of $\alpha$. The results were then compared to the exact values obtained numerically with stochastic simulations using the Gillespie algorithm [31] by calculating the percent error $\epsilon$.

For the protein production rate $J$, with the zeroth order of the perturbative expansion (predicting that $J = \alpha$) we obtain an error of $\epsilon < 5\%$ for only 11% of the genes. Remarkably, that percentage jumps to 80.7% when the first-order coefficients in Eq. (6) are taken into account (Fig. 2). Including the second-order coefficients (computed numerically from Eq. (4)) does not significantly improve results, due to a large value of $\alpha \gtrsim 0.15 \text{ s}^{-1}$ in about 20% of the genes. Since the coefficients in Eq. (7) typically alternate in sign, truncating the series will ultimately lead to a wrong result when the value of $\alpha$ is large enough. On the other hand, the IPA does not suffer from this problem and leads to an error $\epsilon < 5\%$ in 94% of genes, whereby only 1% of genes have $\epsilon > 20\%$. The success of the IPA also suggests that ribosome collisions and traffic jams have a minor effect on the rate of translocation, which is in accordance with recent experimental evidence [26, 32]. This is also apparent from the density

<figure><img src="image_3.png" /><figcaption>FIG. 2: (a): Histogram of the percent error $\epsilon$ measuring the discrepancy between the protein production rate obtained by stochastic simulations and the perturbative expansion including the independent-particle approximation (IPA). (b) and (c): Ribosomal density profile obtained by stochastic simulations (black) compared to Eq. (12) (red) for two values of the initiation rate $\alpha$, one close to the median value $\approx 0.09 \text{ s}^{-1}$ (gene YAL045C, b) and the other that is close to the 90th percentile value $\approx 0.23 \text{ s}^{-1}$ (gene YGL034C, c).</figcaption></figure>

profile $\rho_i$, which is very well approximated by the linear approximation in Eq. (12), even for larger values of $\alpha$ (Fig. 2 b and c).

**Identifying determinants of mRNA translation.** Our analytic prediction allows us to decompose the contributions from initiation and elongation to the rate of translocation $J$, thereby addressing a long-standing question about main determinants of protein production rate. Remarkably, the expressions for the current obtained with both the first-order and the independent particle approximation involve only the first $\ell = 10$ codons. This result therefore strongly indicates that, together with the initiation rate, the first 10 codons of the mRNA are the key determinants of the protein production rate by preventing a new ribosome from binding which effectively decreases the initiation rate.

If we assume that the cell maximizes the rate of protein production, given the above result, we would ex-
<figure><img src="image_4.png" /><figcaption>FIG. 3: Histogram of η<sub>J</sub> (yellow) and η<sub>ρ</sub> (blue) for the *S. cerevisiae* genome.</figcaption></figure>

pect to find a signature in the genome for selecting ef-
ficient *fast* codons at the beginning of each gene. We
test this hypothesis by computing the ribosomal current
for the fastest ($J^F$) and the slowest ($J^S$) set of first
$l$ synonymous codons for each gene, using the IPA. A
score $\eta_J = (J - J^S)/(J^F - J^S)$ is then assigned to each
gene, which is 1 (0) when the sequence corresponds to
the fastest (slowest) codon sequence. On the other hand,
one might assume that the cell not only tries to maximize
protein production rates, but at the same time it tries
to minimize the ribosome density $\rho$ on mRNAs. This as-
sumption is motivated by the fact that ribosomes are lim-
iting [22, 33] and highly costly in terms of cellular energy
resources [34], and therefore ribosome queues are to be
avoided. Hence, we also compute $\eta_\rho = (\rho - \rho^S)/(\rho^F - \rho^S)$
for each gene, where $\rho^F$ and $\rho^S$ denote the ribosome den-
sity for the fastest and slowest set of synonymous codons,
respectively.

Figure 3 shows the histogram of η<sub>J</sub> (yellow) and η<sub>ρ</sub>
(blue) computed for 5836 genes of *S.cerevisiae*. Both his-
tograms show an average of 0.7 suggesting the selection
of fast codons near the start codon to maximize *J*, as well
as an overall selection of fast codons along the mRNA to
minimize ρ. However, the width of the distribution of ηρ
is substantially smaller than the one of ηJ. This might
indicate that the optimization of protein production rate
is strongly dependent on the particular gene, since dif-
ferent proteins are needed at different concentrations. In
contrast, the minimization of the number of ribosomes
on mRNAs could be a more general constraint.

Conclusions. We have presented an analytic method that allows us to quantify sequence determinants of protein production rates, using a model for translation that is based on an inhomogeneous exclusion process. Our results demonstrate that the rate of protein production is largely determined by the initiation rate and the elongation rates of the first 10 codons (the ribosome footprint length on the mRNA), which control how fast ribosomes load onto the mRNA, whereby ribosome collisions and queues have a minor effect under physiological conditions.

**Acknowledgments.** MCR and LC contributed equally to this work. JSN was supported by the Leverhulme Trust Early Career Fellowship. MCR was supported by the Biotechnology and Biological Sciences Research Council (BBSRC) BB/N017161/1 and the Scottish Universities Life Sciences Alliance. LC thanks the CNRS for having granted him a demi-délégation (2017–18).

<footnote>* Electronic address: jszavits@staffmail.ed.ac.uk</footnote>

<footnote>[1] H. Gingold, and Y. Pilpel, Mol. Syst. Biol. **7**, 481 (2011).</footnote>

<footnote>[2] A.J. Kemp, R. Betney, L. Ciandrini, A. C. Schwenger, M. C. Romano, and I. Stansfield, Mol. Microbiol. **87**(2), 284 (2013).</footnote>

<footnote>[3] D. Chu, E. Kazana, N. Bellanger, T. Singh, M. F. Tuite, and T. von der Haar, EMBO J. **33**(1), 21 (2014).</footnote>

<footnote>[4] B. Gorgoni, L. Ciandrini, M. R. McFarland, M. C. Romano, and I. Stansfield, Nucleic Acids Res. **44**(19), 9231 (2016).</footnote>

<footnote>[5] S. Varenne, J. Buc, R. Lloubes, and C. Lazdunski, J. Mol Biol. **180**(3), 549 (1984).</footnote>

<footnote>[6] M. A. Sorensen, C. G. Kurland, and S. Pedersen, J. Mol. Biol. **207**, 365–377 (1989).</footnote>

<footnote>[7] T. E. Gorochowski, I. Avcilar-Kucukgoze, R. A. Bovenberg, J.A. Roubos, and Z. Ignatova, ACS Synth. Biol. **5**(7), 710 (2016).</footnote>

<footnote>[8] G. Boël et al, Nature **529**(7586), 358 (2016).</footnote>

<footnote>[9] C. T. MacDonald, J. H. Gibbs, and A. C. Pipkin, Biopolymers **6**, 1–25 (1968); C. T. MacDonald and J. H. Gibbs, Biopolymers **7**, 707–25 (1969).</footnote>

<footnote>[10] H. Zur, and T. Tuller, Nucleic Acids Res. **44**(19), 9031–9049 (2016).</footnote>

<footnote>[11] B. Derrida, M. R. Evans, V. Hakim, and V. Pasquier, J. Phys. A: Math. Gen. **26**, 1493 (1993); G. M. Schütz, and E. Domany, J. Stat. Phys. **72**, 277–96 (1993).</footnote>

<footnote>[12] R. A. Blythe, and M. R. Evans, J. Phys. A: Math. Theor. **40**, R333–R441 (2007).</footnote>

<footnote>[13] T. Chou, and G. Lakatos, Phys. Rev. Lett. **93**, 198101 (2004).</footnote>

<footnote>[14] L. B. Shaw, J. P. Sethna, and K. H. Lee, Phys. Rev. E **70**, 021901 (2004).</footnote>

<footnote>[15] R. K. P. Zia, J. J. Dong, and B. Schmittmann, J. Stat. Phys. **144**, 405 (2011).</footnote>

<footnote>[16] J. Schmidt, V. Popkov and A. Schadschneider, EPL **110**, 20008 (2015).</footnote>

<footnote>[17] N. T. Ingolia, S. Ghaemmaghami, J. R. S. Newman, and J. S. Weissman, Science **324**, 218–223 (2009).</footnote>

<footnote>[18] M. A. Gilchrist, A. Wagner, J. Theor. Biol., **239**(4), 417-434 (2006).</footnote>

<footnote>[19] N. Mitarai, K. Sneppen, and S. Pedersen, J. Mol. Biol., **382**(1), 236-245 (2008).</footnote>

<footnote>[20] S. Reuveni, I. Meilijson, M. Kupiec, E. Ruppin, and T. Tuller, PLoS Comput. Biol. **7**, e1002127 (2011).</footnote>

<footnote>[21] L. Ciandrini, I. Stansfield, and M. C. Romano, PLoS Comput. Biol. **9**(1), e1002866 (2013).</footnote>

<footnote>[22] P. Shah, Y. Ding, M. Niemczyk, G. Kudla, and J.B. Plotkin, Cell **153**(7), 1589 (2013).</footnote>

<footnote>[23] L. Ciandrini, I. Stansfield, and M. C. Romano. Phys. Rev. E **81**, 051904 (2010).</footnote>

<footnote>[24] S. Klumpp, Y. Chai, and R. Lipowsky, Phys. Rev. E **78**,</footnote>
041909 (2008).

[25] J. Szavits-Nossan, J. Phys. A: Math. Theor. **46**, 315001 (2013).

[26] G. Kudla, A. W. Murray, D. Tollervey, and J. B. Plotkin, Science **324**, 255–258 (2009).

[27] H. M. Salis, E. A. Mirsky, and C. A. Voigt, Nature Biotechnol. **27**, 946–950 (2009).

[28] R. Percudani, A. Pavesi, and S. Ottonello, J. Mol. Biol. **268**(2), 322–330 (1997).

[29] A. Savelsbergh, V. I. Katunin, D. Mohr, F. Peske, M. V. Rodnina, and W. Wintermeyer, Mol. Cell **11**, 1517–1523 (2003).

[30] Y. Arava, Y. Wang, J. D. Storey, C. L. Liu, P. O. Brown, and D. Herschlag, PNAS **100**(7), 3889–3894 (2003).

[31] D. T. Gillespie, J. Phys. Chem. **81**(25), 2340 (1977).

[32] B. S. Laursen, H. P. Sørensen, K. K. Mortensen, and H. U. Sperling-Petersen, Microbiol. Mol. Biol. Rev. **69**, 101 (2005).

[33] P. Greulich, L. Ciandrini, R. J. Allen, and M. C. Romano. Phys. Rev. E **85**, 011142 (2012).

[34] K. C. Keiler, Nature Rev. Microbiol. **13**, 285 (2015).
Supplemental Material to:
Deciphering mRNA Sequence Determinants of
Protein Production Rate

EXACT SOLUTION OF THE STEADY-STATE MASTER EQUATION

The exact solution, Eq. (3) in the main text, pertains to any ergodic Markov jump process with a finite number of states. In order to derive the exact solution, let us write again the steady-state master equation in a matrix form

$$MP = 0, \qquad (S1)$$

where **P** is a column vector whose elements are steady-state probabilities $P_i$, the $N$ states are indexed by $i = 1...N$ and the $N \times N$ matrix $M$ is given by

$$M_{ij} = \begin{cases} W_{jk} & i \neq j \\ -\sum_{k \neq i} W_{ik} & i = j, \end{cases} \qquad (S2)$$

where $W_{ij}$ is a transition rate from a state $i$ to $j$. We note from Eq. (S2) that the sum of all elements in each column of $M$ is zero

$$\sum_{i=1}^{N} M_{ij} = \sum_{\substack{i=1 \\ i \neq j}}^{N} W_{ji} - \sum_{\substack{k=1 \\ k \neq j}}^{N} W_{jk} = 0. \qquad (S3)$$

Consequently, the sum all row vectors of $M$ is zero, which means that the vectors are linearly dependent and thus the determinant of $M$ is equal to zero,

$$\det M = 0. \qquad (S4)$$

Combining Eq. (S4) with the Laplace expansion of a determinant yields

$$0 = \det M = \sum_{i=1}^{N} M_{ij} C_{ij} = \sum_{j=1}^{N} M_{ij} C_{ij}. \qquad (S5)$$

Here $C_{ij}$ is a cofactor of $M$ defined as $C_{ij} = (-1)^{i+j}\det M^{(i,j)}$, where the matrix $M^{(i,j)}$ is obtained from $M$ by removing $i$-th row and $j$-th column. Inserting Eq. (S2) into (S5) gives

$$\begin{aligned} \det M &= \sum_{\substack{i=1 \\ i \neq j}}^{N} M_{ij} C_{ij} + M_{jj} C_{jj} = \sum_{\substack{i=1 \\ i \neq j}}^{N} M_{ij} C_{ij} - \sum_{\substack{k=1 \\ k \neq j}}^{N} M_{kj} C_{jj} \\ &= \sum_{\substack{i=1 \\ i \neq j}}^{N} M_{ij} (C_{ij} - C_{jj}) = 0, \end{aligned} \qquad (S6)$$

from which we conclude that

$$C_{ij} = C_{jj}, \qquad (S7)$$

for any $i$ and $j$, i.e. the cofactor $C_{ij}$ does not depend on the state $i$. Inserting Eq. (S7) back into Eq. (S4) gives

$$\sum_{j=1}^{N} M_{ij} C_{jj} = 0, \qquad (S8)$$

which is precisely the starting steady-state master equation, Eq. (S1). We thus conclude that

$$P_i = \frac{C_{ii}}{\sum_j C_{jj}} = \frac{\det M^{(i,i)}}{\sum_{j=1}^{N} \det M^{(j,j)}}. \qquad (S9)$$
ARGUMENT FOR $f^{(n)}(C) = 0$ WHENEVER THE NUMBER OF PARTICLES IN A CONFIGURATION C IS LARGER THAN $n$

As noted in the main text, a key observation in our analysis is that $f^{(n)}(C) = 0$ whenever the number of particles in $C$ is larger than $n$. We present here the argument for the case of $n = 1$; a similar argument applies to higher orders. For $n = 1$, using Eq. (5) in the main text we obtain

$$f^{(1)}(C) = \frac{1}{e(C)} \left[ \sum_{C'} I_{C' \to C} \delta_{C', \emptyset} + \sum_{C'} (1 - I_{C' \to C}) W_{C' \to C} f^{(1)}(C') - \delta_{C, \emptyset} \right]. \quad (S10)$$

Consider $C$ to be one full lattice configuration with all $N = \lfloor(L - 1 + \ell)/\ell\rfloor$ particles in state 1. The first term on the right hand side of Eq. (S10) is zero because there is no direct transition from the empty configuration $\emptyset$ to $C$. The second term on the right hand side is also zero, because $C$ can only be accessed through the transition rate $\alpha$. Finally, the third term on the right hand side is clearly zero. Hence, $f^{(1)}(C) = 0$. Now, it is easily seen that all full lattice configurations $C'$ that can be accessed only from $C$ will fulfil $f^{(1)}(C') = 0$. Iterating further this procedure to the configurations that can be accessed from those, it becomes clear that all full lattice configurations $C$ have $f^{(1)}(C) = 0$. This argument can be further extended to all configurations with $N - 1, N - 2, N - 3, ..., 2$ particles, resulting in $f^{(1)}(C) = 0$ for all of them.

The situation changes when we consider configurations with only 1 particle on the lattice, since then there is a direct transition from the empty configuration $\emptyset$ to one of those configurations, leading to $f^{(1)}(C) \neq 0$. The rest of the coefficients $f^{(1)}(C)$ for 1-particle configurations will depend on each other, and hence, they do not vanish. Therefore, we can conclude that only 1-particle configurations have non-vanishing $f^{(1)}(C)$. For $n = 1$, this key observation allows us to write Eq. (5) taking into account configurations with only one particle and discarding all the others, which yields Eqs. (6a)-(6c) in the main text.

EQUATIONS FOR THE SECOND-ORDER COEFFICIENTS $f^{(2)}(C)$

As stated in the main text, all second-order coefficients $f^{(2)}(C)$ whereby a configuration $C$ has more than two particles are equal to zero. The equations for the non-zero coefficients $f^{(2)}(C)$ are presented below, where we use the notation $\gamma_j = \gamma$ for $1 \le j \le L - 1$ and $\gamma_L = \beta$.

First, we look at configurations with one particle at site $i = 1$ and the other at site $j = i + l, ..., L$:

$$f^{(2)}(1_1 1_j) = \frac{1}{k_1 + k_j} \left[ f^{(1)}(1_j) + (1 - \delta_{1+l,j}) \gamma_{j-1} f^{(2)}(1_1 2_{j-1}) \right] \quad (S11a)$$

$$f^{(2)}(1_1 2_j) = \frac{1}{k_1 + \gamma_j} \left[ f^{(1)}(2_j) + k_j f^{(2)}(1_1 1_j) \right] \quad (S11b)$$

$$f^{(2)}(2_1 1_j) = \frac{1}{(1 - \delta_{1+l,j})\gamma_1 + k_j} \left[ k_1 f^{(2)}(1_1 1_j) + (1 - \delta_{1+l,j})\gamma_{j-1} f^{(2)}(2_1 2_{j-1}) \right] \quad (S11c)$$

$$f^{(2)}(2_1 2_j) = \frac{1}{(1 - \delta_{1+l,j})\gamma_1 + \gamma_j} \left[ k_1 f^{(2)}(1_1 2_j) + k_j f^{(2)}(2_1 1_j) \right]. \quad (S11d)$$

where we used the Kronecker delta to account for the excluded volume interaction. Next, we look at configurations with particles at sites $i = 2, ..., L - l$ and $j = i + l, ..., L$.

$$f^{(2)}(1_i 1_j) = \frac{1}{k_i + k_j} \left[ \gamma_{i-1} f^{(2)}(2_{i-1} 1_j) + (1 - \delta_{i+l,j}) \gamma_{j-1} f^{(2)}(1_i 2_{j-1}) \right] \quad (S12a)$$

$$f^{(2)}(1_i 2_j) = \frac{1}{k_i + \gamma_j} \left[ \gamma_{i-1} f^{(2)}(2_{i-1} 2_j) + k_j f^{(2)}(1_i 1_j) \right] \quad (S12b)$$

$$f^{(2)}(2_i 1_j) = \frac{1}{(1 - \delta_{i+l,j})\gamma_i + k_j} \left[ k_i f^{(2)}(1_i 1_j) + (1 - \delta_{i+l,j})\gamma_{j-1} f^{(2)}(2_i 2_{j-1}) \right] \quad (S12c)$$

$$f^{(2)}(2_i 2_j) = \frac{1}{(1 - \delta_{i+l,j})\gamma_i + \gamma_j} \left[ k_i f^{(2)}(1_i 2_j) + k_j f^{(2)}(2_i 1_j) \right]. \quad (S12d)$$
Finally, the equations for configurations with only one particle are given by

$$f^{(2)}(1_1) = \frac{1}{k_1} \left[ \gamma_L f^{(2)}(1_1 2_L) + f^{(1)}(\emptyset) \right], \qquad (S13a)$$

$$f^{(2)}(2_1) = \frac{k_1}{\gamma_1} f^{(2)}(1_1) + \frac{\gamma_L}{\gamma_1} f^{(2)}(2_1 2_L), \qquad (S13b)$$

$$f^{(2)}(1_i) = \frac{1}{k_i} \left[ \gamma_{i-1} f^{(2)}(2_{i-1}) + \theta(L-l+1-i) \gamma_L f^{(2)}(1_i 2_L) - f^{(1)}(1_i) \right], \quad i = 2, \dots, L \qquad (S13c)$$

$$f^{(2)}(2_i) = \frac{1}{\gamma_i} \left[ k_i f^{(2)}(1_i) + \theta(L-l+1-i) \gamma_L f^{(2)}(2_i 2_L) - f^{(1)}(2_i) \right], \quad i = 2, \dots, L \qquad (S13d)$$

where $\theta(n) = 0$ for $n < 0$ and $1$ for $n \ge 0$.

The equations (S11)-(S13) can be easily solved numerically by iteration. We first solve Eqs. (S11) for $i = 1, j = l+1$. We then iterate the recursion relation in Eqs. (S11) for fixed $i = 1$ and $j = 2+l, \dots, L$. We then solve equations (S12) for fixed $i = 2$ and $j = l+3, \dots, L$ and repeat this procedure until $i = L-l$ and $j = L$. Finally, we calculate the coefficients $f^{(2)}(1_i)$ and $f^{(2)}(2_i)$ for the single-particle configurations using $f^{(2)}(1_i 2_L)$ and $f^{(2)}(2_i 2_L)$ that we solved in the previous steps. We note that the only unknown that we cannot determine is $f^{(1)}(\emptyset)$, but it turns out that all terms containing $f^{(1)}(\emptyset)$ will cancel out later when we calculate the coefficients in the series expansions (11)-(13) in the main text.

## ESTIMATES OF THE CODON ELONGATION RATES

The total translation rate $\omega_i$ of the codon $i$ is:

$$\frac{1}{k_i} + \frac{1}{\gamma} = \frac{1}{\omega_i}, \qquad (S14)$$

where the translocation rate $\gamma$ is set to $35s^{-1}$ as mentioned in the main text. The codon elongation rate $k_i$, which represents the average arrival and recognition time of the cognate tRNA, is determined by following the procedure introduced in [1], which we report here.

For each of the 41 tRNAs types $j$ we consider their gene copy number (GCN), which allows us to provide a first estimate of the rate $k_j$:

$$k_j = r \frac{GCN_j}{\sum_{j=1}^{41} GCN_j}, \qquad (S15)$$

where $GCN_j$ is the gene copy number of the tRNA of type $j$ with $j = 1, \dots, 41$, and $r$ is a proportionality constant. These rates were then adjusted to take into account experimental evidence suggesting that the translation rates of codons using the G-U wobble are reduced by 39% compared to their G-C counterparts; analogously, codons using the wobble I-C and codons using the wobble I-A are reduced by 36% relative to their I-U counterparts [2]. To calculate the proportionality constant $r$, we used the experimental value of 10 codons/s for the average codon translation rate $\langle \omega_i \rangle$ defined as

$$\langle \omega_i \rangle = \sum_{i=1}^{61} \left( \frac{k_i \gamma}{k_i + \gamma} \right) \frac{n_i}{n} \qquad (S16)$$

where $n_i/n$ is the relative abundance of all different codon types in the cell, $n_i$ is total number of codons in the cell of exactly type $i$, and $n = \sum_{i=1}^{61} n_i$ (the three STOP codons are excluded) [1]. Table I summarises the resulting elongation rates $k_i$ for 61 codons.

The estimates provided are valid in physiological conditions, when amino acids are not limiting. In case of amino acid starvation the rates of the codons affected should be modified (the elongation rates depend on the abundance of tRNAs charged with the correct amino acid).

## DISTRIBUTION OF THE TRANSLATION INITIATION RATES

The translation initiation rates $\alpha$ were estimated for each gene in Ref. [1] by matching the average ribosome density predicted by the model with the experimental value of the ribosome density obtained in Ref. [3]. The resulting values
<table><thead><tr><th>tRNA</th><th>anti-codon</th><th>codon</th><th>GCN</th><th>k<sub>i</sub> [1/s]</th></tr></thead><tbody><tr><td>Ala1</td><td>IGC</td><td>GCU</td><td>11.00</td><td>18.34</td></tr><tr><td>Ala1</td><td>IGC</td><td>GCC</td><td>11.00</td><td>11.74</td></tr><tr><td>Ala2</td><td>UGC</td><td>GCA</td><td>5.00</td><td>8.34</td></tr><tr><td>Ala2</td><td>UGC</td><td>GCG</td><td>5.00</td><td>8.34</td></tr><tr><td>Arg1</td><td>CCU</td><td>AGG</td><td>1.00</td><td>1.67</td></tr><tr><td>Arg2</td><td>ICG</td><td>CGU</td><td>6.00</td><td>10.01</td></tr><tr><td>Arg2</td><td>ICG</td><td>CGC</td><td>6.00</td><td>6.40</td></tr><tr><td>Arg2</td><td>ICG</td><td>CGA</td><td>6.00</td><td>6.40</td></tr><tr><td>Arg3</td><td>UCU</td><td>AGA</td><td>11.00</td><td>18.34</td></tr><tr><td>Arg4</td><td>CCG</td><td>CGG</td><td>1.00</td><td>1.67</td></tr><tr><td>Asn</td><td>GUU</td><td>AAU</td><td>10.00</td><td>10.17</td></tr><tr><td>Asn</td><td>GUU</td><td>AAC</td><td>10.00</td><td>16.68</td></tr><tr><td>Asp</td><td>GUC</td><td>GAU</td><td>15.00</td><td>15.26</td></tr><tr><td>Asp</td><td>GUC</td><td>GAC</td><td>15.00</td><td>25.01</td></tr><tr><td>Cys</td><td>GCA</td><td>UGU</td><td>4.00</td><td>4.07</td></tr><tr><td>Cys</td><td>GCA</td><td>UGC</td><td>4.00</td><td>6.67</td></tr><tr><td>Gln1</td><td>UUG</td><td>CAA</td><td>9.00</td><td>15.01</td></tr><tr><td>Gln2</td><td>CUG</td><td>CAG</td><td>1.00</td><td>1.67</td></tr><tr><td>Glu3</td><td>UUC</td><td>GAA</td><td>14.00</td><td>23.35</td></tr><tr><td>Glu4</td><td>CUC</td><td>GAG</td><td>2.00</td><td>3.34</td></tr><tr><td>Gly1</td><td>GCC</td><td>GGU</td><td>16.00</td><td>16.28</td></tr><tr><td>Gly1</td><td>GCC</td><td>GGC</td><td>16.00</td><td>26.68</td></tr><tr><td>Gly2</td><td>UCC</td><td>GGA</td><td>3.00</td><td>5.00</td></tr><tr><td>Gly3</td><td>CCC</td><td>GGG</td><td>2.00</td><td>3.34</td></tr><tr><td>His</td><td>GUG</td><td>CAU</td><td>7.00</td><td>7.12</td></tr><tr><td>His</td><td>GUG</td><td>CAC</td><td>7.00</td><td>11.67</td></tr><tr><td>Ile1</td><td>UAU</td><td>AUA</td><td>2.00</td><td>3.34</td></tr><tr><td>Ile2</td><td>IAU</td><td>AUU</td><td>13.00</td><td>21.68</td></tr><tr><td>Ile2</td><td>IAU</td><td>AUC</td><td>13.00</td><td>13.87</td></tr><tr><td>Leu1</td><td>UAG</td><td>CUA</td><td>3.00</td><td>5.00</td></tr><tr><td>Leu1</td><td>UAG</td><td>CUG</td><td>3.00</td><td>5.00</td></tr></tbody></table>

<table><caption>TABLE I: Elongation rates $k_i$ considering supply (gene copy number of tRNAs) and wobble base-pairing.</caption><thead><tr><th>tRNA</th><th>anti-codon</th><th>codon</th><th>GCN</th><th>k<sub>i</sub> [1/s]</th></tr></thead><tbody><tr><td>Leu3</td><td>CAA</td><td>UUG</td><td>10.00</td><td>16.68</td></tr><tr><td>Leu4</td><td>UAA</td><td>UUA</td><td>7.00</td><td>11.67</td></tr><tr><td>Leu5</td><td>GAG</td><td>CUU</td><td>1.00</td><td>1.02</td></tr><tr><td>Leu5</td><td>GAG</td><td>CUC</td><td>1.00</td><td>1.67</td></tr><tr><td>Lys1</td><td>CUU</td><td>AAG</td><td>14.00</td><td>23.35</td></tr><tr><td>Lys2</td><td>UUU</td><td>AAA</td><td>7.00</td><td>11.67</td></tr><tr><td>Met</td><td>CAU</td><td>AUG</td><td>5.00</td><td>8.34</td></tr><tr><td>Phe</td><td>GAA</td><td>UUU</td><td>10.00</td><td>10.17</td></tr><tr><td>Phe</td><td>GAA</td><td>UUC</td><td>10.00</td><td>16.68</td></tr><tr><td>Pro1</td><td>UGG</td><td>CCA</td><td>10.00</td><td>16.68</td></tr><tr><td>Pro1</td><td>UGG</td><td>CCG</td><td>10.00</td><td>16.68</td></tr><tr><td>Pro2</td><td>IGG</td><td>CCU</td><td>2.00</td><td>3.34</td></tr><tr><td>Pro2</td><td>IGG</td><td>CCC</td><td>2.00</td><td>2.13</td></tr><tr><td>Ser2</td><td>IGA</td><td>UCU</td><td>11.00</td><td>18.34</td></tr><tr><td>Ser2</td><td>IGA</td><td>UCC</td><td>11.00</td><td>11.74</td></tr><tr><td>Ser3</td><td>GCU</td><td>AGU</td><td>4.00</td><td>4.07</td></tr><tr><td>Ser3</td><td>GCU</td><td>AGC</td><td>4.00</td><td>6.67</td></tr><tr><td>Ser4</td><td>UGA</td><td>UCA</td><td>3.00</td><td>5.00</td></tr><tr><td>Ser5</td><td>CGA</td><td>UCG</td><td>1.00</td><td>1.67</td></tr><tr><td>Thr1</td><td>IGU</td><td>ACU</td><td>11.00</td><td>18.34</td></tr><tr><td>Thr1</td><td>IGU</td><td>ACC</td><td>11.00</td><td>11.74</td></tr><tr><td>Thr2</td><td>CGU</td><td>ACG</td><td>1.00</td><td>1.67</td></tr><tr><td>Thr3</td><td>UGU</td><td>ACA</td><td>4.00</td><td>6.67</td></tr><tr><td>Trp</td><td>CCA</td><td>UGG</td><td>6.00</td><td>10.01</td></tr><tr><td>Tyr</td><td>GIA</td><td>UAU</td><td>8.00</td><td>8.14</td></tr><tr><td>Tyr</td><td>GIA</td><td>UAC</td><td>8.00</td><td>13.34</td></tr><tr><td>Val1</td><td>IAC</td><td>GUU</td><td>14.00</td><td>23.35</td></tr><tr><td>Val1</td><td>IAC</td><td>GUC</td><td>14.00</td><td>14.94</td></tr><tr><td>Val2</td><td>UAC</td><td>GUA</td><td>2.00</td><td>3.34</td></tr><tr><td>Val2b</td><td>CAC</td><td>GUG</td><td>2.00</td><td>3.34</td></tr></tbody></table>

for 5836 genes that we analyzed in the main text are in the range $0.005 - 4.2 \text{ s}^{-1}$ with the median value of $0.09 \text{ s}^{-1}$. The histogram of the rates is presented in Figure S1, showing that the initiation rate $\alpha$ for the majority of genes is indeed much smaller than the smallest elongation rate, which has the value of 1.02 codons/s for the CUU codon (see Table I). The list of all translation initiation rates can be found in the Supplemental Material of Ref. [1].

* Electronic address: jszavits@staffmail.ed.ac.uk

[1] L. Ciandrini, I. Stansfield, and M. C. Romano, PLoS Comput. Biol. 9(1), e1002866 (2013).

[2] M. A. Gilchrist and A. Wagner, J. Theor. Biol. 239 417-434 (2006).

[3] Y. Arava, Y. Wang, J. D. Storey, C. L. Liu, P. O. Brown, and D. Herschlag, PNAS 100(7), 3889-3894 (2003).
<figure><img src="image_5.png" /><figcaption>FIG. S1: Histogram of the estimated initiation rates for 5836 genes from the *S. cerevisiae* genome. The 95th percentile value of α is ≈ 0.3, which is much smaller than the smallest elongation rate, which is 1.02 codons/s for the CUU codon (vertical dashed line). Only 15 genes have α larger than min;{kᵢ} = 1.02 codons/s.</figcaption></figure>