Hunting for the SARS-CoV-2 Spike Glycoprotein Ancestor

Objective: In December 2019 a novel coronavirus (SARS-CoV-2) that is causing the current COVID-19 pandemic was identied in Wuhan, China. Many questions have been raised about its origin and adaptation to humans. In the present work we performed a genetic analysis of the Spike glycoprotein (S) of SARS-CoV-2 and other related coronaviruses (CoVs) isolated from different hosts in order to trace the evolutionary history of this protein and the adaptation of SARS-CoV-2 to humans. Results: Based on the sequence analysis of the S gene, we suggest that the origin of SARS-CoV-2 is the result of recombination events between bat and pangolin CoVs. The hybrid SARS-CoV-2 ancestor jumped to humans and has been maintained by natural selection. Although the S protein of RaTG13 bat CoV has a high nucleotide identity with the S protein of SARS-CoV-2, the phylogenetic tree and the haplotype network suggest a non-direct parental relationship between these CoVs. Moreover, it is likely that the basic function of the receptor-binding domain (RBD) of S protein was acquired by the SARS-CoV-2 from the MP789 pangolin CoV by recombination and it has been highly conserved.

Viruses, like other pathogenic microorganisms, are subject to different evolutionary forces that allow them to adapt and jump from one host to another. Mutation, gene ow and recombination generate genetic variation that is maintained or removed by natural selection and gene drift. One of the key factors in the process of adaptation by CoVs to different species is their ability to infect cells of their new host.
CoVs are able to infect cells through a membrane glycoprotein called Spike (S) [1]. This protein contains an amino-terminal subunit 1 (S1) and a carboxyl-terminal subunit 2 (S2) [9]. The S1 has a Receptor-Binding Domain (RBD) that allows the virus to bind to different receptors on cells of its different hosts.
For example, MERS-CoV binds to dipeptidyl dipeptidase 4 (DPP4) [10], while SARS-CoV and SARS-CoV-2 bind to the receptor for the Angiotensin-Converting Enzyme 2 (ACE2) [4,11]. Although these last two viruses bind to the same receptor, their RBDs present variations in the amino acids involved in the ACE2 recognition site [12].
Previous reports show that protein S is one of the most variable regions of the SARS-CoV-2 genome that is under natural selection and where recombination signals have been recorded [5,7,13], suggesting that this protein changes constantly and plays an important role in adapting to humans. In our study, we performed a genetic analysis of the protein S of SARS-CoV-2 and related CoVs isolated from different hosts in order to trace the evolutionary history of protein S and the adaptation of SARS-CoV-2 to humans.

Methods
The sequences of the Glycoprotein Spike gene from 76 CoVs isolated from different hosts and 148 clinical isolates of SARS-CoV-2 were retrieved from NCBI's GenBank [14] and GISAID [15] databases (Additional le 1). The most representative sequence, the China Wuhan H1 sequence, was used as a reference sequence in the gures of this study. The sequences were aligned and analyzed with MAFFT v7.3 [16] and edited using BioEdit v7.2 [17] and Jalview v2.11 [18]. The nucleotide substitution model of the sequence set was determined with the j Modeltest2 program [19] and a phylogenetic tree was obtained using the Mr. Bayes program [20] under a General Time-Reversible substitution plus proportion of invariable sites and rate of variation across sites (GTR+I+G). We ran 5 Markov Chains Monte Carlo for 500,000 generations and discarded 25% of the initial trees. The consensus tree was edited using the Figtree v1.4.3 program [21]. A haplotype network was built using the PopArt v1.7 program [22]. Nucleotide similarity analysis was carried out using the SimPlot v3.5 program with a 100 bp window at 30 bp step and Kimura 2-parameter model [23]. The DnaSP v5.1 program [24] was used to perform the McDonald-
Since we did not nd any phylogenetic incongruence re ected on the tree, suggesting a lack of recombination between clusters, we focused more speci cally on the cluster where SARS-CoV-2 was found (A7) (Fig. 1a). We analyzed the genealogical relationships between the CoVs that comprise cluster A7. The haplotype network showed the formation of a loop between the group of pangolins and 4 hypothetical ancestral haplotypes (Fig. 1b), suggesting recombination within this cluster. Furthermore, the network suggests that 4 of the isolates analyzed here (SARS-CoV-2, RaTG13, MP789 and RmYN02) diverged from these 4 hypothetical ancestors (Fig. 1b). Despite the fact that SARS-CoV-2 and RaTG13 CoV share a genomic nucleotide identity of 96.2% [4], and an S gene nucleotide identity of 93.15% [7], the divergence showed in the phylogenetic tree and in the haplotype network rules out a direct parental relationship between these two isolates ( Fig. 1a and b).
In 2019, various pangolin CoVs were isolated, among which the isolate MP789 CoV is the most interesting because it shares a nucleotide similarity of 85%-92% with SARS-CoV-2, and 90% with RaTG13 CoV [7]. The similarity analysis of the S nucleotide sequences of cluster A7 shows a mosaic similarity pattern across the S gene between SARS-CoV-2, RaTG13, MP789 and RmYN02, which suggests a probable ancestral genetic exchange between the 4 hypothetical ancestors of these CoVs (Fig. 1c). The most notable differences between SARS-CoV-2 and the rest of the CoVs S gene were found in the RBD, indicating a hybrid zone between RaTG13 and MP789 CoVs in this region (Fig. 1c). This result suggests a probable ancestral cross-species recombination between bat and pangolin CoVs.
S protein is thought to be under natural selection and plays an important role in cross-species transmission [5,[26][27][28]. A recent study reported negative selection in the S gene when SARS-CoV-2 was compared with RaTG13 and a group of pangolin CoVs [26]. We performed an MK test between SARS-CoV-2, RaTG13 and MP789, the results of which showed that between SARS-CoV-2 and RaTG13 CoV there were more synonymous (dS) than non-synonymous (dN) substitutions, indicating negative selection (NI>1). Whereas, between SARS-CoV-2 and MP789 CoV the contrary was found, indicating positive selection (NI<1) ( Table 1). The negative selection predicted for SARS-CoV-2 is due to its high similarity to RaTG13 CoV, therefore, the xation of dN substitutions are not favored. On the other hand, the incongruences found in the pangolin CoV results compared to a previous study [26] could be due to differences in the strategies and methods used.
The S protein RBD plays a key role during the infection process of SARS-CoV-2 to human cells because it contains the six amino acids (L455, F486, Q493, S494, N501 and Y505) that are essential for e cient binding of SARS-CoV-2 to ACE2 [12]. SARS-CoV-2 RBD shows a closer similarity to MP789 CoV RBD (96.8% homology) than to RaTG13 CoV RBD (89.56% homology) [7]. Interestingly, we found that 26 of 33 dN substitutions between SARS-CoV-2 and RaTG13 CoV were located in the RBD, while 7 of 505 dN substitutions between SARS-CoV-2 and MP789 CoV were also located in the RBD (Table 1). This indicates that in this region, MP789 CoV has suffered less dN changes than RaTG13 CoV when compared with SARS-CoV-2. Since only one polymorphism was detected in RBD in the 148 SARS-CoV-2 sequences, the MK test did not determine any value, suggesting that this is a highly conserved region.
The comparison between SARS-CoV-2 and MP789 CoV RBD shows that they share the 6 amino acids that are essential for binding to ACE2 receptor, while in RaTG13 CoV these amino acids are missing (Fig.   2). These results could indicate that both the pangolin and humans have similar ACE2 at the interacting domain with S protein, as reported by others [29,30]. As a consequence, the ACE2 binding sites and the region in general should be conserved (70% homology), being su cient for the interaction to take place.

ND Undetermined
A genetic feature that makes SARS-CoV-2 more infectious is the fact that the S protein harbors an insert of 12-nucleotides between the S1 and S2 subunits that encode for a polybasic cleavage site (RRAR) that is recognized by furins (Fig. 2). This cleavage site is related with an increased e ciency of entry during infection [31,32]. Nevertheless, this insertion is not present in all betacoronaviruses, like in SARS-CoV [13,30,31]. However, the human HKU1 CoV and MERS-CoV have variants of polybasic insertions that are also recognized by furins [33][34][35]. The presence of these polybasic insertions have been seen to increase the pathogenicity of viruses, such as in avian in uenza [36][37][38], MERS-like CoV [39], and in bovine CoV [40].
We also found that RmYN02 CoV has an insert in the same position as that in SARS-CoV-2, but it is not a polybasic cleavage insertion (-AAR). There have been suggestions that instead, it could be the product of recombination between wild bat CoVs [13]. On the other hand, several experiments have shown that this polybasic cleavage site is acquired and xed during the serial passage of CoVs in cell cultures or in animals [36,41]. The aforementioned leads to two possible explanations for the polybasic cleavage insertion in SARS-CoV-2 and the role in its adaptation to humans: 1) the ancestor of SARS-CoV-2 acquired it in a host, went through a recombination process in an unidenti ed intermediary host, and then jumped to humans, or 2) it was acquired in humans during a cycle of human to human transmission that helped its adaptation and virulence process. Rambaut et al., [42] determined that the most recent common ancestor of SARS-CoV-2 appeared in November 2019 and proposes that the virus had enough time to acquire the insert during transmission between humans.
Whether the ancestral CoV that gave rise to SARS-CoV-2 came from a bat or a pangolin is still not yet known, but based on the analysis of the S gene performed in this study, we suggest that it is more likely to have come from a bat. However, the region essential for human ACE2 (RBD) binding is a hybrid between RaTG13 and PM789 CoVs, and it is likely that the basic function of the RBD was acquired by the SARS-CoV-2 from the MP789 CoV by recombination with an ancestral CoV, which had a RaTG13 genomic background. Subsequently, the SARS-CoV-2 ancestor jumped to humans where the S protein has been maintained by positive selection and where the RBD has been highly conserved. This illustrates the complexity of CoV cross-species infection dynamics and the relevance of genetic exchange between CoVs. In the future, molecular epidemiological surveillance studies in wild isolates will be vital for identifying genetic changes in viruses that could result in novel adaptations to humans and thereby, enabling the development of another pandemic, such as the one we are currently experiencing.

Limitations
Here we mentioned two characteristics in S protein as key factors in the adaptation and infectivity of SARS-CoV-2, the presence of six amino acids essentials for binding to ACE2 receptor and the polybasic cleavage insert. Nevertheless, the knowledge about this virus has increased very fast and new possible features involved in virus adaptation to human beings has been reported. To fully understand the origin and adaptation of this virus, it would be necessary to encourage the extensive sampling of wild animals associated with coronaviruses and perform deeper and wider genomic analyses. Availability of data and materials

Abbreviations
The datasets used and/or analyzed during the current study are available in the Additional le 1.

Competing interests
The authors declare that they have no competing interests.

Funding
This work was supported by DGAPA-PAPIIT grant IN213816.
Author's contributions RME, AFA and LSM contributed to the study design; AFA performed the genetic and phylogenetic analyses and the interpretation of results; AFA, LSM and GD wrote the manuscript. RME and AC did the critical review of the manuscript. All authors read and approved the nal version of the manuscript.

Figure 1
Genetic comparison of S gene of SARS-CoV-2 and other coronaviruses a. Phylogenetic tree depicting the relationship between S gene of 76 coronavirus (CoVs) isolates from different hosts and 148 SARS-CoV-2 isolates. CoVs were classi ed by genus (α, β, γ and δCoVs) and each cluster is shown in a different color.
The A7 cluster in red, grouped the SARS-CoV-2 and other related CoVs, each one is marked with a colored circle. The number presented in the nodes indicates posterior probabilities, and the scale bar represents the average number of nucleotide substitutions. All nodes in the A7 cluster had posterior probabilities of 1.0. b. Haplotype network of S gene, the colored circles represent each CoV according to the A7 cluster in the tree. The PcoV_GX_P5L, P5E, P1E, P4L and P2V were named as Pangolin groups. Small black circles represent hypothetical or unsampled haplotypes. The small black circles, marked with a (?), indicate the hypothetical ancestors of SARS-CoV-2_Wuhan_Hu-1, RaTG13, RaTG13, MP789 and YmN02. The numbers indicate the number of mutations that separate each haplotype. c. Patterns of nucleotide sequence identity in the S gene for SARS-CoV-2_Wuhan_Hu-1, RaTG13, MP789, YmN02 and Pangolin group. Similarity was calculated within sliding window of 300 bp moving with steps of 30 bp.