© Copyright 1997 by the American Chemical Society
VOLUME 101, NUMBER 38, SEPTEMBER 18, 1997
LETTERS Experiment vs Force Fields: DNA Conformation from Molecular Dynamics Simulations Michael Feig and B. Montgomery Pettitt* Department of Chemistry, UniVersity of Houston, Houston, Texas 77204-5641 ReceiVed: April 4, 1997; In Final Form: June 26, 1997X
Molecular dynamics simulations of the DNA duplex d(CCCCCTTTTT)2 in ionic solution with AMBER and CHARMM force fields are compared. CHARMM parameters result in an A-DNA conformation with a heterogeneous backbone structure, but in the AMBER force field the DNA is found in a B-like helix with a homogeneous B-DNA backbone during the first 3.5 ns before moving close to the CHARMM structure in the A family at around 4 ns. Both structures are in partial disagreement with experiment, which reports the C/G part in A and the A/T part in B form.
The local structure of DNA serves as an important structural recognition motif in the biological context.1-3 Numerous experimental and theoretical studies have been carried out over many years to understand the influence of base pair sequence and environment on DNA structure,4-6 but many features are not yet fully understood. Molecular dynamics simulations can provide detailed model structures and short time dynamics of DNA in solution. The recent introduction of Ewald summations for electrostatic interactions7,8 has been shown to produce stable trajectories of fully solvated DNA in the nanosecond time range with structures close to crystallographic and NMR data.9-11 However, the convergence of the same DNA dodecamer under similar conditions, but with different force fields (CHARMM12 vs AMBER13), to different structures in A and B form14,15 indicates difficulties in accurate modeling of DNA conformation. In order to evaluate the force field influence on the DNA structure in the absence of other protocol differences with molecular dynamics simulations, we performed simulation runs of the nucleic acid decamer d(CCCCCTTTTT)2 with CHARMM12 and AMBER13 parameters under otherwise identical simulation conditions using a single program and compared them to experiment. X
Abstract published in AdVance ACS Abstracts, August 15, 1997.
S1089-5647(97)01180-2 CCC: $14.00
In intermediate to high-salt solutions, nonalternating C/G sequences have been shown to favor A-DNA,16,17 while nonalternating A/T tracts generally prefer B-DNA conformation.18 The combination of both homopolymeric sequences has been shown to form an A/B junction under high-salt conditions for the sequence d(GGGGGTTTTT)219 and very recently for d(CCCCCTTTTT)2 in water/TFE mixtures.20 Starting from model-built canonical A-DNA,21 d(CCCCCTTTTT)2 was simulated in the NVT ensemble at 300 K with the most recent CHARMM and AMBER all atom nucleic acid force fields12,13 in explicit TIP3P water22 and Na/Cl ions. The simulation program from our group23 employs periodic boundary conditions, the velocity Verlet integration scheme24 and SHAKE25 to constrain chemical bonds. The integration time step was set to 2 fs. Electrostatic interactions were calculated using the Ewald summation technique with a multiple time step algorithm.26 The simulation boxes have a size of 4.0 × 4.0 × 5.0 nm containing 2285 water molecules, 18 Na counterions, and 32 additional Na/Cl ions. The equilibration procedure described in detail previously11,14 consists of an initial 20-step steepest descent minimization followed by alternating runs with either solvent or solute kept fixed while reassigning velocites every 50 steps from a canonical Maxwell distribution at 300 K for © 1997 American Chemical Society
7362 J. Phys. Chem. B, Vol. 101, No. 38, 1997
Figure 1. Averaged pseudo rotation angle of backbone sugar ring excluding edge and junction base pairs in AMBER (s) and CHARMM (‚‚‚) simulations.
150 ps. Simulation runs were carried out for 4 ns with CHARMM and 4.5 ns with AMBER parameters. We expect the results for the converged structures to be independent of the starting structure. Simulations started from initial structures in A and B conformations have been shown to converge to similar structures with the AMBER force field15 and with CHARMM parameters.27 The DNA phosphate-sugar backbone conformation is most easily described by the pseudo rotation angle28 characteristic of the furanose ring pucker. C3′-endo modes are observed for A-DNA and O2′-endo through C2′-endo are typical of B-DNA conformation. Figure 1 shows the time evolvement of the average pseudo rotation for the four different base types excluding edge and C-T junction base pairs in both simulations. With the AMBER force field the sugars of all bases quickly adopt typical B-type C1′-exo to C2′-endo conformations of 120°-140°. Different base specific behavior is found in the CHARMM force field simulation. The adenine and guanine sugars remain in C3′-endo throughout the whole simulation. On the pyrimidine strand, the sugars at the thymine bases fluctuate from C3′-endo to C1′-exo/C2′-endo after 2 ns but the equilibrium lies toward the A conformation, while the cytosine sugars convert completely to C1′-exo/C2′-endo at 2 ns after fluctuating for 1 ns. Most characteristic of A and B-type base geometries are the inclination of a base with respect to the helical axis and the rise between successive base pairs as shown in Figure 2 for the C/G and A/T tracts of both simulations. Average values over the whole simulation time after 500 ps from the CHARMM force field are typical of A-DNA base geometry with 11.6° inclination and 2.8 Å rise. AMBER parameters produced more B-like base stacking with a smaller inclination angle of 4.0° and a larger 3.2 Å rise. However, toward the end of the
Letters
Figure 2. Average base inclination to helical axis and rise between successive base steps calculated with NEWHEL9335 in AMBER (s) and CHARMM (‚‚‚) simulations.
simulation at around 4 ns, both force fields generate similar structures in the A family. Average inclination and rise for the last 500 ps of each simulation are 8.7° and 3.0 Å with CHARMM and 7.6° and 3.1 Å with AMBER. In comparison to experimental results, the CHARMM force field simulates the C/G sequence close to the expected A-DNA conformation with an inclination of 10° and a rise of 2.9 Å29,30 but with the cytosine base backbone in B form. However, the A/T sequence found also in A-DNA conformation in the simulation is in apparent disagreement with available experimental data. AMBER produces on average a B-like helix, with a B-type backbone and a base geometry in the B-DNA family, closer to the experimental data for A/T tracts (0° to -5° inclination and 3.4 Å rise,31-33 but disagreeing with experiment in the C/G part. From this comparison we find very different equilibria of A and B conformations with AMBER and CHARMM force fields in nanosecond molecular dynamics simulations. During the first 3.5 ns CHARMM generates A-DNA base geometries in contrast with B-DNA structures from the AMBER simulations. However, the approach to a more A-like helix close to the CHARMM structures during the last 1 ns in the AMBER simulations might suggest even longer convergence time scales of several nanoseconds for the AMBER force field. Different behavior was also found for the backbone structures. While AMBER exhibits exclusively B-DNA conformations, CHARMM displays basetype sensitive backbone conformations with an A-form backbone at purine bases and a B-form backbone at pyrimidine bases. These results clearly impose difficulties on the interpretation of DNA structures in current molecular dynamics simulations. Previous results on A/B transitions observed with molecular dynamics studies14,15 must be viewed in the context of these results. Improvements in DNA modeling with the discussed force fields are primarily expected from modifications in the backbone
Letters parameters, particularly in the sugar pseudo rotation angle potential, to better reproduce experimental data, as theoretical studies34 have shown how helix geometries are driven by backbone conformations. A detailed analysis of our comparisons from the presented simulations, extended to 10 ns each, is in preparation. Acknowledgment. We thank the Robert A. Welch Foundation, the NIH, and the National Science Foundation for partial support of this research. We further thank Paul E. Smith and Gillian C. Lynch for valuable discussions. The Metacenter is thanked for computational resources, and the Institute for Molecular Design and MSI are thanked for graphics support. References and Notes (1) von Hippel, Peter H. Science 1994, 263, 769. (2) Eisenstein, Miriam; Shakked, Zippora. J. Mol. Biol. 1995, 248, 662. (3) Calladine, C. R.; Drew, H. R. Curvature and flexibility of DNA: Sequence-directed effects seem from a structural mechanics viewpoint. In Mol. Struct. Life 1992, 43. (4) Dickerson, R. E. Methods Enzymol. 1991, 211, 67. (5) Kochoyan, Michel; Louis Leroy, Jean. Curr. Opin. Struct. Biol. 1995, 5, 329. (6) Louise-May, Shirley; Auffinger, Pascal; Westhof, Eric. Curr. Opin. Struct. Biol. 1996, 6, 289. (7) Ewald, P. P. Ann. Phys. 1921, 64, 253. (8) Toukmaji, Abdulnour Y.; Board, John A. Comp. Phys. Commun. 1996, 95, 73. (9) York, D. M.; Yang, W.; Lee, H.; Darden, T. A.; Pedersen, L. J. Am. Chem. Soc. 1995, 117, 5001. (10) Cheatham, T. E.; Miller, J. L.; Fox, T.; Darden, T. A.; Kollman, P. A. J. Am. Chem. Soc. 1995, 117, 4193. (11) Weerasinghe, Samantha; Smith, Paul E.; Mohan, V.; Cheng, Y.K.; Pettitt, B. Montgomery J. Am. Chem. Soc. 1995, 117, 2147. (12) MacKerell, Alexander D., Jr.; Wiorkiewicz-Juczera, Joanna; Karplus, Martin. J. Am. Chem. Soc. 1995, 117, 11946. (13) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1995, 117, 5179.
J. Phys. Chem. B, Vol. 101, No. 38, 1997 7363 (14) Yang, Liqiu; Pettitt, B. Montgomery. J. Phys. Chem. 1996, 100, 2564. (15) Cheatham, T. E.; Kollman, P. A. J. Mol. Biol. 1996, 259, 434. (16) Nishimura, Yoshifumi; Torigoe, Chikako; Tsuboi, Masamichi. Nucleic Acids Res. 1986, 14, 2737. (17) Peticolas, W. L.; Wang, Y.; Thomas, G. A. Proc. Natl. Acad. Sci. U.S.A. 1988, 85, 2579. (18) Taillandier, E.; Ridoux, J.-P.; Liquier, J.; Leupin, W.; Denny, W. A.; Wang, Y.; Thomas, G. A.; Peticolas, W. I. Biochemistry 1987, 26, 3361. (19) Wang, Yang; Thomas, Gerald A.; Peticolas, Warner L. J. Biomol. Struct. Dyn. 1989, 6, 1177. (20) Ivanonv, Valery I.; Minchenkova, Lucy E.; Burckhardt, Gu¨nther; Birch-Hirschfeld, Eckhard; Fritzsche, Hartmut; Zimmer, Christoph. Biophys. J. 1996, 71, 3344. (21) QUANTA, Release 4.0, 1994. (22) Jorgensen, W.; Chandrasekhar, J.; Madura, J.; Impey, R.; Klein, M. J. Chem. Phys. 1983, 79, 926. (23) Smith, Paul E.; Holder, Michael E.; Dang, Liem X.; Feig, Michael; Pettitt, B. Montgomery. ESP; University of Houston, Houston, TX, 1996. (24) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids, 1st ed.; Oxford University: New York, 1987. (25) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (26) Smith, Paul E.; Pettitt, B. M. Comput. Phys. Commun. 1995, 91, 339. (27) Feig, Michael; Pettitt, B. Montgomery. Unpublished results, 1996. (28) Saenger, W. Principles of Nucleic Acid Structure; Springer, New York, 1984. (29) McCall, M.; Brown, T.; Kennard, O. J. Mol. Biol. 1985, 183, 385. (30) Haran, T. E.; Shakked, Z.; Wang, A. H.-J.; Rich, A. J. Biomol. Struct. Dyn. 1987, 5, 199. (31) Alexeev, D. G.; Lipanov, A. A.; Skuratovskii, I. Ya. Nature 1987, 325, 821. (32) Coll, Miquel; Frederick, Christin A.; Wang, Andrew H.-J.; Rich, Alexander. Proc. Natl. Acad. Sci. U.S.A. 1987, 84, 8385. (33) Nelson, H. C. M.; Finch, J. T.; Luisi, B. F.; Klug, A. Nature 1987, 330, 221. (34) Tung, C.-S. J. Biomol. Struct. Dyn. 1992, 9, 1185. (35) Dickerson, R. E. Unpublished work, University of California at Los Angeles, CA, 1992.
