Exploration of the activation pathway of Δα-Chymotrypsin with molecular dynamics simulations and correlation with kinetic experiments

Correlating the experimentally observed kinetics of protein conformational changes with theoretical predictions is a formidable and challenging task, due to the multitude of degrees of freedom (>5,000) in a protein and the huge gap between the timescale of the kinetic event of interest (ms) and the typical timescale of computer simulations (ns). In this study we show that using the targeted molecular dynamics (TMD) method it is possible to simulate conformational changes of the ms time range and to correlate multiple simulations of single pathways with ensemble experiments on both the structural and energetic basis. As a model system we chose to study the conformational change of rat-Δα-chymotrypsin from its inactive to its active conformation. This activation process has been analyzed previously by experimental and theoretical methods, i.e. fluorescence stopped-flow spectroscopy (FSF), molecular dynamics (MD) and TMD. Inspired by the results of these studies on the wild type (WT) enzyme, several mutants were constructed to alter the conformational pathway and studied by FSF measurements. In the present work WT and mutant N18G were subjected to multiple MD and subsequent TMD simulations. We report the existence of two main activation pathways, a feature of chymotrypsin activation that has been hitherto unknown. A method to correlate the energetics of the different pathways calculated by TMD and the kinetic parameters observed by experimental methods such as FSF is presented. Our work is relevant for experimental single molecule studies of enzymes in general. Electronic supplementary material  The online version of this article (doi:) contains supplementary material, which is available to authorized users.     

Structure and dynamics of two β-peptides in solution from molecular dynamics simulations validated against experiment

We have studied two different beta-peptides in methanol using explicit solvent molecular dynamics simulations and the GROMOS 53A6 force field: a heptapeptide (peptide 1) expected to form a left-handed 3(14)-helix, and a hexapeptide (peptide 2) expected to form a beta-hairpin in solution. Our analysis has focused on identifying and analyzing the stability of the dominant secondary structure conformations adopted by the peptides, as well as on comparing the experimental NOE distance upper bounds and 3J-coupling values with their counterparts calculated on the basis of the simulated ensembles. Moreover, we have critically compared the present results with the analogous results obtained with the GROMOS 45A3 (peptide 1) and 43A1 (peptide 2) force fields. We conclude that within the limits of conformational sampling employed here, the GROMOS 53A6 force field satisfactorily reproduces experimental findings regarding the behavior of short beta-peptides, with accuracy that is comparable to but not exceeding that of the previous versions of the force field. GCE legend Conformational clustering analysis of the simulated ensemble of a ss-hexapeptide with two different simulation setups (a and b). The central members of all of the clusters populating more than 5% of all of the structures are shown, together with the most dominant hydrogen bonds and the corresponding percentages of cluster members containing them.     

Influence of the SNPs on the structural stability of CBS protein： Insight from molecular dynamics simulations

Cystathionine β-synthase is an essential enzyme of the trans-suifuration pathway that condenses serine with homocysteine to form cystathionine. Missense mutations in CBS are the major cause of inherited homocystinuria, and the detailed effect of disease associated amino acid substitutions on the structure and stability of human CBS is yet unknown. Here, we apply a unique approach in combining in silico tools and molecular dynamics simulation to provide structural and functional insight into the effect of SNP on the stability and activity of mutant CBS. In addition, principal component analysis and free energy landscape were used to predict the collective motions, thermodynamic stabilities and essential subspace relevant to CBS function. The obtained results indicate that C109R, E176K and D376N mutations have the diverse effect on dynamic behavior of CBS protein. We found that highly conserved D376N mutation, which is present in the active pocket, affects the protein folding mechanism. Our strategy may provide a way in near future to understand and study effects of functional nsSNPs and their role in causing homocystinuria.     

Study of the Alzheimer's Aβ40 peptide in SDS micelles using molecular dynamics simulations

The interaction of the Alzheimer's amyloid beta peptide, Aβ40, with sodium dodecyl sulfate (SDS) micelles, together with the self-assembly of SDS molecules around the peptide from an initial random distribution were studied using atomistic and coarse-grained (CG) molecular dynamics simulations. In atomistic simulations, the peptide structure in the micelle was characterized by two helical regions connected through a short hinge. The initial structure of the system was shown to affect the simulation results. The atomistic self-assembly of SDS molecules resulted in a 38-molecule micelle around the peptide, along with some globules and individual molecules. Coarse-grained simulation results, however, did not show such a difference, and at the end of all CG simulations, a complete 60-molecule micelle was obtained, with the peptide located at the interface of the micelle with water. The obtained CG radial density profiles and SDS micelle size and shape properties were identical for all CG simulations.     

Study of the influence of temperature on the dynamics of the catalytic cleft in 1,3-1,4-β-glucanase by molecular dynamics simulations

The dependence of some molecular motions in the enzyme 1,3-1,4-β-glucanase from Bacillus licheniformis on temperature changes and the role of the calcium ion in them were explored. For this purpose, two molecular dynamics simulated trajectories along 4 ns at low (300 K) and high (325 K) temperatures were generated by the GROMOS96 package. Several structural and thermodynamic parameters were calculated, including entropy values, solvation energies, and essential dynamics (ED). In addition, thermoinactivation experiments to study the influence of the calcium ion and some residues on the activity were conducted. The results showed the release of the calcium ion, which, in turn, significantly affected the movements of loops 1, 2, and 3, as shown by essential dynamics. These movements differ at low and high temperatures and affect dramatically the activity of the enzyme, as observed by thermoinactivation studies. The first two authors contributed equally to this work     

The inhibitors of cyclin-dependent kinases and GSK-3β enhance osteoclastogenesis

Osteoclasts are multinucleated cells with bone resorption activity that is crucial for bone remodeling. RANK‐RANKL (receptor activator of nuclear factor κB ligand) signaling has been shown as a main signal pathway for osteoclast differentiation. However, the molecular mechanism and the factors regulating osteoclastogenesis remain to be fully understood. In this study, we performed a chemical genetic screen, and identified a Cdks/GSK-3β (cyclin-dependent kinases/glycogen synthase kinase 3β) inhibitor, kenpaullone, and two Cdks inhibitors, olomoucine and roscovitine, all of which significantly enhance osteoclastogenesis of RAW264.7 cells by upregulating NFATc1 (nuclear factor of activated T cells, cytoplasmic 1) levels. We also determined that the all three compounds increase the number of osteoclast differentiated from murine bone marrow cells. Furthermore, the three inhibitors, especially kenpaullone, promoted maturation of cathepsin K, suggesting that the resorption activity of the resultant osteoclasts is also activated. Our findings indicate that inhibition of GSK-3β and/or Cdks enhance osteoclastogenesis by modulating the RANK–RANKL signaling pathway.     

Characterization of the internal dynamics and conformational space of zinc-bound amyloid β peptides by replica-exchange molecular dynamics simulations

Amyloid β (Aβ) peptides and metal ions have been associated with the pathogenesis of Alzheimer’s disease. The conformational space of Aβ fragments of different length with and without binding of metal ions has been extensively investigated by replica-exchange molecular dynamics (REMD) simulation. However, only trajectories extracted at relatively low temperatures have been used for this analysis. The capability of REMD simulations to characterize the internal dynamics of such intrinsically disordered proteins (IDPs) as Aβ has been overlooked. In this work, we use an approach recently developed by Xue and Skrynnikov (J Am Chem Soc 133:14614–14628, 2011) to calculate NMR observables, including ¹⁵N relaxation rates and ¹⁵N–¹H nuclear Overhauser enhancement (NOE), from the high-temperature trajectory of REMD simulations for zinc-bound Aβ peptides. The time axis of the trajectory was rescaled to correct for the effect of the high temperature (408 K) compared with the experimental temperature (278 K). Near-quantitative agreement between simulated values and experimental results was obtained. When the structural properties and free-energy surfaces of zinc-bound Aβ_(1–40) and Aβ_(1–42) were compared at the physiological temperature 310 K it was found that zinc-bound Aβ_(1–42) was more rigid than Aβ_(1–40) at the C terminus, and its conformational transitions were also more preferred. The self-consistent results derived from trajectories at high and low temperatures demonstrate the capability of REMD simulations to capture the internal dynamics of IDPs.     

Insights on the acetylated NF-κB transcription factor complex with DNA from molecular dynamics simulations

The nuclear factor-κB (NF-κB) is a DNA sequence-specific regulator of many important biological processes, whose activity is modulated by enzymatic acetylation. In one of the best functionally characterized NF-κB complexes, the p50/p65 heterodimer, acetylation of K221 at p65 causes a decrease of DNA dissociation rate, whilst the acetylation of K122 and K123, also at p65, markedly decreases the binding affinity for DNA. By means of molecular dynamics simulations based on the X-ray structure of the p50/p65 complex with DNA, we provide insights on the structural determinants of the acetylated complexes in aqueous solution. Lysine acetylation involves the loss of favorable electrostatic interactions between DNA and NF-κB, which is partially compensated by the reduction of the desolvation free-energy of the two binding partners. Acetylation at both positions K122 and K123 is associated with a decrease of the electrostatic potential at the p65/DNA interface, which is only partially counterbalanced by an increase of the local Na(+) concentration. It induces the disruption of base-specific and nonspecific interactions between DNA and NF-κB and it is consistent with the observed decrease of binding affinity. In contrast, acetylation at position K221 results in the loss of nonspecific protein-DNA interactions, but the DNA recognition sites are not affected. In addition, the loss of protein-DNA interactions is likely to be counterbalanced by an increase of the configurational entropy of the complex, which provides, at a speculative level, a justification for the observed decrease of NF-κB/DNA dissociation rate.     

Unfolding of the amyloid β-peptide central helix: mechanistic insights from molecular dynamics simulations

Alzheimer's disease (AD) pathogenesis is associated with formation of amyloid fibrils caused by polymerization of the amyloid β-peptide (Aβ), which is a process that requires unfolding of the native helical structure of Aβ. According to recent experimental studies, stabilization of the Aβ central helix is effective in preventing Aβ polymerization into toxic assemblies. To uncover the fundamental mechanism of unfolding of the Aβ central helix, we performed molecular dynamics simulations for wild-type (WT), V18A/F19A/F20A mutant (MA), and V18L/F19L/F20L mutant (ML) models of the Aβ central helix. It was quantitatively demonstrated that the stability of the α-helical conformation of both MA and ML is higher than that of WT, indicating that the α-helical propensity of the three nonpolar residues (18, 19, and 20) is the main factor for the stability of the whole Aβ central helix and that their hydrophobicity plays a secondary role. WT was found to completely unfold by a three-step mechanism: 1) loss of α-helical backbone hydrogen bonds, 2) strong interactions between nonpolar sidechains, and 3) strong interactions between polar sidechains. WT did not completely unfold in cases when any of the three steps was omitted. MA and ML did not completely unfold mainly due to the lack of the first step. This suggests that disturbances in any of the three steps would be effective in inhibiting the unfolding of the Aβ central helix. Our findings would pave the way for design of new drugs to prevent or retard AD.     

A salt-bridge switch in the molecular recognition between RS receptor and RGD ligand from the ABEEM σπ molecular dynamics simulations

Abstract

How the receptor and ligand recognise each other is a challenging subject in explaining the mechanism of recognition at the molecular level. As a starting point, here, a synthesised RS receptor and its RGD ligand were investigated as a proper model to simulate their recognition process in terms of ABEEMσπ/MM polarisable force field. It is found that a switch of forming up a salt bridge in the ligand triggers the recognition of the receptor and ligand. Through the salt-bridge switch that undergoes several cycles from on-state with parallel hydrogen bonds to off-state with bifurcated hydrogen bonds, the active site of ligand can flex easily to interact with the active site of the receptor. In addition, the water molecules form a decisive bridge connecting the active sites of the bound system. The salt-bridge switch and water-mediated movement are cooperative as the important factors for the receptor-ligand recognition. In addition, the properties, such as binding free energy, conformational flexibility and solvent accessible surface area have been calculated to provide adequate evidence for the whole recognition process. According to the simulation, a detailed mechanism was derived involving diffusion, a switch triggered cooperative water-mediated movement, and conformational folding, for the flexible recognition.     

Ligand-escape pathways from the ligand-binding domain of PPARγ receptor as probed by molecular dynamics simulations

Conformational rearrangements of peroxysome proliferator activated receptor (PPARγ) ligand-binding domain (LBD) that accompany the release and binding of ligands are not well understood. To determine the major events associated with the escape of the partial agonist GW0072, molecular dynamic (MD) simulations were performed using two different methods: reversed targeted molecular dynamics (TMD⁻¹) and time-dependent distance restraints (TDR) using as restraints either the root mean square deviation from a reference structure (TMD⁻¹) or the distance between the geometrical centers of the binding pocket and of the ligand (TDR). Both methods do not assume any a priori route for ligand extraction. To avoid artifacts, different initial simulation conditions were used and particular attention was paid for giving time to the protein to relax during the extraction process by running 10–12 ns simulations within explicit water. Two distinct exit gates A and B were found, independently of initial conditions and method. During the exit process no interaction between GW0072 and the transactivation AF-2 helix was observed. Our results suggest that the ligand uses the intrinsic flexibility of the protein to move within the receptor. Paths A and B are very similar to those found for other nuclear receptors, suggesting that these routes are a common characteristics of nuclear receptors that are used by different kinds of ligands. Finally, the knowledge of entry/exit pathways of a receptor should be very useful in discriminating between different ligands that could have been favorably docked in the binding pocket by introducing docking along these pathways into computational drug design protocols.     

Binding modes of phosphotriesterase-like lactonase complexed with δ-nonanoic lactone and paraoxon using molecular dynamics simulations

Phosphotriesterase-like lactonases (PLLs) have received much attention because of their physical and chemical properties. They may have widespread applications in various fields. For example, they show potential for quorum-sensing signaling pathways and organophosphorus (OP) detoxification in agricultural science. However, the mechanism by which PLLs hydrolyze, which involves OP compounds and lactones and a variety of distinct catalytic efficiencies, has only rarely been explored. In the present study, molecular dynamics (MD) simulations were performed to characterize and contrast the structural dynamics of DrPLL, a member of the PLL superfamily in Deinococcus radiodurans, bound to two substrates, δ-nonanoic lactone and paraoxon. It has been observed that there is a 16-fold increase in the catalytic efficiency of the two mutant strains of DrPLL (F26G/C72I) vs. the wild-type enzyme toward the hydrolysis of paraoxon, but an explanation for this behavior is currently lacking. The analysis of the molecular trajectories of DrPLL bound to δ-nonanoic lactone indicated that lactone-induced conformational changes take place in loop 8, which is near the active site. Binding to paraoxon may lead to conformational displacement of loop 1 residues, which could lead to the deformation of the active site and so trigger the entry of the paraoxon into the active site. The efficiency of the F26G/C72I mutant was increased by decreasing the displacement of loop 1 residues and increasing the flexibility of loop 8 residues. These results provide a molecular-level explanation for the experimental behavior.     

Partial activation of α7 nicotinic acetylcholine receptors: insights from molecular dynamics simulations

Nicotinic acetylcholine receptors (nAChRs) are drug targets for neuronal disorders and diseases. Partial agonists for nAChRs are currently being developed as drugs for the treatment of neurological diseases for their relative safety originated from reduced excessive stimulation. In the current study, molecular docking, molecular dynamics simulations and binding energy calculations were performed to theoretically investigate the interactions between the partial agonists, 4-OH-DMXBA and tropisetron with α7-nAChR. The results suggest that the partial agonists 4-OH-DMXBA and tropisetron bind with α7-nAChR in a binding mode similar to that with AChBP. The non-conserved residues in the binding sites contribute to the orientation deviation of these partial agonists from their orientation in AChBP. Energy calculation and decomposition using MM-GB/SA suggests that the van der Waals term (ΔE_VDW) is the main driving force for the binding of the partial agonists to α7-nAChR. The molecular dynamics simulations showed that the opening of the C-loop binding with the partial agonists is in-between the openings for the binding with the full agonist and in the apo state. This conformation difference for the C-loop sheds light on the partial agonism of nAChR.     

Can molecular dynamics simulations improve the structural accuracy and virtual screening performance of GPCR models?

The determination of G protein-coupled receptor (GPCR) structures at atomic resolution has improved understanding of cellular signaling and will accelerate the development of new drug candidates. However, experimental structures still remain unavailable for a majority of the GPCR family. GPCR structures and their interactions with ligands can also be modelled computationally, but such predictions have limited accuracy. In this work, we explored if molecular dynamics (MD) simulations could be used to refine the accuracy of in silico models of receptor-ligand complexes that were submitted to a community-wide assessment of GPCR structure prediction (GPCR Dock). Two simulation protocols were used to refine 30 models of the D₃ dopamine receptor (D₃R) in complex with an antagonist. Close to 60 μs of simulation time was generated and the resulting MD refined models were compared to a D₃R crystal structure. In the MD simulations, the receptor models generally drifted further away from the crystal structure conformation. However, MD refinement was able to improve the accuracy of the ligand binding mode. The best refinement protocol improved agreement with the experimentally observed ligand binding mode for a majority of the models. Receptor structures with improved virtual screening performance, which was assessed by molecular docking of ligands and decoys, could also be identified among the MD refined models. Application of weak restraints to the transmembrane helixes in the MD simulations further improved predictions of the ligand binding mode and second extracellular loop. These results provide guidelines for application of MD refinement in prediction of GPCR-ligand complexes and directions for further method development.     

Oxidation-induced destabilization of the fibrinogen αC-domain dimer investigated by molecular dynamics simulations

Upon activation, fibrinogen is converted to insoluble fibrin, which assembles into long strings called protofibrils. These aggregate laterally to form a fibrin matrix that stabilizes a blood clot. Lateral aggregation of protofibrils is mediated by the αC domain, a partially structured fragment located in a disordered region of fibrinogen. Polymerization of αC domains links multiple fibrin molecules with each other enabling the formation of thick fibrin fibers and a fibrin matrix that is stable but can also be digested by enzymes. However, oxidizing agents produced during the inflammatory response have been shown to cause thinner fibrin fibers resulting in denser clots, which are harder to proteolyze and pose the risk of deep vein thrombosis and lung embolism. Oxidation of Met⁴⁷⁶ located within the αC domain is thought to hinder its ability to polymerize disrupting the lateral aggregation of protofibrils and leading to the observed thinner fibers. How αC domains assemble into polymers is still unclear and yet this knowledge would shed light on the mechanism through which oxidation weakens the lateral aggregation of protofibrils. This study used temperature replica exchange molecular dynamics simulations to investigate the αC-domain dimer and how this is affected by oxidation of Met⁴⁷⁶. Analysis of the trajectories revealed that multiple stable binding modes were sampled between two αC domains while oxidation decreased the likelihood of dimer formation. Furthermore, the side chain of Met⁴⁷⁶ was observed to act as a docking spot for the binding and this function was impaired by its conversion to methionine sulfoxide.     

Effects of ligands on unfolding of the amyloid β-peptide central helix: mechanistic insights from molecular dynamics simulations

Polymerization of the amyloid β-peptide (Aβ), a process which requires that the helical structure of Aβ unfolds beforehand, is suspected to cause neurodegeneration in Alzheimer's disease. According to recent experimental studies, stabilization of the Aβ central helix counteracts Aβ polymerization into toxic assemblies. The effects of two ligands (Dec-DETA and Pep1b), which were designed to bind to and stabilize the Aβ central helix, on unfolding of the Aβ central helix were investigated by molecular dynamics simulations. It was quantitatively demonstrated that the stability of the Aβ central helix is increased by both ligands, and more effectively by Pep1b than by Dec-DETA. In addition, it was shown that Dec-DETA forms parallel conformations with β-strand-like Aβ, whereas Pep1b does not and instead tends to bend unwound Aβ. The molecular dynamics results correlate well with previous experiments for these ligands, which suggest that the simulation method should be useful in predicting the effectiveness of novel ligands in stabilizing the Aβ central helix. Detailed Aβ structural changes upon loss of helicity in the presence of the ligands are also revealed, which gives further insight into which ligand may lead to which path subsequent to unwinding of the Aβ central helix.     

Simulating the antimicrobial mechanism of human β-defensin-3 with coarse-grained molecular dynamics

Human β-defensin-3 (HβD-3) is an endogenous antimicrobial peptide with potent and broad killing activity against various microorganisms, and thus, it is an attractive candidate for the development of novel peptide antibiotics, but its antimicrobial mechanism remains elusive. To characterize the mechanism, we used multi-microsecond coarse-grained simulations with the MARTINI force field. These simulations show HβD-3 peptides can form oligomers on the surface of bacterial membrane and make anionic lipids (POPG) clustered. Furthermore, two kinds of regions (one is composed of pure POPG lipids, and the other is enriched in POPE lipids) are formed in the membrane; on the border of them, there are some obvious defects, which result in the membrane disruption. By contrast, the simulations also reveal that the contacts between the HβD-3 peptides and mammalian membrane are not stable. These results provide biophysical insights into HβD-3 selectivity and suggest a possible antimicrobial mechanism.     

Influencing the direct formation of ‘Janus’ particles in molecular dynamics simulations

The binary nucleation of phase-separated Lennard-Jones clusters was analysed under various system conditions using molecular dynamics simulations. The modified potential model provides a simple gateway to observe non-wetting behaviour and imitates the more complex interactions of non-miscible substances. Thus, not only the transition from ideally mixed clusters to so-called ‘Janus’ particles, but also the structural aspects and dynamic formation processes of nanoscopic droplets are directly observable from the gas phase. Various shapes and sizes of these inhomogeneous clusters have been found via simple tuning of system parameters. From this analysis, we gained further insight into the direct formation of ‘Janus’ particles from the gas phase.