Configurational diffusion for lattice models in proteins

Are these general ideas developed by Bryngelson and Wolynes (BW) valid for quantitatively predicting folding times in model proteins with a realistic energy landscape topography? We show in this section that as long as the glass transition falls after the transition region (top of the barrier in the free energy profile for the collective reaction coordinate) that this is the case. In this limit the single dominant funnel picture is appropriate. The system under study is the designed three-letter code 27-mer lattice model (see figure 3) used in our recent studies.[10, 11, 3]

  
Figure 3: The native state (minimum energy) cube for the sequence studied in this work. The sequence ( ABABBBCBACBABABACACBACAACAB) consist of three monomer types.

Although this 27-mer is far from being a real protein, it has been shown[3] that by developing a law of corresponding states, the results obtained for such a system can be used to describe a small -helical protein. In these simulations the units of temperature and energy are such that . This still leaves an arbitrary scale factor since the only important quantity is the ratio of energy to temperature. We have chosen small, integer values for the coupling energies for convenience and efficiency. The kinetic glass temperature is .

  
Figure 4: The energy and order parameter (Q, the number of correct native contacts) plotted as a function of time for a sample folding simmulation. The temperture is the folding temperature ( ). The time is roughly 30 times the folding time ( ). The plot shows the two-state-like behavior of this system with transition between the native state (E=-84, Q=28) and a molten globule region ( , ). The transition are extremly rapid with respect to the folding time. In addition there are significant fluctations about the two states.

Figure 4 shows the time evolution at for energy and the Q order parameter. This Q parameter is a measure of similarity to the native state and is equal to the number of correct contacts (i.e. contacts that exist in the native state). It ranges from 0, no correct contacts, to 28, the maximum number of possible contacts. Figure 5 plots the energy, entropy and free energy as a function of Q for various temperatures. At temperatures above the glass temperature, there is a significant gradient in both energy and entropy.

 

 

Figure: Energy, entropy and free energy as a function of Q for several temperatures calculated using the Monte Carlo histogram method [11]. For a broad range of temperatures the free energy has a clear double well form: one at the native state (Q=28) and one at the molten globule region ( ). The double well form is consistent with the two state behavior scene in figure 4. Above the glass transition temperature ( ) there is a significant energy and entropy gradient between the molten globule region and the transition region ( ). Both and are less than zero, so the unfavorable reduction in entropy is offset by a loss in energy, leaving a small free energy barrier at the folding temperature (T=1.509). As the temperature approaches the glass temperature the energy and entropy gradients decrease as does the free energy barrier. However the folding time diverges due to the increase in the diffusion constant of the system (i.e., the roughness of the energy landscape).

These profiles are computed from the density of states obtained using the Monte Carlo histogram technique.[11] Starting from a random configuration, collapse occurs at times very short compared to the folding time. Thus in this parameter range, the radius of gyration need not be considered as a separate dynamical reaction coordinate; however, in determining the free energies one must note that the mean radius of gyration does vary with temperature. A molten globule band, where configurations have an average of 20 contacts and Q=7.5 ( similarity to the native state), describes the region where this sequence spends most of its time on its way into the folded state. The shape of the free energy in the vicinity of this mean Q is quasi- harmonic, although in fact large deviations are possible at high temperatures. As shown in our earlier work, since at the local glass transition as a function of Q occurs after the transition region has been overcome, the folding time is determined primarily by the time taken to overcome the free energy barrier, i.e., to cross the transition region. Figure 6 shows a folding trajectory of the Q coordinate superimposed on a plot of the free energy.

  
Figure 6: A folding trajectory ( ) superimposed on the free energy (F[Q]) at the folding temperature. The trajectory is approximately Monte Carlo steps long (approximately 27% of ). The right axis shows the time counted backwards from the folded state. Note, most of the time is spent diffusing in the molten globule region. Once the barrier is surmounted folding proceeds rapidly. Numerous recrossings of the barrier region indicate the need to use a diffusive theory for this reaction coordinate note transition state theory.

The plot shows a time span which is roughly one quarter of the mean folding time ( ) at this temperature. Most of the trajectory consists of diffusive motion about the molten globule region. Once the barrier has been surmounted folding occurs rapidly, taking roughly Monte Carlo steps ( ).

However, to estimate the folding times at a variety of temperatures, knowledge of the free energy barrier alone is not sufficient. Information about the dynamics must be obtained by calculating the configurational diffusion coefficient through the complex energy landscape. As described in the previous section, when a single reaction coordinate is considered, for example Q, can be computed using the double integral give by eq. 2. In general the diffusion coefficient will depend on Q but, one more simplification is assumed here. Only the average value of D, computed for states in the molten globule band, is inferred from simulations. We do this by computing the correlation function of the fluctuations of the reaction coordinate . Within the quasi-harmonic diffusive approximation this correlation function should decay exponentially at long times. The configurational diffusion coefficient will be related to the corresponding correlation time and the mean square instantaneous value of the reaction coordinate fluctuations, . (This calculation is constrained for values of Q before the transition region. ) At a given temperature, the diffusive harmonic model gives . Since , it is a quantity that is much easier to obtain from numerical simulations for more complex systems such as atomistic simulations of proteins. For experimentalists we especially wish to point out this viewpoint and the corresponding approximation allow one to use dynamic probes of fluctuations in the molten globule state at equilibrium to predict the rate of the non-equilibrium folding process.

In figure 7 we show the simulated and diffusion coefficient for various temperatures above the glass transition.

  
Figure 7: Autocorrelation time ( ) and diffusion coefficient ( ) plotted as a function of the inverse temperature ( ). At low temperature (high ) the diffusion coefficient decreases rapidly while the correlation time increases, indicating a slow down of the dynamics due to trapping in local minimum.

The autocorrelation time as been calculated in the molten globule region. To do this, the correlation function <Q(t)Q(0)> was calculated over paths which were confined to the molten globule region of the conformation space. For these calculations we define any conformation with Q<17 as being in the molten globule region. The decays are in fact non-exponential for short times. We have not analyzed the non-exponentiality in detail. The information would provide the dynamics of individual conformation steps and should contain information on the distribution of substate lifetimes. We identify the long time decay with the effective diffusive harmonic value. A more complete model should account for this short time decay via a frequency dependent diffusion coefficient as proposed by BW. Frequency effects may be important for determining the prefactor of barrier crossing as discussed by Grote and Hynes for simple chemical reactions.[19]