Repurposing the natural compounds as potential therapeutic agents for COVID-19 based on the molecular docking study of the main protease and the receptor-binding domain of spike protein.

Severe acute respiratory syndrome coronavirus (SARS-CoV-2) enters the cell by interacting with the human angiotensin-converting enzyme 2 (ACE2) receptor through the receptor-binding domain (RBD) of spike (S) protein. In the cell, the viral 3-chymotrypsin-like cysteine protease (3CLpro) enzyme is essential for its life cycle and controls coronavirus replication. Therefore, the S-RBD and 3CLpro are hot targets for drug discovery against SARS-CoV-2. This study was to identify repurposing drugs using in silico screening, docking, and molecular dynamics simulation. The study identified bentiamine, folic acid, benfotiamine, and vitamin B12 against the RBD of S protein and bentiamine, folic acid, fursultiamine, and riboflavin to 3CLpro. The strong and stable binding of these safe and cheap vitamins at the important residues (R403, K417, Y449, Y453, N501, and Y505) in the S-protein-ACE2 interface and 3CLpro binding site residues especially active site residues (His 41 and Cys 145), indicating that they could be valuable repurpose drugs for inhibiting SARS-CoV-2 entry into the host and replication.

Introduction

The ongoing COVID-19 pandemic has resulted in over 179,255,976 affected individuals and 3,882,080 deaths by the end of 21 June 2021 (World Health Organization) [1]. Therefore, it is vital to control and prevention this spreadable disease.

Coronaviruses (CoVs) belong to a group of positive single-stranded (+ ss) RNA viruses that are classified in the family Coronaviridae [2]. The coronaviridae family includes four genera (α, β, γ, and δ) and are enveloped viruses [2, 3]. The SARS-COV-2 genome comprises approximately 30 kb nucleotides, which contain 10 open reading frames (ORFs). The 5′ terminal region contains two large replicase ORFs, ORF1a, and ORF1b, which encode two viral polyproteins, pp1a and pp1ab. The autocatalytically processing of polyproteins pp1a and pp1ab results in the production of 16 nonstructural proteins (nsp1 to nsp16) (Fig. 1A). Among them, nsp5 (3CLpro enzyme) is indispensable to the viral replication and infection process, therefore is considered an interesting target for the development of potential inhibitors against COVID-19 [4-6].

Fig. 1

Schematic representation of the SARS-CoV-2 structure and its genomic organization. A Schematic representation of the genome sequence of SARS-CoV-2 and proteome showing different polyproteins (pp1a and pp1b) along with the structural and accessory proteins. B Structure of 3CLpro with catalytic residues and S-RBD, which represented the ACE2-binding residues. C Diagram of coronavirus structure showing M (membrane) protein, S (Spike) protein, E (envelope) protein, N (nucleocapsid) protein, and RNA along with the ACE2 receptor [7]

The 3-D structure of the 3CLpro enzyme, also called the main protease (Mpro), contains three domains. Domain I (residues 1–100) and domain II (residues 101–183) have a two-β-barrel fold, which is alike to chymotrypsin, whereas the third domain (residues 198–303) consists of five α-helical structures. The substrate-binding cleft, located between domain I and II, with a catalytic dyad of His41 and Cys145, in which the cysteine thiol allow the nucleophilic attacks and His acts as a proton acceptor (Fig. 1B). Previous studies have shown other residues Thr25, Met49, Phe140. Gly143, His163, Met165, Glu166, His172, and Gln189, which can also interact with ligands [6, 8].

The 3′ terminal of a coronavirus genome encodes structural viral proteins: nucleocapsid (N) protein, membrane (M) protein, envelope (E) protein, and spike (S) protein (Fig. 1C). The spike protein localizes on the virion surface and consists of the ectodomain region (ED), intracellular domain, and TM region [9]. The ED region (S1 + S2) (aa 13–1273) includes the S1 subunit (aa 13–685) mediating SARS-CoV attachment to the host angiotensin-converting enzyme 2 receptor via receptor-binding domain (aa 319–541) and S2 chain (aa 686–1273), which serving host and viral membrane fusion, as a result, potentiate the CoV to release its RNA genome in host cell [9, 10]. Several studies have recently been clarifying the residues involved in the interaction between RBD and human ACE2. The key residues of spike RBD interact with ACE2 receptors are as follows: Arg403, Glu406, Gly446, Lys417, Tyr449, Tyr453, Ala475, Asn487, Gln493, Tyr495, Gly496, Phe497, Gln498, Thr500, Asn501, and Tyr505 [11, 12] Among the above residues, Lys 417, Tyr 449, Tyr 489, Gln 493, Asn 501, and Tyr 505 play a crucial role in anchoring RBD to ACE2 receptor (Fig. 1B) [12]. It therefore appears designing drugs for inhibiting the interaction of the S protein with its receptor to be also another attractive strategy for controlling of SARS-CoV-2.

In this study to find the potent drug molecules to inhibit the main protease and also prevent the spike protein interaction with the host receptor, virtual screening was performed against vitamins available in Selleckchem Inc. (WA, USA) and subsequently, molecular dynamics simulation was carried out on selected ligands. Bentiamine, folic acid, riboflavin, fursultiamine, benfotiamine, calcipotriene, cocarboxylase, and ergosterol are discovered as a potent inhibitor against 3CLpro and vitamin B12 and bentiamine, folic acid, and benfotiamine were found to inhibit RBD of the spike protein.

Methods

Receptor and ligand preparation

The crystal structures of receptor-binding domain/RBD of the spike protein (PDB ID.6M0J_E) and 3CL-protease (PDB ID. 6LU7_A) were obtained from the Protein Data Bank [13] and after cleaning with Discovery Studio 4.1[14], minimized and changed to pdbqt format using MGLTools [15].

Forty-seven FDA-approved vitamins were collected from the Selleckchem Inc. website (https://www.selleckchem.com/). The names of vitamins used as ligands in this study are listed in the Supplementary Information. The 3D structure of the ligands were retrieved from the PubChem database [16] in SDF file format and converted to pdb and pdbqt format using Avogadro [17] and MGLTools, respectively.

Molecular docking and post-docking analysis

To compute the binding score between protein and ligands, AutoDock vina in PyRx 0.8 [18] was used to perform the docking-based virtual screening over 47 candidate compounds against the 3-D structure of 3CLpro and S-RBD proteins. For 3CLpro inhibition calculation, the grid box was set at 26 A° × 26 A° × 26 A° (x, y, and z) and center − 15.518 A° × 21.151 A° × 66.865 A° (x, y, and z) with a grid point spacing of one angstrom. The grid box was set into the His41, Cys145 (catalytic dyad) at the ligand-binding site.

For the RBD domain of S-protein inhibition calculation, the default parameters for the grid box were set to 36 A° × 55 A° × 32 A° (x, y, and z) and center − 33.688 A° × 30.917 A° × 7.143 A° (x, y, and z) to cover the possible binding sites of the SARS-CoV-2 S protein with the ACE-2 human receptor, in agreement with previous results [12].

The ligands with the highest binding scores—i.e., the most negative binding energies—were selected for more analysis. The protein–ligand complexes were visualized by Discovery Studio 4.1. The hydrophobic interactions and hydrogen bonds of the selected docked complexes were analyzed by Discovery Studio and LigPlot + (v 1.4.5) [19].

Molecular dynamic simulation

In order to assess the stability of the protein–ligand complexes, molecular dynamics (MD) simulation was performed on docked complexes using GROMACS package 2019.6 [20]. GROMOS96 54A7 [21] and Prodrg [22] force fields were used to create proper coordinate and topology files for proteins and ligands, respectively. The aqueous environment was created using a simple point charge water model in a cubic box with a distance of 0.8 nm from the box to the surface of the protein. After being neutralized, the steepest descent algorithm was utilized for energy minimization. For each simulation, 50,000 steps of energy minimization were performed. After that, the equilibrations of the system were performed under NVT up to 100 ps at 300 K with restraint forces of 1000 kJ/mol, followed by 100 ps under NPT at the pressure of 1 bar and using restraint forces of 1000 kJ/mol with modified Berendsen thermostat and Parinello–Rahman barostat algorithms [23], respectively. The electrostatic interactions were measured with the particle-mesh Ewald (PME) method [24]. Finally, the MD run was carried out with no restraint for 100 ns on ligand-receptor complexes. All simulations were repeated three times. GROMACS in-built tools were applied to analyze the MD trajectories to calculate the root-mean-square deviation (RMSD), root-mean-square fluctuation (RMSF), and Xmgrace [25] was used for plotting graphs. The molecular dynamic was done on a dimeric form of 3CLpro, which was retrieved using PyMOL [26].

In another attempt for comparing the variation of RMSD, the RMSD standard deviation (RMSD-SD) was calculated for each residue through MD simulation among the three repeated MD trajectories.

Result and discussion

Virtual screening of vitamins against the key viral proteins and analysis

Many of the current pieces of research focused on repurposing FDA-approved drugs for anti-SARS drug development. In this work, 47 approved FDA vitamins were virtually docked against 3CL-protease and S-RBD of coronavirus 2.

As shown in Table 1, the binding energy of folic acid (vitamin B9 or vitamin M), bentiamine (dibenzoyl thiamine), riboflavin (vitamin B2), fursultiamine, ergosterol, calcipotriene (calcipotriol), cocarboxylase (thiamine pyrophosphate hydrochloride), benfotiamine, and vitamin B1 to 3CLpro (PDBid: 6lu7) exhibit the appropriate score.

Table 1

Top screened vitamins against main protease (section I) and RBD domain (section II) of Spike protein

Name	Pubchem ID	Docking score (Kcal/mol)	Protein–ligand interactions
			HP	HB
Section I
Bentiamine	3,036,235	− 7.9	H41, L141, N142, H164, M165, and Q189	G143, S144, C145, and E166
Folic acid	135,398,658	− 7.9	T24, T26, H41, S46, M49, G143, C145, H164, M165, Q189	T25, T45, and E166
Riboflavin	493,570	− 7.7	H41, Met49, F140, E166, M165, and Q189	L141, N142, G143, S144, and C145
Fursultiamine	3,002,119	− 7.0	T26, T25, H41, M49, N142, H164, M165, D187, and Q189	T26, L141, G143, S144, and C145
Ergosterol	444,679	− 7.7	M49, N142, G143, C145, M165, E166, R188, Q189, and T190	T26
Calcipotriene	5,288,783	− 7.2	T25, H41, M49, L141, N142, C145, E166, and Q189	T26 and G143
Cocarboxylase	9068	− 7.1	H41, L141, N142, G143, S144, M165, E166, R188, and Q189	Y54, C145, H164, and D187
Benfotiamine	3,032,771	− 7.1	T25, L27, H41, M49, F140, L141, N142, C145, M165, E166, D184, R188, and Q189	Y54 and G143
Vitamin B1	6042	− 6.0	M49, F140, N142, G143, C145, H163, H164, M165, E166, and Q189	H41 and L141
Section II
Vitamin B12	5,311,498	− 7.6	K417, L455, F456, Tyr473, Tyr489, Gln493, S494, Y495, Q498, and Y505	R403, Y449, and G496
Folic acid	135,398,658	− 7.2	R403, Y453, Y495, N501, and Y505	Q406, Q409, K417 Y449, G496, and Q498
Bentiamine	3,036,235	− 7.0	Y449, Y453, S494, Y495, N501, and Y505	R403, G496, and Q498
Benfotiamine	3,032,771	− 6.2	R403, Y495, Q498, N501, and Y505	Y449, Q493, S494, and G496

HP, hydrophobic interaction; HB, hydrogen bond

The His 41 and Cys 145 are catalytic dyad residues in active site cavity of 3CLpro

For the S-RBD, vitamin B12, folic acid, bentiamine, and benfotiamine displayed the highest affinity (Table 1).

The selected compounds were visually inspected using Discovery Studio to determine their residues, which are involved in interactions with 3CLpro and S-RBD. The results in Table 1 show that all selected ligands interact with the catalytic dyad residues (Cys-145 and His-41) as well as significant interactions with the most ligand-binding residues, id; Thr24, Thr25, Thr26, Met49, Asn142, Gly143, His164, Glu166, and Gln189 in the binding pocket of 3CLpro.

For spike protein, the ligands cover the major part of the RBD-ACE2 interface and interact with important residues in the RBD-ACE-2 interface, such as R403, K417, Y449, Y453, S494, N501, and Y505. Trapping these amino acids is a barrier to interaction between spike RBD and the human ACE2 interface (Table 1).

MD simulations and LigPlot analysis

The stability of the docked structures was examined using the parameters RMSD and RMSF.

The comparison of the standard deviation of backbone RMSD among the three repeated simulations of all the protein–ligand complexes revealed small SD values (< 1A°), indicating that all the simulated protein–ligand complexes were stable during 100 ns simulation (Fig. 2A, B). The RMSD of protein–ligand complexes was also compared with native apoprotein. In the case of 3CLpro, because of the enzymatic nature of the protein (3CLpro), the placement of the ligand in the active site of the enzyme (binding cavity) stabilizes the structure of the protein–ligand complex, as can be seen in Fig. 2A, the complexes are very stable in comparison to apoprotein.

Fig. 2

Analysis of the molecular dynamics simulations for 3CLpro-ligand and S-RBD-ligand complexes. A Root mean square deviation (RMSD) analysis of the protein backbone for the complexes of 3CLpro-with bentiamine, folic acid, riboflavin, and fursultiamine. B S-RBD-with bentiamine, folic acid, benfotiamine, and vitamin B12. C Root mean square fluctuation (RMSF) plot for the protein backbone atoms for 3CLpro-ligand complexes and D S-RBD-ligand complexes

However, in the case of S-RBD protein, there is no significant deviation between the stability of apoprotein and complexes, and both structures have acceptable stability, although in the free state, the protein is somewhat more stable than in complex with ligands (Fig. 2B).

The RMSD plots depicted in Fig. 2A demonstrate 3CLpro-ligand complex structures, except 3CLpro-folic acid reached a steady state during the first few nanoseconds and remained stable throughout the simulation.

The LigPlot analysis of 3CLpro-bentiamine complex in dimeric form (the GROMACS output), exhibits six H bonds in chain A and five H bonds in chain B (Fig. 3 and Table 1). It is evident that key amino acid residues, such as Gly143, Ser144, Cys145, Glu166, and Gln189 produced H-bond interactions with bentiamine in chain A and Gly143, Ser144, Cys145, and Glu166 in chain B. In addition, the amino acids involved in the ligand interaction in both chains and all three repeats are largely the same (complete data not shown) and also match with pre-MD docked interactions (Supplementary Fig. 1). In other words, the ligand was largely constant during the simulation and experienced little fluctuation (Fig. 2C).

Fig. 3

Ligplot analysis of 3CLpro-ligand and S-RBD-ligand interaction. A 2D representation of the hit compounds in the binding pocket of chain A and chain B in the dimeric form of 3CL-protease (6lu7) after molecular dynamics for 3CLpro-bentiamine and 3CLpro-folic acid, 3CLpro-riboflavin, and 3CLpro-fursultiamine. B 2D representation of the hit compounds in the interface area of S-RBD, after molecular dynamics simulation of S-RBD-vitamin B12, folic acid, bentiamine, and benfotiamine. Green lines indicate the hydrogen bonds, and red-dotted lines indicate the hydrophobic interactions (images are drawn for the repeat one of each ligand)

For 3CLpro-folic acid, a very mild variation was observed during simulation. This slight instability may be due to the elongated structure of folic acid, which has to change its conformation to fit into the binding cavity of 3CLpro. Folic acid forms H bonds with Ser144, Cys145, and Gln189 in chain A and Thr24, Asn142, Cys145, Glu166, and Gln189 in chain B of the 3CLpro-folic acid complex after 100 ns MD simulation (Fig. 2B). The residues involved in interaction were almost the same in repeat 1 and repeat 2 and match a lot with residues involved in interaction before MD simulation (Supplementary Fig. 1). The repeat 3 forms two H bonds with Ser144 and Cys145 in chain A and two H bonds with Gly143 and Glu166 in chain B (data not shown). The residues of the 3CLpro receptor, which are involved in the interaction with folic acid in repeat 3 were to some extent different from the ones in repeats 1 and 2.

Riboflavin forms four H bonds with residues Leu141, Gly143, and Glu166 in chain A and five H bonds with Gly143, Ser144, and His164 in chain B of the 3CLpro-riboflavin complex, and fursultiamine forms H bonds with residues Thr26, Gly143, and Asn119 in chain A and Glu166 in chain B of the 3CLpro-fursultiamine complex.

Investigation in SD of complexes 3CLpro-fursultiamine and 3CLpro-riboflavin implied a slight deviation among the three repeats of each complex. This behavior can be attributed to the relatively small size of riboflavin and fursultiamine in comparison to bentiamine and folic acid. The proportion of the size of the ligands with the enzyme cavity leads to very little deviation in the RMSD of repetitions.

Since the structure of S-RBD, which does not have a binding cavity like 3CLpro (at the same time, the S-RBD binding site, is much larger than the 3CLpro cavity), the interaction of S-RBD with ligands illustrated some deviation among the members of repeats in complexes, except for S-RBD-benfotiamine.

In relation to S-RBD-bentiamine, a relatively sharp deviation was seen at about 20 ns for two repeats and the system then reached equilibrium (Fig. 2B). In repeat 1, bentiamine forms two H bonds with Ser494 and Gly496, which is very similar to the docked form applied for MD (Supplementary Fig. 1). Two repeats slightly deviated from the stable state, but there were no considerable deviations.

The simulation of S-RBD-folic acid also revealed the deviation for the three repeats, which could be explained through the special properties of the substrate and the ligand. As mentioned above, on account of the large area of the S-RBD binding site and elongated shape of folic acid, the possibility of this slight deviation thus seems quite logical. Meanwhile, in all three repeats, the ligand form hydrogen and hydrophobic interactions with sensitive and important residues of the receptor (3CLpro). For example in repeat 1, the folic acid forms H bonds with Arg 454 and Asn422 and hydrophobic interactions with Arg403, Glu406, Lys417, Tyr421, Tyr 453, Tyr495, and Tyr505 (Fig. 3), and in repeat 2, H bonds were observed for Gln409 and Tyr453 and hydrophobic interactions for Arg403, Glu406, Gly416, Lys417, Leu455, Tyr495, and Gly496 (data not shown).

Moreover, S-RBD-benfotiamine exhibits a steady increase during the first 20 ns and reached a stable state throughout equilibration. In the case of benfotiamine, due to the small size of the ligand, in all three repeats, the ligand is attached to the upper part of S-RBD where it is a part of the S-RBD-ACE-2 junction. In repeat 1, the benfotiamine forms H bonds with Gln498 and hydrophobic interactions with Gly447, tyr449, Gly496, Phe497, Asn501, and Tyr505 (Fig. 3B) and in repeat3 (which seems more stable than the other two repetitions), H bonds were observed for Tyr449, Gly496, and Gln498 and hydrophobic interactions for Arg403, Try495, Phe497, Asn501, Gly502, and Tyr505 (data not shown).

Additionally, the S-RBD-vitamin B12 complex presented an acceptable equilibration during MD simulation. Vitamin B12 has a large structure, and although it was optimized before MD, its conformational somewhat changed during simulation; hence, its three replications reached an acceptable equilibrium after about 70 ns.

In all three repeats, the ligand had at least one hydrogen bond with the S-RBD and did not leave the protein binding site.

RMSF measures the fluctuations of each residue during simulation. The RMSF plots indicated that each complex demonstrates a similar fluctuation with an apoprotein (Fig. 2C, D). The RMSF plot for each residue of 3CLpro-ligand complexes revealed that the protein did not fluctuat in the 100-ns simulation periods (Fig. 2C). Comparatively, the RMSF plot of S-RBD complexes (Fig. 2D) demonstrated the overall stability of these constructs with the exception of a few residues in C-terminal domains showing greater fluctuations. Investigating S-RBD in 3D structure illustrated that the residues in the loop region fluctuated more during the simulation, but in the presence of bentiamine and folic acid, the fluctuation was even reduced. Finally, for benfotiamine, due to its small size and consequently the interaction with a limited number of receptor residues (Fig. 3B), somewhat more fluctuation was observed compared to other ligands among its three repeats of the RMSF plot (Fig. 2D).

Hence, the present in silico study suggests that the vitamin compounds were able to adapt to an acceptable orientation and H bond and hydrophobic interaction with the important residues of the binding site, viz., Thr25, Met49, Phe140. Gly143, His163, Met165, Glu166, His172, and Gln189 of Mpro and Lys 417, Tyr 449, Tyr 489, Gln 493, Asn 501, and Tyr 505 of S-RBD, respectively.

The above results suggest that vitamins might be the potential 3CLpro and S-RBD inhibitors and could probably be used for treating SARS-CoV-2.

The descriptions of the selected compounds are as follows:

Benfotiamine is a thioester that is a synthetic S-acyl derivative of thiamine (vitamin B1). Benfotiamine may also be beneficial for the treatment and prevention of diabetic nephropathy and type 2 diabetes mellitus [15-Selleckchem].

Bentiamine (dibenzoyl thiamine), a lipophilic derivative of vitamin B (thiamine), is a kind of food additive that can be rapidly absorbed into the body and converted to thiamine (Selleckchem). It appears to break down incompletely to thiamine; therefore, there will always be some of it in our body intact and may be able to inhibit 3CLpro. Bentiamine has been approved as a food additive in Japan [27]. The toxicity, teratogenicity, and mutagenicity studies have shown no side effects for bentiamine [28].

Cocarboxylase is a thiamine (vitamin B1) derivative, which is produced by the enzyme thiamine diphosphokinase (Selleckchem). Thiamine pyrophosphate was found less toxic in comparison to thiamine and thiamine monophosphate. Thiamine pyrophosphate chloride is used in food supplements, as a source of vitamin B1 [27].

Calcipotriene (calcipotriol, MC903) is a synthetic derivative of calcitriol, a form of vitamin D (Selleckchem). Calcipotriene is used as a cream or solution drug in controlling psoriasis and leads to the normalization of epidermal growth [29]. The side effects of this drug are very limited [30].

Folic acid, a B vitamin, plays an important role in cell division and in the synthesis of amino acids and nucleic acids like DNA (Selleckchem). Currently, this vitamin is prescribed as an oral supplement for all women planning, or capable of, pregnancy. Also, due to the function of folate for the detoxification of arsenic through arsenic methylation, folic acid can be used as an effective supplement in promoting health in all age and gender groups [31].

Fursultiamine is a nutritional supplement and vitamin B1 derivative, with potential antineoplastic activity. Fursultiamine can be used for vitamin B1 deficiency [15-Selleckchem].

Ergosterol is a sterol and a biological precursor (a provitamin) to vitamin D2 (Selleckchem). In recent years, the ergosterol-rich extract presented notable antioxidant and antimicrobial properties, besides showing no hepatotoxicity. In fact, ergosterol is a type of plant sterol found in mushrooms. Mushrooms have been widely consumed as food, especially by the Chinese and Japanese [32].

Riboflavin (Vitamin B2) is a heat-stable and water-soluble vitamin found in food and used as a dietary supplement to prevent and treat riboflavin deficiency. This vitamin is essential for healthy skin, nails, and hair [15-Selleckchem].

Vitamin B12 is a water-soluble vitamin with a key role in the normal functioning of the brain and nervous system (Selleckchem). It is better to use natural resources to meet the body’s need for vitamin B12, but doctors sometimes prescribe fortified foods and supplements for people over the age of 50 [33].

Conclusion

This in silico study suggests that FDA-approved vitamins; bentiamine, folic acid, riboflavin, fursultiamine, benfotiamine, ergosterol, cocarboxylase, and calcipotriene could serve as potential inhibitors against the viral Mpro enzyme and vitamin B12, folic acid, bentiamine, and benfotiamine for S-RBD spike protein of COVID-19 with further experimental validation research.

Below is the link to the electronic supplementary material.

Supplementary file1 (DOCX 1108 KB)