key: cord-0911684-fpz114ka authors: Sethi, Aaftaab; Sanam, Swetha; Munagalasetty, Sharon; Jayanthi, Sivaraman; Alvala, Mallika title: Understanding the role of galectin inhibitors as potential candidates for SARS-CoV-2 spike protein: in silico studies date: 2020-08-13 journal: RSC advances DOI: 10.1039/d0ra04795c sha: d564565110cef710dbca8d0ce2150e7cfd0f961a doc_id: 911684 cord_uid: fpz114ka The Severe Acute Respiratory Syndrome Corona Virus 2 (SARS-CoV-2) has been rapidly transmitting and leaving its footprints across the globe. Stringent measures like complete lockdown and extensive testing have been employed by many countries to slow it down in its tracks until a viable treatment is found. Therefore, in the current scenario, prompt solutions need to be uncovered to tackle the virus. In the present study, 330 galectin inhibitors were tested against SARS-CoV-2 spike (S) protein with the aid of molecular docking and molecular dynamics. Finally, the binding free energy and contributing energies were calculated for 2 top scoring ligands by using MM–GBSA method. Many of the galectin inhibitors displayed high binding score against the S protein. They were found to bind to the site of contact of S protein to ACE2. Thus, they show promise of disrupting the ACE2–S protein binding and prevent the virus from invading the host cell. Among the ligands screened, TD-139, a molecule currently in Phase IIb clinical trials, was found to be a potential hit. The present study paves the way for in vitro and in vivo testing of galectin inhibitors against SARS-CoV-2. In addition, it warrants a swift examination of TD-139 for treating COVID-19. Since the beginning of 2020 the world has found itself in a precarious situation with the emergence of 'public enemy 1, COVID-19' as named by the World Health Organization (WHO). The pandemic which started on 12th December 2019 has now infected more than 4.68 million people across 216 countries and territories as of 21st May 2020, resulting in over 0.3 million deaths and is still rapidly spreading. 1 It was learned based on its phylogenetic characteristic and genetic structure that the causal pathogen of COVID-19 belongs to genera of beta-coronavirus (b-CoV). 2 Human pathogenic Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2), is the 7th human corona virus among the other six related species namely 229E, NL63, OC43, HKU1, Middle East Respiratory Syndrome Coronavirus (MERS-CoV) and SARS-CoV. 3 Unlike the related b-CoVs i.e. SARS-CoV and MERS-CoV, which caused a viral outbreak in November 2002 and September 2012 respectively, SARS-CoV-2 has shown rapid transmission and consequential morbidity. 4 Variety of therapies for SARS-CoV and MERS-CoV and other viral diseases that include current and potential treatments with antiviral drugs, steroids, plasma from recovered patients and Chinese medicine are being tested on COVID-19 patients. 5 Intravenous administration with remdesivir, a nucleotide analog and an antiviral drug developed by Gilead initially for the treatment of diseases caused by Ebola, Marburg, MERS and SARS viruses has been found to be efficacious in an American patient with COVID-19. 6 Baricitinib, interferon-a, lopinavir/ ritonavir, and ribavirin have been suggested as potential therapies for patients with acute respiratory symptoms and interactions are monitored carefully. 7 Despite so many experimental and computational studies currently searching for a potential breakthrough, to date, there is no conrmed effective treatment specically available for COVID- 19 . The therapies against SARS-CoV-2 can be divided into two broad categories. One acting on the human immune system and the other acting on the coronavirus itself. The latter can be further sub-divided into two categoriesone which prevents the viral RNA synthesis and replication and the second which blocks the binding to human cell receptors. 8 Two groups of proteins characterize CoVs; structural proteins such as spike (S), membrane (M), nucleocapsid (N), and envelope (E), in addition to the non-structural proteins, like proteases and RNAdependent RNA polymerase (RdRp). 9,10 The S protein is the sole protein responsible for mediating viral entry into the host cell and thus a crucial recognition factor for attachment to human cellular receptor angiotensin-converting enzyme 2 (ACE2) to gain entry. 11, 12 The S protein protrudes out from the viral surface that gives coronaviruses a crownlike appearance by forming spikes on their surface (Fig. 1) . Many computational approaches have been explored to search for small molecular inhibitors against SARS-CoV-2 main protease and RdRp. [13] [14] [15] [16] Analogous screening of potential drugs against the S protein of SARS-CoV-2 provided small molecular compounds with a high binding affinity. However, most of these compounds were not found to attach to the binding interface of the receptor binding domain (RBD) of S protein and ACE2. Fortunately, hesperidin was predicted to lie on the surface of RBD, and was suggested to disrupt the interaction of ACE2-RBD complex. 8 In another development, Bioxytran, a Biotechnology company based in United States, recently concluded that a galectin inhibitor could bind to the coronavirus spikes and reduce viral load. 17 Encouraged by the potential ability of a glycosylated ligand to cause disruption of ACE2 interaction and ndings of Bioxytran, we decided to test other reported galectin inhibitors like TD139 (currently in Phase IIb clinical trial) and others for their ability to disrupt RBD-ACE2 complex (Fig. 2) . 18 In this work, the binding affinity of several monosaccharide and disaccharide galectin inhibitors with SARS-CoV-2 S protein was evaluated using molecular docking, molecular dynamics and molecular mechanics-generalized born surface area (MM-GBSA) calculations. This study provides an insight into the probable repurposing of galectin inhibitors like TD139 against SARS-CoV-2 and related coronaviruses. TD139, thus, could be used as an inhaled therapeutic for topical lung delivery, providing a valuable and swi therapy to combat COVID-19. All the 330 galectin inhibitors which reported in 28 individually published articles were selected for study. The compounds were sketched by using the 2D-sketcher in Maestro (Schrödinger Suite 2018-4). These compounds were then prepared using the "LigPrep" module by generating low energy ionization and tautomeric states with a pH of 7.4 using OPLS-2005 force eld for minimization. Further, the ligand pdb les were converted to pdbqt format using Open Babel. 47 The QikProp analysis was performed for the prepared This journal is © The Royal Society of Chemistry 2020 RSC Adv., 2020, 10, 29873-29884 | 29875 ligands was performed using admetSAR 48 and Toxicity Estimation Soware Tool (TEST). 49 Target protein, SARS-CoV-2 spike protein (PDB ID: 6M0J) with resolution 2.45Å was obtained from Research Collaboratory for Structural Bioinformatics -Protein Data Bank (RCSB -PDB). 50 All bound waters and cofactors were removed, Kolmann charges were computed, polar hydrogen atoms were subsequently added and the AutoDock atom types were dened using Auto-DockTools (ADT), Graphical User Interface (GUI) of AutoDock implemented in Molecular Graphics Laboratory (MGL) Tools. 51 To generate putative binding poses, we used the AutoDock Vina soware. 52 In general, the docking parameters were kept to the default values. The size of the docking grid was xed at 40Å Â 40Å Â 40Å. The X, Y, and Z coordinates of the grid box were xed at À36.126, 32.573 and 3.383 respectively, thus encompassing the active site. The output le generated was analyzed using Discovery Studio Visualizer. 53 Induced Fit Docking (IFD) gives a more accurate picture of binding by providing the exibility to both ligand and the active site residues. Under the Schrödinger's IFD sampling protocol, 54 the protein was optimized further and a grid was generated by choosing the bound ligand. A soen potential docking with van der Waals scaling of 0.5 each for protein and two top scoring ligands was performed. A maximum of 80 poses were generated for each selected ligand. Prime renement was carried out for residues within 5Å of the ligand. Glide redocking was carried out for structures within 30 kcal mol À1 of the best structure. The most favourable binding conformations of each ligand complex, as measured by the IFDScore, were selected for analysis. In order to better understand the mechanism of interaction between the protein and the ligand, molecular dynamics simulations of 100 ns were performed for the two top scoring ligands using Desmond module of Schrödinger suite (2020-2). 55 Aqueous biological system was built by using OPLS_2005 force eld and TIP3P model was used to stimulate the water molecules. Orthorhombic periodic boundary conditions were set up to specify the shape and size of the repeating unit buffered at 10Å distances. The physiological pH was neutralized by adding 0.15 M NaCl. 300 K temperature and 1.01325 bar pressure was maintained by using Nose-Hoover temperature coupling and Martina-Tobias-Klein method for the constant pressure, respectively. Reversible reference system propagation algorithm (REPSA), a time stepping algorithm was used for near nonbonded (2 fs), far non-bonded (6 fs), and bonded interactions (2 fs). In the end, molecular dynamics simulations (100 ns) were performed for the two top scoring ligands. The MM/GBSA method (Schrödinger Release 2020-2: Prime, Schrödinger, LLC, New York, NY, 2020-2) 56 was used to calculate the binding free energy (DG bind ), for receptor-binder complex systems. One thousand snapshots were taken at 5 ns time intervals throughout the MD simulation trajectory to compute the MM/GBSA free energy difference. Molecular docking is an indispensable technique in the path to drug discovery and it attains an increased signicance when there are no leads/hits available to direct a way forward. 57 Docking studies were carried out to ascertain the binding affinity of the inhibitors. The calculated binding affinity of selected galectin inhibitors with SARS-CoV-2 S protein along with their 2D structures has been tabulated in ESI Table 1 (S1). † The pharmacokinetic properties for all the ligands commonly referred to as ADME (Absorption, Distribution, Metabolism and Excretion) properties were also calculated and can be found in ESI Table 2 (S2). † Based on molecular docking analysis some key observations were made. The thiodigalactoside and disaccharide derivatives with rigid substituents/linkers like triazoles and amides showed better binding energy than inhibitors with freely rotating linkers as observed for ligand 2. This is quite probable as they won't quite t in the binding site due to their exibility. The aryl anking groups as substituents led to better scores due to their engagement in stacking interactions with amino acid side chains present in the S protein. When compared to disaccharides, monosaccharides showed signicantly lesser binding affinity towards the SARS-CoV-2 S protein. This is reasonable since the interface of RBD-ACE2 complex encompasses an appreciably large area and the monosaccharide derivatives are unable to t completely in the cavity of RBD. Majority of ligands screened were able to show a better binding score than hesperidin (À7.3 kcal mol À1 ), which was taken as a reference and is reported to bind to the SARS-CoV-2 RBD-ACE2 interface. Top 4 molecules (all disaccharides) with the best binding affinity were picked to study their interactions with the S protein while one monosaccharide was chosen with an appreciably high binding energy in comparison to other monosaccharides and similar or higher to that of many disaccharides. The 2D and 3D interaction diagrams of ligand 3, 125, 152, 213 and 238 along with their respective docking scores can be found in Fig. 3 & 4 respectively. Among the top ligands, ligand 238 is TD139, a galectin inhibitor, is currently in Phase IIb clinical trials to treat idiopathic pulmonary brosis. 18 As mentioned before, small molecules were able to bind strongly to the S protein, however, they were unable to bind to the interface of RBD-ACE2 complex which is crucial to stop the entry of the virus into the human cell. Thus, it was pertinent to check whether the screened ligands were able to bind to the SARS-CoV-2 RBD. The contact residues of SARS-CoV-2 RBD involved in binding to ACE2 include Lys417, Gly446, Tyr449, Tyr453, Leu455, Phe456, Ala475, Phe486, Asn487, Tyr489, Gln493, Gly496, Gln498, Thr500, Asn501, Gly502 and Tyr505. 58 On examination of the interactions of the screened ligands, the following was observedligand 3 displayed H-bond interactions with Gln493, Gly496 ( Leu455 . In all the ligands, a common feature was observed throughout. The free hydroxyls present on the sugar moiety were involved in the polar interactions while aryl groups linked to the core structure via a linker were involved in varied pi stacking interactions. In addition, another key observation was made regarding the linker. In all, the top scoring ligands, the linker was found to be a rigid and considerably longer. An aryl group linked via a triazole was seen in ligand 125, 213 and 238 while an aryl linked via a triazolo amide was sighted in ligand 152. Ligand 3 on the other hand contained an aryl linked via a benzamide derivative. It is also worth noting that monosaccharides substituted at C2 position or disaccharides substituted at positions other than the C3 (188, 189, 200 and 257); ligands where the pyranose core was observed that among all the top scoring ligands, thiodigalactoside core with a C3 substituted unit linked to a hydrophobic moiety via triazole or amide was a common link (1, 3, 109, 110, 151, 152, 210, 213 , 214 and others). The above evidence clearly suggests that many of the galectin inhibitors screened in the study are able to bind strongly to the This journal is © The Royal Society of Chemistry 2020 RSC Adv., 2020, 10, 29873-29884 | 29879 SARS-CoV-2 RBD and more importantly exhibit interaction with the contact residues of RBD which are responsible for binding to ACE2. Consequently, they have phenomenal potential to disrupt the SARS-CoV-2 and ACE2 interaction and thus prevent the virus from infesting the host cell. All the ve ligands showed good partition coefficient (QPlog P o/w : octanol/water partition coefficient) values which were in the acceptable range (À2 to 6.5). The predicted aqueous solubility (QPlog S) of ligand 125, 152 & TD139 were in the acceptable range (À6.5 to 0.5). The PSA has a great inuence on the oral bioavailabity of the molecules. The optimum value of PSA (7-200 A) was exceeded by three ligands (3, 152 & 213) whereas the other two ligands (125 & TD139) has the PSA value in acceptable range. As the molecules taken for the study are polar saccharide molecules, they have high molecular weight and a greater number of OH groups (H-bond donor & acceptors) which led to the 3 violations in the Lipinski's rule of ve in all ligands except ligand 125 (Table 1) . All further studies hereon were performed This journal is © The Royal Society of Chemistry 2020 RSC Adv., 2020, 10, 29873-29884 | 29881 for ligand 213 (the best ranked ligand) and TD139 (currently in clinical trials and thus having the potential to be marketed early). The in silico toxicity proling was carried out for ligand 213 and TD-139. The median lethal dose (LD 50 ) of a compound describes the ability of a compound to kill 50% of the test population. LD 50 was predicted in silico in a rat model. In our study, we found that ligand 213 was predicted to have a LD 50 value of 2.4156 mol kg À1 while for TD139 this value was slightly higher at 2.5355 mol kg À1 . The human Ether-à-go-go Related Gene (hERG) encodes potassium channels, which are responsible for the normal repolarization of the cardiac action potential. Blockage of these channels could lead to fatal cardiac problems. 59 Therefore, it is essential for any potential drug to not inhibit it. Both the molecules were predicted to be non-inhibitors of hERG. At last, it is important to test for toxicity arising from carcinogenic and mutagenic properties of molecules in addition to its reproductive toxicity. Both, TD139 and ligand 213 were predicted to be non-carcinogens and nonmutagens. However, ligand 213 was predicted to have developmental toxicity. Developmental toxicity includes any effect interfering with normal development, both pre and post birth. This includes embryotoxic effects such as abortion, organ toxicity, teratogenic effects, reduced growth, impaired postnatal mental and physical development among others. 49 Thus, this needs to be carefully evaluated when considering ligand 213's therapeutic potential for a disease. The results of the toxicity proling are summarized in Table 2 . To account for movements in the protein structure, side chains and beyond, we employed the IFD protocol, which is crucial for accurate docking of the ligand. IFD protocol of Schrödinger allows incorporation of receptor exibility by adjusting the receptor structure based on the docked ligand. The IFDScore obtained with the IFD protocol detailed in the Material and methods section was found to be in accordance with the docking scores. Ligand 213 was predicted to have IFDScore of À7904.77 while it was predicted to be À7878.48 for TD139. However, an interesting point to note was while both the ligands bound to a similar site as they did in rigid docking, the interactions were signicantly different (Fig. 5 ). Proteins are dynamic macromolecules under physiological conditions, and therefore the results derived from molecular docking is less persuasive. Molecular docking predicts the spatial orientation of a ligand in the active pocket of the protein. However, molecular dynamics simulations, in addition to the t of the binding pocket, takes into account factors like a conformational stability and exibility. Molecular dynamics trajectory was used to examine the equilibration of dynamics over period of time. Conguration root mean square deviation (RMSD) for TD139-6M0J with respect to its initial structure was found to increase and 6M0J complex was found to converge at 45 ns. The RMSD of the C-alphas, heavy atoms, backbone and side chains of TD139-6M0J and 213-6M0J displayed little uctuations in between, but maintained constant values during the major course of the simulation. At 10 ns, RMSD values of Calphas, backbone side chains and heavy atoms were 1.312, 1.327, 2.169 and 1.594Å, respectively, for TD-139-6M0J. At the end of 100 ns, the observed RMSD values of C-alphas (Fig. 6) , heavy atoms, backbone and side chains were 2.082, 2.424, 2.080 and 2.937Å, respectively. TD139 aer initial uctuations maintained a constant RMSD of 1.6Å. At 10 ns, RMSD values of C-alphas, backbone side chains and heavy atoms were 1.329, 1.801, 1.350 and 2.423Å, respectively, for 213-6M0J. At the end of 100 ns, the observed RMSD values of C-alphas (Fig. 6) , heavy atoms, backbone and side chains were 1.757, 2.213, 1.726 and 2.923Å, respectively. 213 maintained a constant RMSD initially but later showed uctuations between 1.7 and 4.3Å throughout the simulation. Protein-ligand interactions of both the complexes were examined during the course of MD simulation. Fig. 7 demonstrates the type of protein-ligand contacts exhibited by the complexes employed in MD simulation. Mainly, there are three types of protein-ligand interactions: Hbonds, hydrophobic and water bridges. TD139 displayed number of H-bond interactions either directly or mediated by water molecule. Interactions with the residues Asn370, Asn481, Gly482 were seen along with strong and consistent interaction with Glu471, Ile472 and Gln474 during the course of MD simulation. However, no hydrophobic interactions could be observed. In the results obtained from the MD studies, it was found that the hydroxyls of pyranose unit and nitrogen atoms belonging to the triazole were responsible for the strong interactions in TD139. This suggests that hydrogen bonding plays a signicant role in accommodating and stabilizing TD-139 at the binding site. Fig. 8 illustrates timeline representation of ligand-protein contacts including both H-bond and hydrophobic interactions for ligands TD139 and 213 with the receptor during 100 ns MD simulation period. Further, ligand-protein interactions of 6M0J-213 complex were also analyzed. 213 displayed hydrophobic interactions with the residues Tyr449, Leu452, Tyr473 and Phe490. In addition, it displayed hydrogen bond interaction with Thr479, Glu471, Ile472 and Gln474. The initial stability of 6M0J can be assigned to its ability to form Hbond which quickly disappears beyond 50 ns along with the stability of ligand. This further strengthens our conclusion that H-bond interactions are crucial for stabilization of ligand within the binding site. The negative values of DG bind indicate that the selected compounds favourably interact with the receptor. The free energy of binding values and parameters contributing towards it are summarized in Table 3 . The DG bind for ligand 213 is À54.11 kcal mol À1 while it is À45.58 kcal mol À1 for TD139. On further investigating the contribution of each energy term it was apparent that van der Waals interaction play a crucial role. In addition, the contribution of H bond and coulombic term was more for TD139 when compared to ligand 213 which is consistent with the MD studies. The coronavirus (SARS-CoV-2) outbreak came to light on December 31, 2019 when China informed the WHO of a cluster of cases of pneumonia of an unknown cause in the city of Wuhan. Subsequently, this viral infection has indeed got viral and spread to the rest of the world. The present study was aimed at discovering small molecular inhibitors that could inhibit SARS-CoV-2 by binding to the spike protein and inevitably prevent the virus from inltrating the host cell. Virtual screening was carried out with reported galectin inhibitors based on evidence gathered from literature. Several ligands showed promising binding affinity and interaction with residues responsible for binding of ACE2 to SARS-CoV-2 RBD, which strongly suggests the potential of the inhibitors to disrupt the RBD-ACE2 interactions. Among the many promising ligands, 213 and TD139 came to light as an effective inhibitor. Luckily TD139 is also currently in phase IIb clinical trials to treat idiopathic pulmonary brosis and has shown no adverse side effect. Thus, the results of our study warrants the testing of TD-139 and related galectin inhibitors for in vitro and in vivo testing against SARS-CoV-2 and related coronaviruses to halt the spread of the virus in present and the future. There are no conicts to declare. World Health Organization Features, evaluation and treatment coronavirus (COVID-19) Bioxytran Releases Details on a Novel Carbohydrate Galectin Inhibitor designed to Eliminate COVID-19 in Patients A study to test the efficacy and safety of inhaled TD139 in subjects with idiopathic pulmonary brosis Crystal structure of SARS-CoV-2 spike receptor-binding domain bound with ACE-2 Schrödinger Release 2020-2: Induced Fit Docking protocol Schrödinger Release 2020-2: Prime, Schrödinger, LLC Authors are thankful to DoP, Ministry of Chemicals & Fertilizers, Govt. of India, New Delhi, for the award of NIPER fellowship and supporting the work. Authors are also thankful to Vinod Devaraji, Schrodinger, Bengaluru for molecular dynamics studies.