Initial virtual screening

Based on the binding affinity obtained through ADV, 69 of the 104 APD peptides bound stronger to the RBD active region than the ACE2 peptide (Table 1 and Table S2, Supporting Information File 1). Three peptides based on lysozyme were also designed for screening to compare with the APD peptides given the antimicrobial role of lysozyme as part of the innate immune system. The ADV results show that most of the APD peptides successfully docked on the active region of the RBD (Figure 1), suggesting that these APD peptides actually bind to the RBD active region, blocking the entry of SARS-CoV-2 to host cells. Additionally, according to Figure 1, peptides are posed in different ways on the RBD, covering different areas of the active surface of RBD. For instance, MVL (74–87), cyanovirin-N (70–80) and dermaseptin-S4 (Figure S2, Supporting Information File 1) are posed laterally to the RBD active surface. Similar results have been reported previously by Qiao & Olvera, who designed a negatively charged EELE tetrapeptide to neutralize the SARS-CoV-2-RBD–ACE2 binding [25]. It is important to note that the analyzed peptides have a nearly neutral charge, thus they have a low probability of unspecific interactions with other molecules, cellular uptake, or macrophage recognition [42-46].

Table 1: Potential peptide candidates against SARS-CoV-2 obtained by molecular docking. The table shows the physical and biochemical properties of the potential peptides.

Peptide PDB/UNIPROT ID Residues Number of amino acids H bonds/residue Affinity (kcal/mol) alpha basrubrin P83186 1–20 20 0.95 −5.2 human beta defensin 3 1KJ6 27–44 19 0.95 −5.0 sesquin P84868 1–10 10 0.80 −5.6 indolicidin 1G89 1–13 13 0.77 −8.0 GF-17 2L5M 1–17 17 0.76 −5.3 cyanovirin-N (70–80) 2EZM 70–80 10 0.73 −5.3 protegrin 5 2NC7 1–18 18 0.72 −7.2 MVL (94–110) 1ZHS 94–110 17 0.71 −5.0 temporin B 6GIL 1–13 13 0.69 −5.6 dermaseptin-S4 2DD6 1–13 13 0.69 −5.5 MVL (74–87) 1ZHS 74–87 14 0.64 −5.9 MVL (16–34) 1ZHS 16–34 19 0.63 −5.8 ACE2 6VYB 21–44 24 0.63 −4.6 lysozyme (1–20) 1REX 20 20 0.35 −4.9 lysozyme (61–80) 1REX 20 20 0.60 −5.7 lysozyme (111–130) 1REX 20 20 0.25 −5.2

It is important to note that RBD residues from Glu484 to Tyr505, Arg403 to Tyr421, and Tyr449 to Ala475 are involved in the docking of the APD peptides, as was previously reported by Othman and co-workers [6]. Figure 2 shows the mapping of these amino acid residues on the active region of RBD, suggesting that APD peptides are docking on five principal regions of the RBD. The center region of the RBD (Figure S4b and Figure S4c, Supporting Information File 1) is of particular interest because it binds directly to the ACE2 cellular receptor. After this analysis, the number of hydrogen bonds and hydrophobic interactions involved in the APD–RBD complexes was determined by the LigPlot+ software (Table 1) [35]. The numbers of hydrogen bonds per residue and hydrophobic interactions go from 0.25 to 0.95 and from 17 to 31, respectively. The binding affinity values are higher for those APD peptides bound to the RBD with a high number of hydrogen bonds and hydrophobic interactions [47]. Table 1 shows that twelve APD and three lysozyme peptides surpass the binding energy calculated for the ACE2 peptide to the RBD (−4.6 kcal/mol), indicating that these peptides bind to the RBD active region more strongly than the ACE2 peptide.

Figure 2: Peptide candidates (blue) docked to the SARS-CoV-2 RBD (white). (a) ACE2 control peptide (red), (b) cyanovirin-N (70–80), (c) lysozyme (1–20), (d) lysozyme (61–80), (e) lysozyme (111–130), (f) MVL (16–34), (g) P1, (h) P2a, (i) P4a, (j) P6a, (k) P7, (l) P8, (m) P9, (n) P10, (o) P11, (p) P12, (q) P13, (r) P15, (s) P17, (t) PH1, (u) PH2, (v) sesquin, and (w) temporin B.

To validate the docking results for the APD peptides, the crystallized ACE2 peptide was tested using the same ADV parameters. The peptide bound on the RBD active region, and superimposing the docked complex onto the crystallized complex showed a low RMSD of 0.31 Å (Figure 3). Generally, an RMSD value of 2 Å or lower is considered a good docking, thus confirming the validity of the protocol [48].

Figure 3: Superimposition of docked ACE2 (blue) onto the crystallized complex (red) in the active site using PyMOL (RMSD = 0.31 Å).

Proposing new peptides based on hydrogen-bond formation

Since hydrogen bonds play an important role in the formation and stabilization of the protein–ligand complex, novel peptides were designed considering the most common amino acid residues from the 13 APD and lysozyme peptides that bind to the RBD active region through hydrogen bonds [35]. This task is easy to develop in comparison to the typical procedures used in standard peptide design, in which complex algorithms are used to generate a large peptide library [49,50]. In contrast, designing peptides based on hydrogen bond interactions allows one to generate peptides that target specific sites while reducing computation time [51]. Figure 4 shows the most frequent amino acid residues binding to the RBD active region, together with the implicated amino acid residues in the formation of the RBD–ACE2 complex (Figure 4a) [34]. 41 theoretical peptides (denominated HB peptides) composed of 20 amino acid residues were designed, and from these, 23 HB peptides docked to the RBD.

Figure 4: Mapping of the number of hydrogen bonds formed between APD and lysozyme peptide candidates to the SARS-CoV-2 RBD. (a) The graph depicts the residue location and frequency of hydrogen bonds formed with the SARS-CoV-2 RBD. (b) The location of the hydrogen bonds concentrates on the active region of SARS-CoV-2 RBD. The color coding differentiates the frequency of hydrogen bonds formed on each residue with blue being the lowest and red being the highest formation frequency.

In agreement with ADV analysis, the HB peptides docked laterally (left and right) and horizontally to the center region of the RBD. Furthermore, the number of hydrogen bonds formed between the amino acid residues from HB peptides and RBD as well as the binding energy between HB peptides and RBD were increased. These results support the success of the chosen strategy for designing peptides based on the hydrogen bond interactions in an easy way.

Immunogenicity

The immune system can recognize external molecules introduced into our body, which, in many cases, leads to the production of antibodies [12]. Specifically, during a viral infection, viral antigens are presented by MHC I to be recognized by T cells, which, in turn, promote cytosine release and the cytotoxic activity of CD8+ T cells [12,52,53]. Due to the efficiency of this system, many biological therapeutics (proteins, peptides, nucleic acids, and even drugs) do not reach their target since they are eliminated by cells of the immune system, which limits their activity. Therefore, proposed peptides with antiviral activity must be evaluated from the immunological point of view. In this context, an immunogenicity prediction of the proposed peptides (APD, lysozyme, and HB peptides) was developed by the binding of the peptide candidates to MHC I. The final peptide selection was carried out using TepiTool to determine which peptides are capable of evading the immune system, in particular MHC I molecules. The TepiTool platform was used to select peptides of MHC I that do not bind to alleles with an IC50 value below 500 nM since, according to Calis and and Adhikari, binding to alleles with an IC50 above 500 nM would present low to zero immunogenic response [37,38]. From 55 peptides tested, including APD, lysozyme, and HB peptides, only 22 peptides (Table S4, Supporting Information File 1) had a low probability of being recognized by MHC I. This suggests that these peptides can be used to neutralize the SARS-CoV-2 virus without activating the immune system. These peptides, denominated here “OAPs”, optimally attached to the RBD, and their interaction with the RDB is discussed in the following section.

Physicochemical parameters and peptide–RBD interaction

Solubility, net charge, and size are important physical parameters that need to be considered in the design of novel drugs since these play a role in the distribution in the human body and in targeting specific cells, bacteria, viruses, or proteins. Therefore, the physicochemical parameters of the peptides and the peptide conformation after binding to RBD were obtained by using WinHydroPro software. The results are given in Table 2. The peptide net charge, isoelectric point, and water solubility for peptides were determined by the on-line software INNOVAGEN’s peptide calculator (PEPCALC). Almost all OAPs are soluble in aqueous media, independently on their isoelectric point, due to the high ratio between hydrophilic and hydrophobic amino acid residues, except the peptides temporin B and lysozyme (61–80). The net charge calculated for the OAPs varies according to the number of negatively and positively charged amino acids present in the primary structure. The net charge values are in the range of −3.9 to 2.9 at pH 7. The OAPs were selected based on an absolute value of the electrical net charge smaller than 3.9 (|charge| ≤ 3.9, which is the net charge of ACE2), to avoid possible cytotoxic effects [54,55].

Table 2: Summary of the physicochemical properties of the final peptide candidates.

Peptide Residues Number of amino acids Affinity (kcal/mol) H bonds/ residue Water solubility Rg (nm) Net charge at pH 7 Isoelectric point Molecular weight (kDa) ACE2 21–44 24 −4.6 0.63 good 1.64 −3.9 4 2890.07 cyanovirin-N 70–80 10 −5.3 0.80 good 0.83 0.9 8.9 1252.40 lysozyme (1–20) 1–20 20 −4.9 0.35 good 0.98 1.9 9.5 2385.81 lysozyme (61–80) 61–80 20 −5.7 0.60 good 1.07 2.9 8.1 2408.69 lysozyme (111–130) 111–130 20 −5.2 0.25 poor 0.83 1.0 9.9 2179.40 MVL (16–34) 16–34 19 −5.8 0.63 good 0.91 0.1 7.9 1938.11 P1 1–20 20 −6.3 0.65 good 0.90 1.0 10 2394.64 P2a 1–20 20 −4.9 0.50 good 0.92 0.1 5.2 2463.58 P4a 1–20 20 −4.8 0.55 good 0.89 0.1 5.2 2463.58 P6a 1–20 20 −4.6 0.55 good 0.86 0.1 5.2 2463.58 P7 1–20 20 −5.0 0.45 good 0.89 0.1 9.5 2463.58 P8 1–20 20 −5.6 0.40 good 0.90 0.1 7.5 2463.58 P9 1–20 20 −5.2 0.50 good 0.87 0.1 7.5 2521.66 P10 1–20 20 −5.5 0.60 good 0.92 −0.9 7.5 2491.59 P11 1–20 20 −5.6 0.55 good 0.86 −0.9 7.5 2491.59 P12 1–20 20 −5.4 0.65 good 0.93 −0.9 7.5 2541.60 P13 1–20 20 −5.2 0.65 good 0.87 1.1 7.5 2532.69 P15 1–20 20 −4.9 0.55 good 0.87 0.1 7.5 2463.58 P17 1–20 20 −4.6 0.40 good 0.91 0.1 7.5 2463.58 PH1 1–20 20 −5.0 0.55 good 0.92 2.0 11.8 2379.51 PH2 1–20 20 −5.3 0.55 good 0.87 0 6.7 2338.41 sesquin 1–10 10 −5.6 0.80 good 0.70 −1.1 3.9 1157.25 temporin B 1–13 13 −5.6 0.69 poor 0.73 1 10.1 1392.77

Given that, at physiological pH, the RBD active region is positively charged, it could be assumed that negatively charged peptides, such as sesquin and MVL (74–87) (Figure S3, panels 4a,b and 14a,b, Supporting Information File 1), would present a stronger binding to the RBD active region than those peptides with slightly positive charge or neutral charge (Figure S3, Supporting Information File 1) due to electrostatic repulsion [56]. However, as seen from the ADV results, cationic peptides, such as lysozyme (61–80) showed a higher binding affinity than the anionic peptides ACE2 and MVL (74–87). The higher binding affinity observed for positively charged peptides can be explained based on the distribution of the electrically charged patches located on the active surface of the RBD. Figure 1b shows the electrostatic surface potential of the RBD active region, in which negatively charged, neutral, and positively charged patches can be identified, depending on the amino acid residues. Therefore, it can be assumed that electrostatic repulsion forces are negligible. This suggests that the intermolecular interactions (hydrogen bonds and hydrophobic interactions) in the APD peptide–RBD complexes are favored. The secondary structure of the APDs changes and adopts a proper conformation to bind to the RBD protein.

An analysis of the secondary structure for free and docked OAPs was carried out using the PEP-FOLD 3.5 (RPBS Web) web server, while their surface area was analyzed using PyMol. Similar to previous results, the secondary structure of the OAPs changes from α-helices to random-coil conformations, when they docked to the RBD protein, as it is shown in Table 3 [57-59]. The secondary structure for free OAPs consists of α-helices and random coils at different fractions, except for P4a, P6a, temporin B, lysozyme (61–80), PH2, and sesquin, which adopt a fully random-coil conformation. Afterwards, the binding energy of the OAP–RBD complexes, as well as the contact area (Ac) with the RBD active region were determined using docking analysis. The secondary structure and net charge of the OAPs were plotted as functions of the binding energy (Figure 5) with the aim of understanding the relationship between these parameters. Figure 5 shows the secondary structure versus the binding energy of the OAP–RBD complexes (Figure 5a) and is divided into three regions. Region I shows the peptides that have binding energy values similar to the ACE2 binding energy (−4.6 kcal/mol), such as P9 (secondary structure composition 0.75 α-helix and 0.25 random coil) and P12 (fully random-coil conformation). In region II, several peptides present different conformations such as random coils (4), high α-helix-to-random-coil ratio (7), or high-random-coil-to-α-helix ratio (2) conformation. The binding energy in this second region is in the range of −4.8 to −5.6 kcal/mol. The OAPs included in region III are characterized by a high fraction of random-coil secondary structures with binding energies between −5.7 and −6.3 kcal/mol. Additionally, the final surface area (Af) of peptides docked to the RBD increased, indicated by positive values of ∆A = Af − A0, where Af and A0 are the final and initial OAP surface area, respectively (Table 3). In contrast, ACE2 and MVL (16–34) show negative values of ∆A. The observed increase of Af suggests that the OAPs have a large contact area (Ac) with the active surface area, blocking key amino acid residues involved in the association of RBD with ACE2 (Table 3 and Figure 5), as will be shown next. Figure 5b shows that the binding energy values vary independently of the net charge of the OAPs. This can be explained based on the electrostatic surface potential of the active region of the RBD (inset in Figure 5b), which can be divided into three characteristic regions: (i) The upper region is characterized by a negative potential (red ellipse). (ii) The middle region has neutral patches, slightly negative and positive patches (white ellipse); and (iii) the bottom region is characterized by a positive potential (blue ellipse). Recently, it has been reported that the residues Phe486, Tyr489 (located in the upper region), Gln493, Gly496 (located in the middle region), Thr500, and Asn501 (located in the bottom region), are involved in the association of RBD with the ACE2 protein [6,9]. Interestingly, OAPs were attached in different configurations around the active regions of the RBD (inset in Figure 4b), and those principally occupied the middle region of the active surface. These OAPs interact with residues Gln493, Gly496, Thr500, and Asn501, among others amino acidic residues located in the upper and bottom regions. Similar results have been reported by Debmalya and co-workers, who analyzed the potential of chimeric peptides to block the RBD using an in silico approach. They found that a peptide with 26 amino acids binds to the Thr500 and Asn501 residues of the RBD, while a peptide with 23 amino acids binds to the Tyr489 and Thr500 residues of the RBD, and a peptide with 20 amino acid binds to the Gln493 and Asn501 residues of the RBD. However, neither of these peptides was able to block the three regions of the RBD [8]. The results reported herein suggest that OAPs have a great potential as drug inhibitors of SAR-CoV2 and can block the entry of viruses to the cell host through the ACE2 cellular receptor.

Table 3: Secondary structure and contact area of peptides. The secondary structure is presented as the number of amino acids in each structure divided by the total number of amino acids of the peptide. A0 represents the initial surface area, Af corresponds to the final surface area, ∆A is the change in area (a positive value indicates an increase in area and a negative value indicates a decrease in peptide area), and Ac is the contact area of the peptide with the SARS-CoV-2 RBD.

OAP α-Helix Random coil A0 (Å2) Af (Å2) ΔA (Å2) Ac (Å2) ACE2 0.25 0.75 3689.47