Fusion of Structure and Ligand Based Methods for Identification of

Article pubs.acs.org/jcim

Fusion of Structure and Ligand Based Methods for Identification of Novel CDK2 Inhibitors Priya Mahajan,†,∥ Gousia Chashoo,‡ Monika Gupta,† Amit Kumar,† Parvinder Pal Singh,§,∥ and Amit Nargotra*,†,∥ †

Discovery Informatics, ‡Cancer Pharmacology, §Medicinal Chemistry, and ∥Academy of Scientific and Innovative Research, CSIR-Indian Institute of Integrative Medicine, Canal Road, Jammu 180001, India

Downloaded via IOWA STATE UNIV on January 21, 2019 at 20:38:18 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Cyclin dependent kinases play a central role in cell cycle regulation which makes them a promising target with multifarious therapeutic potential. CDK2 regulates various events of the eukaryotic cell division cycle, and the pharmacological evidence indicates that overexpression of CDK2 causes abnormal cell-cycle regulation, which is directly associated with hyperproliferation of cancer cells. Therefore, CDK2 is regarded as a potential target molecule for anticancer medication. Thus, to decline CDK2 activity by potential lead compounds has proved to be an effective treatment for cancer. The availability of a large number of X-ray crystal structures and known inhibitors of CDK2 provides a gateway to perform efficient computational studies on this target. With the aim to identify new chemical entities from commercial libraries, with increased inhibitory potency for CDK2, ligand and structure based computational drug designing approaches were applied. A druglike library of 50,000 compounds from ChemDiv and ChemBridge databases was screened against CDK2, and 110 compounds were identified using the parallel application of these models. On in vitro evaluation of 40 compounds, seven compounds were found to have more than 50% inhibition at 10 μM. MD studies of the hits revealed the stability of these inhibitors and pivotal role of Glu81 and Leu83 for binding with CDK2. The overall study resulted in the identification of four new chemical entities possessing CDK2 inhibitory activity.

1. INTRODUCTION Cyclin-dependent kinase (CDKs) belongs to a class of serine/ threonine protein kinases that act as a key regulatory element in cell cycle progression and transcription. CDKs exist in various isoforms, and deregulation of CDKs has been known for its largespectrum therapeutic potential.1 CDK2 is an extensively studied target for its antitumor activity due to its significant contributions in signal transduction pathways involved in cell cycle regulation. It plays an important role in the cell division cycle by regulating multiple events i.e. centrosome duplication, DNA synthesis, G1-S transition, and modulation of G2 progression. Anomalous expression of CDK2 at the restriction point may lead to abnormal cell-cycle regulation, which has been identified in many human cancers.2 Development of efficient drugs for CDK2 to combat the cancer is the need of the hour. CDK2 has been intensively investigated as a therapeutic target for cancer for many decades, and many inhibitors have surfaced out, which belong to the diverse scaffolds viz. indazoles, thiazoles, isothiazoles, acylaminopyrazoles, cabolines, thiazole, pyrimidines, etc.3 To date, many CDK2 inhibitors, i.e. SNS-032 (BMS-387032), Flavopiridol (Alvocidib) HCl, Milciclib (PHA-848125), PHA848125, BAY 10-00394, AT7519, SCH 727965, and Roscovitine (Seliciclib, CYC274), are in clinical trials for anticancer studies.3 Although several CDK2 inhibitors have © 2017 American Chemical Society

been investigated clinically for their potential as anticancer agents, none have been approved for clinical use due to isoform selectivity, solubility, toxicity, etc. issues. Still, many pharmaceutical companies such as Novartis, AstraZeneca, etc. are putting a great deal of effort into designing new and potential inhibitors for CDK2. Day by day the number of crystal structures of CDK2 and its inhibitors is increasing. The availability of such data makes CDK2 a very interesting target for computational studies. In the present study, diverse CDK2 inhibitors and the X-ray crystal structure of protein in complex with inhibitors were explored for designing structure and ligand based in silico models. A ligand based substructure similarity search and pharmacophore based 3D QSAR models were built from known CDK2 inhibitors which identifies compounds that share similar structural properties to known actives. Structure based e-Pharmacophore models and docking studies were performed guided by the crystal structure information on CDK2 in complex with bound inhibitors, which identify compounds that complement with the active site of the target protein. The developed in silico models were applied in parallel to screen a druglike library of Received: May 22, 2017 Published: July 19, 2017 1957

DOI: 10.1021/acs.jcim.7b00293 J. Chem. Inf. Model. 2017, 57, 1957−1969

Article

Journal of Chemical Information and Modeling

biological activity, substructure and similarity search of the reported CDK2 inhibitors was carried out. A data set of 77 CDK2 inhibitors consisting of clinical, preclinical, biological testing compounds was collected from the literature (Table S2). These 77 inhibitors, along with the 27 scaffolds identified from cocrystallized inhibitors, were considered for substructure and similarity search using ChemAxon software with a similarity threshold of 0.5 i.e. identified compounds must have 50% similarity to the query compounds.50 2.2.2. Pharmacophore Based 3D-QSAR Modeling. A pharmacophore model identifies the minimum necessary structure features of inhibitors, that a compound must possess for inhibition of an up regulated enzyme involved in a disease. These structural features are explored to identify potential candidates from large compound databases. Herein, four congeneric series of CDK2 inhibitors were considered to build pharmacophore models for the identification of potent inhibitors from the Institutional compound repository. The reason for building four predictive models, considering SetA-D series, was that the compounds attained via screening must have a maximum number of common pharmacophore features with respect to each series. The number of compounds considered for building pharmacophore models for each series is given in Table 1. Before building the pharmacophore model, inhibitory activity IC50 values of inhibitors were converted into their logarithmic value (pIC50) so as to get the linear corelation while generating QSAR models. The structures of inhibitors were sketched and minimized at the OPLS_2005 force field. Prepared ligands were then used for generating a common pharmacophore hypothesis using the Phase module of the Schrodinger software. To build pharmacophore models the data set of compounds from each series was divided into active, moderately active, and less active compounds on the basis of activity (Tables S3−S6). The chemical structure of each compound conveys a particular set of pharmacophoric features i.e., H-bond acceptor (A), H-bond donor (D), hydrophobic group (H), negatively charged group (N), positively charged group (P), and aromatic ring (R) calculated using Phase.51 Based on these pharmacophoric features, the pharmacophore models were developed. Each model/hypothesis conveys a particular 3D conformation for a set of active ligands in which the compound was going to bind with the receptor. In order to validate the developed pharmacophore hypothesis for the selected series, 3D-QSAR studies were carried out. To build QSAR models, the initial data set was divided into

50,000 compounds procured from the ChemDiv and ChemBridge library. Thereafter, based on the interaction of the in silico hits with the key residues, energy parameter calculation, and structural diversity, a few compounds were selected and subjected for biological activity. The results emphasized on the successful application of virtual screening identify new CDK2 inhibitor scaffolds that can be taken forward for carrying out further medicinal chemistry and optimization studies for designing selective CDK2 inhibitors.

2. MATERIALS AND METHODS 2.1. Data Set. There are more than 300 X-ray crystal structures of CDK2 reported in PDB,4and out of these, nearly 200 are bound with diverse inhibitors. These cocrystallized inhibitors were categorized into 27 diverse scaffolds, and the respective PDB entry formed the data set for structure based computational studies (structures of these scaffolds with their respective PDB Ids are given in Table S1). Twenty-six other CDK2 inhibitors in clinical and preclinical studies were collected from Drug Bank5 and Scifinder databases.6 In addition, CDK2 inhibitors with IC50 less than 100 μM were collected from the ChEMBL database,7 out of which the 51 most potent and diverse inhibitors were identified (Table S2). This data set of 77 known CDK2 inhibitors (26 from Drug Bank and Scifinder databases, 51 from the ChEMBL database)8−39 along with 27 diverse scaffolds identified from cocrystallized inhibitors was considered for ligand based computational studies. In addition, four congeneric series, i.e. indenopyrazole derivative (Set-A series), pyrazolo[4,3-h] qinazoline-3-carboxamides (Set-B series), pyrazoles with imidazole pyrimidine amides derivatives (Set-C series), and 3-aminopyrazoles (Set-D series) reported as CDK2 inhibitors,35−49 were considered to build pharmacophore based 3D-QSAR models (Tables S3−S6). A library of 50,000 druglike compounds procured from ChemDiv and ChemBridge databases and housed in the Institutional compound library of IIIM was used for screening of potent CDK2 inhibitors. This compound library is referred to as the Institutional compound library in the rest of the article. 2.2. Ligand Guided Screening Strategies. 2.2.1. Substructure and Similarity Search. Substructure and similarity search identifies compounds that contain the substructure of the active scaffold/query compound and the compounds which share a similar structure to the query compound, respectively. Based on the widely accepted assumption that compounds which share a similar structure have similar chemical properties and

Table 1. Data Set of Compounds Taken for Building Pharmacophore Based 3D-QSAR sets

series

Set-A Set-B Set-C Set-D

no. of compds in data set

active:moderately active:less active compds

active data set for pharmacophore modeling

training:test set compds

117 65

89:18:10 46:13:6

78 37

88:29 50:15

155 45

101:38:16 33:8: 4

82 28

117:38 34:11

indenopyrazole derivatives pyrazolo[4,3-h]qinazoline-3carboxamides pyrazole-pyrimidine derivatives 3-aminopyrazole series

Table 2. Statistical Parameters of the Four Congeneric Series series

PLSa

hypothesis

R2b

Q2c

stability

RMSEd

no. of hits

Set-A Set-B Set-C Set-D

6 6 6 5

AAADDR AAADDH AAAHHR ADDHR

0.932 0.828 0.804 0.946

0.841 0.702 0.735 0.505

0.744 0.492 0.511 0.345

0.355 0.41 0.22 0.287

166 122 156 289

a

PLS: partial least-squares. bR2: corelation coefficient of the predictive model. cQ2: cross-validated coefficient of the predictive model. dRMSE: root mean square error in the test set prediction. 1958


Article


coefficient of the training set (R2) and the internal predictivity coefficient of the test set (Q2) compounds. 2.3. Structure Guided Screening Strategy. 2.3.1. EnergyOptimized Pharmacophore Mapping (e-Pharmacophore).

training and test sets (Table 1 and Tables S3−S6). The data set of training set compounds was considered for building 3D-QSAR models and validated through test set compounds. The best hypothesis was selected based on the values of the correlation

Figure 1. Best pharmacophore hypothesis from each series aligned on the reference ligand here: light red spheres represent hydrogen bond acceptors (A), orange torus rings are aromatic (R), green spheres are hydrophobic (H), and light blue spheres are the hydrogen bond donor (D) which defines the geometric location of the pharmacophoric site. A) Best pharmacophore hypothesis AAADDR (indenopyrazole derivative) aligned on best-fit compound 81. B) Best pharmacophore hypothesis AAADDH (pyrazolo[4,3-h]qinazoline-3-carboxamides) aligned on best-fit compound 24. C) Best pharmacophore hypothesis AAAHHR (imidazole pyrimidine amides analogs) aligned on best-fit compound 65. D) Best pharmacophore hypothesis ADDHR (3-aminopyrazoles) aligned on best-fit compound 17.

Figure 2. ROC plots of screening a data set of CDK2 inhibitors and decoys at 3S1H using Glide docking: A) HTVS scoring function, B) SP scoring function, and C) XP scoring function. 1959


Article


without disturbing their coordinates, using the XP scoring function of Glide.52 The e-Pharmacophore model generates a receptor based excluded volume with the energy contribution of each atom present in the ligand calculated through the Glide scoring function, and the pharmacophoric feature of the Phase module was utilized to screen the database. The e-Pharmacophore models built on 27 X-ray crystals with their respective PDB Ids are given in Figure S1 generated using e-Pharmacophore mapping of Schrodinger.53 The identified compounds from this study will contribute to important structural features with descriptor which must be important for the affinity as well as the efficacy of the compound for the therapeutic target. 2.3.2. Molecular Docking. In PDB there are nearly 200 X-ray crystal structures of CDK2 reported with cocrystallized inhibitors. In order to select the most appropriate structure for carrying out the molecular docking studies, root mean square deviation (RMSD) calculation was performed on 27 representative crystal structures containing different core moieties (Table S1). The bound inhibitors were prepared and docked at the active site of CDK2 using XP, SP, and HTVS scoring functions of Glide.52 RMSD was calculated between bound and docked conformations of the cocrystallized inhibitor. The lesser the RMSD difference between bound and docked conformations, the greater would be the chance to attain bioactive conformations of the compounds identified from in silico docking studies. The calculated RMSD of all the 27 cocrystallized inhibitors at ligprep and confgen using XP, SP, and HTVS scoring functions is given in the Supporting Information with their respective PDB Ids (Table S7). Eleven PDB Ids, which showed an RMSD of less than 1.5 Å, were selected for enrichment factor (EF) calculation. For EF calculation, a decoy library of 1000 compounds (provided by Schrodinger) and 77 CDK2 inhibitors (from literature) were

Table 3. AUC, ROC, RIE, and BEDROC Values Calculated at HTVS, SP, and XP Docking Results for the 3S1H CDK2 Crystal Structure scoring function

AUCa

ROCc

RIEb

BEDROC (α = 20)d

HTVS SP XP

0.84 0.91 0.94

0.87 0.94 0.97

4.09 5.26 6.27

0.382 0.491 0.585

a

AUAC: area under the curve is the probability that a randomly chosen known active will rank higher than a randomly chosen decoy (the value is bounded between 1 and 0, with 1 being ideal screen performance). bROC: Receiver Operator Characteristic area under the curve (value is bounded between 1 and 0, with 1 being ideal screen performance and 0.5 reflecting random behavior). cRIE: Robust Initial Enhancement. Active ranks are weighted with a continuously decreasing exponential term (large positive RIE values indicate better screen performance). dBEDROC: Boltzmann-enhanced Discrimination Receiver Operator Characteristic area under the curve. The value is bounded between 1 and 0, with 1 being ideal screen performance.

e-Pharmacophore, an amalgamation of molecular docking and pharmacophore modeling, is mainly used to screen large compound libraries. Here the 27 crystal structures of CDK2 were considered for building e-Pharmacophore models. As the bioactive conformations of CDK2 cocrystallized ligands were considered for this study, it increases the chance of identifying potent CDK2 inhibitors from screening the Institutional compound library to having similar pharmacophore features and interactions within the protein−ligand complexes. Initially these 27 crystal structures, bound with diverse inhibitors, were prepared using the protein preparation wizard of Schrodinger. Further, on these 27 prepared crystal structures their respective cocrystallized inhibitors were docked at their reference position

Figure 3. In silico screening strategy applied for screening the ChemDiv and ChemBridge library. 1960


Article


Figure 4. Graph plot of IC50 calculation.

2.4. In Silico Screening Protocol. All three in silico screening filters, i.e. e-Pharmacophore mapping, pharmacophore based 3D QSAR screening, and substructure and similarity search, were applied in parallel to screen the Institutional compound library to identify potent compounds for CDK2 inhibition. Common compounds were identified among these three in silico screening filters. The compounds attained via

prepared and docked at the respective grids of the selected PDB IDs using various scoring functions. PDB Id which identifies the maximum number of actives in the top 1% via screening this decoys + actives library was considered for molecular docking studies.54,55 This approach would increase the chances of attaining bioactive conformations and identification of the maximum number of actives from screening the large compound libraries. 1961


Article


Figure 5. Structure of seven identified hits. The core moieties are shown in bold.

applying these filters were subsequently docked on CDK2 in order to analyze their interaction with the critical amino acid residue i.e. the Leu83 residue reported for CDK2 inhibition.56 A dock score with a cut off value of −8.0 was kept as a benchmark for the selection of compounds. Considering the docking score threshold, interaction with the critical residue, and diversity of identified compounds, certain compounds were selected for their bio evaluation. 2.5. In Vitro Kinase Assay. ADP-Glo Kinase Assay is a lumniscent kinase assay that measures ADP formed from a kinase reaction; ADP is converted into ATP, which is converted into light by Ultra-Glo Luciferase. The lumniscent signal positively correlates with the ADP amount and kinase activity. The assay is well suited for measuring the effects chemical compounds have on the activity of a broad range of purified kinases making it ideal for both primary screening as well as kinase selectivity profiling. Briefly, the assay is performed in white 96-well plates taking both the reaction mixture (kinase reaction in the presence of substrate) and blank control (kinase reaction in the absence of substrate) into consideration. The respective CDK/cyclin reaction is initiated by the addition of 5 μL of 250 μM ATP assay solution (1 mL of ATP assay solution is prepared by adding 25 μL of ATP solution (10 mM) to 500 μL of 2X buffer and 475 μL of dH20) bringing the final volume up to 25 μL, and the reaction mixture is incubated at 30 °C for 15 min. After the incubation period, the reaction termination and the remaining ATP depletion is done by adding 25 μL of ADP-Glo Reagent to each well, and the reaction mixture is further incubated at ambient temperature for another 40 min. After this, 50 μL of Kinase Detection Reagent (prepared by mixing Kinase Detection Buffer with Lyophillized Kinase Detection Substrate) is added to each well, and the plate is incubated again for 30 min. Finally, the 96-well reaction plate is read on a luminescence plate reader, and the ADP produced (nmol) in the presence and absence of the substrate is determined. Percent Kinase Inhibition is calculated as

(Molecular Mechanics Generalized Born Surface Area) DG bind, which is calculated as under ΔG(bind) = E_complex(minimized) − (E_ligand(minimized) + E_receptor(minimized))

where E_complex(minimized) is the absolute free energy of the complex, E_receptor(minimized) is the absolute free energy of the protein, and E_ligand(minimized) is the absolute free energy of the ligand. 2.7. Molecular Dynamics. The binding pocket of CDK2 comprises of the following hydrophobic residues i.e. Ile 10, Val18, Ala31, Val64, Phe80, Phe82, Leu83, Leu134, Ala144; polar residues Thr14, His84, Gln85, Gln131; negatively charged Glu51, Glu81, Asp86, Asp145; positively charged Lys33, Lys89, Gly11, and Gly13 which resides within 4 Å from the centroid of the active site. In order to identify the amino acid residues responsible for providing stability to the compounds found active in in vitro studies for CDK2 inhibition, the docked pose of these actives was considered for MD simulation studies. Initially the protein−ligand complex was merged in the SCP solvent system with cubic boundary conditions. Each system was neutralized by adding Na+ and Cl− counterions for simulation studies. All the complexes were minimized at the default parameters. The simulation was performed in the NPT ensemble for 20 ns. The Maestro-Desmond interoperability tool (version 4.1, Schrodinger, LLC, 2015) was used for performing all the MD studies.57

3. RESULTS AND DISCUSSION 3.1. Ligand Guided Screening Strategies. 3.1.1. Substructure and Similarity Search of Structurally Diverse Known CDK2 Inhibitor. The 27 diverse scaffolds and 77 CDK2 inhibitors resulted in the identification of 8,422 and 897 compounds respectively from substructure and similarity search against the Institutional compound library. 3.1.2. Pharmacophore Based 3D-QSAR Modeling. Four predictive pharmacophore based 3D-QSAR models with R2 close to 1, Q2 more than 0.5, and stability higher than root mean square error (RMSE) were selected to screen the Institutional compound library. The values of statistical parameters identified for the selected predictive pharmacophore models for Set A-D series are given in Table 2. Six point pharmacophoric features were obtained from Set A-C series, which resulted in the identification of 444 compounds, whereas the 5 point

% Kinase Activity Luminescence of Test − Luminescence of Blank × 100 = Luminescence of Control − Luminescence of Blank

% Kinase Inhibition = 100 − % Kinase Activity

2.6. Free Binding Energy Calculation. The compounds which showed higher binding affinity were further considered for binding free energy calculation using Prime MMGBSA 1962


Article

Journal of Chemical Information and Modeling pharmacophore hypothesis of Set-D series could identify 289 compounds upon screening the Institutional compound library (Table 2). All the models showed good prediction of the test set compounds as indicated by the Q2 values (Table 2) which describes the robustness of the QSAR model. The pharmacophoric features of these Set-A to Set-D series with the best hypothesis with the best fit compound are shown in Figure 1. The intersite distances between the pharmacophoric features of the four models were also calculated (Tables S8−11). The plots

of actual vs predicted activities of training and test set compounds of four series are depicted in Figures S2−S5. The pharmacophore screening of the four congeneric series resulted in the identification of the total 733 compounds. 3.2. Structure Guided Screening Strategy. 3.2.1. e-Pharmacophore Screening. The bound coordinates of diverse CDK2 inhibitors were considered to generate 27 e-Pharmacophore models. Screening of the Institutional compound library through these e-Pharmacophore models resulted in the

Table 4. Structure of Actives with Core Moieties, Protein Ligand Interaction Studies before and after MD Simulation

1963


Article

Journal of Chemical Information and Modeling Table 4. continued

by applying ligand and structure based in silico filters in parallel i.e. substructure plus similarity search, e-Pharmacophore mapping, and pharmacophore based 3DQSAR modeling which resulted in the identification of 9253, 11736, and 733 compounds respectively from screening the Institutional compound library. A large number of compounds were identified from each ligand based screening filter. Thus, in order to reduce the number of false positives, the parallel selection approach was applied. These identified compounds were further subjected to molecular docking studies on the crystal structure of CDK2 (3S1H) using the XP scoring function. The compounds having a dock score of

Fusion of Structure and Ligand Based Methods for Identification of

Recommend Documents