Introduction

COVID-19 caused by a member of family Coronaries (CoV) has threatened the survival of human beings on the Earth and it has been declared as global health emergency by the World Health Organization (WHO) [1••]. Coronaviruses (CoV) are known for their ability to cause illness, with severe diseases such as Middle East respiratory syndrome (MERS-CoV) and severe acute respiratory syndrome (SARS-CoV). Though there are other coronaviruses known to infect humans (hCoVsOC43 and 229E), but they are mild pathogens responsible for common cold. In Latin, Corona means “halo” or “crown”; thus, the name represents the structure of the virus which consists of crown like projections on its surface [2]. In 1937, coronavirus was isolated from an infectious bronchitis virus in birds which was responsible to ruin the poultry stocks [3]. Coronaviruses are the type of viruses that directly affect the respiratory tract. These are associated with the common cold, pneumonia, gut, and severe acute respiratory syndrome. Coronaviruses are zoonotic, which means they are transmitted between animals and humans [4]. A new strain, novel coronavirus (nCoV) has come into knowledge since 2019 and has emerged as a threat to mankind.

The first case of Middle East respiratory syndrome (MERS-CoV) was seen in the year 2012, a businessman in Saudi Arabia who died from viral pneumonia [5]. In 2016, a report on 1998 was published by the World Health Organization (WHO) regarding the confirmed cases of MERS-CoV infection and the death rate was approximately 36% (Middle East respiratory coronavirus (MERS-CoV) [5]. The biggest outbreak with first ever confirmed case of this disease came into existence in the year 2015 in South Korea. Including the China, the confirmed cases extend to 186 with total 36 deaths [6, 7]. Cases regarding the novel coronavirus came in to existence among the population of Wuhan, China, on December 8, 2019. Pneumonia was the first symptom of infection and most of the cases were linked to a local fish and animal market. During the research, it was seen that 2019 novel coronavirus was recognized as pathogenic agent responsible for evolution of pneumonia [8]. On January 20, 2020, laboratory in Korea confirmed the first case of coronavirus. On 23 January, 2020, the government of China announced total shutdown of country and advised the people for undergoing personal isolation. In the USA, there are five variants of SARS-Cov-2. B.1.1.7: This variant was discovered for the first time in December 2020 in the USA. It was first discovered in the UK. B.1.351: This variant was discovered for the first time in the USA at the end of January 2021. It was first discovered in December 2020 in South Africa. P.1: In January 2021, this variant was discovered for the first time in the USA. B.1.427 and B.1.429: These two variants were discovered in February 2021 in California (https://www.cdc.gov/coronavirus/2019-ncov/transmission/variant.html).

SARS-CoV-2 consists of four structural proteins: spike (S), membrane (M), envelop (E), and nucleocapsid (N) proteins [9]. Among all, S protein plays an important role in viral attachment, fusion, entry, and also act as a target for development of antibodies, entry inhibitors, and vaccines [10, 11]. The S1 domains of coronaviruses contain receptor-binding domains (RBDs) that directly bind to the cellular receptors [12, 13]. In general, SARS-CoV surface exhibits two components: S1, which contains the receptor binding domain (RBD); and S2, which contains the fusion peptide. SARS-CoV gains entry into cells through interaction of the SARS-SRBD with the cell surface receptor angiotensin-converting enzyme 2 (ACE2) [14, 15]. These interactions are followed by endocytosis, and at the low pH in endosomes, SARS-S is cleaved by a cellular protease called cathepsin L, thereby exposing the S2 domain of the spike protein for membrane fusion [16, 17]. The minimal RBD of SARS-CoV S protein is located in the S1 subunit (AA 318–510) and is responsible for viral binding to host cell receptors [18, 19]. Besides the main receptor for the angiotensin-converting enzyme 2, there are several alternative receptors, such as dendritic cell-specific intercellular adhesion molecule-3-grabbing non-integrin and liver/lymph node-specific intercellular adhesion molecule-3-grabbing integrin [20]. SARS-CoVs recognizes angiotensin-converting enzyme 2 (ACE2) as its receptor, whereas MERS-CoV recognizes dipeptidyl peptidase 4 (DPP4) as its receptor [21, 22]. Two residues (AA 479 and AA 487) in RBD determine SARS progression and tropism, and their mutations may enhance animal-to-human or human-to-human transmission [13]. Some residues (AA 109, 118, 119, 158, 227, 589, and 699) in S protein are critical strategies against this deadly viral agent, especially in high-risk groups, including people of every age group [23]. According to the previous data, the ACE2 receptor expressing cell fused with SARS-S-expressing cells adds to the cell surface by pH-independent mechanism [19]. It enhances the cell stress responses and apoptosis [24]. Binding is very critical for pathogenesis and if the binding of SARS-S with ACE2 receptor is blocked, infection can be stopped. Traditional medicinal plants produce large number of compounds which are used as therapeutics to kill the pathogens [25]. In the recent years, many reports published on antimicrobial activity of the medicinal plants [25,26,27]. It is expected that plant extracts and phytocompounds showing the target site other than antibiotics, a very little information is available on this type of activity of medicinal plants [26, 27]. Extracts of medicinal plants have been used from ancient times and these plants are known for their antiviral properties and less side effects. Traditionally, thyme was acclimated to treat asthma and loosen congestion in the throat and stomach [28]. The pharmacological manuscript of Chailander medical codex (fifteenth and sixteenth centuries) mentions the utilizations of wild thyme for the treatment of headaches caused by cold and laryngitis [29]. During the Renaissance period (sixteenth and seventeenth centuries), wild thyme was utilized internally to treat malaria and epilepsy [30]. Traditionally in many countries, areal part of T. serpyllum is utilized as anthelmintic, a vigorous antiseptic, an antispasmodic, a carminative, deodorant, diaphoretic, disinfectant, expectorant, sedative, and tonic. Thymus serpyllum additionally used to treat respiratory quandaries [29]. In western Balkans, thymus species used to amend blood circulation and as anticholesterolemic, immunostimulant [31]. Carvacrol and thymol are isomers, belonging to the group of monoterpenic phenols with potent antiseptic properties. Chauhan et al. [32] reported thymol (25–200 mg kg−1) as immunomodulatory in cyclosporine A-treated Swiss albino mice by enhancing the expression of cluster of differentiation 4 (CD4), cluster of differentiation 8 (CD8), and Th1 cytokines via upregulation of IFN-4 expression and enhanced secretion of interleukin-12 (IL-12).

Antiviral property of Thymus serpyllum [33] and thymol is already reported [34]. Pilau et al. [35] reported the antiviral activity of carvacrol from Lippia graveolens against human and animal virus (herpes simplex virus, acyclovir-resistant herpes simplex virus 1, bovine herpesvirus 2, respiratory syncytial virus; human rotavirus, bovine viral diarrhea virus). Antiviral nature of Emodin was also reported in several studies [36, 37]. Study from Efferth et al. [38] showed in vitro antiviral properties of artemisinin against hepatitis B virus, hepatitis C virus, and bovine viral diarrhea. Keeping in view the antiviral potential of Himalayan herbs, the current study was focused on the identification of potent phytocompounds from Himalayan herbs (Rheum emodi, Thymus serpyllum, and Artemisia annua) to cure a dangerous COVID-19.

Material and Methods

Bioinformatics Tools

Open Babel GUI [39], UCSF Chimera 1.8.1 [40], Pubchem (www.pubchem.com), RCSB PDB (http://www.rscb.org/pdb), PDBsum (www.ebi.ac.uk/pdbsum), and Autodock/vina software [41, 42] were used in the present investigation.

Ligand Preparation

Four major phytocompounds of three medicinal plants—emodin of Rheum emodi, thymol and carvacrol of Thymus serpyllum, and artemisinin of Artemisia annua—were used for the docking studies. The 3-dimensional structures of all the phytocompounds and chloroquine were obtained from PubChem (www.pubchem.com) in .sdf format. The .sdf file of the phytocompounds was converted into PDB and pdbqt format by using the Open Babel tool [43]. Table 1 shows molecular structure, molecular weight, pharmacological properties, plant source, and percentage of phytocompounds in the respective plants and antimalarial drug chloroquine. Targets of phytocompounds and standard drug chloroquine were predicted by using SwissADME online server.

Table 1 Molecular structure, molecular weight, pharmacological properties, plant source, and percentage of selected phytocompounds and chloroquine
Full size table

Protein Preparation

Two spike proteins of SARS-CoV-2 spike glycoprotein (PDB ID: 6VXX, closed conformation), SARS-CoV-2 spike ectodomain structure (PDB ID: 6VYB, open conformation) [52] and one mutated variant of SARS-CoV-2 B.1.351 (South African variant) variant of Spike glycoprotein (PDB ID: 7NXA) [53] and two receptor of SARS-CoV-2 (Human TMPRSS2 (PDB ID: 7MEQ) [54], Angiotensin-converting enzyme-2 (ACE2 PDB ID: 6M1D)) [55••], and neuropilin-1 (PDB ID: 4DEQ) were used to analyze the interactions of major phytocompounds of R. emodi, T. serpyllum, and A. annua. It has been shown that SARS CoV-2 SB open protein conformation is necessary for binding with ACE2 at host surface; and coronavirus with open surface S-glycoprotein trimers found to be highly pathogenic to human [56••].

The 3-dimensional structures of selected target proteins were retrieved from the Protein Data Bank (PDB) (http://www.rscb.org/pdb). Non-essential water molecules, including heteroatoms, were removed from the target receptor molecule and hydrogen atoms were added to the target receptor molecule. Binding site of both the target proteins of COVID-19 (SARS-CoV-2 spike glycoprotein (PDB ID: 6VXX), SARS-CoV-2 spike ectodomain structure (PDB ID: 6VYB)), SARS-CoV-2 B.1.351 variant Spike glycoprotein (PDB ID: 7NXA), Human TMPRSS2 (7MEQ), Angiotensin-converting enzyme-2, ACE2 (PDB ID: 6M1D), and neuropilin-1 (PDB ID: 4DEQ) was determined by grid box generation. Grid box was generated by adjusting the grid parameter x, y, z coordinate values; grid value for 6VYB and 6VXX was center x: −189.229, y: −255.9, z: 229.87 Å; 7NXA was x: −14.806, y: −19.528, z: −51.972 Å; 7MEQ was x: −1.028, y: −0.352, z: 10.912; and 6MID was x: 126.806, y: 133.196, z:121.533. Size of the grid was same for all the target proteins (i.e., x—40, y—40, z—40 Å) using AutoDock software. The grid values were recorded in the config.txt file format [57••].

Prediction of Drug Likeness of Selected Phytocompounds

The aim of the drug scan was to see whether selected phytochemicals met the drug-likeness criteria. Lipinski’s filters using Molinspiration (http://www.molinspiration.com) were applied for examining drug-likeness attributes, including quantity of hydrogen acceptors (should not be more than 10), quantity of hydrogen donors (should not be more than 5), molecular weight (mass should be more than 500 Daltons), and partition coefficient log P (should not be less than 5). The smiles format of each of the phytochemical was uploaded for the analysis [58].

ADME and Toxicity Prediction of Selected Phytocompounds

Absorption, distribution, metabolism, excretion, and toxicity (ADMET) screening was done to determine the absorption, toxicity, and drug-likeness properties of the selected ligands. The 3-dimensional structures of ligands such as emodin, thymol, carvacrol, artemisinin, and chloroquine were saved in .smiles format and chloroquine was uploaded on SWISSADME (Molecular Modeling Group of the SIB (Swiss Institute of Bioinformatics), Lausanne, Switzerland), admetSAR (Laboratory of Molecular Modeling and Design, Shanghai, China), and PROTOX web servers (Charite University of Medicine, Institute for Physiology, Structural Bioinformatics Group, Berlin, Germany) for ADMET screening. SWISSADME is a web tool used for the prediction of ADME and pharmacokinetic properties of a molecule. The predicted results consist of lipophilicity, water solubility, physicochemical properties, pharmacokinetics, drug likeness, medicinal chemistry, and Brain or Intestinal Estimated permeation method (blood-brain barrier and PGP ± prediction). AdmetSAR provides ADMET profiles for query molecules and can predict about fifty ADMET properties. Toxicity classes are as follows: (i) category I contains compounds with LD50 values ≤ 50 mg kg−1, (ii) category II contains compounds with LD50 values > 50 mg kg−1 but 500 mg kg−1 but 5000 mg kg−1 (Cheng, 2020). PROTOX is a rodent oral toxicity server predicting LD50 value and toxicity class of query molecule. The toxicity classes are as follows: class 1: fatal if swallowed (LD50 ≤5), class 2: fatal if swallowed (55000) [59].

Docking of COVID-19 Receptors and Phytocompounds

The docking of selected ligands to the catalytic triad of protein was performed by using AutoDock/Vina [41]. Docking was performed to obtain populations of conformations and orientation for ligands at binding sites. Docking was performed to study the interactions between SARS-CoV-2 receptors such as 6VXX (closed state), 6VYB (open state), 7NXA (SARS-CoV-2 B.1.351 variant), 7MEQ, 6MID, and neuropilin-1 (4DEQ) with major phytocompounds of R. emodi, T. serpyllum, and A. annua and .pdb file of proteins-ligand complexes was generated. All the bonds of ligands were set to be rotatable. All calculations for protein-ligand docking were performed using the Lamarckian genetic algorithm (LGA) method. The best conformation was chosen with the lowest docked energy, after the completion of docking search. The .pdb complex of protein and ligand was analyzed by using Discovery Studio (https://discover.3ds.com/d) to study the list of interactions between protein and ligand complex. Detailed visualization and comparison of the docked sites of target proteins and ligands were done by using Chimera [60].

Molecular Dynamics (MD) Simulation

In order to further verify the accuracy of docking observations, the two best complexes of S-protein were selected for extensive MD simulation for 50ns. Both the complexes were introduced into Desmond software to study the binding stability of the ligand within the binding site of S-protein [60]. Both complexes were prepared prior to simulation to remove any structural error as described earlier [61, 62]. Both the complexes were solvated in TIP3P water model and 0.15 M NaCl to mimic a physiological ionic concentration. This MD simulation was performed with OPLS3e force field.

Results

Drug Likeness of Selected Phytocompounds

Early preclinical production is supported by drug-likeness filters, which help to prevent expensive late-stage preclinical and clinical failure. The Lipinski rule of 5 was used to evaluate the drug-likeness properties of the emodin, thymol, carvacrol, artemisinin, and chloroquine. The emodin, thymol, carvacrol, and artemisinin followed the Lipinski’s rule of five, whereas chloroquine showed one violation (Table 2).

Table 2 Prediction of drug-likeness activity of selected phytocompounds
Full size table

ADMET Prediction and Toxicity Analysis of Selected Phytocompounds

The comparative ADME properties predicted by SwissADME of emodin, thymol, carvacrol, artemisinin, and chloroquine are summarized in Table 3. Consensus Log Po/w value of ˂5 indicates good aqueous solubility, which means that an adequate amount of drug can reach and be maintained inside the body through oral administration. The emodin, thymol, carvacrol, artemisinin, and chloroquine showed consensus Log Po/w value of < 5 (Table 3). TPSA indicates permeability of compounds into the cells. A TPSA value of < 140 Å2 is required for good permeation of compound into the cell membrane and value < 90 Å2 is required to permeate through blood-brain barrier. All the selected phytocompounds showed TPSA value < 90 Å2, except 94.83 Å2 for emodin, indicating good permeability of selected phytocompounds into the cell as well as through blood-brain barrier. Lipinski’s rule of five helps to determine drug likeness of the compound and orally active drug should not violate the Lipinski’s rule. The emodin, thymol, carvacrol, artemisinin, and chloroquine followed Lipinski’s rule of five (Table 3). The predicted cellular targets of emodin, thymol, carvacrol, artemisinin, and chloroquine are shown in Table 4.

Table 3 ADME properties of selected phytocompounds and chloroquine predicted by SwissADME
Full size table
Table 4 Predicted targets of phytocompounds and standard drug chloroquine
Full size table

Toxicity of emodin, thymol, carvacrol, artemisinin, and chloroquine was predicted by PROTOX-II and admetSAR and results are summarized in Table 5. It was observed that emodin, thymol, carvacrol, artemisinin, and chloroquine are non-carcinogenic and non-cytotoxic in nature and are safe to administer.   However, LD50 value of emodin (5000 mg kg−1) and artemisinin (4228 mg kg−1) calculated from Protox II was higher than that of all other phytocompounds and chloroquine. This suggests that natural phytocompounds are safer even at higher dosage than that of chemically synthesized chloroquine.

Table 5 Toxicity prediction of phytocompounds and chloroquine by admetSAR and PROTOX-II softwares
Full size table

Molecular Docking (MD) Analysis of Phytocompounds and Chloroquine With Spike Protein and Its Variant of SARS-CoV-2

Docking results of phytocompounds of medicinal plants showed good binding affinity and modes of interaction with both the spike protein of SARS-CoV-2 (6VXX closed state and 6VYB open state), SARS-CoV-2 B.1.351 variant, Human TMPRSS2 (7MEQ), Angiotensin-converting enzyme-2, ACE2, and neuropilin-1 as compared to chloroquine. Among all the selected phytocompounds, artemisinin showed the best binding affinity (−10.5 kcal mol−1), followed by thymol (−6.9 kcal mol−1), carvacrol (−6.8 kcal mol−1), emodin (−6.4 kcal mol−1), and chloroquine (−5.6 kcal mol−1) with SARS-CoV-2 spike glycoprotein (6VXX). Artemisinin makes week hydrogen bonds with SER 205, HIS 207, and hydrophobic interactions with ILE 119, VAL 126, ILE 128, PHE 192, PHE 194, ILE 203, LEU 226, and VAL 227 (Table 6, Fig. 1, and Figure S1).

Table 6 E-total of ligands (emodin, thymol, carvacrol, artemisinin, and chloroquine) with targets of SARS-CoV-2using Autodock/vina software
Full size table
Fig. 1
figure 1

Interactions between targeted protein receptor SARS-CoV-2 spike glycoprotein (PDB ID: 6VXX) with artemisinin (A) and emodin (B) using Chimera

Full size image

The binding interaction almost follows the similar pattern of interaction with SARS-CoV-2 spike ectodomain structure (6VYB, open state) with a binding affinity of −10.3 kcal mol−1 (artemisinin), −8.8 kcal mol−1 (emodin), −6.7 kcal mol−1 (thymol), −6.8 kcal mol−1 (carvacrol), and −5.9 kcal mol−1 (chloroquine). Artemisinin makes week hydrogen bonds with SER 730 and THR 778 and showed hydrophobic interactions with TRP 104, ILE 119, ILE 126, VAL 128, PHE 194, and VAL 227 (Table 6 and Fig. 2 and S2). Artemisinin in complex with 6VYB makes week hydrogen bonds with SER 730 and THR 778 and showed hydrophobic interactions with TRP 104, ILE 119, ILE 126, VAL 128, PHE 194, and VAL 227 (Table 5 and Fig. 2 and Figure S2).

Fig. 2
figure 2

Interactions between targeted protein receptor SARS-CoV-2 spike ectodomain structure (PDB ID: 6VYB) with artemisinin (A) and emodin (B) using Chimera

Full size image

In case of SARS-CoV-2 B.1.351 variant, the binding affinity was −6.4, −4.4, −4.7, −5.9, and −4.9 kcal mol−1 for emodin, thymol, carvacrol, artemisinin, and chloroquine respectively. Emodin in complex with SARS-CoV-2 B.1.351 variant (7NXA) makes strong hydrogen bonds with GLU 6, GLN 111, and LYS 207 and showed hydrophobic interactions with SER 7, GLY 9, VAL 92, GLY 112, THR 113, PRO 155, and PRO 208 (Table 6 and Fig. S3). Emodin also showed the best interactions (−7.1 kcal mol−1) with Human TMPRSS2 (7MEQ), followed by artemisinin (−6.9 kcal mol−1), chloroquine (−5.7 kcal mol−1), thymol (−5.5 kcal mol−1), and carvacrol (−5.3 kcal mol−1). Emodin makes hydrogen bonds (moderate strength) with ASN 277 and showed hydrophobic interactions with HIS 274, GLN 276, VAL 278, TRP 306, THR 309, PHE 311, TYR 322, GLN 323, ALA 324, GLY 325, and GLN 327 (Table 6 and Fig. S4). Neuropilin-1 interacts with emodin (−6.6 kcal mol−1), carvacrol and artemisinin (−5.8 kcal mol−1), thymol (−5.6 kcal mol−1), and chloroquine (4.9 kcal mol−1). Emodin makes hydrogen bond (moderate strength) with ILE 147 and hydrophobic interactions with TYR 24, TRP 28, THR 43, ASP 47, THR 76, LYS 78, TYR 80, and GLY 141 (Table 6 and Fig. S5). Similarly, Angiotensin-converting enzyme-2, ACE2 showed best interaction with artemisinin (−7.4 kcal mol−1), followed by emodin (−7.3 kcal mol−1), thymol (−6.9 kcal mol−1), chloroquine (−6.5 kcal mol−1), and carvacrol (−6.1 kcal mol−1). Artemisinin do not make any hydrogen bonds with ACE2, but showed hydrophobic interactions with LEU 545, ALA 273, LEU 495, SER 487, LEU 269, SER 491, ILE 492, ASP 270, TYR 488, VAL 552, PHE 549, and GLU 553 (Table 6 and Figure S6).

The nature of hydrogen bonds was determined on the basis of donor-acceptor distances in protein secondary structure elements. Jeffrey and Jeffrey [63••] categorize hydrogen bonds with donor-acceptor distances of 2.2–2.5 Å as “strong,” 2.5–3.2 Å as “moderate” (mostly electrostatic), and 3.2–4.0 Å as “weak, electrostatic.” Hydrogen bond interaction and hydrophobic interactions of these phytocompounds with target proteins were analyzed through discovery studio and results are summarized in Table 6.

Molecular Dynamics Analysis of Complexes of 6VXX and 6VYB With Artemisinin

Among all the selected phytocompounds, artemisinin showed the best binding affinity with both the spike proteins of SARS-CoV-2 (6VXX closed conformation, 6VYB open confirmation) and these complexes were selected for MD simulation [60, 64, 65]. The binding and conformational stability of the artemisinin complex with the spike receptor protein is a major factor to advocate the inhibitory action of the artemisinin against SARS-CoV-2 infection. MD simulations were carried out for 50ns at 300 Kelvin temperature and 1.01325 bar pressure. Both the complexes of artemisinin with spike receptors (PDB ID: 6VXX, and 6VYB) have 6Na+ ions to neutralize complexes. The artemisinin complexes with 6VXX and 6VYB have water molecules 58601 and 72196, respectively. The Na+ and Cl− ions were added in the environment to these complexes, such that artemisinin +6VXX (Cl− 50.573 nM; Na+ 52.435 nM) and artemisinin+6VYB (Cl− 50.620 nM; Na+ 50.620 nM) could mimic physiological ionic concentration.

The Stereochemical Geometry of S-Protein in Complex With Ligand After Molecular Dynamics Simulation

The Ramachandran mapping of S-protein residues after analyzing the stereochemical geometry of artemisinin+6VXX, and artemisinin+6VYB after MD simulation revealed a very acceptable number of residues in favored region (Fig. 3; Table 7). These complexes possess outlier residues within the acceptable range (less than 1.0%). Overall, spike protein in complex with artemisinin showed a sterically acceptable conformation of molecule after MD simulation, indicating the stability of the complexes.

Fig. 3
figure 3

Ramachandran plots of S-proteins: (A) 6VXX (closed conformation) complex with artemisinin, and (B) 6VYB (open conformation) complex with artemisinin

Full size image
Table 7 The Ramachandran mapping of S-protein residues for analyzing stereochemical geometry of artemisinin complexes with 6VXX and 6VYB after MD simulation
Full size table

Conformational Deviation in Cα of S-Protein (6VXX and 6VYB) in Complex With Artemisinin During Molecular Dynamics Simulation for 50ns

The root mean square deviation (RMSD) of the Cα of the S-protein molecule was analyzed for 50-ns simulation for 6VXX and 6VYB complexed with artemisinin. The Cα RMSD plot for 6VXX (Fig. 4A) gets stabilized at 35ns and remains stabled for the entire simulation period. Ligand RMSD was 4Å at 30ns and then changed to 10Å at 40ns, which is further changed to 14Å at 50-ns time line (Fig. S7). The ligand RMSF plot for artemisinin fit over 6VXX protein showed ligand fluctuation with respect to protein. The high RMSF of atom 7 and 16–18 is mainly due to the exposure to solvent (Fig. 4B). The trajectory analysis reveals conformational shift in artemisinin within binding pocket, which indicates artemisinin is less stable in the binding pocket of 6VXX (closed conformation) complex. It has been proposed that 6VYB (open conformation) is more significant for drug designing. The RMSD plot for 6VYB (Fig. 5C) clearly showed stable interaction between artemisinin and binding pocket residues of 6VYB (open conformation) during the simulation period. Artemisinin showed less fluctuation in open state protein than the closed state protein (Fig. 4C, D).

Fig. 4
figure 4

The RMSD and RMSF plots of S-protein (6VXX and 6VYB) in complex with artemisinin. (A) RMSD of artemisinin + 6VXX complex, (B) RMSF of artemisinin + 6VXX complex, (C) RMSD of artemisinin + 6VYB complex, and (D) RMSF of artemisinin + 6VYB complex as indicated for the backbone and ligand

Full size image
Fig. 5
figure 5

Histogram of ligand contacts with amino acid residues of S-protein. (A) Artemisinin + 6VXX (closed conformation) complex, (B) artemisinin + 6VYB (open conformation) complex. Color codes for hydrogen bonds, hydrophobic, and water bridges interactions are as indicated. X-axis showed amino acid residues and y-axis indicated interaction fraction

Full size image

Furthermore, 6VYBor 6VXXcomplex with artemisinin showed hydrophobic interaction (Trp104, Ile119, Val126, Ile128, Tyr170, Phe192, Phe194, Ile203, Leu226, and Leu229), H-bond interaction (Asn121, Arg192, Ser 205), and water bridges (Asn121, Thr124, Arg192, Ser205, His207) (Fig. 5A, B). There was no salt bridge interaction developed between the 6VYB or 6VXX complex with artemisinin. The observed binding free energy for artimisnin+6VYB complex was −74.54 kcal mol−1, while observed binding free energy for artimisnin+6VXX complex was −55.5 kcal mol−1. The above results clearly advocate the stability of artemisinin in open state protein instead of closed state protein.

Discussion

From the beginning of twenty-first century, three coronaviruses, viz. severe acute respiratory syndrome coronavirus (SARS-CoV) [4], Middle East respiratory syndrome coronavirus (MERS-CoV) [5], and SARS-CoV-2 have crossed the species barrier and resulted in deadly pneumonia in humans [66, 67]. SARS-CoV-2 has caused death of approximately 3.3 million people all around the world. The treatment is symptomatic and oxygen therapy represents the major treatment intervention for patients with severe infection. Mechanical ventilation may be necessary in cases of failure of respiratory, whereas hemodynamic support is essential for managing septic shock. Researchers from Worldwide are continuing to work on developing vaccine againstSARS-CoV-2. Professor Didier Raoult from infectious diseases institute, IHU Méditerranée Infection in Marseille (France), has reported successful results from a new treatment for COVID-19, with early tests suggesting it can stop the virus from being contagious in just 6 days. Chloroquine phosphate and hydroxychloroquine have previously been used to treat coronavirus patients in China, in ongoing COVID-19 clinical trials. Kaletra, a US-based antiviral drug used to treat HIV, is another medicine that is being tested in the fight against the SARS-CoV-2. Emodin has been shown to act as inhibitor of 3a ion channel of SARS-CoV and HCoV-OC43 coronaviruses as well as virus release from HCoV-OC43. It has been reported that emodin is a potent inhibitor of the 3a channel with a K1/2 value of about 20 M. The reduction of extracellular viral RNA copies by emodin reflects inhibition of virus release. At high concentrations of emodin, intracellular levels of viral RNA were reduced suggesting that the high concentrations may also inhibit other stages of the virus life cycle [67]. Ho et al. [68] identified emodin as an effective to block the interaction of the SARS-CoV S protein with the ACE2 and the infection by S protein-pseudo-typed retrovirus. Ahmed et al. [69] reported SARS-CoV-2 spike protein (6VYB) is highly stable protein and it is difficult to un-stabilize the integrity of these proteins by individual drugs. They also reported that inserting of NH2 halogen and vinyl group can increase the binding affinity of coulerpin with 6VYB, while inserting an alkyl group decreases the binding affinity of coulerpin with 6VYB. This work is unique in a way that in silico approach has been utilized to compare the open (6YVB) and closed (6VXX) conformations of spike proteins. It is interesting that open state of the spike protein (6YVB) which is more pathogenic showed more stable interaction with artemisinin as compared to closed state (6VXX) (Data is Fig. 5). Also, artemisinin contacts with amino acid residues of S-protein were different for open and closed conformation (Fig. 5B). Similar to our study, Kumar [60] reported the binding affinity of Nelfinavir (−8.4), Rhein (−8.1), Withanolide D (−7.8), Withaferin A (−7.7), Enoxacin (−7.4), and Aloe-emodin (−7.4) with COVID-19 main protease (6LU7). Rolta et al. [70••] reported the binding affinity of emodin, aloe-emodin, anthrarufin, alizarine, and dantron phytocompound Rheum emodi with three active sites of RNA binding domain of nucleocapsid phospho protein of COVID-19. They reported the binding energies of emodin, aloe-emodin, anthrarufin, alizarine, and dantron were −8.299, −8.508, −8.456, −8.441, and −8.322 Kcal mol−1 respectively with binding site A and −7.714, −6.433, −6.354, −6.598, and −6.99 Kcal mol−1 respectively with binding site B, and −8.299, 8.508, 8.538, 8.841, and 8.322 Kcal mol−1 respectively with binding site C. Similarly, Adem et al. [71] reported khainaoside C, 6-O-Caffeoylarbutin, khainaoside B, khainaoside C, and vitexfolin are potent modulator of open and closed state of SARS-CoV-2 spike S2 proteins. Suravajhala et al. [72] reported the antiviral binding affinity of curcumin with different SARS-CoV-2 Proteins (Spike Glycoprotein-6VYB, nucleocapsid phosphoprotein- 6VYO, membrane glycoprotein-6M17, nsp10-6W4H, and RNA-dependent RNA polymerase- 6M71). Selailia and Chemat [73••] reported hydroxyl chloroquine and artemisinin interact in the same binding pocket of SARS-CoV-2 protein (6LZG); artesunate, artemisinin, and artenimol showed two mode of interaction with LYS 353 and LYS 31; and they also reported the extraction protocol of artemisinin from Artemisia annua. Walls et al. [52] suggested that S-protein is highly pathogenic in human coronaviruses and appears to exist in partially opened states, while S-protein remains largely closed in human coronaviruses that are responsible for common colds. It was also proposed that the S-protein of pathogenic coronaviruses exists in open and close conformation. The current in silico study provides evidence that ligand binding affinity is different for open and closed conformation of S protein and artemisinin interacts more stably with the open conformation of spike protein, than the closed conformation, thus can be used as a potent drug to cure COVID-19. Basu et al. [74••] studied the molecular docking of five phytocompounds (hesperidin, anthraquinone, thein, chrysin, and emodin) with spike protein of SARS-CoV2 and ACE2 receptor. It was shown that hesperidin can bind with ACE2 protein and bound structure of ACE2 protein and spike protein of SARS-CoV2 noncompetitively. The study proposed that the presence of hesperidin, the bound structure of ACE2, and spike protein fragment become unstable. Srivasta et al. [75••] reported the interactions of antimalarial compounds (Mepacrine, Chloroquine, Quinin, Hydroxychloroquine, Artemisinin, Phomarin, and Proguanil) with main protease (PDB ID 6LU7) of SARS-CoV2. They found that mepacrine showed best interactions with 6LU7 (−8.89 kcal mol−1) followed by chloroquine (−8.15 kcal mol−1), quinin (−7.77 kcal mol−1), hydroxychloroquine (−7.62 kcal mol−1), artemisinin (−7.34 kcal mol−1), phomarin (−7.13 kcal mol−1), and proguanil (−6.69 kcal mol−1). Previous studies reported the binding affinity of emodin, artemisinin, and chloroquine against RNA binding domain of nucleocapsid phosphoprotein and main protease of SARS-CoV-2. The current study provides evidence that ligand binding affinity is different for open and closed conformation of S protein. This study also provides the evidence that phytocompounds can inhibit spike protein variant of SARS-CoV-2. Artemisinin interacts more stably with the open conformation of spike protein, than the closed conformation and emodin binds strongly with variant of spike protein SARS-CoV-2; thus, artemisinin and emodin need further attention through in vitro and in vivo studies to be tested to inactivate the SARS-CoV-2.

Conclusions

In this study, we are proposing the artemisinin as a lead phytocompound to inactivate the SARS-CoV-2 virus through inhibiting S-protein, especially in open state conformation. The MD simulation for 50ns showed both the S-protein complexes were stable as there are less than 1.0% outlier amino acid residues in Ramachandran plot developed after MD. The RMSD plot for both complexes and ligand RMSF has shown the stability of artemisinin within the binding pocket of S-protein. The present study proposed a safe and less toxic artemisinin for the treatment for SARS-CoV-2 infection, which can be further validated through in vitro and in vivo studies.