Towards molecular dynamics simulation of large proteins with a hydration shell at constant pressure.

Molecular dynamics simulation of a large protein in explicit water with periodic boundary conditions is extremely demanding in terms of computation time. Consequently, we have sought approximations of the solvent environment that model its important features. Here, we describe our SAPHYR (Shell Approximation for Protein HYdRation) model in which the protein is surrounded by a shell of water molecules maintained at constant pressure. In addition to the usual pairwise interatomic interactions, these water molecules are subjected to forces approximating van der Waals and dipole-dipole interactions with the implicit surrounding bulk solvent. The SAPHYR model is tested for a system of one argon atom in water and for the protein ubiquitin, and then applied to cytochrome P450cam, a protein with over 400 residues. The results demonstrate that structural and dynamic properties of the simulated systems are improved by use of the SAPHYR model, and that this model provides a significant computational saving over simulations with periodic boundary conditions.     

Ab initio molecular dynamics study of the hydration of the formohydroxamate anion

We apply ab initio molecular dynamics (AIMD) to study the hydration structures and electronic properties of the formohydroxamate anion in liquid water. We consider the cis- nitrogen-deprotonated, cis- oxygen-deprotonated, and trans- oxygen-deprotonated formohydroxamate tautomers. They form an average of 6.3, 6.9, and 6.0 hydrogen bonds with water molecules, respectively. The predicted pair correlation functions and time dependence of the hydration numbers suggest that water is highly structured around the nominally negatively charged oxime oxygen in O-deprotonated tautomers but significantly less so around the nitrogen atom in the N-deprotonated species. Wannier function analysis suggests that, in the O-deprotonated anions, the negative charge is concentrated on the oxime oxygen, while in the N-deprotonated case, it is partially delocalized between the nitrogen and the adjoining oxime oxygen atom.     

Structure of hydration water in proteins: a comparison of molecular dynamics simulations and database analysis

Hydration layer water molecules play important structural and functional roles in proteins. Despite being a critical component in biomolecular systems, characterizing the properties of hydration water poses a challenge for both experiments and simulations. In this context we investigate the local structure of hydration water molecules as a function of the distance from the protein and water molecules respectively in 188 high resolution protein structures and compare it with those obtained from molecular dynamics simulations. Tetrahedral order parameter of water in proteins calculated from previous and present simulation studies show that the potential of bulk water overestimates the average tetrahedral order parameter compared to those calculated from crystal structures. Hydration waters are found to be more ordered at a distance between the first and second solvation shell from the protein surface. The values of the order parameter decrease sharply when the water molecules are located very near or far away from the protein surface. At small water-water distance, the values of order parameter of water are very low. The average order parameter records a maximum value at a distance equivalent to the first solvation layer with respect to the water-water radial distribution and asymptotically approaches a constant value at large distances. Results from present analysis will help to get a better insight into structure of hydration water around proteins. The analysis will also help to improve the accuracy of water models on the protein surface.     

A molecular dynamics study of conformational changes and hydration of left-handed d(CGCGCGCGCGCG)2 in a nonsalt solution.

Twelve dinucleotides (one complete turn) of left-handed, flexible, double-helix poly(dG-dC) Z-DNA have been simulated in aqueous solution with K+ counterions for 70 ps. Most of the d(GpC) phosphates have rotated in accordance with a ZI----ZII transition. The ZII conformation was probably partly stabilized by counterions, which coordinate one of the anionic oxygens and the guanine-N7 of the next (5'----3' direction) base. The presence of base-coordinating ions close to the helical axis rotated and pulled about half of the d(CpG) phosphates further into the groove. These ions also gave rise to rather large deviations from the crystal structure (ZI) with their tendency of pulling the bases closer toward the helical axis. A flipping of the orientation about the glycosyl bond from the +sc to the -sc region was observed for one guanosine, also leading to deviations from the crystal structure. Many bridges containing one or two water molecules were found, with a dominance for the latter. They essentially formed a network of intra- and interstrand bridges between anionic and esterified phosphate oxygens. A "spine" of water molecules could be distinguished as a dark zig-zag pattern in the water density map. The lifetime of a bridge containing one water was about twice as long as that of a two-water bridge and it lasted 5-15 times longer than a hydrogen bond in water. The lifetimes were also calculated for a selection of bridge types, in order of decreasing stability: O1P/O2P ... W ... O'4 much greater than O1P/O2P ... W ... guanine-N2 greater than O1P/O2P ... W ... O1P/O2P. The reorientational motion of water molecules in the first hydration shell around selected groups was slowed down considerably compared to bulk water and the decreasing order of correlation times was guanine-N2 greater than O'4 greater than O'3/O'5 greater than O1P/O2P.     

A systematic molecular dynamics simulation study of temperature dependent bilayer structural properties

Although lipid force fields (FFs) used in molecular dynamics (MD) simulations have proved to be accurate, there has not been a systematic study on their accuracy over a range of temperatures. Motivated by the X-ray and neutron scattering measurements of common phosphatidylcholine (PC) bilayers (Ku?erka et al. BBA. 1808: 2761, 2011), the CHARMM36 (C36) FF accuracy is tested in this work with MD simulations of six common PC lipid bilayers over a wide range of temperatures. The calculated scattering form factors and deuterium order parameters from the C36 MD simulations agree well with the X-ray, neutron, and NMR experimental data. There is excellent agreement between MD simulations and experimental estimates for the surface area per lipid, bilayer thickness (D_B), hydrophobic thickness (D_C), and lipid volume (V_L). The only minor discrepancy between simulation and experiment is a measure of (D_B − D_HH) / 2 where D_HH is the distance between the maxima in the electron density profile along the bilayer normal. Additional MD simulations with pure water and heptane over a range of temperatures provide explanations of possible reasons causing the minor deviation. Overall, the C36 FF is accurate for use with liquid crystalline PC bilayers of varying chain types and over biologically relevant temperatures.     

A comparative molecular dynamics study on thermostability of human and chicken prion proteins

To compare the thermostabilities of human and chicken normal cellular prion proteins (HuPrP(C) and CkPrP(C)), molecular dynamics (MD) simulations were performed for both proteins at an ensemble level (10 parallel simulations at 400 K and 5 parallel simulations at 300 K as a control). It is found that the thermostability of HuPrP(C) is comparable with that of CkPrP(C), which implicates that the non-occurrence of prion diseases in non-mammals cannot be completely attributed to the thermodynamic properties of non-mammalian PrP(C).     

A molecular dynamics simulation study of an aminoglycoside/A-site RNA complex: conformational and hydration patterns

Aminoglycoside antibiotics interfere with the translation mechanism by binding to the tRNA decoding site of the 16S ribosomal RNA. Crystallographic structures of aminoglycosides bound to A-site systems clarified many static aspects of RNA-ligand interactions. To gain some insight on the dynamic aspects of recognition phenomena, we conducted molecular dynamics simulations of the aminoglycoside paromomycin bound to a eubacterial ribosomal decoding A-site oligonucleotide. Results from 25 ns of simulation time revealed that: (i) the neamine part of the antibiotic represents the main anchor for binding, (ii) additional sugar rings provide limited and fragile contacts, (iii) long-resident water molecules present at the drug/RNA interface are involved in the recognition phenomena. The combination of MD simulations together with systematic structural information offers striking insights into the molecular recognition processes underlying RNA/aminoglycoside binding. Important methodological considerations related to the use of medium resolution starting structures and associated sampling problems are thoroughly discussed.     

A molecular dynamics analysis of protein structural elements

The relation between protein secondary structure and internal motions was examined by using molecular dynamics to calculate positional fluctuations of individual helix, beta-sheet, and loop structural elements in free and substrate-bound hen egg-white lysozyme. The time development of the fluctuations revealed a general correspondence between structure and dynamics; the fluctuations of the helices and beta-sheets converged within the 101 psec period of the simulation and were lower than average in magnitude, while the fluctuations of the loop regions were not converged and were mostly larger than average in magnitude. Notable exceptions to this pattern occurred in the substrate-bound simulation. A loop region (residues 101-107) of the active site cleft had significantly reduced motion due to interactions with the substrate. Moreover, part of a loop and a 3(10) helix (residues of 67-88) not in contact with the substrate showed a marked increase in fluctuations. That these differences in dynamics of free and substrate-bound lysozyme did not result simply from sampling errors was established by an analysis of the variations in the fluctuations of the two halves of the 101 psec simulation of free lysozyme. Concerted transitions of four to five mainchain phi and psi angles between dihedral wells were shown to be responsible for large coordinate shifts in the loops. These transitions displaced six or fewer residues and took place either abruptly, in 1 psec or less, or with a diffusive character over 5-10 psec. Displacements of rigid secondary structures involved longer timescale motions in bound lysozyme; a 0.5 A rms change in the position of a helix occurred over the 55 psec simulation period. This helix reorientation within the protein appears to be a response to substrate binding. There was little correlation between the solvent accessible surface area and the dynamics of the different structural elements.     

Role of flexibility and polarity as determinants of the hydration of internal cavities and pockets in proteins

Molecular dynamics simulations of Staphylococcal nuclease and of 10 variants with internal polar or ionizable groups were performed to investigate systematically the molecular determinants of hydration of internal cavities and pockets in proteins. In contrast to apolar cavities in rigid carbon structures, such as nanotubes or buckeyballs, internal cavities in proteins that are large enough to house a few water molecules will most likely be dehydrated unless they contain a source of polarity. The water content in the protein interior can be modulated by the flexibility of protein elements that interact with water, which can impart positional disorder to water molecules, or bias the pattern of internal hydration that is stabilized. This might explain differences in the patterns of hydration observed in crystal structures obtained at cryogenic and room temperature conditions. The ability of molecular dynamics simulations to determine the most likely sites of water binding in internal pockets and cavities depends on its efficiency in sampling the hydration of internal sites and alternative protein and water conformations. This can be enhanced significantly by performing multiple molecular dynamics simulations as well as simulations started from different initial hydration states.     

Investigations of peptide hydration using NMR and molecular dynamics simulations: A study of effects of water on the conformation and dynamics of antamanide

Summary The influence of water binding on the conformational dynamics of the cyclic decapeptide antamanide dissolved in the model lipophilic environment chloroform is investigated by NMR relaxation measurements. The water-peptide complex has a lifetime of 35 s at 250 K, which is longer than typical lifetimes of water-peptide complexes reported in aqueous solution. In addition, there is a rapid intracomplex mobility that probably involves librational motions of the bound water or water molecules hopping between different binding sites. Water binding restricts the flexibility of antamanide. The experimental findings are compared with GROMOS molecular dynamics simulations of antamanide with up to eight bound water molecules. Within the simulation time of 600 ps, no water molecule leaves the complex. Additionally, the simulations show a reduced flexibility for the complex in comparison with uncomplexed antamanide. Thus, there is a qualitative agreement between the experimental NMR results and the computer simulations.     

Comparative thermodynamic elucidation of the structural stability of thermophilic proteins

Differential scanning calorimetry, circular dichroism, and visible absorption spectrophotometry were employed to elucidate the structural stability of thermophilic phycocyanin derived from Cyanidium caldarium, a eucaryotic organism which contains a nucleus, grown in acidic conditions (pH 3.4) at 54°C. The obtained results were compared with those previously reported for thermophilic phycocyanin derived from Synechococcus lividus, a procaryote containing no organized nucleus, grown in alkaline conditions (pH 8.5) at 52°C. The temperature of thermal unfolding (t_d) was found to be comparable between C. caldarium (73°C) and S. lividus (74°C) phycocyanins. The apparent free energy of unfolding (ΔG_[urea]=0) at zero denaturant (urea) concentration was also comparable: 9.1 and 8.7 kcal/mole for unfolding the chromophore part of the protein, and 5.0 and 4.3 kcal/mole for unfolding the apoprotein part of the protein, respectively. These values of t_d and ΔG_[urea]=0 were significantly higher than those previously reported for mesophilic Phormidium luridum phycocyanin (grown at 25°C). These findings revealed that relatively higher values of t_d and ΔG_[urea]=0 were characteristics of thermophilic proteins. In contrast, the enthalpies of completed unfolding (ΔH_d) and the half-completed unfolding (ΔH_d)1/2 for C. caldarium phycocyanin were much lower than those for S. lividus protein (89 versus 180 kcal/mole and 62 versus 115 kcal/mole, respectively). Factors contributing to a lower ΔH_d in C. caldarium protein and the role of charged groups in enhancing the stability of thermophilic proteins were discusse.     

Modeling of the structure of bacteriorhodopsin. A molecular dynamics study.

Secondary structure predictions for membrane proteins are relatively reliable and permit the construction of model structures that may serve as initial conformations for molecular dynamics simulations. This might provide a scheme to predict the three-dimensional structures of membrane proteins. The feasibility of such an approach is tested for bacteriorhodopsin. We were not able to fully predict the kidney-shaped structure of bacteriorhodopsin. However, features compatible with this structure developed in a simulation starting from a circular arrangement of the seven predicted helices. When instead we started from the kidney shape, assigning the seven predicted helices in different ways to those on the structure, we could distinguish between the different assignments on the basis of energy and tilt of the helices. In this way we could select the correct assignment from a few others. For the correct assignment, the helices spontaneously adopted a tilt that agrees remarkably well with the experimental model structure derived by others. The root-mean-square deviation between our best molecular dynamics structure and the experimental model structure is 3.8 A, caused mainly by deviations in the internal degrees of freedom of the helices.     

Water and polypeptide conformations in the gramicidin channel. A molecular dynamics study.

Theoretical studies of ion channels address several important questions. The mechanism of ion transport, the role of water structure, the fluctuations of the protein channel itself, and the influence of structural changes are accessible from these studies. In this paper, we have carried out a 70-ps molecular dynamics simulation on a model structure of gramicidin A with channel waters. The backbone of the protein has been analyzed with respect to the orientation of the carbonyl and the amide groups. The results are in conformity with the experimental NMR data. The structure of water and the hydrogen bonding network are also investigated. It is found that the water molecules inside the channel act as a collective chain; whereas the conformation in which all the waters are oriented with the dipoles pointing along the axis of the channel is a preferred one, others are also accessed during the dynamics simulation. A collective coordinate involving the channel waters and some of the hydrogen bonding peptide partners is required to describe the transition of waters from one configuration to the other.     

A study of the structural properties and thermal stability of human Bcl-2 by molecular dynamics simulations

The anti-apoptotic B-cell lymphoma 2 (Bcl-2) protein interacts with several proteins that regulate the apoptotic properties of cells. In this research, we conduct several all-atom molecular dynamics (MD) simulations under high-temperature unfolding conditions, from 400 to 800?K, for 25?ns. These simulations were performed using a model of an engineered Bcl-2 human protein (Bcl-2-Δ22Σ3), which lacks 22 C-terminal residues of the transmembrane domain. The aim of this study is to gain insight into the structural behavior of Bcl-2-Δ22Σ3 by mapping the conformational movements involved in Bcl-2 stability and its biological function. To build a Bcl-2-Δ22Σ3 three-dimensional model, the protein core was built by homology modeling and the flexible loop domain (FLD, residues 33-91) by ab initio methods. Further, the entire protein model was refined by MD simulations. Afterwards, the production MD simulations showed that the FLD at 400 and 500?K has several conformations reaching into the protein core, whereas at 600?K some of the alpha-helices were lost. At 800?K, the Bcl-2 core is destabilized suggesting a possible mechanism for protein unfolding, where the alpha helices 1 and 6 were the most stable, and a reduction in the number of hydrogen bonds initially occurs. In conclusion, the structural changes and the internal protein interactions suggest that the core and the FLD are crucial components of Bcl-2 in its function of regulate ng access to the recognition sites of kinases and caspases.     

Structural insights into the effect of hydration and ions on A-tract DNA: a molecular dynamics study

DNA structure is known to be sensitive to hydration and ionic environment. To explore the dynamics, hydration, and ion binding features of A-tract sequences, a 7-ns Molecular dynamics (MD) study has been performed on the dodecamer d(CGCAAATTTGCG)(2). The results suggest that the intrusion of Na(+) ion into the minor groove is a rare event and the structure of this dodecamer is not very sensitive to the location of the sodium ions. The prolonged MD simulation successfully leads to the formation of sequence dependent hydration patterns in the minor groove, often called spine of hydration near the A-rich region and ribbon of hydration near the GC regions. Such sequence dependent differences in the hydration patterns have been seen earlier in the high resolution crystal structure of the Drew-Dickerson sequence, but not reported for the medium resolution structures (2.0 approximately 3.0 A). Several water molecules are also seen in the major groove of the MD simulated structure, though they are not highly ordered over the extended MD. The characteristic narrowing of the minor groove in the A-tract region is seen to precede the formation of the spine of hydration. Finally, the occurrence of cross-strand C2-H2.O2 hydrogen bonds in the minor groove of A-tract sequences is confirmed. These are found to occur even before the narrowing of the minor groove, indicating that such interactions are an intrinsic feature of A-tract sequences.     

A molecular dynamics study of gating in dioxolane-linked gramicidin A channels.

The gating transition of the RR and SS dioxolane ring-linked gramicidin A channels were studied with molecular dynamics simulations using a detailed atomic model. It was found that the probable reaction path, describing the transition of the ring from the exterior to the interior of the channel where it blocked the permeation pathway, involved several steps including the isomerization of the transpeptide plane dihedral angle of Val1. Reaction coordinates along this pathway were defined, and the transition rates between the stable conformers were calculated. It was found, in good accord with experimental observations, that the calculated blocking rate for the RR-linked channel was 280/s with a mean blocking time of 0.04 ms, whereas such blocking did not occur in the case of the SS-linked channel. An important observation is that the resulting lifetime for the blocked state of the RR-linked channel was in good accord with the experimental observations only when the calculations were performed in the presence of a potassium ion inside the channel.     

A Connected-cluster of hydration around myoglobin: Correlation between molecular dynamics simulations and experiment

An analysis of a molecular dynamics simulation of metmyoglobin in an explicit solvent environment of 3,128 water molecules has been performed. Both statics and dynamics of the protein-solvent interface are addressed in a comparison with experiment. Three-dimensional density distributions, temperature factors, and occupancy weights are computed for the solvent by using the trajectory coordinates. Analysis of the hydration leads to the localization of more than 500 hydration sites distributed into multiple layers of solvation located between 2.6 and 6.8 Å from the atomic protein surface. After locating the local solvent density maxima or hydration sites we conclude that water molecules of hydration positions and hydration sites are distinct concepts. Both global and detailed properties of the hydration cluster around myoglobin are compared with recent neutron and X-ray data on myoglobin. Questions arising from differences between X-ray and neutron data concerning the locations of the protein-bound water are investigated. Analysis of water site differences found from X-ray and neutron experiments compared with our simulation shows that the simulation gives a way to unify the hydration picture given by the two experiments. © 1994 John Wiley & Sons, Inc.     

A molecular dynamics approach to the structural characterization of amyloid aggregation

A novel computational approach to the structural analysis of ordered beta-aggregation is presented and validated on three known amyloidogenic polypeptides. The strategy is based on the decomposition of the sequence into overlapping stretches and equilibrium implicit solvent molecular dynamics (MD) simulations of an oligomeric system for each stretch. The structural stability of the in-register parallel aggregates sampled in the implicit solvent runs is further evaluated using explicit water simulations for a subset of the stretches. The beta-aggregation propensity along the sequence of the Alzheimer's amyloid-beta peptide (Abeta(42)) is found to be highly heterogeneous with a maximum in the segment V(12)HHQKLVFFAE(22) and minima at S(8)G(9), G(25)S(26), G(29)A(30), and G(38)V(39), which are turn-like segments. The simulation results suggest that these sites may play a crucial role in determining the aggregation tendency and the fibrillar structure of Abeta(42). Similar findings are obtained for the human amylin, a 37-residue peptide that displays a maximal beta-aggregation propensity at Q(10)RLANFLVHSSNN(22) and two turn-like sites at G(24)A(25) and G(33)S(34). In the third application, the MD approach is used to identify beta-aggregation "hot-spots" within the N-terminal domain of the yeast prion Ure2p (Ure2p(1-94)) and to design a double-point mutant (Ure2p-N4748S(1-94)) with lower beta-aggregation propensity. The change in the aggregation propensity of Ure2p-N4748S(1-94) is verified in vitro using the thioflavin T binding assay.     

Microtubule dynamics in the spindle. II. A thermodynamic and kinetic description

We have previously presented a model for the assembly and disassembly of mitotic spindle microtubules (MTs) (Pickett-Heaps et al., 1986). In this paper, we describe the thermodynamics of such spindle MT assembly and present equations to describe the polymerization kinetics of different classes of spindle MTs. These equations are used to predict, in terms of kinetics parameters, the magnitude of forces extant on spindle MTs and to define the critical force needed to halt MT assembly. We calculate several of these forces for a hypothetical model cell; our predicted value for the force generated along kinetochore fibers is in close agreement with measured values taken from living cells. The model and its implications are discussed with reference to other recent models of spindle and MT dynamics.