Computational drug repurposing studies on the ACE2-Spike (RBD) interface of SARS-CoV-2
The S-glycoprotein (Spike) of the SARS-CoV-2 forms a complex with the human transmembrane protein Angiotensin-Converting Enzyme 2 (ACE2) during infection. It is one of the hot targets, as it forms the first line of contact with the human cell. Drug repurposing would help in identifying drugs that are safe with no or fewer side effects. The Food and Drug Administration (FDA) approved drugs and the phytochemicals from Indian medicinal plants were explored. Molecular docking and simulations of these molecules targeting the ACE2-Spike complex were performed.Rutin DAB10 and swertiapuniside were obtained as the top-ranked drugs. The molecular dynamics simulations of ligand-free, rutin DAB10-bound, and swertiapuniside-bound ACE2-Spike complex revealed crucial ACE2-Spike interface residues forming strong interactions with these ligand molecules. Inferring that they may affect the ACE2-Spike binding. The conformational flexibility of the drug-binding pocket was captured using the RMSD-based clustering of the ligand-free simulations. An ensemble docking was performed wherein the two databases were docked on each of the cluster representatives of ACE2-Spike. The potential phytochemicals identified belonged to Withania somnifera, Swertia chirayita, Tinospora cordifolia, Andrographis paniculata, Piper longum, and Azadirachta indica. The FDA molecules identified were rutin DAB10, fulvestrant, cefoperazone acid, escin, chlorhexidine diacetate, echinacoside, capreomycin sulfate, and elbasvir.
Introduction<!>Methodology<!>Molecular Dynamics<!>Molecular Docking<!>Clustering and ensemble docking<!>Analysis<!>Results and Discussion<!>Rutin DAB10<!>Cefoperazone acid<!>Pinaverium bromide<!>Abitrexate<!>Direct Docking on crystal structure with phytochemical database<!>Rank<!>Octadecanoate<!>Guineesine<!>Oleic acid<!>3-O-caffeoyl-D-quinic acid<!>Conformational variation in the ACE2-Spike complex<!>Root Mean Square Deviation (RMSD)<!>Principal Component Analysis<!>3.3.1.2a ACE2 conformational variation<!>3.3.1.2b Spike conformational variation<!>Binding efficiency of Rutin DAB10 and Swertiapuniside<!>Rutin DAB10 and Swertiapuniside interactions with ACE2-Spike complex<!>ACE2-Spike active site flexibility and ensemble docking<!>Ensemble docking of FDA database<!>3.3.4.1b Chlorhexidine diacetate<!>3.3.4.1c Echinacoside<!>3.3.4.1d Capreomycin Sulfate<!>3.3.4.1e Elbasvir<!>Ensemble docking of phytochemical database<!>3.3.4.2a Amarogentin<!>3.3.4.2c Cardifolioside B<!>3.3.4.2d Swertiapuniside<!>Conclusion<!>Conflict of Interest Statement
<p>The Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2), also known as novel coronavirus was first identified in Wuhan province of China, owing to its high degree of pathogenicity, it spread all over the world rapidly [1][2][3]. The outbreak due to SARS-CoV-2 has claimed more than 0.47 million lives and hence infected more than 0.9 million people worldwide [4]. The virus is spreading at an alarming rate and has become a serious concern for human health.</p><p>The first sequence of the viral genome from an infected patient was released in January 2020 [5].</p><p>The phylogenetic analysis reported in the literature suggests that, this virus belonged to the SARS family and hence was later on named as SARS-CoV-2 [6,7]. The SARS-CoV-2 virus belongs to the family of positive single-stranded RNA viruses [8,9]. The SARS-CoV-2 belongs to the genus of Betacoronavirus, which also includes the other viruses viz. MERS-CoV, SARS-CoV, and bat SARS related coronavirus to name a few [7][8][9][10]. Phylogenetic analysis showed that SARS-CoV-2 is closely related to the bat SARS and with a sequence identity of 93.1% for the S-glycoprotein (spike) gene [11]. The genome size of the virus is approximately 29Kb. It codes for nearly 30 proteins responsible for various functions. These 30 proteins include structural proteins, nonstructural proteins, and accessory proteins. The virus encodes for four structural proteins viz. spike (S), membrane (M), envelope (E), and nucleocapsid (N) proteins. Among the structural proteins, spike protein forms an important part of the virus as it is known to form the first line of contact with the human cells and controls the downstream viral processes of attachment, fusion, and host cell entry.</p><p>The human Angiotensin-Converting Enzyme-2 receptor (ACE2) was found to be the functional receptor for SARS-CoV through in vivo and in vitro studies [11][12]. Zhou et. al. [1] through his studies showed that ACE2 is also a functional receptor for SARS-CoV-2. ACE2 plays a major role in lowering the blood pressure by catalyzing the maturation of angiotensin (1-7) (a vasodilator, a peptide hormone) from angiotensin II (a vasoconstrictor peptide) [13][14][15]. ACE2 is present on the membranes of cells of most organs viz. lung type II alveolar cells, heart, kidney, endothelium, and intestine. Thus, the presence of ACE2 on different cells explains the damage caused to the different organs due to the SARS-CoV-2 infection [1,[15][16][17]. The ACE2 protein comprises of three domains, the extracellular domain (residue range 18-470), the transmembrane domain (residue range 741-761), and the cytoplasmic domain (residue range 762-805) [18]. The spike protein binds to the extracellular enzymatic domain of ACE2 protein resulting in endocytosis and translocation of the virus [1,19].</p><p>The spike protein comprises two subunits S1 and S2 [20][21][22]. The spike protein exists as a heterotrimer to facilitate the binding to ACE2 protein. It has been reported that spike protein binds to ACE2 protein on the surface of the human cell which is followed by cleavage of spike protein at the junction of S1 and S2 subunit [21]. Based on the studies of SARS-CoV spike protein it was proposed that this cleavage is followed by the release of S1-ACE2 complex [22][23][24][25]. Following the release of the S1-ACE2 complex, S2 transit to a stable post-fusion state from the metastable pre-fusion state, an essential step for membrane fusion. Hence, the binding of spike protein to ACE2 forms a critical step for viral infection [22][23][24][25]. Also, studies on the HeLa cell lines, have shown that the cells which don't express ACE2 are not susceptible to viral inception [1]. In-vitro binding studies suggest, good binding affinity at low nanomolar range between ACE2 and the Receptor Binding Domain (RBD) (residues range 331 to 524 of spike protein) of the SARS-CoV-2 spike protein [28]. All these studies show that for binding to ACE2 receptor on human cell the RBD domain of S1 subunit of spike protein plays a key role. Within the RBD domain, the Receptor Binding Motif (RBM) (residue range 438 to 506) directly interacts with ACE2 [21]. The interacting residues on the surface of the RBD domain of the spike protein are Lys417, Gly446, Tyr449, Tyr453, Leu455, Phe456, Ala475, Phe486, Asn487, Tyr489, Gln493, Gly496, Gln498, Thr500, Asn501, Gly502 and Tyr505 [21,29]. The interacting residues belonging to the ACE2 receptor are Gln24, Thr27, Phe28, Asp30, Lys31, His34, Lys35, Glu37, Asp38, Tyr41, Gln42, Leu79, Met82, Tyr83, Asn330, Lys353, Gly354, Asp355, Arg357 and Arg393. These residues form various interactions between the ACE2 receptor and the RBD of spike, which has been shown in Figure 1. These interactions include 13 hydrogen bonds and 2 salt bridges between the ACE2 receptor and the RBD-spike interface residues (Table 1). Henceforth, the ACE2 receptor and RBD of spike have been referred to as ACE2 and spike respectively.   All these atomistic details of interactions confirm that the spike and ACE2 interaction is a must for virus entry in the cell. Studies have been reported earlier, wherein the ACE2-Spike complex of SARS-CoV-2 has been explored for druggable pockets [30,31]. One of the works by Patil et.al in June 2020 reveals the interface of ACE2-Spike as one of the potential druggable pockets of this complex [30]. There have been in-silico studies which suggest that the ACE2-Spike complex proves a potential drug target. One of such studies, performed by Smith et. al. in February 2020, explores the conformational changes of the complex using molecular dynamics simulations [31].</p><p>Followed by targeting the different ensembles of ACE2-Spike interface. This study enabled them to identify potential drugs that would block the ACE2 and spike interactions by binding at the ACE2-Spike interface [31]. Hence, the current study aims to find the inhibitors which can bind at the interface of ACE2 and RBD domain of spike through docking and simulation methods. The work reported in this article involves the docking of the entire Food and Drug Administrations (FDA) drug database. Besides this, the use of natural products from medicinal plants was also studied as these molecules also pose as a safe option for therapeutics. A database of phytochemicals belonging to the medicinal plants used to treat respiratory disorders was built. The docking of these phytochemicals from various medicinal plants viz. Ocimum sanctum, Withania somnifera, Piper longum, Tinospora cordifolia, Curcuma longa, Terminalia arjuna, Andrographis paniculata, Swertia chirata, Azadirachta indica, Aloe barbadensis to name a few were also carried out. The top-ranked molecules obtained through the direct docking of the ACE2-Spike complex (PDB ID: 6LZG) were rutin DAB10 (FDA) and swertiapuniside (Phytochemicals). However, a detailed analysis of the interactions between the top five ranked molecules and the ACE2-Spike complex (PDB ID: 6LZG) from both the databases was performed. The ligand-free ACE2-Spike, rutin DAB10-bound, and swertiapuniside-bound complexes were further studied through molecular dynamics (MD) simulations. The MD simulations of the ligand-free ACE2-Spike (APO) complex was used as a control to observe changes occurring in the ligand-bound systems.</p><p>The comparative analysis of the conformational parameters revealed the residues that are crucial in binding to these two molecules and they induce an inhibitory effect. All the three systems namely, APO, rutin DAB10-bound ACE2-Spike and swertiapuniside-bound ACE-Spike were simulated for a cumulative of 300 ns. Further, to understand the flexibility of participating residues in ACE2-Spike complex, the simulation data of APO complex was subjected to clustering to obtain instances of various ensembles visited by this protein complex. These cluster representative structures were further used for ensemble docking of the FDA and phytochemical database of natural compounds. Ensemble docking helped in identifying a few other potential molecules from FDA and phytochemical databases that may have a role in inhibiting the ACE2-Spike complex.</p><!><p>The crystal structure of ACE2-Spike with PDB ID: 6LZG was retrieved from the Protein data bank [32]. The PDB was pre-processed by removing the co-ordinates of all the molecules except for the co-ordinates that belonged to the spike protein and ACE2. This PDB was considered as APO form and was used for molecular docking and MD simulation studies. The methodology has been depicted in the Figure 2. The molecular docking was performed for the screening viz. FDA molecules database and phytochemical database of medical plants. The phytochemical database of medicinal plants was created by downloading the phytochemicals from PubChem database. The phytochemicals from medicinal plants which were known to have a role in treating respiratory tract infections were downloaded. This selection resulted in a database consisting of 250 phytochemicals. The top-ranked drug obtained from the docking of each of the two databases on the APO form (PDB ID: 6LZG) of ACE2-Spike complex were further subjection to MD simulations. As a control, the MD simulations of APO form of the ACE2-Spike complex was also carried out. Hence, three ACE2-Spike systems, one APO and two ligand-bound were simulated for a cumulative of 100 ns each. A total of 300 ns simulation data has been presented here and analyzed. In addition to docking of the experimental structure of ACE2-Spike complex (direct docking) and MD simulations of APO and ligand-bound complexes, an ensemble of different conformations of ACE2-Spike complex was generated by clustering of the APO's MD simulation data. The docking of the two mentioned databases was performed considering these different conformations of ACE2-Spike complex obtained from clustering of APO simulation data (ensemble docking). The details about each step has been given below.</p><!><p>MD simulations were carried out using the AMBER 16 simulation package [33]. AMBER14ffSB</p><p>force field [34] was used to parameterize the protein molecules. Na+ ions were added to neutralize the system. TIP3P water model was used to represent water molecules and the box size used was of 14 Å. The solvated system was minimized using steepest descent method for 20000 steps followed by the conjugate gradient method. Following minimization, the system employed the Langevin thermostat for gradually heating it to 300 K. In order to deal with the hydrogen restraints, the SHAKE algorithm was employed. After heating system up to 300 K, it was equilibrated at NPT conditions for 1 ns with pressure and temperature being 1 atm and 300K respectively. The production run was carried out for a cumulative time of 100 ns. The parameters for the ligand in the ligand-bound simulations were derived from the antechamber module of Ambertools 17 and the force field used was General Amber Force Field (gaff) [35,36].</p><!><p>Molecular docking was carried out using DOCK6 [37]. The pre-processing of the ACE2-Spike complex involving addition of hydrogens and charges was carried out using UCSF Chimera [38].</p><p>The standard parameters were used for docking. The binding pocket identification was done using the sphgen, sphere_generation and sphere_selection module of DOCK6. This module identifies largest available cavity in the protein as binding cavity. The top-ranked molecules were chosen based on grid score. The grid scores signify the strength of binding for any small molecule to the receptor protein. It is an energy-based function which sums up all the non-bonded interactions formed by the small molecule with the active site of the receptor molecule. Hence, a more negative value of a grid score indicates better binding of ligand molecule to the receptor protein [37]. In case of direct docking, the receptor molecule was the APO form of ACE2-Spike complex obtained from the PDB ID: 6LZG. The top-ranked molecule obtained by screening of the FDA database was rutin DAB10. The top-ranked molecule obtained by screening of the phytochemical database was swertiapuniside. Rutin DAB10-bound and swertiapuniside-bound ACE2-Spike complexes were further subjected to MD simulations for a cumulative of 100 ns each. The parameters and protocol used for these two simulations were as described in section 2.1.</p><!><p>The MD simulation data of APO form of ACE2-Spike was subjected to clustering to obtain clusters representing different conformations visited by the protein complex throughout the simulations.</p><p>Root Mean Square Deviation (RMSD)-based clustering was carried out using the dbscan method of cpptraj module of AmberTools 17 [39]. The RMSD cut-off used was of 2 Å. This RMSD-based clustering of the simulation data led to the formation of five clusters. Representative structure of ACE2-Spike from each of these five clusters were considered for docking. This approach has been referred as ensemble docking. Five independent docking were performed for screening the FDA and phytochemical databases separately. Hence, a total of 10 independent docking were performed on the five cluster representatives of ACE2-Spike complex. The docking protocol used in each of these cases was similar to that explained in section 2.2.</p><!><p>The cpptraj module of AMBERTOOLS17 was used for the analysis of the MD simulation data.</p><p>The interaction analysis between the docked molecules and the protein were done using PLIP [40],</p><p>UCSF-Chimera and Ligplot [41].</p><!><p>3.1 Direct Docking on crystal structure using FDA Database</p><p>The direct docking of the ACE2-Spike complex for screening the FDA database was performed using the protocol defined in section 2.2 (Figure 2). The top 5 hits from FDA approved database along with their grid scores have been given in Table 2. The top five ranked molecules were rutin DAB10, fulvestrant, cefoperazone acid, pinaverium bromide, and abitrexate. These molecules are known to have diverse roles in terms of their therapeutic properties. The known indication for these molecules has also been listed in Table 2.</p><p>Table 2: Top five ranked molecules from the FDA approved drug database with their therapeutic properties based on earlier use and grid scores. Interaction analysis for these five molecules with the residues of the spike and ACE2 protein was carried out. Figure 3 pictorially depicts the hydrogen bonding interactions formed by each of these five molecules with the residues of both the proteins of the ACE2-Spike complex. The supplementary table S1 also reports all the non-bonded interactions namely, hydrophobic (HP), hydrogen bonding (HB), salt bridge (SB) and π-cation (PC), between the atoms of the ligand and the ACE2-Spike complex. Each of these five molecules and their interactions with the ACE2-Spike complex have been discussed individually below.</p><!><p>Rutin DAB10 is a flavonoid and known to have a variety of biological activities. It has antiallergic, anti-inflammatory, antiproliferative, and anticarcinogenic properties. It also acts as a deliverer of quercetin to the small intestine. The docked pose of rutin DAB10 with ACE2-Spike complex showed hydrophobic interactions with Asn33 and Glu37 of ACE2 protein. It showed hydrogen bonds with His34, Glu37, Lys353, Phe390 of ACE2, and Arg403, Glu406, Gln409 of spike. It also showed salt bridge interaction with His34 of ACE2 and Lys417 of spike. Besides these, it showed π-cation interaction with His34 of ACE2 (Table S1).  S1).</p><!><p>It is an antibiotic agent and used for the treatment of different bacterial infections caused by susceptible organisms in the body, including respiratory tract infections, peritonitis, skin infections, endometritis, and bacterial septicemia. Hydrophobic interactions with Pro389 of ACE2 and Tyr495 of the spike protein were observed. It showed the formation of hydrogen bonds with Asn33, Glu37, Ala386 of ACE2, and Gly504, Tyr505 of the spike protein. It also formed π-cation interaction with His34 of the ACE2 and Arg403 of the spike protein.</p><!><p>It is used for gastrointestinal disorders. It is an anti-spasmodic and helps to restore the normal contraction process of the bowel. It was observed to form hydrophobic interaction with Glu 37, Gln96 of the ACE2, and Tyr495 of the spike. Along with, hydrogen bonds with the Asp30, His34, and Asp38 of the ACE2 protein. It also formed π-cation interaction with Arg403 of the spike (Table S1).</p><!><p>It is used in the treatment of certain types of cancer. It is also used in the control of severe psoriasis or rheumatoid arthritis. It also has implications in the control of juvenile rheumatoid arthritis. It formed hydrogen bonding interactions with Asn33, His34, Glu37, Gln96, Lys335 of the ACE2 and Arg403, Ser494, and Gly496 with the spike protein. It also formed a salt bridge with His34 of the ACE2.</p><p>Interaction analysis suggests that, all the top five ranked ligands formed good interactions either with the residues ACE2 or the spike protein. Also, all these ligands tend to form interactions with the residues which belong to the interface region of the ACE2-Spike complex. This may infer that they may be able to interfere with the ACE2 and spike interactions and thereby block the formation of a stable complex.</p><!><p>The direct docking of the ACE2-Spike complex was performed to screen the molecules belonging to the phytochemical database. The phytochemical database consisted of compounds with their origin lying in the medicinal plants. The grid scores obtained for all these phytochemicals were sorted and the top five ranked molecules were selected for analysis of their interactions with the ACE2-Spike complex. The top five ranked phytochemicals were swertiapuniside (Plant source:</p><p>Swertia chirayita), octadecanoate (Plant source: Azadirachta indica), guineensine (Plant source:</p><p>Piper longum), oleic acid (Plant source: Azadirachta indica) and 3-O-caffeoyl-D-quinic acid (Plant source: Andrographis paniculata). The grid score and the plant source with their popular names of these docked compounds have been listed in Table 3.</p><p>Table 3: Top five ranked molecules from the phytochemical database with the details of their plant source and grid scores. These plants listed in Table 3 are also known to have anti-bacterial properties and help in increasing the immunity. In order, to find out their binding mode and interactions with ACE2-Spike complex, a detailed analysis of these top-ranked five molecules was done. Figure 4 pictorially depicts the hydrogen bonding interactions formed between these five molecules with the residues of ACE2-Spike complex.</p><!><p>Supplementary Table S2 shows all the other non-bonded interactions like, hydrophobic, hydrogen bonding, salt-bridges and π-cation between the molecules and the residues of the ACE2-Spike complex. The detailed description of the interactions made by these five molecules with the ACE2-Spike has been discussed individually.</p><!><p>The second-ranked molecule octadecanoate has its source in the plant Azadirachta indica popularly known as Neem or India lilac. This phytochemical is known to have multiple therapeutic</p><!><p>The third-ranked molecule, guineesine has its source in the plant named Piper longum, it is an alkaloid. It has anti-inflammatory and analgesic properties. It was observed to form hydrophobic interactions with Glu37 of ACE2 and Glu406, Gln409, Tyr769 of the spike protein. It formed hydrogen bonds with His34, Glu37 of ACE2, and Gln409 of the spike. It forms π-cation interaction with His34 of ACE2 and Arg403 of spike protein.</p><!><p>The fourth-ranked phytochemical was oleic acid from the plant Azadirachta indica. It is known to have a beneficial effect in the case of cancer. It is also found to be useful in autoimmune and inflammatory diseases and has the ability to facilitate wound healing. It formed hydrophobic interactions with Asn33, Thr92, and Val93 of ACE2. It formed hydrogen bonds with Lys353, Ser494 of ACE2 protein, and Gly496 of spike protein. It also formed π-cation interaction with Lys353.</p><!><p>The phytochemical database, referred as Swertia). Total three systems were simulated, each system was simulated for 100 ns. Thus, a total 300 ns of simulation data was analyzed.</p><!><p>The conformational variation in the ACE2-Spike complex systems was measured based on three parameters namely, root mean square deviation (RMSD) and principal component analysis (PCA).</p><!><p>The RMSD tells about how much a structure has deviated during the simulations from its initial structure. Figure 5 A shows a histogram plot for RMSD for APO, Rutin, and Swertia systems.</p><p>From the figure, it is clear that the APO system showed a minimum deviation from the starting structure with a maximum population having an RMSD around 2.1 Å. The Rutin system showed an RMSD deviation in the range of 2-3 Å. The Swertia system showed an RMSD deviation in the range of 2.5-4.5 Å. Further, in order to see which of the protein, ACE2, or the spike was responsible for this deviation, histogram plots for RMSD for these individual proteins were plotted. Figure 5 B shows a histogram plot for RMSD of ACE2 protein for APO, Rutin, and Swertia system. The ACE2 protein in the APO system showed RMSD in the range of the 1-2.5 Å. The ACE2 protein in the Rutin system showed RMSD in the range of 1.7-2.5 Å. The ACE2 protein in Swertia system showed RMSD in the range of 2.5-4 Å. Figure 5 C shows the histogram plot for RMSD of Spike protein for APO, Rutin, and Swertia system. The spike protein in the APO system showed RMSD in the range of the 1-2.3 Å. The spike protein in the Rutin system showed RMSD in the range of 1.6-3.1 Å. The spike protein in Swertia system showed RMSD in the range of 1.6-2.6 Å. The RMSD analysis shows that for the entire complex (ACE2-Spike) the maximum deviation was seen in the Swertia system and minimum in the APO system. However, for the ACE2 protein, the maximum deviation was seen in the Swertia system, and in the case of the spike protein, the maximum deviation was seen in the Rutin system. Thus, both the ligands tend to induce significant conformation changes in the ACE2-Spike complex.</p><!><p>In order to identify statistically significant conformational changes, the MD simulation data were subjected to principal component analysis (PCA). The reaction co-ordinates considered for the three systems were all the atoms except hydrogen. The aim behind performing the principal component analysis was to identify the most dominant conformation that is attained by the ACE2-Spike complex in APO and ligand-bound simulations. This would further help in identifying the crucial residues undergoing conformational changes from the ACE-Spike complex. Figure S1 shows the normalized population distribution along the principal component (PC) 1 (A), PC 2 (B), and PC 3 (C). The principal component 1 showed two significant populations in the case of all the three simulations systems (Figure S1 A). The population observed at eigenvalue 50 showed an overlap in the case of all the three simulation systems. This suggests that the conformation captured in this population was observed in all the three simulation systems. However, the ligand-bound systems sampled one more population along the eigenvalue -50±10. The occurrence of this population in the case of APO was comparatively lower. Hence, it may infer that the conformations captured in this population were because of the presence of the ligand in these simulations. In the case of the second and third principal components, the maximum population overlapped in all three systems. Even if two more distinct conformations were observed for Rutin simulations (Figure S1 B), the population was comparatively lower than the most populated conformation.</p><!><p>The residue-wise RMSF captured by the first three principal components for ACE2 and spike has been presented in Figures 6 and 7 respectively. Figure 6 A represents the fluctuations captured by the principal component 1. In case all the three systems the region around 112 to 150 residue number showed significant variation. This variation was maximum in the case of the APO as compared to the ligand-bound systems. Figures 6 D, E, and F show the variance captured by PC1 in the residues of the ACE2-Spike complex for the APO, Rutin and Swertia systems. The color code follows the red-green-blue (RGB) gradation corresponding to increase in variance/fluctuations. Figure 6 D clearly depicts the highly varying residues within the range 112-150 of the ACE2, which can be seen in green color. This helical region lies in the core of the ACE2 receptor. Hence, variation in this region suggests conformational dynamics. It was observed that APO and Swertia showed similar motions for the residue range112 to 150. However, in the case of the Rutin system, it was in the opposite direction. This suggests that the changes induced due to Rutin lead to a conformational change that is not native to the changes observed in the APO.</p><p>The fluctuations over the entire ACE2 receptor appeared to reduce in the second principal component (Figure 6 B). The fluctuations captured by the third component suggest that these occurred but were not present for a longer time in the simulation (Figure 6 C).</p><!><p>The principal component 1 captured fluctuations in different regions of the spike protein. Figure 7A clearly shows the Rutin system fluctuated most within the residue range 375 to 400 and 512 to 525. The residues within the range 375 to 400 fluctuated in the case of the APO and Swertia systems too, however, their magnitude was less as compared to Rutin. The residue stretches from 350 to 400 appeared to be highly fluctuating as the variation was being captured with significant magnitude till the third principal component in all the three systems. The variation in the region Although some variation in this region was captured by PC2 for the other two systems. This PC1 captured variation projected on to the structure of spike protein can be seen in Figure 7 D, E and F. The color gradation from red to green to blue (RGB) corresponds to increase in variance/fluctuations of the residues along the spike protein.</p><!><p>The binding efficiency of both the ligands to the ACE2-Spike complex was calculated using the MMGBSA free energy of binding. The three simulations system, APO, Rutin and Swertia were compared. The free energy contribution made by the interface residues of ACE2 receptor and spike protein were calculated in all three systems. Figure 8 A shows comparison of average free energy contribution in binding between ACE2 and spike by the ACE2 residues for the APO, Rutin and Swertia systems. Lys31 and His34 showed a decrease in free energy contribution in ACE2-Spike binding in the presence of the ligands and Rutin and Swertia. Remaining all the ACE2 interface residues showed the free energy contribution either similar or better than APO as compared the ligand bound systems. Figure 8 B shows the comparison of average free energy contribution in binding between ACE2 and spike by the spike residues for the APO, Rutin and Swertia systems.</p><p>Gln493 showed a decrease in free energy contribution in ACE2-Spike binding in the presence of either of the ligands. Rest of all the spike protein residues were observed to have either similar or better free energy contribution for APO as compared to both the ligand systems. Hence, in order to understand the role of ligand in altering the binding efficiencies of the ACE2 and spike residues, they were compared for contributing in binding between the two protein and the ligand molecule.  for the Tyr453 of spike protein (Figure 9 D), where in the free energy contribution in binding to the rutin DAB10 was better than to the ACE2. In case of Tyr453 too, the free energy difference was observed to be around 2 kcal/mole.</p><p>Figure S3 represents the average free energy contribution in binding obtained for the residues of ACE2 and spike in the Swertia simulation system. Figure S3 A represents the residues of ACE2</p><p>and their average free energy contribution in binding to spike (black) and swertiapuniside (green).</p><p>The His34 and Arg393 of ACE2 were observed to contribute more in binding to the ligand, swertiapuniside as compared to the spike. Tyr453 of spike too, showed more contribution in binding to swertiapuniside as compared to the ACE2 (Figure S3 B). The free energy contribution made by these three residues throughout the 100 ns simulation has been shown in Figure 10. Figure 10 A shows the free energy contribution for His34 of ACE2. The free energy contribution was observed to be fluctuating as multiple peaks were obtained when binding with swertiapuniside.</p><p>However, significant conformations showed better binding with the ligand, swertiapuniside with a maximum free energy difference of around 4 kcal/mole. In case of, Arg393 of ACE2 the free energy difference was observed to be 1 kcal/mole (Figure 10 B). This residue too was observed to have better binding to swertiapuniside. Tyr453 of spike too showed better free energy contribution in binding to swertiapuniside as compared to ACE2 (Figure 10 C).</p><!><p>In order to identify the interaction between two proteins viz. ACE2 and spike, interaction analysis was carried out. The supplementary figure S4 shows the number of hydrogen bonds between ACE2 and spike protein for the three systems. The average number of hydrogen bonds for the APO system was around 10 to 11 throughout the system. For the Rutin system, the average number of hydrogen bonds were around 10 to 11 for throughout the system. For the Swertia system, the average number of hydrogen bonds was around 8 to 9 for throughout the system. At the start of the simulation, the number of hydrogen bonds was around 12. Hence, this number decreased in the case of Swertia system. This shows that the presence of ligand destabilizes the number of hydrogen bonds between ACE2 and spike protein. Further to check which are the hydrogen bonds were present throughout the simulations contact analysis was done.</p><p>As mentioned in the introduction the crystal structure showed the presence of 13 hydrogen bonds and two salt bridges, these bonds were first analyzed during the course of simulations. The histogram plot shown in Figure 11 A, depicts their occupancy during the course of simulation in Figure 11: Percentage occupancy of native(A) and newly identified (B) hydrogen bonds between the ACE2 and spike proteins for the APO (black), Rutin (red) and Swertia (green) systems the three systems. The hydrogen bonds were Gln24-Asn487, Tyr83-Asn487 and Lys353-Gly502 had 100% occupancy and Tyr41-Thr500 had 80% occupancy for all the three systems. Besides these three hydrogen bonds, His34-Tyr453, Asp38-Tyr449, Tyr41-Asn501, Gln42-Gly466, Gln42-Tyr449, Gln42-Tyr489 and Arg393-Tyr505 has higher occupancy in Apo simulations as compared to Rutin and Swertia simulations. The remaining two hydrogen bond Glu35-Gln493 and Glu37-Tyr505 were presented at higher occupancy in the APO and Rutin system as compared to the Swertia system. The salt bridge between Asp30-Lys417 had higher occupancy in Rutin followed by APO and then Swertia system.</p><p>Besides these interactions, few new interactions were formed during the course of simulations and are reported in Figure 11 B. These interactions are hydrogen bonds, hydrophobic interactions, and charged interactions. The histogram plot shown in Figure 11 B shows hydrogen bond occupancy during the course of simulation in the three systems. The hydrogen bond Lys353-Asn501, Lys353-Gly496, Gly354-Gly502, Asp355-Thr500, Arg357-Thr500 show 100% occupancy in all the three system. The hydrogen bond Lys353-Gln498 shows 100% occupancy in APO and Rutin systems, while in the Swertia system it is almost absent. The Phe28-Thr489 and Asn330-Thr500 shows 100% occupancy in APO and Swertia system while Phe28-Thyr489 shows 50% occupancy and Asn330-Thr500 is almost absent in the Rutin system. The His34-Gln493 shows 70% occupancy in the APO system while it is around 2 % in the Rutin and Swertia systems. Besides hydrogen bonds, hydrophobic interactions were also seen. The residues Thr27-Phe456, Tyr41-Gln498, and Lys353-Tyr505 showed 100% occupancy for all the three systems. The Met82-Phe486 and Tyr83-Phe486 showed 100% occupancy for APO and Rutin systems, while Met82-Phe486 showed 40% occupancy and Tyr83-Phe486 showed 5% occupancy in Swertia system. The Thr27-Tyr489 and Gly354-Tyr505 showed 100% occupancy in the APO system, while Thr27-Tyr489 showed a 2 % occupancy in the Rutin system, 30% occupancy in Swertia system and Gly354-Tyr505 showed 2 % occupancy in Rutin system, 50% occupancy in Swertia system.</p><p>Figure 12 A shows the number of hydrogen bonds between the ACE2 protein and the ligands viz.</p><p>rutin DAB10 (Red) and swertiapuniside (Green). Figure 12 B shows the number of hydrogen bonds between the spike protein and the ligands viz. rutin DAB10 (Red) and swertiapuniside (Green). From the figures, it is clear that rutin DAB10 tends to more hydrogen bonds with ACE2 protein as compared to swertiapuniside. The average number of hydrogen bonds between ACE2 and rutin DAB10 is around 5 and that with swertiapuniside is around 3. However, swertiapuniside tends to form more hydrogen bonds with the spike protein as compared to rutin DAB10. The average number of hydrogen bonds between spike and swertiapuniside is around 4 and that with rutin DAB10 is around 3. Figure 13 A shows the percentage occupancy for hydrogen bonds between ACE2 and Ligands. The residues Lys26, Asp30, Glu37, Ala387, Gln388, and Arg393 of ACE2 residues are involved in the hydrogen bonding with rutin DAB10 and swertiapuniside. The percentage occupancy varied for each of the residues for rutin DAB10 and swertiapuniside. Figure 13 B shows the percentage occupancy for hydrogen bonds between spike and ligands. The residues Arg403, Glu406, Arg408, Gln409, Thr415, Asp417, Thr453, and Tyr505 of spike residues are involved in the hydrogen bonding with rutin DAB10 and swertiapuniside. Besides hydrogen bonds, ligands are also involved in forming hydrophobic interaction with ACE2 residues viz.</p><p>Leu29, Asn33, Pro389 and spike residues Tyr453, Tyr495 and Tyr505.</p><!><p>In order to investigate the flexibility of participating residues in the binding cavity of the ACE2-Spike complex and its effect on the possibility of other drugs targeting the ACE2-Spike complex have been carried out. The total cumulative time of 100 ns of ACE2-Spike (APO) complex was used for RMSD based clustering. Clustering was done to obtain structure for ensemble docking where the flexibility of the binding cavity has been considered. Five different clusters were obtained with RMSD cut off of 2 Å. Molecular docking was done on the representative structures of each cluster (Figure 2). The docking protocol as explained in the section 2.2. was followed.</p><p>Both the FDA database and phytochemical database were docked. Based on grid score top tanked molecule for each of representative structure were sorted. Similar approach has been reported for the 3C-like protease simulations where, potential drug candidates were identified using the ensemble docking approach [42]. Ensemble docking approach with molecular dynamics simulations of remdesivir bound to RNA dependent RNA polymerase has also been reported in one of the earlier works [43].</p><!><p>Table 4 gives details about the top ranked molecule from the FDA database obtained for each of the five representative clusters of ACE2-Spike complex.  Detailed interaction analysis was done for each of these five docked complexes. The Figure 14 shows the interaction plot for each of the docked complex. Table S3 lists the residues that were involved in different types of non-bonded interactions with top ranked molecules of the FDA approved database for simulation data. Arg408 of the spike protein.</p><!><p>The chlorhexidine diacetate act as an anti-infective agent and antibacterial agent. It is a bestdocked molecule for cluster 2. It showed hydrophobic interaction with Ala71, Phe72, and Glu75 of the ACE2. It showed hydrogen bonds with Glu35, Asp38, Lys68, and Gln42 of the ACE2. The Chlorhexidine diacetate does not show any interaction with the spike protein.</p><!><p>The echinacoside acts as a neuroprotective and has beneficial cardiovascular effects. It is the bestdocked molecule for cluster 3. It showed hydrophobic interaction Lys417 and Tyr453 of spike protein. It showed hydrogen bonds with Gln96, Gln388, Arg393 of the ACE2 and Arg403, Asp405, Tyr505 of the spike protein.</p><!><p>The capreomycin sulfate is an aminoglycoside and acts as antibiotics.</p><!><p>The elbasvir is a direct-acting antiviral medication to treat chronic hepatitis C. It was the bestdocked molecule for cluster 5. It showed hydrophobic interactions with Leu29, Asp30, Met383, Ala387, and Pro389 of the ACE2. It showed formation of hydrogen bonds with Lys26 of ACE2 and Gly504 of the spike protein. It showed π-cation interaction with Arg408 of the spike protein.</p><!><p>The best docked molecules for phytochemical database of medicinal plants for representative structure for each of the cluster obtained from simulation are given Table 5. The best docked phytochemicals were chosen based on grid score. In the case of cluster 1 the best docked molecule was amarogentin whose plant source is Swertia chirayita. Cluster 2 screened sitoindoside IX as the best docked molecule, which is known to be obtained from Withania somnifera.</p><p>Cardiofolioside B was obtained as the best molecule for cluster 3. Its plant source is known to be Tinospora cordifolia. Clusters 4 and 5 screened swertiapuniside as the best molecule, the plant source being Swertia chirayita.</p><p>Table 5: Top ranked molecules from the phytochemical database for each of the five ACE2-Spike cluster representatives with the details of their plant source and grid scores. Detailed interaction analysis was done for each of the docked complex. Figure 15 shows interaction plot for each of the docked complexes.</p><!><p>The amarogentin phytochemical has the plant source Swertia chirayita, shows antibacterial and</p><!><p>The cardifolioside B having plant source Tinospora cordifolia, shows activity against typhoid, malaria. It is also found effective against filariasis, leprosy and has anthelmintic properties. It is also found useful in gout and rheumatoid arthritis. It is best docked molecule for cluster 3. It showed hydrophobic interaction with Arg408 of the spike protein. It showed formation of hydrogen bonds with Met383 of the ACE2 and Arg408, Gln409, Gln414, Gly502, Gly504 of the spike protein.</p><!><p>The swertiapuniside best docked molecule for cluster 4 and cluster 5. For cluster 4 it showed hydrogen bonds with Lys26, Asp30, His34, Asn90, Gln96 of the ACE2 and Glu406, Gln409, Lys417 of the spike protein. It showed salt bridge interaction with ARG417 of the spike protein.</p><!><p>The molecular docking and simulation studies of the ACE2-Spike complex in the ligand-freeform and in the ligand-bound form revealed several crucial interactions that may help blocking the binding between the human ACE2 receptor and SARS-CoV-2 spike glycoprotein's, receptor binding domain. The conformational parameters used to calculate the stability of the ACE2-Spike complex revealed increased variation in certain regions of the ACE2 and spike protein in case of ligand-bound systems. The interaction analysis of the top ranked drugs rutin DAB10 and swertiapuniside with the ACE2-Spike complex revealed that the residues present on the interface of the two proteins formed weaker interactions. Few hydrogen bonds which are known to present in the experimental structure of ACE2-Spike complex were observed to be maintained throughout the simulations when in the ligand free state. However, presence of these two ligands abrogates the formation of these hydrogen bonds, thereby destabilizing the ACE2-Spike complex. The stability of native hydrogen bonds was affected more when the ACE2-Spike was bound to swertiapuniside as compared to rutin DAB10. The free energy analysis also revealed the residues involved in hydrogen bonding namely, His34 from ACE2 and Tyr453 from spike were involved in forming interactions with the ligand molecules rather than one another. The ensemble docking also revealed FDA drugs viz. escin and capreomycin sulfate which are known to possess antiinflammatory and anti-bacterial activity. Elbasvir, an antiviral drug was also obtained as a top hit.</p><p>In case, of phytochemicals swertiapuniside was obtained as the top hit drug in direct as well as ensemble docking. Few other phytochemicals from the medicinal plants Neem and Ashwagandha were also observed to bind efficiently with the ACE2-Spike complex. These natural compounds were also observed to interact well with the residues present on the interface of the ACE2-Spike complexes. Owing to strong interactions with the drug targets few of the FDA approved and natural compounds of medicinal plant origin may prove to probable candidates for drug repurposing against the COVID-19 infection.</p><!><p>The authors declare no conflict of interest.</p>