Computationally Guided Identification of Novel Mycobacterium

Jan 26, 2016 - Mycobacterium tuberculosis (Mtb) infections are causing serious health concerns worldwide. Antituberculosis drug resistance threatens t...
1 downloads 14 Views 5MB Size
Research Article pubs.acs.org/acscombsci

Computationally Guided Identification of Novel Mycobacterium tuberculosis GlmU Inhibitory Leads, Their Optimization, and in Vitro Validation Rukmankesh Mehra,† Chitra Rani,‡,§ Priya Mahajan,†,§ Ram Ashrey Vishwakarma,† Inshad Ali Khan,*,‡,§ and Amit Nargotra*,†,§ †

Discovery Informatics, CSIR-Indian Institute of Integrative Medicine, Canal Road, Jammu 180001, India Clinical Microbiology, CSIR-Indian Institute of Integrative Medicine, Canal Road, Jammu 180001, India § Academy of Scientific and Innovative Research, CSIR-Indian Institute of Integrative Medicine, Canal Road, Jammu 180001, India ‡

S Supporting Information *

ABSTRACT: Mycobacterium tuberculosis (Mtb) infections are causing serious health concerns worldwide. Antituberculosis drug resistance threatens the current therapies and causes further need to develop effective antituberculosis therapy. GlmU represents an interesting target for developing novel Mtb drug candidates. It is a bifunctional acetyltransferase/uridyltransferase enzyme that catalyzes the biosynthesis of UDP-N-acetyl-glucosamine (UDP-GlcNAc) from glucosamine-1-phosphate (GlcN-1-P). UDPGlcNAc is a substrate for the biosynthesis of lipopolysaccharide and peptidoglycan that are constituents of the bacterial cell wall. In the current study, structure and ligand based computational models were developed and rationally applied to screen a druglike compound repository of 20 000 compounds procured from ChemBridge DIVERSet database for the identification of probable inhibitors of Mtb GlmU. The in vitro evaluation of the in silico identified inhibitor candidates resulted in the identification of 15 inhibitory leads of this target. Literature search of these leads through SciFinder and their similarity analysis with the PubChem training data set (AID 1376) revealed the structural novelty of these hits with respect to Mtb GlmU. IC50 of the most potent identified inhibitory lead (5810599) was found to be 9.018 ± 0.04 μM. Molecular dynamics (MD) simulation of this inhibitory lead (5810599) in complex with protein affirms the stability of the lead within the binding pocket and also emphasizes on the key interactive residues for further designing. Binding site analysis of the acetyltransferase pocket with respect to the identified structural moieties provides a thorough analysis for carrying out the lead optimization studies. KEYWORDS: Mycobacterium tuberculosis, in silico, computational screening, GlmU inhibitors, acetyltransferase inhibitors, QSAR

1. INTRODUCTION Tuberculosis (TB) remains a major public health problem worldwide. Mycobacterium tuberculosis (Mtb) infection is a major cause of morbidity and mortality in large parts of the world and is considered to be one of the most important global health problems.1 The recent rise of multidrug-resistant tuberculosis and its association with HIV infection poses a challenging health concern. Therefore, there exists a pressing requirement to identify novel drug targets and develop new antituberculosis drugs that will be effective against multidrug-resistant-tuberculosis. The bacterial GlmU protein, involved in peptidoglycan and lipopolysaccharide biosynthesis, has recently been identified as an important drug target.2,3 GlmU has been identified as essential for optimal growth of M. tuberculosis.4 The GlmU protein is a © XXXX American Chemical Society

bifunctional acetyltransferase/uridyltransferase enzyme that catalyzes the conversion of glucosamine-1-phosphate (GlcN-1P) and acetyl coenzyme A (acetyl CoA) to N-acetyl-glucosamine-1-phosphate (GlcNAc-1-P) and coenzyme A (CoA) at the C-terminal (acetyltransferase) domain followed by conversion of GlcNAc-1-P and UTP to UDP-N-acetyl-glucosamine (UDP-GlcNAc) at the N-terminal (uridyltransferase) domain.5,6 The second step is present in bacteria as well as humans, whereas the first step is present only in bacteria,2 which makes the first step a suitable target for designing selective inhibitors of Mtb GlmU acetyltransferase site. Received: January 30, 2015 Revised: December 4, 2015

A

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

Figure 1. Structure of Mtb GlmU (PDB ID 3ST8). (A) Monomer unit of Mtb GlmU protein. N-terminal uridyltransferase and C-terminal acetyltransferase sites are represented by bound ligands in CPK representation. (B) Trimer unit of Mtb GlmU representing the predicted binding sites. The sites were predicted using the SiteMap tool. The two binding sites are represented in mesh, where red color represents hydrogen bond acceptor regions, blue color shows hydrogen bond donor regions and yellow color shows hydrophobic areas. The chain names A−C were assigned to the three monomer units of the trimer represented as green, cyan, and gray colors, respectively.

Figure 2. ROC graphs of the docking results for assessing the reliability of the docking protocol. The graphs were plotted between sensitivity (true positive rate) versus 1-specificity (false positive rate). Active compounds having IC50 < 30 μM were mixed with 1000 decoy molecules (Schrodinger decoy library) and 100 confirmed inactives, and docking studies were carried out. (A) ROC plot of HTVS docking showed AUC value of 0.70. (B) ROC plot of SP docking showed AUC value of 0.68. (C) ROC plot of XP docking showed AUC value of 0.67.

domain (residues 263−478), which has a left-handed β-helix (LβH) fold. The trimer is the biologically active unit, where the C-terminal active site is formed by the contribution from all three monomer units. The trimer interface is dominated by contacts

The structure of GlmU monomer shows two distinct domains (as shown in Figure 1): an N-terminal domain (residues 1−262), which has a typical uridyltransferase fold based on a dinucleotidebinding Rossmann fold and a C-terminal acetyltransferase B

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science between the C-terminal β-helical domains of the three monomers.7−10 The data of inhibitors of GlmU has been reported in the literature.11−20 PubChem Bioassay AID 1376 contains 533 Mtb GlmU inhibitors, out of which the IC50 values of 125 compounds are known.11,21 These are acetyltransferase inhibitors of GlmU that exhibit a wide range of activity (IC50 range = ∼1−9999 μM) and structural diversity. Li et al. synthesized the analogs of GlcN6-P and GlcN-1-P and evaluated their inhibitory activity against GlmU and GlmM of Mtb.13 They identified an analog of GlcN-1P as a promising inhibitor of Mtb GlmU.13 Further, Tran et al. carried the inhibition studies of the compounds against the uridyltransferase activity of Mtb GlmU.14 They developed a number of aminoquinazoline-based compounds as potential inhibitors with most active among them exhibited an IC50 value of 74 μM for the uridyltransferase activity of Mtb GlmU.14 Besides these, other GlmU inhibitors are thiol-specific reagents such as iodoacetamide (IDA) and N-substituted maleimides, which include N-ethylmaleimide (NEM), N-phenylmaleimide, N,N′-(1,2-phenylene) dimaleimide (oPDM), and N-(1-pyrenyl)-maleimide (PyrM).12 The screening of large libraries and the identification of inhibitors by applying experimental techniques is an expensive and time-consuming task. Thus, it is useful to develop some computational models that can help in predicting inhibitors for therapeutically relevant targets. In this study, such an attempt has been made to develop computational models and screen a compound library. First, similarity search of the library was performed using GlmU inhibitors, and then based on the data of PubChem AID 1376, two QSAR and two pharmacophore models were developed. These computational models were suitably combined to filter a large compound repository to identify potential inhibitor candidates of Mtb GlmU. The inhibitor candidates obtained from in silico filtering were further screened by in vitro plate assay using cloned and purified Mtb GlmU for the identification of inhibitory leads of acetyltransferase activity. The binding site analysis of the target site was also carried in order to provide a robust strategy for lead optimization. Similar kind of approach involving in silico screening, followed by in vitro screening has been applied previously by us for the identification of E. coli GlmU inhibitors.22

monomer units with contribution from left-handed alpha helices of two subunits, insertion loop of one subunit and C-terminal extension of the third subunit (Figure 1B). Therefore, it is important that the protein should be in trimeric form for docking in the acetyltransferase site. However, no trimer Mtb GlmU crystal structure was present in PDB and all the Mtb GlmU crystallized structures were monomers. Therefore, for docking study the predicted trimeric biological assembly of 3ST8 having a resolution of 1.980 Å was retrieved from PDB where the coenzyme-A was bound at the acetyltransferase site. The protein structure was prepared using Protein Preparation Wizard of Schrodinger.24 Binding pockets on this trimer Mtb GlmU were predicted using SiteMap25 that revealed the presence of two main pockets (Figure 1B). These pockets correspond to the uridyltransferase and acetyltransferase sites. Since the C-terminal acetyltransferase step is present in bacteria and absent in humans, this site was the target of interest for developing selective inhibitors of Mtb GlmU acetyltransferase site. A 3D grid was generated around the bound coenzyme-A present at the acetyltransferase site of the prepared protein trimer of Mtb GlmU (PDB ID 3ST8). In docking analysis, the enrichment parameters that are calculated with most active compounds in the data set makes most sense rather than the less active compounds. These metrics give a way of knowing how well the process has picked the most active compounds from the given set, not compounds that are less active. Hence, in order to select the most active reasonable number of compounds from the data set (PubChem AID 1376), the compounds having IC50 < 30 μM were treated to be the active data set for docking studies (S.No. 1−16 in Table S1). This active data set of 16 compounds was mixed with 1000 drug-like decoy molecules (Schrodinger decoy library) and 100 confirmed inactives (from PubChem Bioassay AID 1376) to form the ligand data set for docking studies. All the ligands were prepared using LigPrep.26 Flexible docking of this data set was performed using high-throughput virtual screening (HTVS), standard precision (SP) and extra precision (XP) scoring functions of Glide.27−30 The analysis of the docking results was performed so as to determine the enrichment of actives in a mixture of actives, inactives and decoys. Various enrichment parameters were calculated using Enrichment Calculator tool of the Schrodinger software.24 Enrichment of the actives in the top percentages of hits were analyzed along with the Receiver Operator Characteristic (ROC) area under curve and BEDROC (α = 160.9) analysis.31 The value of ROC area under curve lies between 1 and 0, where 1 reflects the ideal screening protocol and a value of 0.5 represents the random behavior. Closer the value of ROC to 1 the better is the efficiency of the method in distinguishing actives from decoys. BEDROC (α = 160.9) is similar to the ROC, where a larger α value provides more weight to the beginning of the ranking. The value α = 160.9 depicts that approximately 80% of the score is derived from the top 1% of results. 2.3. Ligand-Based Approach for Developing Filtering Criteria. 2.3.1. Similarity Search. Chemical or molecular similarity refers to the similarity between the chemical elements of a compound and its biological activity. The chemical similarity is one of the most important concepts in Chemoinformatics particularly for drug designing.32 This concept has proven successful and has been widely applied for identification of novel inhibitors of various targets of biological importance.33−38 Using this criterion as the basis, similarity search on the 20 000 compound library was performed using 7 most active inhibitors of PubChem AID 1376 (Table S1). Since the focus in the present

2. EXPERIMENTAL PROCEDURES 2.1. Data Set. A data set of 125 GlmU inhibitors was retrieved from PubChem Bioassay AID 137611,21 with IC50 values against Mtb GlmU (Table S1). These inhibitors were structurally diverse23 and showed a wide range of activity (IC50 = ∼1−9999 μM). These structures were used to develop various methods for in silico filtering of the drug-like compound library procured from ChemBridge for the identification of GlmU inhibitory leads. Structure-based as well as ligand-based approaches were applied for the development of models. The 50% inhibitory concentration of the compounds was expressed in the logarithmic scales so as to get a linear relationship during the development of QSAR models. 2.2. Structure-Based Approach for Developing Filtering Criteria. 2.2.1. Molecular Docking Using Trimer Biological Assembly of Mtb GlmU. GlmU protein is reported to exist as trimer as shown in Figure 2 with two types of the binding sites; N-terminal uridyltransferase and C-terminal acetyltransferase sites. The N-terminal uridyltransferase site is formed within a single monomer whereas the C-terminal acetyltransferase binding site is formed by the involvement of all the three C

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science study was the identification of small drug-like compounds, the PubChem compound 5388496 having a macrocycle ring and a high molecular weight of 782.52 was not considered for similarity search. ChemAxon software39 was used for similarity search using a Tanimoto coefficient of 0.5. 2.3.2. QSAR Modeling. In drug designing, quantitative structure activity relationship (QSAR) is a mathematical relationship between the physicochemical properties of compounds and their biological activities. It is one of the widely accepted approaches in ligand-based drug discovery.40−44 A robust QSAR model can be developed with a reasonable number of congeneric series of compounds having their activity values known. Hence, in the present study, the entire data set of compounds (PubChem AID 1376) was carefully clustered into different groups based on the similarity in their structural moieties. Accordingly, the entire data set was divided into three clusters (clusters A−C). Cluster-A contained 51 compounds that were mainly aromatic and heterocyclic substituted acids (Table S2); cluster-B (16 compounds) contained mainly aromatic and heterocyclic substituted sulfonyls and sulfinyls (Table S3), and cluster-C included 42 compounds that were mainly aromatic and heterocyclic substituted amides (Table S4). These three clusters were divided into training set and test set: cluster-A contained 41 training and 10 test compounds, cluster-B contained 13 training and 3 test compounds, and cluster-C contained 34 training set and 8 test compounds. For dividing the data into training and test sets, compounds were first arranged in descending order of activity and then choosing every fifth compound as the test set compound. For obtaining a linear relationship, the activities of the compounds were converted into pIC50 using the formula (−log(IC50) + 6). QSAR studies were carried out using Discovery Studio software.45 Energy minimization of the compounds was performed using default parameters of Smart Minimizer.45 Descriptors in Discovery Studio were categorized into 1D, 2D, and 3D and all these descriptors were calculated for each cluster, because only 1D or 2D or 3D descriptors cannot describe the wide range of properties of the molecules of the training set. Hence, mixed descriptors were used to explain the significance of identifying potent diverse chemical entities from the drug-like compound library. These descriptors included element counts, electronic, topological, electrotopological, spatial, structural, geometric, and thermodynamics descriptors. From the calculated descriptors, those containing zero or constant values for nearly all the compounds were removed. The selection of the variables to build the QSAR model was carried out by Genetic Function Approximation (GFA) method. In GFA method, first of all, a particular number of equations (set to 100 by default in Discovery Studio) were randomly generated. Then, pairs of “parent” equations were randomly selected from this set of 100 equations and “crossover” operations were performed at random. Only the linear equation terms were used for model building. The predictive qualities of the built models were assessed on the basis of statistical parameters: r2, q2, and prediction of the test set compounds. The external validation of the QSAR models was performed by using the methods described by Golbraikh and Tropsha46 and Roy and Roy.47 The values of the square of correlation coefficient (R2) and the square of correlation coefficient when the intercept is zero (R02) between the actual and the predicted activities of the test set compounds were calculated. These values were calculated by using Regression analysis toolpak of Microsoft Excel sheet. The predictive ability of the models was further determined by calculating the values of rm2 using the equation:47,48

rm 2 = R2(1 − | R2 − R 0 2 |)

The values of k and k′, slopes of the regression line of the predicted activity versus actual activity and vice versa, were calculated by the equations47 k=

∑ yy ̃ ii 2

∑ yĩ

and k′ =

∑ yy ̃ ii ∑ yi2̃

where yi is the actual activity and ỹi is the predicted activity of the test set compounds. 2.3.3. Pharmacophore Modeling. The concept of pharmacophore-based screening, where a set of common features among the inhibitors is used to screen the compound library, has found extensive use in drug discovery.41,44,49−56 Pharmacophore modeling was carried out using Phase module57,58 of the Schrodinger software by applying two strategies. In the first strategy (strategy-1), all the 125 PubChem inhibitors (AID 1376) were defined as active data set where it has been applied to build the hypothesis in which there is a maximum number of ligands that match at least 4 common features. In the second strategy (strategy-2), the data set of 125 inhibitors was divided into 16 actives (IC50 < 30 μM), 77 inactives (IC50 > 80 μM), and 32 moderately actives (IC50 = 30−80 μM). Phase is based on a tree based partitioning algorithm that exhaustively finds the spatial arrangements of functional groups that are common and crucial to the biological activity for a set of high-affinity ligands. The structures of the compounds taken in the data set were converted to 3D and prepared using LigPrep module26 of the Schrodinger software. Various conformers of the LigPrep prepared structures were generated by taking the default parameters of ConfGen. 59,60 Pharmacophore sites were generated on each compound and these sites belong to a builtin set of six pharmacophoric features, namely, aromatic ring (R), hydrophobic group (H), hydrogen bond acceptor (A), hydrogen bond donor (D), positively ionizable (P), and negatively ionizable (N). 2.4. In Silico Screening. For the screening of 20 000 druglike compound library (procured from ChemBridge database), the compounds were first prepared by using LigPrep26 for desalting, cleaning, and addition of the hydrogen atoms. Steroisomeric, tautomeric and ionization states at pH 7 ± 2 were also generated for these compounds. This library was selected for in silico screening because it has already been procured by us for our Institutional drug-like compound repository and hence, any compound required for validation studies would therefore be instantly available. These compounds were procured with the purpose of further modification so as to achieve novelty and specificity for any therapeutic target. Therefore, this screening is for identifying novel structural moieties/scaffolds that could be carried further for lead optimization. The prepared library was screened using filters of similarity search, QSAR models and pharmacophore hypotheses. These are widely used techniques in drug designing with each of these having its own limitations. The main aim of these techniques is to screen large libraries based on some rationale and arrive at some small set of compounds that could act as potential inhibitor candidates for validation in the wet lab. Since the enrichment of actives in the top percentages of docking hits was low, all the ligand-based filters were applied independently on the entire data set. Initially, the compound library was screened independently based on the parallel selection data fusion method. This was carried out by filtering D

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

⎤ ⎡⎛ OD of test sample ⎞ % inhibition = 100 − ⎢⎜ ⎟ × 100⎥ ⎥⎦ ⎢⎣⎝ OD of control sample ⎠

the compounds based on (i) QSAR models, (ii) pharmacophore models, and (iii) similarity search of the active compounds from the data set. Thereafter, a hybrid screening approach was applied to arrive at a reasonable number of hits for in vitro screening. In hybrid-1 approach, common hits (i) between QSAR and similarity search filters, (ii) between QSAR and pharmacophore filters, and (iii) between similarity search (most active compound CID 880974) and pharmacophore filters were selected. In hybrid-2 approach, top-ranking compounds identified through QSAR models were selected as hits to fetch (and not miss) the compounds, from the library, which are structurally dissimilar but share common descriptive properties (1D, 2D, or 3D) with the training set compounds. The unique hits from the two abovementioned hybrid approaches were considered as the final list for the in vitro screening. This combined filtering approach would ensure that no compound is missed, during the in silico screening, from similarity, pharmacophoric feature, and structure−activity relationship point of view. 2.5. In vitro Evaluation. 2.5.1. Protein Expression and Purification. Escherichia coli BL-21(DE3) carrying pQE2-Mtb GlmU plasmid used in this study was obtained as a kind gift from Dr. Vinay Nandicori, NII, India. Expression and purification of Mtb GlmU protein were performed as previously described.9 Briefly, plasmid pQE2-GlmU was transformed in E. coli BL21 (DE3) and induction carried out by using IPTG at 0.1 mM concentration. Cells were ruptured by sonication in phosphate buffer (pH 7.3). Clarified cell lysate was loaded onto a preequilibrated Ni-NTA column (Qiagen, India). The protein was eluted using a linear gradient of imidazole to 200 mM in 25 mM Tris-HCL buffer pH 8.0 containing 140 mM NaCl with 1 mM PMSF and 1 mM DTT. Collected fractions were analyzed using 12% SDS-PAGE. The fractions were dialyzed in dialysis buffer containing 25 mM Tris-HCL, pH 7.4, 15% glycerol. The purified protein was checked for purity on 12% SDS-PAGE. 2.5.2. Activity Assay for Glucosamine-1-phosphate Acetyltransferase of Mtb GlmU. The activity of Mtb GlmU was measured by using the DTNB colorimetric assay as described earlier.61 Briefly, the 50 μL reaction mixture containing 50 mM Tris-HCl (pH 7.5), 5 mM MgCl2, 0.4 mM glucosamine-1phospahte (GlcN-1-P), 0.4 mM Acetyl CoA, and 0.03 μg of purified GlmU was incubated in 96-well microtiter plate at 37 °C for 30 min. A blank reaction, containing all the components of the reaction mixture except the GlmU enzyme was used to correct the errors generated from the substrates. The reaction was terminated by adding 50 μL of stop solution containing 50 mM Tris-HCl (pH7.5), 6 M guanidine hydrochloride, then incubated for 10 min by adding 50 μL of Ellman’s reagent solution containing 0.2 mM DTNB in buffer with 50 mM Tris− HCl (pH 7.5) and 1 mM EDTA.TNB2−, the product generated from the reaction of CoA-SH and DTNB, was monitored at 412 nm. 2.5.3. Screening of Compounds for Inhibitory Activity. The identification of the inhibitors for the acetyltransferase activity of GlmU started with the screening of the compounds identified from in silico filtering. Compounds dissolved in DMSO were added to the reaction mixtures to a final concentration of 100 μM and preincubated with the enzyme prior to the start of the reaction with AcCoA and GlcN-1-P. “Z” factor was also calculated as described earlier,15 it is defined as “a measure of the quality of the high-throughput screen.” In the control samples, the compounds were replaced with the same volume of DMSO. Percent inhibition was calculated using the equation

2.5.4. IC50 Determination. The most potent hit identified in the screening of compounds (at 100 μM) was selected for dose− response assessment (IC50). Reaction mixture containing 50 mM Tris-HCl (pH 7.5), 5 mM MgCl2, GlcN-1-P, and AcCoA concentrations were fixed at 0.2 mM and the inhibitor was added at 8 different concentrations ranging from 0.78125 to 100 μM. The inhibition curve was plotted using Graphpad Prism 6.0 software (MountainView, CA). 2.6. Molecular Dynamics Simulation. Docking of the most potent hit (5810599) identified in this study was carried out using Glide27−30 at the acetyltransferase site of the Mtb GlmU trimer biological assembly (PDB structure 3ST8). The binding of this inhibitor was further analyzed in the presence of solvent using molecular dynamics (MD) simulation. The protein− inhibitor complex was solvated using explicit SPC (simple point charge) water molecules in orthorhombic box. The system was neutralized by adding 26 Na+ counterions to stabilize the net charge of the system. After addition of the solvent and counterions, the prepared system of Mtb GlmU-5810599 complex contained 154250 atoms. MD simulation was performed for a period of 10 ns with a time-step of 2 fs in NPT (constant number of atoms N, pressure P, and temperature T) ensemble at 300 K temperature and 1.01325 bar pressure. Temperature and pressure conditions were maintained by Nose−Hover chain thermostat and Martyna−Tobias−Klein barostat, respectively. MD simulation was performed using multistep protocol present in Desmond software (version 3.8).62 Recording interval of 1.2 ps (ps) was used for energy and 4.8 ps was used for trajectory. The structural stability and proteininhibitor interactions were analyzed over a period of 10 ns MD run. 2.7. Molecular Modeling Analysis for Hit/Lead Optimization. Molecular docking analysis was carried out for the compounds identified from in vitro studies to identify the orientation of the ligands within the acetyltransferase binding site. Further, the binding site analysis of the acetyltransferase binding site was carried out with respect to the hydrophobic, donor and acceptor regions using SiteMap module of Schrodinger suite.25 These sites were generated by the SiteMap based on the collection of site points at or near the surface of the protein that were separated in a solvent exposed region by small gaps where the ligand could probably orient efficiently.25 SiteMap calculates the exposure of the site to solvent, its degree of enclosure by the protein, the contacts between the site points and the protein and the degree of hydrophobicity and hydrophilicity of the site.25 The binding poses of the docked ligands were analyzed so as to find the orientation of the common structural moieties present in these compounds. These structural moieties were analyzed with respect to the vacant SiteMap regions within the pocket so as to find scope of adding hydrophobic, H-bond donor and acceptor groups on these moieties. Thereafter, the addition of the desirable groups on the identified structural moieties was done using Interactive Enumeration tool of the CombiGlide.63 Various smaller and bulkier hydrophobic, H-bond donor and acceptor groups were added on the identified moieties based on the size and the behavior of the spaces available around their docked poses. E

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

test set compounds (Table 5). An equation with five descriptors (Table 4) was derived for cluster-B that showed a good r2 value of 0.989 and q2 value of 0.971. For cluster-C, no relevant QSAR model could be derived, as indicated from low r2 value of 0.376 and low q2 value of −0.048 (Table 4), so this model was not carried forward for screening. Therefore, the two equations of cluster-A and cluster-B were selected for the screening of the compound library. The descriptors reported in the equations are explained in Table S5. Further, for the validation of the QSAR models, prediction of the test set was carried out. The values of the actual and the predicted activities for the test set compounds of both the models (as shown in Table 5) were found to be very close to each other, which reflected that the models could be used with confidence for the prediction of potent GlmU inhibitors. The graphs of the actual and the predicted activities of the training set compounds of cluster-A and -B are shown in Figures 3 and 4 respectively. Structures of all these compounds are shown in Table S1. External validation of the models was also carried out by determining the external statistical parameters of the developed QSAR models. All these parameters were found within the acceptable range as shown in Table 6.46,48 The values of R2 and R02 (test set parameters) for cluster-A model were 0.786 and 0.719, and for cluster-B model were 0.980 and 0.907, respectively. Therefore, the value of [(R2 − R02)/R2] for cluster-A model was 0.085 and cluster-B model was 0.074, which were less than 0.1 (required value).46,48 Also, the values of k and k′ for cluster-A model were 1.024 and 0.975, and for cluster-B model were 1.018 and 0.982, respectively, which were found to be in the acceptable range of 0.85 and 1.15.46,48 The value of rm2 for cluster-A model was 0.583 and for cluster-B model was 0.715, which were well within the acceptable range,36 thereby showing good predicting ability of the QSAR models. 3.2.3. Pharmacophore Models. Two four-point common pharmacophore hypotheses were derived, one from each strategy. From strategy 1, AAHR hypothesis (Figures 5 and S1 and Tables 7, S6, and S7) was built from the conformations of 39 active ligands. This model was built based on the idea that maximum number of inhibitors (out of 125 inhibitors) matches at least four pharmacophoric features. Another hypothesis, AAAD (Figures 6 and S2 and Tables 7, S8, and S9), was developed from strategy-2 that matched the conformations of 9 active ligands. Each of these hypotheses defined a set of identical features with a very similar spatial arrangement. The pharmacophore hypotheses were evaluated on the basis of “survival”, which is the weighted combination of the vector, site, volume, survival score, and a term for the number of matches in the active data set. 3.3. In Silico Screening Protocol. The entire compound library was prepared using LigPrep module26 of Schrodinger suite that generated 38 217 structures. Since docking could not be used effectively for library screening, all the ligand-based filters were applied for the same. To cover the diverse structural data set from the library, all the ligand-based methods, namely, two QSAR models, two pharmacophore models and similarity-based search were combined in such a way that there was a least possibility of missing out any compounds based on these filters. 1724 unique compounds were identified from the library by applying the two robust QSAR models (QSAR1 and QSAR2) (Figure 7 and Figure S3). The cutoff of the predicted activity (pIC50) for the selection of the hits was taken as 4.87 as this value specified the highest actual activity of the training set compound of both the clusters (cluster-A and cluster-B). Similarly, 1365

3. RESULTS AND DISCUSSION 3.1. Structure-Based Approach. 3.1.1. Analysis of the Docking Results. The ROC graph of the docking results was built between true positive rate (sensitivity) and false positive rate (1-specificity). It reflects the probability that a randomly selected active compound will rank higher than a randomly selected decoy molecule. The results of docking using Glide HTVS, SP, and XP yielded the ROC area under curve (AUC) values of 0.70, 0.68, and 0.67, respectively (Figure 2 and Table 1), Table 1. Values of the ROC and BEDROC Calculated for HTVS, SP, and XP Docking Results HTVS docking SP docking XP docking

ROC

BEDROC (α= 160.9)

0.70 0.68 0.67

0.004 0.000 0.156

which reflected that a randomly chosen active had a higher score than a randomly chosen decoy/inactive 7.0 times using HTVS, 6.8 times using SP and 6.7 times using XP docking out of 10. The ROC value of 0.7 might seem to be reasonable for screening, but it was not supported by the optimum values of other important statistical parameters that were essential for determining the enrichment factor. As during screening of the compound library top ranked active compounds are of interest, therefore the enrichment of actives in the top percentages of docking results was also analyzed (Table 2). The top percentages of hits were Table 2. Count and Percentage of Actives in Top N% of HTVS, SP, and XP Docking Results HTVS docking SP docking XP docking

% of results

1%

2%

5%

10%

20%

Active count % of actives Active count % of actives Active count % of actives

0 0.0 0 0.0 1 6.2

0 0.0 0 0.0 2 12.5

1 6.2 1 6.2 3 18.8

2 12.5 2 12.5 4 25.0

9 56.2 6 37.5 6 37.5

showing very low enrichment of actives as shown in Table 2. Besides these, the BEDROC (α = 160.9) value, which focuses on the early recognition of the actives, was analyzed that were also showing low performance in recognizing the actives in the top ranks (Table 1), as indicated by the low values of 0.004, 0.000, and 0.156 using HTVS, SP and XP docking, respectively. Therefore, docking protocol was not included in the filtering criteria for screening of the compound library. The probable reason for this might be the use of predicted trimer assembly of Mtb GlmU as its trimer structure was not available in crystallized form. Therefore, ligand-based drug design approaches were then applied to develop robust computational models for the screening of compound library. 3.2. Ligand-Based Approach. 3.2.1. Similarity Search. Similarity search was performed using 7 active query compounds taken from the data set of the PubChem Bioassay AID 1376. The search was done against the 20 000 compound library that showed 655 compounds (Table 3) were similar to the query compounds at 0.5 similarity threshold. However, no similar structure was found for the query compound 864. 3.2.2. QSAR models. For cluster-A, an equation of six descriptors (Table 4) with r2 = 0.663 and q2 = 0.519 showed a good correlation of the actual and the predicted activity of the F

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science Table 3. Hits Reported Based on Similarity Search on 20 000 Compound Library Using 7 Query Compounds

Table 4. QSAR Models and Statistical Parameters for Clusters A−C Compounds QSAR model cluster-A (QSAR 1) cluster-B (QSAR 2) cluster-Ca (QSAR 3) a

training set

test set

r2

41

10

0.663

0.519

13

3

0.989

0.971

34

8

0.376

−0.048

q2

equation activity = 5.44082 + 1.26336 * CHI_V_3_C − 0.113554 * ES_Count_aasC−0.090106 * IAC_Total +0.205093 * Kappa_1−0.311994 * Num_Rings6 + 0.000281796 * V_DIST_equ activity = 5.87889 + 0.296738 * ALogP + 0.211649 8 CHI_3_P − 0.494766 8 Num_Rings −6.20036 * Shadow_XZfrac +0.0177928 * Shadow_YZ activity = 4.60203 + 0.0161502 * Jurs_PNSA_3−0.392061 * Num_H_Donors +0.231355 * O_Count

Not used for screening of the compound library.

3.4. Enzymatic activity of Mtb GlmU. Enzymatic activity of acetyltransferase domain of Mtb GlmU was determined by using DTNB assay, which is the same as described in PubChem bioassay AID 1376. The specific activity of acetyltransferase domain was found to be 86 μM/min/mg, comparable to the already reported data.61 Km values of the substrates acetyl CoA and glucosamine-1-phosphate were determined to be 230 ± 27.9 and 108 ± 24.2 μM, respectively. 3.5. In Vitro Screening and Similarity Analysis. For the validation of bioassay, “Z” factor was calculated and found to be 0.8166 that confirmed the suitability of the bioassay for in vitro

compounds were identified by applying two pharmacophore filters (PHR1 and PHR2) (Figure 7 and Figure S3). 655 compounds were identified based on the similarity search criteria for seven most active compounds from the data set. Based on the hybrid screening approach, 371 unique hits were identified using hybrid-1 approach and top 250 hits of QSAR were selected using hybrid-2 approach. Hence, in total, 606 potential inhibitor candidates of Mtb GlmU were identified based on in silico screening for wet lab evaluation. The entire in silico filtering approach is depicted in Figure 7, and the detailed methodology adopted for the same is shown in Figure S3. G

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science Table 5. Test Set Prediction Using QSAR Models of Cluster-A and -B test set (PubChem CID)

actual activity (pIC50)

predicted activity (pIC50)

49877982a 2436052a 2163124a 6526229a 2836594a 1102906a 1810505a 9671114a 773112ac 3116472a 6393804b 1486984b 1486985b

4.47 4.29 4.17 4.10 4.04 4.01 3.91 3.86 3.76 3.68 4.31 3.97 3.74

4.37 4.28 4.15 4.05 3.84 4.10 4.03 3.46 5.00 3.28 4.31 3.94 3.54

Figure 5. AAHR hypothesis based on strategy 1 showing distances between different sites.

screening of compounds.64 On the basis of the in vitro screening, 93 compounds were found to have more than 20% inhibition of the acetyltransferase activity of GlmU at 100 μM. Out of these, 15 compounds showed more than 40% inhibition. Similarity analysis of all these compounds with respect to the training data set (PubChem compounds) was also carried out to find the structural novelty of the identified hits. The structural details and the pairwise similarity values of top 15 compounds are summarized in Table 8, and for the entire 93 compounds in Table S10. It was observed that 69 compounds out of 93 were having similarity value less than 0.5 with the training set compounds, whereas only 2 compounds had a maximum similarity of more than 0.8. Thus, this filtering methodology could help in the identification of 15 novel inhibitory leads of this target, directly from the repository. Compound ID 5810599 with 82.5% inhibition was the most active compound found in this assay as shown in Table 8. Further, the dose−response assessment of this compound was carried out and the IC50 of this compound was found to be 9.018 ± 0.04 μM (Figure 8). Besides this, similarity analysis of all the 93 compounds was also carried out with respect to the activity range of the training set data. In total, 17 compounds were found to have similarity values >0.6 with the PubChem compounds. The structure of these compounds with their IDs, similarity values, number of similar compounds and the activity range (from training set) are also included in Table S11. The 15 promising inhibitory lead compounds identified through this approach were also searched for similarity using SciFinder,65 which is an online research information service. It was found that all these compounds were novel for the acetyltransferase activity of Mtb GlmU, and none of these compounds were previously reported for any antibacterial activity. Though there has been a recent report about some substituted pyrimidine compounds showing acetyltransferase inhibitory activity of AAC(6′)-Ib.66 3.6. MD Simulation and Analysis of Protein−Inhibitor Interactions. Binding stability of the inhibitor bound at the acetyltransferase site of Mtb GlmU was validated by MD simulation. Analysis of the RMSD of the protein backbone atoms and individual inhibitor with respect to the time over a period of 10 ns of MD simulation was carried out. The RMSD value of protein backbone was stabilized below 3.2 Å. For individual inhibitor, stable trajectory of RMSD was obtained below 1 Å. The RMSD plot of the protein backbone and the inhibitor is shown in

a

Test set compound of cluster-A model. bTest set compound of the cluster-B model. cCompound 773112 was an outlier in this model and was not taken for any validation tests.

Figure 3. Graph of actual and predicted activity of training set compounds of cluster-A model.

Figure 4. Graph of actual and predicted activity of training set compounds of cluster-B model.

Table 6. External Validation Parameters of the Developed QSAR Models QSAR model

R2

R02

[(R2 − R02)/R2]

rm 2 = R2(1 − | R2 − R 0 2 |)

k

k′

cluster A cluster B

0.786 0.980

0.719 0.907

0.085 0.074

0.583 0.715

1.024 1.018

0.975 0.982

H

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science Table 7. Statistical Parameters of the Developed Pharmacophore Models strategy

hypothesis

survivala

inactivea

matchesa

reference compound (PubChem CID)

activitya (IC50)

strategy 1 strategy 2

AAHR AAAD

3.851 2.356

0.873

39 9

4293405 16195402

28.630 13.480

a “Survival” is the survival score of the actives. “Inactive” is the survival score of the inactive compounds. “Matches” represents the total number of active compounds that matched the hypothesis. “Activity” describes the activity of the reference compound.

Figure 9A. The RMSD value of the inhibitor was lower than the RMSD of the protein, which showed that the inhibitor occupied the same site throughout the simulation. This revealed the stability of the inhibitor with respect to the protein and its acetyltransferase site. Analysis of the RMSF plot of the protein showed that except for N- and C-terminal residues, all other residues were showing fewer fluctuations below 3.6 Å as shown in Figure 9B. Further, the significant interactions involved with the inhibitor binding were analyzed over a period of 10 ns MD run. The interactions that were maintained for more than 20% of the simulation time were examined. Residues A:ALA434, A:VAL449, A:ALA451, B:ARG439, C:LYS464, and C:TRP460 were involved in interactions with the inhibitor 5810599. The type and fraction of interactions are shown in Figure 9C and 9D.

Figure 6. AAAD hypothesis based on strategy 2 showing distances between different sites.

Figure 7. Flow diagram of the screening protocol followed for the filtering of 20 000 compound library. The compound repository was screened by applying parallel selection data fusion method where the filters of QSAR, pharmacophore and similarity search were used. The primary hits obtained were passed through the hybrid screening approach. In hybrid 1 approach, common hits were identified between QSAR and pharmacophore, between QSAR and similarity, and between similarity (most active compound CID 880974) and pharmacophore. In hybrid 2 approach, top 250 hits of the QSAR were taken. These two approaches resulted in 606 unique compounds that were carried forward for in vitro screening. Fifteen promising inhibitory leads were identified from in vitro screening. I

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

Table 8. In vitro Hits Having Percentage of Inhibition More than 40% at 100 μM along with Their Most Similar Compound in PubChem Bioassay (AID 1376) and Their Similarity Values

J

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science Table 8. continued

*

Similarity of these 15 hits was calculated with the 125 Mtb GlmU inhibitors reported (with activity) in the PubChem Bioassay AID 1376.

modes of these hits with the PubChem training set compounds (AID 1376) illustrated that the core of the compounds (in vitro hits as well as the active training set compounds) was bound to the same central hydrophobic space within the binding site. The functional moieties of these compounds were found to be involved in H-bonding interactions with A:SER416, A:ALA434, A:SER450, and A:ALA451 in both the data sets. Residues B:TYR398, B:ARG439, and C:TRP460 were involved in πinteractons. Besides this, the prominent residues of the binding pocket that were involved in close interaction with the inhibitors included A:ASN388, A:ALA391, A:TYR431, A:THR432, A:SER450, B:ALA422, C:LYS464, and C:ARG465. The interaction figures of 5810599 (the most potent identified compound) and 880974 (from the training set) are shown in Figure 10. Further, the orientation of all the 93 hits was analyzed within the binding pocket of acetyltransferase domain to design the lead optimization strategies. On the basis of the structural superimposition, eight common structural moieties were identified, which could act as templates for further modification and optimization of their acetyltransferase inhibiting activity. These structural moieties belonged to six different classes of compounds, namely, pyrimidine, thiazolidine, sulfonyl, dikeone, phenol, and xanthenes. It was observed that these common structural moieties occupied the central hydrophobic core region within the binding pocket as shown in Figure 11. Further, the binding site analysis of the acetyltransferase binding site could help in the identification of the vacant hydrophobic, donor and

Figure 8. Dose−response curve for the inhibition of Mtb GlmU by compound ID-5810599. The IC50 value of the compound was 9.018 ± 0.04 μM. Shown are the mean and deviation for triplicate analysis.

Residues A:ALA434, A:ALA451, and C:LYS464 were forming H-bond interactions with the inhibitors. B:ARG439 was identified as a crucial residue forming water bridge, ionic, hydrophobic and H−bond interactions. C:TRP460 was also recognized as an important residue that was forming hydrophobic contacts with the inhibitor for over 20% of the simulation time. 3.7. Scope of Modification around the Hits. To understand the mode of binding and carry out lead optimization studies, all the 93 hits were docked on to the acetyltransferase binding site of the Mtb GlmU. The comparison of the binding K

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

Figure 9. Analysis of the structural stability and protein-inhibitor interaction. (A) RMSD plot of protein backbone and inhibitor heavy atoms. (B) RMSF plot of protein. (C) Type and fraction of protein-inhibitor interactions over a period of 10 ns MD run. (D) Residues interacting for more than 20% of simulation time.

Figure 10. Ligand interaction diagrams showing residues interacting at the acetyltransferase binding site. (A) Residues interacting with the most active compound 5810599 reported in the present study. (B) Residues interacting with the training set PubChem compound 880974.

acceptor regions within the pocket. Based on the in vitro activity data, docked poses and site map analysis, scope of modifications for improving the binding affinity were identified on these

structural moieties as shown in Figure 12. Structural moiety A was a pyrimidine substituted with benzylidene. The orientation of the moiety A with respect to the empty spaces within the L

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

groups, whereas a restricted scope of modification with acceptor group was possible on the right side of the moiety. The various possible modifications around these structural moieties (A−H) are demonstrated in Figure 12. The binding poses of these structural moieties with respect to the binding surface analysis, for further modifications, are also shown in Table S12. Further, based on the binding surface analysis, various fragments representing the hydrophobic, H-bond acceptor and donor groups were added on the identified positions of the structural moieties. The fragments were divided into bulkier and smaller groups, and added as required to occupy the vacant spaces of the binding site. For example, in the structural moiety A, at first position of the pyrimidine ring, that is, at N atom, a bulkier hydrophobic group was required, and therefore, either a phenyl group or a propyl group was added at this position. Further, various fragments containing H-bond acceptor and donor groups were attached next to the added hydrophobic group. Similarly, small H-bond acceptor groups such as carbonyl, carboxylic, aldehyde and small H-bond donor groups, such as NH2, CH2−NH2, were attached to the phenyl ring of moiety A as there was a limited unoccupied space. On the basis of these additions, various compounds containing moiety A were built. Similar modifications were also performed on rest of the structural moieties (B−H). All these modifications led to the building of a combinatorial library of 13 116 structures. The activities of the designed compounds were predicted using the two QSAR models built in the current study, that is, QSAR 1 (cluster-A) and QSAR 2 (cluster-B). The compounds having ≥6.0 pIC50 value were selected as the QSAR hits. Using this criterion, 331 compounds were identified as the QSAR hits. Further, the identified QSAR hits were passed through the two pharmacophore models (built using strategy-1 and strategy-2). The compounds possessing the pharmacophore Fitness score ≥1.5 were selected as the final hits, which resulted in the identification of 18 compounds. On the basis of the developed in silico models, the IC50 value for these compounds was predicted to be close to 1 μM. These predicted inhibitors belonged to two

Figure 11. Superimposition of the structural moieties within the acetyltransferase binding pocket of Mtb GlmU. Yellow colored surface represents the hydrophobic cavity of the binding site.

binding site suggested that at first position of the pyrimidine ring, that is, at N atom, there was a vacant space to accommodate a phenyl ring or a hydrophobic chain of up to three carbon atoms. Further, the addition of H-bond donors and acceptors were preferable on the added hydrophobic groups for increasing the affinity of the resulting compound. There existed an unoccupied space around the phenyl ring of moiety A, where small acceptor groups like carbonyl and small donor groups like NH2 could be added for increasing the affinity. Similar types of vacant spaces were available around moieties B-G, except for structural moiety H where there existed a large vacant space toward the left of the moiety H for the addition of hydrophobic, donor and acceptor

Figure 12. Scope of modification around the identified structural moieties (A−H). Various sites for modification were recognized based on binding site analysis and the mode of binding of these moieties at the acetyltransferase binding site using Glide XP docking. Dotted circles show scope of adding various groups at these positions. Green colored circles indicate scope of adding both H-bond acceptor as well as donor groups, red colored circles show scope of adding H-bond acceptor groups, blue colored circles represent scope of adding H-bond donor groups, and orange colored circles indicate scope of adding hydrophobic groups. M

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

Figure 13. Designed compounds predicted to have inhibitory activity for the acetyltransferase site of Mtb GlmU. Structure 1 contains moiety D and structures 2−6 belong to the moiety F.

Application of these filters and thereafter the in vitro screening helped in the identification of 15 promising inhibitory leads of Mtb GlmU acetyltransferase activity from the compound repository. IC50 value of the most potent compound was determined and found to be 9.018 ± 0.04 μM. Moreover, close to 100 hits were identified that showed more than 20% inhibition of the acetyltransferase activity. Similarity analysis of these compounds with the training set data revealed the structural novelty of these hits. SciFinder search of the 15 promising inhibitory leads depicts the novelty of these compounds as Mtb GlmU inhibitors. Molecular docking analysis of all these hits were carried out on the available 3D crystal structure of Mtb GlmU, for elucidating the binding pose of these compounds within the acetyltransferase binding pocket. It was observed that the hits identified as well as the active compounds of the training data set bind at the similar location within the binding site of acetyltransferase domain. Molecular dynamics (MD) simulation of the most potent inhibitory lead in complex with protein revealed important interacting residues and structural stability of the protein−inhibitor complex. To elucidate the scope of modifications of the identified hits for improving the binding affinity, thorough analysis of the binding pocket was carried out with respect to the vacant hydrophobic, donor and acceptor regions within the binding pocket. Common structural moieties of the hits were also identified, which could be taken forward for further structural modifications, combined with the binding surface analysis, for designing a novel and potent Mtb GlmU inhibitor. Structural modification of these moieties and their combinatorial design resulted in the identification of 18 compounds predicted to have an IC50 value close to 1 μM. To reduce the number of false-negatives in a virtual screening method, the compounds similar to the most active identified hit were also evaluated in vitro. This resulted in the identification of two more compounds having an IC50 value of ∼27 μM from the 20 000 compounds library. Overall, this robust computational model, with in vitro validation, would be very useful for the identification of novel potent inhibitory leads that can be taken forward for further optimization and development.

structural moieties D and F (Figure 12). The top 6 representative structures of the designed compounds predicted to have inhibitory activity for the acetyltransferase site of Mtb GlmU are depicted in Figure 13. Rest of the 12 predicted inhibitors are shown in Figure S4. 3.7.1. Similarity Search of the Best Identified Inhibitory Lead. It is inevitable to miss out some active compounds (falsenegatives), while screening any compound library using computational methods. One way to reduce the number of false-negatives is to revisit the compound library and identify the compounds similar to the most active in vitro hit. To fish out the possible actives that might have been missed during the in silico screening, the compounds similar to the most active identified inhibitory lead (5810599) were retrieved from the 20 000 compound library. The similarity search in this case was performed using a high Tanimoto coefficient of 0.7. Using this approach, 23 in silico hits were identified that were not screened in vitro earlier. These 23 hits were also evaluated using DTNB bioassay, which helped in the identification of 8 additional inhibitory leads with more than 40% inhibition as shown in Table 9. Out of these, 2 compounds, namely, 5680705 and 5673625 exhibited the IC50 values of 27.65 ± 0.021 and 28.49 ± 0.027 μM, respectively. It was found that both these compounds were highly similar to the query compound 5810599, except that a chlorobenzene group was attached to the thioxopyrimidinedione instead of anisole.

4. CONCLUSIONS As GlmU is an interesting target for developing novel M. tuberculosis drug candidates4 and very few inhibitors are known for this target, a robust computational strategy was developed, with in vitro validation for the identification of novel potent inhibitory leads of Mtb GlmU from a drug-like compound repository of 20 000 compounds. The structure-based approach (docking studies) showed low enrichment of actives in the top ranks when mixed with decoys and confirmed inactives, and hence, could not be used as an effective filtering criterion in this case. Therefore, we combined the ligand-based strategies comprising of QSAR studies, pharmacophore modeling and similarity search for developing a robust filtering methodology. The same methodology could be used for the filtering of other compound repositories available in the public domain. N

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science Table 9. In Vitro Hits Identified by the Similarity Search of the Inhibitory Lead 5810599



ASSOCIATED CONTENT

S Supporting Information *



The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acscombsci.5b00019. Compounds and their structures, dataset of Cluster-A, -B, and -C for QSAR modeling, descriptors used in the QSAR equations, distances and angles between different sites of AAHR hypothesis, distances and angles between different sites of AAAD hypothesis, lists of in vitro hits, analysis of the scope of modification, AAHR hypothesis based on strategy-1 displayed on compound 4293405, AAAD hypothesis based on strategy-2 displayed on compound 5811786, in silico methodology followed for the screening of 20 000 ChemBridge compound library, and designed

compounds predicted to have inhibitory activity for the acetyltransferase site of Mtb GlmU (PDF)

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. Phone: 0191-2569021 EPAB Ext.: 269. Author Contributions

R.M. and A.N. conceived and designed the computational experiments. R.M. carried out the computational experiments. A.N. analyzed the computational data. P.M. carried out the binding site analysis. C.R. and I.A.K. conceived, designed, and carried out the in vitro experiments. R.M. and C.R. drafted the O

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science

GlmU enzyme inhibitors against catheter-associated uropathogens. Antimicrob. Agents Chemother. 2006, 50, 1835−1840. (13) Li, Y.; Zhou, Y.; Ma, Y.; Li, X. Design and synthesis of novel cell wall inhibitors of Mycobacterium tuberculosis GlmM and GlmU. Carbohydr. Res. 2011, 346 (13), 1714−1720. (14) Tran, A. T.; Wen, D.; West, N. P.; Baker, E. N.; Britton, W. J.; Payne, R. J. Inhibition studies on Mycobacterium tuberculosis Nacetylglucosamine-1-phosphate uridyltransferase (GlmU). Org. Biomol. Chem. 2013, 11 (46), 8113−26. (15) Pereira, M. P.; Blanchard, J. E.; Murphy, C.; Roderick, S. L.; Brown, E. D. High-Throughput Screening Identifies Novel Inhibitors of the Acetyltransferase Activity of Escherichia coli GlmU. Antimicrob. Agents Chemother. 2009, 53 (6), 2306−2311. (16) Buurman, E. T.; Andrews, B.; Gao, N.; Hu, J.; Keating, T. A.; Lahiri, S.; Otterbein, L. R.; Patten, A. D.; Stokes, S. S.; Shapiro, A. B. In Vitro Validation of Acetyltransferase Activity of GlmU as an Antibacterial Target in Haemophilus inf luenzae. J. Biol. Chem. 2011, 286 (47), 40734−40742. (17) Green, O. M.; McKenzie, A. R.; Shapiro, A. B.; Otterbein, L.; Ni, H.; Patten, A.; Stokes, S.; Albert, R.; Kawatkar, S.; Breed, J. Inhibitors of acetyltransferase domain of N-acetylglucosamine-1-phosphate-uridyltransferase/glucosamine-1-phosphate-acetyltransferase (GlmU). Part 1: Hit to lead evaluation of a novel arylsulfonamide series. Bioorg. Med. Chem. Lett. 2012, 22 (4), 1510−1519. (18) Stokes, S. S.; Albert, R.; Buurman, E. T.; Andrews, B.; Shapiro, A. B.; Green, O. M.; McKenzie, A. R.; Otterbein, L. R. Inhibitors of the acetyltransferase domain of N-acetylglucosamine-1-phosphate-uridylyltransferase/glucosamine-1-phosphate-acetyltransferase (GlmU). Part 2: Optimization of physical properties leading to antibacterial aryl sulfonamides. Bioorg. Med. Chem. Lett. 2012, 22 (23), 7019−7023. (19) Vertesy, L.; Kurz, M.; Markus-Erb, A.; Toti, L. 2-Phenylbenzofuran derivatives, a process for preparing them, and their use. U.S. Patent 7,148,254, December 12, 2006. (20) Pompeo, F.; van Heijenoort, J.; Mengin-Lecreulx, D. Probing the Role of Cysteine Residues in Glucosamine-1-Phosphate Acetyltransferase Activity of the Bifunctional GlmU Protein from Escherichia coli: SiteDirected Mutagenesis and Characterization of the Mutant Enzymes. J. Bacteriol. 1998, 180 (18), 4799−4803. (21) Wang, Y.; Xiao, J.; Suzek, T. O.; Zhang, J.; Wang, J.; Zhou, Z.; Han, L.; Karapetyan, K.; Dracheva, S.; Shoemaker, B. A.; Bolton, E.; Gindulyte, A.; Bryant, S. H. PubChem’s BioAssay Database. Nucleic Acids Res. 2012, 40, D400−12. (22) Mehra, R.; Sharma, R.; Khan, I. A.; Nargotra, A. Identification and optimization of Escherichia coli GlmU inhibitors: An in silico approach with validation thereof. Eur. J. Med. Chem. 2015, 92, 78−90. (23) Singla, D.; Anurag, M.; Dash, D.; Raghava, G. P. S. A web server for predicting inhibitors against bacterial target GlmU protein. BMC Pharmacol. 2011, 11, 5. (24) Maestro, version 9.3; Schrödinger, LLC: New York, 2012. (25) SiteMap, version 2.6; Schrödinger, LLC: New York, 2012. (26) LigPrep, version 2.5; Schrödinger, LLC: New York, 2012. (27) Glide, version 5.8; Schrödinger, LLC: New York, 2012. (28) Friesner, R. A.; Banks, J. L.; Murphy, R. B.; Halgren, T. A.; Klicic, J. J.; Mainz, D. T.; Repasky, M. P.; Knoll, E. H.; Shaw, D. E.; Shelley, M.; Perry, J. K.; Francis, P.; Shenkin, P. S. Glide: A new approach for rapid, accurate docking and scoring. 1. Method and assessment of docking accuracy. J. Med. Chem. 2004, 47, 1739−1749. (29) Halgren, T. A.; Murphy, R. B.; Friesner, R. A.; Beard, H. S.; Frye, L. L.; Pollard, W. T.; Banks, J. L. Glide: A new approach for rapid, accurate docking and scoring, 2. Enrichment factors in database screening. J. Med. Chem. 2004, 47, 1750−1759. (30) Friesner, R. A.; Murphy, R. B.; Repasky, M. P.; Frye, L. L.; Greenwood, J. R.; Halgren, T. A.; Sanschagrin, P. C.; Mainz, D. T. Extra Precision Glide: Docking and scoring incorporating a model of hydrophobic enclosure for protein-ligand complexes. J. Med. Chem. 2006, 49, 6177−6196. (31) Truchon, J. F.; Bayly, C. I. Evaluating virtual screening methods: good and bad metrics for the ‘early recognition’ problem. J. Chem. Inf. Model. 2007, 47, 488−508.

manuscript. R.A.V., I.A.K., and A.N. supervised the work, implemented resources, and participated in manuscript revisions. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS C.R. (SRF) acknowledges the financial support from University Grant Commission (UGC). The authors acknowledge the financial support from the network project BSC-0205, GAP0141, and GAP-2127.



ABBREVIATIONS Mtb:Mycobacterium tuberculosis QSAR:quantitative structure−activity relationship PDB:Protein Data Bank IC50:half maximal inhibitory concentration ROC:receiver operating characteristic AUC:area under curve DMSO:dimethyl sulfoxide CoA:coenzyme-A SPC:simple point charge



REFERENCES

(1) Van den Boogaard, J.; Kibiki, G. S.; Kisanga, E. R.; Boeree, M. J.; Aarnoutse, R. E. New drugs against tuberculosis: problems, progress, and evaluation of agents in clinical development. Antimicrob. Agents Chemother. 2009, 53, 849−862. (2) Barreteau, H.; Kovac, A.; Boniface, A.; Sova, M.; Gobec, S.; Blanot, D. Cytoplasmic steps of peptidoglycan biosynthesis. FEMS Microbiol. Rev. 2008, 32, 168−207. (3) Zhang, W.; Jones, V. C.; Scherman, M. S.; Mahapatra, S.; Crick, D.; Bhamidi, S.; Xin, Y.; McNeil, M. R.; Ma, Y. Expression, essentiality, and a microtiter plate assay for mycobacterial GlmU, the bifunctional glucosamine-1-phosphate acetyltransferase and N-acetylglucosamine1-phosphate uridyltransferase. Int. J. Biochem. Cell Biol. 2008, 40, 2560− 2571. (4) Sassetti, C. M.; Boyd, D. H.; Rubin, E. J. Genes required for mycobacterial growth defined by high density mutagenesis. Mol. Microbiol. 2003, 48, 77−84. (5) Olsen, L. R.; Roderick, S. L. Structure of the Escherichia coli GlmU pyrophosphorylase and acetyltransferase active sites. Biochemistry 2001, 40, 1913−1921. (6) Olsen, L. R.; Vetting, M. W.; Roderick, S. L. Structure of the E. coli bifunctional GlmU acetyltransferase active site with substrates and products. Protein Sci. 2007, 16, 1230−1235. (7) Zhang, Z.; Bulloch, E. M. M.; Bunker, R. D.; Baker, E. N.; Squire, C. J. Structure and function of GlmU from Mycobacterium tuberculosis. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2009, D65, 275−283. (8) Parikh, A.; Verma, S. K.; Khan, S.; Prakash, B.; Nandicoori, V. K. PknB-mediated phosphorylation of a novel substrate, N-acetylglucosamine-1-phosphate uridyltransferase, modulates its acetyltransferase activity. J. Mol. Biol. 2009, 386 (2), 451−464. (9) Verma, S. K.; Jaiswal, M.; Kumar, N.; Parikh, A.; Nandicoori, V. K.; Prakash, B. Structure of N-acetylglucosamine-1-phosphate uridyltransferase (GlmU) from Mycobacterium tuberculosis in a cubic space group. Acta Crystallogr., Sect. F: Struct. Biol. Cryst. Commun. 2009, 65, 435−439. (10) Jagtap, P. K.; Verma, S. K.; Vithani, N.; Bais, V. S.; Prakash, B. Crystal structures identify an atypical two-metal-ion mechanism for uridyltransfer in GlmU: its significance to sugar nucleotidyl transferases. J. Mol. Biol. 2013, 425 (10), 1745−59. (11) National Center for Biotechnology Information. PubChem BioAssay Database: AID = 1376, Source = Southern Research Molecular Libraries Screening Center (SRMLSC). http://pubchem.ncbi.nlm.nih. gov/assay/assay.cgi?aid=1376. (12) Burton, E.; Gawande, P. V.; Yakandawala, N.; Lovetri, K.; Zhanel, G. G.; Romeo, T.; Friesen, A. D.; Madhyastha, S. Antibiofilm activity of P

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX

Research Article

ACS Combinatorial Science (32) Johnson, M. A.; Maggiora, G. M. Concepts and Applications of Molecular Similarity; John Willey & Sons: New York, 1990. (33) Martin, Y. C.; Kofron, J. L.; Traphagen, L. M. Do structurally similar molecules have similar biological activity? J. Med. Chem. 2002, 45, 4350−4358. (34) Franke, L.; Schwarz, O.; Kuhrt, L. M.; Hoernig, C.; Fischer, L.; George, S.; Tanrikulu, Y.; Schneider, P.; Werz, O.; Steinhilber, D.; Schneider, G. Identification of natural-product-derived inhibitors of 5lipoxygenase activity by ligand-based virtual screening. J. Med. Chem. 2007, 50, 2640−2646. (35) Betzi, S.; Restouin, A.; Opi, S.; Arold, S. T.; Parrot, I.; Guerlesquin, F.; Morelli, X.; Collette, Y. Protein−protein interaction inhibition (2P2I) combining high throughput and virtual screening: Application to the HIV-1 Nef protein. Proc. Natl. Acad. Sci. U. S. A. 2007, 104, 19256− 19261. (36) Willett, P. Similarity-based virtual screening using 2D fingerprints. Drug Discovery Today 2006, 11, 1046−1053. (37) Amaro, R. E.; Schnaufer, A.; Interthal, H.; Hol, W.; Stuart, K. D.; McCammon, J. A. Discovery of drug-like inhibitors of an essential RNAediting ligase in Trypanosoma brucei. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 17278−17283. (38) Langdon, S. R.; Westwood, I. M.; van Montfort, R. L. M.; Brown, N.; Blagg, J. Scaffold-focused virtual screening: Prospective application to the discovery of TTK inhibitors. J. Chem. Inf. Model. 2013, 53, 1100− 1112. (39) Instant JChem, version 5.9.4; ChemAxon, 2012; http://www. chemaxon.com. (40) Yang, S. C.; Chang, S. S.; Chen, H. Y.; Chen, C. Y. C. Identification of potent EGFR inhibitors from TCM database@Taiwan. PLoS Comput. Biol. 2011, 7, e1002189. (41) Valasani, K. R.; Vangavaragu, J. R.; Day, V. W.; Yan, S. S. Structure-based design, synthesis, pharmacophore modeling, virtual screening, and molecular docking studies for identification of novel cyclophilin D inhibitors. J. Chem. Inf. Model. 2014, 54, 902−912. (42) Baurin, N.; Mozziconacci, J. C.; Arnoult, E.; Chavatte, P.; Marot, C.; Allory, L. M. 2D QSAR consensus prediction for high-throughput virtual screening: An application to COX-2 inhibition modeling and screening of the NCI database. J. Chem. Inf. Model. 2004, 44, 276−28. (43) Hoffman, B. T.; Kopajtic, T.; Katz, J. L.; Newman, A. H. 2D QSAR modeling and preliminary database searching for dopamine transporter inhibitors using genetic algorithm variable selection of Molconn Z descriptors. J. Med. Chem. 2000, 43, 4151−4159. (44) Evers, A.; Hessler, G.; Matter, H.; Klabunde, T. Virtual screening of biogenic amine-binding G-protein coupled receptors: Comparative evaluation of protein and ligand-based virtual screening protocols. J. Med. Chem. 2005, 48, 5448−5465. (45) Discovery Studio Modeling Environment, release 2.1; Accelrys Software Inc.; http://www.accelrys.com. (46) Golbraikh, A.; Tropsha, A. Beware of q2. J. Mol. Graphics Modell. 2002, 20, 269−276. (47) Roy, P.; Roy, K. On some aspects of variable selection for partial least squares regression models. QSAR Comb. Sci. 2008, 27, 302−313. (48) Mehra, R.; Nargotra, A.; Shah, B. A.; Taneja, S. C.; Vishwakarma, R. A.; Koul, S. Pro-apoptotic properties of parthenin analogs: a quantitative structure−activity relationship study. Med. Chem. Res. 2013, 22, 2303−2311. (49) Poongavanam, V.; Kongsted, J. Virtual screening models for prediction of HIV-1 RT associated RNase H inhibition. PLoS One 2013, 8, e73478. (50) Krishna, S.; Singh, D. K.; Meena, S.; Datta, D.; Siddiqi, M. I.; Banerjee, D. Pharmacophore-based screening and identification of novel human ligase I inhibitors with potential anticancer activity. J. Chem. Inf. Model. 2014, 54, 781−792. (51) Vuorinen, A.; Engeli, R.; Meyer, A.; Bachmann, F.; Griesser, U. J.; Schuster, D.; Odermatt, A. Ligand-based pharmacophore modeling and virtual screening for the discovery of novel 17β- hydroxysteroid dehydrogenase 2 inhibitors. J. Med. Chem. 2014, 57, 5995. (52) Chiang, Y. K.; Kuo, C. C.; Wu, Y. S.; Chen, C. T.; Coumar, M. S.; Wu, J. S.; Hsieh, H. P.; Chang, C. Y.; Jseng, H. Y.; Wu, M. H.; Leou, J. S.;

Song, J. S.; Chang, J. Y.; Lyu, P. C.; Chao, Y. S.; Wu, S. Y. Generation of ligand-based pharmacophore model and virtual screening for identification of novel tubulin inhibitors with potent anticancer activity. J. Med. Chem. 2009, 52, 4221−4233. (53) Shukla, S.; Kouanda, A.; Silverton, L.; Talele, T. T.; Ambudkar, S. V. Pharmacophore modeling of nilotinib as an inhibitor of ATP-binding cassette drug transporters and BCR-ABL kinase using a threedimensional quantitative structure−activity relationship approach. Mol. Pharmaceutics 2014, 11, 2313−2322. (54) Ritschel, T.; Hermans, S. M. A.; Schreurs, M.; Heuvel, J. J. M. W.; Koenderink, J. B.; Greupink, R.; Russel, F. G. M. In Silico identification and in vitro validation of potential cholestatic compounds through 3D ligand-based pharmacophore modeling of BSEP inhibitors. Chem. Res. Toxicol. 2014, 27, 873−881. (55) Yang, S. Y. Pharmacophore modeling and application in drug discovery: challenges and recent advances. Drug Discovery Today 2010, 15, 444−450. (56) Leach, A. R.; Gillet, V. J.; Lewis, R. A.; Taylor, R. Threedimensional pharmacophore methods in drug discovery. J. Med. Chem. 2010, 53, 539−558. (57) Phase, version 3.4; Schrödinger, LLC: New York, 2012. (58) Dixon, S. L.; Smondyrev, A. M.; Knoll, E. H.; Rao, S. N.; Shaw, D. E.; Friesner, R. A. PHASE: A new engine for pharmacophore perception, 3D QSAR model development, and 3D database screening. 1. Methodology and preliminary results. J. Comput.-Aided Mol. Des. 2006, 20, 647−671. (59) ConfGen, version 2.3; Schrödinger, LLC: New York; 2012. (60) Watts, K. S.; Dalal, P.; Murphy, R. B.; Sherman, W.; Friesner, R. A.; Shelley, J. C. ConfGen: A conformational search method for efficient generation of bioactive conformers. J. Chem. Inf. Model. 2010, 50, 534− 546. (61) Zhou, Y.; Xin, Y.; Sha, S.; Ma, Y. Kinetic properties of Mycobacterium tuberculosis bifunctional GlmU. Arch. Microbiol. 2011, 193, 751−757. (62) Schrödinger Release 2014-2: Desmond Molecular Dynamics System, version 3.8; D. E. Shaw Research: New York, 2014. MaestroDesmond Interoperability Tools, version 3.8; Schrödinger: New York, 2014. (63) Zhou, J. Z. Structure-directed combinatorial library design. Curr. Opin. Chem. Biol. 2008, 12, 379−385. (64) Zhang, J. H.; Chung, T. D.; Oldenburg, K. R. A simple statistical parameter for use in evaluation and validation of high throughput screening assays. J. Biomol. Screening 1999, 4, 67−73. (65) SciFinder, https://scifinder.cas.org/scifinder. (66) Lin, D. L.; Tran, T.; Adams, C.; Alam, J. Y.; Herron, S. R.; Tolmasky, M. E. Inhibitors of the aminoglycoside 60-N-acetyltransferase type Ib [AAC(60)-Ib] identified by in silico molecular docking. Bioorg. Med. Chem. Lett. 2013, 23, 5694−5698.

Q

DOI: 10.1021/acscombsci.5b00019 ACS Comb. Sci. XXXX, XXX, XXX−XXX