key: cord-0976588-coxfvt8t authors: Sharma, Kedar; Morla, Sudhir; Goyal, Arun; Kumar, Sachin title: Computational guided drug repurposing for targeting 2′-O-ribose methyltransferase of SARS-CoV-2 date: 2020-07-29 journal: Life Sci DOI: 10.1016/j.lfs.2020.118169 sha: 8fce9e52062e1bfef430ea53fe8de66ea46a8cba doc_id: 976588 cord_uid: coxfvt8t AIMS: The recent outbreak of pandemic severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2) has led the world towards a global health emergency. Currently, no proper medicine or effective treatment strategies are available; therefore, repurposing of FDA approved drugs may play an important role in overcoming the situation. MATERIALS AND METHODS: The SARS-CoV-2 genome encodes for 2-O-methyltransferase (2′OMTase), which plays a key role in methylation of viral RNA for evading host immune system. In the present study, the protein sequence of 2′OMTase of SARS-CoV-2 was analyzed, and its structure was modeled by a comparative modeling approach and validated. The library of 3000 drugs was screened against the active site of 2′OMTase followed by re-docking analysis. The apo and ligand-bound 2′OMTase were further validated and analyzed by using molecular dynamics simulation. KEY FINDINGS: The modeled structure displayed the conserved characteristic fold of class I MTase family. The quality assessment analysis by SAVES server reveals that the modeled structure follows protein folding rules and of excellent quality. The docking analysis displayed that the active site of 2′OMTase accommodates an array of drugs, which includes alkaloids, antivirals, cardiac glycosides, anticancer, steroids, and other drugs. The redocking and MD simulation analysis of the best 5 FDA approved drugs reveals that these drugs form a stable conformation with the 2′OMTase. SIGNIFICANCE: The results suggested that these drugs may be used as potential inhibitors for 2′OMTase for combating the SARS-CoV-2 infection. The recent outbreak of severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2) in Wuhan, China, in December 2019, has infected over 69, 16,233 and death of over 4,00,135 (Dated June, 7, 2020) worldwide since its outbreak. The number of cases is increasing every hour, and due to this reason, the World Health Organization (WHO) has declared it as a global public health emergency. The SARS-CoV-2 is a positive-stranded RNA virus with a genome size of 29.9 kb encoding 12 open reading frames (ORFs) and two untranslated regions of 254 and 229 bp at 5' and 3' ends, respectively (Chan et al. , 2020) . The molecular organization of the coronavirus genome in 5'→3' direction encodes several proteins. The genome includes nonstructural and structural proteins. The structural proteins include spike (S), envelope (E), matrix (M), and nucleocapsid (N) proteins (Chan et al., 2020) . The ~66% of the genome at the 5'end contains two overlapping ORFs (ORF1a and ORF1b) encoding polyproteins 1a (PP1a) and 1ab (PP1ab). These polyproteins are processed into 16 functional non-structural proteins by papain-like protease and 3C-like protease, which form the replication transcription complex (RTC). These non-structural proteins include proteases and RNA processing enzymes (Cui et al. , 2019 , Song et al. , 2019 , Woo et al. , 2010 . These RNA processing enzymes are uncommon and include exoribonuclease, guanine N7methyltransferase (N7-MTase), and 2'O-methyltransferase. Approximately, 33% of the genome at 3' end encodes four major structural proteins S, E, M, and N. These structural proteins play a crucial role in the formation of the new virus particle and its release (Rota et al. , 2003) . Most viral and eukaryotic cellular mRNA containing polyadenine region at 3' terminal and a 5'terminal cap structure of N7-methyl-guanine moiety connected through 5'-5' triphosphate bridge to the first transcribed nucleotide (Bouvet et al. , 2012, Zeng et al. , J o u r n a l P r e -p r o o f 4 2016). The 5' terminal cap structure of mRNA plays a key role in protecting its degradation by 5'-exonuclease and efficient translation. The capping of mRNA (Cap 0 ) is catalyzed by RNA triphosphatase, a guanylyltransferase and guanine-N7-methyltransferase (N7-MTase) in a sequential reaction. The capping of most viruses and eukaryotic mRNA is further methylated at the 2-OH position of the first nucleotide by ribose 2'-O-methyltransferase (2'-O-MTase) enzyme to form the cap 1 structure. The 2'O-methylation of viral mRNA cap also helps in protecting RNA recognition as "nonself" by host innate immune system, which makes it a potential target for antiviral drug development (Zeng et al., 2016) . The design of antiviral drugs with proven efficacy or vaccines against SARS-CoV-2 is limited due to the poor understanding of its molecular mechanisms of infection. There is a need to design an inhibitor for blocking the SARS-CoV-2 replication within the cell. The repurposing of the approved drugs is an alternate way to find out a solution. Recently, clinicians have used the antimalarial and anti-HIV drugs for the treatment of SARS-CoV-2; however, the possible mechanism involved in its inhibition is not known (Colson et al. , 2020 , Lim et al. , 2020 , Mitja and Clotet, 2020 . In the present study, the 2'-O-methyltransferase (2'OMTase) structure was modeled using a comparative modeling approach and further validated. Besides, the FDA approved drugs were screened against 2'OMTase for selecting the potential candidate for drug repurposing using virtual drug screening followed by redocking analysis using Autodock tool. The selected best FDA approved drugs were further validated by performing the molecular dynamics simulation. J o u r n a l P r e -p r o o f 5 from the NCBI Database (https://www.ncbi.nlm.nih.gov/protein/1809484466). The retrieved sequence was pairwise aligned with the 2'-O-ribose methyltransferase (NCBI RefSeq ID: YP_009725311.1) encoding region of SARS-CoV-2 isolate from Wuhan to identify the 2'-O-2'OMTase. PSI-BLAST webserver was used to search for 2'OMTase domains from the different strains across the SwissProt database. 2'OMTase sequence was aligned with its homologous characterized domains with a query coverage of 100% displaying the sequence identity in the range of 53%-93% by T-Coffee Multiple sequence alignment web server. Furthermore, the results were visualized and analyzed for the identification of conserved and semi-conserved amino acid residues by Esprit 3.0 webserver. The structure modeling of 2'OMTase was performed by using the comparative modeling approach through Modeller v9.23 Program (Eswar et al. , 2006) . The homologous protein structures of 2'OMTase with PDB ID: 3R24 (Chen et al. , 2011 , Decroly et al. , 2011 , 2XYR (Chen et al., 2011 , Decroly et al., 2011 from human SARS-CoV and 2'OMTase with PDB ID 5YN5 from Human betacoronavirus 2c EMC/2012 [Unpublished] were retrieved from Protein Data Bank (PDB). The structures mentioned above were aligned by running 'salign( )' command to generate the structure-based multiple sequence alignment. After that, the generated Multiple Sequence Alignment (MSA) was further aligned with the 2'OMTase query sequence to cover the whole sequence for model building. A total of 20 independent models of 2'OMTase were built by running the 'automodel' command, and the best two models displaying least DOPE scores were selected for energy minimization by YASARA webserver. The energy minimized models of 2'OMTase were evaluated by using several programs viz, ProCheck (Thornton, 1993) , ERRAT (Colovos and Yeates, 1993) and VERIFY3D (Eisenberg et al., 1997) To identify the potential antivirals drugs candidate for repurposing against 2'OMTase, showing the free energy of binding (ΔG) between -10 kcal/mol -8.7 kcal/mol were selected for further analysis. The free energy of binding of -8.3 kcal/mol for sinefungin was taken as reference. The selected best FDA approved drugs obtained from virtual screening were redocked The 2D interaction of protein and ligand was generated by using the LigPlot+ program. The conformational stability and behavioral dynamics of apo and ligand-bound 2'OMTase structures were determined by molecular dynamic (MD) simulation. GROMOS96 43a1 force field compiled with GROMACS v5.1.4 software program installed in the High-Performance Computing facility (Param-Ishan) available at Indian Institute of Technology Guwahati, India was used. The topology of FDA approved drug viz. Sinefungin, dihydroergotamine, digitoxin, irinotecan, and teniposide were generated by using PRODRG server (http://prodrg1.dyndns.org/submit.html). The apo form and ligand-bound complex was placed in a triclinic box with a distance of 1.2 nm and solvated with a single point charge (SPC) waters. The system was neutralized by adding 2 Clions as counterions, and the system was energy minimized by using the Steepest descent method using cut-off up to 239 kcal mol −1 . The system was equilibrated in two steps (NVT ensemble (constant number of particles, volume and temperature) followed by NPT ensemble (constant number of particles, pressure and temperature) for 500 ps each. The equilibrated system was simulated for 50 ns to run the final MD simulation, and the trajectories were recorded at 10 ps interval. The final J o u r n a l P r e -p r o o f 8 MD run was analysed by using gmx rms for the estimation of root-mean-square deviation (RMSD), gmx rmsf for room mean square function (RMSF) calculation. The radius of gyration (Rg) and solvent accessible surface area (SASA) were determined by using gmx gyrate, and gmx sasa commands, respectively. The XMGRACE software (http://plasmagate.weizmann.ac.il/Grace/) was used for data plotting and analysis. The pairwise sequence alignment analysis of ORF1ab polyprotein of Indian isolate of SARS-CoV-2 with isolate from Wuhan revealed that both the sequences of the 2'-O-ribose methyltransferase domain are identical without any mutation. The PSI-BLAST analysis of 2'OMTase displayed the amino acid sequence identity with previously characterized 2'OMTase of CoV isolated from different strains ( Table 1 ). The 2'OMTase sequences from different strains mentioned in Table 1 were retrieved from Uniprot Database and aligned using the T-Coffee program revealed that the catalytic amino acid residues Lys46, Asp130, Lys170, and Glu203 are conserved in all the strains (Fig. 1 ). The comparative modeled structure of 2'OMTase ( Fig. 2A ) displayed the characteristic fold comprising seven stranded β-sheets surrounded by α-helices and loops ( Fig. 2A) . A similar type of structural organization is conserved and reported for class I MTase family (Chen et al., 2011 , Decroly et al., 2011 . The different programs available on the SAVES server were used to assess the quality of 2'OMTase energy minimized modeled structure. The Ramachandran plot analysis by Procheck server displayed that 266 amino acid residues i.e., non-glycine and proline, follow the stereochemical properties. The analysis showed that 91.4% amino acid residues of 266 lies in the favored region, and 8.6% amino acid residues accommodate in the allowed region (Fig. 2B ). The absence of amino acid J o u r n a l P r e -p r o o f 9 residues in the generously allowed and disallowed regions of the Ramachandran plot demonstrated that the modeled structure of 2'OMTase follows the stereochemical property rule and follows all the possible dihedral, phi (ϕ) and psi (ψ) angle values. The modeled structure on further evaluation by the VERIFY3D plot, revealed that 92.28% of amino acid residues have an overall average 3D-1D score ≥ 0.2 (Fig. 2C) . Similarly, the 2'OMTase energy minimized structure assessed by ERRAT Plot displayed the quality factor of 96.20%, which further confirmed that the modeling of 2'OMTase structure is error-free (Fig. 2D) . Similarly, ProSA analysis of 2'OMTase displayed the Z score of -7.4, confirming the modeled structure accommodated in the X-ray zone (Fig. 2E ). PDBSum analysis of the modeled structure of 2'OMTase displayed that it contains 13 α-helices and 13 β-strands as secondary structure elements (Fig. 2F ). The quality assessment parameters of the modeled structure of 2'OMTase were superior to the experimentally determined homologous crystal structure of 2'OMTase, which reveals that the selected modeled structure is the best model for further studies. The 2'OMTase and its homologous proteins were subjected to CastP analysis for the determination of active site volume and to predict the open or close conformation. The active site volume of 2'OMTase from SARS-CoV-2 was found to be 526.70 Å 3 . The active site volume of homologous 2'OMTase for PDB ID 3R24_A was 946.75 Å 3 and for PDB ID 2XYR was 537.38 Å 3 . The Active site volume of 2'OMTase from SARS-CoV-2 matched well with the homologous protein of SARS-CoV and confirmed that the modeled protein has closed conformation. The high throughput virtual screening analysis of FDA approved drugs library against Table 2 . Each ligand screened through virtual screening by PyRx displayed nine different conformations, and the ligands forming interaction with catalytic tetrad were selected for further analysis. These selected top hits included antiviral, alkaloids based drugs, cardiac glycoside, anticancer drugs, steroid-based drugs, and other drugs, etc. These selected FDA approved drugs were further employed for re-docking analysis with the 2'OMTase to discover the potential inhibitor(s) by using Autodock tool. The re-docking analysis of the top twenty selected screened FDA approved drugs was carried out by AutoDock 4.2.6. The molecular docking analysis revealed that the top ten FDA approved drugs viz. sinefungin, dihydroergotamine, digitoxin, irinotecan, and teniposide displayed the free energy of binding in the range of -9.8 kcal/mol to -7.8 kcal/mol and predicted inhibitory constant in the range of 20 nM to 3120 nM (Table 3 ). The amino acid residues involved in the formation of hydrogen bonds, other interactions, the free energy of binding and inhibitory constant are shown in Table 3 . The protein-ligand complex of selected FDA approved drugs and 2'OMTase were further analyzed by using PyMOL v2.3 software and LigPlot+ software. The drug, sinefungin, a natural nucleoside, interacts with 2'OMTase through four hydrogen bonds formed by Gly71, Asp130, Tyr132, and Lys170 (Fig. 3A) . The sinefungin interacts with Tyr47, Gly73, Ser74, Asn99, Leu100, Met131, Asp133, Pro134, and Phe149 amino acid residues and forms hydrophobic/Van der Waals interaction (Fig. 3A) . The amino acid residue Asp130 and Lys170 forming hydrogen bonds with sinefungin, are the part of catalytic tetrad and plays a key role in the formation of the catalytic cleft. Similarly, digitoxin, a cardiac glycoside interacts with 2'OMTase by forming four hydrogen bonds by Lys24, Tyr30, Tyr132 and His174 followed by the formation of J o u r n a l P r e -p r o o f 11 hydrophobic interaction through Asn43, Gly71, Gly73, Ser74, Asp75, Asp99, Leu100, Asp130, Met131, Asp133, Pro134, Thr136, Phe149, Lys170, Glu173, Ser202 and Glu203 (Fig. 3B) . The amino acid residues of catalytic tetrad viz., Asp130, Lys170, and Glu203 involved in hydrophobic interactions revealed that these residues are further involved in stabilizing the ligand in the active site region (Fig. 3B) . The alkaloid drug, dihydroergotamine, forms hydrogen bonds with Tyr47, Gly71, and Lys170 and forms hydrophobic interaction with Asn43, Gly73, Ser74, Asp75, Lys76, Pro80, Asp99, Leu100, Asn101, Asp130, Met131, Tyr132, Pro134 and Glu203 of 2'OMTase ( Fig. 3C ). The three amino acid residues of catalytic tetrad stabilize the drug in the catalytic cleft by hydrogen bond and hydrophobic interactions. Anticancer drug irinotecan (Fig. 3D ) forms hydrogen bonds with Tyr47 and Ser200 and involved in hydrophobic interactions through Tyr30, Gly31, Asp32, Ser33, Met42, Asn43, Lys46, Gly71, Gly73, Ser74, Asp99, Leu100, Asp130, Met131, Tyr132, Lys170, Val197, Asn198 and Ser201 of 2'OMTase. Teniposide interacts with Asp99, Leu100, Asp130, Lys137, and Glu203 and forms hydrogen bond (Fig. 3E) . Besides, the amino acid residues Gly71, Gly73, Ser74, Asn101, Met131, Tyr132, Pro134, Lys170, and Thr172 forms hydrophobic interaction. These drugs are stabilized in the catalytic cleft by catalytic residues Asp130, Lys170, and Glu203. The selected FDA approved drugs completely accommodated the catalytic cleft and interaction stabilized mainly by three amino acid residues viz. Asp130, Lys170, and Glu203 of the catalytic tetrad. This may result in preventing the 5' Cap modification by inhibiting the methyltransferase enzyme activity. The molecular dynamics (MD) simulation analysis of both apo and ligand complex was analyzed for the 50 ns to understand the stability and conformational dynamic behavior. The radius of gyration analysis of apo 2'OMTase for the estimation of global compactness displayed the average Rg values of 1.90 nm up to 50 ns; however, a small increment in the Rg value was observed between 18 ns to 23 ns (Fig. 4D) . The ligand-bound 2'OMTase complexes displayed the fluctuation in Rg value for an initial 3 ns, and thereafter the Rg value remains stable in the range of 1.85-1.89 nm throughout the simulation (Fig. 4D) . The Rg analysis concludes that the protein remains in stable conformation resulted in achieving global compactness. Similarly, the estimated average solvent accessible surface area (SASA) for apo form of 2'OMTase up to 25 ns was 160 nm 2 and thereafter, it reduces to 150 nm 2 and remained unchanged throughout the MD simulation (Fig. 4E) . The ligand-bound complex of 2'OMTase also displayed a similar type of pattern, and the average SASA was found to be in the range of 133-141 nm 2 . The ligand-bound complex of 2'OMTase after the MD simulation was further subjected to the estimation of the interaction energy between protein and ligand molecule. The interaction energy of 2'OMTase for sinefungin was -54.35 kcal/mol, for digitoxin -69.47 kcal/mol, for dihydroergotamine -92.74 kcal/mol, irinotecan -75.34 kcal/mol, and teniposide -97.88 kcal/mol. The interaction energy calculation reveals that the digitoxin, dihydroergotamine, irinotecan, and teniposide has more affinity toward 2'OMTase and forms a stable conformation as compared to the reference ligand i.e. sinefungin. Methylation of mRNA Cap present at 5' terminal is a crucial process, which protects its recognition as a foreign particle from the host immune system. 2'OMTase is a crucial enzyme involved in the methylation of 5'capped mRNA, which makes this enzyme a potential drug target for identification and screening of FDA approved drugs for repurposing. The modeled J o u r n a l P r e -p r o o f 14 structure of the 2'OMTase protein displayed the conserved fold of the Class I MTase enzyme (Chen et al., 2011 , Decroly et al., 2011 . The quality assessment of the modeled structure by SAVES server reveals that the 100% non-glycine and non-proline residues attain the favorable dihedral, phi (ϕ), and psi (ψ) angles. The non-bonded interaction analysis of different types of atoms in the modeled structure when compared with the experimentally determined structure also suggested that the modeled structure is of excellent quality and can be used for further analysis. The multiple sequence alignment analysis of 2'OMTase with other homologous proteins revealed that the residues involved in catalytic tetrad formation i.e., Lys46, Asp130, Lys170, and Glu203, are strictly conserved in all the CoV strains (Chen et al., 2011 , Decroly et al., 2011 . displayed that the regions of the protein ranging from 16-40, 72-80, 99-110, 132-147, 212-219, and 242-268 are mainly composed of β and γ turns and provide flexibility. The interaction energy analysis of the selected drugs revealed that digitoxin, dihydroergotamine, irinotecan, and teniposide have more affinity towards 2'OMTase and forms a stable conformation as compared to the earlier known inhibitor, i.e. sinefungin. The nucleoside inhibitor, sinefungine, suppressed the replication of feline herpesvirus type 1 in the host feline kidney cells (Kuroda et al., 2019) . The alkaloid based drug, dihydroergotamine was not used as antiviral but for the migraine problems in HIV patients (Tfelt-Hansen and Koehler, 2008) . Cardiac glycosides such as digitoxin inhibit the human cytomegalovirus (HCMV) replication by inducing AMP-activated protein kinase (AMPK) activity and autophagy flux through the Na + , K +/ ATPase α1 subunit activation (Mukhopadhyay et al. , 2018) . Similarly, teniposide, the anticancer compound, inhibits the replication of SV40 by blocking the topoisomerase type II activity (Richter et al. , 1987) . To date, all these potential FDA approved drugs are used for different applications and not used for the inhibition of 2'OMTase. The strong predicted binding of these approved drugs at the active site of 2'OMTase displayed that they can be used for inactivating the enzyme. The SARS-CoV-2 infection has created a global pandemic and an urgent need for therapeutics. Repurposing of an existing drug based on the bioinformatics study could be a fast way to understand its antiviral efficacy against SARS-CoV-2. Although the understanding of SARS-CoV-2 molecular biology is limited, and there is a pressing need to design some potential therapeutics. The identification of some novel drugs that can display significant effect in silico can be further explored to validate its use clinically to control the disease and associated mortality. This study may pave the way in understanding the biology of SARS-CoV-2 and finding new possible drugs. Table 1 . The conserved amino acid residues are displayed in the red color background and semi-conserved residues are shown in red color. Genomic characterization of the 2019 novel human-pathogenic coronavirus isolated from a patient with atypical pneumonia after visiting Wuhan Biochemical and structural insights into the mechanisms of SARS coronavirus RNA ribose 2'-O-methylation by nsp16/nsp10 protein complex Verification of protein structures: patterns of nonbonded atomic interactions Chloroquine and hydroxychloroquine as available weapons to fight COVID-19 Origin and evolution of pathogenic coronaviruses Crystal structure and functional analysis of the SARS-coronavirus RNA cap 2'-O-methyltransferase nsp10/nsp16 complex VERIFY3D: assessment of protein models with three-dimensional profiles Comparative protein structure modeling using Modeller Antiviral effect of sinefungin on in vitro growth of feline herpesvirus type 1 Case of the Index Patient Who Caused Tertiary Transmission of COVID-19 Infection in Korea: the Application of Lopinavir/Ritonavir for the Treatment of COVID-19 Infected Pneumonia Monitored by Quantitative RT-PCR Use of antiviral drugs to reduce COVID-19 transmission Digitoxin Suppresses Human Cytomegalovirus Replication via Na(+), K(+)/ATPase alpha1 Subunit-Dependent AMP-Activated Protein Kinase and Autophagy Activation Effects of VM26 (teniposide), a specific inhibitor of type II DNA topoisomerase, on SV40 DNA replication in vivo Characterization of a novel coronavirus associated with severe acute respiratory syndrome From SARS to MERS, Thrusting Coronaviruses into the Spotlight History of the use of ergotamine and dihydroergotamine in migraine from 1906 and onward PROCHECK: a program to check the stereochemical quality of protein structures Coronavirus genomics and bioinformatics analysis Identification and Characterization of a Ribose 2'-O-Methyltransferase Encoded by the Ronivirus Branch of Nidovirales The virus research in our lab is currently supported by the Department of Biotechnology, Government of India (BT/562/NE/U-Excel/2016 , and BT/PR24308/NER/95/644/2017). The authors do not have any potential conflict of interest.