In-Silico Design of a Novel Tridecapeptide Targeting Spike Protein of SARS-CoV-2 Variants of Concern.

Several mutations in severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) have increased the transmission and mortality rate of coronavirus disease-19 (COVID-19) across the globe. Although many vaccines have been developed, a large proportion of the global population remains at high risk of infection. The current study aims to develop an antiviral peptide capable of inhibiting the interaction of SARS-CoV-2 spike protein and its six major variants with the host cell angiotensin-converting enzyme 2 (ACE2) receptor. An in-silico approach was employed to design a therapeutic peptide inhibitor against the receptor-binding domain (RBD) of the spike (S) protein of SARS-CoV-2 and its variants (B.1.1.7, B.1.351, P.1, B.1.617.1, B.1.617.2 and B.1.617.3). The binding specificity and affinity of our designed peptide inhibitor Mod13AApi (YADKYQKQYKDAY) with wild-type S-RBD and its six variants was confirmed by molecular docking using the HPEPDOCK tool, whereas complex stability was determined by the MD simulation study. The physicochemical and ADMET (absorption, distribution, metabolism, excretion, and toxicity) properties of inhibitory peptides were determined using the ExPASy tool and pkCSM server. The docking results and its properties from our in-silico analysis present the Mod13AApi, a promising peptide for the rapid development of anti-coronavirus peptide-based antiviral therapy. Blockage of the binding of the spike protein of SARS-CoV-2 variants with ACE2 in the presence of the therapeutic peptide may prevent deadly SARS-CoV-2 variants entry into host cells. Therefore, the designed inhibitory peptide can be utilized as a promising therapeutic strategy to combat COVID-19, as evident from this in-silico study.

Introduction

Coronavirus disease-19 (COVID-19) is caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), an enveloped, positive-sense single-stranded RNA virus (Machhi et al. 2020). The virus particle is 50–200 nm in diameter (Xu et al. 2020) and contains four structural proteins: the spike (S), envelope (E), membrane (M), and nucleocapsid (N). The N protein encapsulates the viral genome, while S, E, and M create the viral envelope (Bosch et al. 2003). The S-protein is divided into two functional parts, S1 and S2 (V'Kovski et al. 2021). Entry of the virus into the host cell is mediated by the combined action of both subunits, with S1 catalyzing attachment and S2 mediating fusion (Shang et al. 2020; Tang et al. 2020). The S1 subunit is further divided into two functional domains: An N-terminal domain and a C-terminal domain (Shang et al. 2020). A 211-amino acid region (residues 319–529) at the S1 C-terminal comprises a receptor-binding domain (RBD), mediating a pivotal role in the interaction of the virus with the host cell angiotensin-converting enzyme-2 (ACE2) receptor (Bian and Li 2021; Davidson et al. 2020; Ni et al. 2020; Shang et al. 2020).

The ongoing pandemic of SARS-CoV-2 first emerged in late December 2019, is primarily due to its evolving mutations that have increased both the transmissibility and lethality of the virus. The first notable substitution in the S glycoprotein was D614G (Korber et al. 2020; Volz et al. 2021). This mutation increased the infectiousness of the virus and showed higher infectious titers than the wild-type virus (Korber et al. 2020). The World Health Organization classified B.1.1.7 (Alpha), B.1.351 (Beta), P.1 (Gamma), and B.1.617.2 (Delta) as SARS-CoV-2 variants of concern (VOC) strains, and the remaining two, B.1.617.1 (Kappa) and B.1.617.3 (Delta+) as variants of interest (VOI) (Table 1). All of these variants show deletions and substitutions of major residues, but significant substitutions are found in the S-RBD as shown in Table 1, which increases the infectivity of the virus (Shah et al. 2020). Regarding the affinity of the variants with the ACE2, B.1.351 binds ACE2 with a five-fold higher affinity than the original SARS-CoV-2 RBD, while the B.1.1.7 (N501Y) variant binds with a two-fold higher affinity (Ramanathan et al. 2021). Furthermore, it has been shown that individual mutations at K417N and E484K enhance cell fusion more than the N501Y mutant, although there was no change in infectivity (Kim et al. 2021).

Table 1

SARS-CoV-2 sub-lineages (variants) undergone mutation in their Spike protein receptor-binding domain (S-RBD) residues

Sr no.	RBD mutation sites	SARS-CoV-2 Sub-lineages (Variants)
		B.1.1.7	B.1.351	P.1	B.1.617.1	B.1.617.2	B.1.617.3
		(Alpha)	(Beta)	(Gamma)	(Kappa)	(Delta)	(Delta+)
1	K417N		✓				✓
2	K417T			✓
3	L452R				✓	✓	✓
4	T478K					✓
5	E484K		✓	✓
6	E484Q				✓		✓
7	N501Y	✓	✓	✓

Neutralizing antibodies show reduced efficacy against variant strains mainly due to deletion of essential residues in the S-RBD domain, which antibodies recognize as an epitope and substitution of critical residues that enhance the binding of the virus to the receptor (Weisblum et al. 2020; Zhou et al. 2021). Both B.1.1.7 and B.1.351 are resistant to neutralization by most monoclonal antibodies (Wang et al. 2021a, b). Although B.1.1.7 was not found to be resistant to convalescent plasma or sera from vaccinated individuals, B.1.351 was found to be resistant to both (Garcia-Beltran et al. 2021; Wang et al. 2021a, b). Moreover, P.1 is less resistant to a naturally acquired or vaccine-induced antibody response than B.1.351 (Dejnirattisai et al. 2021). The B.1.617.1 variant also evades neutralizing antibodies generated from infection or vaccination (Hoffmann et al. 2021). In terms of individual mutations, N501Y imparts moderate resistance against neutralizing antibodies. At the same time, E484K and K417N are found in B.1.351 and P.1, and L452R and E484Q found in B.1.617 impart a higher resistance to neutralizing antibodies (Wang et al. 2021a, b).

Although food and drug administration (FDA) and world health organization (WHO)-approved vaccines are currently available, they have limitations such as stringent storage requirements, and the question of whether they impart long-term immunization remains unanswered (Garcia-Beltran et al. 2021). Furthermore, neutralizing antibodies generated by the vaccines are evaded by some SARS-CoV-2 variants, as mentioned above (Kuzmina et al. 2021; Li et al. 2021). Another approach is to design antiviral therapy such as small peptides that binds to the RBD of S-protein, inhibiting the binding of the virus to the ACE2 receptor (Baig et al. 2020; Han and Kral 2020). Thus, this study aimed to design a novel peptide that binds to the RBD of the S-protein of variants, inhibiting their interaction with the host cell receptor.

Methods

Structure Retrieval and Modeling of SARS-CoV-2 Variants and Their Interactions with ACE2

The crystal structure of SARS-CoV-2 bound to the human ACE2 complex (Protein Data Bank (PDB) ID: 6M17) was retrieved from the RCSB PDB (https://www.rcsb.org). The SARS-CoV-2 S-protein (referred to as wild-type) and the human ACE2 structure were separately prepared by removing the water molecules, heteroatoms, and other ligand groups, followed by the addition of polar hydrogen atoms in Discovery Studio Visualizer (Dassault Systèmes; 2019). Interacting residues between the ACE2 complex and the S-RBD up to 5 Å were analyzed using UCSF Chimera software (Pettersen et al. 2004). The S-RBD structures for all six variants were modeled using the wild-type S-RBD crystal structure (6M17_Chain E) as a reference, and in silico site-specific mutation of B.1.1.7 (N501Y), B.1.351 (K417N, E484K, and N501Y), P.1 (K417T, E484K, and N501Y), B.1.617.1 (L452R and E484Q), B.1.617.2 (L452R, T478K) and B.1.617.3 (K417N, L452R, and E484Q) were created using the “Rotamers” option in UCSF Chimera software by selecting the conformation of the residues with maximum probability followed by minimization of the structure. All variants of spike protein were then docked against the human ACE2 crystal structure (6M17_Chain B) in pyDockWEB (https://life.bsc.es/pid/pydockweb) (Jimenez-Garcia et al. 2013). As a control, the wild-type strain of the SARS-CoV-2 S-protein was also docked with human ACE2 (6M17_Chain B). To cross-check our results, docking was also performed using the PATCHDOCK (Schneidman-Duhovny et al. 2005), and docked results were analyzed using UCSF Chimera software and Discovery Studio Visualizer (Systèmes 2020).

Designing Peptide Inhibitor (PI) and Molecular Docking Studies

Novel peptide inhibitors (PIs) to block the S-protein of the wild-type strain and six variants were designed with a concept inspired by our previous study (Baig et al. 2020). From the previous study, we took the 18-amino acid PI (18AApi; FLDKFNHEAEDLFYQSSL) and modified it to block the RBD of the wild-type strain and those of variants. To do this, we shortened its length, increased its solubility, and created a more stable PI while keeping the binding specificity intact. The five C-terminal residues of 18AApi identified as having minimal involvement in the binding interaction were sequentially removed, and substitution of residues at positions 8, 9, and 10 resulted in a 13-amino acid PI (13AApi; FLDKFNHNFKDLF) with a similar binding specificity as the wild-type strain. Furthermore, based on its physicochemical properties obtained from ExPASy ProtParam (Gasteiger et al. 2005), the 13AApi was further modified to optimize its hydrophobicity and stability index, and a final modified 13AApi (Mod13AApi) with the sequence “YADKYQKQYKDAY” was obtained with the same binding specificity. The 3D structure of Mod13AApi was modeled in the UCSF Chimera tool using the “build structure” and “energy minimization” options with an Amber ff14SB force field. Prediction of the molecular docking of the 3D structures of these PI with wild-type S-protein was performed in the HPEPDOCK (Zhou et al. 2018) in blind docking mode without providing binding site residues to ensure the highest specificity for the target region. Similarly, this peptide Mod13AApi was docked with all other S-proteins of variants using the same tool. All docking results were analyzed in Chimera and Discovery Studio Visualizer, and final snapshots of their interactions were produced using the Discovery Studio Visualizer tool.

MD Simulation and MM/GBSA Calculation

The docking predicted protein-peptide complex structures of Mod13AApi with S-RBD’s from variants strains were prepared by soaking in a dodecahedron box filled with TIP3P water molecules. Na+ or Cl− ions were added to neutralize the MD simulation system. The parameters for protein, solvent and ions molecules were derived from Amber99SB-ILDN force field (Lindorff-Larsen et al. 2010). Classical MD simulation was performed using GROMACS version 5.0.4 using standard protocol (Abraham et al. 2015). Initially, the simulation system was relaxed using 50,000 steps of steepest-descent energy minimization followed by 1 ns equilibration under NVT and NPT ensemble respectively. Finally, production runs were performed for 500 ns at constant temperature of 300 K. The simulations were performed using time step of 2 fs. Long-range interactions were handled using Particle-Mesh Ewald (PME) method while a cut-off of 10 Å was used for short-range interactions. The MD trajectory was processed and analyzed using MDAnalysis (Michaud-Agrawal et al. 2011), MDTraj (McGibbon et al. 2015), and scikit-learn (Pedregosa et al. 2011) python libraries. MM-PBSA binding energy calculations were performed using the g_mmpbsa program (Rashmi Kumari et al. 2014 ). Graphics were prepared using Matplotlib python library (Hunter 2007) and R (R Core Team 2013).

Physicochemical and Absorption, Distribution, Metabolism, Excretion, and Toxicity Analysis of the PI

Next, we analyzed the physicochemical and absorption, distribution, metabolism, excretion, and toxicity (ADMET) property of Mod13AApi. The physicochemical property of the PI was predicted using the ExPASy ProtParam (ExPASy—ProtParam tool) (Gasteiger et al. 2005), Peptide 2.0 (Peptide Hydrophobicity/Hydrophilicity Analysis Tool (peptide2.com)), and PepCalc (Innovagen) (PepCalc.com—Peptide calculator) (Lear and Cobb 2016) tools. Linear sequences of the Mod13AApi were used for the analysis. Furthermore, we generated the SMILE coordinates of the PIs using the PepSMI tool (PepSMI: Convert Peptide to SMILES string (novoprolabs.com)). The obtained SMILE coordinates were used to predict the ADMET properties using the pkCSM tool (biosig.unimelb.edu.au/pkcsm/) (Pires et al. 2015).

Results and Discussion

In Silico Docking of SARS-CoV-2 S-Protein Variants Show Strong Binding Affinity with the Human ACE2 Receptor

SARS-CoV-2 has now adapted to infect its human host at a high rate and severity. Mutations in its spike (S)-protein, mainly in the RBD at the C-terminal region of the S1 protein (Table 1), are an underlying factor for its lethality (Gómez et al. 2021; Ramanathan et al. 2021; Masaud Shah et al. 2020). We aimed to understand and compare the interaction patterns of these SARS-CoV-2 variant S-proteins, as mentioned in Table 1, with the host cell receptor ACE2. Their binding affinities were analyzed in terms of the change in binding energy as calculated from predictive docking studies. The crystal structure complex of the SARS-CoV-2 S-RBD and the human ACE2 receptor obtained from the PDB database (PDB ID 6M17) was used to determine interacting residues between the complex and each of the individual proteins (S and ACE2) within a 3 Å zone.

The SARS-CoV-2 RBD present on the S1 subunit of the S-protein has a well-characterized interaction with the human ACE2 PD domain. The interacting interface residues of this complex obtained from our analysis within the binding pocket can be used in further docking studies. To comparatively examine the docking results of S-protein variants with ACE2, we used the interaction between the wild-type S-protein and ACE2 protein as a docking control. The top-ranked docking results from two tools were analyzed, and we observed that the docking complex orientation (Fig. 1A) was similar to the solved crystal structure, and the interacting residues in the 3 Å zone were matched, confirming the accuracy of our docking predictions. Then, we examined the S-protein structure of the six variants. The sites of mutation in the RBD region of each variant were studied, and in silico substitution, mutations were created using the wild-type S-protein as a reference. The mutated S-protein of each variant was docked with the human ACE2 protein. Our docking results from pyDockWEB agreed with previous reports suggesting a higher binding affinity of variant S-proteins with ACE2 (Khan et al. 2021; Ramanathan et al. 2021) (Table 2). Recently concluded studies have determined the binding affinity between the human ACE2 and variant S-proteins (B.1.1.7, B.1.351, and P.1) and found that B.1.351 and P.1 binding with ACE2 is significantly stronger, while that of B.1.1.7 is comparable to the wild-type strain (Khan et al. 2021; Ramanathan et al. 2021; Masaud Shah et al. 2020). In addition, these findings have also been validated in in vitro experiments (Ramanathan et al. 2021). Similarly, in our study, the interacting interface between variant S-proteins and ACE2 was determined within the 3 Å zone (Table 3 B), while the docking orientations were observed to be similar to the control (Fig. 1B–G). Interestingly, changes in the binding energy (in kcal/mol) obtained from pyDockWEB were somewhat higher for complexes containing S-protein of the variants when compared to the wild-type complex (Table 2). To avoid technical errors, we cross-checked our docking results using PATCHDOCK (Table 2). Briefly, the increase in the binding energy of complexes containing variant S-proteins docking with the human ACE2 protein suggests that they bind more strongly, and this increase in binding affinity may lead to a higher rate of fusion and entry into host cells, possibly one of the major factors underlying their increased transmission and lethality. In addition, current treatment methods for the wild-type SARS-CoV-2 may not work effectively for these variants, meaning that therapeutic molecules require an adaptive modification. Therefore, with our previous experience of developing an anti-COVID-19 peptide (Baig et al. 2020), we sought to develop a novel PI that specifically blocks the interaction of variant S-proteins with the host ACE2 receptor.

Fig. 1

Docking of Spike protein receptor-binding domain (S-RBD) of SARS-CoV-2 (wild type and six variants) with ACE2 in pyDockWEB: Spike RBD is colored as Blue while ACE2 domain is colored as Red. The residues showing interaction in ACE2 are highlighted in blue, while those showing interaction in Spike RBD are highlighted in Green. All the interacting residues in both domains are named in three amino acid codes followed by residue number. A interaction of Wild-type SARS-CoV-2 S-RBD (PDB ID: 6M17_ChainE) with ACE2 (PDB ID: 6M17_ChainB). B interaction of B.1.1.7 Spike RBD (or alpha) (Mutated at N501Y; indicated by arrow) with ACE-2 (PDB ID: 6M17_ChainB). C interaction of B.1.351 Spike RBD (or beta) (Mutated at K417N, E484K, and N501Y; indicated by arrow) with ACE2 (PDB ID: 6M17_ChainB). D interaction of P.1 Spike RBD (or gamma) (Mutated at K417T, E484K, and N501Y; indicated by arrow) with ACE2 (PDB ID: 6M17_ChainB). E interaction of B.1.617.1 Spike RBD (or kappa) (Mutated at L452R, E484Q; indicated by arrow) with ACE2 (PDB ID: 6M17_ChainB). F interaction of B.1.617.2 Spike RBD (or delta) (Mutated at L452R, T478K; indicated by arrow) with ACE2 (PDB ID: 6M17_ChainB), and G interaction of B.1.617.3 Spike RBD (or delta +) (Mutated at K417N, L452R, E484Q; indicated by arrow) with ACE2 (PDB ID: 6M17_ChainB). ACE2 angiotensin-converting enzyme-2, RBD receptor-binding domain

Table 2

Binding energy of protein–protein complexes of SARS-CoV-2 Spike protein receptor-binding domain (S-RBD) and its variants with human ACE2 receptor

Sr no.	Docking proteins			Docking score
	Host receptor	SARS-CoV-2	Mutation site(s)	PyDock	PATCHDOCK
A	ACE2	WT Spike RBD	Wild type	− 76.090	12,468
B	ACE2	B.1.1.7 Spike RBD	N501Y	− 78.401	12,240
C	ACE2	B.1.351 Spike RBD	K417N, E484K & N501Y	− 88.087	13,834
D	ACE2	P.1 Spike RBD	K417T, E484K & N501Y	− 76.552	12,490
E	ACE2	B.1.617.1 Spike RBD	L452R & E484Q	− 89.996	12,586
F	ACE2	B.1.617.2 Spike RBD	L452R & T478K	− 79.346	12,526
G	ACE2	B.1.617.3 Spike RBD	K417N, L452R & E484Q	− 86.518	12,792

Table 3

Interacting residues between wild-type and variant strains of SARS-CoV-2 Spike protein receptor-binding domain (S-RBD) and human ACE2 analyzed from UCSF Chimera tool. (A) Interacting residues between the crystal structure complex of wild-type (WT) spike protein RBD and ACE2, (B) The 3 Å interacting residues between the pyDockWEB docking complex of Spike-RBD and ACE2

A
Sr no.	Crystal structure complex(PDB Id)	Interacting residues
		SARS-CoV-2 spike	Human ACE2
1	ACE2 & SARS-CoV-2 WT Spike (PDB 6M17)	R403 K417 G446 Y449 Y453 L455 F456 Y473 A475 G476 G485 F486 N487 Y489 Q493 Y495 G496 Q498 T500 N501 G502 V503 Y505	Q24 T27 F28 D30 K31 H34 E35 E37 D38 Y41 Q42 L45 L79 M82 Y83 T324 N330 K353 G354 D355 R357 R393

Novel Peptide Inhibitor (PI) Binds to SARS-CoV-2 S-RBD of Variants and Inhibit Their Interaction with ACE2

A novel method to estimate the effectiveness of a peptide before in vitro and in vivo studies is via the in-silico design of peptides along with the prediction of their efficacy on their target protein using computational tools (Schutz et al. 2020). Designing peptides that can inhibit the binding of the spike receptor-binding domain (S-RBD) to ACE2 has been identified as a novel strategy to inhibit the spread of the virus (Jaiswal and Kumar 2020; Schutz et al. 2020). Other studies have also acknowledged the effectiveness of in silico-designed peptides to inhibit the S-RBD (Baig et al. 2020; Han and Kral 2020), and these have been subsequently verified in in vitro models (Zhang et al. 2020). Peptide-derived inhibitors have several advantages, which include high specificity, tolerability, and safety (Craik et al. 2013; Sharma et al. 2014; Sun 2013). Although the wild-type S-RBD has been targeted by various studies, this study designed a novel peptide inhibiting wild-type SARS-CoV-2 along with the variants, which are the main cause of increased infectivity as these strains have an increased survival advantage (Shah et al. 2020) and high binding affinity with the host. The 18AApi used in this study was taken as a template from our previous study (Baig et al. 2020), and here, we further increased its effectiveness. Although the 18AApi was effective in binding with the S-RBD, the instability index computed using the ProtParam tool on the ExPASy was shown to be 68.65, while the maximum limit of the instability index for a peptide should be 40 to classify the peptide as stable. The stability of a peptide determines its efficacy largely due to the maintenance of its secondary structure (Lee et al. 2019; Zapadka et al. 2017). Unstable peptides lose their secondary conformations and form aggregates, hampering their efficacy and inducing toxicity due to aggregation (Zapadka et al. 2017). The length of the peptide also plays a significant role in maintaining its stability; studies have suggested that short peptides are more stable than longer peptides (Zapadka et al. 2017). Thus, to decrease the length of the peptide while increasing its stability, 18AApi was truncated, removing its last five amino acids, and modifications were added at three positions (E7N, A8F, and E10K), after which the peptide was docked against the S-RBD domain. The resultant 13AApi showed an instability index of -2.21, indicating that it is highly stable. Other vital factors for PIs are their hydrophobicity and half-life. A high hydrophobicity can lead to the self-association of peptides, resulting in decreased solubility in the solvent and decreased efficacy (Zapadka et al. 2017). Half-life is another important parameter that determines the bioavailability of peptides, with a short half-life leading to less bioavailability and less therapeutic efficacy (Sharma et al. 2014). Shortening the peptide from 18 amino acids to 13 amino acids increased the stability of the peptide, but the hydrophobicity of 13AApi increased to 46% compared to 38% for 18AApi. To decrease the hydrophobicity and increase the half-life of the peptide, the amino acids in 13AApi were replaced with a similar group of amino acids possessing similar properties based on their hydrophobicity index. These changes in 13AApi resulted in the design of a novel Mod13AApi with a hydrophobicity score of 15% (compared to 46% for 13AApi) and a half-life of 2.8 h (compared to 1.1 h for 13AApi). Mod13AApi also showed an instability index of 26.8, suggesting that it is more stable than 18AApi. Finally, we performed in silico docking of Mod13AApi with the wild type SARS-CoV-2 S-RBD in HPEPDOCK tool (Fig. 2), and interacting residues within the 3 Å region were analyzed between the protein and the peptide, as shown in Table 4. The binding energy obtained from the tool is provided in Table 5. Our results show that along with better physicochemical properties after modification, Mod13AApi, like 18AApi and 13AApi, maintains specific binding to the S-RBD. Therefore, we next aim to analyze the comparative binding of Mod13AApi with S-RBD of variant strains using in silico docking tools.

Fig. 2

Docking of Wild-type SARS-CoV-2 Spike protein receptor-binding domain (S-RBD) with Mod13AApi. S-RBD is shown in Blue while Mod13AApi is shown in Red. The interacting residues of S1 with peptide are shown in green along with 3-letter code amino acid and residue number

Table 4

Interacting residues between wild-type strain of SARS-CoV-2 Spike protein receptor-binding domain (S-RBD) and designed peptide inhibitor (pi) Mod13AApi

Sr no.	Docking complex	Interacting residues in 3 Å region
		Residues in Spike RBD	Residues in peptide
1	SARS-CoV-2 WT Spike RBD & Mod13AApi	K417, F456, Y473, N487, Y489, F490	Y1, A2, Y5, Q6, Q8, Y9

#The 3 Å interacting residues between docking complex of wild-type S-RBD and peptide inhibitor Mod13AApi are analyzed and obtained from UCSF Chimera

WT wild-type, Mod13AApi modified 13 amino acid peptide inhibitor

Table 5

Docking score of wild-type SARS-CoV-2 and its variants spike protein receptor-binding domain (S-RBD) with Mod13AApi in HPEPDOCK

Sr no.	Docking complex	Docking score (kcal/mol)
		HPEPDOCK
1	SARS-CoV-2 WT Spike RBD & Mod13AApi	− 178.44
2	SARS-CoV-2 B.1.1.7 Spike RBD & Mod13AApi	− 168.19
3	SARS-CoV-2 B.1.351 Spike RBD & Mod13AApi	− 182.47
4	SARS-CoV-2 P.1 Spike RBD & Mod13AApi	− 204.11
5	SARS-CoV-2 B.1.617 .1 Spike RBD & Mod13AApi	− 194.77
6	SARS-CoV-2 B.1.617.2 Spike RBD & Mod13AApi	− 195.18
7	SARS-CoV-2 B.1.617.3 Spike RBD & Mod13AApi	− 194.23

Docking of Mod13AApi to Variants of the SARS-CoV-2 S-RBD

In Fig. 3, we showed the docking of Mod13AApi with S-RBD’s from variants strains. Currently, VOC strains are the main cause of the spread of the disease due to enhanced binding with the ACE2 receptor (Khan et al. 2021; Ramanathan et al. 2021; Shah et al. 2020a) and their ability to resist neutralizing antibodies (Weisblum et al. 2020; D. Zhou et al. 2021). Hence, inhibiting VOC strains may be a potential solution to curb the spread of COVID-19. Like S-RBD of the wild-type SARS-CoV-2, the Mod13AApi was docked with the S-RBD of variant strains using the HPEPDOCK (Zhou et al. 2018) tool. Interestingly, the results indicate that the peptide showed effective binding with all the variant strains (Fig. 3A–F) and bound at the same groove where the S-RBD binds to ACE2 structure (Fig. 1B-G), however, with variable dock scores as indicated in Table 5. Though the B.1.1.7 S-RBD showed less dock score with Mod13AApi, the other five strains S-RBD showed better dock scores with Mod13AApi compared to wild-type SARS-CoV-2 S-RBD. (Table 5). Also, the interacting residues in the 3 Å region between the peptide and the protein were analyzed using Chimera (Table 6). The results indicated that Mod13AApi interacts with many of the residues of the S-RBD involved in the interaction with ACE2. Thus, the peptide might successfully inhibit the interaction of the S-RBD and variants with ACE2. Further, to check the stability of binding of the PI to the wild type and the variants of SARS-CoV-2 S-RBD, we performed the molecular dynamic (MD) simulation study.

Fig. 3

Docking of SARS-CoV-2 Spike receptor-binding domain (S-RBD) of variants with Mod13AApi. A Docking of B.1.1.7 S-RBD with Mod13AApi. B Docking of B.1.351 S-RBD with Mod13AApi. C Docking of P.1 S-RBD with Mod13AApi. D Docking of B.1.617.1 S-RBD with Mod13AApi. E Docking of B.1.617.2 S-RBD with Mod13AApi. F Docking of B.1.617.3 S-RBD with Mod13AApi. B.1.1.7, B.1.351, P.1, B.1.617.1, B.1.617.2, and B.1.617.3 S-RBDs are colored in Violet and Mod13AApi is in Red. The 3 Å interacting residue in S-RBDs with the Mod13AApi are colored in Green along with a 3-letter amino acid code and residue number

Table 6

Interacting residues between wild-type SARS-CoV-2 and its variants Spike protein receptor-binding domain (S-RBD) with Mod13AApi

Sr no.	Docking complex	Interacting residues in 3 Å region
		Residues in Spike RBD	Residues in Mod13AApi
1	SARS-CoV-2 WT Spike RBD & Mod13AApi	K417, F456, Y473, N487, Y489, F490	Y1, A2, Y5, Q6, Q8, Y9
2	SARS-CoV-2 B.1.1.7 Spike RBD & Mod13AApi	R346, Y351, A352, Y449, N450, R466, E484, L492. Q493	Y1, Y5, Q8, Y9, Y13
3	SARS-CoV-2 B.1.351 Spike RBD & Mod13AApi	R346, N354, N450, R466, T470, F490, L492	Y1, Y5, Q8, D11, Y13
4	SARS-CoV-2 P.1 Spike RBD & Mod13AApi	R403, T417, D420, Y421, Y449, Q493, Y495, Q498, Y501, Y505	Y1, Y5, Q8, Y9, A12, Y13
5	SARS-CoV-2 B.1.617.1 Spike RBD & Mod13AApi	L455, F456, Q484, G485, F486, Y489	Y1, Y5, Q6
6	SARS-CoV-2 B.1.617.2 Spike RBD & Mod13AApi	Y453, F456, Y473, N487, Y489, Y505	A2, Y5, Q8, Y9, A12, Y13
7	SARS-CoV-2 B.1.617.3 Spike RBD & Mod13AApi	L455, F456, Q484, G485, F486	Y1, Y5, Q6

3 Å interacting residues between Wild-type (WT) SARS-CoV-2, B.1.1.7, B.1.351, P.1, B.1.617.1, B.1.617.2, and B.1.617.3 S-RBD and Mod13AApi as shown in Fig. 3 are analyzed in UCSF Chimera tool

Molecular Dynamic (MD) Simulation of WT SARS-CoV-2 and Its Variants Spike Protein Receptor-Binding Domain (S-RBD) with Mod13AApi

The docking complexes of WT SARS-CoV-2 and its variants S-RBD with Mod13AApi were next subjected to MD simulation to determine the stability of the protein-peptide complexes. MD simulation of 500 ns was performed for each of the complex using the GROMACS program. Analysis of root mean squared deviations (RMSD) throughout the MD trajectory suggested no significant fluctuations in the ACE2 binding interface of the SARS-CoV-2 and its variants S-RBD (Fig. 4). Moreover, no significant changes in the radius of gyration and fraction of native contacts were observed, suggesting an overall stability of the protein-peptide complexes (Fig. 5). The Mod13AApi peptide initially showed significant fluctuations from the docking predicted binding mode, which more or less stabilized after 400 ns (Fig. 4) for most of the SARS-CoV-2 variants S-RBD. Slight conformational and positional rearrangement might have accounted for most of the deviations and fluctuations as Mod13AApi remained bound to WT SARS-CoV-2 and its variants S-RBD for the entire 500 ns MD trajectory (Fig. 4). To further estimate the energetics of Mod13AApi binding to WT SARS-CoV-2 and its variants S-RBD, MM-PBSA analysis was performed using full 500 ns MD trajectory which revealed favorable binding energies (Table 7). However, contrary to docking studies, MM-PBSA binding energy calculation revealed stronger binding of Mod13AApi to B.1.1.7 S-RBD when compared with WT SARS-CoV-2 S-RBD. Predicted binding affinity of Mod13AApi to B.1.351 S-RBD and B.1.617.3 S-RBD was also more or less similar to WT SARS-CoV-2 S-RBD. Our MD simulation study suggests that further optimization of Mod13AApi may be necessary to improve its binding to SARS-CoV-2 S-RBD variants.

Fig. 4

Molecular dynamics simulation to study the dynamics of Mod13AApi in complex with WT SARS-CoV-2 and its variants S-RBD. Root mean squared deviation (RMSD) of Spike protein, Mod13AApi and protein-peptide interface of A WT with Mod13AApi. B B.1.1.7 S-RBD with Mod13AApi. C B.1.351 S-RBD with Mod13AApi. D P.1 S-RBD with Mod13AApi. E B.1.617.1 S-RBD with Mod13AApi. F B.1.617.2 S-RBD with Mod13AApi. G B.1.617.3 S-RBD with Mod13AApi

Fig. 5

A Radius of gyration and B fraction of native contacts calculated throughout the 500 ns MD trajectory

Table 7

Mod13AApi binding energies calculated using the MM-PBSA approach. The MD trajectory of the full 500 ns was used for MM-PBSA binding energy calculations

Sr no.	Energy type	Energy (kJ/mol)
		WT Spike RBD	B.1.1.7 Spike RBD	B.1.351 Spike RBD	P.1 Spike RBD	B.1.617.1 Spike RBD	B.1.617.2 Spike RBD	B.1.617.3 Spike RBD
1	Van der Waal (∆Evdw)	− 236.81 ± 53.55	− 232.83 ± 23.26	− 149.70 ± 45.15	− 150.83 ± 53.43	− 222.47 ± 32.16	− 233.32 ± 31.66	− 190.19 ± 43.95
2	Electrostatic (∆Eelec)	− 320.22 ± 122.43	− 472.05 ± 57.18	− 405.72 ± 154.06	− 184.51 ± 166.84	− 54.80 ± 109.88	− 260.21 ± 82.75	− 314.11 ± 116.59
3	Polar solvation energy (∆Gpolar)	438.40 ± 119.39	558.54 ± 63.03	452.35 ± 161.0	332.49 ± 172.11	268.52 ± 97.75	445.53 ± 91.45	410.63 ± 110.13
4	Nonpolar energy (∆Gnonpolar)	− 28.0 ± 5.14	− 28.49 ± 2.31	− 20.23 ± 4.95	− 19.13 ± 6.49	− 24.38 ± 3.18	− 27.69 ± 3.27	− 24.34 ± 4.80
5	Binding energy (∆Gbind)	− 146.62 ± 76.73	− 174.84 ± 50.60	− 123.29 ± 77.01	− 21.98 ± 71.55	− 33.13 ± 60.48	− 75.70 ± 50.92	− 118.00 ± 67.70

Physicochemical and ADMET Analysis

Knowledge of the physicochemical properties of any therapeutic molecule is crucial before in vitro, in vivo, or preclinical applications. These properties can aid many steps, from synthesizing the molecule to determining its chemical nature and solubility factor. While designing the peptide inhibitor (PI) for SARS-CoV-2 S-RBD, these factors also played a key role. After analyzing the binding specificity of Mod13AApi with the WT Spike-RBD as well as with the variants, we analyzed the detailed physicochemical properties of this peptide and performed a predictive analysis of their ADMET properties in the host (Table 8). Briefly, protein or peptide solubility depends on several factors, including amino acid composition, hydrophobic and hydrophilic properties, polar and non-polar groups, molecular weight, and pH of the medium. In addition, the isoelectric point (pI) of the peptide also helps in predicting the solubility of the peptide in an aqueous medium. Protein-water interactions increase at pH values significantly above or below the pI value as the protein carries positive and negative charges, whereas minimal interactions are observed when the pH is close to the pI value. Similarly, peptide displays a high solubility and dissolve in an aqueous medium when they have > 25% charged residues and < 25% hydrophobic residues, provided that their polar and non-polar residues are evenly distributed along the sequence. The predicted values of the molecular weight, pI, and net charge at a neutral pH will help determine the condition of the solubility of the aqueous medium. In contrast, in our results, the significantly lower hydrophobicity value of the peptide (15.38%), the higher charge percentage, and the regular distribution of acidic and basic residues ensure the high solubility of the PI as listed in Table 8. In terms of stability and half-life, Mod13AApi is predicted to be stable with a half-life of 2.8 h. The instability index obtained from the ProtParam tool indicates that a numeric value of 40 or below for any peptide is considered highly stable; the instability index of Mod13AApi was 26.8, thus suggesting a stable structure.

Table 8

In-silico analysis of physicochemical and ADMET (absorption, distribution, metabolism, excretion, and toxicity) properties of inhibitory peptide Mod13AApi

Sr no.			Peptide properties	Values
1	Physicochemical analysis	Properties	Molecular weight (g/mol)	1683.84
			Theoretical pI	8.3
			Net charge at pH 7	0.9
			Negative + positive residues	2 + 3
			Molecular formula	C78H110N18O24
			Number of atoms	230
			Instability index	26.8
			Half-life (in h)	2.8
			Hydrophobicity (%)	15.38
			Acidic + basic + neutral ratio (%)	15.38 + 23.08 + 46.15
			GRAVY	− 2.1
2	ADMET Analysis	Absorption	Water solubility (log mol/L)	− 2.892
			Skin perm (log Kp)	− 2.735
			Intestinal absorption (%)	0
		Distribution	Fraction unbound (Human) (Fu)	0.345
			BBB permeability (logBB)	− 3.029
			CNS permeability (logPS)	− 7.622
		Metabolism	Cytochrome P substrate	No
			Cytochrome P inhibitor	No
		Excretion	Total clearance (log ml/min/kg)	− 0.588
		Toxicity	AMES toxicity	No
			Skin sensitivity	No
			MRTD human (log mg/kg/day)	0.438
			Rat oral LD50 (mol/kg)	2.482

pI-Isoelectric point; Instability Index- value below 40 classifies stable protein/peptide; GRAVY- grand average of Hydropathy; Water solubility (logS)- defines solubility in water at 25 °C; Skin permeability- logkp > − 2.5 classifies low skin permeability; Fraction unbound- defines unbound state in plasma protein remaining for pharmacological action; BBB permeability- logBB < − 1 classifies poorly distributed to the brain; CNS permeability- logPS > − 2 classifies CNS penetration and logPS < − 3 classifies no CNS penetration; Total clearance- includes both hepatic and renal clearance; MRTD- maximum recommended tolerated dose should be less than 0.477 mg/kg/day

On the other hand, to further understand how these PI behave in a biological system, we performed an in-silico prediction of its ADMET property, as shown in Table 8. ADMET properties define the absorption in the body in terms of water solubility and skin permeability, and our results suggest that our PI has good water solubility and low skin permeability values. Furthermore, the PI showed a below threshold value for crossing the blood–brain barrier or having central nervous system permeability. In addition, the PI is neither an inhibitor nor substrate of cytochrome P. In terms of toxicity, our prediction suggests that the PI is safe and displays no skin toxicity or carcinogenic toxicity. We also predicted the maximum recommended tolerated dose for humans, which was below the maximum set value of 0.477 mg/kg/day (Table 8).

Conclusion

In this study, we designed a novel peptide that could effectively inhibit the interaction of the wild-type SARS-CoV-2 as well as its variants (B.1.1.7 or Alpha, B.1.351 or Beta, P.1 or Gamma, B.1.617.1 or Kappa, B.1.617.2 or Delta, and B.1.617.3 or Delta +) Spike protein with the host cell ACE2 receptor, thereby blocking the cellular entry of SARS-CoV-2 and its variants. The peptide inhibitor (Mod13AApi) designed in our study displayed high specificity and affinity toward their target spike proteins and bears good physicochemical and ADMET properties. This peptide has several advantages, including economical and ease of synthesis due to its short length, low toxicity, and high target specificity and selectivity. This study indicates that designed PI could be further tested and developed as selective therapeutic candidates for COVID-19.