FormalPara Key Points

1. COVID-19 is an ongoing pandemic caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2).

2. SARS-CoV-2 spike protein interacts with the angiotensin-converting enzyme 2 (ACE2) receptor present on the surface of the host cells.

3. A novel peptide has been designed to inhibit SARS-CoV-2 S-glycoprotein interaction with ACE2, thereby blocking the cellular entry of the virus.

1 Introduction

The recent severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) outbreak has posed a great challenge to human health. During the past 2 decades, we have encountered the outbreak of many deadly viruses, such as Ebola [1], Zika [2] and Nipah [3, 4], as well as the evolution of various strains of coronaviruses (CoV), mainly SARS-CoV [5] and MERS-CoV [6], which resulted in high morbidity and mortality. After almost 100 years of the deadly influenza virus (H1N1, or Spanish flu) pandemic, with millions of deaths worldwide (approximately 40 million) [7, 8], the recent outbreak of a novel CoV, or SARS-CoV-2, has left the entire world in helplessness and misery. The clinical spectrum of COVID-19 ranges from mild fever, cough and shortness of breath, to severe clinical conditions characterized by respiratory failure [9, 10]. Old age, together with pre-existing conditions such as lung or heart disease, diabetes, or a compromised immune system, expedite the infection time and severity [11, 12]. Multiple recent reviews [13, 14] could brief-up the statistical dynamics of COVID-19 cases worldwide.

Structurally, the CoV has the largest known RNA genome (26–32 Kb) among other known viruses, characterized by non-segmented, positive-sense, single-stranded RNA. This genome encodes four major structural proteins of the virus, including the nucleocapsid (N), envelope (E), membrane (M), and spike (S) proteins [15, 16]. The membrane and envelope proteins are associated with virus assembly, while the spike protein plays the main role in facilitating virus entry via mediating its interaction with the transmembrane surface receptor on the host cells [15, 17]. The spike protein directly interacts with the peptidase domain (PD) of the angiotensin-converting enzyme 2 (ACE2) receptor [18, 19], which technically marks the virus entry inside the cells [20]. The SARS-CoV-2 shares around 80% sequence identity with the SARS-CoV genome, suggesting similarity in their host interacting functions [21, 22]. Like SARS-CoV, the spike protein of SARS-CoV-2 contains S1 and S2 subunits, which are jointly responsible for fusion and entry of the virus [23,24,25] inside the host cells. The receptor-binding domain (RBD) in the S1 subunit initiates direct binding with the ACE2 PD, whereas the S2 subunit contains basic elements needed for the membrane fusion [19, 20, 23, 24].

ACE2, a single-pass type I transmembrane metallocarboxypeptidase enzyme is primarily involved in the maturation of peptide hormone angiotensin (Ang), which in turn regulates the vasoconstriction and blood pressure [19, 26, 27]. ACE2 is primarily expressed in alveolar epithelial type II cells and serves as a viral receptor [28, 29]. Besides alveolar epithelial type II cells, it is also expressed in several extrapulmonary tissues, including the heart, kidney and intestine [26, 30]. The full-length structure of ACE2 consists of two main domains—the N-terminal PD and the collectrin-like domain (CLD) at the C-terminal end [19, 30,31,32]. In fact, the spike glycoprotein of SARS-CoV-2 binds to the homodimer of ACE2, which facilitates virus entry into the host cells [19, 33]. In addition to this, studies have been conducted in association with ACE2 and the amino acid transporter B0AT1 (or Slc6a19), and how SARS-CoV-2 may bind to the ACE2-B0AT1 complex [19, 34]. ACE2 interaction with B0AT1 could aid in producing antivirals or a vaccine that can block CoV infection by targeting ACE2 [35,36,37,38].

With the current epidemiology of SARS-CoV-2, a vaccine might be considered a highly anticipated therapy. However, given that vaccine development and production is a highly challenging and time-consuming task, the need of the hour is to develop potent therapeutic agents that could effectively curb the infection in the early stages. Several approaches, such as decoy-soluble ACE2 proteins, antibodies from the serum of infected patients, epitope-based vaccines, repurposing of drugs, and designing blocking peptides are underway [39,40,41,42,43,44,45,46]. Peptides possess several attractive features when compared with small molecules and protein therapeutics, including high structural compatibility with target proteins, the ability to disrupt protein–protein interfaces, etc.

Computational tools ease the way to reach a therapeutic solution for COVID-19. Currently, the computational analysis of structural differences in human ACE2 impact its binding to the SARS‐CoV‐2 spike protein, which thereby lays a foundation for the design and development of ACE2-based peptide inhibitors of SARS-CoV-2 [47,48,49]. In the current study, we used computational biology tools to develop a therapeutic strategy utilizing a novel peptide by exploiting the SARS-CoV-2/ACE2 interaction.

2 Material and Methods

2.1 Three-Dimensional Structural Investigation

The protein sequence of SARS-CoV-2 was retrieved from the RCSB Protein Data Bank (PDB ID: 6M17) for peptide design and development. Initially, a stretch of 23 amino acid residues (Glu23 to Leu45) was taken from ACE PD (PDB 6M17) and saved as pdb using the ‘build structure’ module in Discovery studio visualizer (BIOVIA, DS, 2019) [50]. This sequence was further optimized for geometry using the ‘minimize’ functionality in the UCSF Chimera software [51].

2.2 Alanine Scanning and Peptide Design

We performed alanine scanning to better understand the individual roles of the residues within the 23-amino acid peptide in binding with SARS-CoV-2. Alanine substitution of each residue was performed using BIOVIA Discovery Studio 2020 software [50]. Each peptide was then modeled and energy minimized using UCSF Chimera software [51]. We also truncated five residues from the N-terminal end of the 23-amino acid (23aa) peptide to get a shorter 18-amino acids peptide (Phe28 to Leu45). This peptide was further energy minimized using the method described above.

2.3 Docking and Molecular Dynamics (MD) Simulation Studies

The crystal structure of the SARS-CoV-2 spike protein-ACE2 (PDB ID: 6M17) complex was retrieved from the Protein Data Bank (https://www.rcsb.org/) and used as the starting point for the docking studies. Molecular docking was performed using PyDock (https://life.bsc.es/pid/pydockweb) [52], HADDOCK 2.4 [53], and ZDOCK [54] servers, to validate the results for each program, which in turn should be mutually correlated. The docking results were analyzed using the Chimera [51] and DS Visualizer programs [50]. The results obtained were analyzed for binding energies and peptide conformations in the SARS-CoV-2 spike protein binding interface. Molecular dynamics (MD) simulations were performed at 50 ns for the modeled protein. The coordinates of the structure were stored and examined utilizing the analytical tool in the GROMACS 4.5.5 package [55].

2.4 Physiochemical Properties of the Peptide

The toxicity and other relevant molecular mechanistic values of the peptides, such as electrostatics, desolation, and Van der Waals (VdW) force, were predicted using the ‘ToxinPred’ server (https://crdd.osdd.net/raghava/toxinpred/) [56]. ToxinPred is an in-silico method that has been developed to predict and design toxic/non-toxic peptides. The main dataset used in this method consists of 1805 toxic peptides (≤ 35 residues). The physiochemical characteristics of the peptide sequences were determined using the ProtParam tool (https://web.expasy.org/protparam/) of the ExPASy database server [57] (Fig. 1).

Fig. 1
figure 1

ACE2 is an entry receptor for SARS-CoV-2. SARS-CoV-2 severe acute respiratory syndrome coronavirus 2, ACE2 angiotensin-converting enzyme 2

Full size image

3 Results and Discussion

3.1 The N-Terminal Region of the Angiotensin-Converting Enzyme 2 Peptidase Domain is Critical for Binding to the Severe Acute Respiratory Syndrome Coronavirus 2 Spike Protein

The crystal structure of the SARS-CoV-2 spike protein with the ACE2 PD domain (PDB ID: 6M17) was retrieved from the Protein Data Bank (https://www.rcsb.org/). The interface residues between the SARS-CoV-2 spike protein and the ACE2 PD domain were visualized and interpreted using UCSF Chimera Software [51]. After a detailed analysis of interface residues, a small stretch of the ACE2 PD N-terminal region (23-amino acids: Glu23 to Leu45) was found to be interacting majorly with the SARS-CoV-2 spike protein (Fig. 2 and Table 1).

Fig. 2
figure 2

SARS-CoV-2 interaction with the ACE2 PD domain (6M17). The projected view displays the interacting residues at the interface site. SARS-CoV-2 severe acute respiratory syndrome coronavirus 2, ACE2 angiotensin-converting enzyme 2, PD peptidase domain

Full size image
Table 1 Interacting residues at the interface of the SARS-CoV-2 and ACE2 PD domain
Full size table

Furthermore, we evaluated whether the 23-amino acid chain alone and without the remainder of the ACE2 PD domain could bind to the SARS-CoV-2 spike protein. Hence, the 23aa sequence was taken from the ACE2 PD domain and docked to the SARS-CoV-2 spike protein (6M17) using the PyDock program [52]. To perform a non-biased analysis, we performed a blind docking run whereby we did not specify the binding site during the docking simulations. The obtained results were clustered and analyzed by comparing the binding energies and docked conformations of the peptide within the SARS-CoV-2 spike protein (Table 2). Interestingly, we found that the peptide binds to the SARS-CoV-2 spike protein at exactly the same place where the original ACE2 PD domain interacts, which shows that the 23-amino acid sequence independently has the potential to inhibit the interaction of the SARS-CoV-2 spike protein and ACE2 complex. After the docking analysis, the SARS-CoV-2/ACE2 peptide interface was determined and the critical interacting amino acids were identified in the ACE2 PD domain-derived peptide involved in binding to the SARS-CoV-2 spike protein (Table 2).

3.2 Computational Alanine Scanning Analysis of the Peptide

Computational alanine (A) scanning was performed to identify the critically important amino acids of the 23aa peptide inhibitor involved in binding to the SARS-CoV-2 spike protein. The approach of alanine scanning contributes to estimating the binding energy of each residue to the total binding energy of the original peptide [58]. The process involves the substitution of the amino acid with alanine and records the resulting binding energy changes by molecular docking of the substituted three-dimensional structure of the peptide. If substitution by alanine results in a significant drop in overall free binding energy, those residue(s) are considered the critically important residues for binding [59]. Hence, we performed alanine scanning of the 23aa peptide to determine the significance of each residue having the inhibitory potential of independently interacting with the ACE2 binding region of the SRAS-CoV-2 spike protein. All amino acids in the peptide were independently mutated to alanine using BIOVIA Discovery Studio 2020 software [50]. Alanine-substituted peptides were then modeled and minimized using the UCSF Chimera program [51]. Next, we performed docking studies of each mutated peptide and analyzed their binding energies and three-dimensional conformations in the ACE2 binding interface of the SARS-CoV-2 spike protein. As seen in Fig. 3, we observed that the mutated residues E23A, Q24A, A25A, K26A, T27A, L29A, D30A, K31A, N33A, E35A, A36A, E37A, D38A, L39A, Q42A, S43A, and S44A have no substantial effect on complex stability in terms of the change in total binding energy. The total binding energy of the 23aa peptide with the SARS-CoV-2 spike protein (− 55.873 kcal/mol) is comparable to the binding energy of the mutated residues (E23A, Q24A, A25A, K26A, T27A, L29A, D30A, K31A, N33A, E35A, A36A, E37A, D38A, L39A, Q42A, S43A, S44A), ranging from − 56 to − 53 kcal/mol (Fig. 3). However, the mutated residues F28A, F32A, H34A, F40A, Y41A, and L45A were found to be critically important for binding and stabilizing the peptide-SARS-CoV-2 spike protein complex (Fig. 3). We observed the change from − 55.873 kcal/mol to − 42.996 drop in the total binding score (Fig. 3). Furthermore, we deleted the non-significant (identified through alanine scanning) residues and found the deletion of any residue from F28 to L45 would make an unstable structure. Hence, we removed five residues (E23, Q24, A25, K26, T27) from the N-terminal region of the peptide and retained residues from F28 to L45, resulting in the final 18aa peptide. In addition to the molecular mechanistic values such as electrostatics and desolation, we also calculated the van der Waals (vdW) of all individual peptides (electronic supplementary Table 1). vdW is usually a weak force, a type of intermolecular force enhancing the attraction and binding of a ligand to its receptor. VdW forces are involved in the maintenance and stabilization of a drug-receptor complex.

Fig. 3
figure 3

Identification of critically important residues in the 23aa peptide by alanine scanning, and the binding energies of each peptide docked with SARS-CoV-2 spike protein predicted by pyDockWEB. The alanine scanning results provide the 18aa new peptide inhibitor with similar binding efficiency

Full size image

Replacing the respective amino acids (F28, F32, H34, F40, Y41, L45) with alanine significantly reduces the binding score, either by directly losing the interaction or by changing the peptide conformation, which becomes less compatible with the SARS-CoV-2 spike protein. Hence, based on the results of alanine scanning, the peptide was shortened from the extreme N-terminus due to its dispensability in SARS-CoV-2 binding. The final peptide derivative (F28, L29, D30, K31, F32, N33, H34, E35, A36, E37, D38, L39, F40, Y41, Q42, S43, S44, L45) was further docked and analyzed to validate its binding efficiency. We utilized three state-of-the-art programs (pyDock [37], ZDOCK [39] and HADDOCK 2.4 [38]) to validate our docking studies. The top poses retrieved from each software package resulted in conformations very close to the original 23-amino acid peptide, thereby validating its efficacy as a potent SARS-CoV-2 inhibitor (Fig. 4).

Table 2 Interacting residues (1) at the interface of SARS-CoV-2 and the docked 18 aa peptide; and (2) at the interface of SARS-CoV-2 and the docked 23 aa peptide
Full size table
Fig. 4
figure 4

SARS-CoV-2 (Cyan) interaction with the 18aa-derived peptide (blue) along with the original peptide of 23-amino acid (red). Dockings were performed using different software packages: a pyDock; b ZDOCK; and c HADDOCK 2.4. SARS-CoV-2 severe acute respiratory syndrome coronavirus 2

Full size image

3.3 Physiochemical Properties of the Peptide

The final 18aa peptide, with a molecular weight of 2203.39, has an amino acid composition of Ala (A) 5.6%, Asn (N) 5.6%, Asp (D) 11.1%, Gln (Q) 5.6%, Glu (E) 11.1%, His (H) 5.6%, Leu (L) 16.7%, Lys (K) 5.6%, Phe (F) 16.7%, Ser (S) 11.1% and Tyr (Y) 5.6%. The total number of negatively charged residues was four (Asp and Glu) and the total number of positively charged residues was two (Lys and His). The extinction coefficient calculated in units of M−1 cm−1, at 280 nm measured in water, was found to be 1490 ± 10. It is useful to have an estimation of this coefficient for spectrophotometrically following a protein when purifying it. The half-life is a prediction of the time it takes for half of the amount of protein in a cell to disappear after its synthesis in the cell. The presently used ExPASy tool, ProtParam, relies on the ‘N-end rule’, which relates the half-life of a protein to the identity of its N-terminal residue; the prediction is given for three model organisms (human, yeast and Escherichia coli). Considering the N-terminal of the sequence considered is F (Phe), the estimated half-life is 1.1 h (mammalian reticulocytes, in vitro). The instability index (II) is computed to be 68.65, which classifies the protein as unstable. The physical stability of peptides is influenced by both intrinsic and external factors. Prediction of peptide stability suggests precautionary steps be taken to make the therapeutically potent unstable peptides stable through certain biochemical analysis. The predicted aliphatic index of the peptide was found to be 70.56, with a grand average of hydropathicity (GRAVY) of − 0.522. ToxinPred was used to predict the toxicity of the 18aa peptide with desired toxicity by mutating the minimum number of amino acids [40]. The ‘ToxinPred’ tool allows us to submit query peptide in the FASTA format and to optimize the peptide sequence to obtain the maximum/minimum/desired toxicity based on the quantitative matrix-based position-specific scores by comparing them with the toxic/non-toxic data set. The in-silico predicted ‘ToxinPred’ toxicity results of the 18aa peptide indicate that it is non-toxic compared with the mutated peptides.

3.4 MD Simulations

The respective docked complexes of the 23aa and 18aa peptides with target protein were subjected to MD simulation studies for analysing their stability in terms of root mean square deviation (RSMD) and potential interactions at 50 ns. The stability of these two protein–peptide complexes was analyzed using the RMSD of backbone atoms. The calculated average RMSD of these complexes was in the range of 0.94–1.48 nm and 0.24–0.8 nm for the 23aa and 18aa peptide complexes, respectively (Fig. 5).

Fig. 5
figure 5

a Plots of backbone RMSD of the 23aa and 18aa peptide, shown in magenta and green, respectively. b The 23aa and 18aa peptides complexed with SARS-CoV-2. RMSD root mean square deviation, SARS-CoV-2 severe acute respiratory syndrome coronavirus 2

Full size image

By the end of the 50 ns time-scale simulation, both of the protein–peptide complexes attained their stable conformations (Fig. 5). The lower RMSD values obtained after comparing the initial and final peptide complex conformations (Fig. 5b) revealed that the docked complexes were stable and hence there was not much difference in their conformations during the course of dynamics simulation.

4 Conclusions

Several advantages of peptides, including ease of synthesis and modifications, low toxicity, and high target specificity and selectivity, led us to design a potential therapeutic candidate against COVID-19. The present in-silico study is the first step towards designing inhibitory peptides to block SARS-CoV-2 entry into the host cell. The peptide could inhibit binding of the SARS-CoV-2 spike protein with the ACE2 PD domain, which is the earliest stage of SARS-CoV-2 infection. The findings of this study suggest that a peptide of 18-amino acids could be considered as a selective therapeutic candidate for COVID-19. The study conducted by Han and Kral shows similar aspects towards designing a peptide inhibitor [49]. The identification of 15 interacting residues of an ACE2 PD with an RBD of the SARS-CoV-2 spike protein at a 3 Å region in the crystal structure (PDB 6M17) is in agreement with the findings of Han and Kral [49]. However, in comparison with our structurally stable 18aa peptide inhibitor comprising the ACE2 PD residues ranging from positions 28–45, their inhibitor 1 peptide comprising residues 21–55 is structurally unstable and non-specific; hence, they made multiple modifications to it, which thereby changes its properties as well as its similarity to our peptide. Out of a total of 15 interacting residues of ACE2 PD, 10 residues are in the stretch of Q24 to L45. Utilizing this fact in our study, each residue in the 23–45 stretch (23aa peptide) of the ACE2 PD was estimated for functional significance and structural stability based on their binding energy and conformational stability. In contrast to inhibitor 1, we obtained a stable 18aa peptide inhibitor where non-significant residues were truncated and the non-interacting residues within the 18aa peptide at the interfacial site were found important, both for binding and stability of the peptide. Hence, inhibitor 1 from Han and Kral differs in structural stability, composition and specificity due to the presence of non-significant residues leading to deformity in its helical structure. Provision of experimental investigations of the desired synthesized peptide in combination with our computational studies will be needed to confirm and complete this study. It is recommended to compare our computational results with the designed peptides in vitro and in vivo, and eventually to clinical studies.