A biomolecular force field based on the free enthalpy of hydration and solvation: the GROMOS force-field parameter sets 53A5 and 53A6

Successive parameterizations of the GROMOS force field have been used successfully to simulate biomolecular systems over a long period of time. The continuing expansion of computational power with time makes it possible to compute ever more properties for an increasing variety of molecular systems with greater precision. This has led to recurrent parameterizations of the GROMOS force field all aimed at achieving better agreement with experimental data. Here we report the results of the latest, extensive reparameterization of the GROMOS force field. In contrast to the parameterization of other biomolecular force fields, this parameterization of the GROMOS force field is based primarily on reproducing the free enthalpies of hydration and apolar solvation for a range of compounds. This approach was chosen because the relative free enthalpy of solvation between polar and apolar environments is a key property in many biomolecular processes of interest, such as protein folding, biomolecular association, membrane formation, and transport over membranes. The newest parameter sets, 53A5 and 53A6, were optimized by first fitting to reproduce the thermodynamic properties of pure liquids of a range of small polar molecules and the solvation free enthalpies of amino acid analogs in cyclohexane (53A5). The partial charges were then adjusted to reproduce the hydration free enthalpies in water (53A6). Both parameter sets are fully documented, and the differences between these and previous parameter sets are discussed.     

Assessing the accuracy of improved force‐matched water models derived from Ab initio molecular dynamics simulations

The accuracy of water models derived from ab initio molecular dynamics simulations by means on an improved force‐matching scheme is assessed for various thermodynamic, transport, and structural properties. It is found that although the resulting force‐matched water models are typically less accurate than fully empirical force fields in predicting thermodynamic properties, they are nevertheless much more accurate than generally appreciated in reproducing the structure of liquid water and in fact superseding most of the commonly used empirical water models. This development demonstrates the feasibility to routinely parametrize computationally efficient yet predictive potential energy functions based on accurate ab initio molecular dynamics simulations for a large variety of different systems. © 2016 Wiley Periodicals, Inc.     

Density functional theory for efficient ab initio molecular dynamics simulations in solution

We present a density functional for first-principles molecular dynamics simulations that includes the electrostatic effects of a continuous dielectric medium. It allows for numerical simulations of molecules in solution in a model polar solvent. We propose a smooth dielectric model function to model solvation into water and demonstrate its good numerical properties for total energy calculations and constant energy molecular dynamics.     

An ab initio force field for the cofactors of bacterial photosynthesis

This article presents a new ab initio force field for the cofactors of bacterial photosynthesis, namely quinones and bacteriochlorophylls. The parameters has been designed to be suitable for molecular dynamics simulations of photosynthetic proteins by being compatible with the AMBER force field. To our knowledge, this is the first force field for photosynthetic cofactors based on a reliable set of ab initio density functional reference data for methyl bacteriochlorophyll a, methyl bacteriopheophytin a, and of a derivative of ubiquinone. Indeed, the new molecular mechanics force field is able to reproduce very well not only the experimental and ab initio structural properties and the vibrational spectra of the molecules, but also the eigenvectors of the molecular normal modes. For this reason it might also be helpful to understand vibrational spectroscopy results obtained on reaction center proteins.     

A comparative theoretical study of dipeptide solvation in water

Molecular dynamics studies have been performed on the zwitterionic form of the dipeptide glycine-alanine in water, with focus on the solvation and electrostatic properties using a range of theoretical methods, from purely classical force fields, through mixed quantum mechanical/molecular mechanical simulations, to fully quantum mechanical Car-Parrinello calculations. The results of these studies show that the solvation pattern is similar for all methods used for most atoms in the dipeptide, but can differ substantially for some groups; namely the carboxy and aminoterminii, and the backbone amid NH group. This might have implications in other theoretical studies of peptides and proteins with charged -NH(3) (+) and -CO(2) (-) side chains solvated in water. Hybrid quantum mechanical/molecular mechanical simulations successfully reproduce the solvation patterns from the fully quantum mechanical simulations (PACS numbers: 87.14.Ee, 87.15.Aa, 87.15.He, 71.15.Pd).     

Polarizable model potential function for nucleic acid bases

A polarizable model potential (PMP) function for adenine (A), cytosine (C), guanine (G), thymine (T), and uracil (U) is developed on the basis of ab initio molecular orbital calculations at the MP2/6-31+G* level. The PMP function consists of Coulomb, van der Waals, and polarization terms. The permanent atomic charges of the Coulomb term are determined by using electrostatic potential (ESP) optimization. The multicenter polarizabilities of the polarization term are determined by using polarized one-electron potential (POP) optimization in which the electron density changes induced by a test charge are target. Isotropic and anisotropic polarizabilities are adopted as the multicenter polarizabilities. In the PMP calculations using the optimized parameters, the interaction energies of Watson-Crick type A-T and C-G base pairs were -15.6 and -29.4 kcal/mol, respectively. The interaction energy of Hoogsteen type A-T base pair was -17.8 kcal/mol. These results reproduce well the quantum chemistry calculations at the MP2/6-311++G(3df,2pd) level within the differences of 0.6 kcal/mol. The stacking energies of A-T and C-G were -9.7 and -10.9 kcal/mol. These reproduce well the calculation results at the MP2/6-311++G (2d,2p) level within the differences of 1.3 kcal/mol. The potential energy surfaces of the system in which a sodium ion or a chloride ion is adjacent to the nucleic acid base are calculated. The interaction energies of the PMP function reproduced well the calculation results at the MP2/6-31+G* or MP2/6-311++G(2d,2p) level. The reason why the PMP function reproduces well the high-level quantum mechanical interaction energies is addressed from the viewpoint of each energy terms.     

Surveying implicit solvent models for estimating small molecule absolute hydration free energies

Implicit solvent models are powerful tools in accounting for the aqueous environment at a fraction of the computational expense of explicit solvent representations. Here, we compare the ability of common implicit solvent models (TC, OBC, OBC2, GBMV, GBMV2, GBSW, GBSW/MS, GBSW/MS2 and FACTS) to reproduce experimental absolute hydration free energies for a series of 499 small neutral molecules that are modeled using AMBER/GAFF parameters and AM1-BCC charges. Given optimized surface tension coefficients for scaling the surface area term in the nonpolar contribution, most implicit solvent models demonstrate reasonable agreement with extensive explicit solvent simulations (average difference 1.0-1.7 kcal/mol and R(2)=0.81-0.91) and with experimental hydration free energies (average unsigned errors=1.1-1.4 kcal/mol and R(2)=0.66-0.81). Chemical classes of compounds are identified that need further optimization of their ligand force field parameters and others that require improvement in the physical parameters of the implicit solvent models themselves. More sophisticated nonpolar models are also likely necessary to more effectively represent the underlying physics of solvation and take the quality of hydration free energies estimated from implicit solvent models to the next level.     

OPLS all-atom force field for carbohydrates

The OPLS all-atom (AA) force field for organic and biomolecular systems has been expanded to include carbohydrates. Starting with reported nonbonded parameters of alcohols, ethers, and diols, torsional parameters were fit to reproduce results from ab initio calculations on the hexopyranoses, α,β-d -glucopyranose, α,β-d -mannopyranose, α,β-d -galactopyranose, methyl α,β-d -glucopyranoside, and methyl α,β-d -mannopyranoside. In all, geometry optimizations were carried out for 144 conformers at the restricted Hartree–Fock (RHF)/6–31G* level. For the conformers with a relative energy within 3 kcal/mol of the global minima, the effects of electron correlation and basis-set extension were considered by performing single-point calculations with density functional theory at the B3LYP/6–311+G** level. The torsional parameters for the OPLS-AA force field were parameterized to reproduce the energies and structures of these 44 conformers. The resultant force field reproduces the ab initio calculated energies with an average unsigned error of 0.41 kcal/mol. The α/β ratios as well as the relative energies between the isomeric hexopyranoses are in good accord with the ab initio results. The predictive abilities of the force field were also tested against RHF/6–31G* results for d -allopyranose with excellent success; a surprising discovery is that the lowest energy conformer of d -allopyranose is a β anomer. © 1997 John Wiley & Sons, Inc. J Comput Chem 18 : 1955–1970, 1997     

An efficient hybrid explicit/implicit solvent method for biomolecular simulations

We present a new hybrid explicit/implicit solvent method for dynamics simulations of macromolecular systems. The method models explicitly the hydration of the solute by either a layer or sphere of water molecules, and the generalized Born (GB) theory is used to treat the bulk continuum solvent outside the explicit simulation volume. To reduce the computational cost, we implemented a multigrid method for evaluating the pairwise electrostatic and GB terms. It is shown that for typical ion and protein simulations our method achieves similar equilibrium and dynamical observables as the conventional particle mesh Ewald (PME) method. Simulation timings are reported, which indicate that the hybrid method is much faster than PME, primarily due to a significant reduction in the number of explicit water molecules required to model hydration effects.     

On the Reliability of the AMBER Force Field and its Empirical Dispersion Contribution for the Description of Noncovalent Complexes

The reliability of the AMBER force field is tested by comparing the total interaction energy and dispersion energy with the reference data obtained at the density functional theory–symmetry‐adapted perturbation treatment (DFT–SAPT)/aug‐cc‐pVDZ level. The comparison is made for 194 different geometries of noncovalent complexes (H‐bonded, stacked, mixed, and dispersion‐bound), at the equilibrium distances as well as at longer distances (up to a relative distance of two). The total interaction energies agree very well with the reference data and only the strength of H‐bonded complexes is slightly underestimated. In the case of dispersion energy, the overall agreement is even better, with the exception of the stacked aromatic systems, where the empirical dispersion energy is overestimated. The use of AMBER interaction energy and AMBER dispersion energy for different types of noncovalent complexes at equilibrium as well as at longer distances is thus justified, except for a few cases, such as the water molecule, where the dispersion energy is highly inaccurate.     

Assessment of biomolecular force fields for molecular dynamics simulations in a protein crystal

Different biomolecular force fields (OPLS‐AA, AMBER03, and GROMOS96) in conjunction with SPC, SPC/E and TIP3P water models are assessed for molecular dynamics simulations in a tetragonal lysozyme crystal. The root mean square deviations for the C_a atoms of lysozymes are about 0.1 to 0.2 nm from OPLS‐AA and AMBER03, smaller than 0.4 nm from GROMOS96. All force fields exhibit similar pattern in B‐factors, whereas OPLS‐AA and AMBER03 accurately reproduce experimental measurements. Despite slight variations, the primary secondary structures are well conserved using different force fields. Water diffusion in the crystal is approximately ten‐fold slower than in bulk phase. The directional and average water diffusivities from OPLS‐AA and AMBER03 along with SPC/E model match fairly well with experimental data. Compared to GROMOS96, OPLS‐AA and AMBER03 predict larger hydrophilic solvent‐accessible surface area of lysozyme, more hydrogen bonds between lysozyme and water, and higher percentage of water in hydration shell. SPC, SPC/E and TIP3P water models have similar performance in most energetic and structural properties, but SPC/E outperforms in water diffusion. While all force fields overestimate the mobility and electrical conductivity of NaCl, a combination of OPLS‐AA for lysozyme and the Kirkwood‐Buff model for ions is superior to others. As attributed to the steric restraints and surface interactions, the mobility and conductivity in the crystal are reduced by one to two orders of magnitude from aqueous solution. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Comparative study on formamide–water complex

The hydrogen bonding of 1:1 complexes formed between formamide and water molecule have been investigated systematically using Hartree–Fock (HF), hybrid density functional theory (B3LYP), and post‐Hartree–Fock (MP2 and CCSD(T)) methods with range of basis sets 6‐31G(d), cc‐pVXZ (X = D, T, Q) and aug‐cc‐pVYZ (Y = D, T). Three stable structures are considered on the potential energy surface of formamide and water system. The optimized geometric parameters and interaction energies for various isomers at different levels are estimated. The IR frequencies, intensities, and frequency shifts are reported. This study shows that B3LYP/aug‐cc‐pVDZ method gives better performance for formamide‐water complexes. © 2009 Wiley Periodicals, Inc. Int J Quantum Chem, 2010.     

Extending the treatment of backbone energetics in protein force fields: limitations of gas-phase quantum mechanics in reproducing protein conformational distributions in molecular dynamics simulations

Computational studies of proteins based on empirical force fields represent a powerful tool to obtain structure-function relationships at an atomic level, and are central in current efforts to solve the protein folding problem. The results from studies applying these tools are, however, dependent on the quality of the force fields used. In particular, accurate treatment of the peptide backbone is crucial to achieve representative conformational distributions in simulation studies. To improve the treatment of the peptide backbone, quantum mechanical (QM) and molecular mechanical (MM) calculations were undertaken on the alanine, glycine, and proline dipeptides, and the results from these calculations were combined with molecular dynamics (MD) simulations of proteins in crystal and aqueous environments. QM potential energy maps of the alanine and glycine dipeptides at the LMP2/cc-pVxZ//MP2/6-31G* levels, where x = D, T, and Q, were determined, and are compared to available QM studies on these molecules. The LMP2/cc-pVQZ//MP2/6-31G* energy surfaces for all three dipeptides were then used to improve the MM treatment of the dipeptides. These improvements included additional parameter optimization via Monte Carlo simulated annealing and extension of the potential energy function to contain peptide backbone phi, psi dihedral crossterms or a phi, psi grid-based energy correction term. Simultaneously, MD simulations of up to seven proteins in their crystalline environments were used to validate the force field enhancements. Comparison with QM and crystallographic data showed that an additional optimization of the phi, psi dihedral parameters along with the grid-based energy correction were required to yield significant improvements over the CHARMM22 force field. However, systematic deviations in the treatment of phi and psi in the helical and sheet regions were evident. Accordingly, empirical adjustments were made to the grid-based energy correction for alanine and glycine to account for these systematic differences. These adjustments lead to greater deviations from QM data for the two dipeptides but also yielded improved agreement with experimental crystallographic data. These improvements enhance the quality of the CHARMM force field in treating proteins. This extension of the potential energy function is anticipated to facilitate improved treatment of biological macromolecules via MM approaches in general.     

Liquid properties of dimethyl ether from molecular dynamics simulations using Ab Initio force fields

We have used molecular dynamic simulations to study the structural and dynamical properties of liquid dimethyl ether (DME) with a newly constructed ab initio force field in this article. The ab initio potential energy data were calculated at the second order Møller‐Plesset (MP2) perturbation theory with Dunning's correlation consistent basis sets (up to aug‐cc‐pVQZ). We considered 17 configurations of the DME dime for the orientation sampling. The calculated MP2 potential data were used to construct a 3‐site united atom force field model. The simulation results are compared with those using the empirical force field of Jorgensen and Ibrahim (Jorgensen and Ibrahim, J Am Chem Soc 1981, 103, 3976) and with available experimental measurements. We obtain quantitative agreements for the atom‐wise radial distribution functions, the self‐diffusion coefficients, and the shear viscosities over a wide range of experimental conditions. This force field thus provides a suitable starting point to predict liquid properties of DME from first principles intermolecular interactions with no empirical data input a priori. © 2012 Wiley Periodicals, Inc.     

Evaluation of solvation free energies for small molecules with the AMOEBA polarizable force field

The effects of electronic polarization in biomolecular interactions will differ depending on the local dielectric constant of the environment, such as in solvent, DNA, proteins, and membranes. Here the performance of the AMOEBA polarizable force field is evaluated under nonaqueous conditions by calculating the solvation free energies of small molecules in four common organic solvents. Results are compared with experimental data and equivalent simulations performed with the GAFF pairwise‐additive force field. Although AMOEBA results give mean errors close to “chemical accuracy,” GAFF performs surprisingly well, with statistically significantly more accurate results than AMOEBA in some solvents. However, for both models, free energies calculated in chloroform show worst agreement to experiment and individual solutes are consistently poor performers, suggesting non‐potential‐specific errors also contribute to inaccuracy. Scope for the improvement of both potentials remains limited by the lack of high quality experimental data across multiple solvents, particularly those of high dielectric constant. © 2016 The Authors. Journal of Computational Chemistry Published by Wiley Periodicals, Inc.