Search

Thermodynamic data are often used to calibrate or test atomic-level (AL) force fields for molecular dynamics (MD) simulations. In contrast, the majority of coarse-grained (CG) force fields do not rely extensively on thermodynamic quantities. Recently, a CG force field for lipids, hydrocarbons, ions, and water, [1] in which approximately four non-hydrogen atoms are mapped onto one interaction site, has been proposed and applied to study various aspects of lipid systems. To date, no extensive investigation of its capability to describe salvation thermodynamics has been undertaken. In the present study, a detailed picture of vaporization, solvation, and phase-partitioning thermodynamics for liquid hydrocarbons and water was obtained at CG and AL resolutions, in order to compare the two types of models and evaluate their ability to describe thermodynamic properties in the temperature range between 263 and 343 K. Both CG and AL models capture the experimental dependence of the thermodynamic properties on the temperature, albeit a systematically weaker dependence is found for the CG model. Moreover, deviations are found for solvation thermodynamics and for the corresponding enthalpy-entropy compensation for the CG model. Particularly water/oil repulsion seems to be overestimated. However, the results suggest that the thermodynamic properties considered should be reproducible by a CG model provided it is reparametrized on the basis of these liquid-phase properties.

Binary mixtures of two phosphatidylcholines of different chain lengths are simulated in the bilayer state. We find a phase transition between a liquid state and a gel state at all concentrations. This phase transition is characterized by the area per lipid headgroup, the order parameter, and a change in dynamics. At concentrations with a majority of the longer lipid, we find phase separation into a gel and a liquid state in a small temperature window. This leads to a strong dynamic heterogeneity. Experimental phase transition temperatures are reproduced semiquantitatively. We see a clear shift in the phase transition to higher temperatures with increasing concentration of the longer lipid.

Extensive molecular dynamics (MD) simulations of binary systems of phospholipids and bile salts, a model for human bile, have been performed. Recent progress in hardware and software development allows simulation of the spontaneous aggregation of the constituents into small mixed micelles, in agreement with experimental observations. The MD simulations reveal the structure of these micelles at atomic detail. The phospholipids are packed radially with their headgroups at the surface and the hydrophobic tails pointing toward the micellar center. The bile salts act as wedges between the phospholipid headgroups, with their hydrophilic sides exposed to the aqueous environment. The structure of the micelles strongly resembles the previously proposed radial shell model. Simulations including small fractions of cholesterol reveal how cholesterol is solubilized inside these mixed micelles without changing their overall structure.

Haloalkane dehalogenase (DhlA) was used as a model protein to explore the possibility to use molecular dynamics (MD) simulations as a tool to identify flexible regions in proteins that can serve as a target for stability enhancement by introduction of a disulfide bond. DhlA consists of two domains: an alpha/beta-hydrolase fold main domain and a cap domain composed of five alpha-helices. MD simulations of DhlA showed high mobility in a helix-loop-helix region in the cap domain, involving residues 184-211. A disulfide cross-link was engineered between residue 201 of this flexible region and residue 16 of the main domain. The mutant enzyme showed substantial changes in both thermal and urea denaturation. The oxidized form of the mutant enzyme showed an increase of the apparent transition temperature from 47.5 to 52.5 degrees C, whereas the T(m,app) of the reduced mutant decreased by more than 8 degrees C compared to the wild-type enzyme. Urea denaturation results showed a similar trend. Measurement of the kinetic stability showed that the introduction of the disulfide bond caused a decrease in activation free energy of unfolding of 0.43 kcal mol(-1) compared to the wild-type enzyme and also indicated that the helix-loop-helix region was involved early in the unfolding process. The results show that MD simulations are capable of identifying mobile protein domains that can successfully be used as a target for stability enhancement by the introduction of a disulfide cross-link.

In this article we show how the quasi-Gaussian entropy (QGE) theory can be used to treat systems of polyatomic flexible molecules, where the usual semirigid description is not always appropriate. We describe a completely general derivation of the QGE theory which does not make use of any semirigid approximation, and therefore it is very suited for large and flexible molecules. Using molecular dynamics simulations of flexible molecules in vacuo, we investigated the ability of the theory to describe intramolecular energy fluctuations and conformational equilibria of purely classical molecules in the ideal gas condition. Results show that the gamma state level of the theory and its generalization for treating conformational equilibria (multi-gamma state model) provide excellent theoretical models when applied to three polyatomic molecules of increasing conformational freedom.

This paper presents the first atomistic simulation of a cubic membrane phase. Using the molecular dynamics simulation technique both the global and the local organization of glycerolmonoolein molecules inside the diamond cubic phase are studied. Multinanosecond simulations reveal that the center of the cubic bilayer remains close to the infinite periodic minimal surface that describes the diamond geometry. We further show that the equilibrium structure of the surfactant molecules inside the cubic phase is very similar to their structure inside a simulated lamellar bilayer. The small differences arise from the packing constraints of the surfactants within the cubic phase which has an area per surfactant that increases toward the bilayer center.

We present the newest version of the GROningen MOlecular Simulation program package, GROMOS96. GROMOS96 has been developed for the dynamic modelling of (bio)molecules using the methods of molecular dynamics, stochastic dynamics, and energy minimization as well as the path-integral formalism. An overview of its functionality is given, highlighting methodology not present in the last major release, GROMOS87. The organization of the code is outlined, and reliability, testing, and efficiency issues involved in the design of this large (73 000 lines of FORTRAN77 code) and complex package are discussed. Finally, we present two applications illustrating new functionality: local elevation simulation and molecular dynamics in four spatial dimensions.

Mol. dynamics simulation studies on the folding of beta-peptides H-beta3-HVal-beta3-HAla-beta3-HLeu-(S,S)-beta3-HAla(alphaMe)-beta3-HVal-beta3-HAla-beta3-HLeu-OH and H-beta2-HVal-beta3-HAla-beta2-HLeu-beta3-HVal-beta2-HAla-beta3-HLeu-OH were carried out. Despite the small differences in sequence between the two peptides studied, the simulations correctly predict a left-handed 31-helical fold for the beta-heptapeptide and a right-handed helical fold for the beta-hexapeptide.

The implementation of cutinase from Fusarium solani pisi as a fat-stain removing ingredient in laundry washing is hampered by its unfolding in the presence of anionic surfactants. In this work we present molecular dynamics (MD) computer simulations on cutinase and analysis procedures to distinguish the movements related to its functional behavior (e.g., substrate binding) from those related to the unfolding of the enzyme. Two kinds of MD-simulations were performed: a simulation mimicking the thermal motion at room temperature, and several simulations in which unfolding is induced either by high temperature or by using a modified water-protein interaction potential. Essential dynamics analyses (A. Amadei et al., Proteins 17:412–425, 1993) on the simulations identify distinct regions in the molecular structure of cutinase in which the motions occur for function and initial unfolding. The unfolding in various simulations starts in a similar way, suggesting that mutations in the regions involved might stabilize the enzyme without affecting its functionality. Proteins 33:253–264, 1998. © 1998 Wiley-Liss, Inc.

Nuclear spin relaxation is a powerful method for studying molecular dynamics at atomic resolution. Recent methods development in biomolecular NMR spectroscopy has enabled detailed investigations of molecular dynamics that are critical for biological function, with prominent examples addressing allostery, enzyme catalysis, and protein folding. Dynamic processes with similar correlation times are often detected in multiple locations of the molecule, raising the question of whether the underlying motions are correlated (corresponding to concerted fluctuations involving many atoms distributed across extended regions of the molecule) or uncorrelated (corresponding to independent fluctuations involving few atoms in localized regions). Here, we have used 13 C α ( i - 1)/ 13 C α ( i ) differential multiple-quantum spin relaxation to provide direct evidence for correlated dynamics of consecutive amino acid residues in the protein sequence. By monitoring overlapping pairs of residues ( i - 1 and i , i and i + 1, etc.), we identified correlated motions that extend through continuous segments of the sequence. We detected significant correlated conformational transitions in the native state of the E140Q mutant of the calmodulin C-terminal domain. Previous work has shown that this domain exchanges between two major conformational states that resemble the functionally relevant open and closed states of the WT protein, with a mean correlation time of ≈20 μs. The present results reveal that an entire α-helix undergoes partial unraveling in a transient and cooperative manner.

The design of a transmembrane four-helix bundle is described. We start with an idealized four-helix bundle geometry, then use statistical information to build a plausible transmembrane bundle. Appropriate residues are chosen using database knowledge on the sequences of membrane helices and loops, then the packing of the bundle core is optimized, and favorable side chain rotamers from rotamer libraries are selected. Next, we use explicit physical knowledge from biomolecular simulation force fields and molecular dynamics simulations to test whether the designed structure is physically possible. These procedures test whether the designed protein will indeed be α-helical, well packed and stable over a time scale of several nanoseconds in a realistic lipid bilayer environment. We then test a modeling approach that does not include sophisticated database knowledge about proteins, but rather relies on applying our knowledge of the physics that governs protein motions. This independent validation of the design is based on simulated annealing and restrained molecular dynamics simulation in vacuo, comparable to procedures used to refine NMR and X-ray structures.

To understand the effect of the finite size of simulation cells on the equilibrium properties of bilayers, an extensive series of glycerolmonoolein bilayer molecular dynamics simulations in which the surface area and system size were systematically changed have been conducted. Systems ranging from 200 to 1800 lipids were simulated, covering length scales up to 20 nm. The dependence of the surface tension on the area per lipid is shown, although long simulation times were needed (up to 40 ns) to obtain reliable estimates. As the size of simulated patches increases, long wavelength undulatory modes appear with a concomitant increase in the area compressibility due to coupling of undulation modes to area fluctuations. Both the undulatory intensities and the peristaltic intensities of the bilayer fluctuations can be fitted in the long wavelength limit to continuum model predictions. The effect of system size on surface tension appears to depend on the stress conditions.

Molecular dynamics simulations of protein folding and unfolding are often carried out at temperatures (400–600K) that are much higher than physiological or room temperature to speed up the (un)folding process. Use of such high temperatures changes both the protein and solvent properties considerably, compared to physiological or room temperature. Water models designed for use in conjunction with biomolecules, such as the simple point charge (SPC) model, have generally been calibrated at room temperature and pressure. To determine the distortive effect of high simulation temperatures on the behavior of such “room temperature” water models, the structural, dynamic, and thermodynamic properties of the much-used SPC water model are investigated in the temperature range from 300 to 500K. Both constant pressure and constant volume conditions, as used in protein simulations, were analyzed. We found that all properties analyzed change markedly with increasing temperature, but no phase transition in this temperature range was observed.

The ability to determine the free energy of solvation for a number of small organic molecules with varying sizes and properties from the coordinate trajectory of a single simulation of a given reference state was investigated. The relative free energies were estimated from a single step perturbation using the perturbation formula. The reference state consisted of a cavity surrounded by solvent. To enhance sampling a soft-core interaction was used for the cavity. The effect of the size of the cavity, the effective core height, and the length of simulation on the ability to reproduce results obtained from thermodynamic integration calculations was considered. The results using a single step perturbation from an appropriately chosen initial state were comparable to results from thermodynamic integration calculations for a wide range of compounds. Using a large number of compounds the computational efficiency was potentially increased by 2)3 orders of magnitude over traditional free energy approaches. Factors determining the efficiency of the approach are discussed. Q 1999 John Wiley & Sons, Inc. J Comput Chem 20: 1604)1617, 1999

Several series of site-directed mutations in thermolysin-like proteases are presented that show remarkable nonadditivity in their effect on thermal stability. A simple model is proposed that relates this nonadditivity to the occurrence of independent partial unfolding processes that occur in parallel at elevated temperatures. To prove this model, a thermolysin-like protease was designed in which two mutations located approximately 35 A apart in the structure individually exert small stabilizing effects of 2.3 and 4. 1 degreesC, respectively, but when combined stabilize the protease by 14.6 degreesC. This overadditivity, which follows directly from the model, confirms that unfolding of this engineered protease starts in parallel at two different regions of the protein.

The derivation of the van der Waals parameters for the aliphatic CHn united atoms of the GROMOS96 force field is presented. The parameters have been adjusted to reproduce the experimental enthalpies of vaporization and vapor pressures or densities of a set of nine alkanes in the liquid state at 298 K (or at the boiling point in the case of methane), using a cutoff radius for the van der Waals interactions of 1.6 nm. Force fields to be used in molecular simulations are bound to the conditions chosen for their parametrization, for example, the temperature, the densities of the systems included in the calibration set, or the cutoff radius used for the nonbonded interactions. Van der Waals parameters for the CHn united atoms of earlier GROMOS force fields were developed using a cutoff radius of 0.8 nm for the van der Waals interactions. Because the van der Waals interaction energy between aliphatic groups separated by distances between 0.8 and 1.4 nm is not negligible at liquid densities, the use of these parameters in combination with longer cutoffs leads to an overestimation of the attractive van der Waals interaction energy. The relevance of this excess attraction depends on the size of the groups that are interacting, as well as on their local densities. Free energies of hydration have been calculated for five alkanes. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 535–547, 1998

The attachment protein G of respiratory syncytial virus (RSV) has a modular architecture. The ectodomain of the protein comprises a small folded conserved region which is bounded by two mucin-like regions. In this study, a sequence and structural homology is described between this central conserved region of RSV-G and the fourth subdomain of the 55-kDa tumor necrosis factor receptor (TNFr). The three-dimensional structures of RSV-G and human TNFr were previously determined with NMR spectroscopy and X-ray crystallography, respectively. The C-terminal part of both subdomains fold into a cystine noose connected by two cystine bridges with the same spacing between cysteine residues and the same topology. Although a general structural similarity is observed, there are differences in secondary structure and other structural features. Molecular Dynamics calculations show that the BRSV-G NMR structure of the cystine noose is stable and that the TNFr crystal structure of the cystine noose drifts towards the BRSV-G NMR structure in the simulated solution environment. By homology modelling a model was built for the unresolved N-terminal part of the central conserved region of RSV-G. The functions for both protein domains are not known but the structural similarity of both protein domains suggests a similar function. Although the homology suggests that the cystine noose of RSV-G may interfere with the antiviral and apoptotic effect of TNF, the biological activity remains to be proven.

The impact of an extensive, uniform and hydrophobic protein surface on the behavior of the surrounding solvent is investigated. In particular, focus is placed on the possible enhancement of the structure of water at the interface, one model for the hydrophobic effect. Solvent residence times and radial distribution functions are analyzed around three types of atomic sites (methyl, polar, and positively charged sites) in 1 ns molecular dynamics simulations of the alpha-helical polypeptide SP-C in water, in methanol and in chloroform. For comparison, water residence times at positively and negatively charged sites are obtained from a simulation of a highly charged alpha-helical polypeptide from the protein titin in water. In the simulations the structure of water is not enhanced at the hydrophobic protein surface, but instead is disrupted and devoid of positional correlation beyond the first solvation sphere. Comparing solvents of different polarity, no clear trend toward the most polar solvent being more ordered is found. In addition, comparison of the water residence times at nonpolar, polar, positively charged, or negatively charged sites on the surface of SP-C or titin does not reveal pronounced or definite differences. It is shown, however, that the local environment may considerably affect solvent residence times. The implications of this work for the interpretation of the hydrophobic effect are discussed.

The self-aggregation of lipid molecules to form bilayer membranes is a process fundamental to the organization of life. Although qualitatively explained by the hydrophobic effect, 1 the molecular aggregation itself is a complex phenomenon that has not been possible to study in detail experimentally. Here, we report a series of molecular dynamics computer simulations that for the first time demonstrate the possibility to observe the entire process at atomic detail with realistic lipids. Starting from random solutions, bilayers are formed on time scales of 10-100 ns, with properties matching experimental data. Several key steps and approximate time scales of the aggregation can be identified. The final rate-limiting process is the reduction and disappearance of large hydrophilic transmembrane water pores, of biological relevance for, for example, ion permeation. Singer and Nicholson were the first to recognize the implications of the extreme flexibility of membranes for the structure of cellular walls, leading to the famous fluid-mosaic model 2 with diffusing lipids and proteins. The bilayer formation process is, however, extremely fast and involves subtle rearrangements at the molecular level, making it elusive to current experimental methods. Simplified computer models have been used to mimic aggregation of surfactant-like molecules into monolayers and micelles, 3 bilayerlike structures, 4,5 and even vesicles. 6 These models are theoretically important to extend length and time scales, but they do not include atomic detail like hydrogen bonds and represent the collective entropic effects driving aggregation 1 as pairwise interactions. Detailed molecular dynamics simulations have on the other hand, provided accurate models of up to nanometer and nanosecond scales, but previously only for preassembled bilayers. 7-11 This work demonstrates the first simulations of aggregation of lipids into bilayers with atomic detail of the structure and interactions. Compared to micelle aggregation studies, 12,13 bilayer formation is considerably more challenging due to the balance between hydrophobicity and solvation, and the aggregation involves collective mesoscopic dynamics. The phospholipid dipalmitoylphosphatidylcholine (DPPC) was initially chosen for the study, since it is present in biological

Important biological processes, such as vesicle fusion or budding, require the cell matrix to undergo a transition from a lamellar to a nonlamellar state. Although equilibrium properties of membranes are amenable to detailed theoretical studies, collective rearrangements involved in phase transitions have thus far only been modeled on a qualitative level. Here, for the first time, the complete transition pathway from a multilamellar to an inverted hexagonal phase is elucidated at near-atomic detail using a recently developed coarse-grained molecular dynamics simulation model. Insight is provided into experimentally inaccessible data such as the molecular structure of the intermediates and the kinetics involved. Starting from multilamellar configurations, the spontaneous formation of stalks between the bilayers is observed on a nanosecond timescale at elevated temperatures or reduced hydration levels. The stalks subsequently elongate in a cooperative manner leading to the formation of an inverted hexagonal phase. The rate of stalk elongation is approximately 0.1 nm ns(-1). Within a narrow hydration/temperature/composition range the stalks appear stable and rearrange into the rhombohedral phase.

Structural, dynamic and energetic properties of proteins in solution can be studied in atomic detail by molecular dynamics computer simulation. Protein unfolding can be caused by a variety of driving forces induced in different ways: increased temperature or pressure, change of solvent composition, or protein amino acid mutation. The stability and unfolding of four different proteins (bovine pancreatic trypsin inhibitor, hen egg white lysozyme, the surfactant protein C and the DNA-binding domain of the 434 repressor) have been studied by applying the afore-mentioned driving forces and also to some artificial forces. The results give a picture of protein (in)stability and possible unfolding pathways, and are compared to experimental data where possible.

Understanding the binding and insertion of peptides in lipid bilayers is a prerequisite for understanding phenomena such as antimicrobial activity and membrane-protein folding. We describe molecular dynamics simulations of the antimicrobial peptide alamethicin in lipid/water and octane/water environments, taking into account an external electric field to mimic the membrane potential. At cis-positive potentials, alamethicin does not insert into a phospholipid bilayer in 10ns of simulation, due to the slow dynamics of the peptide and lipids. However, in octane N-terminal insertion occurs at field strengths from 0.33V/nm and higher, in simulations of up to 100ns duration. Insertion of alamethicin occurs in two steps, corresponding to desolvation of the Gln7 side chain, and the backbone of Aib10 and Gly11. The proline induced helix kink angle does not change significantly during insertion. Polyalanine and alamethicin form stable helices both when inserted in octane and at the water/octane interface, where they partition in the same location. In water, both polyalanine and alamethicin partially unfold in multiple simulations. We present a detailed analysis of the insertion of alamethicin into the octane slab and the influence of the external field on the peptide structure. Our findings give new insight into the mechanism of channel formation by alamethicin and the structure and dynamics of membrane-associated helices.

In order to investigate the dependence of the viscosity on the mass of the molecules in a liquid, and thus check the validity of Stokes' law for molecules, several molecular dynamics simulations of `water' molecules with different mass and different molecular mass distributions were performed. The viscosity is shown to be sensitive to the mass but less sensitive to the mass distribution. The product of diffusion coefficient and viscosity, which according to Stokes' law should be independent of the mass, varies. We may therefore conclude that Stokes' law is not valid for small molecules.

Several simulation methods are currently in use to study peptide and protein folding. They can be classified in three groups depending on how the solvent is treated. At the simplest level, the solvent is ignored. At a second level, solvent effects are implicitly represented in the atomic interaction function. At the third level, solvent degrees of freedom are treated explicitly. We have performed molecular dynamics simulations in vacuum, stochastic dynamics simulations, and molecular dynamics simulations in methanol of the folding of a β-heptapeptide to test these different levels of approximation. The results clearly show that the solvent cannot be ignored in folding simulations, and highlight the difficulties of implicitly modelling solvent effects in an atomic interaction function for a solute.