Dynamics of TMAO and urea in the hydration shell of the protein SNase.

Using all-atom molecular dynamics simulations of aqueous solutions of the globular protein SNase, the dynamic behavior of water molecules and cosolvents (trimethylamine-N-oxide (TMAO) and urea) in the hydration shell of the protein was studied for different solvent compositions. TMAO is a potent protein-stabilizing osmolyte, whereas urea is known to destabilize proteins. For molecules that are initially located in successive narrow layers at a given distance from the protein, the mean displacements and the distribution of displacements for short time intervals are calculated. For molecules that are initially located in solvation shells of a given thickness around the protein, the characteristic residence times in these shells are determined to characterize the dynamic behavior of the solvent molecules as a function of the distance to the protein. A combined consideration of these characteristics allows to reveal additional features of the dynamics of the cosolvents. It is shown that TMAO molecules leave the nearest vicinity of the protein faster than urea molecules, despite the fact that the mobility of TMAO molecules, measured by their mean displacements, is lower than that of urea. Moreover, we show that the rate of release of TMAO molecules from the hydration shell is lower in ternary (TMAO + urea + H2O) solvent mixtures than in the binary ones. This is consistent with a recent observation that the fraction of TMAO near the protein decreases in the presence of urea. From the analysis of the decay of the number of particles initially located in the region of the first peak of the distribution function of solvent molecules around the protein, we estimated that about 20 water molecules and 6-7 urea molecules stay near the protein for more than 1000 ps.

TMAO and urea in the hydration shell of the protein SNase.

We performed all-atom MD simulations of the protein SNase in aqueous solution and in the presence of two major osmolytes, trimethylamine-N-oxide (TMAO) and urea, as cosolvents at various concentrations and compositions and at different pressures and temperatures. The distributions of the cosolvent molecules and their orientation in the surroundings of the protein were analyzed in great detail. The distribution of urea is largely conserved near the protein. It varies little with pressure and temperature, and does practically not depend on the addition of TMAO. The slight decrease with temperature of the number of urea molecules that are in contact with the SNase molecule is consistent with the view that the interaction of the protein with urea is mainly of enthalpic nature. Most of the TMAO molecules tend to be oriented to the protein by its methyl groups, a small amount of these molecules contact the protein by its oxygen, forming hydrogen bonds with the protein, only. Unlike urea, the fraction of TMAO in the hydration shell of SNase slightly increases with temperature (a signature of a prevailing hydrophobic interaction between TMAO and SNase), and decreases significantly upon the addition of urea. This behavior reflects the diverse nature of the interaction of the two osmolytes with the protein. Using the Voronoi volume of the atoms of the solvent molecules (water, urea, TMAO), we compared the fraction of the volume occupied by a given type of solvent molecule in the hydration shell and in the bulk solvent. The volume fraction of urea in the hydration shell is more than two times larger than in the bulk, whereas the volume fraction of TMAO in the hydration shell is only slightly larger in the binary solvent (TMAO + water) and becomes even less than in the bulk in the ternary solvent (TMAO + water + urea). Thus, TMAO tends to be excluded from the hydration shell of the protein. The behavior of the two cosolvents in the vicinity of the protein does not change much with pressure (from 1 to 5000 bar) and temperature (from 280 to 330 K). This is also in line with the conception of the "osmophobic effect" of TMAO to protect proteins from denaturation also at harsh environmental conditions. We also calculated the volumetric parameters of SNase and found that the cosolvents have a small but significant effect on the apparent volume and its contributions, i.e. the intrinsic, molecular and thermal volumes.

Exploring volume, compressibility and hydration changes of folded proteins upon compression.

Understanding the physical basis of the structure, stability and function of proteins in solution, including extreme environmental conditions, requires knowledge of their temperature and pressure dependent volumetric properties. One physical-chemical property of proteins that is still little understood is their partial molar volume and its dependence on temperature and pressure. We used molecular dynamics simulations of aqueous solutions of a typical monomeric folded protein, staphylococcal nuclease (SNase), to study and analyze the pressure dependence of the apparent volume, Vapp, and its components by the Voronoi-Delaunay method. We show that the strong decrease of Vapp with pressure (ßT = 0.95 × 10(-5) bar(-1), in very good agreement with the experimental value) is essentially due to the compression of the molecular volume, VM, ultimately, of its internal voids, V. Changes of the intrinsic volume (defined as the Voronoi volume of the molecule), the contribution of the solvent to the apparent volume, and of the contribution of the boundary voids between the protein and the solvent have also been studied and quantified in detail. The pressure dependences of the volumetric characteristics obtained are compared with the temperature dependent behavior of these quantities and with corresponding results for a natively unfolded polypeptide.

Disentangling volumetric and hydrational properties of proteins.

We used molecular dynamics simulations of a typical monomeric protein, SNase, in combination with Voronoi-Delaunay tessellation to study and analyze the temperature dependence of the apparent volume, Vapp, of the solute. We show that the void volume, VB, created in the boundary region between solute and solvent, determines the temperature dependence of Vapp to a major extent. The less pronounced but still significant temperature dependence of the molecular volume of the solute, VM, is essentially the result of the expansivity of its internal voids, as the van der Waals contribution to VM is practically independent of temperature. Results for polypeptides of different chemical nature feature a similar temperature behavior, suggesting that the boundary/hydration contribution seems to be a universal part of the temperature dependence of Vapp. The results presented here shine new light on the discussion surrounding the physical basis for understanding and decomposing the volumetric properties of proteins and biomolecules in general.

Calculation of the volumetric characteristics of biomacromolecules in solution by the Voronoi-Delaunay technique.

Recently a simple formalism was proposed for a quantitative analysis of interatomic voids inside a solute molecule and in the surrounding solvent. It is based on the Voronoi-Delaunay tessellation of structures, obtained in molecular simulations: successive Voronoi shells are constructed, starting from the interface between the solute molecule and the solvent, and continuing to the outside (into the solvent) as well as into the interior of the molecule. Similarly, successive Delaunay shells, consisting of Delaunay simplexes, can also be constructed. This technique can be applied to interpret volumetric data, obtained, for example, in studies of proteins in aqueous solution. In particular, it allows replacing qualitatively and descriptively introduced properties by strictly defined quantities, such as the thermal volume, by the boundary voids. The extension and the temperature behavior of the boundary region, its structure and composition are discussed in detail, using the example of a molecular dynamics model of an aqueous solution of the human amyloid polypeptide, hIAPP. We show that the impact of the solute on the local density of the solvent is short ranged, limited to the first Delaunay and the first Voronoi shell around the solute. The extra void volume, created in the boundary region between solute and solvent, determines the magnitude and the temperature dependence of the apparent volume of the solute molecule.

Volumetric properties of hydrated peptides: Voronoi-Delaunay analysis of molecular simulation runs.

The study of hydration, folding, and interaction of proteins by volumetric measurements has been promoted by recent advances in the development of highly sensitive instrumentations. However, the separation of the measured apparent volumes into contributions from the protein and the hydration water, V(app) = V(int) + ΔV, is still challenging, even with the detailed microscopic structural information from molecular simulations. By the examples of the amyloidogenic polypeptides hIAPP and Aß42 in aqueous solution, we analyze molecular dynamics simulation runs for different temperatures, using the Voronoi-Delaunay tessellation method. This method allows a parameter free determination of the intrinsic volume V(int) of complex solute molecules without any additional assumptions. For comparison, we also use fused sphere calculations, which deliver van der Waals and solute accessible surface volumes as special cases. The apparent volume V(app) of the solute molecules is calculated by different approaches, using either a traditional distance based selection of hydration water or the construction of sequential Voronoi shells. We find an astonishing coincidence with the predictions of a simple empirical approach, which is based on experimentally determined amino acid side chain contributions (Biophys. Chem.1999, 82, 35). The intrinsic volumes of the polypeptides are larger than their apparent volumes and also increase with temperature. This is due to a negative contribution of the hydration water ΔV to the apparent volume. The absolute value of this contribution is less than 10% of the intrinsic volume for both molecules and decreases with temperature. Essential volumetric differences between hydration water and bulk water are observed in the nearest neighborhood of the solute only, practically in the first two Delaunay sublayers of the first Voronoi shell. This also helps to understand the pressure dependence of the partial molar volumes of proteins.

Temperature dependence of the solubility of carbon dioxide in imidazolium-based ionic liquids.

The solubility of carbon dioxide in ionic liquids of type 1-alkyl-3-methylimidazolium bis(trifluoromethylsulfonyl)imide ([Cnmim][NTf2]) with varying chain length n=2, 4, 6, 8 is computed from molecular dynamics simulations. By applying both Bennett's overlapping distribution method and Widom's particle insertion technique, we determine solvation free energies that are in excellent agreement with available experimental solubility data over a large temperature range from 300 to 500 K. We find that the computed solvation free energy of carbon dioxide is remarkably insensitive to the alkane chain length, emphasizing the importance of solvent models with accurate volumetric properties. The simulations suggest that the "anomalous" temperature dependence of the CO2 solvation at infinite dilution is characterized by counter-compensating negative entropies and enthalpies of solvation. By systematically varying the interaction strength of CO2 with the solvent, we show that the negative solvation entropy of CO2 is not caused by solvation cavities, but enforced by Coulomb and van der Waals interactions. We observe that solvation free energies and enthalpies obtained for models with different solute-solvent interaction strengths are subject to a linear correlation, similar to an expression that has been suggested for gases in polymers. Despite the apparent chain length insensitivity of the solvation free energy, significant changes in the solvation shell of a CO2 molecule are observed. The chain length insensitivity is found to be a consequence of two counter-compensating effects: the increasing free energy of cavity formation is balanced by a favorable interaction of CO2 with the alkyl chain of the imidazolium cation.

Temperature and concentration effects on the solvophobic solvation of methane in aqueous salt solutions.

We perform molecular dynamics (MD) simulations of aqueous salt (NaCl) solutions using the TIP4P-Ew water model (Horn et al., J. Chem. Phys. 2004, 120, 9665) covering broad temperature and concentration ranges extending deeply into the supercooled region. In particular we study the effect of temperature and salt concentration on the solvation of methane at infinite dilution. The salt effect on methane's solvation free energy, solvation enthalpy and entropy, as well as their temperature dependence is found to be semi-quantitatively in accordance with the data of Ben-Naim and Yaacobi (J. Phys. Chem. 1974, 78, 170). To distinguish the influence of local (in close proximity to ions) and global effects, we partition the salt solutions into ion influenced hydration shell regions and bulk water. The chemical potential of methane is systematically affected by the presence of salt in both sub volumes, emphasizing the importance of the global volume contraction due to electrostriction effects. This observation is correlated with systematic structural alterations similar to water under pressure. The observed electrostriction effects are found to become increasingly pronounced under cold (supercooled) conditions. We find that the influence of temperature and salt induced global density changes on the solvation properties of methane is well recovered by simple scaling relation based on predictions of the information theory model of Garde et al. (Phys. Rev. Let. 1999, 77, 4966).

Thermodynamic and structural characterization of the transformation from a metastable low-density to a very high-density form of supercooled TIP4P-Ew model water.

We explore the phase diagram of the metastable TIP4P-Ew liquid model water from 360 K down to 150 K at densities ranging from 0.950 to 1.355 g cm(-3). In addition to the low-density/high-density (LDL/HDL) liquid-liquid transition, we observe a structural high-density/very high-density (HDL/VHDL) transformation for the lowest temperatures at 1.30 g cm(-3). The characteristics of the isobars and isotherms suggest the presence of a stepwise HDL/VHDL transition with first-order-like appearance. In addition, we also identify an apparent pretransition at 1.24 g cm(-3), which suggests that the experimentally detected LDA/VHDA transformation might evolve into a multiple-step process with different local structures representing local minima in the free-energy landscape. Such a scenario is supported by a pronounced correlation between the isothermal density dependence of the pressure, with a stepwise increase of the oxygen coordination number, due to the appearance of interstitial water molecules.

The solvent-dependent shift of the amide I band of a fully solvated peptide as a local probe for the solvent composition in the peptide/solvent interface.

We determine the shift and line shape of the amide I band of a model AK peptide from molecular dynamics (MD) simulations of the peptide dissolved in methanol/water mixtures with varying composition. The IR spectra are determined from a transition dipole coupling exciton model. A simplified empirical model Hamiltonian is employed, which takes into account both the effect of hydrogen bonding and the intramolecular vibrational coupling. We consider a single isolated AK peptide in a mostly helical conformation, while the solvent is represented by 2600 methanol or water molecules, simulated for a pressure of 1 bar and a temperature of 300 K. Over the course of the simulations, minor reversible conformational changes at the termini are observed, which are found to only slightly affect the calculated spectral properties. Over the entire composition range, which varies from pure water to the pure methanol solvent, a monotonous shift towards higher frequency of the IR amide I band of about 8 wavenumbers is observed. This shift towards higher frequency is comparable to the shift found in preliminary experimental data also presented here on the amide I' band. The shift is found to be caused by two counter-compensating effects. An intramolecular red shift of about 1.2 wavenumbers occurs, due to stronger intramolecular hydrogen bonding in a methanol-rich environment. Dominating, however, is the intermolecular solvent-dependent shift towards higher frequency of about 10 wavenumbers, which is attributed to the less effective hydrogen-bond-donor capabilities of methanol compared to water. The importance of the solvent contribution to the IR shift, as well as the significantly different hydrogen formation capabilities of water and methanol, makes the amide I band sensitive to composition changes in the local environment close to the peptide/solvent interface. This allows, in principle, an experimental determination of the composition of the solvent in close proximity to the peptide surface. For the AK peptide case, we observe at low methanol concentrations a significantly enhanced methanol concentration at the peptide/solvent interface, supposedly promoted by the partially hydrophobic character of the AK peptide's solvent-accessible surface.

Modeling of aqueous poly(oxyethylene) solutions. 2. Mesoscale simulations.

We extend our work on aqueous solutions of poly(oxyethylene) oligomers H-(CH2-O-CH2)n -H (POEn). On the basis of atomistic simulations of trimer and decamer solutions (first part of this series of papers), different sets of coarse-grained implicit-solvent potentials have been constructed using the iterative Boltzmann inversion technique. The comparison of structures obtained from coarse-grained simulations (gyration radii, end-to-end distances, radial distribution functions) with atomistic reference simulations and experiments shows that the state-specific potentials are transferable both to a wide concentration range, if the same molecule size is considered, and to at least 2 orders of magnitude larger molecules (in terms of molecular mass). Comparing the performance of different mesoscale potentials, we find different applicability ranges in terms of molecule sizes. The experimental gyration radii for chains comprising up to 1500 monomers are reproduced almost quantitatively by the decamer-fitted potentials with dihedral interactions included. The trimer-fitted potentials reproduce experimental chain dimensions of up to some hundred monomers but seem to become metastable beyond a certain chain length, as we evidenced some chain collapses. Relaxation of large-scale features is 1-2 orders of magnitude faster in the mesoscale simulations than in the atomistic simulations. The diffusion behavior in dependence of concentration is captured correctly when the decamer potential is applied to the decamer itself. For all other chain lengths, we find that time mapping from coarse-grained to atomistic trajectories has to be determined separately for each concentration. Overall, diffusion is 1-2 orders of magnitude faster on the mesoscale, depending considerably on the Lowe-Andersen thermostat parameters. The CG simulations provide an overall speed-up of about 3 orders of magnitude.

Modeling of aqueous poly(oxyethylene) solutions: 1. Atomistic simulations.

The performance of different recently proposed force fields in combination with TIP4P-Ewald (TIP4P-Ew) water in reproducing experimental data of liquid 1,2-dimethoxyethane (DME) and its aqueous solutions for conformer populations, densities of solutions, and self-diffusion coefficients was explored. A modified version of the OPLS force field ("engineered") showed best performance in describing the conformer equilibria, but extremely high interconformational barriers reduce its applicability in dynamical simulations. The TraPPE-united atom force field (TraPPE-UA) by Siepmann et al. (J. Phys. Chem. B 2004, 108, 17596) was found to perform best in reproducing thermodynamic properties, but it showed some deficiency in describing the conformer equilibria. We reparameterized the dihedral potentials to match recent ab initio data by Anderson and Wilson (Mol. Phys. 2005, 103, 89) and could improve significantly the performance of description of conformer populations of DME in water. Subsequently, this modified TraPPE-UA was used in extensive simulations of poly(oxyethylene) oligomers H(CH2OCH2)nH (POEn) with n = 3, 5, 10, 12, 20, 30 repeat units at mass fractions between 3% and 80% at 298 K. Density, radii of gyration, and diffusion coefficients are in very good agreement with available experimental data. We conclude that this force field in combination with the TIP4P-Ew water model is very suitable for simulations of poly(oxyethylene) oligomers in aqueous solution. The application to real polymeric systems on the atomistic level is however hindered by very slow decorrelation of large-scale features and by slow diffusion.

Temperature-induced conformational transition of a model elastin-like peptide GVG(VPGVG)(3) in water.

The conformation of a single elastin-like peptide GVG(VPGVG)3 in liquid water is studied by computer simulations in the temperature interval between 280 and 440 K. Two main conformational states of the peptide can be distinguished: a rigid conformational state, dominating at low temperatures, and a flexible conformational state, dominating at high temperatures. A temperature-induced transition between these states occurs at about 310 K, rather close to a transition temperature seen in experiments. This transition is accompanied by the thermal breaking of the hydrogen-bonded spanning network of the hydration water via a percolation transition upon heating. This finding indicates that the H-bond clustering structure of the hydration water plays an important role in the conformational stability of biomolecules. A second important observation is the Gaussian distribution of the end-to-end distance in the high-temperature state, which supports the idea of a rubber-like elasticity of the studied elastin-like peptide. Finally our results challenge the idea of the folding of elastin-like peptides upon heating.

Adding salt to an aqueous solution of t-butanol: is hydrophobic association enhanced or reduced?

Recent neutron scattering experiments on aqueous salt solutions of amphiphilic t-butanol by Bowron and Finney [Phys. Rev. Lett. 89, 215508 (2002); J. Chem. Phys. 118, 8357 (2003)] suggest the formation of t-butanol pairs, bridged by a chloride ion via O-H...Cl- hydrogen bonds, leading to a reduced number of intermolecular hydrophobic butanol-butanol contacts. Here we present a joint experimental/theoretical study on the same system, using a combination of molecular dynamics (MD) simulations and nuclear magnetic relaxation measurements. Both MD simulation and experiment clearly support the more classical scenario of an enhanced number of hydrophobic contacts in the presence of salt, as it would be expected for purely hydrophobic solutes. [T. Ghosh et al., J. Phys. Chem. B 107, 612 (2003)]. Although our conclusions arrive at a structurally completely distinct scenario, the molecular dynamics simulation results are within the experimental error bars of the Bowron and Finney data.

Percolation transition of hydration water: from planar hydrophilic surfaces to proteins.

The formation of a spanning hydrogen-bonded network of hydration water is found to occur via a 2D percolation transition in various systems: smooth hydrophilic surfaces, the surface of a single protein molecule, protein powder, and diluted peptide solution. The average number of water-water hydrogen bonds at the percolation threshold varies from 2.0 to 2.3, depending on temperature, system size, and surface properties. Calculation of nH allows an easy estimation of the percolation threshold of hydration water in various systems, including biomolecules.

Low-temperature and high-pressure induced swelling of a hydrophobic polymer-chain in aqueous solution.

We report molecular dynamics simulations of a hydrophobic polymer-chain in aqueous solution between 260 K and 420 K at pressures of 1 bar, 3000 bar, and 4500 bar. The simulations reveal a hydrophobically collapsed structure at low pressures and high temperatures. At 3000 bar and about 260 K and at 4500 bar and about 260 K, however, an abrupt transition to a swelled state is observed. The transition is driven by a smaller volume and a remarkably strong lower enthalpy of the swelled state, indicating a steep positive slope of the corresponding transition line. The swelling is strongly stabilized by the energetically favorable state of water in the polymer's hydrophobic first hydration shell at low temperatures. This finding is consistent with the observation of a positive heat capacity of hydrophobic solvation. Moreover, the slope and location of the estimated swelling transition line for the collapsed hydrophobic chain coincides remarkably well with the cold denaturation transition of proteins.

Liquid-liquid phase transitions in supercooled water studied by computer simulations of various water models.

Liquid-liquid and liquid-vapor coexistence regions of various water models were determined by Monte Carlo (MC) simulations of isotherms of density fluctuation-restricted systems and by Gibbs ensemble MC simulations. All studied water models show multiple liquid-liquid phase transitions in the supercooled region: we observe two transitions of the TIP4P, TIP5P, and SPCE models and three transitions of the ST2 model. The location of these phase transitions with respect to the liquid-vapor coexistence curve and the glass temperature is highly sensitive to the water model and its implementation. We suggest that the apparent thermodynamic singularity of real liquid water in the supercooled region at about 228 K is caused by an approach to the spinodal of the first (lowest density) liquid-liquid phase transition. The well-known density maximum of liquid water at 277 K is related to the second liquid-liquid phase transition, which is located at positive pressures with a critical point close to the maximum. A possible order parameter and the universality class of liquid-liquid phase transitions in one-component fluids are discussed.

Formation of spanning water networks on protein surfaces via 2D percolation transition.

The formation of spanning hydrogen-bonded water networks on protein surfaces by a percolation transition is closely connected with the onset of their biological activity. To analyze the structure of the hydration water at this important threshold, we performed the first computer simulation study of the percolation transition of water in a model protein powder and on the surface of a single protein molecule. The formation of an infinite water network in the protein powder occurs as a 2D percolation transition at a critical hydration level, which is close to the values observed experimentally. The formation of a spanning 2D water network on a single rigid protein molecule can be described by adapting the cluster analysis of conventional percolation studies to the characterization of the connectivity of the hydration water on the surface of finite objects. Strong fluctuations of the surface water network are observed close to the percolation threshold. Our simulations also furnish a microscopic picture for understanding the specific values of the experimentally observed hydration levels, where different steps of increasing mobility in the hydrated powder are observed.

Properties of spanning water networks at protein surfaces.

The formation of a spanning two-dimensional hydrogen-bonded water network at the surface of proteins via a percolation transition enables their biological function. We show in detail how the spanning (percolating) water network appears at the surfaces of model hydrophilic spheres and at the surface of a single protein (lysozyme) molecule. We have found essential correlations of the linear extension, radius of gyration, and position of the center of mass of the largest water cluster with its size. The specific two-peak structure of the probability distribution of the largest cluster size allowed us to study various properties separately for spanning and nonspanning largest clusters. The radius of gyration of the spanning cluster always exceeds the radii of the spheres or the effective radius of the protein. Any spanning cluster envelops essentially more than half of the surface area. The temporal decay of the spanning networks shows a stretched exponential character. Their average lifetime at the percolation threshold is about the lifetime of a water-water hydrogen bond.