Insights into DFG-in and DFG-out JAK2 binding modes for a rational strategy of type II inhibitors combined computational study

Abstract

A number of type I JAK2 inhibitors binding to the active conformation have emerged so far, while the type II JAK2 inhibitors binding also to the allosteric site adjacent to the ATP-binding pocket in the inactive conformation, have been less studied. In this paper, the individual mechanisms and the differences of the type I and type II inhibitors binding into DFG-in and DFG-out JAK2 were characterized by molecular docking, molecular dynamic simulations, free energy calculations and decomposition. Comparison between the type I thienopyridine inhibitors suggest that the binding affinities are mainly determined by hydrogen bonds with residues E930 and L932 and strong hydrophobic contributions with V863, Y931 and R980 in the ATP-binding pocket. For type II JAK2 inhibitors BBT594 and CHZ868, not only do their type I heads bind in a similar mode in the ATP-binding site, but also the type II tail and the linker can be anchored by hydrogen bonds with E898 and D994 and lead to substantial hydrophobic contributions with residues L902, L893 and F995 in the unique allosteric site. In addition, new type II JAK2 inhibitors are reasonably designed and further evaluated by a comprehensive modeling study. The expected interactions in the ATP-binding and allosteric pocket indicate that the structural and energetic insights into JAK2 inhibitors can facilitate the process of more promising type II JAK2 inhibitors.

---

1. Introduction

The janus kinase 2 (JAK2),¹ belongs to the cytoplasmic non-receptor tyrosine kinases,² and is a member of the Janus family comprising three other members JAK1, JAK3 and TYK2, that play an essential role in the cytokine and growth factor signal transduction pathways.³ JAK2 phosphorylates and activates the signal transducers and activators of transcription (STATs),⁴ for which the constitutive activation of JAK-STAT signaling is associated with hematological neoplasias,⁵ such as thrombocythemia (ET) and polycythemia vera (PV),⁶ primary myelofibrosis (PMF)⁷ and myeloproliferative neoplasms (MPNs).⁸ In 2005, a somatic mutation JAK2V617F in the pseudokinase domain⁹ was identified in MPN-like diseases that the progression of expression of which results in downstream activation of relevant signaling pathways.¹⁰ These findings suggest that JAK2 inhibitors have therapeutic benefit to prevent the dysregulation of JAK2 activity.^11,12

JAK2 kinase exists in an equilibrium of active and inactive states, in which the conformational adjustment induced by the tyrosine phosphorylation of the activation loop occurs in the Asp-Phe-Gly (DFG) motif.^13,14 Thereby JAK2 inhibitors can be classified into two groups: the type I inhibitors conventionally bind exclusively the ATP-binding cleft of the active conformation (DFG-in), and the type II inhibitors are not only accommodated in the ATP-binding site but also can utilize an allosteric site adjacent to the ATP pocket of the inactive conformation (DFG-out).^15,16 To trace the history of the type I JAK2 inhibitors, the tyrphostin AG490 as the first type I JAK2 inhibitor has been developed for two decades.¹⁷ The progress of type I JAK2 inhibitors was advanced rapidly when the JAK2V617F mutation was discovered.¹⁸ These compounds belong to several different structural classes including pyrazines, pyrimidines, azaindoles, aminoindazoles, deazapurines, stilbenes, benzoxazoles and quinoxalines. Recently, a series of thienopyridine derivatives emerged to act as type I inhibitors with excellent bioactivities against JAK2-in (active).¹⁹ In contrast, it may be easier to achieve kinase selectivity for the type II inhibitors attributed to an additional and hydrophobic site with less conservation relative to the ATP-binding pocket.¹⁶ Despite some new developments in targeting the inactive JAK2 conformation (JAK2-out), the majority of efforts in type II inhibitors have failed to produce good leads.²⁰ To date, two type II inhibitors BBT594²¹ and CHZ868²² are characterized to target the JAK2-out (inactive) kinase. BBT594 is limited in potency and selectivity for JAK2 as well as pharmacokinetic properties while CHZ868 is effective in preclinical MPN and B-ALL (B cell acute lymphoblastic leukemia) models.²³ Based on the rational structure–activity relationship (SAR) design, the type II inhibitors can be regarded as composed of three portions, the type I head, the linker and the type II tail, to bind into the ATP-binding and allosteric pocket¹⁶ (Fig. 1). In order to provide a very promising entry to the type II inhibitor field, it is imperative to define the mechanisms of inhibition and further study the direct comparison between the type I and type II compounds.

Fig. 1 Three portions of the type II JAK2 inhibitors: the type I head, the linker and the type II tail. Chemical structures of BBT594 (a) and CHZ868 (b).

Herein, we investigated the individual and different mechanisms of the type I and type II inhibitors binding into active and inactive JAK2 conformations, respectively, by an atomic-level and comprehensive modeling protocol touching upon molecular docking, molecular dynamics simulations, free energy calculations and energy decompositions. The previous work has identified the highly JAK2 selective mechanism of the thienopyridine derivatives as type I inhibitors, of which the comprehensive modeling methods and the results for type I thienopyridine inhibitors to JAK2-in consistently match the present paper.²⁴ In addition, several new type II compounds are rationally designed on the basis of the structural and energetic results and their inhibitory potency were further evaluated by performing MD simulations combined with binding free energy analysis.

2. Materials and methods

2.1. Biological data set

As type I JAK2 inhibitors, a set of 20 thienopyridine derivatives with excellent inhibitory bioactivities against JAK2-in (active) were extracted on the basis of the work of Schenkel et al.¹⁹ The biological activity values IC₅₀ of all type I JAK2 inhibitors were turned into pIC₅₀ (−log IC₅₀) values for easier comparison, and the information of their structures and pIC₅₀ is listed in the ESI (Table 1S†). All 20 thienopyridine JAK2-in inhibitors with a wide range of IC₅₀ from 0.001 to 1.17 μM were considered to perform molecular docking. In order to gain a clear picture of energetic and structural features, inhibitors 22 and 25 with remarkable potency (IC₅₀ = 0.002 and 0.003 μM) and inhibitors 12 and 13 with poor activities (IC₅₀ = 1.17 and 1.17 μM) should be regarded as most suitable to run molecular dynamic simulations. In addition, currently only two potential type II JAK2 inhibitors BBT594 and CHZ868 (Fig. 1) with great activities against JAK2-out (inactive)²³ serve as references to investigate the binding characteristics by molecular docking and molecular dynamic simulations.

2.2. Molecular docking

2.2.1. Protein preparation. The two crystallographic structures of DFG-in JAK2 in complexation with thienopyridine inhibitor 19 (PDB entry: 3TJD; resolution, 2.9 Å)¹⁹ and DFG-out JAK2 in complexation with BBT594 (PDB entry: 3UGC; resolution, 1.34 Å)²¹ were obtained from the Protein Data Bank²⁵ (http://www.rcsb.org/pdb), which were perfected in the protein preparation wizard panel in Schrödinger 9.0 software package.²⁶ In view of the homodimer for 3TJD, the only one ligand–receptor subunit (B chain) was retained to prepare by Glide. It was operated on to remove all crystallographic water molecules, repair missing residues, add all hydrogen atoms and adjust the ligand bond orders and formal charges. With the purpose of reorienting the various groups and alleviating potential steric clashes for the OPLS2005 force field,²⁷ the hydrogen-bonding network was optimized, and then the minimization of the protein structure was restrained with converging heavy atoms to RMSD 0.30 Å.

2.2.2. Ligand preparation. All of type I compounds and two type II compounds were constructed in ChemDraw, and further processed and optimized by using the LigPrep panel in Maestro of Schrödinger.²⁸ The corresponding low-energy 3D structures were produced by generating variations on the ionization state, tautomers, stereochemistry and ring conformations.

2.2.3. Docking strategy. To predict the probable binding features and orientations of the inhibitors against DFG-in and DFG-out JAK2 kinase, molecular docking procedures were performed in Maestro of Schrödinger 9.0.²⁶ First, the mass center of the co-crystal ligand was confined as the centroid of the grid box in Glide, and the volume of the receptor in the active-site region was calculated by the default settings in the Receptor Grid Generation panel. As a result, coordinates of the ligands were respectively deleted from the two crystal structures and the biosystems without ligand were set. Then, the receptor rigid docking (RRD) method was used for the molecular docking protocol, in which the receptor was held rigid and the ligand was free to move. To decrease penalties for close contacts, the default was that the scaling factor for van der Waals radii was set to 0.8 with partial atomic charges less than or equal to 0.25. Finally, on the basis of the docking scores the conformations with the best recognition poses were obtained from molecular docking as initial structures of MD simulations.

2.3. Molecular dynamics simulations

2.3.1. Model systems. To identify the thorough binding mechanisms, MD simulations were performed in the Amber 11 software package.²⁹ According to the experimental bioactivities and molecular docking results, four docked complexes of DFG-in JAK2 with type I inhibitors 12, 13, 22 and 25, as well as two docked complexes of DFG-out JAK2 with type II inhibitors BBT594 and CHZ868, were appointed as the initial structures of MD simulations. The geometry optimization of the six ligands was carried out by the Gaussian 09 program at the HF/6-31G* level.³⁰ The general AMBER force field (GAFF)³¹ was loaded to generate the force field parameters for all ligands, and the HF/6-31G* RESP atomic charges were calculated with the antechamber module.^32–34 For the protein, the parameters could be founded by using the AMBER ff99SB force field.³⁵ Chlorine counter ions were added to neutralize the positive charges in LEaP, and the systems were solvated in a rectangular periodic box of TIP3P water molecules with a margin distance of 12 Å.²⁴

2.3.2. Simulation steps. Initially, the minimization of the solvent complexes was realized by three steps in the sander module: (i) all the protein atoms were restrained by a 10.0 kcal mol⁻¹ Å harmonic constraint, 500 steps of steepest descent minimization followed by 500 steps of conjugate gradient minimization; (ii) the protein backbone atoms were constrained with a force constant of 10.0 kcal mol⁻¹ Å by 1000 steps of steepest descent minimization and a subsequent 1000 steps of conjugate gradient minimization; (iii) all atoms were minimized through 2500 steps of steepest descent method and conjugate gradient method. Subsequently, in the NVT ensemble a gradual heating from 0.0 to 298.15 K for 50 ps was executed in a Langevin thermostat with a collision frequency of 2.0 ps⁻¹. Thereafter, the systems were converged in the NPT ensemble by five density equilibrations of 50 ps with a decreased restraint weight on the complex from 10.0, to 5.0, to 2.0, to 1.0, to 0.5, to 0.1, and 150 ps of equilibration without any position restraints. Finally, in the PMEMD module the production procedure (for 12 ns of four ligand–JAK2-in complexes and for 20 ns of two ligand–JAK2-out complexes) was performed in the NPT ensemble at 298.15 K and 1.0 atm pressure. During the MD simulations, the particle mesh Ewald (PME) method³⁶ with a 10.0 Å cutoff was used to calculate the long-range electrostatics, the bond lengths with hydrogen atoms were constrained by the SHAKE algorithm,³⁷ the time step of 2 fs and periodic boundary conditions (PBC) were set.

2.4. Energy calculations

2.4.1. Total binding free energy calculations. For each system, Molecular Mechanics/Generalized Born Surface Area (MM/GBSA) calculations³⁸ were executed to calculate the binding free energies in AMBER11. For four ligand–JAK2-in complexes, a number of 300 snapshots were extracted from the last stable 3 ns MD trajectories with an interval of 10 ps. For two ligand–JAK2-out complexes, a total number of 500 snapshots at time intervals of 10 ps were chosen from the last 5 ns MD trajectories with an interval of 10 ps. The following equations were used to calculate the binding free energy:³⁹

ΔG_bind = G_complex − G_receptor − G_ligand

ΔG = ΔG_MM + ΔG_sol − TΔS

where G_complex, G_receptor and G_ligand are the free energies of the complex, protein and inhibitor, respectively. It is clear that the binding free energy is equal to the sum of the molecular mechanics free energy (ΔE_gas), the solvation free energy (ΔG_sol) and the entropy term (TΔS).

ΔE_gas = ΔE_ele + ΔE_vdw

The molecular mechanics free energy (ΔE_gas) consists of the electrostatic (ΔE_ele) and van der Waals (ΔE_vdw) interactions.

ΔG_sol = ΔG_GB + ΔG_SA

where ΔG_GB and ΔG_SA are the polar and the nonpolar contributions to the solvation free energy (ΔG_sol). ΔG_GB was evaluated via the implicit solvation Generalized-Born (GB) model (igb = 5) and the dielectric constants of the solute and solvent were 1.0 and 80.0, respectively. ΔG_SA is a simple linear relation as follows:

ΔG_SA = γΔSASA + b

where SASA is the solvent-accessible surface (Ų), and the surface-tension (γ) and the offset (b) were set as 0.0072 kcal mol⁻¹ Å⁻² and 0.⁴⁰

The research mainly was not aimed at the calculations of absolute binding free energy but rather the qualitative analysis of binding mode. Due to the high computational demand and the low prediction accuracy, we overlooked the entropy contribution (−TΔS) in the study.

2.4.2. Free energy decomposition analysis for each residue. To make an incisive analysis for the binding site of the complexes, the binding free energies were decomposed into the effect of each residue by the MM/GBSA calculations, and subsequently the responsible residues were identified.

ΔG = ΔE_gas + ΔG_sol = ΔE_ele + ΔE_vdw + ΔG_GB + ΔG_SA

The total binding contribution for each inhibitor-residue is also defined as the electrostatic (ΔE_ele) and van der Waals (ΔE_vdw) energy terms, and the polar (ΔG_GB) and nonpolar (ΔG_SA) solvation free energy.

3. Results and discussion

3.1. Molecular docking

3.1.1. Validation of the docking protocol. Assessing if the accuracy of the molecular docking method is significant we can plot the docking scores as determined X-ray crystallography vs. the experimental activity.⁴¹ As shown in the linear relationship between docking scores and pIC₅₀ (Fig. 2) and in the residual plot (Fig. 3), the correlation coefficient (R = 0.78) was supported significantly by the residuals. The result suggests that the docking performance has the great capability to predict the bioactivity of the thienopyridine JAK2-in inhibitors. In addition, the root mean square deviation (RMSD) value between the best-scoring conformation and the X-ray crystallographic conformation (equal to or less than 2.0 Å)⁴² was calculated to evaluate the docking behavior in the present study. From Fig. 4a, there is a satisfactory superimposition between the JAK2-in crystal structure with compound 19 and two type I inhibitors 22 and 25, the RMSD values between them are equal to 0.580 and 0.889 Å, respectively. For type II inhibitors, the RMSD values between the JAK2-out crystal structure with compound BBT594 and the docked poses for BBT594 and CHZ868 are 1.179 and 1.446 Å, respectively (Fig. 5a). The overlay results indicate that the docking protocol is validated to bind the inhibitors into JAK2 kinase.

Fig. 2 The correlation between experimental activities and the docking scores for type I JAK2 inhibitors.

Fig. 3 The residual plot based on variables of docking score and pIC₅₀.

Fig. 4 (a) The superposition between the JAK2-in crystal structure with the inhibitor 19 (blue) and the re-docking conformations of JAK2-in with inhibitors 22 (green) and 25 (magenta); the docked structures of inhibitors (b) inhibitor 22 and (c) inhibitor 25 in the active site of JAK2-in.

Fig. 5 (a) The superposition between the JAK2-out crystal structure with the inhibitor BBT594 (blue) and the re-docking conformations of JAK2-out with inhibitors BBT594 (magenta) and CHZ868 (green); the docked structures of inhibitors (b) BBT594 and (c) CHZ868 in the active site of JAK2-out.

3.1.2. Initial orientation from molecular docking. It has been observed that the initial orientation of inhibitors can be determined in the binding pocket from molecular docking results. As shown in Fig. 4a, the scaffold of the thienopyridine JAK2-in inhibitors can be well embedded in the ATP-binding cleft and oriented so that the thienopyridine ring points towards the hinge region, the amide group is exposed to the bulk solvent and the tert-butyl group is located in the glycine loop. The docked alignments of inhibitors 22 and 25 in the active site reveal that compound 22 is fixed to JAK2-in by three hydrogen bonds with residues E930 and L932 (Fig. 4b) and compound 25 is stabilized to JAK2-in via four hydrogen bonds with residues E930, L932 and D994 (Fig. 4c). The residues, especially these in hinge region may play a important role in increasing biological potency. In contrast, the type II inhibitors can bind into not only the ATP site as found in DFG-in binding, but also the allosteric site for DFG-out binding.¹⁶ From Fig. 5a, the type I head points towards the hinge region, the type II tail is buried in the allosteric site and the linker is astride the DFG motif. Three hydrogen bonds between compound BBT594 and residues L932 and D994 are formed in the JAK2-out binding pocket (Fig. 5b). The compound CHZ868 is anchored in the ATP and allosteric site by four hydrogen bonds with residues E898, L932 and D994 (Fig. 5c). The results suggest that the key residues in the hinge region and in the additional site may serve as active players to JAK2-out inhibitor binding.

3.2. MD simulations

3.2.1. Stability of the simulation systems. MD simulations of 12 ns for four ligand–JAK2-in systems and that of 20 ns time duration for two ligand–JAK2-out systems were carried out to further obtain dynamic interactions. In order to examine the equilibration and stability of all systems, the temporal root mean square deviations (RMSDs) relative to the respective initial conformations were monitored, which are used to assess the structural dynamic changes during the MD simulations. The time evolution of RMSD values for all protein backbone atoms, the C_α atoms for the residues in the active site (residues within 5 Å around the ligand) and the ligand heavy atoms were plotted in Fig. 1S and 2S.† For all systems, the small RMSD fluctuations suggest that all simulation systems have reached equilibrium and stabilization during the MD simulations, which indicates that the binding conformations are reasonable to further analyze their binding energies. As shown in Fig. 1S,† the RMSD values of the four ligand–JAK2-in complexes were convergent after 7 ns, so that the last stable 3 ns trajectories were taken to the subsequent energy calculations. From Fig. 2S† the RMSDs against time for the two ligand–JAK2-out complexes reached the stable state during the last 10 ns, and the stable trajectories of the last 5 ns were employed to calculate their binding energies.

3.2.2. Structural and energetic characteristics of the DFG-in binding model. To further investigate the binding affinities of type I inhibitors and JAK2-in, the absolute binding free energy and free energy decomposition into individual components were calculated by using the MM/GBSA method in AMBER (Fig. 6). The predicted binding energies (ΔG_pred) of JAK2-in complexes with the studied type I inhibitors are listed in Table 1; what is noteworthy is that the ΔG_pred values of complexes with high activity inhibitors 22 and 25 (−57.27 and −55.77 kcal mol⁻¹, respectively) are greater than that of the complexes with poorly active inhibitors 12 and 13 (−45.93 and −47.99 kcal mol⁻¹, respectively). The results pinpoint that the ranking orders of the binding free energies are in good agreement with the experimental activity data. For the ligand–JAK2-in complexes the major driving force controlling the inhibitors binding into JAK2-in is the striking nonpolar contributions, which suggests that the most hydrophobic residues of the active pocket may give rise to van der Waals interactions to the inhibitor bindings. Conversely, the polar interactions with positive values fail to stabilize the ligand–JAK2-in complexes. In spite of the favorable electrostatic interactions (ΔE_ele), this is more than counteracted owing to the unfavorable polar solvation energy. Based on the comparison of highly and poorly active compound systems from four energy terms (ΔE_vdw, ΔE_ele, ΔG_GB and ΔG_SA), the significantly different binding affinities may be attributed to both the van der Waals force and the electrostatic force. More informatively, the comparisons between high and low activity inhibitors to JAK2-in are studied as follows:

Fig. 6 Ligand–residue interaction spectra for JAK2-in complexes with type I inhibitors (a) 22, (b) 13, (c) 25 and (d) 12.

Table 1 The predicted binding free energies and the individual energy terms for four complexes of DFG-in JAK2 with type I inhibitors based on MM/GBSA method (kcal mol⁻¹)

Contribution	22–JAK2-in	25–JAK2-in	12–JAK2-in	13–JAK2-in
ΔEvdw	−58.93 ± 2.71	−59.11 ± 2.83	−49.89 ± 2.84	−46.14 ± 3.14
ΔEele	−31.81 ± 4.93	−30.28 ± 5.74	−16.96 ± 3.43	−21.72 ± 3.13
ΔGGB	40.60 ± 3.96	41.10 ± 5.99	27.70 ± 3.23	26.32 ± 3.03
ΔGSA	−7.12 ± 0.19	−7.48 ± 0.21	−6.78 ± 0.20	−6.46 ± 0.23
ΔGMM	−90.74 ± 5.45	−89.39 ± 5.95	−66.85 ± 4.40	−67.86 ± 4.50
ΔGsol	33.48 ± 3.89	33.62 ± 5.98	20.92 ± 3.22	19.86 ± 3.02
ΔEpolar	8.79	10.82	10.73	4.61
ΔEnonpolar	−66.05	−66.59	−56.67	−52.60
ΔGpred	−57.27 ± 3.81	−55.77 ± 4.12	−45.93 ± 3.20	−47.99 ± 4.03
pIC50	8.70	8.52	5.93	5.93

---

Comparison of the structure–affinity relationship between high activity inhibitor 22 and low activity inhibitor 13 into JAK2-in. There are three structural differences for the two inhibitors that result in the inhibitory activity increasing from 1.17 μM (compound 13) to 0.003 μM (compound 22), the predicted binding free energy of compound 22 is lower than that of compound 13 by nearly 10 kcal mol⁻¹ (Table 1). The substituents for compound 13 at the R₁ and R₂ position are methyl groups and for compound 22 are a hydrogen and a propoxy. The 4-morpholinophenyl moiety of compound 13 at the R₃ position is replaced by the tert-butylsulfonamide group of compound 22. Concerning the bulkier tert-butylsulfonamide of the high activity inhibitor 22, it is buried fully in the active cavity, especially in the hinge region and glycine loop, and then can have more stable interactions with the key residues V863, M929, E930, Y931, L932, R980 and N981 than the low activity inhibitor 13 (Fig. 7b). As shown in Fig. 7a and c, a strong σ–π interaction is observed between the propoxy of compound 22 and the aromatic ring of Y931 (−2.41 kcal mol⁻¹). The hydrophilic residues R980 and N981 show stronger interaction with the polar sulfonamide of compound 22 than the 4-morpholinophenyl moiety of compound 13, and a stable nonpolar interaction between V863 and the tert-butyl of the R₃ position for compound 22 is generated (−2.55 kcal mol⁻¹). The greater polar contributions for the residues M929, E930 and L932 to compound 22 (−0.64, −2.25 and −1.99 kcal mol⁻¹, respectively) than compound 13 (0.02, −1.67 and −0.64 kcal mol⁻¹, respectively) correspond to the results of hydrogen bond analysis (Table 2). The complex 22–JAK2-in is anchored by a strong hydrogen bonding network that involved residues E930 (one H-bond, 22(N16–H38)⋯E930(O), with 76.80% occupancy), L932 (two H-bonds, L932(N–H)⋯22(N1) and 22(N18–H40)⋯L932(O), with 54.93 and 27.80% occupancy, respectively) and D994 (one H-bond, 22(N24–H43)⋯D994(OD1), with 24.07% occupancy). By contrast, there are only two weak hydrogen bonding interactions for the complex 13–JAK2-in (13(N10–H31)⋯E930(O) and L932(N–H)⋯13(N1) with 46.90 and 44.80% occupancy, respectively). In all, compound 22 has stronger polar and nonpolar contributions over compound 13 into the binding pocket of JAK2-in.

Fig. 7 (a) Comparison of the averaged structures for the 22–JAK2-in and 13–JAK2-in complexes (carbon atoms are colored in red and green, respectively); (b) energy difference of each residue to the binding of inhibitors 22 and 13; (c) the contributions of the individual energy terms for the key residues (ΔG_polar is marked in red and blue, and ΔG_nonpolar is marked in magenta and navy in 22–JAK2-in and 13–JAK2-in, respectively).

Table 2 Hydrogen bonds analysis for the complexes 22–JAK2-in and 13–JAK2-in from MD simulation

Complex	Donor	Acceptor	Distance (Å)	Angle (°)	Occupancy (%)
22–JAK2-in	22(N16–H38)	Glu930(O)	2.866	18.31	76.80
	Leu932(N–H)	22(N1)	2.915	15.58	54.93
	22(N18–H40)	Leu932(O)	2.895	39.91	27.80
	22(N24–H43)	Asp994(OD1)	2.896	22.31	24.07
13–JAK2-in	13(N10–H31)	Glu930(O)	2.870	20.62	46.90
	Leu932(N–H)	13(N1)	2.922	16.16	44.80

---

Comparison of the structure–affinity relationship between high activity inhibitor 25 and low activity inhibitor 12 into JAK2-in. Compound 25 in which a tert-butyl and a 4-morpholinophenyl of compound 12 at R₂ and R₃ positions, respectively, are replaced by a 4-morpholinoethyl and a bulkier tert-butylsulfonamide, is found to exhibit drastic inhibitory potential. The two inhibitors vary in bioactivity from 1.17 μM (compound 12) to 0.002 μM (compound 25), ascribed to the two structural differences, with the variation of the binding free energies for complexes 25–JAK2-in over 12–JAK2-in being nearly 10 kcal mol⁻¹. From the interactional structure and energy (Fig. 3S,†), the important residues S862, M929, E930, L932, P933, G935, R980 and N981 play a vital role in the binding difference between inhibitor 25 and inhibitor 12 against JAK2-in. The scaffold of compound 25 slightly moves to the hinge region mainly caused by the bulkier tert-butylsulfonamide on the R₃ position, which contributes to the advantage over compound 12 in terms of polar and nonpolar interactions. The numbers and occupancies of hydrogen bonds have a critical influence on electrostatic contribution (Table 2S†). The complex 25–JAK2-in shows four strong hydrogen bonds including to residues E930 (one H-bond, 25(N16–H42)⋯E930(O), with 73.73% occupancy), L932 (two H-bonds, L932(N–H)⋯25(N1) and 25(N18–H44)⋯L932(O), with 52.27 and 16.53% occupancy, respectively) and D994 (one H-bond, 25(N24–H47)⋯D994(OD1), with 59.13% occupancy). In contrast, compound 12 participated in two hydrogen bonds with the residues E930 (with 44.90% occupancy) and L932 (with 44.80% occupancy). In addition, the favorable hydrophobic contributions between the bulkier tert-butylsulfonamide of compound 25 and the important residues S862, R980 and N981 are −0.71, −1.61 and −0.86 kcal mol⁻¹, respectively, but they involve the unfavorable polar contributions with the 4-morpholinophenyl of compound 12 (0.08, 0.96 and 0.24 kcal mol⁻¹, respectively) indicating a bioactivity decrease. The significant difference for the residue G935 interaction mainly reflects in the adverse polar force (0.86 kcal mol⁻¹) because the hydrophilic amide group of compound 12 is nearer to the hydrophobic residue G935. Notably, the strong repulsive force (1.56 kcal mol⁻¹) attributed to the polar 4-morpholinoethyl of compound 25 at the R₂ position around the hydrophobic residue P933 acts against the stabilization of complex 25–JAK2-in (Fig. 3S†). Consequently, the similar behavior difference of the two comparisons suggests that the superior polar and nonpolar contributions for high activity compounds are mainly responsible due to the bulkier substituent at the R₃ position.

3.2.3. Structural and energetic characteristics of the DFG-out binding mode. To gain a thorough binding mechanism between type II inhibitors and JAK2-out, the MM/GBSA method was utilized to calculate the absolute binding free energy and free energy decomposition per residue. As listed in Table 3, the predicted binding energies (ΔG_pred) for JAK2-out complexes with BBT594 and CHZ868 (−65.54 and −67.66 kcal mol⁻¹, respectively) are qualitatively consistent with their experimental activity values. Apparently, the nonpolar contributions serve as a vital factor to stabilize the complexes arising from the most hydrophobic residues in the ATP and allosteric site. The polar interactions fail to be favorable to JAK2-out binding to inhibitors BBT594 and CHZ868. More importantly, the van der Waals interaction of complex BBT594–JAK2-out (−76.25 kcal mol⁻¹) is superior to that of complex CHZ868–JAK2-out (−64.43 kcal mol⁻¹); while the electrostatic force of complex BBT594–JAK2-out (−11.06 kcal mol⁻¹) is weaker than that of complex CHZ868–JAK2-out (−30.21 kcal mol⁻¹) in accord with the hydrogen bonding results (Table 4). For the two active type II inhibitors to JAK2-out each has its own merits so that both these favorable candidates should be further investigated to develop promising type II inhibitors.

Table 3 The predicted binding free energies and the individual energy terms for two complexes of DFG-out JAK2 with type II inhibitors based on MM/GBSA method (kcal mol⁻¹)

Contribution	BBT594–JAK2-out	CHZ868–JAK2-out
ΔEvdw	−76.25 ± 2.95	−64.43 ± 3.11
ΔEele	−11.06 ± 4.05	−30.21 ± 3.65
ΔGGB	31.18 ± 3.53	33.92 ± 2.08
ΔGSA	−9.39 ± 0.12	−6.94 ± 0.11
ΔGMM	−87.33 ± 4.88	−94.64 ± 3.72
ΔGsol	21.79 ± 3.52	26.98 ± 2.07
ΔEpolar	20.12	3.71
ΔEnonpolar	−58.64	−71.38
ΔGpred	−65.54 ± 4.22	−67.66 ± 3.41
IC50 (μM)	0.99	0.11

---

Table 4 Hydrogen bonding analysis for the complexes BBT594–JAK2-out and CHZ868–JAK2-out from MD simulation

Complex	Donor	Acceptor	Distance (Å)	Angle (°)	Occupancy (%)
BBT594–JAK2-out	Asp994(N–H)	BBT594(O21)	2.857	27.06	68.48
	Leu932(N–H)	BBT594(N2)	2.915	24.81	52.12
	BBT594(N7–H44)	Leu932(O)	2.901	32.23	27.12
CHZ868–JAK2-out	CHZ868(N23–H47)	Glu898(OE2)	2.838	13.84	62.08
	Asp994(N–H)	CHZ868(N19)	2.919	14.46	51.12
	Leu932(N–H)	CHZ868(N2)	2.917	24.85	45.28
	CHZ868(N7–H35)	Leu932(O)	2.917	24.85	26.76

---

As shown in Fig. 8, both BBT594 and CHZ868 can be well immersed in the ATP site and allosteric site of JAK2-out with drastically hydrophobic contributions. For complex BBT594–JAK2-out, these residues are significant to provide forceful nonpolar interactions connected with residues L902 of the αC-helix, V911, M929, Y931 and L932 of the hinge region, and L983, F995 of the activation loop. It is stabilized by three hydrogen bonds with L932 (two H-bonds, L932(N–H)⋯BBT594(N2) and BBT594(N7–H44)⋯L932(O), with 52.12 and 27.12% occupancy, respectively) and D994 (one H-bond, D994(N–H)⋯BBT594(O21), with 68.48% occupancy). Note, though, that the hydrophilic residue E898 shows an unfavorable polar interaction (3.52 kcal mol⁻¹) with the hydrophobic moiety of the type II tail. For complex CHZ868–JAK2-out, compound CHZ868 is associated with the key residues E898 and L902 of the αC-helix, V911, M929, Y931 and L932 of the hinge region, L983, G993, D994 and F995 of the activation loop. A stable network composed of four hydrogen bonds (Table 4) is formed via residues E898 (one H-bond, CHZ868(N23–H47)⋯E898(OE2), with 62.08% occupancy), L932 (two H-bonds, L932(N–H)⋯CHZ868(N2) and CHZ868(N7–H35)⋯L932(O), with 45.28 and 26.76% occupancy, respectively) and D994 (one H-bond, D994(N–H)⋯CHZ868 (N19), with 51.12% occupancy), which are in accord with the polar contributions of E898 (−2.56 kcal mol⁻¹), L932 (−0.19 kcal mol⁻¹) and D994 (−0.12 kcal mol⁻¹).

Fig. 8 The modes of the type II inhibitors binding to JAK2-out in the ATP and allosteric pocket (a) BBT594–JAK2-out and (d) CHZ868–JAK2-out; the corresponding ligand–residue interaction spectra for the complexes (b) BBT594–JAK2-out and (e) CHZ868–JAK2-out; the contributions of the individual energy terms for the key residues (ΔG_polar and ΔG_nonpolar are marked in black and red, respectively) (c) BBT594–JAK2-out and (f) CHZ868–JAK2-out.

In contrast, several interaction differences (Fig. 9) arise from the two type II inhibitors with different scaffolds (Fig. 1). First, the type I heads of the two compounds (the N-(pyrimidin-4-yl)acetamide moiety of BBT594 and the N-(pyridin-2-yl)acetamide moiety of CHZ868) are extremely similar so that similar interactions are observed between the hinge region and the type I head. Next, the linkers of type II inhibitor may exist in an orientation to leap over the DFG motif attributed to the switch of the ether and amide or amino at the end of the linkers, and the middle hydrophobic sections of the linkers impact on the binding DFG residues. The key residue D994 taking on an anchor role to fit the linker forms a hydrogen bond with both BBT594 and CHZ868. Conspicuous differences are observed at their type II tails, where the 1-methyl-4-(2-(trifluoromethyl)benzyl)piperazine moiety of BBT594 is substituted by the 1,3-difluorobenzene moiety of CHZ868. The acid residue E898 in the allosteric pocket can form a hydrogen bond with the amino at the linker of CHZ868, whereas it causes an opposite contribution with the (trifluoromethyl)benzyl group on the type II tail of BBT594. Compared with the small and short type II tail for CHZ868, the big and long tail for BBT594 can extend to residues I901 and L902 with greater nonpolar interactions. Regarding the powerful hydrophobic (trifluoromethyl)benzyl moiety for BBT594, it is practically completely screened by the adverse polar interactions, resulting in relatively lower interaction, although the van der Waals interactions with residues G993 and D994 provide a favorable force. Factors such as given above need to be considered when searching for new type II inhibitors.

Fig. 9 (a) Comparison of the averaged structures for BBT594–JAK2-out and CHZ868–JAK2-out (carbon atoms are colored in green and red, respectively); (b) energy difference of each residue to the binding of inhibitors BBT594 and CHZ868; (c) the contributions of the individual energy terms for the key residues (ΔG_polar is marked in red and blue, and ΔG_nonpolar is marked in magenta and navy in BBT594–JAK2-out and CHZ868–JAK2-out, respectively).

3.2.4. Comparison between DFG-in and DFG-out binding conformations. The special binding site for DFG-out conformation occurs not only the ATP site along with DFG-in binding but also the allosteric site for DFG-out binding.^14–16 Type I JAK2 kinase inhibitors, as ATP direct competitive inhibitors, can target the ATP binding site of the active conformation (DFG-in). The thienopyridine scaffold invariably occupies the adenine region and the substituent groups develop substantial contributions around the hydrophobic residues. They can normally create a 2–3 hydrogen bonding network with the residues of hinge region. Type II JAK2 kinase inhibitors are regarded as allosteric inhibitors because of the additional hydrophobic pocket adjacent to the ATP binding site in the inactive conformation (DFG-out). In order to perfectly bury into the unique binding pocket with the ATP site and the allosteric site, the type II JAK2 kinase inhibitors are composed of three components the type I head, the linker and the type II tail.¹⁶ The type I heads of BBT594 and CHZ868 have an analogous binding behavior with the type I JAK2 inhibitors that have a stable 2–3 hydrogen bonding network with the hinge residues and hydrophobic interactions in the ATP site. Their linkers are located on the DFG motif of the inactive conformation and have a hydrogen bond with residue D994 and hydrophobic interactions with the residues in the activation loop. The linker portion via throwing over the DFG motif of the inactive state aims at imbedding the type II tail into the allosteric site. The type II tails of BBT594 and CHZ868 attaching hydrophobic substitutions occupy the additional and hydrophobic site from which the stable nonpolar interactions form. In addition, the highly selective potency of the type II inhibitor is responsible for the unique allosteric pocket with a much less high degree of sequence and structural conservation. In all, the action mechanisms of DFG-in and DFG-out JAK2 inhibitors exist in the communication aspect as well as in the special aspect, and it is important to design promising type II JAK2 inhibitors by making full and rational use of insight into the association and characteristics. Type II JAK2 inhibitors are available to be installed drawing on the clear structure-based design strategy that the three portions, the type I head, the linker and the type II tail can perfectly occupy the ATP site and allosteric site by stable hydrogen bonding interactions and strong hydrophobic contributions.

3.3. Rational design of DFG-out kinase inhibitors

The above unequivocal comparisons provide a predominant basis for constructing novel type II JAK2 kinase inhibitors. The structure-based design strategy is that the head portions of BBT594 and CHZ868 can be preserved because of the excellent interactions with the ATP site, or the head portions are assigned to the scaffold of the type I JAK2 inhibitor due to the similar behavior with the hinge region, and that the linker portions are slightly modified or preserved to steer the orientation pointing to the tail binding site and have strong interactions in the DFG motif, and based on the unique features of residues in the allosteric site, the hydrophobic tail portions can be rationally improved to maximise hydrogen bonding interactions and hydrophobic contributions.^43–45 Gratifyingly, six potential type II JAK2 inhibitors optimized here are listed in Table 5. Their inhibitory abilities were predicted preliminarily relying on molecular docking. High docking scores were obtained from the docking results, which suggest that these new type II molecules may show satisfactory bioactivity. In order to further verify the binding potency and gain structural and energetic characteristics, molecules BC1 and C1 as ligands were assigned randomly to perform 20 ns of molecular dynamic simulations. From Fig. 4S,† the RMSD values of all systems tend to equilibration stability after 10 ns, and the last 5 ns trajectory can be used to analyze the binding energies in AMBER.

Table 5 The structures and corresponding predicted docking scores for new type II molecules

No.	Structure	Docking score
B1		−14.020
B2		−14.305
BC1		−13.112
BC2		−13.239
C1		−14.017
C2		−13.608

---

As can be seen in Table 6, the binding free energies predicted for complexes BC1–JAK2-out and C1–JAK2-out are greater (−69.18 and −69.58 kcal mol⁻¹, respectively) than that of complexes BBT594–JAK2-out and CHZ868–JAK2-out (−65.54 and −67.66 kcal mol⁻¹, respectively). There are nearly as great contributions for the van der Waals energy term as a driving force for the two new complexes as for complexes BBT594–JAK2-out and CHZ868–JAK2-out. Notably, it is gratifying that the electrostatic interactions of the two new compounds have been enhanced remarkably. To have a clearer understanding of the interactions and demonstrate the feasibility of the design concept, the binding energy was broken down into individual residue energy contributions. From Fig. 10, similar binding modes between the new DFG-out inhibitors and BBT594 and CHZ868 to JAK2-out are observed. It is obvious that the new inhibitors are anchored by a high number of hydrogen bonding interactions with high occupancies (Table 3S†), which account for the favorable electrostatic energy term. The head and linker portions fixed into the ATP binding cleft are substantially involved in these important hydrophobic and hydrogen bonding contributions with the key residues V911, M929, Y931 and L932 of the hinge region, and L983, G993, D994 and F995 of the activation loop. Also, importantly, the tail groups allow for the perfect cooperation with the unique allosteric site rather than reject the residue E898 bearing the hydrogen bonding acceptor. The strong polar contributions between E989 and BC1 and C1 (−3.56 and −4.95 kcal mol⁻¹, respectively) may be attributed to the stable hydrogen bonds (BC1(N24–H48)⋯E898(OE2) and BC1(N31–H53)⋯E898(OE2) with 28.52 and 51.72% occupancy, C1(N21–H40)⋯E898(OE2) with 74.55% occupancy). The binding affinities of new inhibitors to JAK2-out suggest that their bioactivity potency is proved by the combination of molecular docking and molecular dynamic simulations and is significant to further design promising DFG-out inhibitors.

Table 6 The predicted binding free energies and the individual energy terms between the new type II molecules and JAK2-out based on MM/GBSA method (kcal mol⁻¹)

Contribution	BC1–JAK2-out	C1–JAK2-out
ΔEvdw	−67.59 ± 3.01	−63.03 ± 2.86
ΔEele	−40.92 ± 4.45	−44.02 ± 4.24
ΔGGB	46.67 ± 2.91	44.40 ± 2.93
ΔGSA	−7.33 ± 0.11	−6.93 ± 0.11
ΔGMM	−108.52 ± 4.43	−107.05 ± 4.31
ΔGsol	39.33 ± 2.91	37.47 ± 2.93
ΔEpolar	5.74	0.38
ΔEnonpolar	−74.92	−69.96
ΔGpred	−69.18 ± 3.73	−69.58 ± 3.43

---

Fig. 10 The modes of the new type II molecules binding to JAK2-out in the ATP and allosteric pocket, the corresponding ligand–residue interaction spectra for the complexes and the contributions of the individual energy terms for the key residues (ΔG_polar and ΔG_nonpolar are marked in black and red, respectively). (a)–(c) BC1–JAK2-out and (d)–(f) C1–JAK2-out.

4. Conclusions

Here, we have investigated the structural mechanisms of type I and type II inhibitors binding to active and inactive JAK2 conformations, respectively, by atomic-level and comprehensive modeling protocol touching upon molecular docking, molecular dynamics simulations, free energy calculations and energy decompositions. In view of the comparative analysis between the thienopyridine JAK2-in inhibitors with the high and low activities, the binding affinities of type I inhibitors to JAK2-in are mainly responsible for hydrogen bonding interactions with residues E930 and L932 in the hinge region and strong hydrophobic contributions with the key residues V863, Y931 and R980 in the ATP-binding pocket. For type II JAK2 inhibitors BBT594 and CHZ868, not only can the type I head can form a stable hydrogen bonding network with the hinge residues E930 and L932 and strong hydrophobic interactions in the ATP-binding site, but the type II tail and the linker can be anchored in the unique allosteric site by hydrogen bonding interactions with residues E898 and D994, and substantial hydrophobic contributions with the important residues L902, L893 and F995 are observed in the allosteric site. Based on the above interaction information of ligands binding into JAK2, type II JAK2 inhibitors are reasonably designed and further evaluated by executing MD simulations combined with binding free energy calculation. Three portions consisting of the type I head, the linker and the type II tail of new type II JAK2 inhibitors have the expected interactions with JAK2- in terms of the ATP-binding and allosteric pocket. Therefore, the structural and energetic analysis can provide effective insight into the binding mechanism against JAK2-out to develop more promising type II JAK2 inhibitors.

Acknowledgements

The present work was supported by the National Natural Science Foundation of China (NSFC, No. 21275067).

References

1. A. G. Harpur, A. C. Andres, A. Ziemiecki, R. R. Aston and A. F. Wilks, Oncogene, 1992, 7, 1347–1353.
2. J. J. O’Shea, M. Gadina and R. D. Schreiber, Cell, 2002, 109, S121–S131.
3. J. J. O’Shea and R. Plenge, Immunity, 2012, 36, 542–550.
4. G. R. Stark and J. E. Darnell, Immunity, 2012, 36, 503–514.
5. A. Tefferi and T. Barbui, Mayo Clin. Proc., 2015, 90, 1283–1293.
6. C. James, V. Ugo, J. P. Le Couedic, J. Staerk, F. Delhommeau, C. Lacout, L. Garcon, H. Raslova, R. Berger, A. Bennaceur-Griscelli, J. L. Villeval, S. N. Constantinescu, N. Casadevall and W. Vainchenker, Nature, 2005, 434, 1144–1148.
7. R. L. Levine, M. Wadleigh, J. Cools, B. L. Ebert, G. Wernig, B. J. P. Huntly, T. J. Boggon, L. Wlodarska, J. J. Clark, S. Moore, J. Adelsperger, S. Koo, J. C. Lee, S. Gabriel, T. Mercher, A. D’Andrea, S. Frohling, K. Dohner, P. Marynen, P. Vandenberghe, R. A. Mesa, A. Tefferi, J. D. Griffin, M. J. Eck, W. R. Sellers, M. Meyerson, T. R. Golub, S. J. Lee and D. G. Gilliland, Cancer Cell, 2005, 7, 387–397.
8. P. Koppikar, N. Bhagwat, O. Kilpivaara, T. Manshouri, M. Adli, T. Hricik, F. Liu, L. M. Saunders, A. Mullally, O. Abdel-Wahab, L. Leung, A. Weinstein, S. Marubayashi, A. Goel, M. Gonen, Z. Estrov, B. L. Ebert, G. Chiosis, S. D. Nimer, B. E. Bernstein, S. Verstovsek and R. L. Levine, Nature, 2012, 489, 155–159.
9. E. J. Baxter, L. M. Scott, P. J. Campbell, C. East, N. Fourouclas, S. Swanton, G. S. Vassiliou, A. J. Bench, E. M. Boyd, N. Curtin, M. A. Scott, W. N. Erber and A. R. Green, Lancet, 2005, 365, 1054–1061.
10. K. Gnanasambandan, A. Magis and P. P. Sayeski, Biochemistry, 2010, 49, 9972–9984.
11. S. Ioannidis, M. L. Lamb, T. Wang, L. Almeida, M. H. Block, A. M. Davies, B. Peng, M. Su, H.-J. Zhang, E. Hoffmann, C. Rivard, I. Green, T. Howard, H. Pollard, J. Read, M. Alimzhanov, G. Bebernitz, K. Bell, M. Ye, D. Huszar and M. Zinda, J. Med. Chem., 2011, 54, 262–276.
12. C. Zhao, S. H. Yang, D. B. Khadka, Y. Jin, K. T. Lee and W. J. Cho, Bioorg. Med. Chem., 2015, 23, 985–995.
13. C. McInnes, M. Mezna and G. Kontopidis, Chem.Biol., 2006, 13, 693–694.
14. P. M. Ung and A. Schlessinger, ACS Chem. Biol., 2015, 10, 269–278.
15. Z. Zhao, H. Wu, L. Wang, Y. Liu, S. Knapp, Q. Liu and N. S. Gray, ACS Chem. Biol., 2014, 9, 1230–1241.
16. Y. Liu and N. S. Gray, Nat. Chem. Biol., 2006, 2, 358–364.
17. N. Meydan, T. Grunberger, H. Dadi, M. Shahar, E. Arpaia, Z. Lapidot, J. S. Leeder, M. Freedman, A. Cohen, A. Gazit, A. Levitzki and C. M. Roifman, Nature, 1996, 379, 645–648.
18. R. Baskin, A. Majumder and P. P. Sayeski, Curr. Med. Chem., 2010, 17, 4551–4558.
19. L. B. Schenkel, X. Huang, A. Cheng, H. L. Deak, E. Doherty, R. Emkey, Y. Gu, H. Gunaydin, J. L. Kim, J. Lee, R. Loberg, P. Olivieri, J. Pistillo, J. Tang, Q. Wan, H. L. Wang, S. W. Wang, M. C. Wells, B. Wu, V. Yu, L. Liu and S. Geuns-Meyer, J. Med. Chem., 2011, 54, 8440–8450.
20. O. Silvennoinen and S. R. Hubbard, Cancer Cell, 2015, 28, 1–2.
21. R. Andraos, Z. Qian, D. Bonenfant, J. Rubert, E. Vangrevelinghe, C. Scheufler, F. Marque, C. H. Regnier, A. De Pover, H. Ryckelynck, N. Bhagwat, P. Koppikar, A. Goel, L. Wyder, G. Tavares, F. Baffert, C. Pissot-Soldermann, P. W. Manley, C. Gaul, H. Voshol, R. L. Levine, W. R. Sellers, F. Hofmann and T. Radimerski, Cancer Discovery, 2012, 2, 512–523.
22. S. C. Meyer, M. D. Keller, S. Chiu, P. Koppikar, O. A. Guryanova, F. Rapaport, K. Xu, K. Manova, D. Pankov, R. J. O’Reilly, M. Kleppe, A. S. McKenney, A. H. Shih, K. Shank, J. Ahn, E. Papalexi, B. Spitzer, N. Socci, A. Viale, E. Mandon, N. Ebel, R. Andraos, J. Rubert, E. Dammassa, V. Romanet, A. Dolemeyer, M. Zender, M. Heinlein, R. Rampal, R. S. Weinberg, R. Hoffman, W. R. Sellers, F. Hofmann, M. Murakami, F. Baffert, C. Gaul, T. Radimerski and R. L. Levine, Cancer Cell, 2015, 28, 15–28.
23. S. C. Wu, L. S. Li, N. Kopp, J. Montero, B. Chapuy, A. Yoda, A. L. Christie, H. Liu, A. Christodoulou, D. van Bodegom, J. van der Zwet, J. V. Layer, T. Tivey, A. A. Lane, J. A. Ryan, S. Y. Ng, D. J. DeAngelo, R. M. Stone, D. Steensma, M. Wadleigh, M. Harris, E. Mandon, N. Ebel, R. Andraos, V. Romanet, A. Dolemeyer, D. Sterker, M. Zender, S. J. Rodig, M. Murakami, F. Hofmann, F. Kuo, M. J. Eck, L. B. Silverman, S. E. Sallan, A. Letai, F. Baffert, E. Vangrevelinghe, T. Radimerski, C. Gaul and D. M. Weinstock, Cancer Cell, 2015, 28, 29–41.
24. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, J. Chem. Phys., 1983, 79, 926–935.
25. H. M. Berman, J. Westbrook, Z. Feng, G. Gilliland, T. N. Bhat, H. Weissig, I. N. Shindyalov and P. E. Bourne, Nucleic Acids Res., 2000, 28, 235–242.
26. Schrödinger, version9.0, Schrödinger, LLC, New York, NY, 2009.
27. G. A. Kaminski, R. A. Friesner, J. Tirado-Rives and W. L. Jorgensen, J. Phys. Chem. B, 2001, 105, 6474–6487.
28. LigPrep, version9.0, Schrödinger, LLC, New York, NY, 2009.
29. D. Case, T. Darden, T. Cheatham III, C. Simmerling, J. Wang and R. Duke, AMBER 11, University of California, San Francisco, 2010.
30. M. Frisch, G. Trucks, H. Schlegel, G. Scuseria, M. Robb, J. Montgomery Jr, T. Vreven, K. Kudin and J. Burant, Gaussian 09, revision A. 02, Gaussian Inc., Wallingford, CT, 2009.
31. J. M. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman and D. A. Case, J. Comput. Chem., 2004, 25, 1157–1174.
32. T. Fox and P. A. Kollman, J. Phys. Chem. B, 1998, 102, 8070–8079.
33. P. Cieplak, W. D. Cornell, C. Bayly and P. A. Kollman, J. Comput. Chem., 1995, 16, 1357–1377.
34. C. I. Bayly, P. Cieplak, W. D. Cornell and P. A. Kollman, J. Phys. Chem., 1993, 97, 10269–10280.
35. V. Hornak, R. Abel, A. Okur, B. Strockbine, A. Roitberg and C. Simmerling, Proteins: Struct., Funct., Bioinf., 2006, 65, 712–725.
36. T. Darden, D. York and L. Pedersen, J. Chem. Phys., 1993, 98, 10089–10092.
37. J. P. Ryckaert, G. Ciccotti and H. J. C. Berendsen, J. Comput. Phys., 1977, 23, 327–341.
38. L. Xu, Y. Y. Li, L. Li, S. Y. Zhou and T. J. Hou, Mol. BioSyst., 2012, 8, 2260–2273.
39. H. Gouda, I. D. Kuntz, D. A. Case and P. A. Kollman, Biopolymers, 2003, 68, 16–34.
40. D. Sitkoff, K. A. Sharp and B. Honig, J. Phys. Chem., 1994, 98, 1978–1988.
41. R. A. Friesner, J. L. Banks, R. B. Murphy, T. A. Halgren, J. J. Klicic, D. T. Mainz, M. P. Repasky, E. H. Knoll, M. Shelley, J. K. Perry, D. E. Shaw, P. Francis and P. S. Shenkin, J. Med. Chem., 2004, 47, 1739–1749.
42. B. Kramer, M. Rarey and T. Lengauer, Proteins: Struct., Funct., Genet., 1999, 37, 228–241.
43. B. Okram, A. Nagle, F. J. Adrian, C. Lee, P. Ren, X. Wang, T. Sim, Y. Xie, X. Wang, G. Xia, G. Spraggon, M. Warmuth, Y. Liu and N. S. Gray, Chem. Biol., 2006, 13, 779–786.
44. S. Wenglowsky, D. Moreno, E. R. Laird, S. L. Gloor, L. Ren, T. Risom, J. Rudolph, H. L. Sturgis and W. C. Voegtli, Bioorg. Med. Chem. Lett., 2012, 22, 6237–6241.
45. A. Kohl, P. Amstutz, P. Parizek, H. K. Binz, C. Briand, G. Capitani, P. Forrer, A. Pluckthun and M. G. Grutter, Structure, 2005, 13, 1131–1141.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c6ra06266k

---