Computational analysis of affinity dynamics between the variants of SARS-CoV-2 spike protein (RBD) and human ACE-2 receptor | Virology Journal

Mutations in RBD and structure stability of SARS-CoV-2 virus VOCs

Prior bioinformatic investigations have highlighted the pivotal role of the spike protein’s Receptor-Binding Domain (RBD) in susceptibility to novel mutations and its critical function in host binding [16]. To comprehensively assess the RBD within the SARS-CoV-2 spike protein, a multiple sequence alignment was employed. This procedure involved the comparison of RBD sequences from SARS-CoV-2 strains originating from the top five nations significantly impacted by the pandemic, as documented by the World Health Organization (WHO), where substantial fatalities were recorded. The reference sequence for this alignment was human SARS-CoV-2 virus reported in Wuhan, China, with accession number YP_009724390.1. The Alpha variant (B.1.1.7), first identified in the United Kingdom, displayed the N501Y mutation within its RBD region. This mutation has raised concerns about potential implications, such as increased viral load, extended viral persistence, and an elevated risk of fatality when compared to the reference sequence [17] (Table S1, Supplementary Data). In South Africa, the Beta variant (B.1.351) presented three distinctive mutations within its RBD region, including K417N, E484K, and N501Y, which have been associated with a notable reduction in antibody neutralization [18] (Table S2, Supplementary Data). Similarly, the Gamma variant (P.1) discovered in Brazil exhibited mutations K417T, E484K, and N501Y in the RBD region, with an estimated 2.6-fold increase in transmissibility [19] (Table S3, Supplementary Data). The Delta variant (B.1.617.2), originating in India, was characterized by L452R, E484Q, and T478K mutations within the RBD, leading to a significant surge in cases and hospitalizations during India’s second wave [19] (Table S4, Supplementary Data). Lastly, the Omicron variant (B.1.1.529) in South Africa presented an array of prominent mutations, including G339D, S371L, S373P, S375F, K417N, N440K, G446S, T478K, E484A, Q498R, G496S, Q498R, N501Y, and Y505H. These mutations have raised concerns regarding the potential for increased reinfection risk compared to other VOCs [20] (Table S5, Supplementary Data). Enhanced virus infectivity can stem from improved host receptor binding stability. The strength of each unique RBD mutation in each variant was assessed (Table 1). Missense mutations within RBDs can impact their binding affinity to hACE2R. Using the DUET webserver, missense mutations were evaluated in the RBD protein and protein stability changes were recorded. The mean ΔΔG values in monomer stability in each full-length RBD-residue position ranged from − 1.229 kcal/mol in G496S to 0.876 kcal/mol in N440K. The results depicted in Fig. 1 highlight specific amino acid (AA) positions that play a significant role in either stabilizing or destabilizing the Receptor-Binding Domain (RBD) of the spike glycoprotein, with corresponding percentages listed in Table 6. Noteworthy glycoprotein mutations, including G339D, S371L, N440K, S477N, T478K, E484A, E484K, Q493R, Q498R, and Y505H, demonstrated pronounced stabilizing effects on the protein structure. The analysis revealed that among the 22 unique AA substitutions, only 47.6% exhibited stabilizing effects, while 52.3% had destabilizing effects on the Spike Protein. A particularly notable observation was the high stabilizing effect exerted by the N440K transformation in all Omicron variants, indicated by a ΔΔG value of 0.876 kcal/mol. Comparisons between RBD mutations in Beta and Gamma (E484K; ΔΔG: 0.654 kcal/mol) and Omicron subvariants (Q483R; ΔΔG: 0.495 kcal/mol) revealed that Beta and Gamma variants showed significant stabilizing effects. Conversely, Delta and B.1.617.2 (L452R) and Beta and Omicron variants (K417N) exhibited notable destabilizing impacts, as indicated by low ΔΔG values when predicting missense mutations’ effects on protein stability (see Fig. 1).

The stability changes of each mutation of different variant of concern and the omicron subvariants and their mean stability values

Ab initio 3D modelling of the RBD region to know conformational changes

The primary objective of our study was to assess the trRosetta algorithm’s capability to consistently predict protein structures, particularly in scenarios where no similar structures are known or available in public databases. To achieve this, we generated predicted models for all the variants of concern using trRosetta, and these models were then superimposed onto the original Wuhan variant using the “super” command in PyMol. Upon examination, it was observed that for the majority of pairs involving trRosetta-generated models and the original Wuhan variant, the Root-Mean-Square Deviation (RMSD) ranged from 0.69 to 0.77 Å. A lower RMSD indicates a more accurate alignment of the Receptor-Binding Domains (RBDs). Consequently, we considered the trRosetta-generated models suitable for further analysis (see Table 2). Qualitative assessment of structural changes, observed through superimposing the trRosetta-generated naturally occurring variants (Fig. 2A-a), revealed that the mutations predominantly affected the interaction interface situated on the bottom section of the RBD (Fig. 2A-b, A-c). The Receptor-Binding Motif (RBM), identified as the main functional motif in the RBD, highlighted specific residues (Fig. 2B-a). Further, structural variations in comparison to Delta and Omicron sequences were noted (Fig. 2B-b, B-c), consistent with findings by Bhoumick et al., 2022 [21]. These 3D modeling outcomes unequivocally demonstrate a substantial conservation of the spike protein’s RBD. Comparing the variants to the wild type (Wuhan sequence), variations exhibit a higher alpha-helix structure, while the secondary structure prediction indicates fewer extended strands and a reduced presence of random coil shapes. The projected increase in alpha-helices suggests that beta strands are less prone to mutations than alpha helices, emphasizing the structural resilience of the spike protein’s RBD.

A-a The trRosetta-generated Receptor-Binding Domain (RBD) structures of all SARS CoV 2 variants were super-imposed on 6M0J. A-b Structural changes between the RBM of 6M0J and the Delta sequence. A-c Structural changes of the RBM between 6M0J and the Omicron sequence. B-a RBM of 6M0J; (1), (2), (3), (4) are the changes observed. B-b Structural changes of the RBM between 6M0J and the Delta sequence. B-c Structural changes of the RBM between 6M0J and the Omicron sequence

Comparative analysis of the binding ability of mutant RBD complexes to hACE2

The binding mechanism, affinities between Receptor-Binding Domains (RBDs) and the hACE2 receptor for various SARS-CoV-2 variants, including the original Wuhan variant, were assessed through the HADDOCK 2.0 server. The docking scores, cluster sizes, and dissociation constants (KD) were presented in the Tables 3, 4, 5, 6 and 7 for the USA, India, France, Germany, and Brazil, with corresponding docking complexes in Figs. 2D and 3 DIMPLOT representations of the interactions in Fig. 4. The original Wuhan variant’s predicted HADDOCK score was − 89.5 kcal/mol, indicating the lowest binding affinity. In the context of the five most affected countries (USA > India > France > Germany > Brazil), the Omicron BA.1.1 variant showed a remarkable 1.5-fold increase in affinity towards hACE2 in the USA (-130.1 kcal/mol), surpassing other Omicron subvariants in USA. Consistent trends emerged in binding affinity values, with the BA.1.1 variant at − 17.4 kcal/mol, indicating enhanced binding compared to the original Wuhan variant at − 12.9 kcal/mol. These findings align with data from the global COVID-19 tracker (ourworldindata.org/covid-19). Figure S4 [A, B and C]( Supplementary Data) shows COVID-19 cases and fatalities by the World Health Organization (WHO). In early 2022 Omicron outbreak (BA.1.1), USA reported the highest cases, hospitalizations, and elevated mortality, supporting the hypothesis of BA.1.1’s heightened transmissibility due to its exceptional binding affinity. In the Indian context, a comprehensive evaluation of docking scores revealed that the BA.1.1 Omicron variant displayed a 1.5-fold increase in binding affinity at -114 kcal/mol towards the human ACE2 receptor (hACE2). This heightened affinity closely mirrored that of the Delta variant, at -112.4 +/- 3.7 kcal/mol. Importantly, the computed docking scores aligned with the corresponding binding affinity with the BA.1.1 variant showing an affinity of -15.2 kcal/mol, not significantly different from the Delta variant B.1.617.2’s -14.8 kcal/mol. Additionally, an analysis of global COVID-19 tracker data reinforced these findings. Figure S5 [A, B and C] ( Supplementary Data) depicted the extensive incidence of confirmed COVID-19 cases and fatalities, sourced from World Health Organization (WHO). During mid-2021 Delta variant B.1.617.2 outbreak, India witnessed a devastating initial COVID-19 wave marked by significant fatalities, a surge in cases, hospitalizations, and elevated mortality. Remarkably, vaccines were not widely accessible during this period, emphasizing the pivotal role of the Delta variant’s exceptional binding affinity in driving its heightened transmissibility and virulence, culminating in a significant impact of SARS-CoV-2 in India. In contrast, in the French context, the Delta variant B.1.617.2 had the highest HADDOCK score, registering at -102.3 +/- 21.4 kcal/mol, closely resembling the Gamma variant P.1’s HADDOCK score at -98.1 +/- 9.5 kcal/mol. These scores directly correlated with binding affinity values of -14.4 kcal/mol and − 14.3 kcal/mol. Furthermore, analysis of global COVID-19 tracker data supported these findings. Figure S6 [A, B and C] ( Supplementary Data) graphically displayed the extensive confirmed COVID-19 cases and fatalities during the late 2021 and early 2022 Delta variant B.1.617.2 outbreak in France, offering compelling evidence for the Delta variant’s substantial transmissibility and virulence during this period in the country. Similar findings were observed in both Germany and Brazil, mirroring those seen in the USA. The BA.1.1 Omicron variant exhibited a 1-fold increase in binding affinity compared to the original Wuhan strain with HADDOCK scores ranging from − 96.7 +/- 4.6 kcal/mol to -116.5 +/- 2.0 kcal/mol and binding affinity from − 16.5 kcal/mol and − 15.4 kcal/mol, respectively. Global COVID-19 tracker data aligns with our study’s findings, as depicted in Figure S7 [A, B and C] ( Supplementary Data) and Figure S8 [A, B and C] ( Supplementary Data). This representation illustrates the comprehensive landscape of confirmed COVID-19 cases and confirmed deaths, sourced from WHO records. During the early 2022 Omicron variant BA.1.1 outbreak, both Germany and Brazil experienced significant COVID-19 outbreaks, characterized by a substantial loss of lives, increased confirmed cases, hospitalizations, and elevated mortality rates. This underscores the notion that the exceptional binding affinity of the Omicron variant significantly contributed to its peak transmissibility and virulence during this period in both Germany and Brazil. The ongoing evolution of SARS-CoV‐2 has resulted in considerable genetic diversity within the viral population. The virus gains entry into host cells by interacting with the hACE‐2 receptor, a process facilitated by the binding of the Receptor-Binding Domain (RBD) of the spike protein (SP). Mutations in the spike protein have the potential to enhance viral entry, and alterations in the RBD can impact viral infectivity and stability [22]. Consequently, the interaction interface between the spike protein and the hACE‐2 receptor has become a key focus in the development of novel drugs [23]. This study aims to investigate how mutations in the RBD influence binding affinity and stability across various Variants of Concern (VOCs), ranging from the wild type (WT) to the Omicron variant. Notably, variants originating in the United Kingdom (Alpha‐ B.1.1.7), South Africa (Beta‐ B.1.351), and Brazil exhibit significant structural differences and altered binding properties [24]. As an example, the study identified a specific mutation (N501Y) in the RBD of the Alpha variant compared to the wild type (WT). This mutation represents one of several alterations that researchers are exploring to understand how changes in the RBD may impact the virus’s ability to bind to and enter host cells. The research seeks to provide insights into the dynamic relationship between RBD mutations and their effects on viral infectivity, shedding light on potential targets for drug development and therapeutic interventions. The alteration introduced in the Alpha variant increased its binding capacity by 7.9% compared to the wild type (WT) [25], potentially contributing to the higher reported cases of COVID-19 associated with the Alpha variant. Additionally, we noted variations in salt bridges and hydrogen bonds among different variants, such as the absence of a salt bridge involving Lys478 in Alpha, Beta, and Gamma complexes but its presence in the BA.5 variant. This suggests a potential link between these mutations and variant transmissibility [26]. In our study, we observed the formation of a salt bridge between K417 (RBD) and D38 (hACE‐2) in the WT RBD-hACE2 complex, confirming the results obtained by Wang et al. [27]. However, the substitution of K417 by N417 and T417 in Beta and Gamma, respectively, resulted in the loss of the salt bridge. In contrast to observations by Han et al. [28] and Li et al. [29], where the absence of salt bridges resulted in reduced binding affinity for hACE‐2, our results contradicted their findings. We observed a high binding energy for the Gamma variant compared to the wild type (WT), even in the absence of salt bridge formation. This underscores the significance of other intermolecular interactions in enhancing binding capacities. The binding affinity (∆G‐kcal/mol) plays a crucial role in determining whether complex formation occurs, holding a vital role in controlling interactions. With the exception of Delta and Omicron, other Variants of Concern (VOCs) displayed higher KD and weaker binding affinity values, aligning with the results reported by Khan et al. [24]. In our study, BA.2.12.1 exhibited the lowest binding affinity (-13.6 kcal/mol), indicating high stability of the complex compared to other variants. BA.2.12.1 also demonstrated the lowest KD value (9.8 × 10 − 11) and high hACE‐2 affinity compared to other VOCs and Omicron subvariants. A reduced KD value signifies a strong binding between the spike protein-RBD and hACE‐2, establishing a greater affinity between the receptor and ligand [30]. It is noteworthy that among Alpha, Beta, Gamma, and Omicron variants, the N501Y mutation in their RBDs is a primary concern, as it is one of the active residues directly interacting with hACE2R.The HADDOCK scores highlight that both the Delta variant (-99 kcal/mol) and the Omicron subvariants (OSVs) exhibit more robust binding abilities compared to the wild type (WT) (-87 kcal/mol). The Receptor Binding Domain (RBD) of Omicron (B.1.1.529) has undergone notable mutations. Particularly, the BA.1 and BA.2 OSVs, aligning with findings by Han et al. [31], are more prevalent and display the lowest score (-130.1 kcal/mol) in the case of BA.2. This data reinforces the higher transmissibility observed in the Omicron variant and its subvariants [32]. Studies indicate that the BA.2 Omicron subvariant is 1.5 times more contagious than BA.1, substantiating its high binding energy [33]. Additionally, Omicron BA.1 is predicted to be significantly more contagious than Wuhan‐Hu‐1 and Delta variants, primarily attributed to RBD mutations N440K, T478K, and N501Y [34]. Omicron demonstrates an overall infectivity level that surpasses the ancestral SARS‐CoV‐2 variation and other subsequent variants, including Delta, mainly attributed to its extensive mutations in the Receptor Binding Domain (RBD) region [35]. The docking scores for Delta and Alpha exhibit similarities to those of BA.5, potentially due to amino acid substitutions at L452R and T478K. Notably, these mutations, along with F486V, are absent in BA.2. Similar to Delta, the 478th residue in the RBD of BA.5 forms a salt bridge and hydrogen bond with the human angiotensin‐converting enzyme receptor (hACE2R). In the context of the Delta variant, the presence of the Arg452 mutation is also observed in BA.5. In BA.5, this mutation contributes to the formation of two hydrogen bonds and one salt bridge, leading to increased binding affinity compared to Delta. Notably, despite BA.5 having lower binding affinity than BA.2, it demonstrates heightened pathogenicity and efficient transmission in humans, suggesting that factors beyond binding energy play a role in viral transmission. When examining the active site residues of the Receptor Binding Domains (RBDs), residue 493 displayed hydrogen bond formations in Other SARS-CoV-2 variants (OSVs) complexes, except in BA.2 and BA.5. The formation of hydrogen bonds is crucial for the docking score value, highlighting their significance in the overall interaction dynamics. Hence the binding stability of hACE2 and RBD complexes were in the order of Omicron BA.1.1 > Omicron BA.5.2.1 > Omicron BA.2.12.1 > Delta > Omicron-B.1.529.1 > Gamma > Alpha > Beta > Original Wuhan strain due to the increase of the KD.

a Represents the binding interface of mutant complexes and a surface representation. b Offers the binding interface and stick model of the fundamental hydrogen bonding interactions of the mutant. (Chain A represents hACE2R, and chain B represents the RBD of each variant). Original Wuhan strain (1), Alpha variant(2), Beta variant(3), Gamma variant(4), Omicron variant BA.1.1 which has the highest binding affinity among all variants in USA (5), Omicron variant BA.1.1 which has the highest binding affinity among all variants in India (6), Delta variant B.1.617.2 which has the highest binding affinity among all variants in France (7), Omicron variant BA.1.1 which has the highest binding affinity among all variants in Germany (8), Omicron variant BA.1.1 which has the highest binding affinity among all variants in Brazil (9)

RBD’s 2D DIMPOT interactions representation, including hydrogen interactions in mutant complexes. Original Wuhan strain (A), Omicron variant BA.1.1 which has the highest binding affinity among all variants in USA (B), Omicron variant BA.1.1 which has the highest binding affinity among all variants in India (C), Delta variant B.1.617.2 which has the highest binding affinity among all variants in France (D), Omicron variant BA.1.1 which has the highest binding affinity among all variants in Germany (E), Omicron variant BA.1.1 which has the highest binding affinity among all variants in Brazil (F)