Cutoff Size Need Not Strongly Influence Molecular Dynamics Results for Solvated Polypeptides † David A. C. Beck, ‡ Roger S. Armen, ‡ and Valerie Daggett* ,‡,§ Biomolecular Structure and Design Program, UniVersity of Washington, Seattle, Washington 98195-7610, and Department of Medicinal Chemistry, UniVersity of Washington, Seattle, Washington 98195-7610 ReceiVed June 29, 2004; ReVised Manuscript ReceiVed October 26, 2004 ABSTRACT: The correct treatment of van der Waals and electrostatic nonbonded interactions in molecular force fields is essential for performing realistic molecular dynamics (MD) simulations of solvated polypeptides. The most computationally tractable treatment of nonbonded interactions in MD utilizes a spherical distance cutoff (typically, 8-12 Å) to reduce the number of pairwise interactions. In this work, we assess three spherical atom-based cutoff approaches for use with all-atom explicit solvent MD: abrupt truncation, a CHARMM-style electrostatic shift truncation, and our own force-shifted truncation. The chosen system for this study is an end-capped 17-residue alanine-based R-helical peptide, selected because of its use in previous computational and experimental studies. We compare the time-averaged helical content calculated from these MD trajectories with experiment. We also examine the effect of varying the cutoff treatment and distance on energy conservation. We find that the abrupt truncation approach is pathological in its inability to conserve energy. The CHARMM-style shift truncation performs quite well but suffers from energetic instability. On the other hand, the force-shifted spherical cutoff method conserves energy, correctly predicts the experimental helical content, and shows convergence in simulation statistics as the cutoff is increased. This work demonstrates that by using proper and rigorous techniques, it is possible to correctly model polypeptide dynamics in solution with a spherical cutoff. The inherent computational advantage of spherical cutoffs over Ewald summation (and related) techniques is essential in accessing longer MD time scales. Accurate molecular dynamics (MD) 1 simulations of polypeptides in solution necessitate robust methodologies. The selection of certain methodologies, such as the use of all-atom molecular mechanics potentials and explicit repre- sentation of fully flexible waters, is important for obtaining correct, experimentally verifiable results (1). However, some aspects, such as the choice of treatments for the nonbonded atomic interactions, are up for debate (2-9). As such, it is important that they be continually re-examined and carefully tested against experiment. Furthermore, MD simulations using force fields originally derived from, and parametrized against, static states (i.e., crystal structures) must be evaluated against experiments that reflect the dynamic properties of a system. As computer power increases, it is necessary to continually test and evaluate MD methods using successively longer simulations (10). These methods should be tested for stability (e.g., conservation of conserved properties), accuracy (e.g., proper treatment of short- and long-range forces), and completeness (i.e., conformational sampling). For these purposes, peptides that are well characterized by circular dichroism (CD) and solution nuclear magnetic resonance (NMR) spectroscopy are useful systems for testing the accuracy of these methods against the respective ensemble- averaged experimental observables. Further, because of their small size, explicitly solvated peptide systems are compu- tationally less expensive than identically treated proteins. The inherently dynamic nature of peptides provides a good test of the sampling capabilities of MD. For example, within the time scale of the molecular dynamics trajectory, are the experimentally averaged data (e.g., from both “folded” and “unfolded” conformations) adequately reproduced? In this work, we use robust MD methods (1) to evaluate the stability, accuracy, and completeness of several spherical cutoff treatments. Previous tests for accurate treatment of long-range electrostatic interactions in alanine-based helical peptides have used the stability of the R-helix as the primary criterion for correct treatment of interactions (4-6). For example, in a highly cited paper, Schreiber and Steinhauser (4) focused on the Y(KAAAA) 3 K-NH 2 peptide (11) and compared simulations computed using an abrupt group-based truncation, with 6, 10, and 14 Å spherical cutoffs for the nonbonded interactions, with simulations calculated using Ewald summations. The simulations they performed were as short as 60 ps and as long as 825 ps; all were conducted at 300 K. They observed a lack of convergence of the helical content (derived primarily from visual inspection of struc- † This work was supported by the National Institutes of Health (Grant GM 50789 to V.D.). D.A.C.B. and R.S.A. were supported by an NIH Molecular Biophysics Training Grant (National Research Service Award 5 T32 GM 08268). * To whom correspondence should be addressed. E-mail: daggett@ u.washington.edu. ‡ Biomolecular Structure and Design Program. § Department of Medicinal Chemistry. 1 Abbreviations: MD, molecular dynamics; CD, circular dichroism; NMR, nuclear magnetic resonance; PME, particle mesh Ewald; NVE, constant number of particles, volume, and energy; ilmm, in lucem Molecular Mechanics; NVEp, constant number of particles, volume, energy, and linear momentum; rmsd, root-mean-square deviation. 609 Biochemistry 2005, 44, 609-616 10.1021/bi0486381 CCC: $30.25 © 2005 American Chemical Society Published on Web 12/18/2004