key: cord-0992286-a79kercb
authors: Greco, Samuele; Gerdol, Marco
title: Independent acquisition of short insertions at the RIR1 site in the spike N-terminal domain of the SARS-CoV-2 BA.2 lineage
date: 2022-04-11
journal: bioRxiv
DOI: 10.1101/2022.04.11.487924
sha: 2ff0a55cedbbf751f19c4dcc770251c445845cbd
doc_id: 992286
cord_uid: a79kercb

Although the SARS-CoV-2 variants BA.1 and BA.2 share over 30 non-synonymous substitutions in the spike glycoprotein, they show several unique mutations that were likely acquired after the split between these two major omicron lineages. One of the most intriguing mutations associated with BA.1 is the presence of the inserted tripeptide Glu-Pro-Glu within the N-terminal domain. While the functional implications of this insertion are still unclear, several other SARS-CoV-2 lineages had previously independently acquired similarly short insertions at the very same site, named RIR1. We have previously identified this site, located approximately between codon 212 and codon 216, as a hotspot of insertions, which usually involve small nucleotide sequences including three or four codons. Here we show that similar insertion events have independently occurred at least 13 times in early 2022 within the BA.2 lineage, being occasionally associated with significant community transmission. One of these omicron sublineages, characterized by a Ser-Gly-Arg insertion in position 212, is responsible of over 2% of all SARS-CoV-2 cases recorded in Denmark, as of early April 2022. Molecular surveillance data highlight a slow but steady growth compared with the parental BA.2 lineage in all Danish regions, suggesting that the RIR1 insertion may confer a selective advantage. We report the identification of other currently circulating BA.2 sublineages showing similar insertions, whose spread should be therefore carefully monitored in the upcoming months.

of vaccine-elicited or convalescent sera, omicron displayed extremely significant immune evasion properties, to the point that some authors have suggested the necessity to consider this VOC as belonging to a distinct SARS-CoV-2 serotype (Simon-Loriere and Schwartz, 2022) . The highly divergent antigenic properties of omicron derive from the presence of a very high number of S1 mutations, several of which have been shown to mediate ACE2 recognition and/or immune escape (Cameroni et al., 2022; Mannar et al., 2022) . Curiously, many of these mutations occur at sites subjected to strong purifying selection both in previous SARS-CoV-2 variants and in bat coronaviruses, suggesting that they may cooperatively interact to mitigate slight fitness deficits associated with individual substitution, with a significant impact on the structural and biological properties of the spike trimer (Martin et al., 2022a) . The highly unusual mutational patter of omicron may, to some extent, also explain the preferential use of the TMPRSS2-independent and cathepsin-mediated endosomal cell entry route (Pia and Rowland-Jones, 2022) , as well as the lower fusogenity and altered cell tropism displayed by this variant (Meng et al., 2022) . The extreme phylogenetic divergence between omicron and all the other previous SARS-CoV-2 variants raised important questions concerning its origins, which might be either sought in complex selective processes occurring in immunocompromised patients, or in an unidentified animal reservoir (Mallapaty, 2022) .

Although omicron includes three lineages, only BA.1 and BA.2 have caused major outbreaks in human, whereas BA.3 has been responsible for less than one thousand documented cases worldwide. Although BA.1 and BA.2 share several S1 mutations, they display a remarkable sequence divergence, which is roughly comparable with the divergence observed among earlier VOCs. One of the most puzzling mutations only found in one out of the two omicron sublineages is a short insertion of three codons, which encode the tripeptide Glu-Pro-Glu at position 215, within the spike NTD of BA.1. Nearly all VOCs and VOIs carry small spike deletions compared with Wuhan-Hu-1, which may either confer enhanced antibody escape or act as permissive mutations, counterbalancing the fitness cost of otherwise deleterious RBD mutations (Meng et al., 2021) . These small deletions independently occurred on multiple occasions during SARS-CoV-2 evolution at Recurrent Deletion Regions (RDRs) located within the NTD, marking one of the most significant signatures of convergent evolution documented during the first two years of the pandemics. On the other hand, the occurrence of spike insertions has been much rarer and consequently subjected to far less study. We have previously shown that insertions of short sequence stretches, usually encoding 3 or 4 amino acids, had been independently acquired over 50 times at the very same NTD site, named Recurrent Insertion Region 1 (RIR1).

Before the emergence of BA.1, RIR1 insertions have been documented in alpha, gamma and delta, where they did not lead to significant community spread, and with two lineages (A.2.5 and B.1.214.2) that became dominant in a few Central America and Central African countries for a few months during 2021 (Gerdol et al., 2022) . Although the functional role of RIR1 insertions is presently unclear, their convergent acquisition by multiple viral lineages most certainly suggests the need for increased monitoring, due to the possibility of associated fitness advantages.

While the frequency of observation of BA.1 is progressively declining as the result of the higher transmissibility of BA.2, we here show that a few BA.2 sublineages have independently acquired RIR1 insertions. Molecular surveillance data suggests that one of these sublineages in particular, hereby named BA.2+ins(L), which carries a Ser-Arg-Gly insertion at position 212 and is relatively widespread in Denmark, may display a moderate growth advantage compared with the parental BA.2 lineage.

All SARS-CoV-2 genome sequence data deposited in GISAID (Shu and McCauley, 2017) up to April 11 th , 2022, was screened, looking for sequences belonging to the BA.2 lineage and showing insertions between positions 210 and 220 in the S gene. Genomes displaying unidentified amino acids (i.e. X) in such insertions were flagged as the likely product of mis-assembly and discarded. Similarly, genomes assigned to the BA.2 lineage but carrying an "EPE" insertion in position 215, previously described as one of the lineage-defining mutations of the sister lineage BA.1, were also flagged as likely misclassified and were therefore discarded. All resulting genomes were grouped based on the inserted nucleotide sequence, whose exact position and phase of insertion were defined based on the alignment with the reference SARS-CoV-2 genome sequence (Wuhan-Hu-1, NCBI accession ID: NC_045512.2) and a reference BA.2 genome (EPI_ISL_8128502). Following the nomenclature scheme proposed in our previous work (Gerdol et al., 2022) , each of these groups was labeled with progressive Roman numerals, following their chronological order of identification, starting from insertion L.

The molecular surveillance data from Denmark was subjected to further analysis, due to the local spread of one of the BA.2 sublineages carrying an insertion at RIR1, i.e. BA.2+ins(L). All genome data and associated metadata, deposited up to April 11 th , 2022, were downloaded from GISAID. The share of BA.2 carrying insertion L was estimated, both at a national and at a regional scale (i.e. by separately taking into account Syddanmark, Sjaelland, Nordjylland, Midtjylland and Hovestaden) up to April 1 st , 2022. The frequencies of observation of BA.2-ins(L) relative to all BA.2 sequenced genomes were reported the 7-day moving average, with 95% confidence intervals.

The phylogenetic tree was prepared as follows. All available SARS-CoV-2 genome sequence data deposited in GISAID up to April 5 th , 2022, paired with metadata, were downloaded. Genomes with low sequencing coverage, those including stretches of undefined nucleotides (i.e. N) longer than 25 and those associated with incomplete collection date were discarded. Separate fasta files were created for viral genomes belonging to the following lineages: alpha; beta; gamma; delta; BA.1; BA.1.1; BA.2 (excluding genomes with RIR1 insertions; BA.2+ins(L) and BA.2+ins(LIII) (see Table 1 and Figure 1) ; sequences not belonging to any VOC from the year 2020 Particular attention was placed in an additional filtering step, which involved the fasta file including the BA.2 sequences bearing RIR1 insertions, aimed at discarding the sequences that may possibly derive from low quality assemblies. Briefly, retained sequences had to display all characterizing BA.2 mutations, defined as those present in > 80% of all sequenced genomes, as of April 5 th , 2022. A few BA.2+ins(L) sequences displaying a small NSP1 deletion (codons 81-87) were also discarded due to its high frequency of observation in some BA.2 sublineages. Each fasta file was subsequently subsampled to include a representative group of sequences for the aforementioned lineages. In detail the final sequence dataset included: 50 alpha, beta, gamma, delta, BA.1 and BA.1.1 genomes; 150 genomes belonging to non-VOC lineages of 2020; 600 BA.2 genomes lacking RIR1 insertions; all available BA.2+ins(L) and BA.2+ins(LIII) passing quality filtering (i.e. 531 and 182, respectively). The Wuhan-Hu-1 strain genome (NCBI accession ID: NC_045512.2) was added as reference sequence

All the sampled sequences were merged and aligned with mafft v7.490 (Katoh and Standley, 2013) , configured for closely related viral genomes as suggested by the authors. Gapped regions corresponding to insertions were not discarded. The best-fitting substitution model for the alignment was estimated with modeltest-ng v0.1.7 (Darriba et al., 2020) . The Augur pipeline (Huddleston et al., 2021) was then used to produce the sampling date-refined phylogenetic tree, with the Wuhan-Hu-1 genome set as a reference for rooting. In order to take into account the monophyletic origin of the insertions a guide tree was used in the tree building step. All scripts and command used in this analysis are available at https://github.com/54mu/SARS-CoV-2_BA.2_RIR1insertions. The tree visualization was obtained with Auspice (Hadfield et al., 2018) . Extended methods are available at https://gist.github.com/54mu/5e858e0f4aa6a7236cf5849e6ec371c7

In total, we identified 13 independent events of insertion at RIR1 associated with the BA.2 lineage. These, based on the naming scheme we have previously proposed (Gerdol et al., 2022) , were chronologically ordered, from the least to the most recent one and assigned Roman numerals, from insertion L to insertion LXII ( Table 1) .

As all previously described RIR1 insertions, those associated with BA.2 involve the acquisition of either three (seven cases), four (five cases), or five (a single case) codons. However, their position was somewhat different compared with those that arose in other variants throughout 2020 and 2021 ( Figure 1A) . Indeed, four insertions (i.e. LI, LII, LIV and LX) occurred within codon 212 at phase II and four others (i.e. L, LIII, LVII and LVIII) occurred within codon 213 at phase 0. None of the 49 RIR1 insertions described in our previous work was observed in these locations (Gerdol et al., 2022) . On the other hand, a single insertion (i.e. LVI) was found within codon 215, at phase II, a placement consistent with previously described RIR1 insertions. The exact location of the four remaining BA.2-associated RIR1 insertions could not be ascertained due to the possibility of multiple alternative sequence alignments with the BA.2 reference.

Due to the unusual placement of such insertions in the spike CDS, the inserted tri-, tetra-or penta-peptides in BA.2 would be usually placed between Leu212 and Val213 (albeit this residue is replaced by Gly in BA.2), and only in two cases between Arg214 and Asp215 ( Figure 1B) . Because of the out-of-frame nature of several insertions, in four cases the residue flanking the insertion at the N-terminal side (i.e. Leu212) was replaced with Phe (insertions LI, LII, LIV and LX) ( Figure 1B) . As of note, due to of the presence of the Val213Gly substitution in BA.2, the placement of most of these insertions by GISAID was often offset by a single amino acid (e.g. insertion L is reported as ins213GRG instead of ins212SGR). Consequently, a non-synonymous substitution was often incorrectly called in position 213 (e.g. Val213Ser instead of Val213Gly for insertion L).

Only three out of the 49 RIR1 insertions described in our previous work were associated with a significant community spread, either locally, or globally: insertion III ( The other 46 RIR1 insertions described in our previous work remained confined to small clusters counting no more than a few dozen sequenced cases, even though on some occasions (e.g. lineage B.1.639) undocumented community transmission likely continued for several months. Overall, these observations pointed out that the association between RIR1 insertions and the acquisition of increased fitness over previous variants could not be generalized. Due to the lack of structural overlap between RIR1 and important NTD antibody epitopes and in light of the frequent association between these insertions and nonsynonymous mutations targeting the RBD, in our interpretation RIR1 insertions could be interpreted as compensatory mutations. As previously hypothesized for example for the H69/V70 deletion found in many VOCs and VOIs (Meng et al., 2021) , small insertions at RIR1 may be able to mitigate slight negative fitness costs associated with the presence of non-synonymous RBD mutations affecting immune escape (Gerdol et al., 2022) . From this perspective, the limited community spread of most SARS-CoV-2 variants bearing RIR1

insertions was not surprising, due to the important role played by intra-host selection in shaping spike mutational patterns, since these selective forces are not necessarily expected to lead to viral lineages with increased transmissibility. The previously documented acquisition of RIR1 insertions through in vitro experiments, both in response to several passages in Vero cell cultures (Shiliaev et al., 2021) and after the treatment with monoclonal antibodies (Committee for Medicinal Products for Human Use, 2021), might be consistent with this view. Hence, the acquisition of a selective advantage by of RIR1 insertions-carrying viral lineages may be conditioned by the presence or absence of other spike non-synonymous mutations associated with negative fitness costs.

In total, 11 out of the 13 BA.2-associated RIR1 insertions were only associated with a very small number of cases and were therefore unlikely to be linked with BA.2 sublineages with enhanced fitness ( Table 1) , even though it cannot be excluded that some of these, identified very recently, may further spread in the coming weeks. Nevertheless, two exceptions were immediately evident, i.e. insertion L and insertion LIII ( Table 1) .

The two BA.2 sublineages carrying such insertions will be hereafter named BA.2+ins(L) and BA.2+ins(LIII).

BA.2+ins(L), first detected on January 7 th , 2022, is characterized by the presence of an in-frame TCCGGCAGA insertion between codon 212 and 213, which determines the acquisition of the tripeptide Ser-Gly-Arg in position 212 (Figure 1) . accounts for approximately 1% of all BA.2 cases recorded in this region (Figure 2F) .

At a national level, the comparative analysis of the incidence of infections due to BA.2 and BA.2+ins(L) allows to draw a preliminary estimate of the weekly growth advantage displayed by the latter, which would be in the range of 20%. Although the absolute number of infections linked with this sublineage presently remains relatively low (i.e. with the number of daily national cases recently dropping below 5,000, those linked with BA.2+ins(L) might be in the range of 50-100 per day), this moderate growth trend has been continuing for over two months, following consistent patterns in multiple Danish regions. It would be therefore important to clarify whether this apparent growth advantage over the parental BA.2 lineage is driven by an alteration of the structure and function of the spike protein or rather by other epidemiological factors, such as founder effects, as previously observed for other SARS-CoV-2 lineages in the past .

The second relevant BA.2-related sublineage with an insertion at RIR1 was BA.2+ins(LIII). This sublineage, first detected on February 5 th , 2022 in England, displays a longer insertion compared with BA-2+ins(L), i.e.

ACAGTAGGAGGA, which is also found in-frame between codon 212 and 213, and encodes the tetrapeptide Thr-Val-Gly-Gly (Figure 1) . Curiously, neither these two insertions, nor the other 11 listed in Table 1 were associated with other relevant S1 non-synonymous substitutions, questioning their previously hypothesized "compensatory" role in stabilizing otherwise slightly deleterious mutations with significant effect on immune escape (Gerdol et al., 2022) . In detail, BA.2+ins(L) did not show any other associated non-synonymous mutations, neither in the S gene, nor in other genomic regions, leaving the Ser-Gly-Arg as the only sublineage-defining mutation present.

The only other relevant mutation shared by all BA.2+ins(L) genomes was synonymous, i.e. C22791T (also located within the S gene), was not phylogenetically informative, as it was found in nearly one third of all BA.2 genomes. On the other hand, BA.2+ins(LIII) also displayed the presence of Thr1213Ser, a S2 nonsynonymous conservative substitution rarely observed in SARS-CoV-2 prior to the emergence of BA.2 (only 62 cases had been recorded prior to December 2021), but occasionally found in BA.2. The only other two spike non-synonymous mutations linked with BA.2 RIR1 insertions were G75D in insertion LI, Y144del in insertion LX and T1136S in insertion LXII, even though several point mutations were found in other SARS-CoV-2 genes ( Table 1) .

In summary, most BA.2 associated RIR1 insertions appear to have not occurred in parallel with other mutations targeting the RBD, the polybasic site or known important epitopes recognized by neutralizing antibodies. The two major omicron lineages share a large fraction of the non-synonymous substitutions found in the S gene (i.e. 21), but each of the two also carry a dozen additional mutations not shared with the sister lineage. Some of the mutations shared by BA.1 and BA.2, if taken individually, would have been predicted to significantly decrease viral fitness. Nevertheless, they are thought to cooperatively interact by positive epistasis, mitigating the negative fitness costs associated with their presence and possibly altering some biological properties of the spike protein (Martin et al., 2022b) . Insertion XLI (i.e. ins215EPE), which is included among the mutations found in BA.1 but not shared with BA.2, undoubtedly emerged after the split between the two lineages and its presence may therefore provide a fitness advantage in combination with any of the other mutations exclusively found in BA.1. Nevertheless, we did not observe any overlap between the presence of RIR1 insertions in BA.2 and the acquisition, by convergent evolution, of any spike mutation typical of BA.1. A different interpretation would therefore be needed to explain the emergence of RIR1 insertions in BA.2, as well as the apparent slight fitness advantage associated with BA.2+ins(L). While, to date, no conclusive explanation has been provided concerning the benefits associated with the presence of ins214EPE in BA.1, it is very likely that this small insertion, in combination with other NTD non-synonymous mutations, had a significant impact on the structure of the spike protein in BA.1. These alterations to the original primary sequence of the SARS-CoV-2 spike protein might have possibly contributed to the enhanced inter-domain and inter-subunit packing, as well as to the higher accessibility for ACE2 binding, driven by its higher predisposition towards an open configuration (Ye et al., 2022) . From this perspective, the acquisition of RIR1 insertions such as L, LIII and others, might be also interpreted as a trait conferring increased stability to the BA.2 spike trimer. Hence, the independent acquisition of different short sequence stretches at RIR1 in the BA.2 lineage might simply depend on the presence of pre-existing BA.2-associated mutations and be completely independent on the acquisition of additional non-synonymous mutations. Undoubtedly, further structural studies would be required to clarify the effects of the acquisition of RIR1 insertions on a BA.2 genetic background in terms of the impact on the structural conformation and stability of spike trimers.

A relevant consequence of the lack of other lineage-defining mutations in most BA.2 sublineages carrying RIR1 insertions is the inherent difficulty in their assignment through phylogenetic inference. As a matter of fact, indels are usually simply treated as missing characters and not appropriately weighted-in by most maximum likelihood approaches, including those currently used by the agur pipeline, i.e. IQ-TREE (Nguyen et al., 2015) and faststree (Price et al., 2010) , and this might lead to significant statistical inconsistencies, as previously noted by other authors (McTavish et al., 2015; Warnow, 2012 showing the full mutational profile of BA.2 were discarded, all BA.2+ins(L) were correctly placed in a monophyletic clade (Figure 3) . On the other hand, this adjustment was not required to achieve a good support of the monophyly of the BA.2+ins(LIII) sublineage, which as mentioned above also carries the spike Thr1213Ser mutation. This highlights the fact that similar difficulties in future lineage assignments could arise for other SARS-CoV-2 lineages whose only distinctive mutations compared with their parental variant are either insertions or deletions. Hovestaden (panel F). The graphs report the 7-day moving average frequencies of observations of BA.2+ins(L)

genomes, relative to all BA.2 genomes, with 95% confidence intervals. 

The authors declare they have no competing interests.

Evolutionary origins of the SARS-CoV-2 sarbecovirus lineage responsible for the COVID-19 pandemic

Broadly neutralizing antibodies overcome SARS-CoV-2 Omicron antigenic shift

GlaxoSmithKline use of sotrovimab (VIR-7831/GSK4182136) for the treatment of COVID-19 (No. EMA/304600/2021)

ModelTest-NG: A New and Scalable Tool for the Selection of DNA and Protein Evolutionary Models

Estimated transmissibility and impact of SARS-CoV-2 lineage B.1.1.7 in England

Emergence of a recurrent insertion in the N-terminal domain of the SARS-CoV-2 spike glycoprotein

Evolutionary Arms Race between Virus and Host Drives Genetic Diversity in Bat Severe Acute Respiratory Syndrome-Related Coronavirus Spike Genes

Nextstrain: real-time tracking of pathogen evolution

SARS-CoV-2 variants, spike mutations and immune escape

Spread of a SARS-CoV-2 variant through Europe in the summer of 2020

SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor

Augur: a bioinformatics toolkit for phylogenetic analyses of human pathogens

Human neutralizing antibodies elicited by SARS-CoV-2 infection

MAFFT multiple sequence alignment software version 7: improvements in performance and usability

Tracking Changes in SARS-CoV-2 Spike: Evidence that D614G Increases Infectivity of the COVID-19 Virus

Angiotensin-converting enzyme 2 is a functional receptor for the SARS coronavirus

Structural basis and functional analysis of the SARS coronavirus nsp14-nsp10 complex

Where did Omicron come from? Three key theories

SARS-CoV-2 Omicron variant: Antibody evasion and cryo-EM structure of spike protein-ACE2 complex

Selection analysis identifies unusual clustered mutational changes in Omicron lineage BA.1 that likely impact Spike function

Selection analysis identifies clusters of unusual mutational changes in Omicron lineage BA.1 that likely impact Spike function

Twisted trees and inconsistency of tree estimation when gaps are treated as missing data -The impact of model mis-specification in distance corrections

Altered TMPRSS2 usage by SARS-CoV-2 Omicron impacts infectivity and fusogenicity

IQ-TREE: A Fast and Effective Stochastic Algorithm for Estimating Maximum-Likelihood Phylogenies

Omicron entry route

Reduced sensitivity of SARS-CoV-2 variant Delta to antibody neutralization

FastTree 2 -Approximately Maximum-Likelihood Trees for Large Alignments

Increased risk of SARS-CoV-2 reinfection associated with emergence of Omicron in South Africa

Resurgence of COVID-19 in Manaus, Brazil, despite high seroprevalence

Enhanced fusogenicity and pathogenicity of SARS-CoV-2 Delta P681R mutation

Natural and Recombinant SARS-CoV-2 Isolates Rapidly Evolve In Vitro to Higher Infectivity through More Efficient Binding to Heparan Sulfate and Reduced S1/S2 Cleavage

GISAID: Global initiative on sharing all influenza data -from vision to reality

Towards SARS-CoV-2 serotypes?

Modelling conformational state dynamics and its role on infection for SARS-CoV-2 Spike protein variants

Structure, Function, and Antigenicity of the SARS-CoV-2 Spike Glycoprotein

Standard maximum likelihood analyses of alignments with gaps can be statistically inconsistent

Cryo-EM structure of a SARS-CoV-2 omicron spike protein ectodomain

Cryo-electron microscopy structures of the N501Y SARS-CoV-2 spike protein in complex with ACE2 and 2 potent neutralizing antibodies

The authors would like to acknowledge the tremendous efforts made by the clinicians, researchers and public health authorities that allowed the collection of SARS-CoV-2 genome data and made sequence data available in a timely manner though GISAID, as well as the great efforts made by the developers of nextstrain to assist researchers in SARS-CoV-2 evolution studies.

This research did not receive any specific grant from funding agencies in the public, commercial, or not-forprofit sectors.