Thomas E. Cheatham

The Amber biomolecular simulation programs.

We describe the development, current features, and some directions for future development of the Amber package of computer programs. This package evolved from a program that was constructed in the late 1970s to do Assisted Model Building with Energy Refinement, and now contains a group of programs embodying a number of powerful tools of modern computational chemistry, focused on molecular dynamics and free energy calculations of proteins, nucleic acids, and carbohydrates.

AMBER, a package of computer programs for applying molecular mechanics, normal mode analysis, molecular dynamics and free energy calculations to simulate the structural and energetic properties of molecules

We describe the development, current features, and some directions for future development of the AMBER package of computer programs. This package has evolved from a program that was constructed to do Assisted Model Building and Energy Refinement to a group of programs embodying a number of the powerful tools of modern computational chemistry-molecular dynamics and free energy calculations.

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 Åqvist 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, TIP4PEW, 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.

PTRAJ and CPPTRAJ: Software for Processing and Analysis of Molecular Dynamics Trajectory Data

We describe PTRAJ and its successor CPPTRAJ, two complementary, portable, and freely available computer programs for the analysis and processing of time series of three-dimensional atomic positions (i.e., coordinate trajectories) and the data therein derived. Common tools include the ability to manipulate the data to convert among trajectory formats, process groups of trajectories generated with ensemble methods (e.g., replica exchange molecular dynamics), image with periodic boundary conditions, create average structures, strip subsets of the system, and perform calculations such as RMS fitting, measuring distances, B-factors, radii of gyration, radial distribution functions, and time correlations, among other actions and analyses. Both the PTRAJ and CPPTRAJ programs and source code are freely available under the GNU General Public License version 3 and are currently distributed within the AmberTools 12 suite of support programs that make up part of the Amber package of computer programs (see http://ambermd.org ). This overview describes the general design, features, and history of these two programs, as well as algorithmic improvements and new features available in CPPTRAJ.

Refinement of the Cornell et al. Nucleic acids force field based on reference quantum chemical calculations of glycosidic torsion profiles

We report a reparameterization of the glycosidic torsion χ of the Cornell et al. AMBER force field for RNA, χOL. The parameters remove destabilization of the anti region found in the ff99 force field and thus prevent formation of spurious ladder-like structural distortions in RNA simulations. They also improve the description of the syn region and the syn–anti balance as well as enhance MD simulations of various RNA structures. Although χOL can be combined with both ff99 and ff99bsc0, we recommend the latter. We do not recommend using χOL for B-DNA because it does not improve upon ff99bsc0 for canonical structures. However, it might be useful in simulations of DNA molecules containing syn nucleotides. Our parametrization is based on high-level QM calculations and differs from conventional parametrization approaches in that it incorporates some previously neglected solvation-related effects (which appear to be essential for obtaining correct anti/high-anti balance). Our χOL force field is compared with several previous glycosidic torsion parametrizations.

Molecular dynamics simulation of nucleic acids: successes, limitations, and promise.

In the last five years we have witnessed a significant increase in the number publications describing accurate and reliable all‐atom molecular dynamics simulations of nucleic acids. This increase has been facilitated by the development of fast and efficient methods for treating the long‐range electrostatic interactions, the availability of faster parallel computers, and the development of well‐validated empirical molecular mechanical force fields. With these technologies, it has been demonstrated that simulation is not only capable of consistently reproducing experimental observations of sequence specific fine structure of DNA, but also can give detailed insight into prevalent problems in nucleic acid structure, ion association and specific hydration of nucleic acids, polyadenine tract bending, and the subtle environmental dependence of the A‐DNA–B‐DNA duplex equilibrium. Despite the advances, there are still issues with the methods that need to be resolved through rigorous controlled testing. In general, these relate to deficiencies of the underlying molecular mechanical potentials or applied methods (such as the imposition of true periodicity in Ewald simulations and the need for energy conservation), and significant limits in effective conformational sampling. In this perspective, we provide an overview of our experiences, provide some cautionary notes, and provide recommendations for further study in molecular dynamics simulation of nucleic acids.

The flying ice cube: Velocity rescaling in molecular dynamics leads to violation of energy equipartition

This article describes an unexpected phenomenon encountered during MD simulations: velocity rescaling using standard protocols can systematically change the proportion of total kinetic energy (KE) found in motions associated with the various degrees of freedom. Under these conditions, the simulation violates the principle of equipartition of energy, which requires a mean kinetic energy of RT/2 in each degree of freedom. A particularly pathological form of this problem occurs if one does not periodically remove the net translation of (and rotation about) the center of mass. In this case, almost all of the kinetic energy is converted into these two kinds of motion, producing a system with almost no kinetic energy associated with the internal degrees of freedom. We call this phenomenon “the flying ice cube.” We present a mathematical analysis of a simple diatomic system with two degrees of freedom, to document the origin of the problem. We then present examples from three kinds of MD simulations, one being an in vacuo simulation on a diatomic system, one involving a low resolution model of DNA in vacuo, and the third using a traditional all‐atom DNA model with full solvation, periodic boundary conditions, and the particle mesh Ewald method for treating long‐range electrostatics. Finally, we discuss methods for avoiding the problem. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 726–740, 1998

A systematic molecular dynamics study of nearest-neighbor effects on base pair and base pair step conformations and fluctuations in B-DNA

It is well recognized that base sequence exerts a significant influence on the properties of DNA and plays a significant role in protein–DNA interactions vital for cellular processes. Understanding and predicting base sequence effects requires an extensive structural and dynamic dataset which is currently unavailable from experiment. A consortium of laboratories was consequently formed to obtain this information using molecular simulations. This article describes results providing information not only on all 10 unique base pair steps, but also on all possible nearest-neighbor effects on these steps. These results are derived from simulations of 50–100 ns on 39 different DNA oligomers in explicit solvent and using a physiological salt concentration. We demonstrate that the simulations are converged in terms of helical and backbone parameters. The results show that nearest-neighbor effects on base pair steps are very significant, implying that dinucleotide models are insufficient for predicting sequence-dependent behavior. Flanking base sequences can notably lead to base pair step parameters in dynamic equilibrium between two conformational sub-states. Although this study only provides limited data on next-nearest-neighbor effects, we suggest that such effects should be analyzed before attempting to predict the sequence-dependent behavior of DNA.

Molecular Dynamics Simulations of the Dynamic and Energetic Properties of Alkali and Halide Ions Using Water-Model-Specific Ion Parameters

The dynamic and energetic properties of the alkali and halide ions were calculated using molecular dynamics (MD) and free energy simulations with various different water and ion force fields including our recently developed water-model-specific ion parameters. The properties calculated were activity coefficients, diffusion coefficients, residence times of atomic pairs, association constants, and solubility. Through calculation of these properties, we can assess the validity and range of applicability of the simple pair potential models and better understand their limitations. Due to extreme computational demands, the activity coefficients were only calculated for a subset of the models. The results qualitatively agree with experiment. Calculated diffusion coefficients and residence times between cation−anion, water−cation, and water−anion showed differences depending on the choice of water and ion force field used. The calculated solubilities of the alkali−halide salts were generally lower than the true solubility of the salts. However, for both the TIP4PEW and SPC/E water-model-specific ion parameters, solubility was reasonably well-reproduced. Finally, the correlations among the various properties led to the following conclusions: (1) The reliability of the ion force fields is significantly affected by the specific choice of water model. (2) Ion−ion interactions are very important to accurately simulate the properties, especially solubility. (3) The SPC/E and TIP4PEW water-model-specific ion force fields are preferred for simulation in high salt environments compared to the other ion force fields.

DNA Basepair Step Deformability Inferred from Molecular Dynamics Simulations

The sequence-dependent DNA deformability at the basepair step level was investigated using large-scale atomic resolution molecular dynamics simulation of two 18-bp DNA oligomers: d(GCCTATAAACGCCTATAA) and d(CTAGGTGGATGACTCATT). From an analysis of the structural fluctuations, the harmonic potential energy functions for all 10 unique steps with respect to the six step parameters have been evaluated. In the case of roll, three distinct groups of steps have been identified: the flexible pyrimidine-purine (YR) steps, intermediate purine-purine (RR), and stiff purine-pyrimidine (RY). The YR steps appear to be the most flexible in tilt and partially in twist. Increasing stiffness from YR through RR to RY was observed for rise, whereas shift and slide lack simple trends. A proposed measure of the relative importance of couplings identifies the slide-rise, twist-roll, and twist-slide couplings to play a major role. The force constants obtained are of similar magnitudes to those based on a crystallographic ensemble. However, the current data have a less complicated and less pronounced sequence dependence. A correlation analysis reveals concerted motions of neighboring steps and thus exposes limitations in the dinucleotide model. The comparison of DNA deformability from this and other studies with recent quantum-chemical stacking energy calculations suggests poor correlation between the stacking and flexibility.