Insight into the inhibition mechanism and structure–activity relationship of 2,6-dipicolinic acid and its analogue to New Delhi metallo-β-lactamase-1

In previous literature, it was found that the activity of New Delhi Metallo-β-lactamase-1 (NDM-1) was inhibited by 2,6-dipicolinic acid (DPA) derivatives. To identify the mechanism of interaction between the inhibitors and NDM-1, molecular dynamics simulations were performed for the complex systems. Via the molecular modelling, inhibitors were found to be able bind to the region of catalytic activity of NDM-1. However, the detailed binding sites of the inhibitors differed with their structures. It was determined that His189, Lys211, Met248, Ser249, His250, and Ser251 are key residues for the binding of inhibitor 36 with NDM-1, and Asp124 is the only critical residue in the NDM-1-DPA complex. Furthermore, because of the interaction of the benzene ring in inhibitor 36 with the side chain of Lys211, inhibitor 36 can form 4 strong hydrogen bonds with protein. For the NDM-1-DPA complex, owing to the absence of the aniline group, DPA can only form a weak interaction with the residues around the binding site of NDM-1, except for Asp124, leading to a weaker inhibitory activity. Therefore, we believe that the strong interaction of the inhibitor with Lys211 results in effective inhibition, and the aniline group is the element required for the inhibitory activity.     

A comparative study of trypsin specificity based on QM/MM molecular dynamics simulation and QM/MM GBSA calculation

Hydrogen bonding and polar interactions play a key role in identification of protein-inhibitor binding specificity. Quantum mechanics/molecular mechanics molecular dynamics (QM/MM MD) simulations combined with DFT and semi-empirical Hamiltonian (AM1d, RM1, PM3, and PM6) methods were performed to study the hydrogen bonding and polar interactions of two inhibitors BEN and BEN1 with trypsin. The results show that the accuracy of treating the hydrogen bonding and polar interactions using QM/MM MD simulation of PM6 can reach the one obtained by the DFT QM/MM MD simulation. Quantum mechanics/molecular mechanics generalized Born surface area (QM/MM-GBSA) method was applied to calculate binding affinities of inhibitors to trypsin and the results suggest that the accuracy of binding affinity prediction can be significantly affected by the accurate treatment of the hydrogen bonding and polar interactions. In addition, the calculated results also reveal the binding specificity of trypsin: (1) the amidinium groups of two inhibitors generate favorable salt bridge interaction with Asp189 and form hydrogen bonding interactions with Ser190 and Gly214, (2) the phenyl of inhibitors can produce favorable van der Waals interactions with the residues His58, Cys191, Gln192, Trp211, Gly212, and Cys215. This systematic and comparative study can provide guidance for the choice of QM/MM MD methods and the designs of new potent inhibitors targeting trypsin.     

Intermolecular insights into allosteric inhibition of histone lysine-specific demethylase 1

BackgroundHistone lysine-specific demethylase 1 (LSD1) has become a potential anticancer target for the novel drug discovery. Recent reports have shown that SP2509 and its derivatives strongly inhibit LSD1 as allosteric inhibitors. However, the binding mechanism of these allosteric inhibitors in the allosteric site of LSD1 is not known yet.MethodsThe stability and binding mechanism of allosteric inhibitors in the binding site of LSD1 were evaluated by molecular docking, ligand-based pharmacophore, molecular dynamics (MD) simulations, molecular mechanics generalized born surface area (MM/GBSA) analysis, quantum mechanics/molecular mechanics (QM/MM) calculation and Hirshfeld surface analysis.ResultsThe conformational geometry and the intermolecular interactions of allosteric inhibitors showed high binding affinity towards allosteric site of LSD1 with the neighboring amino acids (Gly358, Cys360, Leu362, Asp375 and Glu379). Meanwhile, MD simulations and MM/GBSA analysis were performed on selected allosteric inhibitors in complex with LSD1 protein, which confirmed the high stability and binding affinity of these inhibitors in the allosteric site of LSD1.ConclusionThe simulation results revealed the crucial factors accounting for allosteric inhibitors of LSD1, including different protein–ligand interactions, the positions and conformations of key residues, and the ligands flexibilities. Meanwhile, a halogen bond interaction between chlorine atom of ligand and key residues Trp531 and His532 was recurrent in our analysis confirming its importance.General significanceOverall, our research analyzed in depth the binding modes of allosteric inhibitors with LSD1 and could provide useful information for the design of novel allosteric inhibitors.     

Understanding microscopic binding of macrophage migration inhibitory factor with phenolic hydrazones by molecular docking, molecular dynamics simulations and free energy calculations

Macrophage migration inhibitory factor (MIF), an immunoregulatory protein, is a potential target for a number of inflammatory diseases. In the current work, the interactions between MIF and a series of phenolic hydrazones were studied by molecular docking, molecular dynamics (MD) simulations, binding free energy calculations, and binding energy decomposition analysis to determine the structural requirement for achieving favorable biological activity of phenolic hydrazones. First, molecular docking was used to predict the binding modes of inhibitors in the binding site of MIF. The good correlation between the predicted docking scores and the experimental activities shows that the binding conformations of the inhibitors in the active site of MIF are well predicted. Moreover, our results suggest that the flexibility of MIF is essential in ligand binding process. Then, MD simulations and MM/GBSA free energy calculations were employed to determine the dynamic binding process and compare the binding modes of the inhibitors with different activities. The predicted binding free energies given by MM/GBSA are not well correlated with the experimental activities for the two subsets of the inhibitors; however, for each subset, a good correlation between the predicted binding free energies and the experimental activities is achieved. The MM/GBSA free energy decomposition analysis highlights the importance of hydrophobic residues for the MIF binding of the studied inhibitors. Based on the essential factors for MIF-inhibitor interactions derived from the theoretical predictions, some derivatives were designed and the higher inhibitory activities of several candidates were confirmed by molecular docking studies. The structural insights obtained from our study are useful for designing potent inhibitors of MIF.     

Computational studies and peptidomimetic design for the human p53-MDM2 complex

The interaction between human p53 and MDM2 is a key event in controlling cell growth. Many studies have suggested that a p53 mimic would be sufficient to inhibit MDM2 to reduce cell growth in cancerous tissue. In order to design a potent p53 mimic, molecular dynamics (MD) simulations were used to examine the binding interface and the effect of mutating key residues in the human p53-MDM2 complex. The Generalized Born surface area (GBSA) method was used to estimate free energies of binding, and a computational alanine-scanning approach was used to calculate the relative effects in the free energy of binding for key mutations. Our calculations determine the free energy of binding for a model p53-MDM2 complex to be -7.4 kcal/mol, which is in very good agreement with the experimentally determined values (-6.6--8.8 kcal/mol). The alanine-scanning results are in good agreement with experimental data and calculations by other groups. We have used the information from our studies of human p53-MDM2 to design a beta-peptide mimic of p53. MD simulations of the mimic bound to MDM2 estimate a free energy of binding of -8.8 kcal/mol. We have also applied alanine scanning to the mimic-MDM2 complex and reveal which mutations are most likely to alter the binding affinity, possibly giving rise to escape mutants. The mimic was compared to nutlins, a new class of inhibitors that block the formation of the p53-MDM2 complex. There are interesting similarities between the nutlins and our mimic, and the differences point to ways that both inhibitors may be improved. Finally, an additional hydrophobic pocket is noted in the interior of MDM2. It may be possible to design new inhibitors to take advantage of that pocket.     

Determinants of reactivity and selectivity in soluble epoxide hydrolase from quantum mechanics/molecular mechanics modeling

Soluble epoxide hydrolase (sEH) is an enzyme involved in drug metabolism that catalyzes the hydrolysis of epoxides to form their corresponding diols. sEH has a broad substrate range and shows high regio- and enantioselectivity for nucleophilic ring opening by Asp333. Epoxide hydrolases therefore have potential synthetic applications. We have used combined quantum mechanics/molecular mechanics (QM/MM) umbrella sampling molecular dynamics (MD) simulations (at the AM1/CHARMM22 level) and high-level ab initio (SCS-MP2) QM/MM calculations to analyze the reactions, and determinants of selectivity, for two substrates: trans-stilbene oxide (t-SO) and trans-diphenylpropene oxide (t-DPPO). The calculated free energy barriers from the QM/MM (AM1/CHARMM22) umbrella sampling MD simulations show a lower barrier for phenyl attack in t-DPPO, compared with that for benzylic attack, in agreement with experiment. Activation barriers in agreement with experimental rate constants are obtained only with the highest level of QM theory (SCS-MP2) used. Our results show that the selectivity of the ring-opening reaction is influenced by several factors, including proximity to the nucleophile, electronic stabilization of the transition state, and hydrogen bonding to two active site tyrosine residues. The protonation state of His523 during nucleophilic attack has also been investigated, and our results show that the protonated form is most consistent with experimental findings. The work presented here illustrates how determinants of selectivity can be identified from QM/MM simulations. These insights may also provide useful information for the design of novel catalysts for use in the synthesis of enantiopure compounds.     

Discovery of potential sclerostin inhibitors from plants with loop2 region of sclerostin inhibition by interacting with residues outside Pro-Asn-Ala-Ile-Gly motif

Abstract

Sclerostin, an antagonist of the Wnt/β-catenin signaling pathway, was discovered as a potential therapeutic target for stimulating bone formation in osteoporosis. In this study, molecular docking was employed to predict the binding of 29 herbal compounds, which were reported as bone formation stimulators, to the loop2 region of sclerostin. Then, the 50 ns molecular dynamics (MD) simulation of the complexes between sclerostin and the top 10 hits obtained from molecular docking were carried out. Root mean square deviations (RMSDs) analysis of MD trajectories pointed out that all ligands-complexes remain stable throughout the duration of MD simulations. In addition, the molecular mechanics/generalized born surface area (MM/GBSA) binding free energy and energy decomposition analyses were determined. The results here suggested that baicalin is the most promising inhibitor of sclerostin. Interestingly, baicalin binds to sclerostin via the hydrophobic interaction with the amino acid residues on loop2 region but outside the Pro-Asn-Ala-Ile-Gly (PNAIG) motif, particularly the Arg-Gly-Lys-Trp-Trp-Arg (RGKWWR) motif. This finding could be a novel strategy for developing new sclerostin inhibitors in the future.

Communicated by Ramaswamy H. Sarma     

Binding free energies for nicotine analogs inhibiting cytochrome P450 2A6 by a combined use of molecular dynamics simulations and QM/MM-PBSA calculations

Molecular dynamics (MD) simulations and hybrid quantum mechanical/molecular mechanical (QM/MM) calculations have been performed to explore the dynamic behaviors of cytochrome P450 2A6 (CYP2A6) binding with nicotine analogs (that are typical inhibitors) and to calculate their binding free energies in combination with Poisson–Boltzmann surface area (PBSA) calculations. The combined MD simulations and QM/MM-PBSA calculations reveal that the most important structural parameters affecting the CYP2A6-inhibitor binding affinity are two crucial internuclear distances, that is, the distance between the heme iron atom of CYP2A6 and the coordinating atom of the inhibitor, and the hydrogen-bonding distance between the N297 side chain of CYP2A6 and the pyridine nitrogen of the inhibitor. The combined MD simulations and QM/MM-PBSA calculations have led to dynamic CYP2A6-inhibitor binding structures that are consistent with the observed dynamic behaviors and structural features of CYP2A6-inhibitor binding, and led to the binding free energies that are in good agreement with the experimentally-derived binding free energies. The agreement between the calculated binding free energies and the experimentally-derived binding free energies suggests that the combined MD and QM/MM-PBSA approach may be used as a valuable tool to accurately predict the CYP2A6-inhibitor binding affinities in future computational design of new, potent and selective CYP2A6 inhibitors.     

Lys169 of Human Glucokinase Is a Determinant for Glucose Phosphorylation: Implication for the Atomic Mechanism of Glucokinase Catalysis

Glucokinase (GK), a glucose sensor, maintains plasma glucose homeostasis via phosphorylation of glucose and is a potential therapeutic target for treating maturity-onset diabetes of the young (MODY) and persistent hyperinsulinemic hypoglycemia of infancy (PHHI). To characterize the catalytic mechanism of glucose phosphorylation by GK, we combined molecular modeling, molecular dynamics (MD) simulations, quantum mechanics/molecular mechanics (QM/MM) calculations, experimental mutagenesis and enzymatic kinetic analysis on both wild-type and mutated GK. Our three-dimensional (3D) model of the GK-Mg²⁺-ATP-glucose (GMAG) complex, is in agreement with a large number of mutagenesis data, and elucidates atomic information of the catalytic site in GK for glucose phosphorylation. A 10-ns MD simulation of the GMAG complex revealed that Lys169 plays a dominant role in glucose phosphorylation. This prediction was verified by experimental mutagenesis of GK (K169A) and enzymatic kinetic analyses of glucose phosphorylation. QM/MM calculations were further used to study the role of Lys169 in the catalytic mechanism of the glucose phosphorylation and we found that Lys169 enhances the binding of GK with both ATP and glucose by serving as a bridge between ATP and glucose. More importantly, Lys169 directly participates in the glucose phosphorylation as a general acid catalyst. Our findings provide mechanistic details of glucose phorphorylation catalyzed by GK, and are important for understanding the pathogenic mechanism of MODY.     

Probing the structural and energetic basis of kinesin-microtubule binding using computational alanine-scanning mutagenesis

Kinesin-microtubule (MT) binding plays a critical role in facilitating and regulating the motor function of kinesins. To obtain a detailed structural and energetic picture of kinesin-MT binding, we performed large-scale computational alanine-scanning mutagenesis based on long-time molecular dynamics (MD) simulations of the kinesin-MT complex in both ADP and ATP states. First, we built three all-atom kinesin-MT models: human conventional kinesin bound to ADP and mouse KIF1A bound to ADP and ATP. Then, we performed 30 ns MD simulations followed by kinesin-MT binding free energy calculations for both the wild type and mutants obtained after substitution of each charged residue of kinesin with alanine. We found that the kinesin-MT binding free energy is dominated by van der Waals interactions and further enhanced by electrostatic interactions. The calculated mutational changes in kinesin-MT binding free energy are in excellent agreement with results of an experimental alanine-scanning study with a root-mean-square error of ~0.32 kcal/mol [Woehlke, G., et al. (1997) Cell 90, 207-216]. We identified a set of important charged residues involved in the tuning of kinesin-MT binding, which are clustered on several secondary structural elements of kinesin (including well-studied loops L7, L8, L11, and L12, helices α4, α5, and α6, and less-explored loop L2). In particular, we found several key residues that make different contributions to kinesin-MT binding in ADP and ATP states. The mutations of these residues are predicted to fine-tune the motility of kinesin by modulating the conformational transition between the ADP state and the ATP state of kinesin.     

Probing the interaction mechanism of small molecule inhibitors with matriptase based on molecular dynamics simulation and free energy calculations

Matriptase is a serine protease associated with a wide variety of human tumors and carcinoma progression. Up to now, many promising anti-cancer drugs have been developed. However, the detailed structure–function relationship between inhibitors and matriptase remains elusive. In this work, molecular dynamics simulation and binding free energy studies were performed to investigate the biochemistry behaviors of two class inhibitors binding to matriptase. The binding free energies predicted by MM/GBSA methods are in good agreement with the experimental bioactivities, and the analysis of the individual energy terms suggests that the van der Waals interaction is the major driving force for ligand binding. The key residues responsible for achieving strong binding have been identified by the MM/GBSA free energy decomposition analysis. Especially, Trp215 and Phe99 had an important impact on active site architecture and ligand binding. The results clearly identify the two class inhibitors exist different binding modes. Through summarizing the two different modes, we have mastered some important and favorable interaction patterns between matriptase and inhibitors. Our findings would be helpful for understanding the interaction mechanism between the inhibitor and matriptase and afford important guidance for the rational design of potent matriptase inhibitors.     

3D-QSAR and molecular recognition of Klebsiella pneumoniae NDM-1 inhibitors

New Delhi metallo-β-lactamase-1 (NDM-1) as a target for the development of anti-superbug agents, plays an important role in the resistance of β-lactam antibiotics and has received worldwide attention. Sulfhydryl propionic acid derivatives can effectively inhibit the catalytic activity of NDM-1, but the quantitative structure–activity relationship (QSAR) and inhibitor-target recognition mechanism both remain unclear. In this work, CoMFA and CoMSIA models of sulfhydryl propionic acid inhibitors with high predictive ability were obtained, from which the effect of different substituents on the inhibitory activity against NDM-1 were revealed at the molecular level. Then, two 120-nanosecond comparative molecular dynamics (MD) simulations for NDM-1 enzyme and NDM-1-inhibitor complex systems were performed to study the recognition and inhibition mechanism of sulfhydryl propionic acid derivatives. The binding of inhibitors alters the entire H-bond network of the NDM-1 system accompanied by the formation of strong interactions with I35, W93, H120, H122, D124, H189 and H250, further weakens the recognition of NDM-1 by its endogenic substrates. Finally, the results of free energy landscape and conformation cluster analyses show that NDM-1 underwent a significant conformational change after the association with sulfhydryl propionic acid inhibitors. Our findings can provide theoretical support and help for anti-superbug agents design based on the structures of NDM-1 and sulfhydryl propionic acid derivatives.     

Thermodynamics calculation of protein–ligand interactions by QM/MM polarizable charge parameters

The calculation of protein–ligand binding free energy (ΔG) is of great importance for virtual screening and drug design. Molecular dynamics (MD) simulation has been an attractive tool to investigate this scientific problem. However, the reliability of such approach is affected by many factors including electrostatic interaction calculation. Here, we present a practical protocol using quantum mechanics/molecular mechanics (QM/MM) calculations to generate polarizable QM protein charge (QMPC). The calculated QMPC of some atoms in binding pockets was obviously different from that calculated by AMBER ff03, which might significantly affect the calculated ΔG. To evaluate the effect, the MD simulations and MM/GBSA calculation with QMPC for 10 protein–ligand complexes, and the simulation results were then compared to those with the AMBER ff03 force field and experimental results. The correlation coefficient between the calculated ΔΔG using MM/GBSA under QMPC and the experimental data is .92, while that with AMBER ff03 force field is .47 for the complexes formed by streptavidin or its mutants and biotin. Moreover, the calculated ΔΔG with QMPC for the complexes formed by ERβ and five ligands is positively related to experimental result with correlation coefficient of .61, while that with AMBER ff03 charge is negatively related to experimental data with correlation coefficient of .42. The detailed analysis shows that the electrostatic polarization introduced by QMPC affects the electrostatic contribution to the binding affinity and thus, leads to better correlation with experimental data. Therefore, this approach should be useful to virtual screening and drug design.     

Molecular modeling of green fluorescent protein: structural effects of chromophore deprotonation

Molecular dynamics (MD) simulations were carried out to study the conformational rearrangement induced by deprotonation of the fluorescent chromophore in GFP, as well as the associated changes in the hydrogen-bonding network. For both the structures with either a neutral or an anionic chromophore, it was found that the beta-barrel was stable and rigid, and the conformation of the chromophore was consistent with the available x-ray structure. The conformational change in Thr203 due to deprotonation was also found to be consistent with the three-state isomerization model. Although GFP is highly fluorescent, denatured-GFP is nonfluorescent, indicating that the environment of the protein plays an important role in its fluorescence behavior. Our MD simulations, which explore the effect of the protein shell on the conformation of the chromophore, find the flexibility of the central chromophore to be significantly restricted due to the rigid nature of the protein shell. The hydrogen-bonding between the chromophore and neighboring residues was also shown to contribute to the chromophore rigidity. In addition to the MD studies, quantum mechanics/molecular mechanics (QM/MM) ONIOM calculations were carried out to investigate the effect of the beta-barrel on the internal rotation in the chromophore. Along with providing quantitative values for torsional rotation barriers about the bridging bond in the chromophore, the ONIOM calculations also validate our MD force field parameters.     

Dynamics and Binding Modes of Free cdk2 and its Two Complexes with Inhibitors Studied by Computer Simulations

Abstract

This article presents a molecular dynamics (MD) study of the cdk2 enzyme and its two complexes with the inhibitors isopentenyladenine and roscovitine using the Cornell et al. force field from the AMBER software package. The results show that inserting an inhibitor into the enzyme active site does not considerably change enzyme structure but it seemingly changes the distribution of internal motions. The inhibitor causes differences in the domain motions in free cdk2 and in its complexes. It was found out that repulsion of roscovitine N9 substituent causes conformational change on Lys 33 side chain. Isopentenyladenine forms with Lys 33 side chain terminal amino group a hydrogen bond. It implies that the cavity, where N9 substituent of roscovitine is buried, can adopt larger substituent due to Lys 33 side chain flexibility. The composition of electrostatic and van der Waals interactions between the inhibitor and the enzyme were also calculated along both cdk2/inhibitor MD trajectories together with MM-PB/GBSA analysis. These results show that isopentenyladenine-like inhibitors could be more effective after modifications leading to an increase in their van der Waals contact with the enzyme. We suggest that a way leading to better inhibitors occupying isopentenyladenine binding mode could be: to keep N9 and N7 purine positions free, to keep 3,3-dimethylallylamino group at C6 position, and to add, e.g., benzylamino group at C2 position. The results support the idea that the isopentenyladenine binding mode can be used for cdk2 inhibitors design and that all possibilities to improve this binding mode were not uncovered yet.     

Influence of Chirality of Crizotinib on Its MTH1 Protein Inhibitory Activity: Insight from Molecular Dynamics Simulations and Binding Free Energy Calculations

As a promising target for the treatment of lung cancer, the MutT Homolog 1 (MTH1) protein can be inhibited by crizotinib. A recent work shows that the inhibitory potency of (S)-crizotinib against MTH1 is about 20 times over that of (R)-crizotinib. But the detailed molecular mechanism remains unclear. In this study, molecular dynamics (MD) simulations and free energy calculations were used to elucidate the mechanism about the effect of chirality of crizotinib on the inhibitory activity against MTH1. The binding free energy of (S)-crizotinib predicted by the Molecular Mechanics/Generalized Born Surface Area (MM/GBSA) and Adaptive biasing force (ABF) methodologies is much lower than that of (R)-crizotinib, which is consistent with the experimental data. The analysis of the individual energy terms suggests that the van der Waals interactions are important for distinguishing the binding of (S)-crizotinib and (R)-crizotinib. The binding free energy decomposition analysis illustrated that residues Tyr7, Phe27, Phe72 and Trp117 were important for the selective binding of (S)-crizotinib to MTH1. The adaptive biasing force (ABF) method was further employed to elucidate the unbinding process of (S)-crizotinib and (R)-crizotinib from the binding pocket of MTH1. ABF simulation results suggest that the reaction coordinates of the (S)-crizotinib from the binding pocket is different from (R)-crizotinib. The results from our study can reveal the details about the effect of chirality on the inhibition activity of crizotinib to MTH1 and provide valuable information for the design of more potent inhibitors.