Impact of hydrodynamic interactions on protein folding rates depends on temperature

We investigated the impact of hydrodynamic interactions (HI) on protein folding using a coarse-grained model. The extent of the impact of hydrodynamic interactions, whether it accelerates, retards, or has no effect on protein folding, has been controversial. Together with a theoretical framework of the energy landscape theory (ELT) for protein folding that describes the dynamics of the collective motion with a single reaction coordinate across a folding barrier, we compared the kinetic effects of HI on the folding rates of two protein models that use a chain of single beads with distinctive topologies: a 64-residue $\alpha /\beta$ chymotrypsin inhibitor 2 (CI2) protein, and a 57-residue $\beta$-barrel $\alpha$-spectrin Src-homology 3 domain (SH3) protein. When comparing the protein folding kinetics simulated with Brownian dynamics in the presence of HI to that in the absence of HI, we find that the effect of HI on protein folding appears to have a “crossover” behavior about the folding temperature. This means that at a temperature greater than the folding temperature, the enhanced friction from the hydrodynamic solvents between the beads in an unfolded configuration results in lowered folding rate; conversely, at a temperature lower than the folding temperature, HI accelerates folding by the backflow of solvent toward the folded configuration of a protein. Additionally, the extent of acceleration depends on the topology of a protein: for a protein like CI2, where its folding nucleus is rather diffuse in a transition state, HI channels the formation of contacts by favoring a major folding pathway in a complex free energy landscape, thus accelerating folding. For a protein like SH3, where its folding nucleus is already specific and less diffuse, HI matters less at a temperature lower than the folding temperature. Our findings provide further theoretical insight to protein folding kinetic experiments and simulations.

Figure 1

Representations of protein models of chymotrypsin inhibitor 2 (CI2) (top row) and the $\alpha$-spectrin Src-homology 3 (SH3) domain (bottom row). The protein models are in [(a), (b)] a cartoon representation, [(c), (d)] a $C_{\alpha }$-only representation, and [(e), (f)] a protein topology cartoon created with Pro-origami [32]. Arrows are $\beta$ strands and cylinders are helices. Structures from [(a), (b), (c), (d)] were created with VMD [33]. The secondary structures were assigned using DSSP [34]. The key residues for the hydrophobic core (A16, L49, and I57) and the minicore (L32, V38, and F50) are represented with green and orange beads, respectively, for CI2 in (a). Residues of the diverging turn (M20, K21, K22, G23, and D24) and distal loop (N42 and D43) are represented with green and orange for SH3 in (b).

Figure 2

The average folding time $t_{\mathrm{fold}}$ in units of reduced time $\tau$ with respect to temperature for CI2 and SH3 in the presence or in the absence of hydrodynamic interactions (HI). Panels (a) and (b) show $t_{\mathrm{fold}}$ over a broad range of temperatures using the Brownian dynamics without HI (BD) for CI2 and SH3, respectively. The temperature for each protein is expressed in units of their corresponding folding temperature $T_f$. Note the U-shaped dependence of the folding time (non-Arrhenius behavior). $t_{\mathrm{fold}}$ using BD with HI (BDHI) is compared to $t_{\mathrm{fold}}$ using only BD in panel (c) for CI2 and (d) for SH3. The crossover occurs when the two curves intersect. Error bars are calculated using the jackknife method.

Figure 3

Free energy ($F$) in units of $k_{\mathrm{B}}T$ with respect to fraction of native contact formation $Q$ without or with HI (BD and BDHI, respectively) for (a) CI2 and (b) SH3 at a temperature below the folding temperature ($T_f$), at $T_f$, and above $T_f$. The transition state region is indicated by TS. Error bars are included.

Figure 4

The mean square displacement of the fraction of the native contact formation $Q$, ${\langle \mathrm{\Delta }{Q}^2\rangle }_{t_0,\mathrm{\Omega }}$, as a function of lag time in units of reduced time $\tau$ in the absence or presence of HI [BD (solid black line) and BDHI (dash-dot red line), respectively] at $T<T_f$ [(a) $0.95T_f^{CI2}$ and (b) $0.91T_f^{SH3}$ for CI2 and SH3, respectively], and $T>T_f$ [(c) $1.06T_f^{CI2}$ and (d) $1.03T_f^{SH3}$ for CI2 and SH3, respectively]. ${\langle \mathrm{\Delta }{Q}^2\rangle }_{t_0,\mathrm{\Omega }}$ was averaged over all initial times $t_0$ and all trajectories $\mathrm{\Omega }$. Shaded width of solid lines represents the error. The inset zooms in to the range of time used for the linear fit of data. Open circles (BD) and open squares (BDHI) are placed every 100 data points for $T<T_f$ and 200 data points for $T>T_f$. The unit for the effective diffusion coefficient is $1/\tau$ since $Q$ is dimensionless.

Figure 5

Impact of HI on key pairs of secondary structure interactions $\mathrm{\Delta }P\left(t\right)$ at $T<T_f$ [(a) $0.95T_f^{CI2}$ and (b) $0.91T_f^{SH3}$ for CI2 and SH3, respectively], and $T>T_f$ [(c) $1.06T_f^{CI2}$ and (d) $1.03T_f^{SH3}$ for CI2 and SH3, respectively]. Time is normalized with respect to maximum simulation time ($9\times 10^6\tau$) and shown in log scale. The dashed, vertical lines correspond to the time average of ${\langle Q\rangle }_{\mathrm{\Omega }}=0.4$ for CI2 and SH3. ${\langle Q\rangle }_{\mathrm{\Omega }}=0.4$ is the onset of the transition state. The curves were smoothed by running averaging over every 500 data points. Legends on curves are displayed for visual guidance. Errors are included but too small to be visible. At the top of the panels, the location of secondary structure elements and unstructured regions along the sequence of CI2 (left) and SH3 (right) is indicated as visual guidance.

Figure 6

Temporal evolution of the fraction of native contact formation $Q$ averaged over all trajectories $\mathrm{\Omega }$ as function of normalized time $t/t_{\max }$ at $T<T_f$ and $T>T_f$ for (a) CI2 and (b) SH3. Time is normalized with respect to maximum simulation time ($9\times 10^6\tau$). Shaded width of lines represents the error, which is calculated using a jackknife method.

Figure 7

(a) The probability of contact formation of each native pair $Q_{ij}\left(t\right)$ (upper triangle) and the displacement correlation $C_{ij}\left(t\right)$ (lower triangle) at $0.95T_f^{CI2}$ for CI2 at a normalized time where ${\langle Q\rangle }_{\mathrm{\Omega }}=0.4$. (b) The distribution of native pairs from panel (a) with $\left|i-j\right|\ge 10$. (c) The same presentation as (a) for $Q_{ij}\left(t\right)$ and $C_{ij}\left(t\right)$ at $1.06T_f^{CI2}$. (d) The distribution of native pairs from panel (c) with $\left|i-j\right|\ge 10$. The arrows and rectangles along (a) and (c) represent $\beta$ strands and helices, respectively. The orange and purple dashed boxes highlight the long-range contacts between the N terminus and C terminus, and the contacts in neighborhood of the minicore, respectively. The black and green dashed boxes highlight the displacement correlation between the N terminus and C terminus, and between the N terminus and the $\alpha$ helix, respectively. (e) The displacement correlation at $0.95T_f^{CI2}$ and $1.06T_f^{CI2}$ for all pairs are classified into three sets based on the magnitude of positive correlations. If the magnitude of the correlation is similar at both temperatures, the pair is colored with a green edge on the left structure. If the magnitude of the correlation is greater at $0.95T_f^{CI2}$ than that at $1.06T_f^{CI2}$ a pair is colored with a blue edge on the right structure. If the magnitude of the correlation is greater at $1.06T_f^{CI2}$ than that at $0.95T_f^{CI2}$ a pair is colored with a red edge on the right structure. Only the pairs with sequence separation greater than 8 residues and magnitude of displacement correlation above the threshold $\mu =0.061$ are considered for this representation. The key residues for the hydrophobic core (A16, L49, and I57) and the minicore (L32, V38, and F50) are illustrated with green and orange beads, respectively. (f) $\mathrm{\Pi }\left(\right|i-j\left|\right)$ for all pairs whose magnitude of the displacement correlation is above $\mu$ are organized according to the sequence separation $\left|i-j\right|$.

Figure 8

(a) The probability of contact formation of each native pair $Q_{ij}\left(t\right)$ (upper triangle) and the displacement correlation $C_{ij}\left(t\right)$ (lower triangle) at $0.91T_f^{SH3}$ for SH3 a normalized time where ${\langle Q\rangle }_{\mathrm{\Omega }}=0.4$. (b) The distribution of native pairs from panel (a) with $\left|i-j\right|\ge 10$. (c) The same representation as (a) for $Q_{ij}\left(t\right)$ and $C_{ij}\left(t\right)$ at $1.03T_f^{SH3}$. (d) The distribution of native pairs from panel (c) with $\left|i-j\right|\ge 10$. The arrows and rectangles along (a) and (c) represent $\beta$ strands and helices, respectively. The orange and purple dashed boxes highlight the long-range contacts between the N terminus and C terminus, and the contacts between the diverging turn (DT) and the distal loop (DL), respectively. The black and green dashed boxes highlight the displacement correlation between the N terminus and C terminus, and between the N terminus and both DT and $\beta 2$, respectively. (e) The displacement correlation at $0.91T_f^{SH3}$ and $1.03T_f^{SH3}$ for all pairs are classified into three sets based on the magnitude of positive correlations. The coloring rules are the same as Fig. 7 with pairs with sequence separation greater than 7 residues and magnitude of displacement correlation above the threshold $\mu =0.095$. (f) $\mathrm{\Pi }\left(\right|i-j\left|\right)$ for all pairs whose magnitude of the displacement correlation is above $\mu$ are organized according to the sequence separation $\left|i-j\right|$.

Figure 9

The Brownian motion of residues without HI (BD) shows small and random displacement correlation for native and non-native pairs of CI2 and SH3. Upper and lower triangles represent the probability of contact formation for each native pair $Q_{ij}\left(t\right)$ and displacement correlation $C_{ij}\left(t\right)$, respectively. The panels are plotted at $T<T_f$ [(a) $0.95T_f^{CI2}$ and (b) $0.91T_f^{SH3}$ for CI2 and SH3, respectively], and $T>T_f$ [(c) $1.06T_f^{CI2}$ and (d) $1.03T_f^{SH3}$ for CI2 and SH3, respectively] at a normalized time where ${\langle Q\rangle }_{\mathrm{\Omega }}=0.4$. The arrows and rectangles represent $\beta$ strands and helices, respectively.

Figure 10

Two representative kinetic pathways are projected on a two-dimensional free energy landscape of the fraction of native contact formation $Q$ and the fraction of native contact formation in the region of the minicore $Q_T$ for CI2 at $0.95T_f^{CI2}$. Panel (a) shows a major pathway of folding kinetics (route I), representing a fast route, in the presence of HI. Panel (b) shows a minor pathway (route II). The folding free energy was colored in gray scale in units of $k_{\mathrm{B}}T$. The kinetic trajectories were colored by a normalized time (time divided by their respective first passage time) and projected on the folding free energy. Key conformations were selected for visual guidance. The significant residues for the hydrophobic core (Ala16, Leu49, and Ile57) and minicore (Leu32, Val38, and Phe50) are illustrated with green and orange beads, respectively. The minicore forms before the hydrophobic core in panel (a), whereas the opposite occurs in panel (b). Structures were created with VMD [33].

Figure 11

One representative kinetic pathway is projected on a two-dimensional free energy landscape of the fraction of native contact formation $Q$ and the fraction of native contact formation in the region formed by the diverging turn (DT) and the distal loop (DL) $Q_T$ for SH3 at $0.91T_f^{SH3}$. It shows a dominant pathway of folding kinetics. The folding free energy was colored in gray scale in units of $k_{\mathrm{B}}T$. The kinetic trajectory was colored by a normalized time (time divided by their respective first passage time) and projected on the folding free energy. Key conformations were selected for visual guidance. The residues of the diverging turn (M20, K21, K22, G23, and D24) and distal loop (N42 and D43) are illustrated with green and orange beads, respectively. The native contacts between DT and DL are formed before reaching the top of the barrier. Structures were created with VMD [33].