key: cord-1026748-funa609b authors: Luan, Binquan; Huynh, Tien title: Insights into SARS-CoV-2’s Mutations for Evading Human Antibodies: Sacrifice and Survival date: 2021-04-09 journal: J Med Chem DOI: 10.1021/acs.jmedchem.1c00311 sha: bd73a82063ede5718204a4b3afccea42059009e6 doc_id: 1026748 cord_uid: funa609b [Image: see text] The highly infectious SARS-CoV-2 variant B.1.351 that first emerged in South Africa with triple mutations (N501Y, K417N, and E484K) is globally worrisome. It is known that N501Y and E484K can enhance binding between the coronavirus receptor domain (RBD) and human ACE2. However, the K417N mutation appears to be unfavorable as it removes one interfacial salt bridge. Here, we show that despite the decrease in binding affinity (1.48 kcal/mol) between RBD and ACE2, the K417N mutation abolishes a buried interfacial salt bridge between the RBD and neutralizing antibody CB6. This substantially reduces their binding energy by 9.59 kcal/mol, thus facilitating the process by which the variant efficiently eludes CB6 (including many other antibodies). Our theoretical predictions agree with existing experimental findings. Harnessing the revealed molecular mechanisms makes it possible to redesign therapeutic antibodies, thus making them more efficacious. Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) 1 is a positive-sense single-stranded RNA virus that causes the ongoing coronavirus disease 2019 (COVID- 19) pandemic. Because SARS-CoV-2 is a highly transmissible and pathogenic novel coronavirus, the global morbidity and mortality rates of COVID-19 have been quite substantial since it was first detected at the end of 2019. 2 This has brought together researchers from the international scientific community to combat COVID-19 with extraordinary effort. Although we have come a long way in a short period of time, it is still challenging to have the pandemic under full control especially with the recent discovery of several new variants of SARS-CoV-2, such as B.1.1.7 (or 501Y.V1) from the United Kingdom, B.1.351 (or 501Y.V2) from South Africa, and P.1 from Brazil. With the increased transmissibility and hence rapid spreading of these newly identified lineages, serious concerns over whether they could undermine the currently available vaccines or natural immunity of people previously infected with COVID-19 have been raised. 3 Unfortunately, a recent experiment provided evidence that variant B.1.351 eludes natural and vaccine-induced sera. 4 Experimental studies have demonstrated that angiotensin converting enzyme 2 (ACE2), a protein expressed on the surface of human cells in various organs, plays a crucial role in the viral infection of SARS-CoV-2. 5−7 As a prelude to SARS-CoV-2's cell entry, the receptor binding domain (RBD) of the spike glycoprotein (S-protein) 8 on the virion surface binds ACE2 on a host cell. Therefore, the S-protein of SARS-CoV-2 becomes a primary target of the host immune system and is used as the leading antigen for vaccine development. Despite being detected from three different continents, the fact that the three aforementioned SARS-CoV-2 variants share some of the same mutations at the S-protein's RBD suggests that these mutations might have conferred an evolutionary advantage to the virus. While various mutations as well as deletions outside the RBD of S-protein are also important for the enhanced fitness (of variants) for entering host cells, here we focus on mutations in the RBD that harbors the binding site of ACE2 and the epitopes for neutralizing monoclonal antibodies (mAbs). For the U.K. variant (with the N501Y mutation on the RBD), experimental studies have demonstrated that the reproductive number that measures its infectiousness is approximately 0.4− 0.7 higher than those of other strains of the virus 9 and determined recently it is unlikely to escape BNT162b2 vaccinemediated protection. 10 Our previous work 11 showed that this mutation in the U.K. variant can increase the RBD's binding affinity for ACE2 but has no obvious effect on mAbs. In addition to the N501Y mutation, the South Africa and Brazil variants also contain the K417N and E484K mutations. Recent studies showed that the E484K mutation not only can enhance RBD− ACE2 binding 12 but also can help the virus escape the therapeutically relevant mAbs. 13 However, reflected in the paucity of existing work, the significance of the K417N mutation so far is still elusive, which prevents us from fully understanding the infection mechanism of these new variants. As shown in Figure 1a , in the wild-type RBD, K417 forms a salt bridge with D30 in ACE2 (human). Consistently, animals such as pangolin, sheep, Siberian tiger, horseshoe bat, hamster, cat, horse, cow, and ferret having either D30 or E30 (after sequence alignment) in their ACE2 are susceptible to SARS-CoV-2, while other animals such as chicken (A30), rat (N30), and duck (A30) are not. 14 This suggests that the salt bridge formed by K417 (in the RBD) and D30/E30 (in ACE2) is important for RBD−ACE2 recognition. Thus, the K417N mutation resulting in the abolishment of this favorable interfacial interaction (i.e., reducing the RBD−ACE2 binding affinity as verified in experiment 12 ) is highly unintuitive. Here, we are motivated to investigate the molecular mechanism of the K417N mutation, to unveil its key benefits for the virus to evolve through this path. Complementary to experimental efforts, all-atom molecular dynamics (MD) simulations with sophisticated and wellcalibrated force fields have been widely used to image nanoscale events and investigate the molecular mechanism of proteins. 15−18 In this work, we conducted a computational analysis of the K417N mutation in the South Africa variant using MD simulations with explicit solvent, aiming to gain a better understanding of its underlying molecular mechanism. In addition to the RBD−ACE2 interaction, we also investigated the RBD's interaction with mAbs. In particular, we explored the mAb CB6 that recognizes an epitope site in the RBD overlapping the binding site of ACE2 and investigated its binding competition with ACE2 over RBD. Our results may help provide invaluable insights into why K417N has been selected during viral evolution and inspire a better design of more efficacious mAbs for treating COVID-19 patients infected with the new SARS-CoV-2 variants. ■ RESULTS Figure 1b illustrates the simulation system for modeling the interaction between ACE2 and RBD (see the Experimental Section for detailed simulation protocols). Briefly, taken from the crystallographic structure (PDB entry 6M0J), the structure of the RBD−ACE2 complex was solvated in a 0.15 M NaCl electrolyte. Similarly, as shown in Figure 1c , we built the simulation system for the complex of RBD and one Fab of human antibody CB6 with their atomic coordinates taken from the crystallographic structure (PDB entry 7C01). Hereafter, we simply refer to the Fab as CB6. During the ∼350 ns MD simulation, the RBD−ACE2 complex starting from the structure in the crystal environment was equilibrated in the physiology-like environment (a 0.15 M electrolyte). Figure S1 shows the root-mean-square deviation (RMSD) of protein-backbone atoms in the RBD bound with either ACE2 or CB6. After ∼40 ns, the RMSD values saturated at around 1.7 Å for the RBD bound on CB6 and 1.5 Å for the RBD bound on ACE2. These small RMSD values suggest that the secondary structure of the modeled RBD (with four internal disulfide bonds) was very stable and properly equilibrated. We also calculated the interfacial contact areas for the complexes using the solvent accessible surface area method. 19 On average, the contact area is ∼8.6 nm 2 between the RBD and ACE2 while that for the RBD−CB6 contact is larger and is ∼10.3 nm 2 , as shown in Figure 2a . During the entire simulation time, values of contact areas fluctuated around the mean, indicating that these two complex structures were stable and equilibrated in the electrolyte. By examining the simulation trajectories, we found that K417 indeed played an important role in stabilizing the RBD−ACE2 and RBD−CB6 complexes. In Figure 2b , we highlight the timedependent number N of heavy (or non-hydrogen) atoms in ACE2 or CB6 that were within 5 Å of K417 in the RBD. Notably, K417 in RBD interacts with many more atoms in CB6 than in ACE2. The average saturated number N for the RBD−CB6 complex is ∼23 atoms (orange line in Figure 2b ), while that for the RBD−ACE2 complex is only ∼6 atoms (blue line in Figure 2b ). Additionally, the number N for the RBD−ACE2 complex fluctuates much more than that for the RBD−CB6 complex, which indicates that the interaction in the latter complex is more stable. Due to the K417−D104 salt bridge, K417 in the RBD stably contacted nearby residues (including D104) in the heavy chain. However, the distance between K417 in the RBD and residues (such as T94 and P95) in the light chain increased slightly to accommodate water molecules that entered into the hydrophilic interfacial area, resulting in a decrease in the contact number N at the beginning of the MD simulation (see the orange line in Figure 2b ). Among those residues (in ACE2 or CB6) in contact with K417 in the RBD, there are two key salt bridges. One is formed by K417 in the RBD and D30 in ACE2 18 (as shown in Figure 1a ). We illustrate in Figure 2c that this salt bridge is on the surface of the RBD−ACE2 complex. Thus, this salt bridge is Proteins are shown in cartoon representation, with RBD colored gray, ACE2 colored blue, and the heavy (V H and C H ) and light (V L and C L ) chains of the CB6's Fab colored purple and green, respectively. Na + and Cl − are shown as yellow and cyan spheres, respectively. Water is shown transparently. exposed to water and constantly disrupted by polar water molecules, accounting for the observed fluctuations in N (blue line in Figure 2b ). The other one is formed by K417 in the RBD and D104 in CB6 20 as shown in Figure 2d . Remarkably, this salt bridge is buried among the heterotrimer composed of the RBD and two variable domains of heavy chain V H and light chain V L of CB6. This salt bridge buried inside the protein complex was very stable during the simulation, consistent with the nearly constant number N in Figure 2b (orange line). To further demonstrate the stability of these two salt bridges, we define distance D between the NZ atom in the lysine (K) and the CG atom in the aspartate (D), as shown in the inset of Figure 3a . The histograms in Figure 3a show the probability distributions of distance D. For the salt bridge in the RBD− CB6 complex, the most probable distance in the single-peak distribution is ∼3.2 Å, corroborating the idea that this buried salt bridge is very stable. For the salt bridge in the RBD−ACE2 complex, the most probable distance is in the first peak and is ∼3.5 Å. This larger distance (as compared with that for the buried salt bridge) reflects the weakened electrostatic interaction inside the water-exposed salt bridge. Furthermore, an additional peak at larger distances (∼5.5 Å) is an indication of a state of the broken salt bridge. Indeed, we observed in the simulation trajectory that this salt bridge is susceptible to polar water molecules and can break and re-form from time to time. Therefore, in the wild-type RBD, K417 interacts more strongly with D104 in CB6 than with D30 in ACE2. The calculations of the binding free energy change induced by the K417N mutation were accomplished using the free energy perturbation (FEP) method. 21 As required in FEP calculations, we performed 177 ns MD simulations of the RBD in a 0.15 M NaCl electrolyte (a free state), as shown in Figure S2 . With protein structures for both bound (or complex) and free states in respective MD simulations, we employed the FEP alchemy method to calculate the binding free energy difference for the K417N mutation on the RBD, using the thermodynamic cycle shown in panels b−e of Figure 3 . By definition, binding free energy changes for the RBD with either ACE2 or CB6 (due to the K417N mutation) can be obtained as ΔΔG = ΔG 2 − ΔG 1 (see Figure 3 ). In practice, it is not easy to directly calculate ΔG 1 and ΔG 2 , which can be circumvented by computing ΔG A and ΔG B instead using the thermodynamic cycle (see Figure 3 ). Therefore, ΔΔG = ΔG A − ΔG B . Through the ensemble average, 21 ΔG A and ΔG B can be calculated theoretically (see the Experimental Section) as where k B is the Boltzmann constant, T is the temperature, and H i and H f are the Hamiltonians for the initial (i) stage with K417 (see Figure 3b ,d) and the final (f) stage with N417 (see Figure 3c ,e), respectively. The results from the FEP calculations are summarized in Table 1 . In the free state (Figure 3d,e) , the K417N mutation yielded a free energy change ΔG B of −32.08 kcal/mol. In the (Figure 3b,c) , free energy change ΔG A = −22.49 kcal/mol. Altogether, the obtained value of ΔΔG is 9.59 kcal/mol, suggesting that the K417N mutation significantly reduced the binding affinity between the RBD and CB6. Similarly, for the RBD−ACE2 complex, we have a ΔG A of −30.60 kcal/mol, and consequently, ΔΔG is 1.48 kcal/mol. Thus, the K417N mutation also reduced the binding affinity between the RBD and ACE2, but the reduction is ∼6.5 times less than that between the RBD and CB6. Previously, it was demonstrated in an experiment that the K417N mutation weakened the binding affinity between the RBD and ACE2, 12 which is consistent with our simulation results. The noticeably small errors listed in Table 1 manifest the accuracy and convergence afforded by the FEP methodology. The free energy loss due to the K417N mutation is mainly caused by the change in electrostatic energy. Dielectric constant ε p inside a protein is approximately 4, 22 while dielectric constant ε w of water is ∼90 for the TIP3P model. 23 If we approximate the protein−water interface with a planar surface separating two dielectric materials (ε p and ε w ), the effective dielectric constant on the interface can be calculated as (ε p + ε w )/2, from a simple continuum electrostatics calculation. Thus, the ratio γ of electrostatic energy changes for salt bridges on the protein surface and in water can be roughly estimated as (ε p + ε w )/2ε w ∼ 11.8 (larger than our simulation result γ = 6.5). However, the more realistic explicit model with atomic details 24 predicts that the average dielectric constant on the protein−water interface is approximately 20−30, yielding γ values of ∼5−7.5, which are in excellent agreement with our numerical result. More importantly, the phenomenon of weakening the binding affinity of the RBD−CB6 complex observed in the K417N mutation is not unique. After comparison of some crystal structures for the complex of the RBD (of the SARS-CoV-2's spike protein) and human antibodies, surprisingly we noticed that K417 in the RBD can form a buried salt bridge with either a glutamate or an aspartate in four other human antibodies as highlighted in IGHV-53 (Figure 4a ), C1AB3 (Figure 4b ), BD-236 (Figure 4c) , and CC12.1 (Figure 4d ). Thus, it is likely that all of these antibodies may not be able to neutralize the 501Y.V2 variant. It is worth noting that during the preparation of this paper, two experimental preprints posted on bioRxiv demonstrating that the South Africa variant can evade human antibodies CB6 (Figure 1c , also known as LyCoV016) and CC12.1 (Figure 4d ) that are known to neutralize wild-type SARS-CoV-2. 13, 25 These experimental results further confirm our predictions. Because K417 in RBD can form a buried salt bridge with a glutamate/aspartate in a class of human antibodies, we speculate that K417 is the Achilles' heel of the virus and can be easily targeted by many human antibodies. Therefore, the K417N mutation could be evidence of viral adaptation to the human immune system or natural selection under the "pressure" of human neutralizing mAbs. Although it is counterintuitive to learn that the K417N mutation actually decreases the RBD's binding affinity with ACE2 by 1.48 kcal/mol, it is in fact beneficial for RBD to escape human antibody CB6 by dramatically reducing the binding affinity by 9.59 kcal/mol. Because of the loss of binding affinity, the variant with this single mutation is less competitive than the wild-type SARS-CoV-2 virus, in terms of the ability to enter host cells. So far, there is no evidence that this particular mutation occurs alone in any known variants of SARS-CoV-2. In the South Africa variant (501Y.V2), actually there are two more mutations, N501Y and E484K. To explore their contributions, we carried out MD simulation for the RBD−ACE2 complex with all three mutations (K417N, N501Y, and E484K), as shown in Figure 5 . During the entire 340 ns MD simulation with the triple mutations (Figure 5a ), Y501 interacted stably with residues Y41 and K353 in ACE2 (see Figure 5b and Movie S1). Namely, Y501 in the RBD formed a T-shaped contact with Y41 in ACE2, and simultaneously, Y501 formed a hydrophobic interaction with the alkane chain of the amphipathic K353. We did not observe any coordination formed by N417 in the RBD and any residue in ACE2. Before the E484K mutation, E484 forms a hydrogen bond with F490 ( Figure S5) , which attaches the loop (residues 477−486) to the remainder of the RBD. After the mutation, this hydrogen bond was abolished and the loop was able to detach from the remainder of the RBD. Interestingly, K484 was gradually pulled toward E75 in ACE2 due to the electrostatic attraction (see Figure 5b and Movie S1). In Figure 5c , we show that the initial distance D between K484 in the RBD and E75 in ACE2 was on average ∼11 Å. After ∼190 ns, the distance was decreased to ∼3.5 Å, suggesting the formation of a salt bridge. Due to the exposure to water, this salt bridge can break and reform frequently. This new state persisted for the remaining 150 ns (Figure 5c ). These dynamic interfacial interactions are illustrated in Movie S1, as well. Therefore, the gain from the improved RBD−ACE2 binding that resulted from both N501Y and E484K mutations might be enough to compensate for the loss caused by the K417N mutation and is likely to yield an interaction with ACE2 even stronger than that of the wild-type virus. In summary, we have shown that K417 plays an important role in the binding between the RBD and ACE2 or CB6 by forming interfacial salt bridges. The salt bridge between K417 in the RBD and D104 in CB6 is buried inside the complex and is therefore much more stable than the water-exposed one formed between K417 in the RBD and D30 in ACE2. Thus, the K417N mutation can weaken the binding of the RBD with CB6 much more than with ACE2. More dramatically, the K417N mutation allows the variant to escape from many human antibodies other than CB6 by removing a salt bridge buried in the RBD−antibody interface. Interestingly, the virus with the K417N mutation seems to sacrifice its binding affinity for ACE2 to survive the attack of the antibodies. In fact, after the sequence alignment, the residue in SARS-CoV corresponding to K417 in SARS-CoV-2 is valine (V404). 26 On the contrary, a recent study has demonstrated that the human immune systems may also quickly respond to these new variants and produce corresponding neutralizing antibodies that contain more somatic hypermutation, increased potency, and resistance to RBD mutations. 27 The continued evolution of the humoral response to recent SARS-CoV-2 variants warrants further studies. Additionally, with the revealed mechanism underlying the K417N mutation, it is possible to design more efficacious antibody cocktails to treat COVID-19 patients infected with variant 501Y.V2 as well as the recently discovered variant P.1 in Brazil. ■ EXPERIMENTAL SECTION MD Simulations. We carried out all-atom MD simulations for the bound (the RBD−ACE2 complex) and free (stand-alone RBD) states using the NAMD2.13 package 28 running on the IBM Power Cluster. To model the RBD−ACE2 and RBD−CB6 complexes, we first obtained the crystal structures with PDB entries 6M0J 29 and 7C01, 30 respectively, and then solvated them in a rectangular water box that measured ∼80 Å × ∼80 Å × ∼135 Å. Bound Zn 2+ ions were included in the RBD−ACE2 complex. Na + and Cl − were added to neutralize the entire simulation system and set the ion concentration to 0.15 M (Figure 1b,c) . Each built system was first minimized for 10 ps and further equilibrated for 1000 ps in the NPT ensemble (P ∼ 1 bar, and T ∼ 300 K), with atoms in the backbones harmonically constrained (spring constant k = 1 kcal mol −1 Å −2 ). The production run was performed in the NVT ensemble, and only atoms in the backbones of ACE2 that are far from the RBD (residues 110−290, 430−510, and 580−615 in ACE2) were constrained, preventing the whole complex from rotating out of the water box. The similar approach was applied in the production run for the RBD−CB6 complex, with atoms in the backbones of the C H (residues 122−218) and C L (residues 109−215) domains constrained. We also performed MD simulation for the RBD alone in the 0.15 M NaCl electrolyte (a free state) using the same protocol ( Figure S1 ). To simulate the binding between the RBD of the South Africa variant and ACE2, we obtained the crystal structure (PDB entry 6M0J) for the RBD−ACE2 complex and introduced three mutations (N501Y, K417N, and E484K) when preparing the .psf and .pdb files for the complex. We retained all glycosylations in both the RBD and ACE2 as well as one Zn 2+ bound on ACE2 resolved in the crystal structure. The complex was solvated in a water box that measures ∼131.8 Å × ∼131.8 Å × ∼131.8 Å; 231 Na + and 208 Cl − ions were added to neutralize the whole system and set the ion concentration to 0.15 M. The final simulation system (Figure 5a ) comprises 232337 atoms. After equilibrating the structures in bound and free states, we carried out free energy perturbation (FEP) calculations. 21 In the perturbation method, many intermediate stages (denoted by λ) whose Hamiltonian H(λ) = λH f + (1 − λ)H i are inserted between the initial and final states to yield a high accuracy. With the softcore potential enabled, λ in each FEP calculation for ΔG A or ΔG B varies from 0 to 1.0 (forward) in 20 perturbation windows (lasting 300 ps in each window), yielding gradual (and progressive) annihilation and exnihilation processes for K417 and N417, respectively. We tested the convergence of the forward FEP runs with the backward ones and with the forward FEP runs with a longer simulation time of 1.5 ns per window; the obtained results are similar ( Figure S3 ). The same approach has been applied in our previous work in which the obtained FEP results are consistent with experimental ones. 31 We followed our previous FEP protocol 20 to obtain the means and errors for ΔG A and ΔG B from the lowest five energy changes among ten data in total. In Figure S4 , we also show that the means and errors for ΔG A and ΔG B obtained from all ten data are similar. In FEP runs for the K417N mutation, the net charge of the MD system changed from 0 to −1 e (where e is the elementary charge). It is important to have similar sizes of the simulation systems for the free and bound states, 32, 33 so that the energy shifts from the Ewald summation approximately cancel out when calculating ΔΔG. We applied the CHARMM36m force field 34 for proteins, the TIP3P model 35, 36 for water, and the standard force field 37 for Na + , Zn 2+ , and Cl − . The periodic boundary conditions (PBC) were used in all three dimensions. Long-range Coulomb interactions were calculated using particle-mesh Ewald (PME) full electrostatics with a grid size of ∼1 Å in Figure 5 . MD simulation of the RBD−ACE2 complex with the N501Y, K417N, and E484K mutations in the RBD. (a) Simulation system. ACE2 (blue) and the RBD (gray) are shown in cartoon representation. The glycosylated residues are not shown. Na + and Cl − ions are shown as yellow and cyan spheres, respectively. Water is shown transparently. (b) Equilibrated atomic coordinations around the three mutated residues. Protein residues are shown in stick representation. (c) Time-dependent distance D between E75 in ACE2 and K484 in the RBD. D is defined as the distance between the CD atom in E75 and the NZ atom in K484. Journal of Medicinal Chemistry pubs.acs.org/jmc Article each dimension. The pairwise van der Waals (vdW) energies were calculated using a smooth (10−12 Å) cutoff. The temperature T was maintained at 300 K by applying the Langevin thermostat, 38 while the pressure was kept constant at 1 bar using the Nose−Hoover method. 39 While the SETTLE algorithm 40 SARS-CoV-2, severe acute respiratory syndrome coronavirus 2 COVID-19, coronavirus disease 2019; ACE2, angiotensin converting enzyme 2; mAb Protein Data Bank; FEP, free energy perturbation; TIP3P, transferable intermolecular potential 3 point; PME, particle-mesh Ewald; vdW, van der Waals ■ REFERENCES Coronaviridae Study Group of the International Committee on Taxonomy of Viruses. The Species Severe Acute Respiratory Syndrome-Related Coronavirus: Classifying 2019-nCoV and Naming it SARS-CoV-2 Comparison of the Characteristics, Morbidity, and Mortality of COVID-19 and Seasonal Influenza: a Nationwide, Population-Based Retrospective Cohort Study Could New COVID Variants Undermine Vaccines? Labs Scramble to Find Out SARS-CoV-2 Cell Entry Depends on ACE2 and TMPRSS2 and Is Blocked by a Clinically Proven Protease Inhibitor Structural Basis for the Recognition of SARS-CoV-2 by Full-Length Human ACE2 Structural and Functional Properties of SARS-CoV-2 Spike Protein: Potential Antivirus Drug Development for COVID-19 Veesler, D. Structure, Function, and Antigenicity of the SARS-CoV-2 Spike Glycoprotein Sahin, U. Neutralization of SARS-CoV-2 Lineage B. 1.1. 7 Pseudovirus by BNT162b2 Vaccine−Elicited Human Sera. Science Molecular Mechanism of the N501Y Mutation for Enhanced Binding between SARS-CoV-2as Spike Protein and Human ACE2 Receptor Deep Mutational Scanning of SARS-CoV-2 Receptor Binding Domain Reveals Constraints on Folding and ACE2 Binding SARS-CoV-2 501Y. V2 Escapes Neutralization by South African COVID-19 Donor Plasma Insights on Cross-Species Transmission of SARS-CoV-2 from Structural Modeling Molecular Dynamics Simulations of Biomolecules Challenges in Protein-Folding Simulations in silico Exploration of Molecular Mechanism of Clinically Oriented Drugs for Inhibiting SARS-CoV-2's Main Protease Enhanced Receptor Binding of SARS-CoV-2 through Networks of Hydrogen-Bonding and Hydrophobic Interactions Solvent-Accessible Surfaces of Proteins and Nucleic Acids Silico Antibody Mutagenesis for Optimizing Its Binding to Spike Protein of Severe Acute Respiratory Syndrome Coronavirus 2 Free Energy Calculations Calculation of the Total Electrostatic Energy of a Macromolecular System: Solvation Energies, Binding Energies, and Conformational Analysis A Simple Polarizable Model of Water Based on Classical Drude Oscillators On the Dielectric Constant of Proteins: Smooth Dielectric Function for Macromolecular Modeling and its Implementation in DelPhi 1.351 and B.1.1.7 Identifying SARS-CoV-2-related Coronaviruses in Malayan Pangolins Scalable Molecular Dynamics with NAMD Structure of the SARS-CoV-2 Spike Receptor-Binding Domain Bound to the ACE2 Receptor A Human Neutralizing Antibody Targets the Receptor Binding Site of SARS-CoV-2 Combined Computational-Experimental Approach to Explore the Molecular Mechanism of SaCas9 with a Broadened DNA Targeting Range Mechanism of Divalent-Ion-Induced Charge Inversion of Bacterial Membranes Free Energy of Ionic Hydration CHARMM36m: an Improved Force Field for Folded and Intrinsically Disordered Proteins Comparison of Simple Potential Functions for Simulating Liquid Water Simulation of Activation Free Energies in Molecular Systems Finite Representation of an Infinite Bulk System: Solvent Boundary Potential for Computer Simulations Constant Pressure Molecular DynamicsAlgorithms An Analytical Version of the SHAKE and RATTLE Algorithm for Rigid Water Molecules Reversible Multiple Time Scale Molecular Dynamics