Comparison of multiple Amber force fields and development of improved protein backbone parameters

Abstract

The ff94 force field that is commonly associated with the Amber simulation package is one of the most widely used parameter sets for biomolecular simulation. After a decade of extensive use and testing, limitations in this force field, such as over-stabilization of alpha-helices, were reported by us and other researchers. This led to a number of attempts to improve these parameters, resulting in a variety of "Amber" force fields and significant difficulty in determining which should be used for a particular application. We show that several of these continue to suffer from inadequate balance between different secondary structure elements. In addition, the approach used in most of these studies neglected to account for the existence in Amber of two sets of backbone phi/psi dihedral terms. This led to parameter sets that provide unreasonable conformational preferences for glycine. We report here an effort to improve the phi/psi dihedral terms in the ff99 energy function. Dihedral term parameters are based on fitting the energies of multiple conformations of glycine and alanine tetrapeptides from high level ab initio quantum mechanical calculations. The new parameters for backbone dihedrals replace those in the existing ff99 force field. This parameter set, which we denote ff99SB, achieves a better balance of secondary structure elements as judged by improved distribution of backbone dihedrals for glycine and alanine with respect to PDB survey data. It also accomplishes improved agreement with published experimental data for conformational preferences of short alanine peptides and better accord with experimental NMR relaxation data of test protein systems.

Figure 1

Blocked alanine tetrapeptide (three alanine residues but four peptide bonds) used in the optimization procedure. Note the definition of dihedral angles: φ =C-N-C^α-C, ψ =N-C^α-C-N, φ′ = C-N-C^α-C^β, ψ′ = C^β-C^α-C-N. There are three φ/ψ pairs and three additional φ′/ψ′ pairs in this tetrapeptide. Glycine tetrapeptide would only have the φ/ψ dihedrals due to the absence of C^β carbon. As an example, ψ₁ and φ′₃ are shown in bold.

Figure 2

Ramachandran plot of 28 Gly₃ (left) and 51 Ala₃ (right) conformers as obtained from QM geometry optimization of structures that were acquired using stochastic search with MM energies. All three dihedral pairs are shown for each conformer (black squares). The gray points in the background are φ/ψ values of glycine and alanine residues collected from a subset of the PDB. Typical secondary structure regions are outlined by boxes in the Ala₃ plot (see Methods for additional details). Ala₃ conformer number 16 (outlier) in the lower left of Ala₃ plot is designated by an arrow.

Figure 3

Free energy φ/ψ maps for Gly₃ (top row) and Ala₃ (bottom row) from 80ns simulations with explicit TIP3P water. The energies are color coded from 0 up to 5 kcal/mol. PDB survey data is represented by a simple Ramachandran plot. Individual force fields are designated as: ff99SB (this work), ff03, ff94, ff99 and ff94gs.

Figure 4

“Lowest energy” profiles for three prototypical decoy systems, each tested with six AMBER force field variants (ff94, ff99, ff99SB, ff03, ff94gs, ff99ϕ): (A) trpzip2, which is known to adopt a β-hairpin, (B) a Baldwin-type sequence as a representative of an α-helix, and (C) trpcage as a representative of mixed secondary structures. RMSD values are calculated using the experimentally determined structure as a reference. Ideally, a force field should show lowest energies for the lowest RMSD values.

Figure 5

Order parameters (S²) derived from experiment (black line) and calculated from MD simulations in explicit water with different force field parameter sets: ff99 (green), ff99SB (red), ff94 (blue) and ff03 (magenta). Error bars reflect the convergence of calculated S² values. (A) Lysozyme. Secondary structures of lysozyme are labeled: helix A (HA: residues 4–15), loop 1 (L1: 16–23), helix B (HB: 24–36), strand 1 (S1: 41–45), turn 1 (T1: 46–49), strand 2 (S2: 50–53), strand 3 (S3: 58–60), long loop 2 (L2: 61–78), 3₁₀ helix 1 (H1: 80–84), loop 3 (L3: 85–89), helix C (HC: 89–99), loop 4 (L4: 100–107), helix D (HD: 108–115), loop 5 (L5: 116–119) and 3₁₀ helix 2 (H2: 120–124). (B) Ubiquitin. Secondary structure elements in ubiquitin: strand 1 (S1: residues 1–7), turn 1 (7–10), strand 2 (S2: 10–17), turn 2 (T2: 18–21), helix 1 (H1: 23–34), turn 3 (T3: 37–40), strand 3 (S3: 41–44), turn 4 (T4: 45–48), turn 5 (T5: 51–54), helix 2 (H2: 56–59), turn 6 (T6: 62–65) and strand 4 (S4: 66–71). Differences are described in the text.