key: cord-0839696-h62ii7ir authors: Kalhor, Hourieh; Sadeghi, Solmaz; Abolhasani, Hoda; Kalhor, Reyhaneh; Rahimi, Hamzeh title: Repurposing of the approved small molecule drugs in order to inhibit SARS-CoV-2 S protein and human ACE2 interaction through virtual screening approaches date: 2020-09-24 journal: Journal of biomolecular structure & dynamics DOI: 10.1080/07391102.2020.1824816 sha: a3b819849a46798b82fde784fa15e0c817673c3d doc_id: 839696 cord_uid: h62ii7ir Most recently, the new coronavirus (SARS-CoV-2) has been recognized as a pandemic by the World Health Organization (WHO) while this virus shares substantial similarity with SARS-CoV. So far, no definitive vaccine or drug has been developed to cure Covid-19 disease, since many important aspects about Covid-19 such as pathogenesis and proliferation pathways are still unclear. It was proven that human ACE2 is the main receptor for the entry of Covid-19 into lower respiratory tract epithelial cells through interaction with SARS-CoV-2 S protein. Based on this observation, it is expected that the virus infection can be inhibited if protein-protein interaction is prevented. In this study, using structure-based virtual screening of FDA databases, several lead drugs were discovered based on the ACE2-binding pocket of SARS-CoV-2 S protein. Then, binding affinity, binding modes, critical interactions, and pharmaceutical properties of the lead drugs were evaluated. Among the previously approved drugs, Diammonium Glycyrrhizinate, Digitoxin, Ivermectin, Rapamycin, Rifaximin, and Amphotericin B represented the most desirable features, and can be possible candidates for Covid-19 therapies. Furthermore, molecular dynamics (MD) simulation was accomplished for three S protein/drug complexes with the highest binding affinity and best conformation and binding free energies were also computed with the Molecular Mechanics/Poisson–Boltzmann Surface Area (MM/PBSA) method. Results demonstrated the stable binding of these compounds to the S protein; however, in order to confirm the curative effect of these drugs, clinical trials must be done. Coronaviruses are a genus of the Coronaviridae family; belong to the Coronavirinae subfamily, and the order of Nidovirales. They are categorized into four genera including Alphacoronavirus, Betacoronavirus, Gammacoronavirus, and Deltacoronavirus (Shanmugaraj et al., 2020; Siddell et al., 1983) . They are enveloped viruses with a large plus-strand RNA genome which are typically present among several species of animals such as cows, bats, camels, cats, and avian. They may transmit from animals to humans, a process termed "spill over" (Mukhtar & Mukhtar, 2020; Shanmugaraj et al., 2020) . More recently, a new Betacoronavirus has been found out provisionally named 2019 novel coronavirus (2019-nCoV) (Elfiky, 2020b; Zhu et al., 2020) . This virus is now officially known as severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) which causes the COVID-19 disease. This virus is probably originated from an animal repository and has recently triggered the epidemic in humans because of rapid transmission from human to human as well as its high mortality rate (Elfiky, 2020a; Mukhtar & Mukhtar, 2020; Wrapp et al., 2020) . Therefore, inhibiting the SARS-CoV-2 virus has been a serious challenge for researchers and clinicians and they have become motivated to introduce and develop vaccines and therapeutic antibodies as well as drugs against this virus. Hence, the first genome sequencing of SARS-CoV-2 was published by Fan Wu et al. from china. They performed Phylogenetic analysis of the whole-genome sequence, containing 29,903 nucleotides, and reported that the virus has 89.1% nucleotide similarity to a group of SARSlike coronaviruses. Comparison of their conserved domains revealed that the receptor-binding domain (RBD) of the spike protein of SARS-CoV-2 was closely related to those of .9% amino acid identity), which makes it capable to use the human ACE2 (Angiotensin-Converting Enzyme 2) receptor for cell entry . Additionally, recent studies have shown that the ACE2 is also the receptor for SARS-CoV-2's entry into lower respiratory tract epithelial cells, (Agostini et al., 2018; Chang et al., 2020; Morgenstern et al., 2005; Shang et al., 2020; Wrapp et al., 2020; Xu et al., 2020) . Designing novel drugs against a new virus through experimental techniques is very time-consuming; however, it is required to find an effective drug immediately to treat the infection and decrease death cases. Therefore, it seems to be logical to search for potential therapeutics among previously approved drugs. Based on the above statement, in this study, potential agents were identified to inhibit the interaction of RBD domain of the SARS-CoV-2 with ACE2 receptor by using virtual screening approaches. Our results showed that among the studied drugs, Diammonium Glycyrrhizinate, Digitoxin, Ivermectin, Rapamycin, Rifaximin, and Amphotericin B might be effective therapeutics for the treatment of Covid-19 infection due to their better binding affinities and conformations. Lastly, MD simulation analysis and binding free energy calculations were accomplished for SARS-CoV-2-RBD/Diammonium Glycyrrhizinate, SARS-CoV-2-RBD/Digitoxin, and SARS-CoV-2-RBD/Ivermectin complexes, which had the highest binding affinity and the best conformations. Results of this study indicated that in silico approaches can be effectively used to develop a drug discovery pipeline using FDA approved drug databases, and it may lead to introduce novel potentials for the old drugs. Virtual screening approaches are extensively being applied in designing and development of new drugs. In this regard, one of the most common virtual screening techniques is structure-based virtual screening (SBVS) which only needs the three-dimensional structure of the interested protein and identifying its potential binding pockets to choose drugs, which interact strongly with these binding pockets, from large databases Kalhor, Sadeghi, et al., 2020; Shiri et al., 2018 Shiri et al., , 2019 . In the SBVS method, identification and preparation of the target receptor is an essential step. Hence, the crystallographic structure of SARS-CoV-2 S protein (RBD domain) in complex with ACE2 (PDB entry: 6VW1) was applied for molecular docking studies (Shang et al., 2020) . Also, other complexes of the SARS-CoV/ACE2 (PDB IDs: 6acj, 6acg, 6ack, 2ajf) ( F. Li et al., 2005; Song et al., 2018 ) were used for structural alignment analysis. Subsequently, the selected SARS-CoV-2-RBD was prepared as known receptor using AutoDock Tools 4.2. First, water molecules were removed, and then atoms were adjusted to the AutoDock atom types. In the next step, bond orders were assigned, and hydrogen atoms and Gasteiger-Marsili charges were added to the crystal structure. In order to conduct molecular docking, the receptor structure was prepared in the PDBQT format (Morris et al., 2009 ). In order to perform virtual screening, the binding pocket of SARS-CoV-2-RBD was selected based on the binding pocket inferred from the crystallographic structure of SARS-CoV-2/ ACE2 complex. To determine the main residues involved in the binding regions, the crystal structure was studied using PDBsum web tool and LigPlot þ software (Laskowski, 2009; Laskowski & Swindells, 2011) . Accordingly, Tyr449, Tyr453, Lus455, Phe456, Ala475, Gly476, Phe486, Asn487, Tyr489, GLn493, Gly496, Gln498, Tyr500, Asn501, Gly502 and Tyr505 residues were chosen as the reference amino acids for virtual screening. By using AutoDock Tools4.2, Grid box for docking studies was predicted into X ¼ 26 Å, Y ¼ 42 Å, Z ¼ 26 Å grid points, and the grid spacing was 1 Å. In this work, two different FDA approved drug databases were downloaded from the zinc database15 (1681 compounds) and used for virtual screening, including FDA Approved drugs in per DrugBank(https://zinc.docking.org/ substances/subsets/fda/) and FDA-approved-Drug-Library (1364 compounds). All compounds in this FDA-approved-Drug-library have well-characterized biological activity, clear targets, safety, and bioavailability properties which could significantly make drug development and optimization faster and easier. Also, since it covers various research areas, it is an comprehensive and ideal library for drug repurposing (https://www.targetmol.com/compound-library/FDA-approve d-Drug-Library) In this study, by using AutoDock version 4.2, each drug of the databases was prepared as follows: non-polar hydrogen bonds were merged, Gasteiger-Marsili charges were added, atoms were adjusted to the AutoDock atom types, and the rotatable bonds were assigned, then were saved in the PDBQT format. Finally, drugs were converted to the SDF using Open Babel software (O'Boyle et al., 2011) to be prepared for docking with the receptor. In this step, approximately all of the FDA approved drugs from two databases were docked in the selected binding pocket using AutoDock Vina (Koes et al., 2013) . (http:// smina.sf.net). Subsequently, based on the highest binding affinity, the result of docking was sorted and then were clustered based on the structural similarities and physicochemical features by ChemMine Web Tools (http://chemmine.ucr. edu/). Afterward, the selected drugs' physiochemical properties were evaluated using OSIRIS Data Warrior (Sander et al., 2015) . Ultimately, the potential drugs were chosen against the binding pocket of SARS-CoV-2-RBD. In order to evaluate the conformational changes of the selected drugs in complex with SARS-CoV-2-RBD, Molecular Dynamics (MD) simulations were performed by Gromacs version 5.1.2 for a period of 100 ns. Firstly, topology file and force field parameters of the selected drugs were produced for Gromacs utilizing the PRODRG server (Sch€ uttelkopf & Van Aalten, 2004) . After that, the topology file of the receptor was generated using pdb2gmx and the united-atom GROMOS 96 43A1 force field (Abraham et al., 2015) . In the next step, the receptor/ligand complexes were solvated in a rectangular box with SPC (simple point charge) water molecules with a 10 Å marginal radius (Van Tilbeurgh et al., 1999) . In order to make the simulation system electrically neutral, we replaced solvent molecules with Cl À or Na þ ions. Then, Particle Mesh Ewald (PME) summation was considered to obtain the long range electrostatic interactions (Darden et al., 1993) . The Linear Constraint Solver (LINCS) algorithm was used for covalent bond constraints (Hess et al., 1997) . The system was relaxed with several energy minimization steps. Subsequently, the system was balanced for 500 ps using NVT (constant Number of particles, Volume, and Temperature), following another 500 ps using NPT (constant Number of particles, Pressure, and Temperature) ensembles at 300 K. Thermostat and barostat of the system were adjusted through the V-rescale22 and Parrinello-Rahman barostat algorithms, respectively (Bussi et al., 2007; Parrinello & Rahman, 1981) . Finally, the balanced systems were simulated for a period of 100 ns with 2 fs time steps. The resulted trajectories of the MD simulations were used for further processing. The interaction free energies of each receptor-ligand complex were evaluated using Molecular Mechanics/ Poisson-Boltzmann Surface Area (MM/PBSA) method. In this study, the last 1000 ps of the MD trajectories were used for the calculation of binding free energies, van der Waals energy, electrostatic energy, surface accessible surface area, and polar solvation energy, by using the equations as described by Kumari et al. (2014) . Three dimensional structures of the protein and receptor/ligand complexes were visualized using PyMOL software, also the hydrophobic interactions and hydrogen bonds between protein and compounds were analyzed using discovery studio software (DeLano, 2002; Visualizer, 2017 It is well proved that the active form of SARS-CoV-2 spike (S) glycoprotein contains important structural domains; NTD (Nterminal domain), RBD (receptor binding domain), SD1/SD2 (subdomain1/2), and S2 domain. The RBD of SARS-CoV-2 S has efficient role in the binding to human ACE2 enzyme. The SARS-CoV-2 S is closely related to SARS-CoV S which also binds to ACE2 through RBD (Shang et al., 2020) . According to the previous studies, the ACE2 binding site of SARS-RBD including; Tyr436, Tyr440, Tyr442, Lus443, Lus472, N473, Tyr 475, Gly482, Tyr484, Thr486, Thr487, Gly488, Ile489, Tyr491, and Gln492 residues and the SARS-RDB binding site of ACE2 including; Gln24, Thr27, Phe28, Asp30, Lys31, His34, Asp38, Tyr41, Gln42, Lus45, Lus79, Met82, Tyr83, Gln325, Gly326, Asn330, Lys353, A355 were key residues that were involved in SARS-RBD/ACE2 interactions (Song et al., 2018) . Interestingly, the ACE2 binding site of SARS-CoV-2 -RBD including; Tyr449, Tyr453, Lus455, Phe456, Ala475, Gly476, Phe486, Asn487, Tyr489, Gln493, Gly496, Gln498, Thr500, Asn501, Gly502 and Tyr505 residues and the SARS-CoV-2-RBD binding site of ACE2 including; Ser19, Gln24, Thr27, Phe28, Lys31, His34, Glu35, Glu37, Asp38, Tyr41, Gln42, Lus45, Lus79, Met82, Tyr83, Asn330, Lys353, Gly354, Asp355, Arg357, and Arg393 were recognized as the key residues which were involved in SARS-CoV-2-RBD/ACE2 interactions (Shang et al., 2020) . Consistent with this findings, superimposition of SARS-CoV-2-RBD/ACE2 (PDB ID: 6vw1) and SARS-RBD/ACE2 (PDB ID: 6cag) complexes demonstrated that the secondary structure of residues involved in SARS-CoV-2-RBD binding site were highly similar to SARS-RBD binding site, so that both SARS-CoV-2 and SARS similarly interact with human ACE2 ( Figure 1A ). Taken together, the aforementioned residues of SARS-CoV-2 were selected as the binding pocket which was located in the RBD, so this selected binding pocket was applied for further virtual screening. In order to conduct structure-based virtual screening (SBVS), approved drugs libraries (FDA databases) was used to select potential drugs which may bind to the RBD of SARS-CoV-2. Molecular docking between the binding sites of SARS-CoV-2-RBD and selected libraries was performed by applying AutoDock Vina. Results were sorted based on the binding affinity. The binding affinity threshold was set as À10 kcal/ mol. After analyzing the result, 20 drugs were found against the binding pocket with the binding affinity higher than À10 kcal/mol. To remove drugs with identical structure, they were grouped based on the structure similarity using ChemMine Web Tools. Eventually, a library with 12 hit drugs with highest binding affinity and best conformations were chosen as hit drugs and recruited for further analysis. To select the best lead drugs, interactions between selected hit drugs with SARS-CoV-2-RBD, binding affinity and number of hydrogen bonds as well as interacting residues were studied. Subsequently, drugs with the highest binding affinity and hydrogen bonds were chosen from the hit drugs as the potential lead drugs (Table 1) . As it can be deduced from the results Diammonium Glycyrrhizinate showed the highest binding affinity (-11.7545 kcal/mol), and, Digitoxin, Ivermectin, Rapamycin, Rifaximin, and Amphotericin B followed average ranged from À11.2534 to À10.5021 (kcal/ mol), respectively. Moreover, the two-dimensional structure of the selected lead drugs as well as their predicted binding modes in complex with SARS-CoV-2-RBD were analyzed using pyMOL and discovery studio softwares, respectively. The Diammonium Glycyrrhizinate small molecule consists of two parts: (1) a Glycoside segment, and (2) a triterpene segment ( Figure 2A ). Analysis of the docking results revealed that seven hydrogen bonds were formed between hydroxyl (OH) and oxygen (-O-and C ¼ O) of Glycoside segment of Diammonium Glycyrrhizinate and Gln493 (angle O-H- residues. Also, one hydrogen bond was formed between triterpene segment of Diammonium Glycyrrhizinate and Tyr489 (angle O-H-O ¼ 78.84 , distance ¼ 5.38 Å) residue in the binding pocket of SARS-CoV-2-RBD ( Figure 2B ). In addition, the triterpene segment of the docked compound was embedded in a hydrophobic pocket formed by Tyr489, Asn487 and Phe486 residues ( Figure 2C ). The structure of Digitoxin small molecule also consists of two parts: (1) a Glycoside segment and (2) a Steroid segment ( Figure 3A ). Detailed study of docking results of Digitoxin revealed the presence of six hydrogen bonds formed between hydroxyl (OH) and oxygen of Glycoside segment from Digitoxin and Tyr453 Additionally, one hydrogen bond was formed between hydroxyl of Steroid segment from Digitoxin and Gly485 (angle O-H-O ¼ 83.41 , distance ¼ 2.86 Å) residue in the binding pocket of SARS-CoV-2-RBD ( Figure 3B ). In addition, Steroid segment of the docked compound was embedded in a hydrophobic pocket formed by Tyr489, Glu484, Gly485 and Phe486 residues ( Figure 3C ). The structure of Ivermectin small molecule consists of two parts: (1) a Glycoside segment, and (2) a Macrolide segment ( Figure 4A ). The analysis of Ivermectin demonstrated that five hydrogen bonds were formed between hydroxyl (OH) and oxygen of Ivermectin and Tyr453 (angle HO-H- Figure 4B ). In addition, the docked compound was embedded in a hydrophobic pocket formed by Gln493, Asn501, Tyr505 and Tyr449 residues ( Figure 4C ). The structure of Rapamycin small molecule consists of two parts: (1) a Cyclohexanol segment and (2) a Polyene macrolide segment ( Figure 5A ). Based on the docking results, eight hydrogen bonds were formed between hydroxyl (OH) and oxygen (-O-and C ¼ O) from Rapamycin and Arg408 98 Å residues, and one hydrogen bond was formed between N of piperidine from Rapamycin and Tyr505 (angle N-H-O ¼ 111.26 , distance ¼ 5.60 Å) residue in the binding site of SARS-CoV-2-RBD ( Figure 5B ). In addition, the docked compound was embedded in a hydrophobic pocket formed by Tyr453, Gly416, Val417, Leu455, Phe456 and Gln409 residues ( Figure 5C ). The structure of Rifaximin small molecule consists of two parts: (1) Macrolide segment, and (2) a Naphto Imidazopyridine segment ( Figure 6A ). Four hydrogen bonds were formed between hydroxyl (OH) and oxygen (-O-and residues, and two hydrogen bonds were formed between N of amine and pyridine from Rifaximin and Ser494 (angle N-H-O ¼ 161.38 , distance ¼ 4.27 Å) and Tyr449 (angle N-H-OH ¼129.54 , distance ¼ 4.35 Å) residues in the binding pocket of SARS-CoV-2-RBD ( Figure 6B ). Moreover, the docked compound was embedded in a hydrophobic pocket formed by Gln493, Tyr505, Tyr453 and Tyr449 residues ( Figure 6C ). Like all the previously mentioned drugs, the structure of Amphotericin B small molecule consists of two parts: (1) Figure 7B ). In addition, the docked compound was embedded in a hydrophobic pocket formed by Leu492, Leu455, Phe456 and Tyr449 residues ( Figure 7C ). In this study, several important physiochemical parameters of the selected lead drugs were evaluated using OSIRIS Data Warrior; including molecular weight (MW), LogP, LogS, number of hydrogen bond donors (HBD), number of hydrogen bond acceptors (HBA), polar surface area (PSA), and number of rotatable bonds (Rb). The results showed that the majority of lead drugs possessed following physicochemical properties; MW< 1000, À9