Development of a new physics-based internal coordinate mechanics force field and its application to protein loop modeling

Abstract

We report the development of internal coordinate mechanics force field (ICMFF), new force field parameterized using a combination of experimental data for crystals of small molecules and quantum mechanics calculations. The main features of ICMFF include: (a) parameterization for the dielectric constant relevant to the condensed state (ε = 2) instead of vacuum, (b) an improved description of hydrogen-bond interactions using duplicate sets of van der Waals parameters for heavy atom-hydrogen interactions, and (c) improved backbone covalent geometry and energetics achieved using novel backbone torsional potentials and inclusion of the bond angles at the C(α) atoms into the internal variable set. The performance of ICMFF was evaluated through loop modeling simulations for 4-13 residue loops. ICMFF was combined with a solvent-accessible surface area solvation model optimized using a large set of loop decoys. Conformational sampling was carried out using the biased probability Monte Carlo method. Average/median backbone root-mean-square deviations of the lowest energy conformations from the native structures were 0.25/0.21 Å for four residues loops, 0.84/0.46 Å for eight residue loops, and 1.16/0.73 Å for 12 residue loops. To our knowledge, these results are significantly better than or comparable with those reported to date for any loop modeling method that does not take crystal packing into account. Moreover, the accuracy of our method is on par with the best previously reported results obtained considering the crystal environment. We attribute this success to the high accuracy of the new ICM force field achieved by meticulous parameterization, to the optimized solvent model, and the efficiency of the search method.

Figure 1

Junction of the flexible and static segments of the polypeptide chain in loop simulations. The arrows show virtual bonds that are parts of the internal coordinate trees of the two segments. Virtual C and C^α atoms at the termini of the two segments are constrained to their physical counterparts. These constraints, in conjunction with rigid covalent geometry within the two segments, maintain (near) ideal geometry of the physical C^α-C bond.

Figure 2

(a) nonbonded plus electrostatics φ/ψ energy map for Ace-Ala-NMe; (b) QM φ/ψ energy map for Ace-Ala-NMe. The color code from purple to red corresponds to the 0–8 kcal/mol range; (c) heatmap of the deviations between the QM and MM (computed without torsion potential) energies for Ace-Ala-NMe. The size of the squares indicates the frequency (on a logarithmic scale) of occurrence of a particular φ/ψ value pair (within a 10°×10° bin); (d) Ramachandran plot for a set of 21 diverse ultra-high resolution structures (resolution between 0.5 and 0.8Å, PDBs 1ejg, 1ucs, 2vb1, 1us0, 2dsx, 1r6j, 2b97, 1x6z, 1gci, 1pq7, 1iua, 2ixt, 1w0n, 2h5c, 1nwz, 1n55, 2o9s, 2jfr, 2pwa, 2o7a and 2hs1), blue points. Background is colored by φ/ψ frequencies for a much broader protein set, calculated within each 10°×10° degrees square on the φ/ψ plane. (e) contour map of the φ/ψ torsion potential. Color code from purple to red corresponds to the 0–4 kcal/mol range; (f) final energy map including the torsion potential. The color code from purple to red corresponds to the 0–8 kcal/mol range.; (g) heatmap of the residual deviations between the QM and final MM energies for Ace-Ala-NMe. Size of the squares indicates the frequency (on a logarithmic scale) of occurrence of a particular φ/ψ value pair (within a 10°×10° bin). Contours in Fig. 2a,b,e,g are drawn with 1kcal/mol step.

Figure 3

(a) QM and (b) MM maps for terminally blocked glycine. The color code from purple to red corresponds to the 0–20 kcal/mol range. Contours are drawn with 1kcal/mol step.

Figure 4

Histogram of the distribution of proline φ angle in 21 diverse ultra-high resolution protein structures (161 proline residues, average φ angle is −68.8°, see Legend to Fig. 2d).

Figure 5

φ/ψ map for Ace-Pro-Nme: (a) QM energy map for trans-proline; (b) QM energy map for cis-proline; (c) Ramachandran map for trans- and cis-proline; (d) final MM energy map for trans-proline; (e) final MM energy map for cis-proline. The color code from purple to redof the energy maps corresponds to the 0–20 kcal/mol range. Contours are drawn with 1kcal/mol step

Figure 6

Overlay of the native (orange) structure and the lowest-energy (green) conformation predicted using ICMFF: (a) for the 12-residue loop in 1oth (residues 69-80); (b) for the 13-residue loop in 1p1m (residues 327-339); (c) for the 7-residue loop in 2rn2. The native loop conformation is stabilized by interactions with symmetric molecules (shown in gray); (d) for the 11-residue loop in 2eng. Red spheres represent water molecules stabilizing the native loop conformation.

Figure 7

Overlay of the native (orange) structure and the lowest-energy (green) conformation predicted using ICMFF for the 9-residue loop containing cis-proline (2ixt, residues 69-80).

Figure 8

Comparison of loop prediction results obtained with ICM (solid line) and reported for HLP (dashed line). The y-axis is the percentage of loops for which the backbone RMSD of the lowest-energy conformation is at or below the RMSD on the x-axis.