A re-parameterization of the standard TIP4P water model for use with Ewald techniques is introduced, providing an overall global improvement in water properties relative to several popular nonpolarizable and polarizable water potentials. Using high precision simulations, and careful application of standard analytical corrections, we show that the new TIP4P-Ew potential has a density maximum at approximately 1 degrees C, and reproduces experimental bulk-densities and the enthalpy of vaporization, DeltaH(vap), from -37.5 to 127 degrees C at 1 atm with an absolute average error of less than 1%. Structural properties are in very good agreement with x-ray scattering intensities at temperatures between 0 and 77 degrees C and dynamical properties such as self-diffusion coefficient are in excellent agreement with experiment. The parameterization approach used can be easily generalized to rehabilitate any water force field using available experimental data over a range of thermodynamic points.

Quantitative free energy computation involves both using a model that is sufficiently faithful to the experimental system under study ͑accuracy͒ and establishing statistically meaningful measures of the uncertainties resulting from finite sampling ͑precision͒. We use large-scale distributed computing to access sufficient computational resources to extensively sample molecular systems and thus reduce statistical uncertainty of measured free energies. In order to examine the accuracy of a range of common models used for protein simulation, we calculate the free energy of hydration of 15 amino acid side chain analogs derived from recent versions of the OPLS-AA, CHARMM, and AMBER parameter sets in TIP3P water using thermodynamic integration. We achieve a high degree of statistical precision in our simulations, obtaining uncertainties for the free energy of hydration of 0.02-0.05 kcal/mol, which are in general an order of magnitude smaller than those found in other studies. Notably, this level of precision is comparable to that obtained in experimental hydration free energy measurements of the same molecules. Root mean square differences from experiment over the set of molecules examined using AMBER-, CHARMM-, and OPLS-AA-derived parameters were 1.35 kcal/mol, 1.31 kcal/mol, and 0.85 kcal/mol, respectively. Under the simulation conditions used, these force fields tend to uniformly underestimate solubility of all the side chain analogs. The relative free energies of hydration between amino acid side chain analogs were closer to experiment but still exhibited significant deviations. Although extensive computational resources may be needed for large numbers of molecules, sufficient computational resources to calculate precise free energy calculations for small molecules are accessible to most researchers.

A rigorous formalism for the extraction of state-to-state transition functions from a Boltzmann-weighted ensemble of microcanonical molecular dynamics simulations has been developed as a way to study the kinetics of protein folding in the context of a Markov chain. Analysis of these transition functions for signatures of Markovian behavior is described. The method has been applied to an example problem that is based on an underlying Markov process. The example problem shows that when an instance of the process is analyzed under the assumption that the underlying states have been aggregated into macrostates, a procedure known as lumping, the resulting chain appears to have been produced by a non-Markovian process when viewed at high temporal resolution. However, when viewed on longer time scales, and for appropriately lumped macrostates, Markovian behavior can be recovered. The potential for extracting the long time scale behavior of the folding process from a large number of short, independent molecular dynamics simulations is also explored.

Abstract:The growing adoption of generalized-ensemble algorithms for biomolecular simulation has resulted in a resurgence in the use of the weighted histogram analysis method (WHAM) to make use of all data generated by these simulations. Unfortunately, the original presentation of WHAM by Kumar et al. is not directly applicable to data generated by these methods. WHAM was originally formulated to combine data from independent samplings of the canonical ensemble, whereas many generalized-ensemble algorithms sample from mixtures of canonical ensembles at different temperatures. Sorting configurations generated from a parallel tempering simulation by temperature obscures the temporal correlation in the data and results in an improper treatment of the statistical uncertainties used in constructing the estimate of the density of states. Here we present variants of WHAM, STWHAM and PTWHAM, derived with the same set of assumptions, that can be directly applied to several generalized ensemble algorithms, including simulated tempering, parallel tempering (better known as replica-exchange among temperatures), and replica-exchange simulated tempering. We present methods that explicitly capture the considerable temporal correlation in sequentially generated configurations using autocorrelation analysis. This allows estimation of the statistical uncertainty in WHAM estimates of expectations for the canonical ensemble. We test the method with a one-dimensional model system and then apply it to the estimation of potentials of mean force from parallel tempering simulations of the alanine dipeptide in both implicit and explicit solvent.

Replica-exchange molecular dynamics simulations in implicit solvent have been carried out to study the folding thermodynamics of a designed 20-residue peptide, or ''miniprotein.'' The simulations in this study used the AMBER (parm94) force field along with the generalized Born͞solvent-accessible surface area implicit solvent model, and they spanned a range of temperatures from 273 to 630 K. Starting from a completely extended initial conformation, simulations of one peptide sequence sample conformations that are <1.0 Å C ␣ rms positional deviation from structures in the corresponding NMR ensemble. These folded states are thermodynamically stable with a simulated melting temperature of Ϸ400 K, and they satisfy the majority of experimentally observed NMR restraints. Simulations of a related mutant peptide show a degenerate ensemble of states at low temperature, in agreement with experimental results.A recent paper by Neidigh et al. (1) describes the design of a stably folded, 20-residue ''miniprotein.'' This protein was derived from C-terminal fragments of the 39-residue exendin-4 peptide (2). Several constructs of increasing stability were made, gradually introducing stabilizing features like helical N-capping residues and a solvent-exposed salt bridge. The final 20-residue peptide showed a cooperative melting transition with a midpoint of Ϸ315 K in aqueous solution at pH 7. The NMR structure [Protein Data Bank ID 1L2Y (3)] of this peptide, denoted TC5b by Neidigh et al. (1), was determined based on 169 interproton distance restraints, and it shows a well structured hydrophobic core where the indole side chain of a Trp residue is buried between the rings of two Pro residues. The small size and stability of this protein make it an ideal choice for simulation studies of protein folding.There are three main types of questions that arise in computational simulations of the protein folding process: questions of structure, thermodynamics, and kinetics. The first and most basic question is whether the simulation can reproduce the folded conformation of the protein. Once it is clear that the simulations can sample between the folded and unfolded states of the protein, thermodynamic questions arise. Specifically, the simulation's ability to reproduce the stability of the folded state as a function of temperature or denaturant should be assessed. If the simulation reproduces the melting curve of the protein, it becomes reasonable to look at further thermodynamic details of the folding process, including its entropic and enthalpic contributions as well as the nature of the heat-capacity changes on protein folding or unfolding. The final area of interest is the question of kinetics: whether the simulation predicts the correct overall rate of protein folding, and the nature of that process. Kinetics is an area of interest where the microscopic insights from simulation are a valuable complement to experiment, helping to establish whether folding is a simple process with several clear intermediates and barriers or one with a b...

Historically, experimental measurements have been used to bias biomolecular simulations toward structures compatible with those observations via the addition of ad hoc restraint terms. We show how the maximum entropy formalism provides a principled approach to enforce concordance with these measurements in a minimally biased way, yielding restraints that are linear functions of the target observables and specifying a straightforward scheme to determine the biasing weights. These restraints are compared with instantaneous and ensemble-averaged harmonic restraint schemes, illustrating their similarities and limitations.

Abstract. Protein folding involves physical timescales-microseconds to seconds-that are too long to be studied directly by straightforward molecular dynamics simulation, where the fundamental timestep is constrained to femtoseconds. Here we show how the long-time statistical dynamics of a simple solvated biomolecular system can be well described by a discrete-state Markov chain model constructed from trajectories that are an order of magnitude shorter than the longest relaxation times of the system. This suggests that such models, appropriately constructed from short molecular dynamics simulations, may have utility in the study of long-time conformational dynamics.

No abstract