Pharmacophore-Based Virtual Screening and Experimental Validation of Pyrazolone-Derived Inhibitors toward Janus Kinases.

Janus kinases (JAKs) are nonreceptor protein tyrosine kinases that play a role in a broad range of cell signaling. JAK2 and JAK3 have been involved in the pathogenesis of common lymphoid-derived diseases and leukemia cancer. Thus, inhibition of both JAK2 and JAK3 can be a potent strategy to reduce the risk of these diseases. In the present study, the pharmacophore models built based on the commercial drug tofacitinib and the JAK2/3 proteins derived from molecular dynamics (MD) trajectories were employed to search for a dual potent JAK2/3 inhibitor by a pharmacophore-based virtual screening of 54 synthesized pyrazolone derivatives from an in-house data set. Twelve selected compounds from the virtual screening procedure were then tested for their inhibitory potency against both JAKs in the kinase assay. The in vitro kinase inhibition experiment indicated that compounds 3h, TK4g, and TK4b can inhibit both JAKs in the low nanomolar range. Among them, the compound TK4g showed the highest protein kinase inhibition with the half-maximal inhibitory concentration (IC50) value of 12.61 nM for JAK2 and 15.80 nM for JAK3. From the MD simulations study, it could be found that the sulfonamide group of TK4g can form hydrogen bonds in the hinge region at residues E930 and L932 of JAK2 and E903 and L905 of JAK3, while van der Waals interaction also plays a dominant role in ligand binding. Altogether, TK4g, found by virtual screening and biological tests, could serve as a novel therapeutical lead candidate.

Introduction

Janus tyrosine kinases (JAKs) are described as nonreceptor protein tyrosine kinases that play an important role in a broad range of cell signaling, such as cell growth, survival, development, and differentiation of cells.[1] The group of JAKs consists of four members, including JAK1, JAK2, JAK3, and TYK2 (tyrosine kinase 2).[2] These enzymes are part of the signal transducers and activators of the transcription (STAT) pathway that is activated by cytokines and induces a cascade of signals for development or homeostasis.[3] JAK1, JAK2, and TYK2 are recognized to be ubiquitously expressed, whereas JAK3 is chiefly expressed in lymphoid tissues and appears to be a selective regulator of lymphoid development and control of function within the immune system.[4] Among them, JAK2 and JAK3 have been associated with myeloid leukemia and inflammatory and autoimmune diseases.[5,6] JAK2 is the only member of the JAK family that pairs with itself. It associates with multicytokine pathways along with other JAK members for hormones like cytokines such as erythropoietin (EPO), thrombopoietin (TPO), and cytokine receptor ligands, involved in hematopoietic cell development such as interleukin-3 (IL-3), interleukin-6 (IL-6), and granulocyte-macrophage colony-stimulating factor (GM-CSF).[7−9] Imbalance in IL-3 and GM-CSF signaling can increase differentiation, leading to lymphoid-derived diseases and myelofibrosis (bone marrow/blood cancers).[10−12] JAK2 mutations such as V617F have been identified in people with myelofibrosis, and aberrant IL-6/JAK2/STAT3 signaling has an essential role in solid tumors of colorectal cancer.[13] JAK3 also plays a crucial role in lymphoid development associated only with the common gamma-chain via the IL-2, IL-4, and IL-7 pathways. The imbalance of IL-2 signaling also contributes to the pathogenesis of lymphoid-derived diseases.[14−16] JAK3 mutations (e.g., A572 V and A573 V) could lead to continual activation of the JAK-STAT signaling, resulting in various leukemias and lymphomas, including monomorphic epitheliotropic intestinal T-cell lymphoma, T-cell acute lymphoblastic leukemia, and hepatosplenic T-cell lymphoma.[17,18] According to the aforementioned statements, JAK2 and JAK3 are considered to be proven therapeutic targets due to their roles in the signaling transduction pathways involved in immune function and cancerous conditions.[19] Therefore, the development of dual JAK2/3 inhibitors could be an effective way to treat a variety of diseases such as lymphoid-derived diseases, myelofibrosis, and cancers.

Several drugs targeting JAKs have been developed to treat rheumatoid arthritis (RA), myelofibrosis, psoriasis, leukemia, and lymphoma.[20−22] Ruxolitinib is the first Food and Drug Administration (FDA)-approved JAK inhibitor (half-maximal inhibitory concentration (IC50) was 3.3 nM for JAK1 and 2.8 nM for JAK2)[23,24] for the treatment of myelofibrosis. Tofacitinib (CP-690550) is a JAK1/2/3 inhibitor with the IC50 from enzyme assay of 15.1, 77.4, and 55.0 nM, respectively, which decreases lymphocyte activation and proliferation for autoimmune diseases.[24,25] In addition, small molecules from natural compounds such as pyridines,[26] flavonoids,[27] and pyridazines[28] have been described as JAKs inhibitors. Also, inhibitors of JAK2 and JAK3 derived from pyrazolone-containing compounds have been reported. Zak et al. found that, from an enzyme-based assay, the methyl derivative of pyrazolone linked with imidazopyridine showed the IC50 value of 250 nM for JAK2.[29] A pyrazol-3-ylamino pyrazine compound containing methyl and cyanide substitutions was effective against JAK2 (IC50 was 3 nM) and JAK3 (IC50 was 11 nM).[30] In addition, the aminopyrazole analogue containing a chloride atom showed potent inhibition on both JAKs (2.2 nM for JAK2 and 3.5 nM for JAK3).[31] Moreover, pyrazolo-nicotinonitrile (AZ960) has been reported as a potent compound toward JAK2 expressing SET-2 cell line (IC50 of <3 nM).[32]

Virtual screening has been recognized as an important in silico technique for finding novel compounds in drug development.[33] Several JAKs inhibitors derived from pharmacophore-based screening techniques have been reported.[13,34,35] The pharmacophore model obtained from the docked JAK3/tofacitinib complex showed four relevant sites—one hydrophobic center, one hydrogen-bond donor, and two hydrogen-bond acceptors.[34] The reference ligand phenylaminopyrimidine for dual inhibitions of JAK2/3 showed the following pharmacophore sites: two aromatic features, one hydrogen-bond donor, and hydrogen-bond acceptor predicted by PHASE.[35] The chloroquinoline-4-amine derivatives (NSC13626) for JAK2 inhibitors (Kd = 6.6 μM) derived from the NCI database were obtained by molecular docking. There, hydrogen-bond interactions with the E930 and L932 at the hinge region as well as D994 and F995 at the activation loop (A-loop) of JAK2 were found.[13]

Herein, the combination of pharmacophore-based virtual screening (PBVS) with molecular docking was applied to identify a novel dual inhibitor toward JAK2 and JAK3 from the pyrazolone derivatives in our in-house data set (Figure ). First, the molecular dynamics (MD) trajectories of tofacitinib/JAKs complexes taken from the previous study[36] were used to build the pharmacophore models. The best pyrazolone derivatives obtained from screening compounds were selected for in vitro JAK2/3 kinase inhibitory activity studies. Finally, the molecular interactions between the potent compound and JAKs were studied by molecular docking followed by MD simulations.

Figure 1

Virtual screening workflow uses a combination of in silico PBVS and molecular docking of JAK2/3 inhibitors derived from pyrazolone derivatives (Figures S1), followed by testing by enzyme-based assay.

Results and Discussion

Pharmacophore Models

The MD snapshots of tofacitinib and JAK2/3 complexes in aqueous solution from the last 150 ns simulation times (1500 MD snapshots) obtained from our previous study[36] were used as a template to generate pharmacophore models for PBVS. To reduce the computational cost and time, the representative pharmacophore models (RPMs) were selected by clustering from the pharmacophore models of each system (Figure ). Several chemical features such as hydrophobicity (yellow spheres) and hydrogen-bond donor (HBD, green arrow) and acceptor (HBA, red arrow) are defined according to the binding interaction of tofacitinib on residues in the adenosine triphosphate (ATP)-binding pocket of JAK2/3.

Figure 2

(A) Combination of 2D and 3D pharmacophores in one model derived from the last 150 ns selected from key interactions and (B) merge of representative RPMs derived from the last 150 ns of the simulation times of tofacitinib in the binding pocket of JAK2/3. The yellow color sphere and green and red color arrows represent pharmacophores with hydrophobic interactions, hydrogen-bond donor (HBD) and acceptor (HBA) abilities, respectively.

The results revealed that crucial pharmacophore features of both JAKs and tofacitinib were slightly similar, consisting of (i) hydrophobic interactions with V863 and L855 in Glycine loop (G-loop), A880, V911, and M929 in the hinge region and L983 in the catalytic loop for JAK2; V836 (G-loop), A853, V884, and L956 in the catalytic loop as well as M902 in the hinge region for JAK3, (ii) HBDs in the hinge region with E930 for JAK2; E903 for JAK3, and (iii) HBAs in the hinge region with L932 and L905 for JAK2 and JAK3, respectively (Figure A). The methyl group of the piperidine ring occupies a small hydrophobic pocket in the C-lobe, in agreement with decomposition free energy for both JAKs.[36] Each RPM result from pharmacophores clustering from the last 150 ns of both systems leads to the conclusions (Figure S2). The JAK2 found (i) 7 features consisting of hydrophobic interactions (16), HBDs (2), and HBAs (1), whereas JAK3 found (ii) 8 features including hydrophobic interactions (20), HBDs (4), and HBAs (4). As a result, it was found that the JAK3 system gave the hydrophobic interaction, and HBD and HBA more than the JAK2 system. The results revealed that tofacitinib fitted well within the active site of JAK3, resulting in a reasonably higher representative pharmacophore model than the JAK2 system. Then, the RPMs of both systems were summarized in Figure B. These pharmacophore patterns could be used to identify the vital contributing features within the binding site of both JAKs. Therefore, these representative pharmacophore models of both JAKs were used as a template for JAKs inhibitor screening. Furthermore, Figure shows the ratio of pharmacophore occurrences at greater than 80%. The pharmacophore feature of pyrrolo[2,3-d]pyrimidine scaffold of tofacitinib forms a high appearance of two hydrogen bonds in the hinge region residues E930 (99%) and L932 (97%) for JAK2 as well as E903 (99%) and L905 (97%) for JAK3, which agrees well with the previous report,[34,36] indicating that hydrogen bonds are important interactions between tofacitinib and both JAKs. In addition, tofacitinib/JAKs complexes were also stabilized by hydrophobic interactions.

Figure 3

Interaction map of tofacitinib in the binding pocket of JAK2/3 derived from the last 150 ns of the simulation times. The H, HBA, and HBD abbreviations represent the hydrophobic interaction, hydrogen-bond acceptor, and hydrogen-bond donor pharmacophore features.

Virtual Screening

The RPMs obtained from JAK2/3 and tofacitinib complexes were used as templates for the screening of novel inhibitors against two JAKs by using PBVS (Figure ). Structures of selected pyrazolone derivatives (54) (given in Figure S1) and tofacitinib were included in an in-house library database. The hits from the in-house library were filtered by the pharmacophore fit score of the tofacitinib template for JAK2 (34.37) and JAK3 (36.39). The 23 hit compounds for JAK2 inhibitors with higher pharmacophore fit scores than tofacitinib were found (scores between 35.71 and 46.77) as well as 28 hit compounds for JAK3 inhibitors (scores between 36.44 and 38.72).

To validate the structure-based three-dimensional (3D) pharmacophore models, the involved virtual screening of 23 active compounds, tofacitinib, and 6588 decoy compounds for JAK2 and 28 active compounds, tofacitinib, and 4673 decoy compounds for JAK3 were used to analyze the receiver operating characteristics (ROC, Figure ). The ROC curve is a graphical representation of the sensitivity (proportion of true positives) as a function of specificity (proportion of false positives).[37,38] We found that the AUC100% (area under the curve) value was 0.86 and 0.89 for JAK2 and JAK3, respectively. The AUC value is close to 1.00, indicating that the pharmacophore model is able to identify the true active compounds from the decoy JAKs compounds. Therefore, hit compounds toward both JAKs derived from PBVS are acceptable for this study.

Figure 4

ROC curve validation of the 3D structure-based pharmacophore model of JAK2/3 inhibitors.

Subsequently, hit compounds with JAK2/3 as previously mentioned were screened by molecular docking using two different programs, namely, GOLD and FlexX (Figure S3). Note that both types of docking programs use different algorithms. GOLD uses a genetic algorithm to explore the ligand conformational flexibility with partial flexibility of the protein,[39] while FlexX is based on an incremental construction algorithm including three phases, (i) base selection, (ii) base placement, and (iii) complex construction.[40] If the results are in similar patterns; therefore, the results should be reliable. Docking results from GOLD and FlexX showed a similar pattern (Figure S3A). The compounds with fitness score and docking energy better than tofacitinib (55.50 by GOLD and −17.11 kcal/mol by FlexX for JAK2; 63.45 by GOLD and −24.90 kcal/mol by FlexX for JAK3) were selected (summarized in Figure S3B). Fifteen compounds obtained from GOLD and FlexX could be inhibitors against JAK2, whereas 20 and 12 compounds obtained from GOLD and FlexX, docking, respectively, could be inhibitors for JAK3. Altogether, overlapped screened compounds from both docking programs for JAK2/3 were 12 compounds.

The two-dimensional (2D) structures of these 12 screened pyrazolone derivatives are shown in Figure . The obtained results suggested that the 12 screened pyrazolone derivatives consisted of diverse functional groups, including aromatic rings with nitro (-NO2), halogen (Cl), carbonyl (C=O), and sulfonamide (-S(=O)2-NH2) groups, could be inhibitors of JAK2/3. These functional groups could interact within conserved regions of JAKs, including the hinge region, which is the region that accommodates the adenine ring of ATP or G loop, which controls the mobility of ligand binding.[41] For example, the previously reported −NH and =N moieties of pyrazole in 4-amino-(1H)-pyrazole compounds could form hydrogen bonds with the hinge region of JAK2, including E930 and L932,[42] and pyrazolo-pyrimidine analogues with fluorine atoms showed potent JAK2 inhibition (27 nM).[43] A halogen-substituted compound has also been considered in drug discovery and design.[44−46] The pharmacological properties of these compounds were also predicted to verify whether these screened pyrazolone derivatives can be used as a good candidate for JAKs inhibitors. All the screened compounds obey Lipinski’s rule[47] (Table S1), suggesting that they could be candidates for novel JAK2/3 inhibitors. Altogether, these screened compounds were selected to investigate the JAK 2/3 kinase inhibitory activity.

Figure 5

2D structure of screened compounds. R3 is equal to −CH3 for all compounds.

Janus Kinase Inhibitory Activity

JAK2 and JAK3 inhibitory activities of screened pyrazolone derivatives at concentrations of 1 μM were evaluated, and tofacitinib was used as a positive control (Figure ). Almost all screened compounds showed considerable JAK3 inhibition activity at 1 μM. Among them, 3h (65.86 ± 2.13%), TK4b (57.40 ± 1.92%), and TK4g (79.35 ± 0.20%) displayed the best inhibitory activity toward JAK2. Therefore, these three compounds were selected to determine the IC50 values for JAK2 and JAK3. The IC50 values of all studied compounds against JAKs are summarized in Table and Figure S4. The IC50 values of 3h, TK4b, and TK4g toward JAK2 were 23.85 ± 0.35, 19.40 ± 1.14, and 12.61 ± 1.25 nM, respectively. These results correspond to the IC50 values of pyrazolone derivatives consisting of chloro-methoxyphenyl-amino groups (19 nM) and pyrrole carboxamide (20 nM) toward JAK2 from previous reports.[48,49] The IC50 values of the three mentioned compounds toward JAK3 were 18.90 ± 1.13 nM for 3h, 18.42 ± 0.34 nM for TK4b, and 15.80 ± 0.81 nM for TK4g. Previous reports showed that an aminopyrazole analogue with three chlorine atoms, which is to some extent similar to compound 3h compound, gave IC50 values toward JAK2 and JAK3 of 2.2 and 3.5 nM, respectively.[31] These values are lower than our results (23.85 and 18.90 nM). On the other hand, the reported IC50 of an aminopyrazole analogue consisting of heterocycles as derivatives against JAK2 and JAK3 was 98 and 39 nM, respectively, which are higher than our findings.[50] The previous report showed that the IC50 value of tofacitinib against JAK2/3 was 77.4 nM for JAK2 and 55.0 nM for JAK3.[51] These values are in a similar range to ours. However, they showed the IC50 values in a similar trend. Compound TK4g against JAK2 showed an IC50 value significantly different from that of tofacitinib (p ≤ 0.05), indicating that this compound could bind with JAK2 better than tofacitinib[36] (Table ), whereas the other two compounds (3h and TK4b) bind to JAK2 in a manner similar to tofacitinib. In addition, three compounds and tofacitinib could inhibit JAK3 at a similar level. Due to the high ability of binding between TK4g and JAK2, this compound was selected for further study to understand the molecular interactions toward both JAKs.

Figure 6

Relative inhibition of Janus kinases activity derived from hit pyrazolone derivatives at 1 μM concentrations. Bars represent as mean ± SEM. * p ≤ 0.05, ** p ≤ 0.01, and *** p ≤ 0.001 vs tofacitinib.

Table 1

IC50 Value of Enzymatic Activities and Docking Energies between Potent Pyrazolone Derivatives and JAK2/3

	JAK2			JAK3
		GOLD score	binding energy (kcal/mol)		GOLD score	binding energy (kcal/mol)
compounds	IC50 (nM)a	GOLD	FlexX	IC50 (nM)a	GOLD	FlexX
3h	23.85 ± 0.35	59.58	–17.22	18.90 ± 1.13	66.53	–26.18
TK4b	19.40 ± 1.41	59.33	–17.87	18.42 ± 0.34	65.48	–28.21
TK4g	12.61 ± 1.25b	66.88	–24.01	15.80 ± 0.81	71.24	–27.37
Tofacitinibc	26.90 ± 0.57	55.50	–17.11	20.69 ± 0.42	63.45	–24.90

Results are presented as mean ± SEM.

p ≤ 0.05 vs tofacitinib.

Results are obtained from previously reported.[36]

Molecular Interactions of Pyrazolone Derivatives from Molecular Docking

The docking score of pyrazolone derivatives, including 3h, TK4b, and TK4g, which showed good inhibitory activity against JAK2/3, are summarized in Table . GOLD and FlexX docking results demonstrate that all potent compounds show higher binding ability toward JAK2/3 than tofacitinib. In the case of the GOLD, a higher score value means better binding; the docking score of JAK2/3 with TK4g (66.88 and 71.24) was higher than that of 3h (59.58 and 66.53), TK4b (59.33 and 65.48), and tofacitinib (55.50 and 63.45), while, with FlexX docking results, a lower value means better binding. It was shown that the complexation of JAK2/3 with TK4g (−24.01 and −27.37 kcal/mol) exhibits a lower energy value than with 3h (−17.22 and −26.18 kcal/mol), TK4b (−17.87 and −28.21 kcal/mol), and tofacitinib (−17.11 and −24.90 kcal/mol). Interestingly, docking results correspond with IC50 values between three interested compounds and JAK2/3. As a result, the docking score of the TK4g compound toward JAK2/3 showed the highest binding energy (most negative energy value). Therefore, this compound was selected for the study of molecular interactions. Superimpositions between TK4g and tofacitinib from both programs were done (Figure S5). A similar orientation of both JAKs with TK4g was generated from GOLD and FlexX, where the sulfonamide group -S(=O)2-NH2 interacts with the hinge region. However, the structures of TK4g and JAK2/3 generated from GOLD were selected to study the dynamic of the binding mode. Note that the binding affinity of TK4g against JAK1 was investigated. We found that the complexation of JAK1 with TK4g (−24.43 kcal/mol from FlexX and 67.84 from GOLD) exhibits a higher binding affinity than tofacitinib (−17.51 kcal/mol and 54.46 from FlexX and GOLD, respectively). Therefore, JAK1 could be inhibited by TK4g.

Dynamic Interactions of Pyrazolone Derivatives

The binding mechanisms of TK4g against JAK2/3 along 150 ns simulations were studied by (i) root-mean-square deviation (RMSD), (ii) #atom contacts, and (iii) #H-bonds (# = “number of”; Figure S6). It was found that TK4g was stable within the active site along the simulation times, as the RMSD result showed low fluctuation. In addition, the #atom contacts (15.46 ± 4.39 for JAK 2 and 16.04 ± 5.09 for JAK 3) and #H-bonds (1.83 ± 0.85 for JAK 2 and 1.88 ± 1.16 for JAK 3) of TK4g binding against both JAKs was obtained. However, the last 50 ns of each system, which showed stability within the active site (Figure ), was selected to perform the binding pattern by MM/GBSA per-residue decomposition energy method. The interactions between TK4g and JAK2/3 are illustrated in Figure .

Figure 7

Structural snapshot per time of last 50 ns simulations of TK4g against JAK2/3.

Figure 8

Binding patterns between TK4g (black ball-and-stick model) and JAK2/3 derived from the last 50 ns MD simulation are described as follows. (A) Per-residue decomposition free energy and (B) percentage of hydrogen-bond occupation. The binding orientation of TK4g within the binding pocket obtained from the representative MD snapshot. The lowest and highest energies are ranged from green to light gray, respectively.

The pyrazolone core structure and the aromatic rings of TK4g were stabilized by crucial residues of JAK2/3 in the four regions as follows (Figure A): (i) hinge region; M929, E930, Y931, L932 of JAK2, and E903, Y904, L905 of JAK3, (ii) G-loop; L855, G856, K857, G858, S862, V863 of JAK2, and L828, G829, V836 of JAK3, (iii) catalytic loop; L983 of JAK2, and I955, L956 of JAK3, and (iv) DFG activation loop; D994 of JAK2 and near DFG activation loop, A966 of JAK3. The hydrophobic interaction between TK4g and the hinge region residues M929 and Y931 of JAK2 was similar to the tetrazole ring interaction.[52] Our binding model reveals that compound TK4g is located in the binding pockets of both JAKs mainly by the van der Waals interactions (vdW) (Figure S7). In addition, the methyl group of the pyrazolone ring forming van der Waals contacts pointed toward the C-terminal lobe, with (i) catalytic loop residue L983 for JAK2 and (ii) I955 and L956 at the catalytic loop and C909 for JAK3. These findings are consistent with a previous study reporting that hydrophobic region pyrrolepyrimidine and piperidine of tofacitinib can interact with the L983 residue of JAK2.[53] Furthermore, the TK4g was stabilized by vdW in the catalytic loop of JAK3 (I955 and L956), similar to tofacitinib/JAK3 binding, which corresponds with the previous report.[54] Therefore, the vdW interaction of the hydrophobic aromatic ring in the binding pocket of JAKs plays a vital role in the potency of the inhibitors.

Moreover, hydrogen bonds are also crucial for TK4g binding. The sulfonamide group of TK4g forms two hydrogen bonds with the hinge region of both JAKs in the residue as follows; (i) oxygen atom in −SO2 group with leucine residue in position 932 (73.42%) of JAK2 and 905 (91.86%) of JAK3, which corresponds with an oxygen atom in the carbonyl group (C=O) of pyrazolone derivative compound showed hydrogen bonding with L932,[50] (ii) −NH group with E930 (70.70%) of JAK2 and E903 (99.32%) of JAK3 (Figure B). The hydrogen bonding result in these residues in the hinge region is similar to interactions between the pyrrolopyrimidine ring of tofacitinib and JAK2/3 in the previous report.[55] Therefore, the hydrogen interactions are crucial for the JAK2/3 inhibition. Due to their functions in signaling transduction pathways implicated in immune regulation and cancerous disorders, JAK2 and JAK3 are considered therapeutic targets. Others drugs approved by the FDA for JAKs that target the ATP-binding site have been found, including baricitinib (JAK1/2)[24] and fedratinib (JAK2).[56] Furthermore, pyrazolone analogues for JAK2 inhibitors in investigational status have been reported, including AZD1480 (phase 1)[57] and momelotinib (phase 2).[58] Therefore, pyrazolone derivatives in this work might be used as an inhibitor toward JAK2/3. Note that, the aforementioned drugs, along with other drugs, ruxolitinib (JAK1/2); gandotinib and pacritinib (JAK2); decernotinib and peficitinib (JAK3)[59] were studied for JAK2/3 binding in comparison to our TK4g and tofacitinib by GOLD molecular docking (Figure S8). Among drugs, tofacitinib showed the lowest JAK2/3 inhibitions. Furthermore, it was found that the TK4g could be able to inhibit JAK2 better than ruxolitinib, momelotinib, and pacritinib, whereas TK4g showed JAK3 inhibition better than decernotinib and peficitinib.

Conclusions

The development of inhibitors against human kinases has emerged as a promising strategy for decreasing the risk of diseases. In this work, pharmacophore-based virtual screening, molecular docking, and drug-like prediction followed by experimental testing were successfully applied to investigate novel dual JAK2/3 inhibitors from 54 in-house synthesized pyrazolone derivatives. The combination of pharmacophore-based virtual screening with molecular docking was used to identify 12 inhibitors toward JAK2 and JAK3. Among these tested compounds, 3h, TK4g, and TK4b showed both JAKs inhibition at the low nanomolar level. Molecular interactions between the most potent compound, TK4g, with JAK2/3 were observed by molecular docking followed by MD simulations. The TK4g binding within JAK2/3 mainly occurs by van der Waals interaction. The complexes with TK4g are stabilized by the sulfonamide group via hydrogen bonding at the hinge region of JAK2 (E930 and L932) and JAK3 (E903 and L905). These findings provide a combination of pharmacophore screening, molecular docking, and experimental testing that seems to be useful for the achievements of potential compounds in the development of potent JAK2/3 inhibitors. This fruitful strategy, a combination process, to search for JAK2/3 inhibitors may lead to higher efficiency in searching for other types of drugs or other enzyme inhibitors.

Materials and Methods

Chemicals and Reagents

Janus kinase 2 (JAK2; SRP0171), Janus kinase 3 (JAK3; SRP0173), and dimethyl sulfoxide (DMSO) were purchased from Sigma-Aldrich (Darmstadt, Germany). The ADP-Glo Kinase Assay kit was purchased from Promega (Madison, WI, USA). The Poly(4:1 Glu, Tyr) Peptide (P61–58) was purchased from SignalChem Biotech (Canada). The series of pyrazolone derivatives was synthesized at the Department of Chemistry, Faculty of Science, Mahidol University.[60−62]

Computational Methods

Pharmacophore-Based Virtual Screening

The geometries of the complexes between tofacitinib and JAK2/3 were obtained by all-atom molecular dynamics (MD) simulations in three independent simulations for 500 ns at 310 K using AMBER16 in the periodic boundary condition as in a previous report.[36] Pharmacophore models of the ligands were used to determine the type and geometric constraints of the chemical features derived from the structure-based method. This is an essential technique to reveal the interaction pattern of the active compound in the binding pocket according to chemical properties and biological interactions. To enlarge the likelihood of accomplishment, pharmacophore features of molecular interaction can be obtained from a collection of MD trajectories. The ligand extraction and pharmacophore feature identification were used to build pharmacophore models from the sets of the trajectory of tofacitinib/JAK(s) complexes from the production phase of the three independent simulations (1500 frames in total). LigandScout 4.2.6 Advance program via the KNIME analysis platform[36,38,63−65] was used for the extended investigation. The pharmacophore features between tofacitinib and JAK(s) were generated using the “MD pharmacophore creator” node in the KNIME program with default parameters. All water molecules and ions for all systems were removed in this step. The pharmacophore models from each system were clustered by features using the “CHA filter” node to reduce the computational time. The representative pharmacophore models (RPMs) were analyzed and visualized again using LigandScout 4.2.6 Advance.[63]

The in-house library of pyrazolone-containing compounds was prepared for virtual screening using the idbgen option in the LigandScout 4.2.6 Advance program.[63] The conformers were produced as a maximum number of 200 conformations for each molecule processed. The pyrazolone derivatives were screened by pharmacophore features using the “IScreen” node of LigandScout 4.2.6 Advance via the KNIME analysis platform.[64] The resulting 3D-interaction feature model was validated for true active compounds from decoys by screening a set of hit compounds and a set of decoys (inactive compounds) obtained from DecoyFinder.[66] Decoy libraries were converted to the 3D multiconformational databases for virtual screening, computing conformations and annotating each conformation with pharmacophore features.

Molecular Docking

The hit compounds resulting from pharmacophore-based virtual screening were studied by molecular docking using the crystal structures of JAK2 (PDB:3FUP) and JAK3 (PDB:3LXK) with tofacitinib bound. The protonation states of all ionizable pyrazolone derivatives were fixed using ChemAxon.[67] Hydrogen atoms were added to the protein. Note that the tofacitinib was redocked into the ATP-binding pocket of the protein for validation of the docking study using FlexX[68] and GOLD.[69] The reference ligand (tofacitinib) was defined as the docking center, and the sphere of a 12 Å radius around tofacitinib was generated for docking. FlexX considers ligand conformational flexibility by an incremental fragment placing technique,[68] and Gasteiger atomic charges[70] were applied. The ligand’s poses were ranked by binding energy. GOLD uses a genetic algorithm technique[69] and Chemscore (rescore) parameters. Poses of all the ligands were sorted based on the GOLD fitness score. Each pyrazolone derivative was docked, collecting 100 docking poses. The docking results were visualized by the UCSF Chimera package[71] and Accelrys Discovery Studio 2.5 (Accelrys Inc.).[72] Then, the pharmacological characteristics of the screened compounds were predicted using the SwissADME program.[73]

MD Simulations

The potent compound (TK4g) against JAK2/3 was performed by all-atom MD simulations for 150 ns in the periodic boundary condition using AMBER20.[74] The ligand was optimized using the Gaussian09 program at the HF/6-31g(d) level. The parmchk module was used to generate the ligand’s restrained ESP (RESP) charges, which were converted from the electrostatic potential (ESP) charges. The protein and ligand were treated with the AMBER ff19SB force field[75] and GAFF2,[76] respectively. All missing hydrogen atoms were added using the tleap module and subsequently were minimized by the 1000 steps of steepest descent (SD) followed by 4000 steps of conjugated gradient (CG). The TIP3P model was used to solvate the system in the 10 Å octa box. The water molecules were minimized using 500 SD steps followed by 1000 CG steps. Then, all complexity was fully minimized using a similar procedure.

Nonbonded interactions were considered using the short-range cutoff of 12 Å, whereas long-range electrostatic interactions were studied using Ewald’s approach.[77] The Berendsen method was used to regulate the pressure.[78] All covalent bonds involving hydrogen atoms were constrained using the SHAKE algorithm.[79] The simulated models were heated to 310 K for 100 ps of relaxation. The temperature was controlled by a Langevin thermostat with a collision frequency of 2.0 ps, and the time step was set to 2 fs.[80−83] Finally, at 310 K, 150 ns NPT simulation of the TK4g/JAK(s) complex was performed. The root-mean-square deviation (RMSD), intermolecular hydrogen bonding, and the number of contact atoms of the TK4g/JAK(s) complex were calculated by the CPPTRAJ module.[84] The criteria between the hydrogen-bond donor (HBD) and hydrogen-acceptor (HBA) were considered as follows: (i) at most 3.5 Å for distance and (ii) at least 120° for the angle. The protein–ligand binding pattern was analyzed using the MM-GBSA per-residue decomposition free energy calculation (ΔGbind) with the MMPBSA.py[85] in AMBER20 with a set of 100 snapshots obtained from the last 50 ns.

Data and Software Availability

The crystal structures of JAK2 (PDB ID:3FUP) and JAK3 (PDB ID:3LXK) with tofacitinib bound are available in the RCSB protein data bank (https://www.rcsb.org/). The free software SwissADME (http://www.swissadme.ch/) was used for predicting the pharmacological properties. The Chimera USCF (https://www.cgl.ucsf.edu/chimera/) and MarvinSketch (https://chemaxon.com/products/marvin) software, which are all free for academic users, were used for molecular visualization and compound structure generation.

Experimental Methods

Janus Kinase Inhibitory Activity Assay

The ability to inhibit kinase activity toward JAK2/3 of the screened pyrazolone derivatives derived from virtual screening at 1 μM was determined by an ADP-Glo Kinase Assay Kit. Tofacitinib was used as a positive control. The reaction system contained 2.5 ng/μL of JAK2 or JAK3 with 5 μM ATP and 2 ng/μL poly(Glu, Tyr) in a buffer containing 40 mM Tris-HCI pH 7.5, 20 mM MgCl2, and 0.1 mg/mL bovine serum albumin and was incubated for 1 h at room temperature. Then, 5 μL of ADP-Glo reagent was added, and the reaction was incubated for 40 min. Kinase detection reagent at 10 μL was added and incubated at room temperature for 30 min to convert adenosine diphosphate (ADP) to ATP. The ATP product was measured by luminescence using a microplate spectrophotometer (Synergy HTX Multi-Mode reader, BioTek, Winooski, VT, USA). The relative inhibition (%) of tofacitinib and the pyrazolone derivatives was calculated according to eq , where positive and negative are, respectively, the reactions with and without JAK2/3, and samples are reactions with inhibitor. The IC50 values were determined using the GraphPad Prism 7.0 software.

Statistical Analysis

Data are expressed as mean ± standard error of the mean (SEM) from triplicate experiments. Statistical analysis between groups of kinase inhibitions was performed using a one-way analysis of variance (ANOVA) with Tukey’s post hoc test. A p-value of at most 0.05 was considered to be significant.