A molecular dynamics investigation of CDK8/CycC and ligand binding: conformational flexibility and implication in drug discovery

Abnormal activity of cyclin-dependent kinase 8 (CDK8) along with its partner protein cyclin C (CycC) is a common feature of many diseases including colorectal cancer. Using molecular dynamics (MD) simulations, this study determined the dynamics of the CDK8-CycC system and we obtained detailed breakdowns of binding energy contributions for four type-I and five type-II CDK8 inhibitors. We revealed system motions and conformational changes that will affect ligand binding, confirmed the essentialness of CycC for inclusion in future computational studies, and provide guidance in development of CDK8 binders. We employed unbiased all-atom MD simulations for 500 ns on twelve CDK8-CycC systems, including apoproteins and protein–ligand complexes, then performed principal component analysis (PCA) and measured the RMSF of key regions to identify protein dynamics. Binding pocket volume analysis identified conformational changes that accompany ligand binding. Next, H-bond analysis, residue-wise interaction calculations, and MM/PBSA were performed to characterize protein–ligand interactions and find the binding energy. We discovered that CycC is vital for maintaining a proper conformation of CDK8 to facilitate ligand binding and that the system exhibits motion that should be carefully considered in future computational work. Surprisingly, we found that motion of the activation loop did not affect ligand binding. Type-I and type-II ligand binding is driven by van der Waals interactions, but electrostatic energy and entropic penalties affect type-II binding as well. Binding of both ligand types affects protein flexibility. Based on this we provide suggestions for development of tighter-binding CDK8 inhibitors and offer insight that can aid future computational studies.     

(Thermo)dynamic Role of Receptor Flexibility,Entropy, and Motional Correlation in Protein–Ligand Binding

The binding of 2‐amino‐5‐methylthiazole to the W191G cavity mutant of cytochrome c peroxidase is an ideal test case to investigate the entropic contribution to the binding free energy due to changes in receptor flexibility. The dynamic and thermodynamic role of receptor flexibility are studied by 50 ns‐long explicit‐solvent molecular dynamics simulations of three separate receptor ensembles: W191G binding a K⁺ ion, W191G–2a5mt complex with a closed 190–195 gating loop, and apo with an open loop. We employ a method recently proposed to estimate accurate absolute single‐molecule configurational entropies and their differences for systems undergoing conformational transitions. We find that receptor flexibility plays a generally underestimated role in protein–ligand binding (thermo)dynamics and that changes of receptor motional correlation determine such large entropy contributions.     

Labeling substrates of protein arginine methyltransferase with engineered enzymes and matched S-adenosyl-L-methionine analogues

Elucidating physiological and pathogenic functions of protein methyltransferases (PMTs) relies on knowing their substrate profiles. S-adenosyl-L-methionine (SAM) is the sole methyl-donor cofactor of PMTs. Recently, SAM analogues have emerged as novel small-molecule tools to efficiently label PMT substrates. Here we reported the development of a clickable SAM analogue cofactor, 4-propargyloxy-but-2-enyl SAM, and its implementation to label substrates of human protein arginine methyltransferase 1 (PRMT1). In the system, the SAM analogue cofactor, coupled with matched PRMT1 mutants rather than native PRMT1, was shown to label PRMT1 substrates. The transferable 4-propargyloxy-but-2-enyl moiety of the SAM analogue further allowed corresponding modified substrates to be characterized through a subsequent click chemical ligation with an azido-based probe. The SAM analogue, in combination with a rational protein-engineering approach, thus shows potential to label and identify PMT targets in the context of a complex cellular mixture.     

Effect of the R119G mutation on human P5CR structure and its interactions with NAD: Insights derived from molecular dynamics simulation and free energy analysis

Pyrroline-5-carboxylate reductase (P5CR), an enzyme with conserved housekeeping roles, is involved in the etiology of cutis laxa. While previous work has shown that the R119G point mutation in the P5CR protein is involved, the structural mechanism behind the pathology remains to be elucidated. In order to probe the role of the R119G mutation in cutis laxa, we performed molecular dynamics (MD) simulations, essential dynamics (ED) analysis, and Molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) binding free energy calculations on wild type (WT) and mutant P5CR-NAD complex. These MD simulations and ED analyses suggest that the R119G mutation decreases the flexibility of P5CR, specifically in the substrate binding pocket, which could decrease the kinetics of the cofactor entrance and egress. Furthermore, the MM-PBSA calculations suggest the R119G mutant has a lower cofactor binding affinity for NAD than WT. Our study provides insight into the possible role of the R119G mutation during interactions between P5CR and NAD, thus bettering our understanding of how the mutation promotes cutis laxa.     

A molecular dynamics and computational study of human KAT3 involved in KYN pathway

Kynurenine aminotransferases (KATs) catalyze the transamination of kynurenine (KYN) pathway and endogenous KYNs have been suggested to highly correlate to abnormal brain diseases. HKAT3 is a key member of KAT family, while the binding mechanism of KYN and cofactor with HKAT3 has not been determined yet. In this study, we focus on the structure-function relationship among KYN, cofactor and HKAT3. The binding models of KYN complex and KYN&cofactor complex were obtained and were studied by molecular dynamics (MD) simulations. We identified several critical residues and influence of conformational changes in human kynurenine aminotransferase 3 (HKAT3) complexes. The cofactor may contribute largely not only to the catalysis, but also to the binding. In addition, a hypothesis is proposed that a strong hydrophobic interaction between Tyr159 and Lys280 may influence the binding mode and the binding region of the substrate and the cofactor. Our results will be a good starting point for further determination of the biological role.     

Prediction of protein–protein complexes using replica exchange with repulsive scaling

The realistic prediction of protein–protein complex structures is import to ultimately model the interaction of all proteins in a cell and for the design of new protein–protein interactions. In principle, molecular dynamics (MD) simulations allow one to follow the association process under realistic conditions including full partner flexibility and surrounding solvent. However, due to the many local binding energy minima at the surface of protein partners, MD simulations are frequently trapped for long times in transient association states. We have designed a replica-exchange based scheme employing different levels of a repulsive biasing between partners in each replica simulation. The bias acts only on intermolecular interactions based on an increase in effective pairwise van der Waals radii (repulsive scaling (RS)-REMD) without affecting interactions within each protein or with the solvent. For a set of five protein test cases (out of six) the RS-REMD technique allowed the sampling of near-native complex structures even when starting from the opposide site with respect to the native binding site for one partner. Using the same start structures and same computational demand regular MD simulations sampled near native complex structures only for one case. The method showed also improved results for the refinement of docked structures in the vicinity of the native binding geometry compared to regular MD refinement.     

Ligand Binding and Release Investigated by Contact-Guided Iterative Multiple Independent Molecular Dynamics Simulations

Binding and releasing ligands are critical for the biological functions of many proteins, so it is important to determine these highly dynamic processes. Although there are experimental techniques to determine the structure of a protein-ligand complex, it only provides a static picture of the system. With the rapid increase of computing power and improved algorithms, molecular dynamics (MD) simulations have diverse of superiority in probing the binding and release process. However, it remains a great challenge to overcome the time and length scales when the system becomes large. This work presents an enhanced sampling tool for ligand binding and release, which is based on iterative multiple independent MD simulations guided by contacts formed between the ligand and the protein. From the simulation results on adenylate kinase, we observe the process of ligand binding and release while the conventional MD simulations at the same time scale cannot.     

Calculation of binding free energies of inhibitors to plasmepsin II

An understanding at the atomic level of the driving forces of inhibitor binding to the protein plasmepsin (PM) II would be of interest to the development of drugs against malaria. To this end, three state of the art computational techniques to compute relative free energies-thermodynamic integration (TI), Hamiltonian replica-exchange (H-RE) TI, and comparison of bound versus unbound ligand energy and entropy-were applied to a protein-ligand system of PM II and several exo-3-amino-7-azabicyclo[2.2.1]heptanes and the resulting relative free energies were compared with values derived from experimental IC(50) values. For this large and flexible protein-ligand system, the simulations could not properly sample the relevant parts of the conformational space of the bound ligand, resulting in failure to reproduce the experimental data. Yet, the use of Hamiltonian replica exchange in conjunction with thermodynamic integration resulted in enhanced convergence and computational efficiency compared to standard thermodynamic integration calculations. The more approximate method of calculating only energetic and entropic contributions of the ligand in its bound and unbound states from conventional molecular dynamics (MD) simulations reproduced the major trends in the experimental binding free energies, which could be rationalized in terms of energetic and entropic characteristics of the different structural and physico-chemical properties of the protein and ligands.     

Insights into Conformational Dynamics and Allostery in DNMT1-H3Ub/USP7 Interactions

DNA methyltransferases (DNMTs) including DNMT1 are a conserved family of cytosine methylases that play crucial roles in epigenetic regulation. The versatile functions of DNMT1 rely on allosteric networks between its different interacting partners, emerging as novel therapeutic targets. In this work, based on the modeling structures of DNMT1-ubiquitylated H3 (H3Ub)/ubiquitin specific peptidase 7 (USP7) complexes, we have used a combination of elastic network models, molecular dynamics simulations, structural residue perturbation, network modeling, and pocket pathway analysis to examine their molecular mechanisms of allosteric regulation. The comparative intrinsic and conformational dynamics analysis of three DNMT1 systems has highlighted the pivotal role of the RFTS domain as the dynamics hub in both intra- and inter-molecular interactions. The site perturbation and network modeling approaches have revealed the different and more complex allosteric interaction landscape in both DNMT1 complexes, involving the events caused by mutational hotspots and post-translation modification sites through protein-protein interactions (PPIs). Furthermore, communication pathway analysis and pocket detection have provided new mechanistic insights into molecular mechanisms underlying quaternary structures of DNMT1 complexes, suggesting potential targeting pockets for PPI-based allosteric drug design.     

Differences between G‐Protein‐Stabilized Agonist–GPCR Complexes and their Nanobody‐Stabilized Equivalents

Protein nanobodies have been used successfully as surrogates for unstable G‐proteins in order to crystallize G‐protein‐coupled receptors (GPCRs) in their active states. We used molecular dynamics (MD) simulations, including metadynamics enhanced sampling, to investigate the similarities and differences between GPCR–agonist ternary complexes with the α‐subunits of the appropriate G‐proteins and those with the protein nanobodies (intracellular binding partners, IBPs) used for crystallization. In two of the three receptors considered, the agonist‐binding mode differs significantly between the two alternative ternary complexes. The ternary‐complex model of GPCR activation entails enhancement of ligand binding by bound IBPs: Our results show that IBP‐specific changes can alter the agonist binding modes and thus also the criteria for designing GPCR agonists.     

Differences between G-Protein-Stabilized Agonist–GPCR Complexes and their Nanobody-Stabilized Equivalents

Protein nanobodies have been used successfully as surrogates for unstable G-proteins in order to crystallize G-protein-coupled receptors (GPCRs) in their active states. We used molecular dynamics (MD) simulations, including metadynamics enhanced sampling, to investigate the similarities and differences between GPCR–agonist ternary complexes with the α-subunits of the appropriate G-proteins and those with the protein nanobodies (intracellular binding partners, IBPs) used for crystallization. In two of the three receptors considered, the agonist-binding mode differs significantly between the two alternative ternary complexes. The ternary-complex model of GPCR activation entails enhancement of ligand binding by bound IBPs: Our results show that IBP-specific changes can alter the agonist binding modes and thus also the criteria for designing GPCR agonists.     

α‐Cyclodextrin Host–Guest Binding: A Computational Study of the Different Driving Forces

Free‐energy differences govern the equilibrium between bound and unbound states of a host and its guest molecules. The understanding of the underlying entropic and enthalpic contributions, and their complex interplay are crucial for the design of new drugs and inhibitors. In this study, molecular dynamics (MD) simulations were performed with inclusion complexes of α‐cyclodextrin (αCD) and three monosubstituted benzene derivatives to investigate host–guest binding. αCD Complexes are an ideal model system, which is experimentally and computationally well‐known. Thermodynamic integration (TI) simulations were carried out under various conditions for the free ligands in solution and bound to αCD. The two possible orientations of the ligand inside the cavity were investigated. Agreement with experimental data was only found for the more stable orientation, where the substituent resides inside the cavity. The better stability of this conformation results from stronger Van der Waals interactions and a favorable antiparallel host–guest dipole–dipole alignment. To estimate the entropic contributions, simulations were performed at three different temperatures (250, 300, and 350 K) and using positional restraints for the host. The system was found to be insensitive to both factors, due to the large and symmetric cavity of αCD, and the nondirectional nature of the host–guest interactions.     

Se-Adenosyl-l-selenomethionine Cofactor Analogue as a Reporter of Protein Methylation

Posttranslational methylation by S-adenosyl-l-methionine(SAM)-dependent methyltransferases plays essential roles in modulating protein function in both normal and disease states. As such, there is a growing need to develop chemical reporters to examine the physiological and pathological roles of protein methyltransferases. Several sterically bulky SAM analogues have previously been used to label substrates of specific protein methyltransferases. However, broad application of these compounds has been limited by their general incompatibility with native enzymes. Here we report a SAM surrogate, ProSeAM (propargylic Se-adenosyl-l-selenomethionine), as a reporter of methyltransferases. ProSeAM can be processed by multiple protein methyltransferases for substrate labeling. In contrast, sulfur-based propargylic SAM undergoes rapid decomposition at physiological pH, likely via an allene intermediate. In conjunction with fluorescent/affinity-based azide probes, copper-catalyzed azide-alkyne cycloaddition chemistry, in-gel fluorescence visualization and proteomic analysis, we further demonstrated ProSeAM's utility to profile substrates of endogenous methyltransferases in diverse cellular contexts. These results thus feature ProSeAM as a convenient probe to study the activities of endogenous protein methyltransferases.     

The large scale conformational change of the human DPP III-substrate prefers the "closed" form

Human dipeptidyl peptidase III (DPP III) is a two domain metallo-peptidase from the M49 family. The wide interdomain cleft and broad substrate specificity suggest that this enzyme could experience significant conformational change. Long (>100 ns) molecular dynamics (MD) simulations of DPP III revealed large range conformational changes of the protein, suggesting the pre-existing equilibrium model for a substrate binding. The binding free energy calculations revealed tighter binding of the preferred synthetic substrate Arg-Arg-2-naphtylamide to the "closed" than to the "open" DPP III conformation. Our assumption that Asp372 plays a crucial role in the large scale interdomain closure was proved by the MD simulations of the Asp372Ala variant. During the same simulation time, the variant remained more "open" than the wild type protein. Apparently, Ala was not as efficient as Asp in establishing the interdomain interactions. According to the MM-PBSA calculations, the electrostatic component of the free energy of solvation turned out to be higher for the "closed" protein than for its less compact form. However, the gain in entropy due to water released from the interdomain cleft nicely balanced this negative effect.     

Effect of ligand binding on the intraminimum dynamics of proteins

Effects of ligand binding on protein dynamics are studied via molecular dynamics (MD) simulations on two different enzymes, dihydrofolate reductase (DHFR) and triosephosphate isomerase (TIM), in their unliganded (free) and liganded states. Domain motions in MD trajectories are analyzed by collectivities and rotation angles along the principal components (PCs). DHFR in the free state has well‐defined domain rotations, whereas rotations are slightly damped in the binary complex with nicotinamide adenine dinucleotide phosphate (NADPH), and remarkably distorted in the presence of NADP⁺, showing that NADP⁺ is solely responsible for the loss of correlation of the domains in DHFR. Although mean square fluctuations of MD simulations in the same PC subspaces are similar for different ligation states, linear stochastic time series models show that backbone flexibility along the first five PCs is decreased upon NADPH and NADP⁺ binding in subpicosecond scale. This shows that mobility of the protein along the PCs is closely related with intraminimum dynamics, and alterations in ligation states may change the intraminimum dynamics significantly. Low vibrational frequencies of the alpha‐carbon atoms of DHFR are determined from the time series models of a larger number of low indexed PCs, and it is found that number of modes in the lowest frequencies is reduced upon ligand binding. A similar result is obtained for TIM in the unliganded and dihydroxyacetone phosphate bound states. We suggest that stochastic time series modeling is a promising method to be used in determining subtle perturbations in protein dynamics. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2011     

Stabilization of a transition-state analogue at the active site of yeast cytosine deaminase: importance of proton transfers

It is believed that the binding of pyrimidin-2-one to cytosine deaminase (CD) leads to the formation of 4-[R]-hydroxyl-3,4-dihydropyrimidine (DHP). Here the formation of transition-state analogue (TSA) at the active site of yeast cytosine deaminase (yCD) is investigated by quantum mechanical/molecular mechanical (QM/MM) molecular dynamics (MD) and free energy simulations. It is shown that DHP may in fact be unstable in the active site and a proton transfer from the Zn hydroxide group to Glu-64 may occur during the nucleophilic attack, leading to an alkoxide-like TSA complex instead. The free energy simulations for the nucleophilic attack process show that the proton transfer from the Zn hydroxide to Glu-64 may play an important role in stabilizing the TSA complex.     

分子动力学研究V82A和L90M变异对HIV PR-IDV复合物的影响

为了说明V82A和L90M变异对蛋白酶(PR)和茚地那韦(IDV)复合物的影响,进行了5.5ns的MD模拟.用MM-PBSA方法计算了体系的结合自由能,计算和实验结果一致.分解自由能为不同能量项说明,这两个变异引起熵的贡献变化大于焓的贡献变化.分解自由能到每个残基说明Wild,V82A和L90M具有相似的结合模式,结合能的贡献主要来源于A28/A28',I50/I50'和I84/I84'这六个残基组,详细分析了Wild和IDV的结合模式,对比分析了V82A和L90M变异引起结合模式的细小变化.V82A变异引起结合模式的变化是由于变异后位阻减小导致的.L90M变异引起D25和L90间的作用增强并引起结合模式的细小变化.研究结果有助于更好地理解变异对抑制剂和HIV-1PR结合模式的影响,并可以用来帮助设计更高效的PR抑制剂.     

Molecular dynamics simulation of the association of nonpolar spheres in water

We apply molecular dynamics (MD) simulations to the study of the association of nonpolar spheres of effective radii between 1.6 and 6.1 A dissolved in water. The constrained MD method is used to calculate the potential of mean force (PMF) of the interaction between spheres. The depth of the potential of mean force increases with increasing radius of the nonpolar sphere. Our results suggest that the PMF is largely governed by size or entropic effects, and that energetic effects associated with the breaking or distortion of hydrogen bonds are of minor importance.     

Insights into the dynamics of HIV-1 protease: a kinetic network model constructed from atomistic simulations

The conformational dynamics in the flaps of HIV-1 protease plays a crucial role in the mechanism of substrate binding. We develop a kinetic network model, constructed from detailed atomistic simulations, to determine the kinetic mechanisms of the conformational transitions in HIV-1 PR. To overcome the time scale limitation of conventional molecular dynamics (MD) simulations, our method combines replica exchange MD with transition path theory (TPT) to study the diversity and temperature dependence of the pathways connecting functionally important states of the protease. At low temperatures the large-scale flap opening is dominated by a small number of paths; at elevated temperatures the transition occurs through many structurally heterogeneous routes. The expanded conformation in the crystal structure 1TW7 is found to closely mimic a key intermediate in the flap-opening pathways at low temperature. We investigated the different transition mechanisms between the semi-open and closed forms. The calculated relaxation times reveal fast semi-open ? closed transitions, and infrequently the flaps fully open. The ligand binding rate predicted from this kinetic model increases by 38-fold from 285 to 309 K, which is in general agreement with experiments. To our knowledge, this is the first application of a network model constructed from atomistic simulations together with TPT to analyze conformational changes between different functional states of a natively folded protein.     

Significant enhancement of docking sensitivity using implicit ligand sampling

The efficient and accurate quantification of protein-ligand interactions using computational methods is still a challenging task. Two factors strongly contribute to the failure of docking methods to predict free energies of binding accurately: the insufficient incorporation of protein flexibility coupled to ligand binding and the neglected dynamics of the protein-ligand complex in current scoring schemes. We have developed a new methodology, named the 'ligand-model' concept, to sample protein conformations that are relevant for binding structurally diverse sets of ligands. In the ligand-model concept, molecular-dynamics (MD) simulations are performed with a virtual ligand, represented by a collection of functional groups that binds to the protein and dynamically changes its shape and properties during the simulation. The ligand model essentially represents a large ensemble of different chemical species binding to the same target protein. Representative protein structures were obtained from the MD simulation, and docking was performed into this ensemble of protein conformation. Similar binding poses were clustered, and the averaged score was utilized to rerank the poses. We demonstrate that the ligand-model approach yields significant improvements in predicting native-like binding poses and quantifying binding affinities compared to static docking and ensemble docking simulations into protein structures generated from an apo MD simulation.