Ion solvation thermodynamics from simulation with a polarizable force field

Thermodynamic measurements of the solvation of salts and electrolytes are relatively straightforward, but it is not possible to separate total solvation free energies into distinct cation and anion contributions without reference to an additional extrathermodynamic assumption. The present work attempts to resolve this difficulty using molecular dynamics simulations with the AMOEBA polarizable force field and perturbation techniques to directly compute absolute solvation free energies for potassium, sodium, and chloride ions in liquid water and formamide. Corresponding calculations are also performed with two widely used nonpolarizable force fields. The simulations with the polarizable force field accurately reproduce in vacuo quantum mechanical results, experimental ion-cluster solvation enthalpies, and experimental solvation free energies for whole salts, while the other force fields do not. The results indicate that calculations with a polarizable force field can capture the thermodynamics of ion solvation and that the solvation free energies of the individual ions differ by several kilocalories from commonly cited values.     

The effects of charge transfer on the aqueous solvation of ions

Ab initio-based charge partitioning of ionic systems results in ions with non-integer charges. This charge-transfer (CT) effect alters both short- and long-range interactions. Until recently, the effects of CT have been mostly neglected in molecular dynamics (MD) simulations. The method presented in this paper for including charge transfer between ions and water is consistent with ab initio charge partitioning and does not add significant time to the simulation. The ions of sodium, potassium, and chloride are parameterized to reproduce dimer properties and aqueous structures. The average charges of the ions from MD simulations (0.900, 0.919, and -0.775 for Na(+), K(+), and Cl(-), respectively) are consistent with quantum calculations. The hydration free energies calculated for these ions are in agreement with experimental estimates, which shows that the interactions are described accurately. The ions also have diffusion constants in good agreement with experiment. Inclusion of CT results in interesting properties for the waters in the first solvation shell of the ions. For all ions studied, the first shell waters acquire a partial negative charge, due to the difference between water-water and water-ion charge-transfer amounts. CT also reduces asymmetry in the solvation shell of the chloride anion, which could have important consequences for the behavior of chloride near the air-water interface.     

The effects of dissolved halide anions on hydrogen bonding in liquid water

It is widely believed that the addition of salts to water engenders structural changes in the hydrogen-bond network well beyond the adjacent shell of solvating molecules. Classification of many ions as "structure makers" and "structure breakers" has been based in part on corresponding changes in the vibrational spectra (Raman and IR). Here we show that changes in O-H vibrational spectra induced by the alkali halides in liquid water result instead from the actions of ions' electric fields on adjacent water molecules. Computer simulations that accurately reproduce our experimental measurements suggest that the statistics of hydrogen-bond strengths are only weakly modified beyond this first solvation shell.     

Determination of alkali and halide monovalent ion parameters for use in explicitly solvated biomolecular simulations

Alkali (Li(+), Na(+), K(+), Rb(+), and Cs(+)) and halide (F(-), Cl(-), Br(-), and I(-)) ions play an important role in many biological phenomena, roles that range from stabilization of biomolecular structure, to influence on biomolecular dynamics, to key physiological influence on homeostasis and signaling. To properly model ionic interaction and stability in atomistic simulations of biomolecular structure, dynamics, folding, catalysis, and function, an accurate model or representation of the monovalent ions is critically necessary. A good model needs to simultaneously reproduce many properties of ions, including their structure, dynamics, solvation, and moreover both the interactions of these ions with each other in the crystal and in solution and the interactions of ions with other molecules. At present, the best force fields for biomolecules employ a simple additive, nonpolarizable, and pairwise potential for atomic interaction. In this work, we describe our efforts to build better models of the monovalent ions within the pairwise Coulombic and 6-12 Lennard-Jones framework, where the models are tuned to balance crystal and solution properties in Ewald simulations with specific choices of well-known water models. Although it has been clearly demonstrated that truly accurate treatments of ions will require inclusion of nonadditivity and polarizability (particularly with the anions) and ultimately even a quantum mechanical treatment, our goal was to simply push the limits of the additive treatments to see if a balanced model could be created. The applied methodology is general and can be extended to other ions and to polarizable force-field models. Our starting point centered on observations from long simulations of biomolecules in salt solution with the AMBER force fields where salt crystals formed well below their solubility limit. The likely cause of the artifact in the AMBER parameters relates to the naive mixing of the Smith and Dang chloride parameters with AMBER-adapted Aqvist cation parameters. To provide a more appropriate balance, we reoptimized the parameters of the Lennard-Jones potential for the ions and specific choices of water models. To validate and optimize the parameters, we calculated hydration free energies of the solvated ions and also lattice energies (LE) and lattice constants (LC) of alkali halide salt crystals. This is the first effort that systematically scans across the Lennard-Jones space (well depth and radius) while balancing ion properties like LE and LC across all pair combinations of the alkali ions and halide ions. The optimization across the entire monovalent series avoids systematic deviations. The ion parameters developed, optimized, and characterized were targeted for use with some of the most commonly used rigid and nonpolarizable water models, specifically TIP3P, TIP4P EW, and SPC/E. In addition to well reproducing the solution and crystal properties, the new ion parameters well reproduce binding energies of the ions to water and the radii of the first hydration shells.     

New evaluation of reconstructed spatial distribution function from radial distribution functions

Although three dimensional (3D) solvation structure is much more informative than one dimensional structure, its evaluation is difficult experimentally and theoretically. In our previous Communication [Yokogawa et al., J. Chem. Phys. 123, 211102 (2005)], we proposed a new method to present reconstructed spatial distribution function (RC-SDF) from a set of radial distribution functions (RDFs). In this article, we successfully extended the method more accurately with new basis sets. This new method was applied to two liquid solvation structures, methanol and dimethyl sulfoxide, as examples. Their RC-SDFs evaluated here clearly show that the former solvation structure is well defined while the latter one is broad, which agrees well with the SDFs calculated directly from molecular dynamics simulations. These results indicate that the method can reproduce well these 3D solvation structures in reasonable computational cost.     

Continuum solvation models in the linear interaction energy method

The linear interaction energy (LIE) method in combination with two different continuum solvent models has been applied to calculate protein-ligand binding free energies for a set of inhibitors against the malarial aspartic protease plasmepsin II. Ligand-water interaction energies are calculated from both Poisson-Boltzmann (PB) and Generalized Born (GB) continuum models using snapshots from explicit solvent simulations of the ligand and protein-ligand complex. These are compared to explicit solvent calculations, and we find close agreement between the explicit water and PB solvation models. The GB model overestimates the change in solvation energy, and this is caused by consistent underestimation of the effective Born radii in the protein-ligand complex. The explicit solvent LIE calculations and LIE-PB, with our standard parametrization, reproduce absolute experimental binding free energies with an average unsigned error of 0.5 and 0.7 kcal/mol, respectively. The LIE-GB method, however, requires a constant offset to approach the same level of accuracy.     

Shared solvation of sodium ions in alcohol-water solutions explains the non-ideality of free energy of solvation

In order to explain the discrepancies between theories and experiments regarding the non-ideality in the free energy of solvation, here we present a microscopic picture of sodium ions dissolved in water-alcohol mixed solvents. We used X-ray absorption spectroscopy to probe the K-edge of sodium ions in mixed solvents of water and alcohols (methanol, ethanol) and in the respective pure solvents. In the mixed solvents a shared solvation of the sodium ions is observed. We find that specifically the water component plays a key role in stabilizing the solvation shell in mixed solvents, which was revealed by a selective photochemical process occurring only in the pure alcohol solvents.     

MST continuum study of the hydration free energies of monovalent ionic species

In this study, we revisit the protocol previously proposed within the framework of the Miertus-Scrocco-Tomasi (MST) continuum model to define the cavity between the solute and solvent for predicting hydration free energies of univalent ions. The protocol relies on the use of a reduced cavity (around 10-15% smaller than the cavity used for neutral compounds) around the atom(s) bearing the formal charge. The suitability of this approach is examined here for a series of 47 univalent ions for which accurate experimental hydration free energies are available. Attention is also paid to the effect of the charge renormalization protocol used to correct uncertainties arising from the electron density located outside the solute cavity. The method presented here provides, with a minimum number of fitted parameters, reasonable estimates within the experimental error of the hydration free energy of ions (average relative error of 4.7%) and is able to reproduce solvation in water of both small and large ions.     

Comparison of two simulation methods to compute solvation free energies and partition coefficients

The thermodynamic integration (TI) and expanded ensemble (EE) methods are used here to calculate the hydration free energy in water, the solvation free energy in 1‐octanol, and the octanol‐water partition coefficient for a six compounds of varying functionality using the optimized potentials for liquid simulations (OPLS) all‐atom (AA) force field parameters and atomic charges. Both methods use the molecular dynamics algorithm as a primary component of the simulation protocol, and both have found wide applications in fields such as the calculation of activity coefficients, phase behavior, and partition coefficients. Both methods result in solvation free energies and 1‐octanol/water partition coefficients with average absolute deviations (AAD) from experimental data to within 4 kJ/mol and 0.5 log units, respectively. Here, we find that in simulations the OPLS‐AA force field parameters (with fixed charges) can reproduce solvation free energies of solutes in 1‐octanol with AAD of about half that for the solute hydration free energies using a extended simple point charge (SPC/E) model of water. The computational efficiency of the two simulation methods are compared based on the time (in nanoseconds) required to obtain similar standard deviations in the solvation free energies and 1‐octanol/water partition coefficients. By this analysis, the EE method is found to be a factor of nine more efficient than the TI algorithm. For both methods, solvation free energy calculations in 1‐octanol consume roughly an order of magnitude more CPU hours than the hydration free energy calculations. © 2012 Wiley Periodicals, Inc.     

Absolute solvation free energy of Li+ and Na+ ions in dimethyl sulfoxide solution: a theoretical ab initio and cluster-continuum model study

The solvation of the lithium and sodium ions in dimethyl sulfoxide solution was theoretically investigated using ab initio calculations coupled with the hybrid cluster-continuum model, a quasichemical theory of solvation. We have investigated clusters of ions with up to five dimethyl sulfoxide (DMSO) molecules, and the bulk solvent was described by a dielectric continuum model. Our results show that the lithium and sodium ions have four and five DMSO molecules into the first coordination shell, and the calculated solvation free energies are -135.5 and -108.6 kcal mol(-1), respectively. These data suggest a solvation free energy value of -273.2 kcal mol(-1) for the proton in dimethyl sulfoxide solution, a value that is more negative than the present uncertain experimental value. This and previous studies on the solvation of ions in water solution indicate that the tetraphenylarsonium tetraphenylborate assumption is flawed and the absolute value of the free energy of transfer of ions from water to DMSO solution is higher than the present experimental values.     

Free energy of solvation of simple ions: molecular-dynamics study of solvation of Cl- and Na+ in the ice/water interface

Molecular-dynamics simulations of Cl(-) and Na(+) ions are performed to calculate ionic solvation free energies in both bulk simple point-charge/extended water and ice 1 h at several different temperatures, and at the basal ice 1 h/water interface. For the interface we calculate the free energy of "transfer" of the ions across the ice/water interface. For the ions in bulk water in the NPT ensemble at 298 K and 1 atm, results are found to be in good agreement with experiments, and with other simulation results. Simulations performed in the NVT ensemble are shown to give equivalent solvation free energies, and this ensemble is used for the interfacial simulations. Solvation free energies of Cl(-) and Na(+) ions in ice at 150 K are found to be approximately 30 and approximately 20 kcal mol(-1), respectively, less favorable than for water at room temperature. Near the melting point of the model the solvation of the ions in water is the same (within statistical error) as that measured at room temperature, and in the ice is equivalent and approximately 10 kcal mol(-1) less favorable than the liquid. The free energy of transfer for each ion across ice/water interface is calculated and is in good agreement with the bulk observations for the Cl(-) ion. However, for the model of Na(+) the long-range electrostatic contribution to the free energy was more negative in the ice than the liquid, in contrast with the results observed in the bulk calculations.     

Structural properties of water: comparison of the SPC, SPCE, TIP4P, and TIP5P models of water

Molecular-dynamics simulations were carried out for the SPC, SPCE, TIP4P, and TIP5P models of water at 298 K. From these results we determine the following quantities: the absolute entropy using the two-particle approximation, the mean lifetime of the hydrogen bond, the mean number of hydrogen bonds per molecule, and the mean energy of the hydrogen bond. From the entropy calculations we find that nearly all contributions to the total entropy originates from the orientation effects. Moreover, we determine the contributions to the total entropy which originate from the first, second, and higher solvation shells. It is interesting that the limits between solvation shells are clearly visible. The first solvation shell (0.22 < r < 0.36 nm) contributes approximately 43 J mol K to the total entropy; the second solvation shell (0.36 < r < 0.60 nm) contributes approximately 12 J mol K, while contributions from the third and other solvation shells are very small, approximately 2 J mol K in summary. This indicates that water molecules are strongly ordered up to 0.55-0.6 nm around the central water molecule, and beyond this limit the ordering diminishes. The results of calculations (entropy and hydrogen bonds) are compared with the experimental data for the choosing of the best water model. We find that the SPC and TIP4P models reproduce the best experimental values, and we recommend these models for computer simulations of the aqueous solution of biomolecules.     

Dependence of ion hydration on the sign of the ion's charge

The solvation of simple ions in water is studied using molecular dynamics simulations with a polarizable force field. Previous simulations using this potential demonstrated that anions are more favorably solvated in water than cations. The present work is an attempt to explain this result by examining the effects of ions on the surrounding water structure, with particular focus on the first solvation shell and its interactions with the surrounding water. We conclude that while the first solvation shell surrounding cations is frustrated by competition between ion-water and water-water interactions, solvation of anions is compatible with good water-water interactions.     

How ions affect the structure of water

We model ion solvation in water. We use the MB model of water, a simple two-dimensional statistical mechanical model in which waters are represented as Lennard-Jones disks having Gaussian hydrogen-bonding arms. We introduce a charge dipole into MB waters. We perform (NPT) Monte Carlo simulations to explore how water molecules are organized around ions and around nonpolar solutes in salt solutions. The model gives good qualitative agreement with experiments, including Jones-Dole viscosity B coefficients, Samoilov and Hirata ion hydration activation energies, ion solvation thermodynamics, and Setschenow coefficients for Hofmeister series ions, which describe the salt concentration dependence of the solubilities of hydrophobic solutes. The two main ideas captured here are (1) that charge densities govern the interactions of ions with water, and (2) that a balance of forces determines water structure: electrostatics (water's dipole interacting with ions) and hydrogen bonding (water interacting with neighboring waters). Small ions (kosmotropes) have high charge densities so they cause strong electrostatic ordering of nearby waters, breaking hydrogen bonds. In contrast, large ions (chaotropes) have low charge densities, and surrounding water molecules are largely hydrogen bonded.     

Predicting hydration free energies of polychlorinated aromatic compounds from the SAMPL-3 data set with FiSH and LIE models

Next-generation solvation models are devised to mimic the accuracy and generality of explicit solvation models at the speed of current popular implicit solvation models. One such method is the first-shell of hydration (FiSH) continuum model that was trained on hydration energetics from LIE calculations and molecular dynamics simulations in explicit solvent. Here we tested prospectively the FiSH model on the SAMPL-3 hydration data set that zooms in the effect of chlorination on solvation. We compare these FiSH predictions with those from retrospective LIE calculations. We find that neither FiSH nor LIE can reproduce well the absolute values and the trend of hydration free energies in the biphenyl and dioxin aromatic chlorination series. Some of the hypotheses behind this performance are discussed and tested. The LIE explicit-solvent model shows some improvement relative to the FiSH continuum model, and we correct a systematic deviation in the continuum van der Waals term of FiSH associated with aromatic Cl atom type.     

A new set of atomic radii for accurate estimation of solvation free energy by Poisson–Boltzmann solvent model

The Poisson–Boltzmann implicit solvent (PB) is widely used to estimate the solvation free energies of biomolecules in molecular simulations. An optimized set of atomic radii (PB radii) is an important parameter for PB calculations, which determines the distribution of dielectric constants around the solute. We here present new PB radii for the AMBER protein force field to accurately reproduce the solvation free energies obtained from explicit solvent simulations. The presented PB radii were optimized using results from explicit solvent simulations of the large systems. In addition, we discriminated PB radii for N‐ and C‐terminal residues from those for nonterminal residues. The performances using our PB radii showed high accuracy for the estimation of solvation free energies at the level of the molecular fragment. The obtained PB radii are effective for the detailed analysis of the solvation effects of biomolecules. © 2014 The Authors Journal of Computational Chemistry Published by Wiley Periodicals, Inc.     

Structural transitions in ion coordination driven by changes in competition for ligand binding

Transferring Na(+) and K(+) ions from their preferred coordination states in water to states having different coordination numbers incurs a free energy cost. In several examples in nature, however, these ions readily partition from aqueous-phase coordination states into spatial regions having much higher coordination numbers. Here we utilize statistical theory of solutions, quantum chemical simulations, classical mechanics simulations, and structural informatics to understand this aspect of ion partitioning. Our studies lead to the identification of a specific role of the solvation environment in driving transitions in ion coordination structures. Although ion solvation in liquid media is an exergonic reaction overall, we find it is also associated with considerable free energy penalties for extracting ligands from their solvation environments to form coordinated ion complexes. Reducing these penalties increases the stabilities of higher-order coordinations and brings down the energetic cost to partition ions from water into overcoordinated binding sites in biomolecules. These penalties can be lowered via a reduction in direct favorable interactions of the coordinating ligands with all atoms other than the ions themselves. A significant reduction in these penalties can, in fact, also drive up ion coordination preferences. Similarly, an increase in these penalties can lower ion coordination preferences, akin to a Hofmeister effect. Since such structural transitions are effected by the properties of the solvation phase, we anticipate that they will also occur for other ions. The influence of other factors, including ligand density, ligand chemistry, and temperature, on the stabilities of ion coordination structures are also explored.