Alexander G Lucaci1, Jordan D Zehr1, Stephen D Shank1, Steven Weaver1, Dave Bouvier2, Han Mei2, Anton Nekrutenko2, Sergei L Kosakovsky Pond1
- Institute for Genomics and Evolutionary Medicine, Temple University, Philadelphia, Pennsylvania, USA
- Department of Biochemistry and Molecular Biology, The Pennsylvania State University, University Park, PA, USA
In order to monitor the continued evolution of the SARS-CoV-2 genome, we developed a rapid assessment tool (RASCL) designed to investigate the nature and extent of selective forces acting on viral genes in clade sequences. We describe our method and its application to the Delta and Kappa lineages in a recent Virologial post. Here, we present an early application of our method to the emerging Lambda variant of interest, which was designated on June 14, 2021 (VOI, WHO Int). We find continued evolution in Spike outside of clade-defining mutations. As well, we report on extensive evolution of the Nucleocapsid gene within this variant of interest.
We investigated the nature and extent of selective forces acting on viral genes in C.37 clade sequences. As sequence collection, sampling, and deposition into public data repositories is an ongoing process, we present an early analysis of the emerging lambda variant. By using complete linkage distance clustering (method described in TN93) we selected a median of 80 C.37 (Lambda) sequences from available sequences to represent genomic diversity. A global reference dataset from publicly available sequences via Virus Pathogen Database and Analysis Resource (ViPR) (Pickett et al., 2012) is included as background.
For C.37 clades sequences, we find strong statistical support for multiple modes of selection: diversifying positive (FEL/MEME), negative (FEL), and directional (FADE). There is also evidence of coevolution among sites within individual gene segments (BGM), and at a number of sites clade evolution occurs with selective forces that are different from those acting on background sequences (Contrast-FEL). These results are summarized in Figure 1 and can be fully explored in the following interactive notebook: Link to full RASCL results for C.37 clade sequences.
As we described in our previous Virological post on the Delta and Kappa variants, we again notice that a key observation is that selection appears to be operating on many sites that are not a part of the clade-defining signature mutation set, implying that there may be ongoing diversification and adaptation in this clade. Additionally, as temporal changes in the frequency of amino-acid residues at sites outside of clade defining mutations occur, these may become seeds for future variant lineages.
Figure 1. Schematic for genome-level results for the C.37/Lambda clade. The color key is defined as follows: Green = MEME (episodic diversifying selection) on All branches (Track 1, or the outermost track) and MEME results on Internal branches (Track 2). Purple = Clade defining sites from Pangolin. Orange = S1 Functional sites shown to influence antibody binding. Yellow = Contrast-FEL results on Internal branches. Blue = FADE results, a test for directional selection. Grey box = FEL results for negatively selected sites. Grey lines = BGM results for coevolving sites.
We focus our attention on positively selected sites located in the Spike gene. We identify 8 sites subject to episodic positive selection in the C.37 clade (Table 1), with only 2 of those appearing on the canonical clade signature list (452, 490). The trends in most common haplotypes at the remaining (non-clade-defining) six positions are shown in Figure 2.
|#||Coordinate (SARS-CoV-2)||Gene/ORF||Codon (in gene/ORF)||# of selected branches||p-value||q-value|
Table 1. Sites in the Spike gene from C.37 sequences which show evidence for episodic selection. Sites which are bolded indicate clade-defining mutations from Pangolin.
Figure 2. Temporal trends of the substitution combinations at all sites represented in Table 1 in the Spike gene for C.37 sequences in 2021. The symbol . denotes the reference residue at that site . Figures like this can be generated using Trends in mutational patterns across SARS-CoV-2 Spike enabled by data from / Sergei Pond / Observable
We notice the high level of selection activity in the Nucleocapsid gene within C.37 sequences (Figure 1 and Table 2). Therefore, we focus our attention on positively selected sites located in the Nucleocapsid gene. We identify 19 sites subject to episodic positive selection in the C.37 clade (Table 2), with none of those appearing on the canonical clade signature list. Of note, the Nucleocapsid 366I mutation seems to occur at high frequency (~49% of C.37 sequences in 2021) at a number of branches on the phylogenetic tree (Table 2). This mutation have played a role in the early evolution of C.37 (see Figure 3.) and allowed other mutations in Nucleocapsid to be tolerated. The trends in most common haplotypes at the nineteen positions are shown in Figure 3.
|#||Coordinate (SARS-CoV-2)||Gene/ORF||Codon (in gene/ORF)||# of selected branches||p-value||q-value||Physiochemical Properties|
|14||29371||N||367||2||0.000329097||0.00846248||bipolar, volume, composition, charge|
Table 2. Sites in the Nucleocapsid gene in C.37 sequences which show evidence for episodic selection. We also note (where appropriate) on physiochemical properties influencing selection on specific sites (per the PRIME method).
Figure 3. Temporal trends of the substitution combinations at all sites reported in Table 2 in the Nucleocapsid gene in C.37 sequences. . denotes the reference residue. Figures like this can be generated using Trends in mutational patterns across SARS-CoV-2 Spike enabled by data from / Sergei Pond / Observable
Ongoing monitoring of emergent VOI and VOC can detect adaptive mutations before they rise to high frequency, and help establish their relationship to key clinical parameters including pathogenicity and transmissibility. Additionally, continued evolution within a particular clade may form the foundation for a subclade with further functional sites of interest. Based on current information for the Spike gene from Stanford Coronavirus Antiviral & Resistance Database (CoVDB) we included several annotations with clinical relevance.
SARS-CoV-2 C.37 Spike gene, sites of interest from Table 1.
- Spike / 12, 20, 72: No annotation is currently available.
- Spike / 246, located in Spike subdomain NTD, and is involved in epitope binding to mAbs
- Spike / 253, located in Spike subdomain NTD, and an epitope binding to mAbs
- Spike / 452, located in Spike subdomain RBD, and is involved in epitope binding to mAbs
- Spike / 490, located in Spike subdomain RBD, and is involved in ACE2 binding
- Spike / 826, No annotation is currently available.
The RASCL application is designed to accept two inputs: a “query” whole genome sequence dataset, and a “background” whole genome sequence dataset. Users are requested to input a whole genome dataset (fasta, unaligned), which will be defined as the “query” sequence dataset, and can be typically retrieved from data repositories such as GISAID. Typically, a background sequence dataset is assembled from publicly available sequences via the ViPR database, and is provided by default unless otherwise specified. However, our application also tolerates the usage of another clade of interest to be used as “background”, which allows one to directly interrogate between-clade evolutionary dynamics. Whole genome datasets are broken down into their respective coding sequences from the NCBI SARS-CoV-2 reference annotation. The gene set includes structural (Spike, Matrix, Nucleocapsid, Envelope) and non-structural (leader, nsp2, nsp3, nsp4, 3C, nsp6, nsp7, nsp8, nsp9, nsp10, helicase, exonuclease, endornase, ORF3a, ORF6, ORF7a, ORF8, RdRp, methyltransferase) genes. Individual genes from query and background datasets are processed from here onwards. Our procedure includes several quality control steps including sanitizing fasta input headers and striking ambiguous codons from an alignment. In order for our analyses to complete rapidly, alignments are subsampled using genetic distances provided by the TN93 (Tamura and Nei 1993) algorithm. A combined (query + background) alignment is created and with sequences kept from the background dataset divergent enough to be useful for subsequent selection analysis. Inference of a maximum likelihood phylogenetic tree (RAxML-NG, Kozlov et al., 2019) is performed on the combined dataset. Trees are labelled for “query” and “background” branches. Analyses SLAC (Pond and Frost, 2005), BGM (Poon et al., 2008), FEL (Pond and Frost, 2005), MEME (Murrell et al., 2012), aBSREL (Smith et al., 2015), BUSTEDS (Wisotsky et al., 2020), RELAX (Wertheim et al., 2015), CFEL (Pond et al., 2021) are performed with state of the art molecular evolution models from the HyPhy (Pond et al., 2020) suite of bioinformatics tools. We report on sites of interest which together may represent regions of ongoing evolution within a specific clade. Results are combined into easily readable and web-accessible JSON files used for web processing. Visualization of results is done with full-feature interactive notebooks through ObservableHQ.
Data retrieval was performed through GISAID on July 7th 2021 with 2,127 sequences downloaded for the C.37 variant of interest. An additional parameter for the search included: ‘low coverage excl.’.
RASCL is freely available via a dedicated Github repository:
Another implementation is also available via a user-friendly Galaxy workflow.
We thank the global community of health-care workers and scientists who have worked tirelessly to face the pandemic head-on. We are also thankful for the GISAID (Elbe et al., 2017) submitters from across the world. Additionally, we thank members of the Datamonkey / HyPhy and Galaxy teams for their continued assistance in the development of this application.
- Elbe, S., and Buckland-Merrett, G. (2017) Data, disease and diplomacy: GISAID’s innovative contribution to global health. Global Challenges, 1:33-46. DOI:10.1002/gch2.1018 PMCID: 31565258
- Pickett BE, Sadat EL, Zhang Y, Noronha JM, Squires RB, Hunt V, Liu M, Kumar S, Zaremba S, Gu Z, Zhou L, Larson CN, Dietrich J, Klem EB, Scheuermann RH. ViPR: an open bioinformatics database and analysis resource for virology research. Nucleic Acids Res. 2012 Jan;40(Database issue):D593-8. doi: 10.1093/nar/gkr859. Epub 2011 Oct 17. PMID: 22006842; PMCID: PMC3245011.
- Alexey M Kozlov, Diego Darriba, Tomáš Flouri, Benoit Morel, Alexandros Stamatakis, RAxML-NG: a fast, scalable and user-friendly tool for maximum likelihood phylogenetic inference, Bioinformatics, Volume 35, Issue 21, 1 November 2019, Pages 4453–4455, https://doi.org/10.1093/bioinformatics/btz305
- Sergei L Kosakovsky Pond, Art FY Poon, Ryan Velazquez, Steven Weaver, N Lance Hepler, Ben Murrell, Stephen D Shank, Brittany Rife Magalis, Dave Bouvier, Anton Nekrutenko, Sadie Wisotsky, Stephanie J Spielman, Simon DW Frost, Spencer V Muse (2020) HyPhy 2.5—A Customizable Platform for Evolutionary Hypothesis Testing Using Phylogenies. 1 Molecular Biology and Evolution 37.1 (2020): 295-299.
- Sergei L. Kosakovsky Pond and Simon D. W. Frost (2005) Not So Different After All: A Comparison of Methods for Detecting Amino Acid Sites Under Selection Molecular Biology and Evolution 22(5): 1208-1222
- Art F. Y. Poon, Fraser I. Lewis, Simon D. W. Frost, Sergei L. Kosakovsky Pond, Spidermonkey: rapid detection of co-evolving sites using Bayesian graphical models, Bioinformatics, Volume 24, Issue 17, 1 September 2008, Pages 1949–1950, https://doi.org/10.1093/bioinformatics/btn313
- Ben Murrell, Joel O. Wertheim, Sasha Moola, Thomas Weighill, Konrad Scheffler and Sergei L. Kosakovsky Pond (2012) Detecting Individual Sites Subject to Episodic Diversifying Selection PLoS Genetics 8(7): e1002764
- M. D. Smith, J. O. Wertheim, S. Weaver, B. Murrell, K. Scheffler and S. L. Kosakovsky Pond Less Is More: An Adaptive Branch-Site Random Effects Model for Efficient Detection of Episodic Diversifying Selection Molecular Biology and Evolution 32: 1342-1353
- Wisotsky SR, Kosakovsky Pond SL, Shank SD, Muse SV. Synonymous Site-to-Site Substitution Rate Variation Dramatically Inflates False Positive Rates of Selection Analyses: Ignore at Your Own Peril. Mol Biol Evol. 2020 Aug 1;37(8):2430-2439. doi: 10.1093/molbev/msaa037. PMID: 32068869; PMCID: PMC7403620.
- Wertheim JO, Murrell B, Smith MD, Kosakovsky Pond SL, Scheffler K. RELAX: detecting relaxed selection in a phylogenetic framework. Mol Biol Evol. 2015 Mar;32(3):820-32. doi: 10.1093/molbev/msu400. Epub 2014 Dec 23. PMID: 25540451; PMCID: PMC4327161.
- Kosakovsky Pond SL, Wisotsky SR, Escalante A, Magalis BR, Weaver S. Contrast-FEL-A Test for Differences in Selective Pressures at Individual Sites among Clades and Sets of Branches. Mol Biol Evol. 2021 Mar 9;38(3):1184-1198. doi: 10.1093/molbev/msaa263. PMID: 33064823; PMCID: PMC7947784.