Experimental Hydration Free Energies Research Articles

Development of polarizable and hydration-focused water models for the Martini 3 force field

Coarse-grained simulations have been increasingly employed in many fields of molecular chemistry. These models enable the evaluation of complex molecular structures by simplifying their description compared to all-atoms models. Specifically, the Martini 3 force field has been developed to describe biomolecular structures such as proteins, lipids, and carbohydrates. Their “building block” approach has been successfully used in modeling complex systems ranging from carbohydrate solutions, biological membranes, and even deep eutectic solvents. However, one main issue arises when considering interfacial systems: water is represented as a single particle with no effect of electrostatic interactions. In this work, we developed a polarizable water model for the Martini 3 force field, introducing charged particles to model the electrostatic interactions present in water systems accurately. Gaussian Processes were trained with simulation data for both bulk water and cross-interaction parameterization to consistently explore the parameter space and fine-tune our model to the target properties. The developed model showed good performance in reproducing the density and dielectric constant of bulk water. Regarding cross-interaction with other beads, two different sets of parameters were developed. The first, called the polarizable model, is parameterized to reproduce the hydration free energies of the original Martini 3 organic beads. This allows for a good description of free energies of transfer between water and organic phases, the main target property for the Martini 3 force field. The second, called the polarizable-hydration model, is parameterized to match the experimental hydration free energy of selected molecules. This allows for more accurate evaluation of hydration properties, which can be useful for water solutions. In particular, we show that the partitioning and interfacial behavior of nonionic surfactants is improved by the polarizable-hydration model when compared to the standard Martini 3 water model.

Inclusion of Water Multipoles into the Implicit Solvation Framework Leads to Accuracy Gains.

The current practical "workhorses" of the atomistic implicit solvation─the Poisson-Boltzmann (PB) and generalized Born (GB) models─face fundamental accuracy limitations. Here, we propose a computationally efficient implicit solvation framework, the Implicit Water Multipole GB (IWM-GB) model, that systematically incorporates the effects of multipole moments of water molecules in the first hydration shell of a solute, beyond the dipole water polarization already present at the PB/GB level. The framework explicitly accounts for coupling between polar and nonpolar contributions to the total solvation energy, which is missing from many implicit solvation models. An implementation of the framework, utilizing the GAFF force field and AM1-BCC atomic partial charges model, is parametrized and tested against the experimental hydration free energies of small molecules from the FreeSolv database. The resulting accuracy on the test set (RMSE ∼ 0.9 kcal/mol) is 12% better than that of the explicit solvation (TIP3P) treatment, which is orders of magnitude slower. We also find that the coupling between polar and nonpolar parts of the solvation free energy is essential to ensuring that several features of the IWM-GB model are physically meaningful, including the sign of the nonpolar contributions.

Fast calculation of hydrogen-bond strengths and free energy of hydration of small molecules

Hydrogen bonding is an interaction of great importance in drug discovery and development as it may significantly affect chemical and biological processes including the interaction of small molecules with other molecules, proteins, and membranes. In particular, hydrogen bonding can impact drug-like properties such as target affinity and oral availability which are critical to developing effective pharmaceuticals, and therefore, numerous methods for the calculation of properties such as hydrogen-bond strengths, free energy of hydration, or water solubility have been proposed over time. However, the accessibility to efficient methods for the predictions of such properties is still limited. Here, we present the development of Jazzy, an open-source tool for the prediction of hydrogen-bond strengths and free energies of hydration of small molecules. Jazzy also allows the visualisation of hydrogen-bond strengths with atomistic resolution to support the design of compounds with desired properties and the interpretation of existing data. The tool is described in its implementation, parameter fitting, and validation against two data sets of experimental hydration free energies. Jazzy is also applied against two chemical series of bioactive compounds to show that hydrogen-bond strengths can be used to understand their structure–activity relationships. Results from the validations highlight the strengths and limitations of Jazzy, and suggest its suitability for interactive design, screening, and machine-learning featurisation.

Machine Learning Prediction of Hydration Free Energy with Physically Inspired Descriptors

We present machine learning models for predicting experimental hydration free energies of molecules without any atom-, bond-, or geometry-specific input feature. Four types of physically inspired descriptors are adopted for predictions. The first type is composed of the total dipole moment, anisotropic polarizability, and vibrational analysis results of the solute molecule. The second and third types are derived from the electrostatic potential distribution of the solute. The last type includes the solvent accessible surface area and shape similarities. Several machine learning regression models are built on the basis of the FreeSolv database with ∼600 samples, showing a better performance in comparison with that of most traditional approaches and other prediction methods based on molecular fingerprints. In particular, the present descriptors are capable of predicting hydration free energies of new compounds with elements or fragments that are never seen in the training set. The importance of these descriptors, the impact of dissociation energies of specific covalent bonds, and the outliers with relatively large prediction errors are also discussed.

Development of AMOEBA Polarizable Force Field for Rare-Earth La3+ Interaction with Bioinspired Ligands.

Rare-earth metals (REMs) are crucial for many important industries, such as power generation and storage, in addition to cancer treatment and medical imaging. One promising new REM refinement approach involves mimicking the highly selective and efficient binding of REMs observed in relatively recently discovered proteins. However, realizing any such bioinspired approach requires an understanding of the biological recognition mechanisms. Here, we developed a new classical polarizable force field based on the AMOEBA framework for modeling a lanthanum ion (La3+) interacting with water, acetate, and acetamide, which have been found to coordinate the ion in proteins. The parameters were derived by comparing to high-level ab initio quantum mechanical (QM) calculations that include relativistic effects. The AMOEBA model, with advanced atomic multipoles and electronic polarization, is successful in capturing both the QM distance-dependent La3+-ligand interaction energies and experimental hydration free energy. A new scheme for pairwise polarization damping (POLPAIR) was developed to describe the polarization energy in La3+ interactions with both charged and neutral ligands. Simulations of La3+ in water showed water coordination numbers and ion-water distances consistent with previous experimental and theoretical findings. Water residence time analysis revealed both fast and slow kinetics in water exchange around the ion. This new model will allow investigation of fully solvated lanthanum ion-protein systems using GPU-accelerated dynamics simulations to gain insights on binding selectivity, which may be applied to the design of synthetic analogues.

Multitask Deep Ensemble Prediction of Molecular Energetics in Solution: From Quantum Mechanics to Experimental Properties.

The past few years have witnessed significant advances in developing machine learning methods for molecular energetics predictions, including calculated electronic energies with high-level quantum mechanical methods and experimental properties, such as solvation free energy and logP. Typically, task-specific machine learning models are developed for distinct prediction tasks. In this work, we present a multitask deep ensemble model, sPhysNet-MT-ens5, which can simultaneously and accurately predict electronic energies of molecules in gas, water, and octanol phases, as well as transfer free energies at both calculated and experimental levels. On the calculated data set Frag20-solv-678k, which is developed in this work and contains 678,916 molecular conformations, up to 20 heavy atoms, and their properties calculated at B3LYP/6-31G* level of theory with continuum solvent models, sPhysNet-MT-ens5 predicts density functional theory (DFT)-level electronic energies directly from force field-optimized geometry within chemical accuracy. On the experimental data sets, sPhysNet-MT-ens5 achieves state-of-the-art performances, which predict both experimental hydration free energy with a RMSE of 0.620 kcal/mol on the FreeSolv data set and experimental logP with a RMSE of 0.393 on the PHYSPROP data set. Furthermore, sPhysNet-MT-ens5 also provides a reasonable estimation of model uncertainty which shows correlations with prediction error. Finally, by analyzing the atomic contributions of its predictions, we find that the developed deep learning model is aware of the chemical environment of each atom by assigning reasonable atomic contributions consistent with our chemical knowledge.

A Polarizable Cationic Dummy Metal Ion Model.

A novel locally polarizable multisite model based on the original cation dummy atom (CDA) model is described for molecular dynamics simulations of ions in condensed phases. Polarization effects are introduced by the electronegativity equalization model (EEM) method where charges on the metal ion and its dummy atoms can fluctuate to respond to the environment. This model includes explicit polarization and ion-induced interactions and can be coupled with nonpolarizable or polarizable water models, making it more transferable to simpler force fields. This approach allows us to enhance the original fixed charge CDA model where the charge distribution cannot adapt to the local solvent structure. To illustrate the new CDApol model, we examined properties of the Zn2+, Al3+, and Zr4+ ions in aqueous solution. The polarizable model and Lennard-Jones parameters were refined for octahedrally coordinated Zn2+, Al3+, and Zr4+ CDAs to reproduce thermodynamic and geometrical properties. Using this locally polarizable model, we were able to obtain the experimental hydration free energy, ion-oxygen distance, and coordination number coupled with the standard 12-6 Lennard-Jones model. This model can be used in myriad additional applications where local polarization and charge transfer is important.

Harnessing Deep Learning for Optimization of Lennard-Jones Parameters for the Polarizable Classical Drude Oscillator Force Field.

The outcomes of computational chemistry and biology research, including drug design, are significantly influenced by the underlying force field (FF) used in molecular simulations. While improved FF accuracy may be achieved via inclusion of explicit treatment of electronic polarization, such an extension must be accompanied by optimization of van der Waals (vdW) interactions, in the context of the Lennard-Jones (LJ) formalism in the present study. This is particularly challenging due to the extensive nature of chemical space combined with the correlated nature of LJ parameters. To address this challenge, a deep learning (DL)-based parametrization framework is developed, allowing for sampling of wide ranges of LJ parameters targeting experimental condensed phase thermodynamic properties. The present work utilizes this framework to develop the LJ parameters for atoms associated with four distinct groups covering 10 different atom types. Final parameter selection was facilitated by quantum mechanical data on rare-gas interactions with the training set molecules. The chosen parameters were then validated through experimental hydration free energies and condensed phase thermodynamic properties of validation set molecules to confirm transferability. The ultimate outcome of utilizing this framework is a set of LJ parameters in the context of the polarizable Drude FF, which demonstrated improvement in the reproduction of both experimental pure solvent and crystal properties and hydration free energies of the molecules compared to the additive CHARMM General FF (CGenFF) including the ability of the Drude FF to accurately reproduce both experimental pure solvent properties and hydration free energies. The study also shows how correlations between difference in the reproduction of condensed phase data between model compounds may be used to direct the selection of new atom types and training set molecules during FF development.

Rational Design of Nonbonded Point Charge Models for Divalent Metal Cations with Lennard-Jones 12-6 Potential.

Exploring a metal-involved biochemical process at a molecular level often requires a reliable description of metal properties in aqueous solution by classical nonbonded models. An additional C4 term for considering ion-induced dipole interactions was previously proposed to supplement the widely used Lennard-Jones 12-6 potential (known as the 12-6-4 LJ-type model) with good accuracy. Here, we demonstrate an alternative to modeling divalent metal cations (M2+) with the traditional 12-6 LJ potential by developing nonbonded point charge models for use with 11 water models: TIP3P, SPC/E, SPC/Eb, TIP4P-Ew, TIP4P-D, and TIP4P/2005 and the more recent OPC3, TIP3P-FB, OPC, TIP4P-FB, and a99SB-disp. Our designed models simultaneously reproduce the experimental hydration free energy, ion-oxygen distance, and coordination number in the first hydration shell accurately for most of the metal cations, an accuracy equivalent to that of the complex 12-6-4 LJ-type and double exponential potential models. A systematic comparison with the existing M2+ models is presented as well in terms of effective ion radii, diffusion constants, water exchange rates, and ion-water interactions. Molecular dynamics simulations of metal substitution in Escherichia coli glyoxalase I variants show the great potential of our new models for metalloproteins.

Bridging the 12-6-4 Model and the Fluctuating Charge Model.

Metal ions play important roles in various biological systems. Molecular dynamics (MD) using classical force field has become a popular research tool to study biological systems at the atomic level. However, meaningful MD simulations require reliable models and parameters. Previously we showed that the 12-6 Lennard-Jones nonbonded model for ions could not reproduce the experimental hydration free energy (HFE) and ion-oxygen distance (IOD) values simultaneously when ion has a charge of +2 or higher. We discussed that this deficiency arises from the overlook of the ion-induced dipole interaction in the 12-6 model, and this term is proportional to 1/r 4 based on theory. Hence, we developed the 12-6-4 model and showed it could solve this deficiency in a physically meaningful way. However, our previous research also found that the 12-6-4 model overestimated the coordination numbers (CNs) for some highly charged metal ions. And we attributed this artifact to that the current 12-6-4 scheme lacks a correction for the interactions among the first solvation shell water molecules. In the present study, we considered the ion-included dipole interaction by using the 12-6 model with adjusting the atomic charges of the first solvation shell water molecules. This strategy not only considers the ion-induced dipole interaction between ion and the first solvation shell water molecules but also well accounts for the increased repulsion among these water molecules compared to the bulk water molecules. We showed this strategy could well reproduce the experimental HFE and IOD values for Mg2+, Zn2+, Al3+, Fe3+, and In3+ and solve the CN overestimation issue of the 12-6-4 model for Fe3+ and In3+. Moreover, our simulation results showed good agreement with previous ab initio MD simulations. In addition, we derived the physical relationship between the C 4 parameter and induced dipole moment, which agreed well with our simulation results. Finally, we discussed the implications of the present work for simulating metalloproteins. Due to the fluctuating charge model uses a similar concept to the 12-6 model with adjusting atomic charges, we believe the present study builds a bridge between the 12-6-4 model and the fluctuating charge model.

Automation of some macromolecular properties using a machine learning approach

In this study, we employed a newly developed method to predict macromolecular properties using a swarm artificial neural network (ANN) method as a machine learning approach. In this method, the molecular structures are represented by the feature description vectors used as training input data for a neural network. This study aims to develop an efficient approach for training an ANN using either experimental or quantum mechanics data. We aim to introduce an error model controlling the reliability of the prediction confidence interval using a bootstrapping swarm approach. We created different datasets of selected experimental or quantum mechanics results. Using this optimized ANN, we hope to predict properties and their statistical errors for new molecules. There are four datasets used in this study. That includes the dataset of 642 small organic molecules with known experimental hydration free energies, the dataset of 1475 experimental pKa values of ionizable groups in 192 proteins, the dataset of 2693 mutants in 14 proteins with given experimental values of changes in the Gibbs free energy, and a dataset of 7101 quantum mechanics heat of formation calculations. All the data are prepared and optimized using the AMBER force field in the CHARMM macromolecular computer simulation program. The bootstrapping swarm ANN code for performing the optimization and prediction is written in Python computer programming language. The descriptor vectors of the small molecules are based on the Coulomb matrix and sum over bond properties. For the macromolecular systems, they consider the chemical-physical fingerprints of the region in the vicinity of each amino acid.

Advanced Electrostatic Model for Monovalent Ions Based on Ab Initio Energy Decomposition.

Ions play important roles in the structures and functions of biomolecules. In biomolecular simulations, ions either directly interact with biomolecules or provide an ionic environment that influences electrostatic interactions of solutes. The AMOEBA+ water model has demonstrated significant advancement of the classical force field for describing molecular interactions due to its improvements on the functional forms to account for essential physics. This work expands the applicability of the AMOEBA+ model toward alkali metal (Li, Na, K, Rb, and Cs) and halogen (F, Cl, Br, and I) ions. Various quantum chemical data on ion-ion and ion-water interactions, experimental ion hydration free energies, and lattice energies of salt crystals are used in the parametrization. The final parameters are verified with other properties outside of the parametrization data, including lattice energies of additional salt crystals and ionic activity coefficients in solution. The new model captures a wide range of ion properties from the gas phase to solution phase and crystals. More importantly, AMOEBA+ provides energy components that are consistent with ab initio energy decomposition. Thus, we expect AMOEBA+ to be more general, transferable, and valuable for the interpretation of intermolecular forces in efficient classical simulations.

Machine learning of free energies in chemical compound space using ensemble representations: Reaching experimental uncertainty for solvation.

Free energies govern the behavior of soft and liquid matter, and improving their predictions could have a large impact on the development of drugs, electrolytes, or homogeneous catalysts. Unfortunately, it is challenging to devise an accurate description of effects governing solvation such as hydrogen-bonding, van der Waals interactions, or conformational sampling. We present a Free energy Machine Learning (FML) model applicable throughout chemical compound space and based on a representation that employs Boltzmann averages to account for an approximated sampling of configurational space. Using the FreeSolv database, FML's out-of-sample prediction errors of experimental hydration free energies decay systematically with training set size, and experimental uncertainty (0.6 kcal/mol) is reached after training on 490 molecules (80% of FreeSolv). Corresponding FML model errors are on par with state-of-the art physics based approaches. To generate the input representation for a new query compound, FML requires approximate and short molecular dynamics runs. We showcase its usefulness through analysis of solvation free energies for 116k organic molecules (all force-field compatible molecules in the QM9 database), identifying the most and least solvated systems and rediscovering quasi-linear structure-property relationships in terms of simple descriptors such as hydrogen-bond donors, number of NH or OH groups, number of oxygen atoms in hydrocarbons, and number of heavy atoms. FML's accuracy is maximal when the temperature used for the molecular dynamics simulation to generate averaged input representation samples in training is the same as for the query compounds. The sampling time for the representation converges rapidly with respect to the prediction error.

Parameterization of Monovalent Ions for the OPC3, OPC, TIP3P-FB, and TIP4P-FB Water Models.

Monovalent ions play significant roles in various biological and material systems. Recently, four new water models (OPC3, OPC, TIP3P-FB, and TIP4P-FB), with significantly improved descriptions of condensed phase water, have been developed. The pairwise interaction between the metal ion and water necessitates the development of ion parameters specifically for these water models. Herein, we parameterized the 12-6 and the 12-6-4 nonbonded models for 12 monovalent ions with the respective four new water models. These monovalent ions contain eight cations including alkali metal ions (Li+, Na+, K+, Rb+, Cs+), transition-metal ions (Cu+ and Ag+), and Tl+ from the boron family, along with four halide anions (F-, Cl-, Br-, I-). Our parameters were designed to reproduce the target hydration free energies (the 12-6 hydration free energy (HFE) set), the ion-oxygen distances (the 12-6 ion-oxygen distance (IOD) set), or both of them (the 12-6-4 set). The 12-6-4 parameter set provides highly accurate structural features overcoming the limitations of the routinely used 12-6 nonbonded model for ions. Specifically, we note that the 12-6-4 parameter set is able to reproduce experimental hydration free energies within 1 kcal/mol and experimental ion-oxygen distances within 0.01 Å simultaneously. We further reproduced the experimentally determined activity derivatives for salt solutions, validating the ion parameters for simulations of ion pairs. The improved performance of the present water models over our previous parameter sets for the TIP3P, TIP4P, and SPC/E water models (Li, P. et al J. Chem. Theory Comput. 2015 11 1645 1657) highlights the importance of the choice of water model in conjunction with the metal ion parameter set.

Predicting Hydration Free Energies of the FreeSolv Database of Drug-like Molecules with Molecular Density Functional Theory.

We assess the performance of molecular density functional theory (MDFT) to predict hydration free energies of the small drug-like molecules benchmark, FreeSolv. The MDFT in the hypernetted chain approximation (HNC) coupled with a pressure correction predicts experimental hydration free energies of the FreeSolv database within 1 kcal/mol with an average computation time of 2 cpu·min per molecule. This is the same accuracy as for simulation-based free energy calculations that typically require hundreds of cpu·h or tens of gpu·h per molecule.

Systematic Parametrization of Divalent Metal Ions for the OPC3, OPC, TIP3P-FB, and TIP4P-FB Water Models.

Divalent metal ions play important roles in biological and materials systems. Molecular dynamics simulation is an efficient tool to investigate these systems at the microscopic level. Recently, four new water models (OPC3, OPC, TIP3P-FB, and TIP4P-FB) have been developed and better represent the physical properties of water than previous models. Metal ion parameters are dependent on the water model employed, making it necessary to develop metal ion parameters for select new water models. In the present work, we performed parameter scanning for the 12-6 Lennard-Jones nonbonded model of divalent metal ions in conjunction with the four new water models as well as four previous water models (TIP3P, SPC/E, TIP4P, and TIP4P-Ew). We found that these new three-point and four-point water models provide comparable or significantly improved performance for the simulation of divalent metal ions when compared to previous water models in the same category. Among all eight water models, the OPC3 water model yields the best performance for the simulation of divalent metal ions in the aqueous phase when using the 12-6 model. On the basis of the scanning results, we independently parametrized the 12-6 model for 24 divalent metal ions with each of the four new water models. As noted previously, the 12-6 model still fails to simultaneously reproduce the experimental hydration free energy (HFE) and ion-oxygen distance (IOD) values even with these new water models. To solve this problem, we parametrized the 12-6-4 model for the 16 divalent metal ions for which we have both experimental HFE and IOD values for each of the four new water models. The final parameters are able to reproduce both the experimental HFE and IOD values accurately. To validate the transferability of our parameters, we carried out benchmark calculations to predict the energies and geometries of ion-water clusters as well as the ion diffusivity coefficient of Mg2+. By comparison to quantum chemical calculations and experimental data, these results show that our parameters are well designed and have excellent transferability. The metal ion parameters for the 12-6 and 12-6-4 models reported herein can be employed in simulations of various biological and materials systems when using the OPC3, OPC, TIP3P-FB, or TIP4P-FB water model.

Balanced polarizable Drude force field parameters for molecular anions: phosphates, sulfates, sulfamates, and oxides.

Polarizable force fields are emerging as a more accurate alternative to additive force fields in terms of modeling and simulations of a variety of chemicals including biomolecules. Explicit treatment of induced polarization in charged species such as phosphates and sulfates offers the potential for achieving an improved atomistic understanding of the physical forces driving their interactions with their environments. To help achieve this, in this study we present balanced Drude polarizable force field parameters for molecular ions including phosphates, sulfates, sulfamates, and oxides. Better balance was primarily achieved in the relative values of minimum interaction energies and distances of the anionic model compounds with water at the Drude and quantum mechanical (QM) model chemistries. Parametrization involved reoptimizing available parameters as well as extending the force field to new molecules with the goal of achieving self-consistency with respect to the Lennard-Jones and electrostatic parameters targeting QM and experimental hydration free energies. The resulting force field parameters achieve consistent treatment across the studied anions, facilitating more balanced simulations of biomolecules and small organic molecules in the context of the classical Drude polarizable force field. Graphical abstract.

Statistical efficiency of methods for computing free energy of hydration.

The hydration free energy (HFE) is a critical property for predicting and understanding chemical and biological processes in aqueous solution. There are a number of computational methods to derive HFE, generally classified into the equilibrium or non-equilibrium methods, based on the type of calculations used. In the present study, we compute the hydration free energies of 34 small, neutral, organic molecules with experimental HFE between +2 and -16 kcal/mol. The one-sided non-equilibrium methods Jarzynski Forward (JF) and Backward (JB), the two-sided non-equilibrium methods Jarzynski mean based on the average of JF and JB, Crooks Gaussian Intersection (CGI), and the Bennett Acceptance Ratio (BAR) are compared to the estimates from the two-sided equilibrium method Multistate Bennett Acceptance Ratio (MBAR), which is considered as the reference method for HFE calculations, and experimental data from the literature. Our results show that the estimated hydration free energies from all the methods are consistent with MBAR results, and all methods provide a mean absolute error of ∼0.8 kcal/mol and root mean square error of ∼1 kcal for the 34 organic molecules studied. In addition, the results show that one-sided methods JF and JB result in systematic deviations that cannot be corrected entirely. The statistical efficiency ε of the different methods can be expressed as the one over the simulation time times the average variance in the HFE. From such an analysis, we conclude that ε(MBAR) > ε(BAR) ≈ ε(CGI) > ε(JX), where JX is any of the Jarzynski methods. In other words, the non-equilibrium methods tested here for the prediction of HFE have lower computational efficiency than the MBAR method.

Hydration Free Energies in the FreeSolv Database Calculated with Polarized Iterative Hirshfeld Charges.

Computer simulations of biomolecular systems often use force fields, which are combinations of simple empirical atom-based functions to describe the molecular interactions. Even though polarizable force fields give a more detailed description of intermolecular interactions, nonpolarizable force fields, developed several decades ago, are often still preferred because of their reduced computation cost. Electrostatic interactions play a major role in biomolecular systems and are therein described by atomic point charges. In this work, we address the performance of different atomic charges to reproduce experimental hydration free energies in the FreeSolv database in combination with the GAFF force field. Atomic charges were calculated by two atoms-in-molecules approaches, Hirshfeld-I and Minimal Basis Iterative Stockholder (MBIS). To account for polarization effects, the charges were derived from the solute's electron density computed with an implicit solvent model, and the energy required to polarize the solute was added to the free energy cycle. The calculated hydration free energies were analyzed with an error model, revealing systematic errors associated with specific functional groups or chemical elements. The best agreement with the experimental data is observed for the AM1-BCC and the MBIS atomic charge methods. The latter includes the solvent polarization and presents a root-mean-square error of 2.0 kcal mol-1 for the 613 organic molecules studied. The largest deviation was observed for phosphorus-containing molecules and the molecules with amide, ester and amine functional groups.

Classical Drude Polarizable Force Field Model for Methyl Phosphate and Its Interactions with Mg2.

Phosphate groups are essential components of nucleic acids and proteins, whose interactions with solvent, metal ions, and ionic side chains help control folding and binding. Methyl phosphate (MP) represents a simple analog of phosphate moieties that are post-translation modifications in proteins and present at the termini of nucleic acids, among other environments. In the present study, we optimized parameters for use in polarizable molecular dynamics simulations of MP in its mono- and dianionic forms, MP- ≡ CH3HPO4- and MP2- ≡ CH3PO42-, along with P i2- ≡ HPO42-, in the context of the classical Drude oscillator model. Parameter optimization was done in a manner consistent with the remainder of the Drude molecular mechanics force field, choosing atomic charges and polarizabilities to reproduce molecular properties from quantum mechanics as well as experimental hydration free energies. Optimized parameters were similar to existing dimethyl phosphate parameters, with a few significant differences. The developed parameters were then used to compute magnesium binding affinities in aqueous solution, using alchemical molecular dynamics free energy simulations. Good agreement with experiment was obtained, and outer sphere binding was shown to be predominant for MP- and MP2-.