Evolutionary shortcuts via multi-nucleotide substitutions and their impact on natural selection analyses

Abstract

Inference and interpretation of evolutionary processes - in particular of the types and targets of natural selection affecting coding sequences, are critically influenced by the assumptions built into statistical models for such analyses. If certain aspects of the substitution process (even when they are not of direct interest) are presumed absent or are modeled with too crude of a simplification, estimates of key model parameters can become biased - often systematically, and lead to poor statistical performance. Here, we performed a detailed characterization of how modeling instantaneous multi-nucleotide (or multi-hit, MH) substitutions impacts dN/dS based inference of episodic diversifying selection at the level of the entire alignment. The inclusion of MH reduces the rate (1.37-fold or 26.8%) at which positive selection is called based on the analysis of N = 9,861 empirical data-sets, while offering significantly better statistical fit to sequence data in 8.37% of cases. Through additional simulation studies, we show that this reduction is not simply due to loss of power because of additional model complexity. After a detailed examination of 21 benchmark alignments and a new high-resolution analysis showing which parts of the alignment provide support for positive selection, we reveal that MH substitutions occurring along shorter branches in the tree are largely responsible for discrepant results in selection detection. Our results add to the growing body of literature which examines decades-old modeling assumptions and finds them to be problematic for biological data analysis. Because multi-nucleotide substitutions have a significant impact on natural selection detection even at the level of an entire gene, we recommend that routine selection analysis of this type consider their inclusion. To facilitate this procedure, we developed a simple model testing selection detection framework able to screen an alignment for positive selection with two biologically important confounding processes: synonymous rate variation, and multi-nucleotide instantaneous substitutions.

Rapid evolution of genes with anti-cancer functions during the origins of large bodies and cancer resistance in elephants

“…These data are consistent with previous observations synonymous substitution rate variation is common Wisotsky et al, 2020) and that shorter branches tend to be more impacted by multiple hits than longer branches (Lucaci et al, 2023b(Lucaci et al, , 2021 such as the elephant lineage. The best-fitting ABSREL model identified 393 genes with a class of sites with ω>1 at P≤0.05 in the elephant lineage (Figure 1C), whereas the best-fitting BUSTED model identified 427 genes with a class of sites with ω>1 at P≤0.05 in the elephant lineage (Figure 1C); 324 genes were identified as positively selected by both methods (Figure 1D).…”

“…The ABSREL and BUSTED models are variants of the branch-sites random effect likelihood (BSREL) model of coding sequence evolution (Pond et al, 2011a;Smith et al, 2015), which allows for variation in d N /d S (ω) rates across lineages and sites; however, while ABSREL selects the optimal number of rate categories for each gene, BUSTED imposes three rate categories such that ω 1 ≤ω 2 ≤ω 3 (Figure 1B). Both models also accommodate synonymous rate variation across sites [S] , double and triple multi-nucleotide mutations per codon [MH] (Lucaci et al, 2021), and both synonymous rate variation and multi-nucleotide mutations [SMH]; failure to account for synonymous rate variation (Dunn et al, 2019;Rahman et al, 2021;Wisotsky et al, 2020) and multi-nucleotide mutations can lead to widespread false positive inferences of positive selection (Dunn et al, 2019;Lucaci et al, 2023a;Nozawa et al, 2009;Suzuki, 2008;Venkat et al, 2018;Yang and Reis, 2011). We used small sample corrected Akaike Information Criterion values Including synonymous rate variation across sites improved model fit for all genes; however, only 31 genes had significant evidence for multi-nucleotide mutations (Figure 1C).…”

“…Goodman et al (2009), for example, used the free ratio model in PAML (Yang, 2007) to identify positively selected genes in the elephant lineage using a dataset of 11 species and 5,714 protein-coding genes; they identified 67 genes with d N /d S ≥1.01, which were enriched in mitochondrial functions but none of the GO terms we identified (Supplementary dataset 1). There are several important differences between Goodman et al (2009) and our study: 1) The free ratio model does not test if the d N /d S ratio in a particular lineage is significantly different than one, rather it tests if d N /d S rate variation is different across lineages (Yang, 1998); 2) Unlike ABSREL and BUSTED, the free ratio model does not allow for rate variation across sites, d N /d S >1 is only inferred if the average d N /d S across all sites is greater than 1; and 3) The method implemented in PAML does not account for synonymous rate variation across sites or multiple simultaneous codon substitutions, which can bias inferences of both d S (Dunn et al, 2019;Rahman et al, 2021;Wisotsky et al, 2020) and d N (Dunn et al, 2019;Lucaci et al, 2023a;Nozawa et al, 2009;Suzuki, 2008;Venkat et al, 2018;Yang and Reis, 2011). While we tested 81.6% (4,663) of the genes tested by Goodman et al (2009), only one (CD3G) was inferred to have been positively selected in our analyses; therefore, it is unsurprising that there was little overlap in our GO enrichment results.…”

“…The copyright holder for this preprint this version posted February 29, 2024. ; https://doi.org/10.1101/2024.02.27.582135 doi: bioRxiv preprint does not allow for rate variation across sites, d N /d S >1 is only inferred if the average d N /d S across all sites is greater than 1; and 3) The method implemented in PAML does not account for synonymous rate variation across sites or multiple simultaneous codon substitutions, which can bias inferences of both d S (Dunn et al, 2019;Rahman et al, 2021;Wisotsky et al, 2020) and d N (Dunn et al, 2019;Lucaci et al, 2023a;Nozawa et al, 2009;Suzuki, 2008;Venkat et al, 2018;Yang and Reis, 2011). While we tested 81.6% (4,663) of the genes tested by Goodman et al (2009), only one (CD3G) was inferred to have been positively selected in our analyses; therefore, it is unsurprising that there was little overlap in our GO enrichment results.…”

Evolution towards simplicity in bacterial small heat shock protein system

“…This indicates that this change possibly had an adaptive character in the ancestral background. It is important to note, however, that models used for the analysis do not account for variation of synonymous substitution rate and multinucleotide substitution events, which in some cases might lead to false positive results (Lucaci et al , 2023). Because of that, the alternative hypothesis that substitutions Q66H and G109D occurred in the common ancestor of Erwiniaceae due to genetic drift, enabling the subsequent loss of ibpB gene, cannot be fully discounted.…”