Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Open Access

Peer-reviewed

Research Article

Molecular dynamics study on the strengthening behavior of Delta and Omicron SARS-CoV-2 spike RBD improved receptor-binding affinity

  • Kanchanok Kodchakorn,

    Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Project administration, Validation, Visualization, Writing – original draft, Writing – review & editing

    Affiliation Thailand Excellence Center for Tissue Engineering and Stem Cells, Department of Biochemistry, Faculty of Medicine, Chiang Mai University, Chiang Mai, Thailand

    ⨯
  • Prachya Kongtawelert

    Roles Conceptualization, Data curation, Funding acquisition, Resources, Supervision, Writing – review & editing

    prachya.k@cmu.ac.th

    Affiliation Thailand Excellence Center for Tissue Engineering and Stem Cells, Department of Biochemistry, Faculty of Medicine, Chiang Mai University, Chiang Mai, Thailand

    ⨯

Abstract

The COVID-19 pandemic caused by a virus that can be transmitted from human to human via air droplets has changed the quality of life and economic systems all over the world. The viral DNA has mutated naturally over time leading to the diversity of coronavirus victims which has posed a serious threat to human security on a massive scale. The current variants have developed in a dominant way and are considered “Variants of Concern” by the World Health Organization (WHO). In this work, Kappa (B.1.617.1), Delta (B.1.617.2), and Omicron (B.1.1.529) variants were obtained to evaluate whether naturally occurring mutations have strengthened viral infectivity. We apply reliable in silico structural dynamics and energetic frameworks of the mutated S-RBD protein for ACE2-binding to analyze and compare the structural information related to the wild-type. In particular, the hotspot residues at Q493, Q498, and N501 on the S-RBD protein were determined as contributing factors to the employment stability of the relevant binding interface. The L452R mutation induces an increment of the hydrogen bonds formed by changing the Q493 environment for ACE2 binding. Moreover, the Q493K exchange in Omicron enables the formation of two additional salt bridges, leading to a strong binding affinity by increased electrostatic interaction energy. These results could be used in proposing concrete informative data for a structure-based design engaged in finding better therapeutics against novel variants.

Citation: Kodchakorn K, Kongtawelert P (2022) Molecular dynamics study on the strengthening behavior of Delta and Omicron SARS-CoV-2 spike RBD improved receptor-binding affinity. PLoS ONE 17(11): e0277745. https://doi.org/10.1371/journal.pone.0277745

Editor: Nicolas N. Nassar, Cincinnati Children’s Hospital Medical Center, UNITED STATES

Received: June 14, 2022; Accepted: November 2, 2022; Published: November 17, 2022

Copyright: © 2022 Kodchakorn, Kongtawelert. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: All relevant data are within the paper and its Supporting Information files.

Funding: This work was supported by grants from the Post-Doctoral Fellowship, grant number: R000030567, by CMU Presidential Scholarship, Chiang Mai University, Chiang Mai, Thailand. However, the founder had no role in study design, data collection and analysis, decision to publish, or preparation of this manuscript.

Competing interests: The authors have declared that no competing interests exist.

Introduction

The severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) is driven by newly emerging variants nearly two years into the pandemic worldwide that has deeply affected the entire population, pushing young health professionals into new careers outside the medical field. A series of coronavirus variants has emerged in different countries, some of which have caused significant outbreaks. The variants of Alpha [1], and more recently Delta [2], have had the greatest global reach. The variants of Beta [3], Gamma [4], and Lambda [5], did not become dominant mutations in most regions of the world, although they did cause substantial outbreaks in Southern Africa (SA). Nowadays, the growing presence of a novel dominant strain of the virus is a reminder the pandemic is far from over. Meanwhile, the personal choices made by people during the pandemic, especially concerning vaccinations, have shaken the foundations of health systems.

On November 26, 2021, the WHO designated one recent variant B.1.1.529, called “Omicron”, a “Variant of Concern (VOC)” with high potential to spread rapidly across the world [6, 7]. This variant is spreading significantly faster than the Delta variant and has recently appeared in SA [8]. At the same time, Omicron is quickly overtaking Delta to be the dominant strain of SARS-CoV-2 in the US and many other countries [8–10]. Both variants are complicating the work of the medical professionals who are responsible for assisting patients with breathing challenges. Thus, it is a time of year when numerous respiratory illnesses emerge. The rapid emergence of Omicron in the background of the immune response implies that the virus may have evolved to escape neutralizing antibody responses by Beta-specific serum [11]. Omicron is a highly divergent variant with a high number of mutations (thirty mutations in the spike (S) protein, 15 of which occur in the receptor-binding domain (RBD), as well as three small deletions and one minor insertion), some of which may be associated with humoral immune escape potential and higher transmissibility. Particular potential hotspot residues for the mutations are the S-RBD region as the site of ACE2 binding that dramatically improves RBD-ACE2 binding affinity, causing fast dissemination in human populations [12].

The S-RBD protein is an essential part of SARS-CoV-2 as it mediates interaction with human cells and is the target for therapeutic neutralizing antibodies [13, 14]. Several studies have proposed the mechanism of the ACE2-bound S-RBD binding in the viral entry process [15, 16]. The mutation’s effects on S-RBD binding affinity and viral expression are found to correlate with an amino acid residue as the key hotspots clustering at the binding interface and has been well documented [17–19]. Indeed, a bacterial surface-displayed variety of the mutated S-RBD selected for tighter ACE2 affinity presents in the variants of concern N501Y (Alpha, Beta, and Gamma) and E484K (Beta and Delta), and multiple other variants, such as S477N [20]. The variant harbors two mutations within the RBD region (L452R and E484Q), the region responsible for the viral entry [21]. Furthermore, the Q498R exchange has shown a major contribution to ACE2 binding (in epistasis with N501Y). A modest change in affinity might cause a significant increase in the infection rate. As mentioned above, Omicron is found to harbor a common feature of N501Y mutation in the S-RBD protein that is not present in Delta, indicating that this new variant of Omicron provides different affinity enhancing mutations that may arise [20, 22], and is also more transmissible or severe than the Delta variant. The side chain of the Q498R mutation forms two additional hydrogen bonds with Q42 and D38 from ACE2, and the N501Y mutation forms extensive packing interactions with the ACE2 residues Y41 and K353 [23]. Compared to Omicron, the Delta variant carries a characteristic L452R mutation, which is the strong binding interaction of the S-RBD protein described for ACE2 binding [24]; however, this mutation would not clash with the Y102 epitope from the heavy chain of the antibody, leading to its loss of potency against the Delta variant [25]. Both the Delta and Omicron variants share a common T478K exchange, while a low affinity for neutralizing antibodies at this position was observed for Delta [26, 27]. Thus, further surveillance, diagnosis, evaluation, and treatment of the mutated SARS-CoV-2 strains are necessary.

Based on the field of computational drug discovery of protein-protein and protein-ligand interactions, we proposed detailed structural behaviors at the potential binding interface that would allow us to gain a better understanding of the mutated molecular effects on the spike protein of SARS-CoV-2 with its cellular receptor [28–33]. The main purpose of the present work was to compare the structural information and energetic affinity between the Wuhan Hu 1 wild-type and the mutated (Kappa, Delta, and Omicron) variants in complex with ACE2 protein by using a variety of computational techniques. Despite the computational demands, the binding free energy approach, such as MM-PBSA method, is one of the most successful and precise in silico techniques for accurate prediction of the ligand selectivity [34], protein-protein interactions [35], and protein stability [36]. It has been shown that the major changes in the binding affinity between the S-RBD and ACE2 protein occurred at various residue positions that almost mutated in the S-RBD region. Moreover, most of the remaining residue exchanges at the binding interface did not generally affect binding with the receptor but may instead change epitopes for the neutralizing antibody response [37, 38].

Materials and methods

Mutation variants and structure preparation

The FASTA sequence of the SARS-CoV-2 of Wuhan-Hu-1 (wild-type) was obtained from Uniport 13. (GenBank: QHD43415.1) [39]. The Kappa (GenBank: MZ157006.1), Delta (GenBank: QWK65230.1) and Omicron (R40B60 BHP 3321001247/2021) variants were obtained from GSAID database. The translated sequence was performed to collect the omicron spike protein. Clustal-X2, a bioinformatics program, was used to align the wild-type sequence with variants of delta and omicron sequences. The box shade application in each sequence alignment of the SARS-CoV-2 S-RBD variants: Kappa (B.1.617.1), Delta (B.1.617.2), and Omicron (B.1.1.529), in comparison with the wild-type (WT) was showed in Fig 1. The starting structure for the wild-type SARS-CoV-2 S-RBD protein in complex with ACE2 were obtained from the RCSB Protein Data Bank (PDB ID: 6M0J) [40]. We construct a various type of the mutated S-RBD protein according to Fig 1 by direct replacing in silico mutated amino acids of each mutation in the S-RBD region (438–505, red color) of the wild-type S-RBD in an experimental PDB structure before proceeding to the next MD simulations step.

thumbnail
Download:
Fig 1. Amino acid sequence alignment in each SARS-CoV-2 S-RBD variant; wild-type, Kappa (B.1.617.1), Delta (B.1.617.2), and Omicron (B.1.1.529).

The ACE2 receptor binding domain (RBD, or spike RBD region) are located in residues 438–505.

https://doi.org/10.1371/journal.pone.0277745.g001

Molecular dynamics simulation and binding free energy analysis

All molecular dynamics (MD) simulations were performed by PMEMD.CUDA [41, 42] from AMBER 18 suite of programs [43] on NVIDIA Geforce GTX-1070 Ti graphic card for speeding up the simulation times. Each mutation complex under periodic boundary condition was solvated in a cubic box of TIP3P water molecules extending to 10 Å along each direction from the complex model, and Na+ ions were added as neutralizing counterions. The cutoff distance was kept to 12 Å in order to compute the non-bonded interactions. The AMBER ff14SB force field parameters were used to apply the description of the complex characterization. The long-range electrostatic were treated using the particle mesh Ewald (PME) method [44]. The SHAKE algorithm and Langevin dynamics were applied to constrain the bonds that involved hydrogen atoms and to control the temperature.

The mutated three-dimensional models were first minimized using single point energy calculation to relax and remove overlapping atoms. The first step was to allow, out of 10,000 iterations, only water molecules to move. In the second step, each of the 10,000 iterations, hydrogen and protein side chains were relaxed, in a fixed order. Finally, 20,000 iterations were calculated with the restriction-free system. After that, the optimized MD systems were carried out gradually heating (H) phase NVT (constant number of atoms, volume, and temperature) ensemble with the fixed protein atoms over 100 ps from 0 K to 310.15 K by using a force constant of 10 kcal mol-1 Å-2. This was followed by 1,000 ps of equilibration-1 (Eq 1) phase NVT-MD at 310.15 K at a force constant of 5.0 kcal mol-1 Å-2. Then, 10,000 ps without any constrained forces of equilibration-2 (Eq 2) phase NPT (NPT; constant number of atoms, pressure, and temperature) were performed for each fully flexible equilibrium system at the same temperature and 1 atm pressure. The density of each system was about 1.0 g cm3. At the end of the production period, all memory the initial configuration was intentionally lost. This was done so that the final results would not depend on the initial configuration. Finally, 50 ns of unrestrained production (Prod) phase MD simulation were applied in the NVT ensemble at a constant temperature of 310.15 K. The time step of 2 fs was set and the trajectory was recorded every 0.2 ps. Trajectory analyses (root mean square deviation and fluctuation, dynamic cross-correlation, hydrogen bond) were carried out using CPPTRAJ module [45] from Amber 18 program from all 50 ns trajectories data. The structural images were presented using DS software.

Amber molecular mechanics Poisson–Boltzmann surface area (MM-PBSA) [46] approach was used to calculate the binding free energy of S-RBD to ACE2 using the snapshots extracted from MD simulation trajectories. The 2500 numbers of structural frequencies were used to extract the structures from last 50 ns trajectories data of the Prod phase MD simulations. The 1000 snapshots were collected from the trajectory data to calculate the binding free energy by the MM-PBSA procedure with the grid size of 0.5 Å. The binding free energy (ΔGbinding, kcal mol-1) was calculated for each molecular species (each ACE2-bound S-RBD variant complex) in the energetic framework of the MM-PBSA method, and the ΔGbinding was computed using Eq 1: (1) in which ΔEvdW and ΔEEEL refer to the van der Waals and electrostatic contribution as calculated by the molecular mechanics force field, and ΔGSOL represent the solvation free energy. The internal dielectric constant was set to the values of 2, 3 and 9 for nonpolar, polar, and charged residues, respectively.

Computational alanine-scanning mutagenesis analysis.

For alanine-scanning mutagenesis, the role of the protein-protein interface key residues was studied by performing computational alanine scanning calculation [47]. Each mutant structure was prepared by Amber 18 software. The mutated amino acid of S-RBD mutant system was directly replaced in the wild-type model of SARS-CoV-2 S-RBD protein (pdb file). Here, this mutation involves all side chain atoms except the beta-carbon (CB) were removed and also their corresponding information (name, number, and coordinates) from the pdb files. Finally, change the residue name to “ALA” for all remaining atoms of the mutant residue. Other mutations will be able to follow similar procedures where the group of atoms after CB but before the carbonyl-corbon (C) may be removed from the pdb file. It can be noted that only one mutation can be measured during a single calculation. Then, the corresponding topology files in each mutant structure was built and calculated the binding free energy by MM-PBSA method. The ΔGbinding of each mutant structural system, in which each key residue was replaced by alanine (A) by truncating the mutated residue at the Cγ-atom, was calculated with the MMPBSA method. The difference in the binding free energy between the wild-type (ΔGwild-type) protein and its alanine mutant (ΔGalanine) counterpart, ΔΔGbinding (kcal mol-1), is given by Eq 2: (2)

A negative ΔΔG values indicate unfavorable substitution for alanine in the relevant position, meanwhile, indicating a favorable contribution for the wild-type residue in that position.

Hydrogen bond analysis.

To search for hydrogen bonds formation, the number of hydrogen bonds present at each frame will be determined, hydrogen bond donors and acceptors of heavy atoms, will be counted using hbond command in CPPTRAJ module which were printed with the avgout format [48]. The distance cutoff for hydrogen bonds (acceptor to donor heavy atom) was set at 3.5 Å, while angel cutoff for hydrogen bonds was performed at 135°. Hydrogen bond occupancy for the residue pairs which are calculated by averaging over production phase NPT-MD simulations were written in term of Eq 3: (3) where Frames is the number of each coordinate structure that the observed hydrogen bond is present.

Dynamic cross-correlation matrix analysis.

Dynamic movement and correlations between the Cα–atoms positions of SARS-CoV-2 spike RBD and ACE2 protein over the simulation period were calculated by dynamic cross-correlation matrix (DCCM) analysis. DCCM was performed using CPPTRAJ module of the AMBER 18 suite. DCCM diagrams are displayed as a color-coded matrix of Pearson correlation coefficients. Where there is a highly-correlated motion, the movement towards the same direction between the residue pairs shows a positive value in the color range from light green to deep red (+1); while in anti-correlated motions, the movement shows a negative value in the color range from gray to royal blue (-1) [49]. The diagonal square relates to the correlation of a residue with itself, i.e., only a region remarked to have highly-positive values (red). The potential amino acids in binding site of S-RBD protein were clarified at residue number of 438–505.

Principal component analysis.

To evaluate the displacement of atoms and conformational dynamics of a protein complex, principal component analysis (PCA) was performed and analyzed using a covariance-matrix-based approach [50, 51]. The elements of the positional covariance matrix C were obtained based on the following Eq 4: (4) where xi and xj are the instant coordinates of the ith and jth Cα-atoms of the systems for using in building matrix C, while 〈xi〉 and 〈xj〉 refers to an ensemble average. The averaged values are computed over the Prod phase MD simulations after superimposition on a reference structure using the CPPTRAJ module of the AMBER 18 package, solvent water molecules and neutralizing ions added by the Leap module are stripped prior to MD trajectory generation. PCA was performed for Cα-atoms on 10,000 snapshots each. PC1 and PC2, which represent the first two principal components, are created from the trajectories averaged from the wild-type and the variant system. The trajectories were analyzed for the relative motions about their center of masses.

Results and discussion

Multiple alignment and physical parameters of spike-RBD variant proteins

Many mutations in the S-RBD region were found in the Omicron variant compared to the others according to the multiple sequence alignment, as illustrated in Fig 1. This indicates that Omicron may be immunologically resistant to the neutralizing response. The mutated spike protein of Omicron includes more than 30 mutations, half of which are located in the RBD regions. The S-RBD loop at the K440-Y505 positions was observed as a viral determinant for ACE2 binding [52, 53]. Furthermore, the residue of T478 is a common mutation seen in both the Delta and Omicron variants.

Initially, the primary structure analysis shared by all SARS-CoV-2 variants reveals that there is an increase in the following amino acid compositions in the S-RBD of the Omicron variant compared to the Delta (Table 1): arginine (R), lysine (K), histidine (H), and glutamic acid (E), suggesting Omicron has more amino acids with electrostatically charged side chains that can contribute to the salt bridge interaction of the ACE2 protein [54]. The higher percentage residue composition of phenylalanine (F), leucine (L), and alanine (A) was observed in the whole protein of Omicron compared to other variants. These non-polar residues are located inside the protein backbone core and are inaccessible to the environmental solvent [55]. Additionally, the residue composition of Omicron is low in polar uncharged side chains, such as asparagine (N) and glutamine (Q), when compared to the Delta mutation.

thumbnail
Download:
Table 1. Amino acid percentage composition of the mutated SARS-CoV-2 variants with reference to the wild-type in the residue range of 333–526 positions and contact to the RBD protein.

https://doi.org/10.1371/journal.pone.0277745.t001

Mutational scanning and energetic analysis of the wild-type S-RBD residues at the binding interface with ACE2

According to recent X-ray results [40, 16, 56, 57], we found that the simulated structure of the ACE2-bound SARS-CoV-2 S-RBD protein preserves the same binding interface information, which does not induce any conformational change in the receptor-binding site [24, 58–61]. One essential feature at the binding interface hotspot between S-RBD and ACE2 proteins is the number of hydrophilic residue interactions [16, 40], which are conserved in the corresponding MD simulation calculation. Certainly, the hydrogen-bonding hotspots and two salt bridges (K31-E35 and K353-D38) in the present work stably populate S-RBD binding with the ACE2 receptor [58, 62–64]. Here, the energetic contributions of wild-type S-RBD residues were firstly analyzed over the production phase MD simulations for understanding the mapping key population hotspot in the SARS-CoV-2 S-RBD wild-type protein. That is to say, the atomistic structural characterisation of the binding interface residues could better quantify the binding energetic differences between the spike protein RBD variants and its cellular ACE2 receptor, providing fundamental information about the design of structure-based vaccines and/or neutralizing antibodies development.

To provide the qualification of the binding energetic characterization and why they are causing concern in the S-RBD variants, we employed the combined description for all topical protein-protein interface interactions along with the corresponding energetic quantification discussed in detail by an alanine scanning calculation through the MM-PBSA method on the binding free energy of the mutated S-RBD complex with ACE2. Each mutant structure in the pdb file was prepared by Amber 18 software as listed in S1 Table. For instance, mutated A498 residue of the S-RBD mutant system was directly replaced in Q498 of the wild-type model of the SARS-CoV-2 S-RBD protein (pdb file). Other mutations will be able to follow similar procedures in each mutation during a single calculation. Then, the corresponding topology files in each mutant structure were built and the binding free energy calculated by MM-PBSA method for comparison to the wild-type MD simulation as the trajectories to be analyzed. As a result, the significant relative binding free energy (ΔΔGbinding, kcal mol-1) of each mutant was reported, demonstrating the relative effect the mutation has on the ΔGbinding of binding for the complex as shown in Fig 2. It can be noted that the negative values of the ΔΔGbinding obtained from the computational alanine scanning mutagenesis indicate unfavorable substitutions for alanine in the relevant position. Meanwhile, it also indicates a favorable contribution for the wild-type residue in that position and vice versa.

thumbnail
Download:
Fig 2. An alanine scanning of spike-RBD binding region in the wild-type SARS-CoV-2 structure over the production phase MD simulations.

The binding energy changes (ΔΔGbinding, kcal mol-1) are reported that demonstrates the relative affect the mutation in each the ΔG of binding complex, ΔΔGbinding = ΔGwild-type−ΔGalanine. The binding energy for each alanine mutant based on the MM-PBSA method is plotted.

https://doi.org/10.1371/journal.pone.0277745.g002

The strongest binding energies of hotspot residues are located in the central segment of the interface, i.e., N487, Q493, and Q498. In other words, a large and comparable loss of binding energy is observed, which forms favorable hydrophilic contacts with ACE2. N487 contributes to ACE2 binding with two interface-anchoring Q24 and Y83 hydrogen-bonding interactions (Fig 3A). Alanine substitution at the N487 position is accompanied by a ~6.0 kcal mol-1 loss in binding free energy (ΔΔGN487A = –5.92 ± 2.55 kcal mol-1). At the Q493 position, one intermolecular hydrogen bond with K31 and E35 of the ACE2 protein and one intramolecular interaction with the S494 side chain (occupancy 22.35%, distanceavg = 3.082 Å) were found to form at the binding interface (Fig 3B). Likewise, Q498 of the S-RBD protein establishes an extensive network of favorable interactions across the ACE2 binding interface with K353-D38 salt bridges and further stabilizes with the Q42 side chain on the receptor (Fig 3C). A significant loss in the binding free energy is predicted as ΔΔGQ498A = –13.64 ± 2.53 kcal mol-1 by replacing these Q493 and Q498 positions with alanine, which make the Q498 key residue of a viral binding hotspot with respect to the less effective Q493, for which ΔΔGQ493A = –8.77 ± 2.88 kcal mol-1, reflecting the relevant roles played by these two residues at the binding interface.

thumbnail
Download:
Fig 3. Structural details of the protein-protein interaction involving the S-RBD and ACE2 residues at the binding interface obtained from the Prod phase MD simulations.

The secondary structures of S-RBD and ACE2 are portrayed as white ribbons and light-blue ribbons, respectively. Each interacting amino acid residue is highlighted in black-color labeled. The hydrogen bonds and salt bridges are represented in black and red dashed line, respectively.

https://doi.org/10.1371/journal.pone.0277745.g003

Other important binding residues showing significant energy loss upon alanine modifications included the L455, F456, F486, T500, and Y505 positions (with a threshold of –2.5 kcal mol-1). In these positions, the binding free energy losses caused by alanine substitutions can decrease in the energy range of 2.6 ~ 3.0 kcal mol-1 (Fig 2). Furthermore, the side chain of the L455 residue points to a charged pocket sealed with D30, K31, N33, and H34 of the ACE2 binding receptor through a moderate stabilizing vdW interaction, while F456 provides the intermolecular π-cation-stabilizing interaction with T27, D30, and K31 of ACE2 (Fig 3D) [65]. F486 presents to stabilize the hydrophobic side chain of Y83 on ACE2 binding, and also covers the neighboring side chain of L79 and M82 residues (Fig 3E). When all of these residues are replaced with alanine, a loss of the corresponding binding free energy is reflected as ΔΔGL455A = –2.67 ± 2.14 kcal mol-1, ΔΔGF456A = –2.62 ± 1.29 kcal mol-1, and ΔΔGF486A = –2.95 ± 1.31 kcal mol-1.

Similarly, Y505 interacts in the central segment of the binding interface with the side chain of the E37-R393 salt bridge on the receptor, while exchanging an internal π-cation interaction between an aromatic group of this residue and the guanidinium moiety of R403, making this interaction pattern reduce the binding affinity for the ACE2 receptor by ~3.0 kcal mol-1 (ΔΔGY505A = –2.97 ± 2.52 kcal mol-1) (Fig 3F). Besides, R403 provides faintly stabilizing polar or electrostatic interactions (occupancy 5%, distanceavg = 3.228 Å) with the E37 side chain of ACE2 at the binding interface. The relative importance of these key residues agrees well with its interaction formation and clearly supports results on the interface structure and energetic framework in further experiments.

Structural analysis and stability of the wild and mutant variants SARS-CoV-2 spike protein complexes

The crystal structure of the wild-type SARS-CoV-2 spike receptor-binding domain bound with ACE2 (PDB ID: 6M0J) was taken as the reference [40]. The structure of the wild-type complex along with the variants of the Kappa (L452R and E484Q), Delta (L452R and T478K), and Omicron (N440K, G446S, S477N, T478K, E484Q, Q493K, G496S, Q498R, N501Y, and Y505H) models were subjected to the MD simulation calculations. Various structural analyses were performed to study the structural impact of these variants. The final MD structures of each variant complex at the binding interface are illustrated in S1–S4 Figs.

Inspecting the dynamics and complex stability of the ACE2-bound S-RBD variant complexes, Root Mean Square Deviations (RMSD) with respect to their optimized structures were monitored for the simulated trajectories of all phase MD simulations, as shown in Fig 4A. Steady oscillations and small fluctuations of RMSD were observed in each complex model compared to the wild-type model. The RMSD plots for the wild-type complex were found to be the highest as compared with the values for the observed variants (S5 Fig). The wild-type model showed increases in the RMSD value at around 39 ns; however, after 45 ns, the RMSD gradually decreased and showed variable fluctuations until the end of the simulation. The slight conformational shift was seen in the Kappa and Omicron variants between 35–40 ns, which later stabilized [59]. Comparatively, Delta showed stable RMSD with inconsiderable fluctuations and endured lesser conformational changes during the whole MD simulations [24, 59, 60]. The average RMSD undulation amplitude was below 0.5 Å throughout most of the simulation time.

thumbnail
Download:
Fig 4. The molecular dynamics results of the ACE2-bound S-RBD complex with the wild-type, Kappa, Delta, and Omicron variants over 50 ns simulation times.

(a) RMSD distribution, and (b) RMSF of Cα-atom of spike RBD; wild-type, Kappa, Delta, and Omicron, in complex with the ACE2. The positions of spike RBD-ACE2 binding region are denoted the fragment of residues 438–505.

https://doi.org/10.1371/journal.pone.0277745.g004

Furthermore, to explore the structural analysis of the potential model in sufficient detail, we calculated the dynamics of the residues in terms of the analysis of root mean square fluctuations (RMSF), which is useful to situate the flexible and disordered regions as well as the heterogeneity of the system [66, 67]. The plot showed a similar pattern of fluctuations with different magnitudes for all simulated systems (Fig 4B). The higher degree of evaluation was found in the native complex, while the Kappa variant witnessed the lowest level of fluctuations within the amino acid residues, indicating the better stability of the Kappa variant complex, followed by the Delta and Omicron variants [24, 59]. As the binding interface residues array in a random coil conformation lacks structural rigidity, a structural arrangement in this region should be necessary to sustain their configurations of the binding surface, which may facilitate the binding affinity. Three loop regions (residues 474–486, 488–490, and 494–505) of SARS-CoV-2 demonstrated limited fluctuation compared to the corresponding region in the SARS-CoV wild-type [60]. Lower fluctuations, observed in the case of SARS-CoV-2 that bind to the K353 of ACE2, could be another reason for the better affinity of the S-RBD protein to ACE2 [40, 68]. The residues 480–488 (a random coil near the binding interface), showed are markable decrease in RMSF value for the variant models (S6 Fig). That is to say, these strong affinity positions stabilized the complex as a whole, indicating the reduced fluctuation in these interfacial residues for the ACE2 binding. In addition to terminal loops, the β-strand conformation (residues 510–524) displayed greater structural rigidity with less fluctuation than the wild-type, which mostly consisted of loops and might account for the binding affinity in each variant [58].

Effect of SARS-CoV-2 variants on the hydrogen bond formation with ACE2

To evaluate the interfacial binding behavior on the ACE2-bound S-RBD variant complexes, the simulated trajectories were analyzed with the cut-off of 3.5 Å distances from acceptor to donor atomic pairs that play critical roles in the stronger binding tendency towards the ACE2 receptor compared to the wild-type. Here, three mutagenic viruses in the S-RBD protein of SARS-CoV-2, currently designated as VOC by WHO, i.e., B.1.617.1 (Kappa), B.1.617.2 (Delta) and B.1.1.529 (Omicron), were simulated. The evolution and stability of the hydrogen bond formations between the binding interface residues were monitored as a function of MD simulation trajectory to trace the occupancy and hydrogen bond formation, in which more than a 10% occupancy rate was shown, as listed in Table 2. Recently, the relative strength and hydrogen bond variation between each atomic pair have been evaluated computationally with the findings substantially in line with our measurements [21, 24, 59, 63, 69]. That is to say, the high-affinity protein-protein complex of S-RBD toward ACE2 is attributed to the explanation of why the binding affinities are distinctly different from each other.

thumbnail
Download:
Table 2. Hydrogen bond occupancy for the residue pairs at the S-RBD binding interface with ACE2 (cut-off 3.5 Å) over the production phase MD simulations.

https://doi.org/10.1371/journal.pone.0277745.t002

Based on the analysis of hydrogen bonding at the S-RBD binding interface with ACE2 (cut-off 3.5 Å), we observed that the mutant residue was found to have a cascading effect on the interfacial surface of the S-RBD region. A similar pattern of an amino acid pairing on the S-RBD regions, i.e., Y449, N487, Q493, Q498, T500, and Y505, was monitored as key hotspots for all complex systems (bold in Table 2), which are oriented in the potential site through negligible hydrogen bonding with two-salt bridges at the ACE2 binding [63, 64, 70–72]. In addition, the high average number of hydrogen bonds at the interfacial binding surface of the Kappa and Delta variants can maintain the stronger interaction of the ACE2-bound S-RBD complex, which is relatively higher compared with the wild-type [24]. Apart from these residues, the Kappa variant has two additional residues (S477 and T478), the Delta variant has five additional residues (Y473, Y489, F490, Y495, and N501), and the Omicron variant showed more two additional residues (Y453 and S496), that were found to be involved in forming the intermolecular hydrogen-bonded interactions [24, 73]. Moreover, the Q493K mutation in the Omicron RBD protein resulted in the loss of the salt bridge with K31 and the hydrogen bond with Q42 of ACE2 [59], while the contact with D38 was newly formed. The likelihood of an interaction between the Q493 of SRAS-CoV-2 and K31 residues of ACE2 was reported recently [60, 68]. In the case of the G446S mutant, no interaction was observed during the simulation as supported by the previous finding. The Q493 position in all complex systems interacts with E35 of the ACE2 receptor with significant occupancy. T500 of the wild-type, Kappa, and Omicron variants mediates a hydrogen bond with D355 of the ACE2 with high occupancy, while the interaction between T500 of Delta and Y41 of ACE2 was observed. Our finding was consistent with previous studies [59, 74, 75]. However, some contradictorty observations were also found. Notably, residues L455, F456, F486, and S494, which were reported to form hydrogen bonds and electrostatic interactions leading to an enhanced binding affinity of SARS-CoV-2, were not observed in these simulations [39, 56, 68].

The structural environment changes in the Q493 beta-sheet can be observed according to the L452R mutation in the Delta variant (Fig 5A) [73]. We found that the number of hydrogen bond occupancies of Q493 with ACE2 at the binding interface increased significantly for the variants that express the L452R mutation (Fig 5B), while the Kappa variant changed slightly. Although the Omicron variant did not mutate at the L452, we found an increased amount of intermolecular hydrogen bonding of Q493 for ACE2 binding, which might be due to this variant carrying the Q493K mutation. That is to say, the changing of the structural environment neighbor 493 position and/or by itself can enhance binding affinity through hydrogen bond formation with the ACE2 protein relative to the Delta. The basic side chain of mutated R452 orients itself for stabilizing and increases the hydrogen bond interaction with the salt-bridge residues (K31 and E35) by Q493 [76, 77]. The oxygen atom of the Q493 side chain creates hydrogen-bonding interaction with the K31 side chain, while one of the hydrogens of the amine group forms the intermolecular hydrogen bonds with the carboxyl group of E35 from ACE2 (Fig 5C). In addition, we also mentioned that the Q493K mutation in the Omicron RBD protein presented several intermolecular hydrogen bonds with E35, and the contact with D38 was newly formed, which was not observed in the wild-type residue expressing a glutamine 493 (Fig 5D and 5E, and S7 Fig). However, the stable salt-bridge with K31 was not observed, which might be too short of a lysine side chain for ACE2 binding [58]. Thus, the Omicron variant can form stable salt-bridges between Q493 and E35 and D38 of ACE2 interaction.

thumbnail
Download:
Fig 5. Hydrogen bonding evolution in the Delta and Omicron variants.

(a) Superimposition between S-RBD of the Delta and the wild-type that represents the Q493 (blue), L452 in the wild-type (yellow), and R452 in the Delta (green) as expressed on two neighboring β-strands forming a small antiparallel β-sheet (red). (b) Combined analysis of the wild-type and three variants that carry on the L452R exchange with regarding to hydrogen bond occupancies (%) formed by Q493 with ACE2 (* p < 0.05, two-tailed Student’s T-test). (c) Hydrogen bond formation represents between Q493 in S-RBD of the Delta and K31 and E35 salt bridges at binding interface with ACE2. (d) The number of the intermolecular hydrogen bond between Q493 (or K493 in the Omicron) from S-RBD and ACE2 at K31, E35, and D38 (* p < 0.05, ** p < 0.005, two-tailed Student’s T-test). (e) Structural representation of the two newly formed salt bridges between K493 on the S-RBD in the Omicron and E35 and D38 that expressed on ACE2. The hydrogen bonds are indicated by black dashed line.

https://doi.org/10.1371/journal.pone.0277745.g005

Energetic affinity of the binding interfaces with ACE2

Supporting the previously stated results, we also calculated the values of binding free energy (ΔGbinding, kcal mol-1) obtained from the MM-PBSA method to compare the binding interface behavior of the mutated S-RBD protein relative to the wild-type [24, 58, 64, 70–74, 78]. The estimated binding free energies (ΔGbinding) calculated for the ACE2-bound S-RBD complex for each variant are listed in Table 3. It was interesting to note that the heavily mutated Delta protein (Delta, Omicron, and Kappa variants) showed the strongest binding affinity of –73.46 kcal mol-1 followed by the Omicron (–61.75 kcal mol-1) and Kappa (–58.90 kcal mol-1) variants. Although, our results are a little different from the previous literature that proposed the highest binding affinity for the Omicron variant [59], the ΔGbinding magnitude of the mutated S-RBD proteins was also found to be lower than that of the wild-type complex by –5.79 - –20.35 kcal mol-1. Furthermore, the calculated ΔGbinding value is in a good agreement with the findings of previous work that showed greater cell-surface binding affinity between the mutated S-RBD (N510Y, K417N, and E484K) and the ACE2 protein than the wild-type by a competition binding assay [73]. The enhanced interaction of the Omicron S-RBD protein with ACE2 is consistent with previously published data [79], and may contribute to the increased infectivity of the Omicron variant. Our result again emphasized the role of the multimutated S-RBD protein in increasing binding affinity toward the ACE2 receptor as supported by the findings in previous studies [58, 73, 74, 78].

thumbnail
Download:
Table 3. Free energy terms (kcal mol–1) for SARS-CoV-2 variants on the S-RBD protein binding to ACE2 estimated by MM-PBSA method.

https://doi.org/10.1371/journal.pone.0277745.t003

Furthermore, the individual energy term (vdW, polar, and nonpolar solvation free energies) was also calculated from the total binding free energy to understand the contribution of individual energy components in the binding. The results revealed that the main driving force for increasing the binding affinity was compensated by favorable change in the electrostatic contributions (EEL) and van der Waals (vdW) interactions. This energy term was found to relate to the analysis of the electrostatic interaction energy at the 493 key residue that showed a strong increase for Omicron when compared to others (S8 Fig), supporting the structural environment changes of Q493 as the result of L452R mutation. The sum of EEL and the polar component of solvation (EPS) contributions is more positive for the mutations, indicating that stabilization is mostly due to the hydrophobic interactions. Taken as a whole, the structural analysis could show the binding characteristics for the ACE2-bound S-RBD variant interaction compared to the wild-type SARS-CoV-2. Specifically, the high mutation rate in the Omicron S-RBD corresponding with the electrically charged side-chains leads to an improved binding affinity for the ACE2 interaction compared to the Delta variant, which may allow the Omicron variant to facilitate higher transmissibility in contrast to the wild-type.

With the help of the structural information provided by the analysis of the binding affinity, the binding free energies were decomposed into the energy contributions of the individual residues at the binding interface, which hopefully will allow researchers to detail the cross-reactivity of the neutralizing antibody response. Fig 6 illustrates the decomposed per-residue free energy (kcal mol-1) calculated by the MM-PBSA method upon the interfacial residues of the S-RBD regions (at position 438–505) with ACE2 receptor (S8 Fig). The negative and positive energy values represent favorable and unfavorable contributions, respectively.

thumbnail
Download:
Fig 6. Per-residue free energy decomposition of the key hotspot at the S-RBD domain (residues 438–505) of the binding interface with ACE2 protein.

All values were given in kcal mol-1.

https://doi.org/10.1371/journal.pone.0277745.g006

In addressing the highest mutagenic rate in the Omicron variant, the dominant mutational Q493K and N501Y hotspots highlighted more favorable interactions with strong binding affinity for the ACE2 interaction. Consequently, we believe that N501Y is a critical mutation that affects the transmission of SARS-CoV-2 by strengthening the interaction between S-RBD and the ACE2 protein [80]. The strong interaction of the N501Y mutant leads to the tighter binding of SARS-CoV-2 to the host cell, allowing complete membrane fusion or the internalization of the receptor with the virus. Several other S-RBD variants, i.e., Kappa and Delta that have very high transmission rates even though they do not contain the N501Y mutation, might increase the binding affinity for the ACE2 receptor, leading to similar effects [81]. In addition, the other mutated residues (N440K, S477N, T478K, and Y505H) were found to significantly promote favorable interaction to the binding of the Omicron S-RBD complex with ACE2 protein [59], while the mutated G446S and E484Q positions contributed less to the interaction increment than those of the mutated residues in the Omicron variant, which is related to the previous experimental and theoretical analysis [73]. On the contrary, the decomposed energy forms unfavorable contact with ACE2 binding even though the mutual residues of G496S and Q498R slightly shift position. This interesting result implied the impact of these mutagenic positions for higher infectivity and transmission of the Omicron variant.

Similarly, L452 and T478, which are mutated to L452R in Kappa, and T478K in the Delta variant, showed high potential individual decomposed energy by mediating a consistent hydrogen bond [59]. R454, L455, and F456 showed significantly less contribution difference in the binding energy for all systems. The Delta mutant also specifically promoted the dominant positions, N487, Y489, F490, and L492, with stronger interaction at the binding interface of the complex compared to the others. Although these positions together with F486, N487, and Y505 presented the lowest binding free energy values in agreement with the hydrogen bonding analysis (Table 2), they did not present mutagenicity in the Delta and Kappa variants.

Previous studies have proposed that different neutralizing antibodies target different regions on the S-RBD protein [24, 82]. The epitome mapping analysis on the monoclonal antibody specific to the spike protein revealed that 20% of more than 400 epitopes were derived from the RBD region. The mutated S-RBD protein in the Kappa, Delta, and Omicron variants of SARS-CoV-2, which contains the ACE2 binding site as well as the epitopes over 80% of antibodies induced by infection or vaccination, decrease the binding affinity of the antibody [23]. For instance, the L452R mutation in the Delta and Kappa variants was found to be the center position of the JMB2002-binding epitope with Y102 from the heavy chain of the Fab, thus providing an explanation for its loss of potency against the mutation [25]. It should be noted that the structural inspection alone may not be sufficient for identifying the key contributions to binding affinity, where the effects of the solvation term are considered. These coevolving residues can mediate protein recognition in multiprotein complexes and are often spatially close to each other, forming clusters of interacting residues that are located near functionally important sites [83, 84].

Dynamic cross-correlation matrix and conformational flexibility analysis

To further investigate the conformational dynamics in detail, the dynamic cross-correlation matrix (DCCM) was calculated to measure the occurrence of correlated motions according to the simulation period based on the positions of backbone carbon–atoms of the S-RBD protein. The phenomenon of the dynamic motions during the simulations of the Cα–atoms of the SARS-CoV-2 spike RBD protein in each bound system are shown in Fig 7. Highly correlated motions, also mentioned as positive correlation, range from the color green (low), through yellow (medium), to deep red (high, +1), while anti-correlated motions, also represented as negative correlation, range from gray (low), through cyan (medium), to royal blue (high, –1). The diagonal square corresponds to the relationship of a residue with itself, i.e., the only region observed to present highly positive values (deep red).

thumbnail
Download:
Fig 7.

Dynamic cross-correlation diagrams of the fluctuations of Cα–atoms of SARS-CoV-2 spike RBD protein in the various variants of (a) Kappa, (b) Delta, (c) Omicron in comparison to the (d) wild-type RBD protein. Positive and negative values are represented in range of color deep red to royal blue, respectively. The diagonal square relates to the correlation of a residue with itself, i.e., only a region remarked to have highly-positive values (red). Black arrows point to the mutation position of each variant that observed for ACE2 binding. The ACE2 binding regions of the SARS-CoV-2 spike RBD protein, at residues number of 438–505, are demarcated with dashed lines.

https://doi.org/10.1371/journal.pone.0277745.g007

The DCCM map results revealed that the S-RBD variants can affect the structural remodeling for the ACE2-bound S-RBD complex, as illustrated by the change in the dynamic patterns and correlations [51, 85]. As shown in dashed lines in Fig 7, three S-RBD variants dominated primarily by correlated motion (lesser gray region) in comparison to the wild-type system. Overall, the Kappa and Omicron conformations display strong correlated motions (red and yellow), while the Delta conformations are embodied by low correlated motion (green). The black arrows point to the mutation position of each variant that was observed for ACE2 binding. The biggest differences observed from the correlated motion maps are that the Omicron variant triggered a correlated dynamic on the binding site as opposed to the Delta and Kappa variants. Since the terminal loops at residues 510–524 accounted for the ACE2 binding based on the structural rigidity (Fig 4B), a more widely correlated motion is observed in the S-RBD variant protein. The above regions relating with obvious changes of motion mode agree with the RMSF changes of the S-RBD variant proteins. These changes of the internal dynamics may reflect different alternation of relative positions between key residues, which may play an important role in the conformational diversity of the S-RBD variants. Furthermore, this finding also explains the lower decomposition binding free energy for the binding site observed during the MM-PBSA calculation and how they were compensated by the additional contacts from the mutation region.

Considering the fact that the biological function of a protein-protein interaction is influenced by its conformational dynamics, we also employed Principal Component Analysis (PCA) as the computational method to investigate the conformational transitions of each complex system using the diagonalization of the covariance matrix of the Cα–atoms fluctuations over the Prod-phase MD simulation. Fig 8 shows the PCA scatter plot generated for the ACE2-bound SARS-CoV-2 spike RBD wild-type and their spike RBD variant complexes along the first two principal components (PC). Most of the simulation ensembles in the variant complex related with the strong binding affinity are concentrated to a narrow range of the conformational space compared to the wild-type [86–88]. This suggests that the variant complex offers more stability and compact packing than the wild-type model, which corroborates the previous reports on the compact cluster of stable states in the variants [24]. Overall, the less correlated motion of the wild-type system confers with the observed higher residue flexibility in Fig 6D, implying that the mutational behaviors of the S-RBD variants triggered conformational dynamics as conferred by the conformational flexibility and dynamic cross-correlation analysis. The variation in structural motions allowed us to assess the binding phenomena of each S-RBD variant for ACE2 binding.

thumbnail
Download:
Fig 8. PCA projection of the motion of Cα–atoms of the RBD-bound S-RBD complex.

The PCA profiles were obtained by plotting the first two principal components (PC1 and PC2) in the various conformation of the Kappa (blue), Delta (green), Omicron (red) and the wild-type (black) systems. PC1 and PC2, respectively, represent a covariance matrix after elimination of eigenvectors. Each point between the single-directional motions represents a unique conformation during the simulation, therefore, similar structural conformations overlap in the graph.

https://doi.org/10.1371/journal.pone.0277745.g008

Conclusions

In the present work, the S-RBD variants in both Omicron and Delta were investigated by focusing on the effectiveness of the structural analysis using several computational tools and a computational saturation mutagenesis structure, as well as the dynamic changes that are likely to contribute to both the protein stabilization and binding affinity of the S-RBD protein for ACE2 interaction. Based on the analysis of the amino acid composition, the amino acids with electrically charged side chains (arginine (R), lysine (K), histidine (H), and glutamic acid (E)) were increased in Omicron, which may increase the contribution of the salt bridge of ACE2 at the binding interface. Consequently, we performed systematic alanine scanning on all different residues in the S-RBD region of the wild-type SARS-CoV-2 structures that form most of the protein binding interfaces, and evaluated the variation in the corresponding binding free energy. These mutagenesis studies provide a clearer picture of the conserved molecular hotspots at the binding interface in the S-RBD protein of SARS-CoV‑2 (wild-type) and highlight residues Q498, Q493, and N487 on the spike protein of the SARS-CoV-2 receptor-binding domain as potential residues contributing to the shaping and determination of stability in the relevant protein-protein binding interface.

Another challenge in this study was providing structural-based information and an energetic framework at the binding interface of the spike-ACE2 complexes. A detailed analysis of hydrogen bonding formation, the estimated binding affinity of the binding complex, and the per-residues decomposition free energies indicated that the S-RBD variants (Delta and Omicron) interacted closely with the ACE2 protein at the binding interface. Nowadays, most vaccines are designed on the basis of the S-RBD protein, leading to the application of the neutralizing antibodies for targeting the S-RBD to weaken its ACE2 binding. We found that the potential hotspots in the S-RBD mutation had the strongest binding affinity compared to the wild-type for ACE2. Although the Omicron S-RBD disrupted some favourable hotspot residues for ACE2 binding, such variants can promote new advisable binding interactions. Additionally, some hydrogen bond formation and π-stacking interactions are formed, which can further encourage the S-RBD protein for ACE2 binding. Additional mutations within the receptor-binding site (RBD region) may affect changes in the amino acid sequence of the epitope and may contribute to the escape from the neutralizing antibody binding. Optimistically, the information and results provided by this work will lead to a concrete informative article for viral community scientists occupied in the context of structure-based design and therapeutic neutralizing antibodies as well as vaccine development, since this information is urgently required in the current fight against future pandemics.

Supporting information

S1 Fig. Structural mapping at the binding interface structure of the wild-type S-RBD protein complex with ACE2.

The structure of ACE2 is shown in gray ribbon, while S-RBD is in green ribbon. The close-up of the binding interfaces is shown the crucial hotspots responsible for binding interaction (right panel). The interfacial residues of SARS-CoV-2 S-RBD residues are annotated and shown in green sticks, while ACE2 residues are in gray sticks. The figure was drawn by Discovery Studio 2019 Client (Biovia software).

https://doi.org/10.1371/journal.pone.0277745.s001

(TIF)

S2 Fig. Structural mapping at the binding interface structure of the Kappa S-RBD protein complex with ACE2.

The structure of ACE2 is shown in gray ribbon, while S-RBD is in green ribbon. The close-up of the binding interfaces is shown the crucial hotspots responsible for binding interaction (right panel). The interfacial residues of SARS-CoV-2 S-RBD residues are annotated and shown in green sticks, while ACE2 residues are in gray sticks. The figure was drawn by Discovery Studio 2019 Client (Biovia software).

https://doi.org/10.1371/journal.pone.0277745.s002

(TIF)

S3 Fig. Structural mapping at the binding interface structure of the Delta S-RBD protein complex with ACE2.

The structure of ACE2 is shown in gray ribbon, while S-RBD is in green ribbon. The close-up of the binding interfaces is shown the crucial hotspots responsible for binding interaction (right panel). The interfacial residues of SARS-CoV-2 S-RBD residues are annotated and shown in green sticks, while ACE2 residues are in gray sticks. The figure was drawn by Discovery Studio 2019 Client (Biovia software).

https://doi.org/10.1371/journal.pone.0277745.s003

(TIF)

S4 Fig. Structural mapping at the binding interface structure of the Omicron S-RBD protein complex with ACE2.

The structure of ACE2 is shown in gray ribbon, while S-RBD is in green ribbon. The close-up of the binding interfaces is shown the crucial hotspots responsible for binding interaction (right panel). The interfacial residues of SARS-CoV-2 S-RBD residues are annotated and shown in green sticks, while ACE2 residues are in gray sticks. The figure was drawn by Discovery Studio 2019 Client (Biovia software).

https://doi.org/10.1371/journal.pone.0277745.s004

(TIF)

S5 Fig. RMSD profiles of MD simulations.

RMSD plot of each S-RBD variant complex with ACE2 protein over the production (Prod) phase MD simulations.

https://doi.org/10.1371/journal.pone.0277745.s005

(TIF)

S6 Fig. Superimposed structure of wild-type and mutated S-RBD protein.

The superimposition structure of the Cα-atom trace of four different structures between the wild-type and mutated S-RBD protein of SARS-CoV-2. The proteins superimpose almost exactly in most parts of the protein expect few regions which exhibit conformation variability (highlighted in the close-up), especially, at the RBD-ACE2 binding interfaces.

https://doi.org/10.1371/journal.pone.0277745.s006

(TIF)

S7 Fig. Hydrogen bonding evolution in the Delta and Omicron variants.

The intermolecular hydrogen bond occupancy between Q493 (or K493 in the Omicron) from S-RBD protein and ACE2 at K31, E35, and D38 for binding complex (* p < 0.05, ** p < 0.005, two-tailed Student’s T-test).

https://doi.org/10.1371/journal.pone.0277745.s007

(TIF)

S8 Fig. Energy decomposition analysis of wild-type and mutated S-RBD protein in complex with ACE2.

Per-residue decomposed energy of (a) electrostatic, and (b) vdW interaction on the key hotspots of the S-RBD regions (at position 438–505) for ACE2 binding. All values were given in kcal mol-1.

https://doi.org/10.1371/journal.pone.0277745.s008

(TIF)

S1 Table. Computational alanine scanning mutation analysis.

Relative binding free energy terms (kcal mol–1) calculated by the computational alanine scanning mutagenesis using MM-PBSA method approach for the S-RBD regions (at position 438–505) of SARS-CoV-2 residues effectively involved in the binding interface with the ACE2.

https://doi.org/10.1371/journal.pone.0277745.s009

(DOCX)

Acknowledgments

Authors would like to acknowledge the Erawan HPC Project, Information Technology Service Center (ITSC), Chiang Mai University, and Thailand Excellence Center for Tissue Engineering and Stem Cells (PK), Department of Biochemistry, Faculty of Medicine, Chiang Mai University, Chiang Mai, Thailand for supporting and providing access to their computing resources.

References

  1. 1. Supasa P, Zhou D, Dejnirattisai W, Liu C, Mentzer AJ, Ginn HM, et al. Reduced neutralization of SARS-CoV-2 B.1.1.7 variant by convalescent and vaccine sera. Cell. 2021;184(8):2201–11.e7. pmid:33743891
  2. 2. Liu C, Ginn HM, Dejnirattisai W, Supasa P, Wang B, Tuekprakhon A, et al. Reduced neutralization of SARS-CoV-2 B.1.617 by vaccine and convalescent serum. Cell. 2021;184(16):4220–36.e13. pmid:34242578.
  3. 3. Zhou D, Dejnirattisai W, Supasa P, Liu C, Mentzer AJ, Ginn HM, et al. Evidence of escape of SARS-CoV-2 variant B.1.351 from natural and vaccine-induced sera. Cell. 2021;184(9):2348–61.e6. pmid:33730597.
  4. 4. Dejnirattisai W, Zhou D, Supasa P, Liu C, Mentzer AJ, Ginn HM, et al. Antibody evasion by the P.1 strain of SARS-CoV-2. Cell. 2021;184(11):2939–54.e9. pmid:33852911
  5. 5. Colmenares-Mejía CC, Serrano-Díaz N, Quintero-Lesmes DC, Meneses L, Salazar Acosta I, Idrovo ÁJ, et al. Seroprevalence of SARS-CoV-2 Infection among Occupational Groups from the Bucaramanga Metropolitan Area, Colombia. Int J Environ Res Public Health. 2021;18(8):4172. pmid:33920843.
  6. 6. Viana R, Moyo S, Amoako DG, Tegally H, Scheepers C, Althaus CL, et al. Rapid epidemic expansion of the SARS-CoV-2 Omicron variant in southern Africa. Nature. 2022. pmid:35042229
  7. 7. Mannar D, Saville James W, Zhu X, Srivastava Shanti S, Berezuk Alison M, Tuttle Katharine S, et al. SARS-CoV-2 Omicron variant: Antibody evasion and cryo-EM structure of spike protein–ACE2 complex. Science. 2022;375(6582):760–4. pmid:35050643
  8. 8. McCallum M, Walls Alexandra C, Sprouse Kaitlin R, Bowen John E, Rosen Laura E, Dang Ha V, et al. Molecular basis of immune evasion by the Delta and Kappa SARS-CoV-2 variants. Science. 2021;374(6575):1621–6. pmid:34751595
  9. 9. McCallum M, Bassi J, De Marco A, Chen A, Walls Alexandra C, Di Iulio J, et al. SARS-CoV-2 immune evasion by the B.1.427/B.1.429 variant of concern. Science. 2021;373(6555):648–54. pmid:34210893
  10. 10. McCallum M, Czudnochowski N, Rosen Laura E, Zepeda Samantha K, Bowen John E, Walls Alexandra C, et al. Structural basis of SARS-CoV-2 Omicron immune evasion and receptor engagement. Science. 2022;375(6583):864–8. pmid:35076256
  11. 11. Liu C, Zhou D, Nutalai R, Duyvesteyn HME, Tuekprakhon A, Ginn HM, et al. The antibody response to SARS-CoV-2 Beta underscores the antigenic distance to other variants. Cell Host Microbe. 2022;30(1):53–68.e12. pmid:34921776.
  12. 12. Dudas G, Hong SL, Potter BI, Calvignac-Spencer S, Niatou-Singa FS, Tombolomako TB, et al. Emergence and spread of SARS-CoV-2 lineage B.1.620 with variant of concern-like mutations and deletions. Nature Communications. 2021;12(1):5769. pmid:34599175
  13. 13. Huang Y, Yang C, Xu X-f, Xu W, Liu S-w. Structural and functional properties of SARS-CoV-2 spike protein: potential antivirus drug development for COVID-19. Acta Pharmacologica Sinica. 2020;41(9):1141–9. pmid:32747721
  14. 14. Samrat SK, Tharappel AM, Li Z, Li H. Prospect of SARS-CoV-2 spike protein: Potential role in vaccine and therapeutic development. Virus Research. 2020;288:198141. pmid:32846196
  15. 15. Structure Li F., Function, and Evolution of Coronavirus Spike Proteins. Annual Review of Virology. 2016;3(1):237–61. https://doi.org/10.1146/annurev-virology-110615-042301.
  16. 16. Shang J, Ye G, Shi K, Wan Y, Luo C, Aihara H, et al. Structural basis of receptor recognition by SARS-CoV-2. Nature. 2020;581(7807):221–4. pmid:32225175
  17. 17. Starr TN, Greaney AJ, Hilton SK, Ellis D, Crawford KHD, Dingens AS, et al. Deep Mutational Scanning of SARS-CoV-2 Receptor Binding Domain Reveals Constraints on Folding and ACE2 Binding. Cell. 2020;182(5):1295–310.e20. pmid:32841599.
  18. 18. Li Q, Nie J, Wu J, Zhang L, Ding R, Wang H, et al. SARS-CoV-2 501Y.V2 variants lack higher infectivity but do have immune escape. Cell. 2021;184(9):2362–71.e9. pmid:33735608
  19. 19. Weisblum Y, Schmidt F, Zhang F, DaSilva J, Poston D, Lorenzi JCC, et al. Escape from neutralizing antibodies by SARS-CoV-2 spike protein variants. eLife. 2020;9:e61312. pmid:33112236
  20. 20. Zahradník J, Marciano S, Shemesh M, Zoler E, Harari D, Chiaravalli J, et al. SARS-CoV-2 variant prediction and antiviral drug design are enabled by RBD in vitro evolution. Nature microbiology. 2021;6(9):1188–98. pmid:34400835.
  21. 21. Hoffmann M, Hofmann-Winkler H, Krüger N, Kempf A, Nehlmeier I, Graichen L, et al. SARS-CoV-2 variant B.1.617 is resistant to bamlanivimab and evades antibodies induced by infection and vaccination. Cell Reports. 2021;36(3). pmid:34270919
  22. 22. ECDC. https://www.ecdc.europa.eu/en/covid-19/variants-concern. 2021.
  23. 23. Yin W, Xu Y, Xu P, Cao X, Wu C, Gu C, et al. Structures of the Omicron spike trimer with ACE2 and an anti-Omicron antibody. Science. 2022;375(6584):1048–53. pmid:35133176.
  24. 24. Khan MI, Baig MH, Mondal T, Alorabi M, Sharma T, Dong J-J, et al. Impact of the Double Mutants on Spike Protein of SARS-CoV-2 B.1.617 Lineage on the Human ACE2 Receptor Binding: A Structural Insight. Viruses. 2021;13(11): 2295. https://doi.org/10.3390/v13112295.
  25. 25. Kumar S, Chandele A, Sharma A. Current status of therapeutic monoclonal antibodies against SARS-CoV-2. PLOS Pathogens. 2021;17(9):e1009885. https://doi.org/10.1371/journal.ppat.1009885.
  26. 26. Wilhelm A, Toptan T, Pallas C, Wolf T, Goetsch U, Gottschalk R, et al. Antibody-Mediated Neutralization of Authentic SARS-CoV-2 B.1.617 Variants Harboring L452R and T478K/E484Q. Viruses. 2021;13(9):1693. pmid:34578275.
  27. 27. Focosi D, Maggi F. Neutralising antibody escape of SARS-CoV-2 spike protein: Risk assessment for antibody-based Covid-19 therapeutics and vaccines. Reviews in Medical Virology. 2021;31(6):e2231. pmid:33724631
  28. 28. Sakkhachornphop S, Jiranusornkul S, Kodchakorn K, Nangola S, Sirisanthana T, Tayapiwatana C. Designed zinc finger protein interacting with the HIV-1 integrase recognition sequence at 2-LTR-circle junctions. Protein Science. 2009;18(11):2219–30. pmid:19701937
  29. 29. Tue-ngeun P, Kodchakorn K, Nimmanpipug P, Lawan N, Nangola S, Tayapiwatana C, et al. Improved scFv Anti-HIV-1 p17 Binding Affinity Guided from the Theoretical Calculation of Pairwise Decomposition Energies and Computational Alanine Scanning. BioMed Research International. 2013;2013:713585. pmid:24308004
  30. 30. Kodchakorn K, Dokmaisrijan S, Chong WL, Payaka A, Wisitponchai T, Nimmanpipug P, et al. GPU Accelerated Molecular Dynamics Simulations for Protein-Protein Interaction of Ankyrin Complex. Integrated Ferroelectrics. 2014;156(1):137–46. https://doi.org/10.1080/10584587.2014.906894.
  31. 31. Fanhchaksai K, Kodchakorn K, Pothacharoen P, Kongtawelert P. Effect of sesamin against cytokine production from influenza type A H1N1-induced peripheral blood mononuclear cells: computational and experimental studies. In Vitro Cellular & Developmental Biology—Animal. 2016;52(1):107–19. pmid:26424131
  32. 32. Kodchakorn K, Poovorawan Y, Suwannakarn K, Kongtawelert P. Molecular modelling investigation for drugs and nutraceuticals against protease of SARS-CoV-2. Journal of Molecular Graphics and Modelling. 2020;101:107717. pmid:32861974
  33. 33. Kodchakorn K, Viriyakhasem N, Wongwichai T, Kongtawelert P. Structural Determination, Biological Function, and Molecular Modelling Studies of Sulfoaildenafil Adulterated in Herbal Dietary Supplement. Molecules. 2021;26(4). https://doi.org/10.3390/molecules26040949.
  34. 34. Abel R, Wang L, Harder ED, Berne BJ, Friesner RA. Advancing Drug Discovery through Enhanced Free Energy Calculations. Accounts of Chemical Research. 2017;50(7):1625–32. pmid:28677954
  35. 35. Clark AJ, Negron C, Hauser K, Sun M, Wang L, Abel R, et al. Relative Binding Affinity Prediction of Charge-Changing Sequence Mutations with FEP in Protein–Protein Interfaces. Journal of Molecular Biology. 2019;431(7):1481–93. pmid:30776430
  36. 36. Ford MC, Babaoglu K. Examining the Feasibility of Using Free Energy Perturbation (FEP+) in Predicting Protein Stability. Journal of Chemical Information and Modeling. 2017;57(6):1276–85. pmid:28520421
  37. 37. Dejnirattisai W, Zhou D, Ginn HM, Duyvesteyn HME, Supasa P, Case JB, et al. The antigenic anatomy of SARS-CoV-2 receptor binding domain. Cell. 2021;184(8):2183–200.e22. pmid:33756110
  38. 38. Eileen S, Lukas H, Friedrich P, Friederike Z, Philipp A. Molecular dynamics simulations of the delta and omicron SARS-CoV-2 spike–ACE2 complexes reveal distinct changes between both variants. Research Square. 2022.
  39. 39. Wu F, Zhao S, Yu B, Chen Y-M, Wang W, Song Z-G, et al. A new coronavirus associated with human respiratory disease in China. Nature. 2020;579(7798):265–9. pmid:32015508
  40. 40. Lan J, Ge J, Yu J, Shan S, Zhou H, Fan S, et al. Structure of the SARS-CoV-2 spike receptor-binding domain bound to the ACE2 receptor. Nature. 2020;581(7807):215–20. pmid:32225176
  41. 41. Götz AW, Williamson MJ, Xu D, Poole D, Le Grand S, Walker RC. Routine Microsecond Molecular Dynamics Simulations with AMBER on GPUs. 1. Generalized Born. Journal of Chemical Theory and Computation. 2012;8(5):1542–55. pmid:22582031.
  42. 42. Salomon-Ferrer R, Götz AW, Poole D, Le Grand S, Walker RC. Routine Microsecond Molecular Dynamics Simulations with AMBER on GPUs. 2. Explicit Solvent Particle Mesh Ewald. Journal of Chemical Theory and Computation. 2013;9(9):3878–88. pmid:26592383
  43. 43. Case D.A., Ben-Shalom I.Y., Brozell S.R., Cerutti D.S., Cheatham I T.E., Cruzeiro V.W.D., et al. AMBER 2018, University of California, San Francisco. 2018.
    • 44. Essmann U, Perera L, Berkowitz ML, Darden T, Lee H, Pedersen LG. A smooth particle mesh Ewald method. The Journal of Chemical Physics. 1995;103(19):8577–93. https://doi.org/10.1063/1.470117.
    • 45. Roe DR, Cheatham TE. PTRAJ and CPPTRAJ: Software for Processing and Analysis of Molecular Dynamics Trajectory Data. Journal of Chemical Theory and Computation. 2013;9(7):3084–95. pmid:26583988
    • 46. Wang E, Sun H, Wang J, Wang Z, Liu H, Zhang JZH, et al. End-Point Binding Free Energy Calculation with MM/PBSA and MM/GBSA: Strategies and Applications in Drug Design. Chemical Reviews. 2019;119(16):9478–508. pmid:31244000
    • 47. Moreira IS, Fernandes PA, Ramos MJ. Computational alanine scanning mutagenesis—an improved methodological approach. Journal of computational chemistry. 2007;28(3):644–54. pmid:17195156.
    • 48. Simões IC, Costa IP, Coimbra JT, Ramos MJ, Fernandes PA. New Parameters for Higher Accuracy in the Computation of Binding Free Energy Differences upon Alanine Scanning Mutagenesis on Protein-Protein Interfaces. Journal of Chemical Information and Modeling. 2017;57(1):60–72. pmid:27936711.
    • 49. Hünenberger PH, Mark AE, van Gunsteren WF. Fluctuation and Cross-correlation Analysis of Protein Motions Observed in Nanosecond Molecular Dynamics Simulations. Journal of Molecular Biology. 1995;252(4):492–503. pmid:7563068
    • 50. Martinez AM, Kak AC. PCA versus LDA. IEEE Transactions on Pattern Analysis and Machine Intelligence. 2001;23(2):228–33. https://doi.org/10.1109/34.908974.
    • 51. Agoni C, Ramharack P, Soliman MES. Allosteric inhibition induces an open WPD-loop: a new avenue towards glioblastoma therapy. RSC Advances. 2018;8(70):40187–97. pmid:35558220
    • 52. Xu C, Wang Y, Liu C, Zhang C, Han W, Hong X, et al. Conformational dynamics of SARS-CoV-2 trimeric spike glycoprotein in complex with receptor ACE2 revealed by cryo-EM. Science Advances. 7(1):eabe5575. pmid:33277323
    • 53. Wu K, Chen L, Peng G, Zhou W, Pennell Christopher A, Mansky Louis M, et al. A Virus-Binding Hot Spot on Human Angiotensin-Converting Enzyme 2 Is Critical for Binding of Two Different Coronaviruses. Journal of Virology. 2011;85(11):5331–7. pmid:21411533
    • 54. Li W, Moore MJ, Vasilieva N, Sui J, Wong SK, Berne MA, et al. Angiotensin-converting enzyme 2 is a functional receptor for the SARS coronavirus. Nature. 2003;426(6965):450–4. pmid:14647384
    • 55. Ortega JT, Serrano ML, Pujol FH, Rangel HR. Role of changes in SARS-CoV-2 spike protein in the interaction with the human ACE2 receptor: An in silico analysis. EXCLI J. 2020;19:410–7. pmid:32210742.
    • 56. Wang Q, Zhang Y, Wu L, Niu S, Song C, Zhang Z, et al. Structural and Functional Basis of SARS-CoV-2 Entry by Using Human ACE2. Cell. 2020;181(4):894–904.e9. pmid:32275855
    • 57. Wrapp D, Wang N, Corbett Kizzmekia S, Goldsmith Jory A, Hsieh C-L, Abiona O, et al. Cryo-EM structure of the 2019-nCoV spike in the prefusion conformation. Science. 2020;367(6483):1260–3. pmid:32075877
    • 58. Kodchakorn K, Chokepaichitkool T, Kongtawelert P. Mutational scanning of spike RBD protein for enhanced ACE2 affinity emerging Southeast Asia in the late transmission phase. Scientific Reports. 2022;12(1):5896. pmid:35393512
    • 59. Khan A, Khan SA, Zia K, Altowyan MS, Barakat A, Ul-Haq Z. Deciphering the Impact of Mutations on the Binding Efficacy of SARS-CoV-2 Omicron and Delta Variants with Human ACE2 Receptor. Frontiers in Chemistry. 2022;10. pmid:35755247
    • 60. Ali A, Vijayan R. Dynamics of the ACE2–SARS-CoV-2/SARS-CoV spike protein interface reveal unique mechanisms. Scientific Reports. 2020;10(1):14214. pmid:32848162
    • 61. Li Q, Wu J, Nie J, Zhang L, Hao H, Liu S, et al. The Impact of Mutations in SARS-CoV-2 Spike on Viral Infectivity and Antigenicity. Cell. 2020;182(5):1284–94.e9. pmid:32730807
    • 62. Yan R, Zhang Y, Li Y, Xia L, Guo Y, Zhou Q. Structural basis for the recognition of SARS-CoV-2 by full-length human ACE2. Science. 2020;367(6485):1444–8. pmid:32132184
    • 63. Shah M, Ahmad B, Choi S, Woo HG. Mutations in the SARS-CoV-2 spike RBD are responsible for stronger ACE2 binding and poor anti-SARS-CoV mAbs cross-neutralization. Computational and Structural Biotechnology Journal. 2020;18:3402–14. pmid:33200028
    • 64. Coevolution Verkhivker G., Dynamics and Allostery Conspire in Shaping Cooperative Binding and Signal Transmission of the SARS-CoV-2 Spike Protein with Human Angiotensin-Converting Enzyme 2. International Journal of Molecular Sciences. 2020;21(21). https://doi.org/10.3390/ijms21218268.
    • 65. Laurini E, Marson D, Aulic S, Fermeglia M, Pricl S. Computational Alanine Scanning and Structural Analysis of the SARS-CoV-2 Spike Protein/Angiotensin-Converting Enzyme 2 Complex. ACS Nano. 2020;14(9):11821–30. pmid:32833435
    • 66. Sousa SF, Fernandes PA, Ramos MJ. Molecular dynamics simulations on the critical states of the farnesyltransferase enzyme. Bioorganic & Medicinal Chemistry. 2009;17(9):3369–78. pmid:19369081
    • 67. Król M, Roterman I, Piekarska B, Konieczny L, Rybarska J, Stopa B, et al. Analysis of correlated domain motions in IgG light chain reveals possible mechanisms of immunological signal transduction. Proteins: Structure, Function, and Bioinformatics. 2005;59(3):545–54. pmid:15778960
    • 68. Wan Y, Shang J, Graham R, Baric Ralph S, Li F. Receptor Recognition by the Novel Coronavirus from Wuhan: an Analysis Based on Decade-Long Structural Studies of SARS Coronavirus. Journal of Virology. 2020;94(7):e00127–20. pmid:31996437
    • 69. Wang Y, Liu M, Gao J. Enhanced receptor binding of SARS-CoV-2 through networks of hydrogen-bonding and hydrophobic interactions. Proceedings of the National Academy of Sciences. 2020;117(25):13967–74. pmid:32503918
    • 70. Teng S, Sobitan A, Rhoades R, Liu D, Tang Q. Systemic effects of missense mutations on SARS-CoV-2 spike glycoprotein stability and receptor-binding affinity. Briefings in bioinformatics. 2021;22(2):1239–53. pmid:33006605.
    • 71. Wang J, Xu X, Zhou X, Chen P, Liang H, Li X, et al. Molecular simulation of SARS-CoV-2 spike protein binding to pangolin ACE2 or human ACE2 natural variants reveals altered susceptibility to infection. The Journal of general virology. 2020;101(9):921–4. pmid:32538738.
    • 72. Chen J, Wang R, Wang M, Wei G-W. Mutations Strengthened SARS-CoV-2 Infectivity. Journal of Molecular Biology. 2020;432(19):5212–26. pmid:32710986
    • 73. Tian F, Tong B, Sun L, Shi S, Zheng B, Wang Z, et al. N501Y mutation of spike protein in SARS-CoV-2 strengthens its binding to receptor ACE2. eLife. 2021;10:e69091. pmid:34414884
    • 74. Antony P, Vijayan R. Molecular Dynamics Simulation Study of the Interaction between Human Angiotensin Converting Enzyme 2 and Spike Protein Receptor Binding Domain of the SARS-CoV-2 B.1.617 Variant. Biomolecules. 2021; 11(8). pmid:34439910
    • 75. Kumar V, Singh J, Hasnain SE, Sundar D. Possible Link between Higher Transmissibility of Alpha, Kappa and Delta Variants of SARS-CoV-2 and Increased Structural Stability of Its Spike Protein and hACE2 Affinity. International Journal of Molecular Sciences. 2021; 22(17). pmid:34502041
    • 76. Wu K, Peng G, Wilken M, Geraghty RJ, Li F. Mechanisms of Host Receptor Adaptation by Severe Acute Respiratory Syndrome Coronavirus *. Journal of Biological Chemistry. 2012;287(12):8904–11. pmid:22291007
    • 77. Li F, Li W, Farzan M, Harrison SC. Structure of SARS Coronavirus Spike Receptor-Binding Domain Complexed with Receptor. Science. 2005;309(5742):1864–8. pmid:16166518
    • 78. Fratev F. N501Y and K417N Mutations in the Spike Protein of SARS-CoV-2 Alter the Interactions with Both hACE2 and Human-Derived Antibody: A Free Energy of Perturbation Retrospective Study. Journal of Chemical Information and Modeling. 2021;61(12):6079–84. pmid:34806876
    • 79. Cameroni E, Bowen JE, Rosen LE, Saliba C, Zepeda SK, Culap K, et al. Broadly neutralizing antibodies overcome SARS-CoV-2 Omicron antigenic shift. Nature. 2022;602(7898):664–70. pmid:35016195
    • 80. Hodcroft EB, Zuber M, Nadeau S, Vaughan TG, Crawford KHD, Althaus CL, et al. Spread of a SARS-CoV-2 variant through Europe in the summer of 2020. Nature. 2021;595(7869):707–12. pmid:34098568
    • 81. Kim S, Liu Y, Lei Z, Dicker J, Cao Y, Zhang XF, et al. Differential Interactions Between Human ACE2 and Spike RBD of SARS-CoV-2 Variants of Concern. bioRxiv. 2021:2021.07.23.453598. pmid:34341794
    • 82. Cho H, Gonzales-Wartz KK, Huang D, Yuan M, Peterson M, Liang J, et al. Bispecific antibodies targeting distinct regions of the spike protein potently neutralize SARS-CoV-2 variants of concern. Science Translational Medicine. 13(616):eabj5413. pmid:34519517
    • 83. Chakrabarti S, Panchenko AR. Coevolution in defining the functional specificity. Proteins: Structure, Function, and Bioinformatics. 2009;75(1):231–40. https://doi.org/10.1002/prot.22239.
    • 84. Lee B-C, Park K, Kim D. Analysis of the residue–residue coevolution network and the functionally important residues in proteins. Proteins: Structure, Function, and Bioinformatics. 2008;72(3):863–72. https://doi.org/10.1002/prot.21972.
    • 85. Isa DM, Chin SP, Chong WL, Zain SM, Rahman NA, Lee VS. Dynamics and binding interactions of peptide inhibitors of dengue virus entry. Journal of Biological Physics. 2019;45(1):63–76. pmid:30680580
    • 86. Moeller NH, Shi K, Demir Ö, Belica C, Banerjee S, Yin L, et al. Structure and dynamics of SARS-CoV-2 proofreading exoribonuclease ExoN. Proceedings of the National Academy of Sciences. 2022;119(9):e2106379119. pmid:35165203
    • 87. Manjula S, Sivanandam M, Kumaradhas P. Probing the “fingers” domain binding pocket of Hepatitis C virus NS5B RdRp and D559G resistance mutation via molecular docking, molecular dynamics simulation and binding free energy calculations. Journal of Biomolecular Structure and Dynamics. 2019;37(9):2440–56. pmid:30047829
    • 88. Chen J, Wang J, Zhu W. Molecular Mechanism and Energy Basis of Conformational Diversity of Antibody SPE7 Revealed by Molecular Dynamics Simulation and Principal Component Analysis. Scientific Reports. 2016;6(1):36900. pmid:27830740