Physics-based potentials for the coupling between backbone- and side-chain-local conformational states in the UNited RESidue (UNRES) force field for protein simulations

Abstract

The UNited RESidue (UNRES) model of polypeptide chains is a coarse-grained model in which each amino-acid residue is reduced to two interaction sites, namely, a united peptide group (p) located halfway between the two neighboring α-carbon atoms (Cαs), which serve only as geometrical points, and a united side chain (SC) attached to the respective Cα. Owing to this simplification, millisecond molecular dynamics simulations of large systems can be performed. While UNRES predicts overall folds well, it reproduces the details of local chain conformation with lower accuracy. Recently, we implemented new knowledge-based torsional potentials (Krupa et al. J. Chem. Theory Comput. 2013, 9, 4620–4632) that depend on the virtual-bond dihedral angles involving side chains: Cα···Cα···Cα···SC (τ(1)), SC···Cα···Cα···Cα (τ(2)), and SC···Cα···Cα···SC (τ(3)) in the UNRES force field. These potentials resulted in significant improvement of the simulated structures, especially in the loop regions. In this work, we introduce the physics-based counterparts of these potentials, which we derived from the all-atom energy surfaces of terminally blocked amino-acid residues by Boltzmann integration over the angles λ(1) and λ(2) for rotation about the Cα···Cα virtual-bond angles and over the side-chain angles χ. The energy surfaces were, in turn, calculated by using the semiempirical AM1 method of molecular quantum mechanics. Entropy contribution was evaluated with use of the harmonic approximation from Hessian matrices. One-dimensional Fourier series in the respective virtual-bond-dihedral angles were fitted to the calculated potentials, and these expressions have been implemented in the UNRES force field. Basic calibration of the UNRES force field with the new potentials was carried out with eight training proteins, by selecting the optimal weight of the new energy terms and reducing the weight of the regular torsional terms. The force field was subsequently benchmarked with a set of 22 proteins not used in the calibration. The new potentials result in a decrease of the root-mean-square deviation of the average conformation from the respective experimental structure by 0.86 Å on average; however, improvement of up to 5 Å was observed for some proteins.

Fig. 1

The UNRES model of polypeptide chains. The interaction sites are peptide-group centers (p), and side-chain centers (SC) attached to the corresponding α-carbons with different C^α ⋯ SC bond lengths, d_SC. The peptide groups are represented as dark gray circles and the side chains are represented as light gray ellipsoids of different size. The α-carbon atoms are represented by small open circles. The geometry of the chain can be described by the virtual-bond lengths, backbone virtual-bond angles θ_i, i = 1, 2,…, n−2, backbone virtual-bond-dihedral angles γ_i, i = 1, 2,…, n−3, and the angles α_i and β_i, i = 2, 3,…, n − 1 that describe the location of a side chain with respect to the coordinate frame defined by ${\mathrm{C}}_{i-1}^{\alpha },{\mathrm{C}}_i^{\alpha }$, and ${\mathrm{C}}_{i+1}^{\alpha }$.

Fig. 2

Illustration of backbone torsional angle γ (A, red) and side-chain backbone torsional angles τ⁽¹⁾ (A, green), τ⁽²⁾ (B, red), τ⁽³⁾ (B, green).

Fig. 3

Illustration of variables used for calculation of potential energy surfaces with the example of the lysine residue.

Fig. 4

Cartoon representation of truncated experimental structures of training set proteins: A: 1BDD, B: 1CLB, C: 1E0G, D: 1E0L, E: 1GAB, F: 1KOY, G: 1POU, H: 1PRU. The chains are colored from blue (N-terminus) to red (C-terminus).

Fig. 5

Plots of average side-chain-backbone potentials (A–C), [$\overline{U_{SC-\mathit{\text{corr}}}^{XY}}({\tau }^{(m)})$, m = 1, 2, 3; eq 7]; and sample side-chain backbone correlation potentials: U_AlaAla(τ^(m)), m = 1, 2, 3 (D–F), U_MetMet(τ^(m)), m = 1, 2, 3, (G–I), and U_ProPro(τ^(m)), m = 1, 2, 3 (J–L). Black circles represent the values of the dimensionless PMFs calculated from histograms (eq 3) and red lines represent one-dimension Fourier series fits (eq 5) to the PMF values.

Fig. 6

Color plots of the correlation coefficients of the respective $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$ potentials with the respective $\overline{U_{SC-\mathit{\text{corr}}}^{XY}}({\tau }^{(m)})$ potentials averaged over all residue types at positions X and Y (eq 6). A: $r_{XY}^{(1)}$ (for the C^α ⋯ C^α ⋯ C^α ⋯ SC angles, τ⁽¹⁾ B: $r_{XY}^{(2)}$ (for the SC ⋯ C^α ⋯ C^α ⋯ C^α angles, τ⁽²⁾), C: $r_{XY}^{(3)}$ (for the SC ⋯ C^α ⋯ C^α ⋯ SC angles, τ⁽³⁾). The color scales are on the right bars of each panel. Types of residues X are on the abscissae and types of residues Y are on the ordinates. The color scale is shown on the right side of each panel.

Fig. 7

Sample convergence plots of the heat capacity of 1E0L for w_SC−corr = 0.0 (standard value), w_tor = 1.84316 (standard value) and w_tord = 1.26571 (standard value) (A), w_SC−corr = 0.25, w_tor = 1.84316 (standard value) and w_tord = 1.26571 (standard value) (B), w_SC−corr = 0.25, w_tor = 1.34316 and w_tord = 1.26571 (standard value) (C), and w_SC−corr = 1.0, w_tor = 0.0, w_tord = 1.26571 (standard value) (D). Different colors denote heat-capacity curves for consecutive windows of the MREMD simulation, for the range from 10,000 to 50,000,000 MD steps divided into 8 equal windows. Windows are colored in order: red, orange, yellow, green, cyan, blue, purple, black.

Fig. 8

Bar diagrams of various RMSDs averaged over 8 proteins. White bars show the lowest RMSD obtained during the corresponding MREMD run [ρ^min of eq 12], light-grey bars show the minimum of the cluster-averaged RMSD [${\langle \rho \rangle }_{\mathit{\text{clust}}}^{\mathit{\text{min}}}(T_a)$ of eq 11] from five temperatures of clustering (210, 240, 270, 290 and 310 K), dark-grey bars show the minimum of the cluster-averaged RMSD [${\langle \rho \rangle }_{\mathit{\text{clust}}}^{\mathit{\text{min}}}(T_a)$ of eq 11] from the temperature of clustering 10 K lower than the heat capacity peak, and the black bars show the RMSD averaged over the conformational ensemble generated during the MD run by WHAM in the last part of the simulation. S stands for w_SC−corr, T stands for w_tor, D stands for w_tord.

Fig. 9

Plot of distances between C^α atoms of the average structures of the most native-like cluster of simulated 1CLB structures from the respective C^α atoms of the experimental 1CLB structure after optimal superposition. Solid black line: calculations with standard force field parameters (w_SC−corr = 0.25, w_tor = 1.34316, w_tord = 1.26571), dashed black line: the force field that includes the $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$, m = 1, 2, 3 potentials derived in this work (w_SC−corr = 0.25, w_tor = 1.34316, w_tord = 1.26571). Red horizontal lines on the abscissa mark α-helices in the experimental structure.

Fig. 10

Plot of distances between C^α atoms of the average structure of the most native-like cluster of simulated 1KOY structure from the respective C^α atoms of the experimental 1KOY structure after optimal superposition. Solid black line: calculations with standard force field parameters (w_SC−corr = 0.25, w_tor = 1.34316, w_tord = 1.26571), dashed black line: the force field that includes the $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$, m = 1, 2, 3 potentials derived in this work (w_SC−corr = 0.25, w_tor = 1.34316, w_tord = 1.26571). Red horizontal lines on the abscissa mark α-helices in the experimental structure.

Fig. 11

A: Superposition of the C^α-trace of the average structure of the most probable cluster of conformations of 1CLB obtained in MREMD simulations with inclusion of the $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$ potentials derived in this work (green lines) on that of the experimental structure of 1CLB (red lines). B: the average structure of the most probable cluster of conformations of 1CLB obtained in MREMD simulations with w_SC−corr = 0.0, w_tor = 1.84316, w_tord = 1.26571 (without new potentials). The RMSDs from the experimental structures are 6.04 Å for panel A and 7.90 Å for panel B, respectively

Fig. 12

A: Superposition of the C^α-trace of the average structure of the most probable cluster of conformations of 1KOY obtained in MREMD simulations with inclusion of the $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$ potentials derived in this work (green lines) on that of the experimental structure of 1KOY (red lines). B: the average structure of the most probable cluster of conformations of 1KOY obtained in MREMD simulations with w_SC−corr = 0.0, w_tor = 1.84316, w_tord = 1.26571 (without new potentials). The RMSDs from the experimental structure are 5.25 Å for panel A, and 9.23 Å for panel B, respectively.

Fig. 13

Bar diagrams of the lowest RMSDs from the experimental structures for 22 proteins from testing set. White bars: the lowest RMSD obtained in MREMD simulations [ρ^min of eq 12] with the UNRES force field without new terms. Light-grey bars: ρ^min with the best set of new terms.

Fig. 14

Bar diagrams of the RMSDs of the mean structures of the most native-like clusters [${\langle \rho \rangle }_{\mathit{\text{clust}}}^{\mathit{\text{min}}}(T_a)$ of eq 11] for the 22 test proteins corresponding to clustering at T = 290 K. White bars: the UNRES from force field without new terms. Light-grey bars: the UNRES force field with the new terms and optimal energy-term weights (w_SC−corr = 0.25, w_tor = 1.34316, w_tord = 1.26571). For each protein and each force field, the error bar represents the standard deviations of the RMSDs of the structures of the most native-like cluster from the RMSD of the average structure of that cluster.

Fig. 15

A: Superposition of the C^α-trace of the average structure of the most probable cluster of conformations of 1FEX obtained in MREMD simulations with inclusion of the $U_{SC-\mathit{\text{corr}}}^{XY}({\tau }^{(m)})$ potentials derived in this work (green lines) on that of the experimental structure of 1FEX (red lines). B: the average structure of the most probable cluster of conformations of 1FEX obtained in MREMD simulations with w_SC−corr = 0.0, w_tor = 1.84316, w_tord = 1.26571 (without new potentials). The RMSDs from the experimental structure are 7.93 Å for panel A, and 9.15 Å for panel B, respectively.