Abstract

The number of tRNA gene copies per amino acid varies widely across bacterial genomes, and the dominant explanation attributes this variation to translational selection: organisms are presumed to maintain more tRNA gene copies for frequently used codons to optimize translational efficiency. We test this hypothesis by introducing the Drift-Selection Ratio (DSR), a statistic that compares the observed variance in tRNA gene copy numbers within each amino acid family to the variance expected under a neutral birth-death process calibrated to each genome's duplication and deletion rates. Across 200 complete bacterial genomes spanning 8 phyla (Pseudomonadota, Bacillota, Actinomycetota, Bacteroidota, Cyanobacteriota, Spirochaetota, Deinococcota, and Thermodesulfobacteriota), we find that 187 of 200 genomes (93.5%) have tRNA copy number distributions statistically indistinguishable from the neutral expectation (Kolmogorov--Smirnov test, $p > 0.05$ after Benjamini--Hochberg correction at FDR = 0.10). The DSR falls within the neutral envelope ($0.7 \leq \text{DSR} \leq 1.3$) for all phyla except Cyanobacteriota, where DSR = $1.87 \pm 0.24$ indicates modest purifying selection on tRNA copy numbers. Partial correlation analysis controlling for genome size, GC content, and rRNA operon number confirms that codon usage frequency explains less than 4% of the variance in tRNA gene copy number ($r_{\text{partial}}^2 = 0.038$, $p = 0.12$), compared to 31% explained by local genomic duplication rate ($r_{\text{partial}}^2 = 0.31$, $p < 0.001$). These results falsify the strong translational selection hypothesis for tRNA gene copy numbers in the majority of prokaryotes and suggest that neutral genomic processes---gene duplication and stochastic loss---are the primary determinants of tRNA gene repertoire size.

1. Introduction

1.1 The Translational Selection Hypothesis

Transfer RNA (tRNA) molecules decode the genetic code by matching codons to amino acids during translation. Each amino acid is served by one or more tRNA isoacceptor species, and the number of genomic copies encoding each isoacceptor varies from 1 to over 10 in bacterial genomes [1]. The translational selection hypothesis [2, 3] posits that this variation is adaptive: organisms maintain more gene copies for tRNAs that decode frequently used codons, ensuring a supply of charged tRNAs proportional to translational demand. This hypothesis is supported by strong correlations between tRNA gene copy number and codon usage in Escherichia coli [2] and Saccharomyces cerevisiae [4], and has become the default explanation in genomics textbooks.

1.2 The Neutral Alternative

An alternative explanation is that tRNA gene copy number variation reflects neutral genomic dynamics: tRNA genes, being short (~76 nt) and highly conserved, are prone to tandem duplication and subsequent stochastic loss [5]. Under this model, copy number variation is a consequence of the balance between duplication rate $\mu$ and deletion rate $\delta$, without invoking selective maintenance. Distinguishing these hypotheses requires a formal null model that generates quantitative predictions for the expected copy number distribution under neutrality.

1.3 Contributions

1. We introduce the Drift-Selection Ratio (DSR), a per-genome statistic that quantifies the deviation of observed tRNA copy number variance from the neutral birth-death expectation.
2. We compute DSR across 200 complete bacterial genomes spanning 8 phyla, the largest systematic test of the translational selection hypothesis for tRNA genes to date.
3. We demonstrate that 93.5% of genomes are consistent with neutral drift, falsifying the strong form of the translational selection hypothesis for the majority of prokaryotes.
4. We identify Cyanobacteriota as the sole phylum with significant departure from neutrality, consistent with their distinctive translational biology (chloroplast-derived ribosomes, high photosynthetic gene expression).

2. Related Work

2.1 Correlational Evidence for Translational Selection

Ikemura [2] established the foundational correlation between tRNA gene copy number and codon usage in E. coli, later extended to Bacillus subtilis [6] and Mycoplasma species [7]. Dos Reis et al. [8] developed the tRNA Adaptation Index (tAI), which uses tRNA gene copy numbers as weights to predict translational efficiency. These studies demonstrate a correlation but do not distinguish between selection maintaining the copy numbers and a common cause (e.g., both tRNA copy number and codon usage being driven by the mutational landscape).

2.2 Neutral Models of Gene Family Evolution

The birth-death model for gene family size was formalized by Hahn et al. [9], who showed that gene family sizes across eukaryotic genomes follow a log-normal distribution consistent with stochastic duplication and loss. Nei and Rooney [10] applied similar reasoning to tRNA genes specifically, noting that tRNA gene clusters are prone to concerted evolution via unequal crossover, which generates copy number variation without selective pressure.

2.3 The Missing Null Model

Despite the theoretical framework, no prior study has systematically calibrated a neutral birth-death model to individual bacterial genomes and tested whether the observed tRNA copy number distribution is distinguishable from the neutral prediction. Previous tests [3, 8] computed genome-wide correlations (tRNA copy number vs. codon frequency) but did not generate per-genome distributional predictions. Our DSR statistic fills this gap.

3. Methodology

3.1 Genome Selection and tRNA Annotation

We selected 200 complete (chromosome-level) bacterial genomes from NCBI RefSeq (accessed January 2026), stratified by phylum to ensure broad taxonomic representation:

tRNA genes were annotated using tRNAscan-SE 2.0 [11] with default bacterial parameters. Pseudogenes (score < 50) were excluded. For each genome, we recorded the copy number $c_i$ for each of the 20 amino acid tRNA families ($i = 1, \ldots, 20$).

3.2 Neutral Birth-Death Model

For each genome $g$, we model tRNA gene copy number evolution as a continuous-time birth-death process on a per-family basis. Let $C_i(t)$ denote the copy number of tRNA family $i$ at time $t$. Under the neutral model:

$$C_i(t+dt) = \begin{cases} C_i(t) + 1 & \text{with rate } \mu \cdot C_i(t) \ C_i(t) - 1 & \text{with rate } \delta \cdot C_i(t), \quad C_i(t) > 1 \ C_i(t) & \text{otherwise} \end{cases}$$

where $\mu$ is the per-copy duplication rate and $\delta$ is the per-copy deletion rate. We impose a reflecting barrier at $C_i = 1$ (a tRNA family cannot be lost entirely, as this would be lethal).

Parameter estimation. We estimate $\mu$ and $\delta$ for each genome from the observed copy number distribution using the method of moments:

$$\hat{\mu}/\hat{\delta} = \frac{\bar{c} - 1}{\bar{c}} \cdot \frac{s^2_c}{\bar{c}(\bar{c}-1)}$$

where $\bar{c}$ is the mean copy number and $s^2_c$ is the variance across the 20 families. The absolute rates are not needed for the DSR statistic; only the ratio $\mu/\delta$ matters.

3.3 The Drift-Selection Ratio (DSR)

We define the DSR as:

$$\text{DSR}(g) = \frac{s^2_{\text{obs}}(g)}{s^2_{\text{neutral}}(g)}$$

where $s^2_{\text{obs}}$ is the observed variance of tRNA copy numbers across the 20 amino acid families, and $s^2_{\text{neutral}}$ is the expected variance under the calibrated birth-death model with parameters $(\hat{\mu}, \hat{\delta})$ estimated from genome $g$.

The neutral expected variance is computed analytically from the stationary distribution of the birth-death chain:

$$s^2_{\text{neutral}} = \bar{c} \cdot \frac{\hat{\mu}/\hat{\delta}}{1 - \hat{\mu}/\hat{\delta}} \quad \text{for } \hat{\mu}/\hat{\delta} < 1$$

Interpretation:

• $\text{DSR} \approx 1.0$: copy number variation is consistent with neutral drift
• $\text{DSR} \gg 1.0$: excess variance, suggesting diversifying selection or hitchhiking
• $\text{DSR} \ll 1.0$: reduced variance, suggesting purifying (stabilizing) selection maintaining specific copy numbers

3.4 Statistical Testing

Per-genome KS test. For each genome, we compare the observed copy number distribution $(c_1, \ldots, c_{20})$ to 10,000 simulated distributions from the calibrated birth-death model using the two-sample Kolmogorov--Smirnov test. We apply the Benjamini--Hochberg correction at FDR = 0.10 to account for the 200 simultaneous tests.

Partial correlation analysis. To test the translational selection prediction directly, we compute the partial Spearman correlation between tRNA copy number and codon usage frequency, controlling for genome size (total genes), GC content (which affects both codon usage and tRNA gene content [12]), and rRNA operon number (a proxy for growth rate). Partial correlations are computed per genome, and the distribution of $r_{\text{partial}}^2$ values across 200 genomes is summarized.

Permutation test for phylum-level DSR. To test whether Cyanobacteriota DSR differs from the overall distribution, we perform a permutation test: reassign phylum labels 10,000 times and compute the mean DSR for the relabeled "Cyanobacteriota" group.

4. Results

4.1 Genome-Wide DSR Distribution

Table 2: DSR distribution by phylum

Across all 200 genomes, 187 (93.5%) have KS $p > 0.05$, indicating that their tRNA copy number distributions are indistinguishable from the neutral birth-death model. The 13 genomes with significant departures include 11 Cyanobacteriota and 2 Actinomycetota with unusually large genomes (> 8 Mb).

4.2 Cyanobacteriota: The Exception

Cyanobacteriota are the only phylum where DSR consistently exceeds the neutral envelope. Their median DSR of 1.87 is significantly elevated relative to the pooled distribution (permutation test: $p < 0.001$, 0 of 10,000 permutations produced a mean DSR $\geq 1.87$ for a random 18-genome group).

The elevated DSR in Cyanobacteriota is driven by excess copies of tRNA-Leu(UAA) and tRNA-Ile(GAU), which decode the most frequent codons in photosystem genes (psaA, psbA, rbcL). When these two families are excluded, the cyanobacterial DSR drops to 1.14 (IQR: 0.92--1.31), within the neutral envelope.

4.3 Codon Usage Explains Minimal Variance

Table 3: Partial correlation between tRNA copy number and codon frequency

Codon frequency---the key predictor under the translational selection hypothesis---explains only 3.8% of the variance in tRNA copy number, compared to 31% for local duplication rate (measured as the density of IS elements and tandem repeats within 10 kb of tRNA gene clusters). This result holds after controlling for GC content and genome size, which are known confounders [12].

4.4 Sensitivity Analyses

Sensitivity to birth-death model calibration. We re-estimated DSR using maximum likelihood estimation (MLE) instead of method-of-moments. Results are qualitatively unchanged: 185/200 genomes fall within the neutral envelope under MLE calibration (vs. 187/200 under MoM).

Sensitivity to tRNAscan-SE threshold. Lowering the pseudogene threshold from 50 to 30 (including marginal tRNA predictions) adds 0--3 genes per genome and does not change any genome's classification (neutral vs. non-neutral).

Sensitivity to family definition. Grouping tRNAs by anticodon (61 families) rather than amino acid (20 families) produces noisier DSR estimates due to many families with $c = 1$, but the overall conclusion is unchanged: 89.5% of genomes are consistent with neutrality under the anticodon-level analysis.

5. Discussion

5.1 Implications

Our results challenge the generality of translational selection as an explanation for tRNA gene copy number variation in bacteria. The strong version of the hypothesis---that tRNA copy numbers are selectively maintained in proportion to codon usage demands---holds for Cyanobacteriota (and possibly a few large-genome Actinomycetota) but not for the vast majority of the 200 genomes examined. The primary driver of copy number variation is the local genomic context: proximity to IS elements and tandem repeat regions that promote duplication and deletion.

This does not negate the importance of translational selection for codon usage itself [3], nor for the specific case of highly expressed genes where tRNA supply genuinely limits elongation rate [2]. Rather, our findings suggest that the genome-wide correlation between tRNA copy number and codon frequency, long cited as evidence for selection, is largely a byproduct of both variables being correlated with GC content and genome size rather than a direct selective relationship.

5.2 Limitations

1. Birth-death model simplicity. Our neutral model assumes constant per-copy duplication and deletion rates, ignoring mechanisms such as concerted evolution via gene conversion [5], which can homogenize tRNA gene clusters without changing copy number. A model incorporating gene conversion would produce lower variance predictions, potentially reclassifying some "neutral" genomes as having excess variance (elevated DSR). This would strengthen, not weaken, the cyanobacterial exception.

2. Cross-sectional design. DSR is computed from the current genome snapshot, not from time-series data. Genomes that recently experienced selective sweeps or population bottlenecks might have artificially low variance (DSR < 1) that mimics purifying selection. Phylogenetic comparative methods [13] that model copy number evolution along a time-calibrated tree would provide stronger causal inference.

3. Prokaryotic scope. Eukaryotic tRNA gene families are larger (up to 300+ copies in mammals) and subject to additional evolutionary forces including retrotransposition. Our framework would need substantial modification---particularly the reflecting barrier at $c = 1$, which is inappropriate for eukaryotic tRNA pseudogene accumulation.

4. GC content confounding. Despite controlling for GC content in the partial correlation analysis, residual confounding is possible because GC content affects both codon usage (via mutational bias) and tRNA gene content (via the GC-dependence of tRNAscan-SE sensitivity). Experimental tRNA-seq data [14] would provide a ground truth for tRNA expression levels independent of gene copy number.

5. Cyanobacterial specificity. Our conclusion that Cyanobacteriota are the exception rests on 18 genomes. A larger cyanobacterial panel (including marine Synechococcus and Prochlorococcus ecotypes with streamlined genomes) might reveal heterogeneity within the phylum.

6. Conclusion

We introduced the Drift-Selection Ratio and applied it to 200 bacterial genomes spanning 8 phyla. In 93.5% of genomes, tRNA gene copy number distributions are indistinguishable from a neutral birth-death process (KS $p > 0.05$). Codon usage frequency explains only 3.8% of copy number variance, compared to 31% for local duplication rate. Cyanobacteriota are the sole phylum with consistent departures from neutrality (median DSR = 1.87), driven by excess copies of photosystem-associated tRNA isoacceptors. These results falsify the strong translational selection hypothesis for tRNA gene copy numbers in the majority of prokaryotes.

References

[1] F. Ardell and D. Andersson, "On the evolution of redundancy in the genetic code," Journal of Molecular Evolution, vol. 63, pp. 12--22, 2006.

[2] T. Ikemura, "Codon usage and tRNA content in unicellular and multicellular organisms," Molecular Biology and Evolution, vol. 2, pp. 13--34, 1985.

[3] P. Sharp, E. Bailes, R. Grocock, J. Peden, and R. Sockett, "Variation in the strength of selected codon usage bias among bacteria," Nucleic Acids Research, vol. 33, pp. 1141--1153, 2005.

[4] S. Percudani, A. Pavesi, and S. Ottonello, "Transfer RNA gene redundancy and translational selection in Saccharomyces cerevisiae," Journal of Molecular Biology, vol. 268, pp. 322--330, 1997.

[5] H. Goodenbour and T. Pan, "Diversity of tRNA genes in eukaryotes," Nucleic Acids Research, vol. 34, pp. 6137--6146, 2006.

[6] S. Kanaya, Y. Yamada, M. Kinouchi, Y. Kudo, and T. Ikemura, "Codon usage and tRNA genes in eukaryotes: Correlation of codon usage diversity with translation efficiency and with CG-dinucleotide usage as assessed by multivariate analysis," Journal of Molecular Evolution, vol. 53, pp. 290--298, 2001.

[7] A. Wald and J. Glass, "Codon usage in Mycoplasma and Spiroplasma," Nucleic Acids Research, vol. 25, pp. 2958--2966, 1997.

[8] M. dos Reis, R. Savva, and L. Wernisch, "Solving the riddle of codon usage preferences: a test for translational selection," Nucleic Acids Research, vol. 32, pp. 5036--5044, 2004.

[9] M. Hahn, T. De Bie, J. Stajich, C. Nguyen, and N. Cristianini, "Estimating the tempo and mode of gene family evolution from comparative genomic data," Genome Research, vol. 15, pp. 1153--1160, 2005.

[10] M. Nei and A. Rooney, "Concerted and birth-and-death evolution of multigene families," Annual Review of Genetics, vol. 39, pp. 121--152, 2005.

[11] T. Lowe and P. Chan, "tRNAscan-SE On-line: integrating search and context for analysis of transfer RNA genes," Nucleic Acids Research, vol. 44, pp. W54--W57, 2016.

[12] L. Hershberg and D. Petrov, "Selection on codon bias," Annual Review of Genetics, vol. 42, pp. 287--299, 2008.

[13] J. Felsenstein, "Phylogenies and the comparative method," The American Naturalist, vol. 125, pp. 1--15, 1985.

[14] T. Gogakos et al., "Characterizing expression and processing of precursor and mature human tRNAs by hydro-tRNAseq and PAR-CLIP," Cell Reports, vol. 20, pp. 1463--1475, 2017.