ACE-2-Derived Biomimetic Peptides for the Inhibition of Spike Protein of SARS-CoV-2.

SARS-CoV-2, a novel coronavirus causing overwhelming death and infection worldwide, has emerged as a pandemic. Compared to its predecessor SARS-CoV, SARS-CoV-2 is more infective for being highly contagious and exhibiting tighter binding with host angiotensin-converting enzyme 2 (hACE-2). The entry of the virus into host cells is mediated by the interaction of its spike protein with hACE-2. Thus, a peptide that has a resemblance to hACE-2 but can overpower the spike protein-hACE-2 interaction will be a potential therapeutic to contain this virus. The non-interacting residues in the receptor-binding domain of hACE-2 have been mutated to generate a library of 136 new peptides. Out of this library, docking and virtual screening discover seven peptides that can exert a stronger interaction with the spike protein than hACE-2. A peptide derived from simultaneous mutation of all the non-interacting residues of hACE-2 yields almost three-fold stronger interaction than hACE-2 and thus turns out here to be the best peptide inhibitor of the novel coronavirus. The binding of the best peptide inhibitor with the spike protein is explored further by molecular dynamics, free energy, and principal component analysis, which demonstrate its efficacy compared to hACE-2. The delivery of the screened inhibitors with nanocarriers like metal-organic frameworks will be worthy of further consideration to boost their efficacy.

Introduction

The current death toll and infected cases surged into millions are caused by the novel coronavirus COVID-19, engendering an unprecedented global containment. While no approved cure is available so far for the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), the COVID-19 pandemic has put health, lives, economies, societal relations, and humanity at risk.[1] Based on epidemiological data, unlike its predecessors, SARS-CoV-2 is highly contagious in humans, spreading through close contact and respiratory droplets, reasoned for the occurrence of this pandemic.[2,3] It targets the lower respiratory system to induce viral pneumonia similar to earlier SARS-CoV and MERS-CoV.[4,5]

SARS-CoV-2 is a β-coronavirus and belongs to Coronavirinae subfamily of Coronaviridae genera, consisting of a single-stranded positive sense RNA genome with nucleocapsid protein, covered by a lipid bilayer membrane containing hemagglutinin-esterase dimer, envelope, and membrane proteins and spangled with spikes of glycoproteins, which give a coronary appearance.[6] Based on the sequence alignment and homology analysis, SARS-CoV and SARS-CoV-2 share a strongly conserved receptor-binding domain (RBD) with a sequence identity of 72.4% (Figure ). There are 14 residues in the RBD region of SARS-CoV-2 which actively participate in the interaction with host angiotensin-converting enzyme 2 (hACE-2). Out of 14 residues, 8 amino acid residues are common between both viruses, Tyr449/Tyr436, Tyr453/Tyr440, Asn487/Asn473, Tyr489/Tyr475, Gly496/Gly482, Thr500/Thr486, Gly502/Gly488, and Tyr505/Tyr491 of SARS-CoV-2/SARS-CoV. Five residues, albeit not common structurally, render similar biochemical properties in both, Leu455/Tyr442, Phe456/Leu443, Phe486/Leu472, Gln493/Asn479, and Asn501/Thr487 of SARS-CoV-2/SARS-CoV. The remaining one amino acid differs—Gln498/Tyr484 (SARS-CoV-2/SARS-CoV).[7−9] The surface plasmon resonance showed that the binding affinity of SARS-CoV-2 spike protein to hACE-2 is around 10–20 times greater than that of SARS-CoV spike protein,[10,11] making it more harmful than the latter. This is attributed to several mutations in the SARS-CoV-2 RBD and a more compact conformation of hACE-2-binding ridge, which stabilize two virus-binding hotspots of hACE-2 more.[8] In passing, the structural analysis clearly indicates that a convergent evolution between the SARS-CoV-2 and SARS-CoV RBDs has boosted the ACE-2 binding.[5,12,13] Thus, the interaction between the viral spike protein and hACE-2 on the host cell surface is of great importance to initiate infection. The host ACE-2 interacts with the RBD of the spike protein, mainly via the α-1 helix, as most of the interactions are from α-1 helix of hACE-2.[9] Therefore, the design of new peptide inhibitors by mutations to the α-1 helix having significantly higher binding affinity to the RBD than α-1 helix/hACE-2 is a very important and an unexplored avenue to discover potent therapeutics of the novel coronavirus to prevent hACE-2–RBD interactions. Despite cumulative efforts at all corners of the globe, no confirmed effective treatment is yet ready for COVID-19. Here, we report for the first time the design of biomimetic peptide inhibitors of SARS-CoV-2 by mutation.

Figure 1

Sequence alignment between the RBD regions of SARS-CoV and SARS-CoV-2 (red and blue colors, respectively, indicate identical and similar residues).

Therapeutic peptides have numerous major advantages over proteins or antibodies because of small size, easy synthesis, and penetration through the membranes of cells, with high specificity and affinity. Earlier, anti-microbial peptides such as HR1P and HR2P had shown antiviral activity against MERS-CoV whose cell fusion is mediated by the spike protein.[14] Likewise, synthesized peptides SARSWW-III, IV, and MHVWW-IV were effective against SARS-CoV and murine hepatitis virus (MHV) (S2 subunits) having IC50 = 2–4 μM.[15]

In this study, a library of 136 biomimetic peptide inhibitors against SARS-CoV-2 was designed, docked, and simulated, which mimic hACE-2 but incorporate necessary mutated components for the non-interacting residues in hACE-2’s virus-binding domain, utilizing the recently published crystal structure (PDB code: 6M0J),[16] to enhance interaction with the spike protein. The α-1 helix of hACE-2 was taken as a reference peptide for designing the peptide library, which contains key amino acids that take part in the interaction with the RBD of the spike protein. The α-1-helix-derived designed biomimetic peptide inhibitors each of 24 residues long have shown greater interaction than α-1 helix of hACE-2. All the inhibitors are structurally stable and have relatively small molecular weights. Out of them, seven peptide inhibitors manifest greater binding affinity to the spike protein than the α-1 helix. Stability and interaction of the best-screened peptide inhibitor complexed with the spike protein were validated using molecular dynamics (MD) simulation and molecular mechanics Poisson–Boltzmann surface area (MMPBSA) binding free energy calculation. Compared to the α-1 helix-spike protein, there is a greater stability of the complex of the designed best peptide implying it is a potent inhibitor of the spike protein which can resist the hACE-2–spike interaction.

Methods

Molecular Docking of Peptide Inhibitors

The all non-interacting residues of α-1 helix were mutated using visual MD.[17] The mutated peptides thus formed are energy minimized applying the amber99sb force field.[18] The peptide docking was performed by using the high ambiguity-driven protein–protein docking (HADDOCK) server.[19] HADDOCK facilitates information-driven, flexible docking for modeling of bimolecular complexes. With a wide capability to deal with protein–protein, protein–nucleic acids, protein–ligand, and so forth, docking, HADDOCK distinctly differs from ab initio docking methods by encoding information from identified or predicted protein interfaces in ambiguous interaction restraints (AIRs) to drive the docking process. The HADDOCK score was the sum of van der Waals energy (EvdW), electrostatic energy (Eelec), desolvation energy (Edesol), energy from restraint violations (EAIR), and the buried surface area. The HADDOCK score was again determined after simulation and water refining that is HADDOCK score = 1.0EvdW + 0.2Eelec + 1.0Edesol + 0.1EAIR.[19]

MD Simulations

MD simulations were run for the spike protein bound to the α-1 helix and the mutated peptide no. 13 using GROMACS 5.1.4 for 150 ns.[20−24] The amber99sb force field was used for the spike protein and peptides.[18] The respective docked complexes were engrossed in a cubical box of simple point charge water molecules. The complexes were neutralized by adding an appropriate number of Na+ and Cl– ions. To get rid of short-range bad contacts, energy minimization was performed using the steepest descent method for 50,000 steps. Then, the systems were subjected to 150 ns MD simulation runs at 298 K temperature and 1 bar pressure, using 0.002 ps time step and the Berendsen thermostat. Finally, the root-mean-square-deviation (RMSD), root-mean-square-fluctuation (RMSF), radius of gyration (Rg), essential dynamics (ED), and the average structure after MD simulations were evaluated to investigate conformational changes and the stability of peptide 13 inhibitor and α-1 helix bound to the spike protein.

Binding Free Energy Calculation

Subsequently, binding free energies of the spike bound best peptide inhibitor and α-1 helix complexes were calculated on the 150 ns of MD simulation trajectories taking 15,000 configuration snapshots. The total binding free energy and its components: electrostatic energy, polar and apolar solvation energies, and the van der Waals energy were calculated using the MM-PBSA method executed by the g_mmpbsa tool.[25,26] The per-residue energy contribution of the complexes was also estimated.

Results and Discussion

Designing Peptide Inhibitor Library

To construct a peptide library, the crystal structure of hACE-2 bound to the spike protein with PDB ID: 6M0J was considered as a reference in the present study.[16] Between the RBD region of the spike protein and hACE-2, there exist a total of fifteen interactions, which include 13 hydrogen bonds (h-bond) and 2 salt–bridge interactions. The interacting residues pertaining to the host receptor hACE-2 are GLN24, ASP30 (2), HIS34, GLU35, GLU37, ASP38, TYR41 (2), GLN42 (2), TYR83 (2), LYS353, and ARG393. Most of the interactions are from the hACE-2 α-1 helix domain because of which it was taken here as a template to design new peptide inhibitors. Having a greater affinity of the designed biomimetic peptides to the RBD may warrant to prevent hACE-2–RBD interactions. To make highly specific and stable peptide binders, the α-1 helix of hACE-2 consisting of 24 amino acids located between 21 and 44 residues was taken. Out of these 24 amino acid residues, GLN24, ASP30, LYS31, GLU35, GLU37, ASP38, TYR41, and GLN42 are seen presently interacting with the spike protein RBD. To design a peptide library, 12 non-interacting residues of α-1 helix were mutated with all possible combinations to generate a total of 136 biomimetic peptides.

A single-point mutation and simultaneous mutation of non-interacting residues of hACE-2 in the active region, namely, for GLU22, GLU23, ALA25, LYS26, THR27, PHE28, LEU29, PHE32, ASN33, ALA36, LEU39, and PHE40, at multiple locations have been carried out to form a 136-peptide library (Table S1). The exhaustive list of the top 13 peptide sequences filtered out is presented in Table , where all these mutations are marked in red color.

Table 1

Sequences of the Top 13 Peptide Inhibitors Derived from the Mutation of α-1 Helix of hACE-2 Are Presenteda

Mutated residues are marked in red color.

Molecular Docking and Interaction Analysis

All the designed 136 mutated peptides and α-1 helix (reference peptide) of hACE-2 were docked in the active site of the spike protein. The active-site residues 417, 446, 449, 487, 489, 493, 500, 501, 502, and 505 from the receptor-binding region of the spike protein were selected owing to participating in the interaction with hACE-2. The best docking poses (Figure ) which have the highest negative HADDOCK[19] scores were filtered out. As a consequence, the best 13 peptides resulted out of a mutation of 12 different residues were identified among the top HADDOCK scores, the complexes of these 13 peptides with spike RBD are shown in Figure .

Figure 2

Binding poses of the top 13 designed biomimetic peptide inhibitors in the RBD of the spike protein. Peptides, spike protein, and hACE-2 are shown in red, orange, and pink colors, respectively.

For the first 12 high scoring mutated peptides reported in the table, the best docking poses inevitably show good interaction with the spike protein. Finally, the best peptide (no. 13) was obtained by mutating all 12 non-interacting residues of the α-1 helix simultaneously with the best scoring residues of 12 initially screened out peptides. The HADDOCK scores are given for all mutations in Table S1. Based on these scores of all probable mutations at every non-interacting site, Figure S1 gives the rationale behind the selection of mutant residues in the best peptide 13 inhibitor owing to their lowest energy. Table presents the docking scores and spike’s interacting residues (via h-bond formation) for the reference and the top 13 mutated peptides screened out pertaining to the mutation of the non-interacting residue (s). The spike protein complexes of the peptide inhibitor 13 having the best docking score and α-1 helix were selected for MD simulation and binding free energy calculation.

Table 2

Docking Scores and H-Bond Interactions of the Best-Screened 13 Mutated Peptide Inhibitors along with the Reference Peptide α-1 Helix

sl. no.	non-interacting residues	mutation	docking score (kcal/mol)	hydrogen bond interaction within 4 Å
α-1 helix (reference)	GLU 22, GLU 23, ALA 25, LYS 26, THR 27, PHE 28, LEU 29, PHE32, ASN33, ALA36, LEU39, PHE40		–118.5 ± 5.5	LYS417, GLY446, TYR449(2), ASN487, THR500, TYR505
1	GLU 22	ASP	–111.4 ± 2.4	GLU484, PHE486, GLN493, GLN498(2), THR500, ASN501 (2), TYR505, GLY496, GLY446 (2), GLY485
2	GLU 23	TRP	–123.1 ± 7.7	ARG403(2), TYR449, ASN487(2), TYR489, TYR505
3	ALA 25	PHE	–104.5 ± 1.4	GLU484, GLN498, ASN501, TYR505
4	LYS 26	TRP	–132.4 ± 6.3	ARG403(2), LYS417(2), TYR449, THR500, TYR505,
5	THR 27	PHE	–125.2 ± 2.8	LYS417, TYR473, ASN487, TYR489, GLN493, THR500(2), TYR505(2)
6	PHE 28	HIS	–99.8 ± 1.8	GLU484, TYR453, TYR505(2), GLN493, PHE490
7	LEU 29	TYR	–127.6 ± 23.6	ARG403, LYS417(2), ARG408, GLN414, THR415, LYS417, TYR505, THR415, TYR421, GLU484
8	PHE32	TRP	–119.4 ± 3.3	GLU484(2), TYR453, PHE486, GLN493, ASN501, GLY496, GLY485, CYS488
9	ASN33	ASP	–102.3 ± 2.9	ARG403(2), PHE486, ASN487, TYR489, TYR505, TYR449, GLN498
10	ALA36	TRP	–116.8 ± 6.0	ARG403, ARG408, TYR505, GLN493, GLY485, GLU484, TYR453
11	LEU39	TRP	–118.9 ± 3.1	ARG403, TYR453, ASN487, TYR489, GLY496, TYR453
12	PHE40	TRP	–111.4 ± 0.8	THR415, LYS417, TYR489, GLN493, TYR505, GLY485
13 (the best designed peptide)	GLU 22, GLU 23, ALA 25, LYS 26, THR 27, PHE 28, LEU 29, PHE32, ASN33, ALA36, LEU39, PHE40	above 12 mutations combined	–150 ± 3.7	ARG403 (2), ARG417 (2), TYR421, TYR449, TYR453, GLN493, TYR505, THR500 and ALA475

There are 7 peptides which have higher docking scores, ranging between −118.9 ± 3.1 to −150 ± 3.7 kcal/mol, than the reference peptide α-1 helix of hACE-2 (−118.5 ± 5.5 kcal/mol). Among all, the peptide inhibitor 13 exhibits the highest docking score of −150 ± 3.7 kcal/mol along with a greater number of interactions with the spike RBD (Table ). The interacting residues of peptide 13 within a 3.10 Å distance of the spike RBD are ARG403 (2), ARG417 (2), TYR421, TYR449, TYR453, GLN493, TYR505, THR500, and ALA475 (Figure ).

Figure 3

Interaction of peptide inhibitor 13 with the spike RBD of SARS-CoV-2 (interacting residues are labelled in the inset for clarity).

Apart from a significantly lower binding energy than that for α-1 helix of hACE-2, the designed peptide 13 has 11 h-bonds and 4 salt–bridges [ARG403(2) and LYS417(2)] in comparison to the 7 h-bonds and 2 salt–bridges (LYS417 and ARG403) of the former. The closer proximity of 3.10 Å indeed reflects a very strong interaction and inhibition capability of the peptide 13 toward the virus spike protein.

MD Simulation and Binding Free Energy Calculation

MD simulations were performed to obtain deeper insights into the stability and dynamic properties of the spike-bound α-1 helix and the designed best-screened peptide inhibitor. This is quite useful to assess how the naturally occurring motions of proteins or peptides can affect the protein–peptide interactions, the stability of the bound complexes, and their conformations due to binding. Superimposition of the initial structure and the final structure (averaged over the MD trajectory) of the spike protein in the α-1 helix and the best peptide inhibitor complexes displays in Figure that the active site residues in the latter have a small variation than that of the former, while in other regions of the protein, it is the reverse.

Figure 4

Superimposition of the initial (green) and final (cyan/pink) structures of the spike protein bound to (A) α-1 helix and (B) designed peptide inhibitor 13 in MD simulations.

Consistent with the higher docking score, the bar plot, Figure shows that the average RMSD of the key amino acids in the RBD of the spike protein in the peptide inhibitor 13 complex is lower (0.190 ± 0.020 nm) in comparison to that of the α-1 helix complex (0.203 ± 0.018 nm). Because of fluctuations in the non-RBD regions, the whole spike protein has a larger RMSD value in the peptide inhibitor 13 complex than the other complex; see Figure S2. Such fluctuations can be ignored as it does not affect the binding, demonstrated by the significantly lower binding free energy of the designed peptide inhibitor than that for α-1 helix in post-dynamics analysis, −252.3 ± 19.6 and −87.5 ± 10.7 kJ/mol, respectively (Table ). Together with Figure S2, the RMSD versus time plot (Figure S3) also indicates conformational changes in the protein of the best inhibitor complex, except for the RBD (Figures and 5) that renders stability at the bonded state.

Figure 5

Average RMSD and binding free energy of the key amino acids of the spike protein in the complexes of α-1 helix and the peptide inhibitor 13 are shown as bar plots with standard deviations as error bars.

Table 3

Decomposition of the Binding Free Energy into the Van der Waals, Electrostatic, Polar Solvation, and Apolar Solvation Energies (kJ/mol) in the α-1 Helix and Peptide Inhibitor 13 Bound Spike Protein Complexes

peptide inhibitors bound spike protein complexes	van der Waals energy	electrostatic energy	polar solvation energy	apolar solvation energy	binding free energy
α-1 helix	–28.0 ± 5.5	–155.0 ± 18.2	100 ± 16.1	–4.1 ± 0.8	–87.5 ± 10.7
peptide inhibitor 13	–116.8 ± 11.1	–478.7 ± 37.0	358.7 ± 32.5	–15.2 ± 1.3	–252.3 ± 19.6

In a strong correlation with the RMSD, the average RMSF value in Figure for the binding-site residues, that is, for the key amino acids of the spike protein is also lower in the best peptide inhibitor complex (0.125 ± 0.067 nm) compared to that of the α-1 helix (0.127 ± 0.070 nm), of course for the entire protein, the trend is reversing as previously, 0.159 ± 0.087 and 0.124 ± 0.056 nm, respectively. Between two complexes, to explain the higher RMSF value of the whole protein in the best inhibitor complex, Figure S4 is very percipient where the RMSF of the Cα atom of each residue has been plotted for both. While most of the residues of the spike protein in both the complexes exhibit similar fluctuations, residues lying between 380 and 390 show larger peaks in the best inhibitor complex. Such residues belong to a loop region of the protein giving rise to larger fluctuations and contributing to the overall higher RMSF. To our expectation, the active region extending from 400 to 505 residues remains relatively intact with the peptide 13 that suggests a stronger binding and larger stabilization with this inhibitor. Undoubtedly, here, the origins of fluctuations are the loops, not the active region.

Figure 6

Average RMSF and Rg of the key amino acids of the spike protein in the complexes of the α-1 helix and the peptide inhibitor 13 are shown as bar plots with error bars.

The mass-weighted root-mean-square distance of the residues determined from the common center of mass is the radius of gyration (Rg) accounting for the level of compaction or the overall dimension of protein. Compaction increases with interaction strength. The lower average value of Rg for the key amino acids in the inhibitor complex (1.348 ± 0.011 nm) than that of the α-1 helix (1.583 ± 0.010 nm) shown in Figure is an indication of that, emphasizing again a stronger interaction and greater stability prevailing between the spike RBD and the best peptide inhibitor. However, the region minus the RBD of the spike protein results in larger values of Rg, see Figures S2 and S5, for the peptide inhibitor than the α-1 helix.

The h-bond interactions turn out to be a major driving force for inhibition. The best peptide complex has a larger number of h-bonds (7.88 ± 1.87) than that for the α-1 helix (5.39 ± 2.13). During dynamics, the h-bonds of the latter undergo large oscillations as obvious in Figure S6. While the RBD contributes only 1.29 ± 0.85 h-bonds to the α-1 helix complex, this anchors all the h-bonds in the case of the best peptide, ascertaining a robust binding that could potentially inhibit the virus. Consequently, the stability of the corresponding complex is higher in comparison to that of the reference peptide.

The binding free energies reported earlier (Figure ) are decomposed into several components in Table . The peptide inhibitor 13 has almost three times higher the binding strength of the α-1 helix in which the van der Waals and electrostatic interactions make a significant difference. These interactions contribute, respectively, −116.8 ± 11.1 and −478.7 ± 37.0 kJ/mol for the designed inhibitor and, respectively, −28.0 ± 5.5 and −155.0 ± 18.2 kJ/mol for the α-1 helix to binding. Although small, another favorable predominant contribution is the apolar solvation energy. The reinforcing negative interactions make the mutant peptide 13 indeed a powerful inhibitor of the spike protein to shield the attack of the virus and consequently its entry into host cells. For a better insight, the contributions of every residue of the α-1 helix and the best peptide inhibitor 13 in complex formation with the spike protein are provided in Table S2.

Projecting these contributions (the per-residue binding free energies) in a bar plot in Figure immediately reflects that most of the mutated residues (red) have favorable or lower binding energy in comparison to the wild-types (black). The E22D, K26W, N33D, L39E, and F40W mutations display larger contributions to binding to the spike protein. Apart from these, some amino acid residues contributing higher energy than α-1 helix are D30, E35, E37, D38, and S44.

Figure 7

Per-residue contribution plot for the designed peptide inhibitor 13 (red) and α-1 helix (black) bound to the spike protein. The major contributors are highlighted in yellow color.

Figure S7 presents the 2D plot for principal component analysis or ED of the whole spike protein to portray the directions of principal motions. Taken the first two PCs into consideration, simulation results reveal a comparable subspace dimension for the protein in the complexes of α-1 helix (black) and peptide inhibitor 13 (red) in the figure. The protein bound to the inhibitor 13 covers subspace similar to that of the α-1 helix; accordingly, the trace values of covariance matrixes extracting the spatial data of the protein are similar.

The overall analysis suggests that the peptide inhibitor 13-spike protein complex has greater stability/interaction than the α-1 helix-spike protein complex.

Conclusions

To hinder the fusion of the spike protein with hACE-2, 136 peptide-inhibitor analogues of hACE-2 have been designed by sequential and simultaneous mutations of the non-interacting residues in α-1 helix of hACE-2. Seven peptides that have a stronger interaction with the spike protein than the α-1 helix of hACE-2 are discovered. The affinity of the best screened peptide inhibitor, derived from mutating all the non-interacting residues, to the spike protein is almost thrice that of hACE-2. Such a strong affinity can misguide the spike protein from attaching to hACE-2 preventing viral entry and infection into host cells. Besides, the designed peptides having docking scores lower than −123 kcal/mol could be also a promising inhibitor. Post-MD analyses, for example, the RMSD, RMSF, Rg, binding free energy, number of h-bonds and salt bridges, and so forth, clearly demonstrate tighter binding and larger stability of the best-screened peptide with the active site of the spike protein compared to the α-1 helix. Therefore, the best peptide inhibitor identified herein definitely has a greater chance to serve as a potential therapeutic of SARS-CoV-2 and can be worthy of further experimental investigation.

Owing to the alterable aperture, controlled drug uptake/release, biocompatibility, biodegradation, and so forth, the metal–organic frameworks (MOFs) have emerged as an effective nanocarrier for target-specific delivery of therapeutics for different biologically more complex diseases.[27−30] The delivery of the identified potential peptide inhibitors with suitable MOF carriers in pleural empyema will produce more robust action against SARS-CoV-2 because of the availability of more inhibitors and prolonged interactions with the target.[30]