Dynamics and unfolding pathways of a hyperthermophilic and a mesophilic rubredoxin.

Molecular dynamics simulations in solution are performed for a rubredoxin from the hyperthermophilic archaeon Pyrococcus furiosus (RdPf) and one from the mesophilic organism Desulfovibrio vulgaris (RdDv). The two proteins are simulated at four temperatures: 300 K, 373 K, 473 K (two sets), and 500 K; the various simulations extended from 200 ps to 1,020 ps. At room temperature, the two proteins are stable, remain close to the crystal structure, and exhibit similar dynamic behavior; the RMS residue fluctuations are slightly smaller in the hyperthermophilic protein. An analysis of the average energy contributions in the two proteins is made; the results suggest that the intraprotein energy stabilizes RdPf relative to RdDv. At 373 K, the mesophilic protein unfolds rapidly (it begins to unfold at 300 ps), whereas the hyperthermophilic does not unfold over the simulation of 600 ps. This is in accord with the expected stability of the two proteins. At 473 K, where both proteins are expected to be unstable, unfolding behavior is observed within 200 ps and the mesophilic protein unfolds faster than the hyperthermophilic one. At 500 K, both proteins unfold; the hyperthermophilic protein does so faster than the mesophilic protein. The unfolding behavior for the two proteins is found to be very similar. Although the exact order of events differs from one trajectory to another, both proteins unfold first by opening of the loop region to expose the hydrophobic core. This is followed by unzipping of the beta-sheet. The results obtained in the simulation are discussed in terms of the factors involved in flexibility and thermostability.     

Effects of environment on the structure of Pyrococcus furiosus rubredoxin: a molecular dynamics study

Molecular dynamics simulations based on a 0.95-A resolution crystal structure of Pyrococcus furiosus have been performed to elucidate the effects of the environment on the structure of rubredoxin, and proteins in general. Three 1-ns simulations are reported here: two crystalline state simulations at 123 and 300 K, and a solution state simulation at 300 K. These simulations show that temperature has a greater impact on the protein structure than the close molecular contacts of the crystal matrix in rubredoxin, although both have an effect on its dynamic properties. These results indicate that differences between NMR solution structures and X-ray crystal structures will be relatively minor if they are done at similar temperatures. In addition, the crystal simulations appears to mimic previous crystallographic experiments on the effects of cryo-temperature on temperature factors, and might provide a useful tool in the structural analysis of protein structures solved at cryo-temperatures.     

Modeling the structure of Pyrococcus furiosus rubredoxin by homology to other X-ray structures.

The three-dimensional structure of rubredoxin from the hyperthermophilic archaebacterium, Pyrococcus furiosus, has been modeled from the X-ray crystal structures of three homologous proteins from Clostridium pasteurianum, Desulfovibrio gigas, and Desulfovibrio vulgaris. All three homology models are similar. When comparing the positions of all heavy atoms and essential hydrogen atoms to the recently solved crystal structure (Day, M. W., et al., 1992, Protein Sci. 1, 1494-1507) of the same protein, the homology model differ from the X-ray structure by 2.09 A root mean square (RMS). The X-ray and the zinc-substituted NMR structures (Blake, P. R., et al., 1992b, Protein Sci. 1, 1508-1521) show a similar level of difference (2.05 A RMS). On average, the homology models are closer to the X-ray structure than to the NMR structures (2.09 vs. 2.42 A RMS).     

Refinement of comparative models of protein structure by using multicanonical molecular dynamics simulations

Comparative modelling is a powerful method that easily predicts a considerably accurate structure of a protein by using a template structure having a similar amino-acid sequence to the target protein. However, in the region where the amino-acid sequence is different between the target and the template, the predicted structure remains unreliable. In such a case, the model has to be refined. In the present study, we explored the possibility of a molecular dynamics-based method, using the human SAP Src Homology 2 (SH2) domain as the modelling target. The multicanonical method was used to alleviate the multiple-minima problem and the generalised Born/surface area model was used to reduce the computational cost. In addition, position restraints were imposed on the atoms in the reliable regions to avoid unnecessary conformational sampling. We analyzed the conformational distribution of the ligand-recognition loop of the domain and found that the most populated conformational clusters in the ensemble of the model agreed well with one of the two major clusters in the ensemble of the reference simulation starting from the crystal structure. This demonstrates that the current refinement method can significantly improve the accuracy of an unreliable region in a comparative model.     

Exploring single-domain antibody thermostability by molecular dynamics simulation

Abstract

Single-domain antibodies also known as nanobodies are recombinant antigen-binding domains that correspond to the heavy-chain variable region of camelid antibodies. Previous experimental studies showed that the nanobodies have stable and active structures at high temperatures. In this study, the thermal stability and dynamics of nanobodies have been studied by employing molecular dynamics simulation at different temperatures. Variations in root mean square deviation, native contacts, and solvent-accessible surface area of the nanobodies during the simulation were calculated to analyze the effect of different temperatures on the overall conformation of the nanobody. Then, the thermostability mechanism of this protein was studied through calculation of dynamic cross-correlation matrix, principal component analyses, native contact analyses, and root mean square fluctuation. Our results manifest that the side chain conformation of some residues in the complementarity-determining region 3 (CDR3) and also the interaction between α-helix region of CDR3 and framework2 play a critical role to stabilize the protein at a high temperature.

Communicated by Ramaswamy H. Sarma     

Reproducible polypeptide folding and structure prediction using molecular dynamics simulations

The folding of a polypeptide from an extended state to a well-defined conformation is studied using microsecond classical molecular dynamics (MD) simulations and replica exchange molecular dynamics (REMD) simulations in explicit solvent and in vacuo. It is shown that the solvated peptide folds many times in the REMD simulations but only a few times in the conventional simulations. From the folding events in the classical simulations we estimate an approximate folding time of 1-2 micros. The REMD simulations allow enough sampling to deduce a detailed Gibbs free energy landscape in three dimensions. The global minimum of the energy landscape corresponds to the native state of the peptide as determined previously by nuclear magnetic resonance (NMR) experiments. Starting from an extended state it takes about 50 ns before the native structure appears in the REMD simulations, about an order of magnitude faster than conventional MD. The calculated melting curve is in good qualitative agreement with experiment. In vacuo, the peptide collapses rapidly to a conformation that is substantially different from the native state in solvent.     

Understanding the effect of tert-butanol on Candida antarctica lipase B using molecular dynamics simulations

The use of organic solvents as reaction media for enzymatic reactions has many advantages. Several organic solvents have been proposed as reaction media, especially for transesterifications using Candida antarctica lipase B (CalB). Among organic solvents, tert-butanol is associated with an enhanced conversion rate in bio-diesel production. Thus, it is necessary to understand the effect of tert-butanol on CalB to explain the high-catalytic efficiency compared with the reaction in other hydrophilic organic solvents. In this study, the effects of tert-butanol on the structure of CalB were investigated by MD simulations. The overall flexibility was increased in the presence of tert-butanol. The substrate entrance and the binding pocket size of CalB in tert-butanol were maintained as in TIP3P water. The distance between the catalytic residues of CalB in tert-butanol indicated a higher likelihood of forming hydrogen bonds. These structural analyses could be useful for understanding the effect of tert-butanol on lipase transesterification.     

Leucine 41 is a gate for water entry in the reduction of Clostridium pasteurianum rubredoxin

Biological electron transfer is an efficient process even though the distances between the redox moieties are often quite large. It is therefore of great interest to gain an understanding of the physical basis of the rates and driving forces of these reactions. The structural relaxation of the protein that occurs upon change in redox state gives rise to the reorganizational energy, which is important in the rates and the driving forces of the proteins involved. To determine the structural relaxation in a redox protein, we have developed methods to hold a redox protein in its final oxidation state during crystallization while maintaining the same pH and salt conditions of the crystallization of the protein in its initial oxidation state. Based on 1.5 A resolution crystal structures and molecular dynamics simulations of oxidized and reduced rubredoxins (Rd) from Clostridium pasteurianum (Cp), the structural rearrangements upon reduction suggest specific mechanisms by which electron transfer reactions of rubredoxin should be facilitated. First, expansion of the [Fe-S] cluster and concomitant contraction of the NH...S hydrogen bonds lead to greater electrostatic stabilization of the extra negative charge. Second, a gating mechanism caused by the conformational change of Leucine 41, a nonpolar side chain, allows transient penetration of water molecules, which greatly increases the polarity of the redox site environment and also provides a source of protons. Our method of producing crystals of Cp Rd from a reducing solution leads to a distribution of water molecules not observed in the crystal structure of the reduced Rd from Pyrococcus furiosus. How general this correlation is among redox proteins must be determined in future work. The combination of our high-resolution crystal structures and molecular dynamics simulations provides a molecular picture of the structural rearrangement that occurs upon reduction in Cp rubredoxin.     

Driven to near-experimental accuracy by refinement via molecular dynamics simulations

Protein model refinement has been an essential part of successful protein structure prediction. Molecular dynamics simulation-based refinement methods have shown consistent improvement of protein models. There had been progress in the extent of refinement for a few years since the idea of ensemble averaging of sampled conformations emerged. There was little progress in CASP12 because conformational sampling was not sufficiently diverse due to harmonic restraints. During CASP13, a new refinement method was tested that achieved significant improvements over CASP12. The new method intended to address previous bottlenecks in the refinement problem by introducing new features. Flat-bottom harmonic restraints replaced harmonic restraints, sampling was performed iteratively, and a new scoring function and selection criteria were used. The new protocol expanded conformational sampling at reduced computational costs. In addition to overall improvements, some models were refined significantly to near-experimental accuracy.     

Effect of temperature on welding of metallic nanowires investigated using molecular dynamics simulations

The mechanism of welding of Au–Au, Ag–Ag and Au–Ag nanowires (NWs) with head-to-head contact is studied using molecular dynamics simulations based on the second-moment approximation of the many-body tight-binding potential. The effect of temperature in the range of 300–900 K is investigated. Simulation results show that at the initial welding, an incomplete jointing area forms through the interactions of the van der Waals attractive force, and that the jointing area increases with increasing the extent of contact between the two NWs during the welding process and temperature. Few defects form along the (1 1 1) close-packed plane during the welding process because the acting stress is quite low. Among the three NW pairs, the Au–Au NWs have the best cold-welding quality, whereas the Au–Ag NWs have the worst cold-welding quality due to the welding of different materials. With an increase in temperature, the weld stress and the mechanical strength of the NWs significantly decrease, and the number of disordered structures increases. The welding fails when the temperature exceeds the molten temperature of the NWs.     

Exploration of the chlorpyrifos escape pathway from acylpeptide hydrolases using steered molecular dynamics simulations

Acylpeptide hydrolases (APH) catalyze the removal of an N-acylated amino acid from blocked peptides. APH is significantly more sensitive than acetylcholinesterase, a target of Alzheimer’s disease, to inhibition by organophosphorus (OP) compounds. Thus, OP compounds can be used as a tool to probe the physiological functions of APH. Here, we report the results of a computational study of molecular dynamics simulations of APH bound to the OP compounds and an exploration of the chlorpyrifos escape pathway using steered molecular dynamics (SMD) simulations. In addition, we apply SMD simulations to identify potential escape routes of chlorpyrifos from hydrolase hydrophobic cavities in the APH-inhibitor complex. Two previously proposed APH pathways were reliably identified by CAVER 3.0, with the estimated relative importance of P1 > P2 for its size. We identify the major pathway, P2, using SMD simulations, and Arg526, Glu88, Gly86, and Asn65 are identified as important residues for the ligand leaving via P2. These results may help in the design of APH-targeting drugs with improved efficacy, as well as in understanding APH selectivity of the inhibitor binding in the prolyl oligopeptidase family.     

Trypsinogen activation as observed in accelerated molecular dynamics simulations

Serine proteases are involved in many fundamental physiological processes, and control of their activity mainly results from the fact that they are synthetized in an inactive form that becomes active upon cleavage. Three decades ago Martin Karplus's group performed the first molecular dynamics simulations of trypsin, the most studied member of the serine protease family, to address the transition from the zymogen to its active form. Based on the computational power available at the time, only high frequency fluctuations, but not the transition steps, could be observed. By performing accelerated molecular dynamics (aMD) simulations, an interesting approach that increases the configurational sampling of atomistic simulations, we were able to observe the N‐terminal tail insertion, a crucial step of the transition mechanism. Our results also support the hypothesis that the hydrophobic effect is the main force guiding the insertion step, although substantial enthalpic contributions are important in the activation mechanism. As the N‐terminal tail insertion is a conserved step in the activation of serine proteases, these results afford new perspective on the underlying thermodynamics of the transition from the zymogen to the active enzyme.     

Ligand unbinding pathways from the vitamin D receptor studied by molecular dynamics simulations

Molecular dynamics simulation techniques have been used to study the unbinding pathways of 1α,25-dihydroxyvitamin D₃ from the ligand-binding pocket of the vitamin D receptor (VDR). The pathways observed in a large number of relatively short (<200 ps) random acceleration molecular dynamics (RAMD) trajectories were found to be in fair agreement, both in terms of pathway locations and deduced relative preferences, compared to targeted molecular dynamics (TMD) and streered molecular dynamics simulations (SMD). However, the high-velocity ligand expulsions of RAMD tend to favor straight expulsion trajectories and the observed relative frequencies of different pathways were biased towards the probability of entering a particular exit channel. Simulations indicated that for VDR the unbinding pathway between the H1–H2 loop and the β-sheet between H5 and H6 is more favorable than the pathway located between the H1–H2 loop and H3. The latter pathway has been suggested to be the most likely unbinding path for thyroid hormone receptors (TRs) and a likely path for retinoic acid receptor. Ligand entry/exit through these two pathways would not require displacement of H12 from its agonistic position. Differences in the packing of the H1, H2, H3 and β-sheet region explain the changed relative preference of the two unbinding pathways in VDR and TRs. Based on the crystal structures of the ligand binding domains of class 2 nuclear receptors, whose members are VDR and TRs, this receptor class can be divided in two groups according to the packing of the H1, H2, H3 and β-sheet region. Electronic supplementary material  The online version of this article (doi:) contains supplementary material, which is available to authorized users.     

Using molecular dynamics for the refinement of atomistic models of GPCRs by homology modeling

Despite GPCRs sharing a common seven helix bundle, analysis of the diverse crystallographic structures available reveal specific features that might be relevant for ligand design. Despite the number of crystallographic structures of GPCRs steadily increasing, there are still challenges that hamper the availability of new structures. In the absence of a crystallographic structure, homology modeling remains one of the important techniques for constructing 3D models of proteins. In the present study we investigated the use of molecular dynamics simulations for the refinement of GPCRs models constructed by homology modeling. Specifically, we investigated the relevance of template selection, ligand inclusion as well as the length of the simulation on the quality of the GPCRs models constructed. For this purpose we chose the crystallographic structure of the rat muscarinic M3 receptor as reference and constructed diverse atomistic models by homology modeling, using different templates. Specifically, templates used in the present work include the human muscarinic M2; the more distant human histamine H1 and the even more distant bovine rhodopsin as shown in the GPCRs phylogenetic tree. We also investigated the use or not of a ligand in the refinement process. Hence, we conducted the refinement process of the M3 model using the M2 muscarinic as template with tiotropium or NMS docked in the orthosteric site and compared with the results obtained with a model refined without any ligand bound.     

Unraveling the role of flexible coil near calcium binding site of levansucrase on thermostability and product profile via proline substitution and molecular dynamics simulations

Due to its bioactivity and versatile applications, levan has appeared as a promising biomaterial. Levansucrase is responsible for the conversion of sucrose into levan. With the goal of enhancing levan production, the strategy for enhancing the stability of levansucrase is being intensively studied. To make proteins more stable under high temperatures, proline, the most rigid residue, can be introduced into previously flexible regions. Herein, G249, D250, N251, and H252 on the flexible coil close to the calcium binding site of Bacillus licheniformis levansucrase were replaced with proline. Mutations at G249P greatly enhance both the enzyme's thermodynamic and kinetic stability, while those at H252P improve solely the enzyme's kinetic stability. GPC analysis revealed that G249P synthesize more levan, but H252P generate primarily oligosaccharides. Molecular dynamics simulations (MD) and MM/GBSA analysis revealed that G249P mutation increased not only the stability of levansucrase, but also affinity toward fructan.     

Sensitivity of molecular dynamics simulations to the choice of the X-ray structure used to model an enzymatic reaction

A subject of great practical importance that has not received much attention is the question of the sensitivity of molecular dynamics simulations to the initial X-ray structure used to set up the calculation. We have found two cases in which seemingly similar structures lead to quite different results, and in this article we present a detailed analysis of these cases. The first case is acyl-CoA dehydrogenase, and the chief difference of the two structures is attributed to a slight shift in a backbone carbonyl that causes a key residue (the proton-abstracting base) to be in a bad conformation for reaction. The second case is xylose isomerase, and the chief difference of the two structures appears to be the ligand sphere of a Mg2+ metal cofactor that plays an active role in catalysis.     

Investigation of ligand selectivity in CYP3A7 by molecular dynamics simulations

Cytochrome P450 (CYP) 3A7 plays a crucial role in the biotransformation of the metabolized endogenous and exogenous steroids. To compare the metabolic capabilities of CYP3A7–ligands complexes, three endogenous ligands were selected, namely dehydroepiandrosterone (DHEA), estrone, and estradiol. In this study, a three-dimensional model of CYP3A7 was constructed by homology modeling using the crystal structure of CYP3A4 as the template and refined by molecular dynamics simulation (MD). The docking method was adopted, combined with MD simulation and the molecular mechanics generalized born surface area method, to probe the ligand selectivity of CYP3A7. These results demonstrate that DHEA has the highest binding affinity, and the results of the binding free energy were in accordance with the experimental conclusion that estrone is better than estradiol. Moreover, several key residues responsible for substrate specificity were identified on the enzyme. Arg372 may be the most important residue due to the low interaction energies and the existence of hydrogen bond with DHEA throughout simulation. In addition, a cluster of Phe residues provides a hydrophobic environment to stabilize ligands. This study provides insights into the structural features of CYP3A7, which could contribute to further understanding of related protein structures and dynamics.