key: cord-0803761-ejen6hpc authors: Ghorbani, Mahdi; Brooks, Bernard R.; Klauda, Jeffery B. title: Critical Sequence Hotspots for Binding of Novel Coronavirus to Angiotensin Converter Enzyme as Evaluated by Molecular Simulations date: 2020-10-28 journal: J Phys Chem B DOI: 10.1021/acs.jpcb.0c05994 sha: 5c80d824721cf9bb82d1c3b9f091840dd37aac42 doc_id: 803761 cord_uid: ejen6hpc [Image: see text] The novel coronavirus (nCOV-2019) outbreak has put the world on edge, causing millions of cases and hundreds of thousands of deaths all around the world, as of June 2020, let alone the societal and economic impacts of the crisis. The spike protein of nCOV-2019 resides on the virion’s surface mediating coronavirus entry into host cells by binding its receptor binding domain (RBD) to the host cell surface receptor protein, angiotensin converter enzyme (ACE2). Our goal is to provide a detailed structural mechanism of how nCOV-2019 recognizes and establishes contacts with ACE2 and its difference with an earlier severe acute respiratory syndrome coronavirus (SARS-COV) in 2002 via extensive molecular dynamics (MD) simulations. Numerous mutations have been identified in the RBD of nCOV-2019 strains isolated from humans in different parts of the world. In this study, we investigated the effect of these mutations as well as other Ala-scanning mutations on the stability of the RBD/ACE2 complex. It is found that most of the naturally occurring mutations to the RBD either slightly strengthen or have the same binding affinity to ACE2 as the wild-type nCOV-2019. This means that the virus had sufficient binding affinity to its receptor at the beginning of the crisis. This also has implications for any vaccine design endeavors since these mutations could act as antibody escape mutants. Furthermore, in silico Ala-scanning and long-timescale MD simulations highlight the crucial role of the residues at the interface of RBD and ACE2 that may be used as potential pharmacophores for any drug development endeavors. From an evolutional perspective, this study also identifies how the virus has evolved from its predecessor SARS-COV and how it could further evolve to become even more infectious. The novel coronavirus (nCOV-2019) outbreak emerging from China has become a global pandemic and a major threat for human public health. According to the World Health Organization as of August 28th 2020, there has been about 25 million confirmed cases and approaching 900,000 deaths because of coronavirus in the world. 1, 2 Much of the human population including the United States of America were under lockdown or official stay-at-home orders to minimize the continued spread of the virus. Coronaviruses are a family of single-stranded enveloped RNA viruses. Phylogenetic analysis of coronavirus genome has shown that nCOV-2019 belongs to the beta-coronavirus family, which also includes Middle East respiratory syndrome coronavirus (MERS-COV), severe acute respiratory syndrome coronavirus (SARS-COV), and bat-SARS-related coronaviruses. 3, 4 It is worth mentioning that SARS-COV, which was widespread in 2002 caused more than 8000 cases and about 800 deaths and MERS-COV in 2012 also spread in more than 25 countries, causing about 2500 cases and more than 850 deaths. (www.who.int/health-topics/coronavirus). 5 In all coronaviruses, a homotrimeric spike glycoprotein on the virion's envelope mediates coronavirus entry into host cells through a mechanism of receptor binding followed by fusion of viral and host membranes. 3, 6 Coronavirus spike protein contains two functional subunits S1 and S2. The S1 subunit is responsible for binding to host cell receptor, and the S2 subunit is responsible for fusion of viral and host cell membranes. 3, 7 The spike protein in nCOV-2019 exists in a meta-stable pre-fusion conformation that undergoes a substantial conformational rearrangement to fuse the viral membrane with the host cell membrane. 7, 8 nCOV-2019 is closely related to bat coronavirus RaTG13 with about 93.1% sequence similarity in the spike protein gene. The sequence similarity of nCOV-2019 and SARS-COV is less than 80% in the spike sequence. 2 The S1 subunit in the spike protein includes a receptor binding domain (RBD) that recognizes and binds to the host cells receptor. The RBD of nCOV-2019 shares 72.8% sequence identity to SARS-COV RBD and the root mean squared deviation (RMSD) for the structure between the two proteins is 1.2Å, which shows the high structural similarity. 4, 8, 9 Experimental binding affinity measurements using surface plasmon resonance (SPR) have shown that nCOV-2019 spike protein binds its receptor human angiotensin converter enzyme (ACE2) with 10 to 20 fold higher affinity than SARS-COV binding to ACE2. 7 Based on the sequence similarity between RBD of nCOV-2019 and SARS-COV and also the tight binding between the RBD of nCOV-2019 and ACE2, it is most probable that nCOV-2019 uses this receptor on human cells to gain entry into the body. 3, 6, 7, 10 The spike protein and specifically the RBD in coronaviruses have been a major target for therapeutic antibodies. However, no monoclonal antibodies targeted to RBD have been able to bind efficiently and neutralize nCOV-2019. 7, 11 The core of nCOV-2019 RBD is a five-stranded antiparallel β-sheet with connected short α-helices and loops ( Figure 1 ). The binding interface of nCOV-2019 and SARS-COV with ACE2 is very similar with less than 1.3ÅRMSD. An extended insertion inside the core containing short strands, α-helices, and loops called the receptor binding motif (RBM) makes all the contacts with ACE2. In nCOV-2019 RBD, the RBM forms a concave surface with a ridge loop on one side and it binds to a convex exposed surface of ACE2. The overlay of SARS and nCOV-2019 RBD proteins is shown in Figure 1A . The binding interface in nCOV-2019 contains loops L 1 to L 4 and short βstrands β 5 and β 6 and a short helix α 5 . The location of RBM in nCOV-2019 RBD as well as different helices, strands, and loops is shown in Figure 1B . The sequence alignment between SARS-COV in human, SARS civet, Bat RaTG13 coronavirus, and nCOV-2019 in the RBM is shown in Figure 2 . There is a 50% sequence similarity between the RBM of nCOV-2019 and SARS-COV. RBM mutations played an important role in the SARS epidemic in 2002. 3, 12 Two mutations in the RBM of SARS-2002 from SARS-Civet were observed from strains of these viruses. These two mutations were K479N and S487T. These two residues are close to the virus binding hotspots in ACE2 including hotspot-31 and hotspot-353. Hotspot-31 centers on the saltbridge between K31-E35 and hotspot-353 are centered on the salt-bridge between K353-E358 on ACE2. Residues K479 and S487 in SARS-Civet are in close proximity with these hotspots and mutations at these residues caused SARS to bind ACE2 with significantly higher affinity than SARS-civet and played a major role in civet-to-human and human-to-human transmission of SARS coronavirus in 2002. 3,13−15 Numerous mutations in the interface of SARS-COV RBD and ACE2 from different strains of SARS isolated from humans in 2002 have been identified and the effect of these mutations on binding ACE2 has been investigated by SPR. 14,16 Two identified RBD mutations (Y442F and L472F) increased the binding affinity of SARS-COV to ACE2 and two mutations (N479K, T487S) decreased the binding affinity. It was demonstrated that these mutations were viral adaptations to either human or civet ACE2. 14, 16 A pseudotyped viral infection assay of the interaction between different spike proteins and ACE2 confirmed the correlation between high affinity mutants and their high infection. 16 Further investigation of RBD residues in binding of SARS-COV and ACE2 was performed through ala-scanning mutagenesis, which resulted in identification of residues that reduce binding affinity to ACE2 upon mutation to alanine. 17 RBD mutations have also been identified in MERS-COV, which affected their affinity to receptor (DPP4) on human cells. 14 Multiple monoclonal antibodies have been developed for SARS since 2002 that neutralized the spike glycoprotein on the SARS-COV surface. 18−22 However, multiple escape mutations exist in the RBD of SARS-COV that affect neutralization with antibodies, which led to the use of a cocktail of antibodies as a robust treatment. 23 Full genome analysis of nCOV-2019 in different countries and the receptor binding surveillance have shown multiple mutations in the RBD of glycosylated spike. The GISAID database 24 (www.gisaid.org/) contains genomes on nCOV-2019 from researchers across the world since December 2019. The latest report by the GISAID database on June 2020 has shown 25 different variants of RBD from strains of nCOV-2019 collected from different countries along with the number of occurrences in these regions which is listed below for the seven most occurring mutations: 213x N439K (211 Scotland, England, and Romania), 65x T478I in England, 30x V483A (26 USA/WA, 2 USA/UN, USA/CT, and England), 10x G476S (8 USA/WA, USA/OR, and Belgium), 7x S494P (3 USA/MI, England, Spain, India, and Sweden), 5x V483F (4x Spain and England), and 4x A475V (2 USA/AZ, USA/NY, and Australia/NSW). It is not known whether these mutations are linked to the severity of coronavirus in these regions. Starr and co-workers 25 performed a deep mutational scanning of nCOV-2019 RBD and used flow cytometry to measure the effect of single mutations on the expression of the folded protein as well as its binding affinity to ACE2. They showed that RBD is very tolerant to these mutations to maintain its expression level as well as binding affinity to ACE2 in most cases. According to their results, most natural mutations exert similar binding affinities to ACE2 as wild-type nCOV-2019. Furthermore, they showed that mutations at critical positions at the RBD-ACE2 interface at nCOV-2019 such as residues Q493 and Q498 do not reduce the binding affinity to ACE2 which shows the substantial plasticity of the interface. 25 Different groups have computationally studied the binding of nCOV-2019 RBD with ACE2. 25−29 All these studies point to the higher binding affinity of nCOV-2019 RBD than SARS-COV RBD to ACE2. Interestingly, the role of water-mediated interactions has been pointed out to be a driving force which is shown to be similar for both SARS-COV and nCOV-2019 RBD. 27 Spinello and co-workers 30 studied the binding of nCOV-2019 and SARS-COV RBD to ACE2 and found that the former binds its receptor with 30 kcal/mol higher affinity than SARS-COV RBD. Gao et al. 31 used free energy perturbation (FEP) and showed that most amino acid mutations at the RBM from SARS-COV to nCOV-2019 increase the affinity of RBD to ACE2. The focus of this article is to elucidate the differences between the interface of SARS-COV and nCOV-2019 with ACE2 to understand with atomic resolution the interaction mechanism and hotspot residues at the RBD/ACE2 interface using long-timescale molecular dynamics (MD) simulation. An alanine-scanning mutagenesis in the RBM of nCOV-2019 helped to identify the key residues in the interaction, which could be used as potential pharmacophores for future drug development. Furthermore, we performed molecular simulations on the seven most common mutations found from the surveillance of RBD mutations N439K, T478I, V483A, G476S, S494P, V483F, and A475V. From an evolutionary perspective this study shows the residues in which the virus might further evolve to be even more dangerous to human health. nCOV-2019 shares 76% sequence similarity with SARS-2002 spike protein, 73% sequence identity for RBD and 50% for the RBM. 1 Bat coronavirus RaTG13 seems to be the closest relative of nCOV-2019 sharing about 93% sequence identity in the spike protein. 6 The sequence alignment of SARS-2002 (accession number: AFR58742), SARS-civet (accession number: AY304486.1), Bat RaTG13 (accession number: MN996532.1), and nCOV-2019 (accession number: MN908947.1) is shown in Figure 2 . 6 To investigate the roles of critical mutations on the complex stability of nCOV-2019 with ACE2, mCSM-PPI2 webserver 32 was used to find the residues in nCOV-2019 that are at the interface with ACE2. This method uses a graph-based signature framework and predicts the effect of alanine substitution at interface residues on the binding energy of the complex, and 21 different residues were identified to be in contact with ACE2 and were chosen for further MD simulation. On the other hand, mutations are also observed in the RBD domain from full genome analysis of different nCOV-2019 variants collected from different countries compiled in the GISAID database. 24 The mutations selected are listed in Table S1 along with their location in the RBD. MD Simulations. The crystal structure of nCOV-2019 in the complex with hACE2 (pdb id:6M0J) 17 as well as the SARS-COV complex with human ACE2 (pdb id: 6ACJ) 33 were obtained from RCSB (www.rcsb.org). 34 The RBD domain of nCOV-2019 comprises 194 residues (333−526) and SARS-COV includes 190 residues (323−512). ACE2 protein contains 597 residues (19−615) in the complex structure for both systems. All the structures including nCOV-2019, SARS-COV, and all 21 alanine substitutions of nCOV-2019 were prepared and solvated in GROMACS. 35 A TIP3P water model was used for the solvent and Param99SB-ILDN AMBER force field 36, 37 was used for all the complexes. A few counter ions were added to each system to neutralize the charges on the RBD and ACE2 as the PBSA method for binding energy calculation is known to be problematic with high ionic strength. Each system contained about 260,000 atoms. It is important to note that, none of the RBD/ACE2 complexes studied here were glycosylated. The glycosylation sites on RBD are far from the binding interface and do not interfere with the binding of RBD to ACE2. Moreover, there are nine Cys residues at the RBD of nCOV-2019 and eight of them form four pairs of disulfide bonds (Cys336-Cys361, Cys379-Cys432, Cys391-Cys525, and Cys480-Cys488). In total, 5000 steps of energy minimizations were done using the steepest descent algorithm. In all steps, the LINCS algorithm was used to constrain all bonds containing hydrogen atoms and a time step of 2 fs was used as the integration time step. Equilibration of all systems was performed in three steps. In the first step, 100,000 steps of simulation were performed using a velocity-rescaling thermostat to maintain the temperature at 310 K with a 0.1 ps coupling constant in an NVT ensemble under periodic boundary conditions and harmonic restraints on the backbone and side-chain atoms of the complex. 38 A velocity rescaling thermostat was used in all other steps of simulation. In the next step, we performed 300,000 steps in the isothermal-isobaric NPT ensemble at a temperature of 310 K and pressure of 1 bar using a Berendsen barostat. 39 This was done through decreasing the force constant of the restraint on the backbone and side-chain atoms of the complex from 1000 to 100 and finally to 10 kJ mol nm 2 . The Berendsen barostat was only used for the equilibration step because of its usefulness in rapidly correcting density. In the next step, the restraints were removed, and the systems were subjected to 1,000,000 steps of MD simulation under the NPT ensemble. In the production run, harmonic restraints were removed and all the systems were simulated using a NPT ensemble where the pressure was maintained at 1 bar using the Parrinello-Rahman barostat 40 with a compressibility of 4.5 × 10 −5 bar −1 and a coupling constant of 0.5 ps. It is important to note that the Berendsen barostat was only used for the equilibration step as it was shown that this barostat can cause unrealistic temperature gradients. 41 The production run lasted for 500 ns for SARS-COV and nCOV-2019 complexes and 300 ns for all the mutants with a 2 fs timestep and the particlemesh Ewald 42 for long range electrostatic interactions using the GROMACS 2018.3 package. 43 All mutant systems were constructed as described before and ran for 300 ns of production run. In addition, the simulation time for a few mutants (Y449A, T478I, Y489A, and S494P) was extended to 500 ns. Gibbs Free Energy and Correlated Motions. The last 400 ns of simulation was used to explore the dominant motions in SARS-COV, nCOV-2019 and the mutations with extended simulation, and last 200 ns for all other mutants using principal component analysis (PCA) as part of the quasiharmonic analysis method. 44 For this method the rotational and translational motions of RBD of all systems were eliminated by fitting to a reference (crystal) structure. Next, 4000 snapshots from the last 400 ns of SARS-COV, nCOV-2019 and mutations with extended simulation time, and 2000 snapshots from last 200 ns of all other mutant systems were taken to generate the covariance matrix between C α atoms of RBD. In the mutant systems with production run, the last 400 ns was used for this analysis. Diagonalization of this matrix resulted in a diagonal matrix of eigenvalues and their corresponding eigenvectors. 43, 45 The first eigenvector which indicate the first principal component was used to visualize the dominant global motions of all complexes through porcupine plots using the (PorcupinePlot.tcl) script in visual MD (VMD). The principal components were used to calculate and plot the approximate free energy landscape (aFEL). We refer to the FEL produced by this approach to be approximate in that the ensemble with respect to the first few PC's (lowest frequency quasiharmonic modes) is not close to convergence, but the analysis can still provide valuable information and insight. g_sham, g_covar, and g_anaeig functions in GROMACS 35 were used to obtain principal components and aFEL. The dynamic cross-correlation maps (DCCM) were obtained using the MD_TASK package to identify the correlated motions of RBD residues. 46 In DCCM, the crosscorrelation matrix C ij is obtained from displacement of backbone C α atoms at a time interval Δt. The DCCM was constructed using the last 400 ns of SARS-COV and nCOV-2019 and the extended mutant systems and last 200 ns of all other mutant systems with a 100 ps time interval. Hydrogen bonds were analyzed in VMD where the distance cutoff was 3.2Åand the angle cutoff between the donor and acceptor was 30°. Binding Free Energy from Molecular Mechanics Poisson-Boltzmann Surface Area Method. The molecular mechanics Poisson-Boltzmann surface area (MMPBSA) method was applied to calculate the binding energy between RBD and ACE2 in all complexes. 47, 48 For SARS-COV and nCOV-2019, 200 snapshots of the last 400 ns and for the mutant systems and 100 snapshots of the last 200 ns simulation were used for the calculation of binding free energies with an interval of 2 ns. Simulation for a few mutant systems (Y449A, T478I, Y489A, and S494P) was extended to 400 ns, and the binding energies were calculated for the last 400 ns to assess the convergence of free energies. The binding free energy of a ligand−receptor complex can be calculated as: In these equations, ΔE MM , ΔG bind, solv , and −TΔS are calculated in the gas phase. ΔE MM is the gas phase molecular mechanical energy changes which includes covalent, electrostatic, and van der Waals energies. Based on previous studies, the entropy change during binding is small and neglected in these calculations. 48−50 ΔG bind, solv is the solvation free energy which comprises the polar and nonpolar components. The polar solvation is calculated using the MMPBSA method by setting a value of 80 and 2 for solvent and solute dielectric ■ RESULTS Structural Dynamics. To compute the RMSD of systems, the rotational and translational movements were removed by first fitting the C α atoms of the RBD to the crystal structure and then computing the RMSD with respect to the C α atoms of RBD in each system. Figure 3 shows the RMSD plot in the RBD of SARS-COV, nCOV-2019, and some of its variants. Comparison of the RMSD of SARS and nCOV-2019 RBD shows that SARS-COV has a larger RMSD throughout the 500 ns simulation. In nCOV-2019, the RMSD is very stable with a value of about 1.5 Å, whereas in SARS-COV, the RMSD increases up to ∼4 Å after 100 ns and then fluctuates between 3 and 4 Å. The change in RMSD of SARS is partially related to the motion in the C-terminal which is a flexible loop. The RMSD plots for the nCOV-2019 mutants show similar behaviors to nCOV-2019-wt. In most of the variants, the RMSD is very stable during the 300 ns simulation which shows the great tolerance of the interface for mutations. However, a few mutations showed some RMSD variance. In mutation Y489A, the RMSD increases from 1.37 ± 0.21 Å to1.88 ± 0.16 Å after 2000 ns. Mutation Y505A resulted in an increase in RMSD up to 100 ns to a value to 1.98 ± 0.20 Å and decreased afterward. The RMSD for mutation N487A shows an increasing behavior with a value of 2.10 ± 0.23 Å. Mutations N439K, V483A, and V483F showed a stable RMSD of ∼1.5 Å. For mutations T478I, G476S, S494P, and A475V, the RMSD increases up to ∼2 Å. These variations in the backbone RMSD show the involvement of these residues in the complex stability. RMSD plots for other mutations are shown in Figure S1 . Since the extended loop (residues 449 to 510 shown in Figure 1B from α 4 to α 5 ) of the RBD makes all contacts with ACE2, the RMSD was computed by also fitting to the C α atoms of this region ( Figure S2 ). The extended loop in nCOV-2019 is very stable with less deviation (RMSD = 0.86 ± 0.017 Å) from the crystal structure compared to SARS-COV having an RMSD of 2.79 ± 0.05 Å. Few of the mutants show an increase in the loop RMSD during the simulation. In mutants N487A and Y449A the loop RMSD jumps to a value of about 1.95 ± 0.12 Å and 1.94 ± 0.24 Å, respectively, after about 200 ns of simulation. Mutants G447A and E484A show a loop RMSD values of 2.22 ± 0.03 Å and 1.96 ± 0.02 Å during the last 100 ns of simulation. Other mutants showed a stable extended loop (observed in loop RMSD) during simulation. The stability of extended loop for mutant systems confirms the high tolerance of this region of RBD. To characterize the dynamic behavior for each amino acid in the RBD, we analyzed the root mean square fluctuation Figure S3 . PCA and aFEL. To identify the dominant motions in the nCOV-2019, SARS-COV, and all the mutants, PCA was performed. Most of the combined motions were captured by the first ten eigenvectors generated from the last 400 ns for SARS-COV, nCOV-2019, and extended mutant systems and the last 200 ns for other nCOV-2019 mutants. The percentage of the motions captured by the first three eigenvectors was 51% for nCOV-2019 and 68% for SARS-COV. In all mutations, more than 50% of the motions were captured by the first three eigenvectors. The first few PC's describe the highest motions in a protein which are related to a functional motion such as binding or unbinding of protein from the receptor. The first three eigenvectors were used to calculate the aFEL using the last 400 ns of simulation for nCOV-2019 and SARS-COV as shown in Figure 5 , which displays the variance in conformational motion. SARS-COV showed two distinct low free energy states shown as blue separated by a metastable state. There is a clear separation between the two regions by a free energy barrier of about 6−7.5 kcal/mol. These two states correspond to the loop motions in the L 3 as well as the motion in C-terminal residues of SARS-COV. The L 3 motion in nCOV-2019 is stabilized by the H-bond between N487 on The Journal of Physical Chemistry B pubs.acs.org/JPCB Article RBD and Y83 on ACE2 as well as a π-stacking interaction between F486 and Y83. It is evident that the nCOV-2019 RBD is more stable than SASR-COV RBD and exists in one conformation whereas the SARS-COV interface fluctuates and the aFEL is separated into two different regions. The first two eigenvectors were used to calculate and plot the aFEL as a function of first two principal components using the last 200 ns of the simulation for mutant systems. aFEL for other systems is shown in Figure S4 . In each system, the first eigenvector was used to construct the porcupine plots to visualize the most dominant movements ( Figure S5 ). nCOV-2019 showed a small motion in L 3 and the core region and the extended loop region are very rigid showing small cones in the porcupine plot. The core structure of the RBD remains dormant as the cones are blue in most of the regions ( Figure S5-A) . In SARS-COV, the C-terminal region shows large motions ( Figure S5-B) . Mutation N487A showed a large motion in L 1 ( Figure S5 -C). Mutations Y449A, G447A, and E484A demonstrate large motions in L 3 ( Figure S5 ). Overall, these plots show the involvement of these residues in the dynamic stability of the RBD/ACE2 complex. DCCM. The correlated motions of RBD atoms were also analyzed with the DCCM based on the C α atoms of RBD from the last 400 ns of simulation for nCOV-2019, SARS-COV, and extended mutant systems and the last 200 ns for the other mutant systems (Figure 6 ). The DCCM for nCOV-2019 showed a correlation between residues 490−505 (containing α 5 , L 4 and β 5 regions) and residues 440−455 (containing α 4 , L 1 and β 5 regions) shown in the red rectangle in Figure 6 . This correlation showed the coordination of these regions for binding ACE2 effectively. Another important correlation that appears in the DCCM of nCOV-2019 is between residues 473−481 with residues 482−491. These residues are in L 3 and β 6 regions and their correlation in nCOV-2019 is stronger than SARS-COV. This is due to the presence of the β 6 strand in nCOV-2019, whereas in SARS-COV these residues all belong to L 3 . This indicates that L 3 in nCOV-2019 has evolved from SARS-COV to adopt a new secondary structure, which causes strong correlation and makes the loop act as a recognition region for binding. The correlation in L 3 is shown as a blue rectangle in Figure 6 . Some of the mutations disrupted the patterns of correlation and anticorrelation in nCOV-2019-wt. Mutation N487A showed a stronger correlation in L 3 and β 6 strand than the wild-type RBD. In mutation E484A, correlation in L 3 is stronger than the nCOV-2019-wt. DCCM for other mutants are shown in Figure S6 . It is worth mentioning that mutation F486A disrupts the DCCM of nCOV-2019 by introducing strong correlations in the core region of RBD as well as the extended loop region. Residue F486 resides in L 3 and plays a crucial role in stabilizing the recognition loop by making a π-stacking interaction with residue Y83A on ACE2. Binding Free Energies. The binding energetics between ACE2 and the RBD of SARS-COV, nCOV-2019, and all its mutant complexes were investigated by the MMPBSA method. 48 The binding energy was partitioned into its individual components including: electrostatic, van der Waals, polar solvation, and SASA to identify important factors affecting the interface of RBD and ACE2 in all complexes. nCOV-2019 has a total binding energy of −50.22 ± 1.93 kcal/ mol, whereas SARS-COV has a much higher binding energy of −18.79 ± 1.53 kcal/mol. Decomposition of binding energy to its components show that the most striking difference between nCOV-2019 and SARS-COV is the electrostatic contribution which is −746.69 ± 2.66 kcal/mol for nCOV-2019 and −600.14 ± 7.65 kcal/mol for SARS-COV. This high electrostatic contribution is compensated by a large polar solvation free energy which is 797.30 ± 3.12kcal/mol for nCOV-2019 and 659.61 ± 8.98 kcal/mol for SARS-COV. nCOV-2019 also possess a higher van der Waals (vdw) contribution (−89.93 ± 0.46 kcal/mol than SARS-COV (−70.07 ± 1.22 kcal/mol. Furthermore, the SASA contribution to binding for SARS-COV was −8.30 ± 0.15 kcal/mol and −10.58 kcal/mol for nCOV-2019. Both hydrophobic and electrostatic interactions play major roles in the higher affinity of nCOV-2019 RBD than SARS-COV RBD to ACE2. The binding free energies for nCOV-2019 and SARS-COV were decomposed into a per-residue based binding energy to find the residues that contribute strongly to the binding and are responsible for higher binding affinity of nCOV-2019 than SARS-COV (Figure 7) . Most of the residues in the RBM of Residue Q498 contributed −6.72 ± 0.18 kcal/mol and its corresponding residue in SARS-COV is a Y484 that contributed to total binding by −1.83 ± 0.06 kcal/mol. Other important residues Y505 and N501 have more negative contribution to total binding energy than their counterparts in SARS-COV residues Y491 and T487, respectively (Figure 7) . Residue D480 in SARS-COV contributed positively to binding energy by 6.2 ± 0.15 kcal/mol and the corresponding residue in nCOV-2019 which is a S494 residue lowered this positive contribution to only 1.17 ± 0.06 kcal/mol. Mutation D480A/ G appeared to be a dominant mutation in SARS-COV in 2002−2003. 51 This mutation was reported to escape neutralization by antibody 80R. 52 To investigate the effect of this point mutation on binding of SARS-COV RBD to ACE2 we performed an additional simulation and calculated the binding affinity for this mutant in SARS-COV RBD with the same approach for other mutation in this study (Figure 8) . D480A mutation showed a binding affinity of 23.46 ± 3.07 kcal/mol which is about 5 kcal/mol higher than the wild-type SARS-COV RBD. In SARS-COV, residue R426 had the highest contribution to the total binding energy (−6.27 ± 0.22 kcal/mol although the corresponding residue in nCOV-2019 is N439 with a contribution of −0.32 ± 0.02 kcal/mol. These important mutations on RBM of nCOV-2019 from SARS-COV caused RBD of nCOV-2019 to bind ACE2 with much stronger (about 30 kcal/mol) affinity. Binding free energy decomposition to its individual components for all mutants is represented in Table S2 . In all complexes, a large positive polar solvation free energy disfavors the binding and complex formation, which is compensated by a large negative electrostatic free energy of binding. All variants had similar SASA energies. The vdw free energy of binding ranged from −84.68 ± 0.68 kcal/mol for Q493A to −103.85 ± 0.66 kcal/mol for Y489A. Mutant K417A had the lowest electrostatic contribution to binding −415.67 ± 5.07 kcal/mol and mutants N439K and E484A had the highest electrostatic binding contribution of −989.80 ± 5.6 kcal/mol and −941.20 ± 3.95 kcal/mol, respectively. Most alanine substitutions exhibited similar or lower total binding affinities to nCOV-2019, however a few mutants had higher binding affinity than the wild type. Mutant Y489A had a total binding energy of −61.78 ± 2.59 kcal/mol which was about 11 kcal/mol lower than wild type binding energy. Mutants G446A, G447A, and T478I also demonstrated higher total binding affinities than nCOV-2019. Other alanine substitutions had similar or lower total binding energy than nCOV-2019. Mutant G502A has the lowest binding affinity among all the mutants with a binding energy of −24.31 ± 2.98 kcal/mol. Mutant systems K417A, L455A, T500A, and N501A are the other mutants with total binding affinities significantly lower than the wild type complex. The electrostatic component of binding contributes the most to the low binding affinities for these mutants. The contribution of RBM residues to binding with ACE2 for nCOV-2019 was mapped to the RBD structure and is shown in Figure 9B . Most natural mutants exhibited similar binding affinities compared to wild-type nCOV-2019 with a few exceptions. Mutation T478I which is one of the most frequent mutations based on the GISAID database has a binding affinity which is about 6 kcal/mol higher than that of the wild-type. S494P and A475V showed a slightly lower binding affinity than the wildtype complex. Other natural mutants showed binding affinities similar to wild-type RBD. N439K demonstrated a high electrostatic energy which is compensated by large polar solvation energy and this mutant has a total binding energy of −48.27 ± 3.07 kcal/mol which is similar to nCOV-2019. Hydrogen Bond, Salt-Bridge, and Hydrophobic Contact Analysis. Important hydrogen bonds (H-bonds) and salt bridges between nCOV-2019 RBD or SARS-COV RBD and ACE2 for the last 400 ns of trajectory are shown in Table 1 . nCOV-2019 RBD makes 10 H-bonds/1 salt bridge with ACE2, whereas SARS-COV makes only 5 H-bonds/1 salt bridge with ACE2 with more than 30% persistence. The evolution of the coronavirus from SARS-COV to nCOV-2019 has reshaped the interfacial hydrogen bonds with ACE2. G502 in nCOV-2019 has a persistent H-bond with residue K353 on ACE2. This residue was G488 in SARS-COV, which also makes the H-bond with K353 on ACE2. Q493 in nCOV-2019 makes H-bond with E35 and another H-bond with K31 on ACE2. This residue was an N479 in SARS-COV, which only makes one H-bond with K31 on ACE2. An important mutation from SARS-COV to nCOV-2019 is residue Q498, which was Y484 in SARS-COV. Q498 makes two H-bonds with residues D38 and K353 on ACE2, whereas Y484 in SARS-COV does not make any H-bonds. Importantly, a salt bridge between K417 and D30 in the nCOV-2019/ACE2 complex contributes to the total binding energy by −12.34 ± 0.23 kcal/mol. This residue is V404 in SARS-COV which is not able to make any salt-bridge and does not make H-bond with ACE2. Gao et al. 27 used a FEP approach and showed that mutation V404 to K417 lowers the binding energy of nCOV-2019 RBD to ACE2 by −2.2 ± 0.9 kcal/mol. A salt bridge between R426 on RBD and E329 on ACE2 stabilizes the complex in SARS-COV/ACE2. This residue is N439 in nCOV-2019 which is unable to make salt-bridge with ACE2 residue E329. One of the most observed mutations in nCOV-2019 according to the GISAID database is N439K which recovers some of the electrostatic interactions with ACE2 at this position. Y436 in SARS-COV and Y449 in nCOV-2019 both make H-bonds with D38 on ACE2. The unchanged T486 in SARS-COV corresponds to T500 in nCOV-2019, both of which make consistent H-bonds with ACE2 residue D355. Hydrophobic interactions also play an important role in stabilizing the RBD/ACE2 complex in nCOV-2019. An important interaction between nCOV-2019 RBD and ACE2 is the π-stacking interaction between F486 (RBD) and Y83 (ACE2). This interaction helps in stabilizing L 3 in nCOV-2019 compared to SARS-COV where this residue is L472. It was observed by Gao et al. 26 that mutation L472 to F486 in nCOV-2019 results in a net change in the binding free energy of −1.2 ± 0.2 kcal/mol. Other interfacial residues in nCOV-2019 RBD that participate in the hydrophobic interaction with ACE2 are L455, F456, Y473, A475, and Y489. It is interesting to note that all these residues except Y489 have mutated from SARS-COV. Spinello and co-workers 30 performed longtimescale (1μs) simulation of nCOV-2019/ACE2 and SARS-COV ACE2 and found that L 3 in nCOV-2019 is more stable The Journal of Physical Chemistry B pubs.acs.org/JPCB Article due to presence of the β 6 strand and existence of two H-bonds in L 3 (H-bonds G485-C488 and Q474-G476). Importantly, an amino acid insertion in L 3 makes this loop longer than L 3 in SARS-COV and enables it to act like a recognition loop and make more persistent H-bonds with ACE2. L455 in nCOV-2019 RBD is important for hydrophobic interaction with ACE2 and mutation L455A lowers the vdw contribution of binding affinity by about 5 kcal/mol. The H-bonds between RBD of nCOV-2019 and SARS-COV are shown in Figure 9A . The structural details discussed here are in agreement with other structural studies of the nCOV-2019 RBD/ACE2 complex. 4, 53 H-bond analysis was also performed for the mutant systems and the results for H-bonds with more than 40% consistency are shown in Table S3 . Few of the alanine substitutions increase the number of interfacial H-bonds between nCOV-2019 RBD and ACE2. Interestingly, the ala-substitution at Y489A increased the number of H-bonds in the wild-type complex. Mutation in some of the residues having consistent H-bonds in the wild type complex such as Q498A and Q493A, stunningly maintain the number of H-bonds in the wild-type complex. This indicates that the plasticity in the network of Hbonds in RBM of nCOV-2019 can reshape the network and strengthen other H-bonds upon mutation in these locations. However, few mutations decrease the number of H-bonds from the wild-type complex. Alanine substitution at residue G502 has a significant effect on the network of H-bonds between nCOV-2019 and SARS-COV. This residue locates at the end of L 4 loop near two other important residues Q498 and T500. This mutation breaks the H-bonds at these residues. Mutation K417A decreases the number of H-bonds to only 5 where the H-bond at residue Q498 is broken. This indicates the delicate nature of the H-bond from residue Q498 which can easily be broken upon ala-substitution at other residues. Furthermore, mutation N487 also decreases the number of Hbonds by breaking the H-bond at Q498. In this work, we preformed MD simulations to unveil the detailed molecular mechanism for the receptor binding of nCOV-2019 and compared it with SARS-COV. The role of key residues at the interface of nCOV-2019 with ACE2 was investigated by computational ala-scanning. A rigorous 500 ns MD simulation was performed for nCOV-2019, SARS-COV, and few mutants (Y449, T478I, Y489A, and S494P) as well as 300 ns MD simulation on each mutant. These simulations aid in understanding the dynamic role of RBD/ACE2 interface residues and estimating the binding free energy of these variants, which shed light on crucial residues for the RBD/ ACE2 complex stability. Moreover, numerous mutations have been identified in the RBD of different nCOV-2019 strains from all over the world not known to be critical for infection. 54 The effect of these mutations on the stability of the RBD/ ACE2 complex was investigated to shed light on their role in the viral infection of coronavirus. Changes in the RBD structure of nCOV-2019, SARS-COV, and mutants from their crystal structure were analyzed by RMSD and RMSF. nCOV-2019 showed a stable structure with a RMSD =1.5 Å, whereas SARS-COV had a larger RMSD value between 3−4 Å during the simulation. Most mutations of nCOV-2019 maintained similar stability to the wild-type. However, a few nCOV-2019 mutations resulted in larger deviations (>2 Å), i.e., Y489A, F456A, Y505A, N487A, K417A, Y473A, and Y449A. We further investigated the structure of the extended loop domain ( Figure 1B ) and discovered that nCOV-2019 is stable with an RMSD of less than 1 Å, whereas the extended loop in SARS-COV shows an RMSD of about 3 Å during simulation ( Figure S2 ). Some mutants showed high RMSD in this region. Alanine-substitution at residue N487 increased the extended loop RMSD to 2.5 Å. Other mutations that increased the extended loop RMSD (>2 Å) include Y449A, G477A, and E484A. The dynamic behavior of RBD was further investigated by analyzing the RMSF of all systems. As shown in Figure 4 , nCOV-2019 shows less fluctuation in L 3 than SARS-COV. This is due to the presence of a four-residue motif (GQTQ) in nCOV-2019 L 3 , which forces the loop to adopt a compact structure by making two H-bonds (G485-C488 and Q474-G476) and thereby reducing the fluctuations in the loop. Residues F486 and N487 play major roles in stabilizing the recognition loop by making π-stacking and Hbond interactions with residue Y83 on ACE2. Alanine substitution at N487 introduced a large RMSF to L 1 . Mutation L472 to F486 in SARS-COV was shown to favor binding by −1.2 ± 0.2 kcal/mol using FEP. 26 In addition, this mutation was shown to be among the five mutations that produce a super affinity ACE2 binder based on SARS-COV RBD. 6 Alanine mutations at residues Y449, G447, and E484 increased the motion in L 3 characterized by a large RMSF in this region. Using PCA, the aFEL for nCOV-2019 and SARS-COV demonstrated that the former occupies only one low energy state whereas the latter forms two distinct low energy basins separated by a metastable state with a barrier of about 6−7.5 kcal/mol. This confirms that the level of binding for the RBD domain is weaker in SARS-COV due to the presence of two basins. Similarly, alanine-substitution for a few residues caused the FEL to degenerate into separate multiple low energy regions ( Figure S4 ). Dominant motions in the RBD are visualized in Figure S5 using the first eigenvector of the PCA. nCOV-2019 and SARS-COV did not show any strong motion in the extended loop region. Porcupine plots of alaninemutants demonstrated that mutant N487A shows large motion in the L 1 region and Y449A, G447A, and E484A showed large motions in L 3 ( Figure S5 ). To better characterize the functional motions of RBD, DCCM for all systems are constructed and shown in Figure 6 and Figure S6 . nCOV-2019 showed a large correlation between the α 4 -L 1 -β 5 and α 5 -L 4 -β 5 regions. This correlation was stronger in SARS-COV and few mutants such as Y449A, G447A, and E484A. Another important correlation in nCOV-2019 is inside L 3 and β 6 . This correlation is stronger in nCOV-2019 than SARS-COV due to the presence of β 6 which makes the loop to adopt correlated motions. Few mutants impact the correlation in this region such as N487A. Interestingly, mutant F486A which is in L 3 and participates in binding by π-stacking interaction with Y83 on ACE2, disrupts the DCCM of wildtype nCOV-2019 and introduces strong correlation in the extended loop region as well as the core structure of RBD. The details of hydrogen bond and salt-bridge pattern in nCOV-2019 and SARS-COV to ACE2 (Table 1) respectively) are in good agreement with the binding energies calculated using the generalized Born method by Spinello et al. 30 Moreover, the patterns of H-bonds between nCOV-2019 and ACE2 were also already characterized by other groups 26, 30 which agrees with our work. An important H-bond between nCOV-2019 and ACE2 is between G502 on RBD and K353 of ACE2. G502 is in the L 4 region, which is populated by 5 Hbonds between RBD and ACE2. The contribution of this residue to the total binding energy is −2.03 ± 0.04 kcal/mol and the Ala-substitution at G502 has the highest effect on the binding energy among all the residues by lowering the total binding affinity to 24.31 ± 2.98 kcal/mol, which is the lowest among all mutants. This mutation breaks the other H-bonds in L 4 such as H-bonds from residues Q498 and T500. This residue is preserved and corresponds to residue G488 in SARS-COV, which also makes a H-bond with residue K353 on ACE2. Residue Q493 in nCOV-2019 participates in binding ACE2 by making two H-bonds with residues E35 and K31 on ACE2. Q493 corresponds to residue N479 in SARS-COV, which only makes one H-bond with residue Lys31 on ACE2. This caused Q493 to have more contribution to total binding than its counterpart N479. However, alanine substitution at Q493 did not affect the total binding energy and this mutant had a total binding energy similar to the wild-type complex as it maintains the number H-bonds in the wild-type complex. Residues Q498 and T500 in nCOV-2019 are crucial for binding by making H-bonds with ACE2 residues D38, D355, and K353. Residue Q498 corresponds to residue Y484 in SARS-COV which does not make any H-bond in the SARS-COV/ACE2 complex. Q498 contributes to binding by −6.72 ± 0.18 kcal/mol which is more than the contribution of Y484 in SARS-COV (−1.83 ± 0.06 kcal/mol). Ala-substitution at Q498 did not show large impact on the total binding energy. Residue T500 is conserved and corresponds to residue T486 which also makes a H-bond with Asp355 on ACE2. Mutation of T500 to Alanine lowers the binding affinity by about 10 kcal/mol. Residue N487 in nCOV-2019 locates in L 3 and plays a crucial role in stabilizing the recognition loop by making a Hbond with Y83 on ACE2. This residue contributes to the total binding energy of nCOV-2019 by −1.52 ± 0.06 kcal/mol, whereas its corresponding residue in SARS-COV does not show any contribution to the binding energy (−0.02 ± 0.05 kcal/mol). This demonstrates that L 3 in SARS-COV has evolved to be an important recognition loop in nCOV-2019, which participates in binding with ACE2. Residue K417 in nCOV-2019 has the most contribution to the total binding energy (−12.34 ± 0.23 kcal/mol by making a salt-bridge with residue D30 on ACE2. This residue is crucial for the binding of RBD and ACE2 and alanine substitution lowers the total binding energy to −29.56 ± 2.95 kcal/mol. This salt-bridge is found to be important for the stability of the crystal structure of the RBD/ACE2 complex in nCOV-2019. 4 K417 is Val404 in SARS-COV which does not participate in binding ACE2. Another important residue in nCOV-2019 is L455 which contributes to binding by −1.86 ± 0.03 kcal/mol. This residue is important for hydrophobic interaction with ACE2 and mutating it to alanine lowers the total binding affinity by about 17 kcal/mol. The hydrophobic residue F456 in nCOV-2019 also has a favorable contribution to the binding energy and F456A lowers the binding affinity by 5 kcal/mol. These results are in fair agreement with experimental binding measurements with deep mutational scanning of RBD in nCOV-2019 where they used flow cytometry for different ACE2 concentrations to measure the dissociation constant K D . 25 It was shown that mutations at K417, N487, T500, and G502 are detrimental for binding to ACE2, which agrees with the results here. These experiments showed that mutations at Q493 and Q498 do not impact the binding affinity of RBD to ACE2 which demonstrates the high plasticity of the network of H-bonds at the interface where upon mutation at these residues the network can reshape to form new H-bonds. Mutations at hydrophobic residues L455 and F456 are shown to reduce the binding affinity in these experiments. The total binding energy calculation of all the variants showed that mutation Y489A has the highest binding affinity among all systems which is about 11 kcal/mol stronger than that of the nCOV-2019 complex. This residue is located in β 6 , which is part of the recognition region of RBD for binding to ACE2. Removal of this bulky hydrophobic residue at the interface with ACE2 caused the extended loop to move closer to the ACE2 interface and make more H-bonds with ACE2 (Table S3) . A high electrostatic interaction energy is the reason for the higher binding energy of mutant Y489A than the wild-type complex. It is interesting to note that among the five residues L455, F456, Y473, A475, and Y489 that make hydrophobic interactions with ACE2, Y489 is the only residue that is conserved from SARS-COV. However, the experimental binding affinity measurements using deep mutational scanning showed that mutations at this position lower the binding affinity to ACE2. Other alanine substitutions that increase the binding energy are G446A and G447A. Residues G446 and G447 reside in L 1 and mutation to alanine can make L 1 take a more rigid form as shown in the RMSF plot in Figure S3 . However, experiment showed that these mutations have similar or lower binding affinities to ACE2 than the wild-type RBD and care must be taken when interpreting these results. 25 This discrepancy could be due to force-field inaccuracy and the deficiencies in the PBSA method for the treatment of solvent in binding energy calculation. Further studies are needed to investigate whether these mutations will increase the binding affinity to ACE2. Deep mutational scanning using flow cytometry is a qualitative method to measure the impact of a large number of mutations of protein−protein interactions and further experiments such as SPR or isothermal titration calorimetry which are conventional methods for measuring binding affinities needed to study the effect of these mutations in detail. Important mutations found in naturally occurring nCOV-2019 appear to influence to some extent the binding to ACE2. Mutation T478I which is one of the most frequent mutations according to GISAID database, increases the binding affinity of nCOV-2019 to ACE2 by about 6 kcal/mol. Mutation N439K has the highest occurrence among all strains of coronavirus in the GISAID database which demonstrated the highest electrostatic interaction among all studied systems. This residue corresponds to R426 in SARS-COV which exerts a salt-bridge interaction with E329 on ACE2. Mutation N439K recovers some of this ACE2 interaction; however, it exerts a binding affinity similar to that of wild-type RBD. Contribution of important interface residues to binding affinity was compared for mutations T478I, N439K, and wild-type nCOV-2019 ( Figure S7 ). The most striking differences between wild-type RBD and mutation T478I are residues Y449 and Q498 which have significantly higher contribution to binding than the wild type residue. Most other residues at the interface have similar binding affinities to the nCOV-2019. A The Journal of Physical Chemistry B pubs.acs.org/JPCB Article higher H-bond persistence is also seen for these two residues Y449 and Q498 compared to the wild type RBD which is the reason for their higher contribution to the total binding energy. Mutation N439K has a slightly lower binding affinity to ACE2 than the wild type RBD. Per residue binding energy decomposition showed that K439 in this system has a favorable contribution of −1.80 ± 0.15 kcal/mol to the total binding energy which is balanced by a lower contribution of K417 which resulted in a binding affinity similar to wild-type RBD. Mutant E484A, which is also one of the observed mutations based on GISAID database, demonstrates a high electrostatic interaction with ACE2. E484 contributes to binding by 3.56 ± 0.15 kcal/mol whereas the corresponding residue in SARS-COV, P469 contributes to binding of SARS-COV to ACE2 by −0.27 ± 0.01 kcal/mol. This residue is close to D30 on ACE2 and has electrostatic repulsion with this residue. Most natural mutants including N439K, A475V, G476S, V483A, V483F, E484A, and S494P showed similar or slightly lower binding affinities to ACE2 compared to wild-type complex which agrees with experimental binding measurements. 25 However, the experimental binding affinity for T478I also showed similar binding affinity to the wild-type complex. This difference could be due to the use of MMPBSA approach for calculation of polar solvation and further studies are needed to study the effect of this mutation on viral infectivity of coronavirus. Additional sequence differences between nCOV-2019 and SARS-COV influence RBD/ACE2 binding. Residue D480 in SARS-COV contributes negatively to total binding energy (6.25 ± 0.14 kcal/mol) and mutating this residue to S494 in nCOV-2019 lowers this negative contribution to 1.17 ± 0.06 kcal/mol. D480 in SARS-COV is located in a region of high negative charge from residues E35, E37, and D38 on ACE2. Electrostatic repulsion between D480 on SARS-COV and the acidic residues on ACE2 is the reason for highly negative contribution of this residue to binding of SARS-COV to ACE2. Mutation to S494 in this location removes this highly negative contribution. Gao and co-workers 26 computed the relative free energies of binding because of mutations from the RBD-ACE2 of SARS-COV to the corresponding residues in nCOV-2019. They used a FEP approach and showed that mutation D480S in SARS-COV changed the binding free energy by −1.9 ± 0.8 kcal/mol which is consistent with our study. Furthermore, we performed an additional simulation on D480A mutant in SARS-COV and found that this mutation has a binding affinity of 23.46 ± 3.07 kcal/mol which is about 5 kcal/mol higher than the wild-type SARS-COV RBD. In addition, experimental binding affinity measurements showed that mutations of S494 to an acidic residue highly reduce the binding affinity to ACE2 which confirms the hypothesis here. To our knowledge this is first detailed molecular simulation study on the effect of mutations on binding of nCOV-2019 to ACE2. Previous computational studies have found that nCOV-2019 binds to ACE2 with a total binding affinity which was about 30 kcal/mol stronger than SARS-COV and is in fair agreement with the results here. 56 The critical role of interface residues is computationally investigated here and in other articles and the results of all the studies indicate the importance of these residues for the stability of the complex and finding hotspot residues for the interaction with receptor ACE2. 26, 30, 55, 56 It is interesting to note the role of L 3 in the stability of the RBD/ACE2 complex. The amino acid insertions in L 3 for nCOV-2019 have converted an unessential part of RBD in SARS to a functional domain of the RBD. This loop participates in binding ACE2 by making H-bond as well as π-stacking interactions with ACE2, which makes this region to act as a recognition loop. Previous studies on SARS-COV have shown that there is a correlation between the higher binding affinity to the receptor and higher infection rate by coronavirus. 6, 13, 57 The higher binding affinity of nCOV-2019 for ACE2 than SARS-COV to ACE2 is suggested to be the reason for its high infection rate. Most natural mutations showed similar binding affinities to wild-type ACE2 which indicates that the virus was already effective at the beginning of the crisis for binding ACE2. A few mutations such as N489A and T478I are shown to increase the binding affinity to ACE2. However, more studies are needed to investigate the effect of these mutations in detail. Mutations of nCOV-2019 RBD that do not change the binding affinity and complex stability could have implications for antibody design purposes since they could act as antibody escape mutants. Escape from monoclonal antibodies is observed for mutations of SARS-COV in 2002 and these mutations should be considered for any antibody design endeavors against these escape mutations. In conclusion, this study unraveled key molecular traits underlying the higher affinity of nCOV-2019 for ACE2 compared to SARS-COV and unveiled critical residues for the interaction by in silico alanine scanning mutations and binding free energy calculations. The higher affinity of nCOV-2019 to binding with ACE2 correlates with higher human-tohuman transmissibility of nCOV-2019 compared to SARS-COV. Ala-scanning mutagenesis of the interface residues of nCOV-2019 RBM has shed light on the crucial interface residues and helped obtain an atomic-level understanding of the interaction between coronavirus and the receptor ACE2 on the host cell. MD simulations on RBD mutations found in strains of nCOV-2019 from different countries aid in the understanding of how these mutations can play an important role in viral infection with ACE2 attachment. In addition to previously reported residues, it was found that residue F486 locating in L 3 plays a crucial role in the dynamic stability of the complex by a π-stacking interaction with ACE2. Per-residue free energy decomposition pinpoints the critical role of residues K417, Y505, Q498, and Q493 in binding ACE2. Alanine scanning of interface residues in nCOV-2019 RBD showed that alanine substitution at some residues such as G502, K417, and L455 can significantly decrease the binding affinity of the complex. Moreover, mutation T478I, which is one of the most probable mutations in RBD of nCOV-2019 is found to bind ACE2 with about 7 kcal/mol higher affinity than wild-type. It is also alerting that some of the alanine substitutions at residues G446, G447, and Y489 substantially increased the binding affinity that may lead to a strong RBD attachment to ACE2 and influence the infection virulence. However, details of interaction between these mutants and ACE2 should be carefully studied using experimental techniques. On the other hand, most mutations are found not to impact the binding affinity of RBD with ACE2 in nCOV-2019 which could have implications for vaccine design endeavors as these mutations could act as antibody escape mutants. Receptor recognition is the first line of attack for coronavirus and this study gives novel insights into key structural features of interface residues for the advancement of A novel coronavirus from patients with pneumonia in China A pneumonia outbreak associated with a new coronavirus of probable bat origin Structure, function, and evolution of coronavirus spike proteins Structural and functional basis of SARS-CoV-2 entry by using human ACE2 Receptor recognition by the novel coronavirus from Wuhan: an analysis based on decade-long structural studies of SARS coronavirus Cryo-EM structure of the 2019-nCoV spike in the prefusion conformation The coronavirus spike protein is a class I virus fusion protein: structural and functional characterization of the fusion core complex Structural insights into the Middle East respiratory syndrome coronavirus 4a protein and its dsRNA binding mechanism Structure, function, and antigenicity of the SARS-CoV-2 spike glycoprotein Potent binding of 2019 novel coronavirus spike protein by a SARS coronavirus-specific human monoclonal antibody Receptor recognition mechanisms of coronaviruses: a decade of structural studies Efficient replication of severe acute respiratory syndrome coronavirus in mouse cells is limited by murine angiotensin-converting enzyme 2 Structural analysis of major species barriers between humans and palm civets for severe acute respiratory syndrome coronavirus infections Receptor recognition and cross-species infections of SARS coronavirus Mechanisms of host receptor adaptation by severe acute respiratory syndrome coronavirus The SARS coronavirus S glycoprotein receptor binding domain: fine mapping and functional characterization Development and characterization of a severe acute respiratory syndromeassociated coronavirusneutralizing human monoclonal antibody that provides effective immunoprophylaxis in mice Neutralizing human monoclonal antibodies to severe acute respiratory syndrome coronavirus: target, mechanism of action, and therapeutic potential Human monoclonal antibody combination against SARS coronavirus: synergy and coverage of escape mutants Structural basis for potent crossneutralizing human monoclonal antibody protection against lethal human and zoonotic severe acute respiratory syndrome coronavirus challenge Escape from human monoclonal antibody neutralization affects in vitro and in vivo fitness of severe acute respiratory syndrome coronavirus Data, disease and diplomacy: GISAID's innovative contribution to global health Deep mutational scanning of SARS-CoV-2 receptor binding domain reveals constraints on folding and ACE2 binding Enhanced receptor binding of SARS-CoV-2 through networks of hydrogen-bonding and hydrophobic interactions De novo design of high-affinity antibody variable regions (scFv) against the SARS-CoV-2 spike protein No. 050534. (30) Spinello, A.; Saltalamacchia, A.; Magistrato, A. Is the Rigidity of SARS-CoV-2 Spike Receptor-Binding Motif the Hallmark for Its Enhanced Infectivity? Insights from All-Atoms Simulations Structural basis for receptor recognition by the novel coronavirus from Wuhan mCSM-PPI2: predicting the effects of mutations on protein−protein interactions Structure of the SARS-CoV-2 spike receptor-binding domain bound to the ACE2 receptor The protein data bank High performance molecular simulations through multi-level parallelism from laptops to supercomputers Improved side-chain torsion potentials for the Amber ff99SB protein force field Optimization of the additive CHARMM all-atom protein force field targeting improved sampling of the backbone ϕ, ψ and side-chain χ1 and χ2 dihedral angles Canonical sampling through velocity rescaling Polymorphic transitions in single crystals: A new molecular dynamics method Constant pressure molecular dynamics simulation: the Langevin piston method Particle mesh Ewald: An N· log (N) method for Ewald sums in large systems Van Der Spoel, D. GROMACS 3.0: a package for molecular simulation and trajectory analysis. Mol. modeling ann Harmonic analysis of large systems Principal component analysis: a method for determining the essential dynamics of proteins MD-TASK: a software suite for analyzing molecular dynamics trajectories Assessing the performance of the g_mmpbsa tools to simulate the inhibition of oseltamivir to influenza virus neuraminidase by molecular mechanics Poisson−Boltzmann surface area methods Recent developments and applications of the MMPBSA method Virtual screening using molecular simulations Insight into the interactive residues between two domains of human somatic Angiotensin-converting enzyme and Angiotensin II by MM-PBSA calculation and steered molecular dynamics simulation Tracking changes in SARS-CoV-2 Spike: evidence that D614G increases infectivity of the COVID-19 virus Broadening of neutralization activity to directly block a dominant antibody-driven SARS-coronavirus evolution pathway Characterization of the receptor-binding domain (RBD) of 2019 novel coronavirus: implication for development of RBD protein as a viral attachment inhibitor and vaccine Phylogenetic network analysis of SARS-CoV-2 genomes The SARS-CoV-2 exerts a distinctive strategy for interacting with the ACE2 human receptor Computational simulations reveal the binding dynamics between human ACE2 and the receptor binding domain of SARS-CoV-2 spike protein SARS-CoV-2 Spike Glycoprotein Receptor Binding Domain is Subject to Negative Selection with Predicted Positive Selection Mutations