Open Access Article
This Open Access Article is licensed under a
Creative Commons Attribution 3.0 Unported Licence

Thermodynamics and kinetics in antibody resistance of the 501Y.V2 SARS-CoV-2 variant

Son Tung Ngo*ab, Trung Hai Nguyenab, Duc-Hung Phamc, Nguyen Thanh Tungde and Pham Cam Nam*f
aLaboratory of Theoretical and Computational Biophysics, Ton Duc Thang University, Ho Chi Minh City, Vietnam. E-mail: ngosontung@tdtu.edu.vn
bFaculty of Applied Sciences, Ton Duc Thang University, Ho Chi Minh City, Vietnam
cDivision of Immunobiology, Cincinnati Children's Hospital Medical Center, Cincinnati 45229, OH, USA
dInstitute of Materials Science, Vietnam Academy of Science and Technology, Hanoi, Vietnam
eGraduate University of Science and Technology, Vietnam Academy of Science and Technology, Hanoi, Vietnam
fDepartment of Chemical Engineering, The University of Da Nang, University of Science and Technology, Da Nang City, Vietnam. E-mail: pcnam@dut.udn.vn

Received 27th May 2021 , Accepted 6th October 2021

First published on 13th October 2021


Abstract

Understanding the thermodynamics and kinetics of the binding process of an antibody to the SARS-CoV-2 receptor-binding domain (RBD) of the spike protein is very important for the development of COVID-19 vaccines. In particular, it is essential to understand how the binding mechanism may change under the effects of RBD mutations. In this context, we have demonstrated that the South African variant (B1.351 or 501Y.V2) can resist the neutralizing antibody (NAb). Three substitutions in the RBD including K417N, E484K, and N501Y alter the free energy landscape, binding pose, binding free energy, binding kinetics, hydrogen bonding, nonbonded contacts, and unbinding pathway of RBD + NAb complexes. The low binding affinity of NAb to 501Y.V2 RBD confirms the antibody resistance of the South African variant. Moreover, the fragment of NAb + RBD can be used as an affordable model to investigate changes in the binding process between the mutated RBD and antibodies.


Introduction

The novel β-coronavirus, SARS-CoV-2, whose sequence is similar to SARS-CoV-1 and MERS-CoV, which induced the human respiratory epidemic at the beginning of this century, is the cause of the human respiratory disease (COVID-19) pandemic worldwide.1,2 This virus has infected more than 160 million people and is associated with more than 3 million deaths.3 SARS-CoV-2 is a single-positive-strand RNA virus, whose genome encodes for four main components: spike, envelope, membrane, and nucleocapsid.4,5 The spike protein (S protein) of SARS-CoV-2 which is used by the virus to bind to human angiotensin-converting-enzyme 2 (ACE2), has been researched thoroughly. ACE2 is present in different tissues in the body, including the lung, heart and liver,6 and is employed by SARS-CoV-2 as a receptor to bind and infect human cells. The S trimer comprises three copies of S1 and S2 subunits. The S1 subunit contains 4 domains: S1A, S1B, S1C, and S1D, in which the S1B domain is also called the receptor-binding domain (RBD), which mediates the attachment of the spike protein to the target cell via binding to the ACE2 receptor.7 Once the RBD is in the ‘up’ conformation, it can recognize and bind to ACE2, which leads to conformational changes of the S2 subunit and enables SARS-CoV-2 to fuse with the cell membrane and to enter host cells.1,7

RBD is the main target of neutralizing antibodies (NAbs) which can be isolated from plasma of COVID-19 patients, immunoglobulin libraries, or immunized laboratory animal models.1 These NAbs can be roughly divided into four main classes, of which class 1s′ and class 2s′ RBD epitopes overlap with the ACE2-binding site, suggesting a neutralization mechanism that involves direct competition with ACE2. Class 1 antibodies, which are encoded by the immunoglobulin V-gene (VH3-53) segment with complementarity-determining regions 1 and 2 (CDRH1 and CDRH2) and a short CDRH3, are mostly elicited by SARS-CoV-2 infection. On the other hand, when class 2 antibodies also target site I10,15 which is also target epitopes of class 1 antibodies, they bind to RBD in both ‘up’ and down’ conformations of S protein.1,8 Additionally, class 3 antibodies bind outside ACE2 and recognize both up and down RBD, while class 4 antibodies comprise previously described antibodies that cannot block ACE2 and target only to RBD in ‘up’ conformation.1 Besides RBD, the N-terminal domain (NTD) of protein S is also a popular target for NAbs and many potent monoclonal antibodies directed against this region show great potential in clinical trials for COVID-19 treatment.8 The majority of these antibodies target a single immunodominant site on NTD, including the N1-loop (NTD N-terminus), N3-loop (supersite b-hairpin), and N5 loop (supersite loop). Subsets of these antibodies and NAbs in class 1 and class 3 form multi-donor classes, with a different set of VH germline restricted mode of spike recognition.8

Due to many reasons, including high transmissibility, the longevity of the pandemic, and encountering with immunocompromised hosts, SARS-CoV-2 undergoes different rounds of mutations, which has altered the structures of the virus, modulated its infectivity, and changed the antigenicity of the surface proteins.9 The variants, including United Kingdom (B1.1.7) and South African (B1.351 or 501Y.V2) variants have associated with increased transmissibility and possibly increased mortality.8 Especially, the SARS-CoV-2 lineage in South Africa, included nine mutations in the spike protein, seems to decrease the efficacy of NAb as well as Covid-19 vaccine efficacy of some vaccines currently being used.10,11 The mutations in B1.351 can be divided into two groups, one concentrates in NTD, including four substitutions and a deletion (L18F, D80A, D215G, Δ242–244, and R246I), and the other involves three substitutions in RBD (K417N, E484K, and N501Y).12 These changes induce S protein biological and structural alterations. Especially, mutation E484K is very critical, which can reduce the effect of NAb.13

Evaluating antibody resistance of the 501Y.V2 SARS-CoV-2 variant is greatly attractive to scientists.8,10,11,14 Understanding the physical insights into the process probably enhances the vaccine developments, but the knowledge is still limited. Therefore, in this context, atomistic simulations were carried out to reveal the insights at the atomic level of the binding process of NAb to 501Y.V2 SARS-CoV-2 RBD. For the first step, structural changes of the 501Y.V2 and wildtype (WT) SARS-CoV-2 RBD + NAb complexes were characterized via unbiased MD simulations. Thermodynamics and kinetics of the binding process were then revealed via biased MD simulations. Moreover, fragment of NAb (fNAb) is often used to study the binding of S protein/RBD to antibody.15,16 In this work we also investigated binding of fNAb (cf. Fig. 1E and S1 of the ESI file) to RBD to evaluate the possibility of using fNAb in studying the influence of RBD mutations on the binding affinity instead of using NAb which is a larger molecule and costs more computing resources. Furthermore, it should be noted that glycosylation of RBD was neglected to clarify the interaction nature between RBD + antibodies, although glycans play an important role in the modulation of the spike conformational dynamics.17–19 Details of simulations were described in Fig. 1 and the ESI file.


image file: d1ra04134g-f1.tif
Fig. 1 Starting structures of RBD + antibody systems. (A) 501Y.V2 RBD + NAb in MD simulations; (B) 501Y.V2 RBD + fNAb in MD simulations; (C) + (D) 501Y.V2 RBD + NAb/fNAb in biased MD simulations; (E) free energy scheme. The WT RBD + antibody complexes were formed an initial conformation similar to the 501Y.V2 one. The fNAb was mobilized from bound to unbound states via the fast pulling of ligand (FPL) calculations, then the free energy profile was calculated via US simulations.

Computational methods

Structure of SARS-CoV-2 antibodies and RBD

The three-dimensional structure of SARS-CoV-2 RBD and their antibodies NAb was found from the Protein Data Bank (PDB) with the identities of 7BWJ.20 The fNAb was extracted from NAb as showed in Fig. 1E. The resolution of 7BWJ and 2.85 Å. Moreover, as mentioned above, it should be noted that the new variant of the SARS-CoV-2 in South Africa, B1.351 or 501Y.V2, forms eight changes in the spike protein. There are four substitutions and a deletion in the N-terminal domain (NTD) including L18F, D80A, D215G, Δ242–244, and R246I. Consequently, three substitutions were found in RBD involving K417N, E484K, and N501Y. The RBD structure with three substitutions was thus prepared via changing three residues of 7BWJ using the PyMOL mutagen tools.21 Besides, 50Y.V2 RBD with glycan was obtained from PDB ID 7LYQ,22 in which glycan linked with the residue Asn343.

Molecular dynamics simulations

The atomistic simulation was performed using the GROMACS version 5.1.5 with general-purpose computing on graphics processing units.23 The protein, antibody, and neutralized ions were parameterized via the Amber99SB-iLDN force field24 since it is a suitable force field for free energy calculation.25 The TIP3P water model was chosen to simulate the water molecule.26 The system's configurations were shown and reported in detail in Fig. 1A and B and Table S1 of the ESI.

The MD simulation parameters were taken from the previous works.27,28 However, in particular, the integration time step was taken to be 3 femtoseconds. A non-bonded pair interaction were cutoff at a radius of 0.9 nm, in which the electrostatic interaction was calculated using the fast Particle-Mesh Ewald electrostatics approach29 as well as the van der Waals (vdW) interaction was computed using the cut-off scheme. The solvated complex was initially optimized using the energy minimization via the steepest descent method. The minimized system was then relaxed in NVT and NPT ensembles with a length of 100 ps each simulation. During NVT and NPT simulations, the integral was attempted every 1 femtosecond. The equilibrium snapshots obtained via NPT simulations were used as starting conformations of MD simulations. The conventional MD simulations were performed with interval 100 ns and repeated 4 times independently. These independent trajectories thus have the same initial conformation but different generated velocities.

Biased molecular dynamics simulations

Steered-MD (SMD) simulation. Representative structures of RBD + antibody systems, which were obtained via MD simulations, were employed as initial structures of FPL simulations. The complexes were reinserted into the rectangular PBC box for saving the computing resources. The configuration information was described in Fig. 1C, D and Table S1. The FPL simulations were carried out to generate unbinding conformations of the systems, which were used as starting shapes of US simulations. From the beginning, the antibodies were forced to dissociate from the binding mode with the WT/501Y.V2 RBD using SMD simulations. In particular, eight SMD trajectories were carried out to probe the most optimal-unbinding pathway. The trajectory, in which the rupture force, FMax, and pulling work, W, formed the smallest deviation in comparison with the median values, was used for generating US windows. In FPL, the antibody was pulled along Z-axis via an external force using cantilever k = 1000 kJ mol−1 nm−2 and constant velocity v = 0.001 nm ps−1. During the simulation, the RBD was softly fixed via Cα restraint. The pulling force, antibody displacement, and systemic coordinates were recorded every 33 integrated steps.
Umbrella sampling simulation. The systemic snapshots, which were extracted from the FPL trajectory since the antibody displaced every ca. 1.0 Å along the unbinding pathway ξ, were used as starting shapes of US simulations. Ca. 25 US windows each complex were simulated with a length of 10 ns of MD simulation to calculate the potential of mean force (PMF) curve. It should be noted that a short NPT simulation was executed to reduce initial fluctuations.30,31 The PMF values were calculated via the weighted histogram analysis method (WHAM).32 The free energy barriers, ΔG++on and ΔG++off, and binding free energy, ΔGb, of the binding process between RBD and NAb were estimated as described as Fig. 1E.

Analyzed tools

The free energy landscape (FEL) of the complex was constructed using the principal component analysis (PCA) method,33 in which coordinates PC1, first eigenvector, and PC2, second eigenvector, were calculated using GROMACS tools “gmx anaeig”. In particular, the PCs of backbone protein were calculated over the entire conformational ensemble. A non-bonded (NB) contact was counted when the pair between two heavy atoms is smaller than 4.5 Å. A hydrogen bond (HB) contact was counted when the angle ∠ between donor (D) – hydrogen (H) – acceptor (A) is larger than 135° and the distance between D and A is smaller than 3.5 Å. The PMF value was estimated via the Weighted Histogram Analysis Method (WHAM) with the execution of auto-correlated time. The computed error was calculated using the bootstrapping method.34 The solvent accessible surface area (SASA) was computed using GROMACS tool “gmx sasa”.

Results and discussion

It should be noted that investigating structures of protein–protein complexes and understanding how they bind together are fundamental issues.35 Moreover, structures of several complexes remain difficult to solve experimentally.36,37 Furthermore, in order to characterize the protein–protein binding mechanisms, powerful experimental approaches are required,38,39 but the obtained data are normally limited or indirect. Obtaining direct data at an atomic level about binding pathways and physical insights into the binding mechanisms are still open issues.35 Atomistic MD simulations emerge as potential approaches for investigating both dynamics and structural change of protein–protein complexes.17,40,41 Using MD simulations, we can easily monitor the associate and dissociate processes of a monomer to the others.28,42 However, in fact, the association of two proteins may take place at a much longer time scale than unbiased simulations can usually reach. Normally, the enhanced sampling methods, which may combine several short simulation trajectories, are used to modeling the unbinding process of two proteins.30,43 The association of protein–protein is thus predicted.43 Therefore, in this work, atomistic simulations will be performed to reveal the insights at the atomic level of the binding process of an antibody to various SARS-CoV-2 variants RBD. Structural changes of the SARS-CoV-2 RBD + antibody complexes were characterized via MD simulations. Besides, because the glycan linked with the residue Asn343,22 which is far from the binding surface of RBD, the effect of the glycan on the binding of RBD and antibody can be thus neglected.44 Moreover, the obtained superposition of RBD with and without glycan confirmed the argument that the neglection of glycan probably adopts small effects on the interacted picture between RBD and antibody (Fig. S2 and S3 of the ESI file). Thermodynamics and kinetics of the binding process were then revealed via a combination of SMD and US simulations.

Unbiased MD simulations were carried out to understand the structural change at the atomistic level of 501Y.V2 RBD + antibodies since the binding affinity of the antibodies to 501Y.V2 RBD was altered according to the recent report.8,10,11,14 The stabilized conformations of the RBD + antibody complexes were investigated over the equilibrium trajectories (cf. Fig. S4 of the ESI file). Moreover, the obtained superposition of calculated metrics in the different intervals confirmed the convergence of unbiased MD simulations (Fig. S5–S8 of the ESI file). Furthermore, the structural changes of the complexes were reported in Fig. S9 and S10 of the ESI file. In particular, the backbone root-mean-square deviation (RMSD) and SASA of the complexes were enlarged when the 501Y.V2 variant were induced. However, the HB and NB contacts between RBD and NAb/fNAb were significantly reduced due to mutations. The obtained results imply that the protein–protein binding affinity between 501Y.V2 RBD to NAb/fNAb was decreased in comparison with WT one.

In order to estimate the representative conformations of the complexes, the two-dimensional FEL was generated using “gmx sham” tool.33,45 Two coordinates constructing FEL were first and second eigenvectors, which were computed using the PCA method.33 The obtained results were described in Fig. 2 (FEL for individual trajectories were presented in Fig. S11–S14 of the ESI). Clearly, the 501Y.V2 variant increases the number of the FEL local minima implying that the 501Y.V2 complex is more flexible than the WT one. It also suggests that the binding free energy ΔGb between 501Y.V2 RBD and NAbs is probably reduced. The obtained results are in good consistent with NB and HB contacts analyses above.


image file: d1ra04134g-f2.tif
Fig. 2 Free energy landscape of RBD + antibody complexes was constructed using PCA method. In particular, (A) presents the FEL of the WT RBD + fNAb over 4 independent MD trajectories; (B) mentions the FEL of the WT RBD + NAb over 4 independent MD trajectories; (C) describes the FEL of the 501Y.V2 RBD + fNAb over 4 independent MD trajectories; and (D) denotes the FEL of the 501Y.V2 RBD + NAb over 4 independent MD trajectories.

The WT RBD + fNAb only formed one minimum noted as w1 in Fig. 2A, which is located at (CV1; CV2) coordinates of (0.40; 0.40). The representative structures of the complex corresponding to w1 was calculated using the clustering method with a backbone RMSD cutoff of 0.2 nm. In particular, the antibody adopted HBs to 4 residues of the WT RBD including G447, Y449, N450, and E484 (cf. Fig. 3). These results suggest that a mutation E484K will significantly alter the binding affinity/mechanism of the RBD + fNAb. Two minima were observed in FEL of 501Y.V2 RBD + fNAb, which are located at (CV1; CV2) coordinates of (1.60; −1.40) and (−3.60; −1.00) denoted as m1 and m2, respectively. The corresponding population of m1 and m2 is 75 and 25%, respectively. Analyzing the representative structure m1, the antibody was found to be able to form HBs to the residues K444, G447, Y449, and N450 of the 501Y.V2 RBD. The corresponding residues of m2, which formed HBs to RBD 2–4, are G447, Y449, N450, and K484 (cf. Fig. 3). The observed structural changes imply that the binding affinity and kinetics between RBD and fNAb probably change. A similar story of RBD + NAb, which is mentioned in detail below, was obtained and confirmed the results.


image file: d1ra04134g-f3.tif
Fig. 3 The representative structures of WT and 501Y.V2 RBD + fNAb in different perspective. The structures corresponds to the minima w1, m1, and m2, which were mentioned in Fig. 2. The structures were obtained using the clustering method with a backbone RMSD cut-off of 0.2 nm. The interaction diagram between RBD and fNAb was estimated using PyMOL tool.

The WT/501Y.V2 RBD + NAb systems were also investigated. FEL of the complexes was significantly altered when the mutations were induced. The WT RBD + NAb formed two minima, which were shown in Fig. 2B. These minima located at (CV1; CV2) coordinates of (0.63; 0.75) and (0.63; 3.38) denoting as W1 and W2, respectively. The corresponding population of W1 and W2 is 63 and 37%, respectively. Besides that, the 501Y.V2 RBD + NAb FEL (Fig. 2D) adopted three minima, which located at (CV1; CV2) coordinates of (7.50; 4.13), (5.63; −4.50), and (−18.8; 1.50) labelling as M1, M2, and M3, respectively. The corresponding population of M1, M2, and M3 is 45, 38, and 17%, respectively. Analyzing the complex W1, the HBs were observed between antibody and residues G447, Y449, N450, and E484 of the RBD that is in good consistency to the w1 case. However, HBs were only found between the NAb and residue E484 of the WT RBD in the complex W2 (Fig. 4). The obtained results indicate that residue E484 plays an important role in the binding process of the antibody to the RBD that is in good consistent with the recent work.13 Replacing the E484 with another residue probably modifies the binding mechanism of the antibodies to RBD rather than substitutions at the different positions. Moreover, it should be noted that in the 501Y.V2 variant induced, a lysine residue substitutes the glutamate residue at the sequence 484. The replacement probably terminates the HBs and weakening the attracted force between the NAb and the RBD. The argument was confirmed via evaluations of the representative structures of 501Y.V2 RBD + NAb complexes. In conformation M1, the HBs between NAb and the residues G447, Y449, and N450 of RBD were found. The residues G447, Y449, N450, and T470 of 501Y.V2 RBD procedure HBs to NAb in conformation M2. Furthermore, the NAb only found two HBs to the residue E471 and N481 of the 501Y.V2 RBD. The free energy approach should be carried out to clarify the change of binding affinity upon the structural changes of the 501Y.V2 RBD + NAb complexes.


image file: d1ra04134g-f4.tif
Fig. 4 The representative structures of WT and 501Y.V2 RBD + NAb corresponding to the minima W1, W2, M1, M2, and M3, which were mentioned in Fig. 2. The structures were obtained using the clustering method with a backbone RMSD cut-off of 0.2 nm. The interaction diagram between RBD and NAb were obtained using PyMOL tool.

As discussion above, the RBD + fNAb structure is more flexible when the 501Y.V2 variant was induced. The binding affinity/mechanism of the complex are thus altered. In this work, a combination of SMD/US simulations were carried out to probe the change in RBD + NAbs association. The SMD was used to generate US windows (cf. the ESI file). The free energy profile was then calculated using the WHAM.32 The binding free energy ΔGb between RBD and NAbs is able to calculate via PMF curve as mentioned in Fig. 1E.30,31 Moreover, the free energy barriers ΔG++on and ΔG++off, which were associated with the binding kinetic rate constant kon and the unbinding kinetic rate constant koff can be also estimated, respectively.

The calculated results for free energy barriers (cf. Table 1) indicated that the NAbs will bind to 501Y.V2 RBD more difficult than WT one because of the larger ΔG++on. NAbs are much easier to bind to than to unbind from RBD, because the ΔG++off is larger than the ΔG++on. However, in the M3 case, the ΔG++off = 0.14 ± 0.18 kcal mol−1 is significantly smaller than the ΔG++on = 2.83 ± 0.65 kcal mol−1 indicating that it takes more time for NAb to bind to 501Y.V2 RBD for them unbind. Moreover, the observations were also supported by the binding free energy, ΔGb, calculations, in which the thermodynamic metric corresponding to the association between NAbs and RBD is significantly decreased when the 501Y.V2 variant was induced (Table 1). The NAb is thus resisted to bind to 501Y.V2 RBD. Therefore, it may be argued that the 501Y.V2 variant could reduce the vaccine efficiency. The observation is in good agreement with the experimental data that the neutralizing antibody is weaker bind to 501Y.V2 spike protein than to WT one.8,10,11,14 However, besides, the free energy value over the population of minima and the predicted kinetic coefficients were reported in Table S2 of the ESI file. The predicted values should only use for qualitative comparisons between the RBD + antibody complexes, because the absolute value of the metrics is quite different from experimental values.46

Table 1 The calculated results using SMD and US simulationsa
No. System FMax W ΔG++on ΔG++off ΔGb
a The calculated results over SMD and US simulations. The details of free energy profile and histograms over US simulations were reported in Fig. S15–S18 of the ESI file.
1 WT RBD + fNAb (w1) 1388.0 ± 18.6 139.9 ± 3.3 0.24 ± 0.20 18.31 ± 0.82 −18.07 ± 0.84
2 501Y.V2 RBD + fNAb (m1) 859.3 ± 40.7 72.8 ± 3.4 0.81 ± 0.26 12.00 ± 0.82 −11.19 ± 0.77
3 501Y.V2 RBD + fNAb (m2) 1007.4 ± 31.2 86.8 ± 3.1 0.73 ± 0.11 11.62 ± 0.54 −10.89 ± 0.55
4 WT RBD + NAb (W1) 1133.6 ± 39.0 178.8 ± 8.5 0.36 ± 0.75 39.82 ± 1.31 −39.46 ± 1.08
5 WT RBD + NAb (W2) 1137.6 ± 25.8 181.5 ± 6.5 0.76 ± 0.29 43.36 ± 0.73 −42.60 ± 0.67
6 501Y.V2 RBD + NAb (M1) 745.5 ± 25.4 96.1 ± 3.6 0.62 ± 0.15 21.64 ± 0.66 −20.93 ± 0.68
7 501Y.V2 RBD + NAb (M2) 748.0 ± 33.1 81.6 ± 6.5 0.42 ± 0.24 16.16 ± 0.88 −15.74 ± 0.91
8 501Y.V2 RBD + NAb (M3) 470.2 ± 25.7 50.5 ± 4.5 0.14 ± 0.18 2.83 ± 0.65 −2.70 ± 0.68


The collective-variable FEL,47 was constructed by number of contacts between two proteins within 0.45 nm and the displacement of the antibody, revealed the unbinding pathway of NAbs. The obtained FEL was shown in Fig. 5 and S19 of the ESI file. The representative structures of the complexes within a backbone RMSD of 0.2 nm were then estimated using clustering method.33 The unbinding pathways were significantly altered under effects of the 501Y.V2 variant. A larger number of transition states of the WT RBD + fNAb complex implies that it is hard to unbind the antibody from WT system than 501Y.V2 variant. Moreover, the representative structures B, b, and b′ correspond to the binding model of the RBD + fNAb complexes. The structures D7, d6, and d4′ respond to the minima where the fNAb completely detached from RBD. The other conformations correspond to dissociated structures along unbinding pathways. The similar picture was also observed when the RBD + NAb complexes were investigated (Fig. S19 of the ESI file).


image file: d1ra04134g-f5.tif
Fig. 5 The collective-variable FEL revealed the unbinding pathways of fNAb from the binding mode with WT/501Y.V2 RBD over US simulations. The representative structures of complexes were also estimated using the clustering method with a backbone RMSD cutoff of 0.2 nm.

Conclusions

In this work, the NAb resistance of 501Y.V2 variant was investigated using atomistic simulations. In particular, the binding pose of NAb/fNAb to WT/501Y.V2 RBD was revealed using atomistic simulations. Reducing number of HB and NB contacts between RBD and antibodies were observed when the 501Y.V2 variant was induced. Increasing FEL minima of 501Y.V2 RBD + NAb/fNAb in comparison with the WT RBD systems infer that the complex 501Y.V2 RBD + NAb/fNAb is more unstable than the WT one. Moreover, thermodynamics and kinetics of the binding process between RBD and NAb were also determined using SMD/US simulations. Interestingly, the binding free energy ΔGb of WT RBD + NAb/fNAb is significantly smaller than that of 501Y.V2 RBD + NAb/fNAb. It is consistent with the results of the binding kinetic rate constant kon and the unbinding kinetic rate constant koff. Poorly binding affinity of NAb/fNAb to 501Y.V2 RBD confirms the antibody resistance of the South African variant.8,10,11,14 Furthermore, the RBD + fNAb system can be used as an affordable model to investigate the change of the binding process between mutations RBD and antibodies. The required computing resources are thus reduced significantly.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was supported by Vietnam National Foundation for Science & Technology Development (NAFOSTED) grant #104.99-2019.57.

References

  1. C. O. Barnes, A. P. West, K. E. Huey-Tubman, M. A. G. Hoffmann, N. G. Sharaf, P. R. Hoffman, N. Koranda, H. B. Gristick, C. Gaebler, F. Muecksch, J. C. C. Lorenzi, S. Finkin, T. Hägglöf, A. Hurley, K. G. Millard, Y. Weisblum, F. Schmidt, T. Hatziioannou, P. D. Bieniasz, M. Caskey, D. F. Robbiani, M. C. Nussenzweig and P. J. Bjorkman, Cell, 2020, 182, 828–842 CrossRef CAS PubMed.
  2. E. de Wit, N. van Doremalen, D. Falzarano and V. J. Munster, Nat. Rev. Microbiol., 2016, 14, 523–534 CrossRef CAS PubMed.
  3. Worldometrics, COVID-19 Coronavirus Pandemic, https://www.worldometers.info/coronavirus/ Search PubMed.
  4. S. T. Ngo, N. Quynh Anh Pham, L. Thi Le, D.-H. Pham and V. V. Vu, J. Chem. Inf. Model., 2020, 60, 5771–5780 CrossRef CAS PubMed.
  5. D. Schoeman and B. C. Fielding, Virol. J., 2019, 16, 69 CrossRef PubMed.
  6. M. Hoffmann, H. Kleine-Weber, S. Schroeder, N. Kruger, T. Herrler, S. Erichsen, T. S. Schiergens, G. Herrler, N. H. Wu, A. Nitsche, M. A. Muller, C. Drosten and S. Pohlmann, Cell, 2020, 181, 271–280.e8 CrossRef CAS PubMed.
  7. J. Lan, J. Ge, J. Yu, S. Shan, H. Zhou, S. Fan, Q. Zhang, X. Shi, Q. Wang, L. Zhang and X. Wang, Nature, 2020, 581, 215–220 CrossRef CAS PubMed.
  8. C. K. Wibmer, F. Ayres, T. Hermanus, M. Madzivhandila, P. Kgagudi, B. Oosthuysen, B. E. Lambson, T. de Oliveira, M. Vermeulen, K. van der Berg, T. Rossouw, M. Boswell, V. Ueckermann, S. Meiring, A. von Gottberg, C. Cohen, L. Morris, J. N. Bhiman and P. L. Moore, Nat. Med., 2021, 27, 622–625 CrossRef CAS PubMed.
  9. WHO, Coronavirus disease 2019 (COVID-19) Situation Report – 52, 2020, https://www.who.int/docs/default-source/coronaviruse/20200312-sitrep-52-covid-19.pdf?sfvrsn=e2bfc9c0_2 Search PubMed.
  10. P. Wang, M. S. Nair, L. Liu, S. Iketani, Y. Luo, Y. Guo, M. Wang, J. Yu, B. Zhang, P. D. Kwong, B. S. Graham, J. R. Mascola, J. Y. Chang, M. T. Yin, M. Sobieszczyk, C. A. Kyratsous, L. Shapiro, Z. Sheng, Y. Huang and D. D. Ho, Nature, 2021, 593, 130–135 CrossRef CAS PubMed.
  11. M. Hoffmann, P. Arora, R. Groß, A. Seidel, B. F. Hörnich, A. S. Hahn, N. Krüger, L. Graichen, H. Hofmann-Winkler, A. Kempf, M. S. Winkler, S. Schulz, H.-M. Jäck, B. Jahrsdörfer, H. Schrezenmeier, M. Müller, A. Kleger, J. Münch and S. Pöhlmann, Cell, 2021, 184, 2384–2393 CrossRef CAS PubMed.
  12. Y. Weisblum, F. Schmidt, F. Zhang, J. DaSilva, D. Poston, J. C. C. Lorenzi, F. Muecksch, M. Rutkowska, H.-H. Hoffmann, E. Michailidis, C. Gaebler, M. Agudelo, A. Cho, Z. Wang, A. Gazumyan, M. Cipolla, L. Luchsinger, C. D. Hillyer, M. Caskey, D. F. Robbiani, C. M. Rice, M. C. Nussenzweig, T. Hatziioannou and P. D. Bieniasz, eLife, 2020, 9, e61312 CrossRef CAS PubMed.
  13. S. Jangra, C. Ye, R. Rathnasinghe, D. Stadlbauer, H. Alshammary, A. A. Amoako, M. H. Awawda, K. F. Beach, M. C. Bermúdez-González, R. L. Chernet, L. Q. Eaker, E. D. Ferreri, D. L. Floda, C. R. Gleason, G. Kleiner, D. Jurczyszak, J. C. Matthews, W. A. Mendez, L. C. F. Mulder, K. T. Russo, A.-B. T. Salimbangon, M. Saksena, A. S. Shin, L. A. Sominsky, K. Srivastava, F. Krammer, V. Simon, L. Martinez-Sobrido, A. García-Sastre and M. Schotsaert, Lancet Microbe, 2021, 2, E283–E284 CrossRef CAS PubMed.
  14. S. Cele, I. Gazy, L. Jackson, S.-H. Hwa, H. Tegally, G. Lustig, J. Giandhari, S. Pillay, E. Wilkinson, Y. Naidoo, F. Karim, Y. Ganga, K. Khan, M. Bernstein, A. B. Balazs, B. I. Gosnell, W. Hanekom, M.-Y. S. Moosa, R. J. Lessells, T. de Oliveira, A. Sigal, S. A. Ngs and C.-K. Team, Nature, 2021, 593, 142–146 CrossRef CAS PubMed.
  15. M. A. Tortorici, M. Beltramello, F. A. Lempp, D. Pinto, H. V. Dang, L. E. Rosen, M. McCallum, J. Bowen, A. Minola, S. Jaconi, F. Zatta, A. De Marco, B. Guarino, S. Bianchi, E. J. Lauron, H. Tucker, J. Zhou, A. Peter, C. Havenar-Daughton, J. A. Wojcechowskyj, J. B. Case, R. E. Chen, H. Kaiser, M. Montiel-Ruiz, M. Meury, N. Czudnochowski, R. Spreafico, J. Dillen, C. Ng, N. Sprugasci, K. Culap, F. Benigni, R. Abdelnabi, S.-Y. C. Foo, M. A. Schmid, E. Cameroni, A. Riva, A. Gabrieli, M. Galli, M. S. Pizzuto, J. Neyts, M. S. Diamond, H. W. Virgin, G. Snell, D. Corti, K. Fink and D. Veesler, Science, 2020, 370, 950–957 CrossRef CAS PubMed.
  16. L. Piccoli, Y.-J. Park, M. A. Tortorici, N. Czudnochowski, A. C. Walls, M. Beltramello, C. Silacci-Fregni, D. Pinto, L. E. Rosen, J. E. Bowen, O. J. Acton, S. Jaconi, B. Guarino, A. Minola, F. Zatta, N. Sprugasci, J. Bassi, A. Peter, A. De Marco, J. C. Nix, F. Mele, S. Jovic, B. F. Rodriguez, S. V. Gupta, F. Jin, G. Piumatti, G. Lo Presti, A. F. Pellanda, M. Biggiogero, M. Tarkowski, M. S. Pizzuto, E. Cameroni, C. Havenar-Daughton, M. Smithey, D. Hong, V. Lepori, E. Albanese, A. Ceschi, E. Bernasconi, L. Elzi, P. Ferrari, C. Garzoni, A. Riva, G. Snell, F. Sallusto, K. Fink, H. W. Virgin, A. Lanzavecchia, D. Corti and D. Veesler, Cell, 2020, 183, 1024–1042 CrossRef CAS PubMed.
  17. L. Casalino, Z. Gaieb, J. A. Goldsmith, C. K. Hjorth, A. C. Dommer, A. M. Harbison, C. A. Fogarty, E. P. Barros, B. C. Taylor, J. S. McLellan, E. Fadda and R. E. Amaro, ACS Cent. Sci., 2020, 6, 1722–1734 CrossRef CAS PubMed.
  18. H. Woo, S.-J. Park, Y. K. Choi, T. Park, M. Tanveer, Y. Cao, N. R. Kern, J. Lee, M. S. Yeom, T. I. Croll, C. Seok and W. Im, J. Phys. Chem. B, 2020, 124, 7128–7137 CrossRef CAS PubMed.
  19. B. Turoňová, M. Sikora, C. Schürmann, W. J. H. Hagen, S. Welsch, F. E. C. Blanc, S. von Bülow, M. Gecht, K. Bagola, C. Hörner, G. van Zandbergen, J. Landry, N. T. D. de Azevedo, S. Mosalaganti, A. Schwarz, R. Covino, M. D. Mühlebach, G. Hummer, J. Krijnse Locker and M. Beck, Science, 2020, 370, 203 CrossRef PubMed.
  20. B. Ju, Q. Zhang, J. Ge, R. Wang, J. Sun, X. Ge, J. Yu, S. Shan, B. Zhou, S. Song, X. Tang, J. Yu, J. Lan, J. Yuan, H. Wang, J. Zhao, S. Zhang, Y. Wang, X. Shi, L. Liu, J. Zhao, X. Wang, Z. Zhang and L. Zhang, Nature, 2020, 584, 115–119 CrossRef CAS PubMed.
  21. L. Schrödinger, The PyMOL Molecular Graphics System, Version 1.3r1, 2010 Search PubMed.
  22. S. M.-C. Gobeil, K. Janowska, S. McDowell, K. Mansouri, R. Parks, V. Stalls, M. F. Kopp, K. Manne, D. Li, K. Wiehe, K. O. Saunders, R. J. Edwards, B. Korber, B. F. Haynes, R. Henderson and P. Acharya, Science, 2021, eabi6226 CrossRef CAS PubMed.
  23. M. J. Abraham, T. Murtola, R. Schulz, S. Páll, J. C. Smith, B. Hess and E. Lindahl, SoftwareX, 2015, 1–2, 19–25 CrossRef.
  24. A. E. Aliev, M. Kulke, H. S. Khaneja, V. Chudasama, T. D. Sheppard and R. M. Lanigan, Proteins: Struct., Funct., Bioinf., 2014, 82, 195–215 CrossRef CAS PubMed.
  25. D. Patel, J. S. Patel and F. M. Ytreberg, J. Chem. Theory Comput., 2021, 17, 2457–2464 CrossRef CAS PubMed.
  26. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, J. Chem. Phys., 1983, 79, 926–935 CrossRef CAS.
  27. S. T. Ngo, H. M. Hung, D. T. Truong and M. T. Nguyen, Phys. Chem. Chem. Phys., 2017, 19, 1909–1919 RSC.
  28. N. T. Tung, P. Derreumaux, V. V. Vu, P. C. Nam and S. T. Ngo, ACS Omega, 2019, 4, 11066–11073 CrossRef CAS PubMed.
  29. T. Darden, D. York and L. Pedersen, J. Chem. Phys., 1993, 98, 10089–10092 CrossRef CAS.
  30. S. T. Ngo, K. B. Vu, L. M. Bui and V. V. Vu, ACS Omega, 2019, 4, 3887–3893 CrossRef CAS PubMed.
  31. S. T. Ngo, J. Comput. Chem., 2021, 42, 117–123 CrossRef CAS PubMed.
  32. J. S. Hub, B. L. de Groot and D. van der Spoel, J. Chem. Theory Comput., 2010, 6, 3713–3720 CrossRef CAS.
  33. E. Papaleo, P. Mereghetti, P. Fantucci, R. Grandori and L. De Gioia, J. Mol. Graph. Model., 2009, 27, 889–899 CrossRef CAS PubMed.
  34. B. Efron, Ann. Stat., 1979, 7, 1–26 Search PubMed.
  35. A. C. Pan, D. Jacobson, K. Yatsenko, D. Sritharan, T. M. Weinreich and D. E. Shaw, Proc. Natl. Acad. Sci. U. S. A., 2019, 116, 4244–4249 CrossRef CAS PubMed.
  36. P. J. Lupardus, M. Ultsch, H. Wallweber, P. Bir Kohli, A. R. Johnson and C. Eigenbrot, Proc. Natl. Acad. Sci. U.S.A., 2014, 111, 8025–8030 CrossRef PubMed.
  37. Y. Shan, K. Gnanasambandan, D. Ungureanu, E. T. Kim, H. Hammarén, K. Yamashita, O. Silvennoinen, D. E. Shaw and S. R. Hubbard, Nat. Struct. Mol. Biol., 2014, 21, 579–584 CrossRef CAS PubMed.
  38. C. Tang, J. Iwahara and G. M. Clore, Nature, 2006, 444, 383–386 CrossRef CAS PubMed.
  39. C. Frisch, A. R. Fersht and G. Schreiber, J. Mol. Biol., 2001, 308, 69–77 CrossRef CAS PubMed.
  40. F. Chen, H. Liu, H. Y. Sun, P. C. Pan, Y. Y. Li, D. Li and T. J. Hou, Phys. Chem. Chem. Phys., 2016, 18, 22129–22139 RSC.
  41. M. Sikora, S. von Bülow, F. E. C. Blanc, M. Gecht, R. Covino and G. Hummer, PLoS Comput. Biol., 2021, 17, e1008790 CrossRef CAS PubMed.
  42. S. T. Ngo, D. T. Truong, N. M. Tam and M. T. Nguyen, J. Mol. Graph. Model., 2017, 76, 1–10 CrossRef CAS PubMed.
  43. J. A. Lemkul and D. R. Bevan, J. Phys. Chem. B, 2010, 114, 1652–1660 CrossRef CAS PubMed.
  44. D. Hao, H. Wang, Y. Zang, L. Zhang, Z. Yang and S. Zhang, J. Chem. Inf. Model., 2021 DOI:10.1021/acs.jcim.1c00233.
  45. Y. Mu, P. H. Nguyen and G. Stock, Proteins: Struct., Funct., Bioinf., 2005, 58, 45–52 CrossRef CAS PubMed.
  46. X. Chi, X. Liu, C. Wang, X. Zhang, X. Li, J. Hou, L. Ren, Q. Jin, J. Wang and W. Yang, Nat. Commun., 2020, 11, 4528 CrossRef CAS PubMed.
  47. T. H. Nguyen, G. Rossetti, F. Arnesano, E. Ippoliti, G. Natile and P. Carloni, J. Chem. Theory Comput., 2014, 10, 3578–3584 CrossRef CAS PubMed.

Footnote

Electronic supplementary information (ESI) available: Systemic configuration information of investigated systems; structure of neutralizing antibody; backbone RMSD of RBD and RBD + glycan; superposition of representative structure of RBD and RBD + glycan over MD simulations; backbone RMSD of RBD + antibody over MD simulations; superposition between computed metrics of WT RBD + fNAb in different intervals; superposition between computed metrics of 501Y.V2 RBD + fNAb in different intervals; superposition between computed metrics of WT RBD + NAb in different intervals; superposition between computed metrics of 501Y.V2 RBD + NAb in different intervals; the difference between computed metrics of WT and 501Y.V2 RBD + fNAb; the difference between computed metrics of WT and 501Y.V2 RBD + NAb; FEL of WT/501Y.V2 RBD + fNAb/NAb in different MD trajectories; free energy profile over the unbinding process of fNAb out of the WT/501Y.V2 RBD + fNAb complexes. The results were obtained using the WHAM calculation; the histograms of US simulations over the unbinding process of fNAb out of the WT/501Y.V2 RBD + fNAb complexes; free energy profile over the unbinding process of NAb out of the WT/501Y.V2 RBD + NAb complexes. The results were obtained using the WHAM calculation; the histograms of US simulations over the unbinding process of NAb out of the WT/501Y.V2 RBD + NAb complexes; and the collective-variable FEL revealed the unbinding pathways of NAb from the WT/501Y.V2 RBD + NAb complexes. See DOI: 10.1039/d1ra04134g

This journal is © The Royal Society of Chemistry 2021
Click here to see how this site uses Cookies. View our privacy policy here.