Observation of complete pressure-jump protein refolding in molecular dynamics simulation and experiment.

Density is an easily adjusted variable in molecular dynamics (MD) simulations. Thus, pressure-jump (P-jump)-induced protein refolding, if it could be made fast enough, would be ideally suited for comparison with MD. Although pressure denaturation perturbs secondary structure less than temperature denaturation, protein refolding after a fast P-jump is not necessarily faster than that after a temperature jump. Recent P-jump refolding experiments on the helix bundle λ-repressor have shown evidence of a <3 μs burst phase, but also of a ~1.5 ms "slow" phase of refolding, attributed to non-native helical structure frustrating microsecond refolding. Here we show that a λ-repressor mutant is nonetheless capable of refolding in a single explicit solvent MD trajectory in about 19 μs, indicating that the burst phase observed in experiments on the same mutant could produce native protein. The simulation reveals that after about 18.5 μs of conformational sampling, the productive structural rearrangement to the native state does not occur in a single swift step but is spread out over a brief series of helix and loop rearrangements that take about 0.9 μs. Our results support the molecular time scale inferred for λ-repressor from near-downhill folding experiments, where transition-state population can be seen experimentally, and also agrees with the transition-state transit time observed in slower folding proteins by single-molecule spectroscopy.

Introduction

Pressure is a fundamental thermodynamic variable that can modulate protein structure, dynamics, and function.[2−5] The effect of hydrostatic pressure on protein stability has been studied extensively,[6−8] including many studies that have explored how protein denaturation occurs under high-pressure conditions.[9−12] An ellipse-shaped zone of stability in the temperature–pressure plane has been proposed on the basis of theoretical models, molecular dynamics (MD) simulations, and experiments.[13−15] In addition to studying the stability of proteins in equilibrium under high-pressure conditions, pressure perturbation can also be used to study protein folding kinetics.[16,17] Microsecond pressure-jump (P-jump) techniques have recently extended the range of such experiments to near the folding “speed limit”.[18,19]

Here we ask whether such fast protein refolding from the pressure-denatured state can be seen from beginning to end in a full-atom explicit solvent MD simulation, and whether the speed limit from the pressure-denatured state has the same time scale as the ∼1 μs speed limit observed in temperature-jump (T-jump) experiments and single-molecule experiments.[20,21] The speed limit, “molecular time scale” or “transition state crossing time”, is slower than backbone torsional transitions or end-to-end contacts of short loops, and attributed to either reduced diffusion on a rough energy landscape, or equivalently to extremely short-lived high-free-energy folding intermediates.[22−24]

MD simulation has been used for studying the effect of pressure on protein thermodynamics[6,14,25] and protein denaturation,[26−28] but kinetic refolding and the molecular time scale have proved elusive so far. Recent developments in both hardware and software have pushed MD simulations up to the millisecond time domain.[25−29] At the same time, a new generation of microsecond kinetics P-jump experiments can access the MD time scale,[18] and would benefit greatly from interpretation based on the atomistic detail available from simulations. So far simulations have only observed the trapping process associated with non-native helix formation after a P-jump,[19] and it remains to be seen whether the burst phase observed in experiments translates into fast and complete refolding in silico.

We previously applied the ultrafast P-jump technology to study the refolding kinetics of an alanine-rich mutant of five-helix bundle λ-repressor, λ*YA. In addition to a <3 μs burst phase, λ*YA refolding exhibited a ∼1.5 ms “slow” phase attributed to off-pathway helix-rich traps. We surmised that two glycines, mutated in place of alanines in helices 2 and 3 of λ*YA, would make the protein more flexible and help it escape from traps with excess helical structure, making it a good test case for complete refolding upon P-jump in silico. Following a T-jump near its melting temperature, this new glycine-rich mutant of λ-repressor, λ*YG, is known to fold with a rate coefficient kf ≈ (22 μs)−1.[30]

We computed over 33 μs of explicit solvent dyanmics in several long trajectories to simulate a P-jump experiment. High-pressure-denatured states, generated through a high-temperature unfolding and high-pressure equilibration simulation procedure,[19] were found to contain significant residual helical structure. Nonetheless, λ*YG refolded into the native state in under 20 μs following a pressure drop.

However, follow-up P-jump experiments on λ*YG are in agreement with previous P-jump experiments on the alanine-rich λ repressor mutant λ*YA. We still see a small <3 μs burst phase and a large 1.7 ms slow phase, the latter previously attributed to non-native helix formation.[19] Thus glycine substitution in helices 2 and 3 did not eliminate the slow phase. The microsecond folding observed here in silico therefore suggests that a fraction of the proteins with native-like residual helix in the unfolded state refolds very rapidly during the experimentally observed microsecond burst phase, while the remaining population with non-native helix content in the turns is trapped for >1 ms.

After 18.6 μs of conformational search, the simulation revealed a 0.9 μs stretch of productive structural assembly, bracketed by three-helix alignment and loop formation motions that were almost, but not quite, concerted. This 0.9 μs time scale agrees with the molecular time scale of 1–2 μs measured by T-jump experiments on near-downhill folding λ-repressor mutants.[20,31]

Materials and Methods

Protein

λ*YG consists of the wild-type λ-repressor sequence residues 6–85, with mutations Tyr22Trp, Gln33Tyr, Ala37Gly, and Ala49Gly.[32] It has a melting temperature of Tm = 55 °C and a folding time constant of 22 μs near Tm in aqueous buffer.[30] For MD simulations the sequence was used without C-terminal amidation or N-terminal acetylation. The initial structure of the λ-repressor fragment was taken from the Protein Data Bank (PDB code 1LMB).[33] From a crystal structure of the similar λ*YA mutant, we know that the mutant and wild-type native structures are quite similar.[1] For P-jump experiments, the protein was grown in E. coli and purified and lyophilized as described previously.[32]

Molecular Dynamics Simulations

MD simulations were performed in explicit solvent using the TIP3P water model[34] and the CHARMM22 force field with CMAP corrections for protein and ions.[35−37] The force field has excess helix propensity,[26,38,39] which may accelerate trapping (from non-native helix) and folding (from native helix). The initial protein structure was placed in a cubic box of 24 282 water molecules at 55 mM NaCl salinity, neutralized with extra ions employing VMD.[40] The simulated system, including protein, water molecules, and ions, measured 91.1 Å in each dimension at T = 325 K and P = 1 bar and contained 74 197 atoms. All simulations were carried out with periodic boundary conditions in a constant particle number, temperature, and pressure ensemble (NPT), in five steps. Step (1): in a P-jump simulation, pressure was increased from 1 bar to 5 kbar in 0.15 μs at a rate of 1 bar/30 ps while maintaining the temperature at T = 325 K. Step (2): in a high-temperature and high-pressure unfolding simulation, temperature was increased to 525 K while maintaining pressure at P = 5 kbar, running the simulation for 0.15 μs. Step (3): after the protein unfolded in step (2), the temperature was dropped back to 325 K and pressure was kept at P = 5 kbar while the denatured protein was equilibrated at high pressure for 1 μs in a high-pressure equilibrium simulation. Step (4): in a pressure-drop simulation, pressure was jumped downward from 5 kbar to 1 bar in 0.15 μs at a rate of −1 bar/30 ps while maintaining the temperature at T = 325 K. Steps (1)–(4) were carried out on general purpose supercomputers using NAMD 2.9.[41] Step (5): the resulting pressure-denatured state under refolding conditions was employed as the initial state for a refolding simulation, carried out on the special purpose supercomputer Anton[42,43] for 32 μs. Constant temperature (T = 325 K) and constant pressure (P = 1 bar) were maintained during the refolding simulation.

MD Simulations Using NAMD

The simulation algorithm and features of the NAMD program are described in ref (41). The system to be simulated was first subjected to 6000 steps of conjugate gradient minimization and equilibrated for 300 ps with harmonic restraints applied to all the heavy atoms of the protein. The simulation was then continued for 3 ns without restraints at a constant pressure of 1 bar using Nosé–Hoover Langevin piston barostat and at a constant temperature of 325 K with Langevin damping constant of 5.0 ps–1. In the subsequent simulations of steps (1)–(4), constant temperature was maintained using Langevin dynamics with a damping constant of 1.0 ps–1 and multiple time stepping employed with an integration time step of 2.0 fs, short-range forces being evaluated every time step and long-range electrostatics evaluated every three time steps. Cutoff for short-range nonbonded interactions was 8.0 Å; long-range electrostatics was calculated using the particle-mesh Ewald method.[44] All bonds involving hydrogen in the protein were constrained using RATTLE,[45] while the geometries of water molecules were maintained using SETTLE.[46]

MD Simulations on Anton

The refolding simulation in step (5) was carried out on the Anton platform.[42,43] Multiple time stepping was employed, with an integration time step of 2.0 fs. Short-range forces were evaluated every time step and long-range electrostatics every three time steps. Cutoff for the short-range nonbonded interactions was 9.28 Å; long-range electrostatics was calculated using the k-Gaussian Split Ewald method[47] with a 64 × 64 × 64 grid. All bonds involving hydrogen atoms were constrained using SHAKE.[48]

Pressure-Jump Experiment

Refolding kinetics experiments were performed on a home-built P-jump apparatus as described previously.[18,19] Briefly, an 8–10 μL dimple was machined into a sapphire cube with a side length of 3/8-in. (Esco Products, Oak Ridge, NJ). The sample consisting of 300 μM protein in 50 mM phosphate buffer at pH 7 with either 0 or 1 M guanidine hydrochloride (GuHCl) was then pipetted into the dimple and sealed with a double-layer of Mylar-coated aluminum foil to prevent mixing between the sample and the pressurization fluid (water). The sealant foil lay 2 mm below a 0.007-in. stainless steel burst membrane, to which it was connected by a pressurization channel. The sample and burst membrane were pressurized hydrostatically to 1.2 kbar using a pressure pump (High Pressure Equipment Company, Erie, PA). The burst membrane was ruptured by passing ∼10 kA of current (95 V) through it, releasing the sample pressure back to 1 bar within 2–3 μs.

The sample was optically excited with a frequency-tripled Ti:sapphire laser (KMLabs, Boulder, CO), which generated femtosecond pulses of 285 ± 3 nm light separated by 12.5 ns. Fluorescence was collected and the photons were delivered to a photomultiplier (R7400U-03, Hamamatsu Corp., Bridgewater, NJ) using an optical waveguide (Oriel Instruments, Stratford, CT). We used a band-pass filter (B370, Hoya, Santa Clara, CA) to avoid interference from the excitation light. The signal was recorded and digitized at 100 ps time resolution using an oscilloscope with a 2.5 GHz bandwidth (DPO7254, Tektronix, Beaverton, OR). The amplitude of the signal was usually on the order of 100–250 mV. The data were analyzed as described previously.[19]

Results

We chose λ-repressor mutant λ*YG (Y22W/Q33Y/A37G/A49G) as a model system to study complete fast protein refolding after a pressure drop. λ-repressor is a five-helix bundle protein with 80 amino acids. It is the largest fast-folding protein folded in silico to date by all-atom MD simulations.[25,28,29] Fast folding of various λ-repressor mutants has been studied previously using temperature-jump,[1,20,31,32,49] pressure-jump,[19,50] and rapid microfluidic mixing techniques.[51]

MD Simulation of the Native State

In order to later compare with the denatured simulation and refolding simulation, we first performed a 0.3 μs MD simulation of λ*YG at T = 325 K and P = 1 bar, starting from the crystal structure.[33] This 0.3 μs simulation will be referred to hereafter as the native simulation. The average values of several structural characteristics, such as radius of gyration (Rgyr), were determined from the native simulation and defined as the protein’s native values, shown as red solid lines in Figure 1.

Figure 1

Structural characterization of the λ*YG unfolding trajectory. Cα-RMSD values have been calculated relative to the crystal structure 3KZ3.[1] α-content is the fraction of residues that are in the α-helical conformation, and Rgyr is the radius of gyration. The native values, calculated from a 0.3 μs equilibrium simulation of the native structure at T = 325 K and P = 1 bar, are shown as red solid lines. The pressure applied through the simulation, shown as the color background, varies from 1 bar (white) to 5 kbar (blue). The temperature is kept at 325 K, except for the time window between 0.15 and 0.3 μs, where 525 K is used to unfold the protein. See Materials and Methods for a description of the unfolding procedure. Representative structures at different time points are shown at the top. Protein coloring runs blue to red from the N-terminus to the C-terminus. The unfolding trajectory is rendered as Movie S1 in the SI.

Upward Pressure-Jump Simulation

The denaturation simulation followed a procedure described in a previous P-jump MD simulation study that did not observe refolding of the λ*YA mutant.[19] Briefly, we started with the native state of λ*YG, shown in the t = 0 μs conformation in Figure 1. The pressure was gradually increased from 1 bar to 5 kbar over 0.15 μs, while temperature was held constant at T = 325 K. In Figure 1, the value of the pressure is depicted by the background color changes from white (1 bar) to blue (5 kbar). The protein remains in its native conformation in the first 0.15 μs of upward P-jump simulation. High pressure can unfold a protein,[9,11] but such high-pressure denaturation is a slow process that takes place on a time scale of seconds or even longer. Therefore, 0.15 μs of pressurizing is too short for observing any discernible conformational change.

Extensive Denaturation at High Pressure and High Temperature

To accelerate the protein unfolding process, we heated the system to T = 525 K and simulated the system for another 0.15 μs, while keeping pressure at P = 5 kbar. As shown in Figure 1, the protein rapidly unfolds as evidenced by the increase of Cα-root-mean-squared deviation (RMSD) relative to the crystal structure (>20 Å). The content of secondary structure, α-helix in particular, drops from the native value of 65.5% to a value in the 10–30% range. During the unfolding, the protein also assumes some extended conformations with Rgyr of more than 30 Å. The high-T-P denatured state, obtained after the high-temperature and high-pressure unfolding simulation, is shown in Figure 1 as the conformation at t = 0.3 μs.

Equilibration at High Pressure and Room Temperature

The high temperature used in the simulation unfolds the protein, but also likely disrupts the protein more than when only high pressure is used for denaturation. To obtain a state more representative of the pressure-denatured ensemble, the high-T-P denatured state was equilibrated for 1 μs at P = 5 kbar and T = 325 K. The most striking observation in the equilibration, shown in Figure 1, is that the α-helix content recovers from ∼30.0% to ∼60.0%, which is already close to the native value of ∼65.5%. The existence of high α-helix content at high pressure indicates that pressure denaturation is mainly breaking the tertiary contacts, but does not perturb the secondary structure considerably, as proposed previously.[19] The result is consistent with the finding that pressure does not affect the helix–coil equilibrium significantly, based on replica exchange MD simulations of α-helical peptide using a different force field.[25] Recent experiments by Neumaier et al. have shown that high pressure can slightly stabilize a helix, which explains the frequently observed helical structures in pressure-denatured proteins.[52] Notably, the pressure-denatured state after the high-pressure equilibration, shown in Figure 1 as the conformation at t = 1.3 μs, already contains helices 1 and 4.

Downward Pressure Jump Simulation

We performed a 0.15 μs downward P-jump simulation that initiated the refolding process. Significant amount of helical structure, resulting from the high-pressure equilibration, did not change substantially during the downward jump. The entire 1.45 μs denaturation simulation and downward P-jump simulation rendered as one trajectory are shown in Movie S1 in the Supporting Information (SI). The time evolution of the Cα-displacement per residue relative to the crystal structure and the secondary structure per residue are shown in Figure S1.

Protein Refolding Simulation

Following the downward P-jump simulation, we carried out a 32 μs MD simulation at P = 1 bar and T = 325 K to investigate fast refolding. The protein, except for the last helix (helix 5), folded into the native state after 19 μs. Representative snapshots along the folding pathway are shown in Figure 2A. The Cα-RMSD of the protein relative to the crystal structure (PDB code 3KZ3)[1] is shown in Figure 2B, and the secondary structure is shown in Figure 2C; see Figure S2 for additional quantities of interest). A particular mechanism by which λ*YG mutant folds in our trajectory is punctuated by two fast events separated by a longer conformational search, as shown in Figure 2D.

Figure 2

Protein refolding trajectory from the simulation at T = 325 K and P = 1 bar after pressure jump. The refolding trajectory is rendered as Movie S2 in SI. (A) Distinct molecular rearrangements observed at the bottleneck (transition-state passage). At each time point, the folded residues (Cα displacement relative to the crystal structure ≤2 Å) are colored blue. (B) Cα-RMSD values for the protein core (residues 7–80), calculated relative to the crystal structure 3KZ3.[1] The native range is defined by the mean value (red solid line) ± standard deviation (green dashed line) from a 0.3 μs equilibrium simulation of the native structure at T = 325 K and P = 1 bar. (C) Time evolution of the secondary structure throughout the trajectory. The secondary structure of the crystal structure is shown on the left side of the panel. (D) Time evolution of per-residue Cα displacements from the crystal structure throughout the refolding trajectory. (E) The time window between 18 and 20 μs is enlarged to reveal the sequence of rearrangements at the bottleneck (see Movie S3 in SI). The color bar runs from blue (close to crystal structure) to red (far from crystal structure).

In a first fast step within 2 μs of the P-jump, helices 1 (residues 7–27) and 4 (residues 58–70) adjusted their orientation and registration to reach a near-native conformation. This conformation remained fairly stable for the next 16 μs of conformational search. Helix 2 also formed individually very quickly, within 1 μs.

Subsequently, a number of factors prolonged the conformational search: helix 3 did not form individually until ∼7 μs (see Figure 2C), and was originally shifted toward the C-terminus by one full helical turn; neither helix 2 nor helix 3 acquired a native orientation; helical overshoots were observed both in helices 1 and 2, where they eroded loop 1 between the two helices as shown in Figure 2C, preventing the correct helix orientation from locking in.

The key bottleneck was crossed in the second fast step between 18.6 and 19.5 μs (Figure 2E). In 0.9 μs, three key rearrangements brought the protein from a compact denatured state through the transition-state region into the native basin (Figure 2A). First, between 18.6 and 18.8 μs, helices 2 and 3, along with loop 2 that connects them, reoriented themselves and assumed the native packing conformation relative to helices 1 and 4. The one-helical-turn shift in helix 3 and the helical overshoot in loop 1 disappeared during this structural transition. Within the next 0.4 μs, loop 1 assumed its native conformation, reached at t = 19.2 μs. Finally, loop 3 between helices 3 and 4 rearranged into its native structure at t = 19.5 μs. The whole refolding trajectory is visualized in Movie S2 in the SI. The trajectory for the time window between 18 and 20 μs is provided in Movie S3 in the SI.

Helix 5 Does Not Form in the Simulation

The formation of helix 5 is presumably the last step of the folding process; this step is not observed in the simulation. It is possible that the λ*YG in solution adopts a native state that is different from the crystal structure.[1] Indeed, an unstable helix 5 in the native state is consistent with the high B factors in the crystal structure[1] and observed in implicit solvent MD simulations.[53] An unstructured helix 5 has also been observed before for another λ-repressor mutant (D14A) and attributed to the absence of C-terminus residues from the wild-type λ-repressor.[28]

Protein Refolding in Pressure-Jump Experiments

Previous fast P-jump experiments up to 0.5 ms on λ*YG revealed a microsecond burst phase.[18] We recently extended the capabilities of our P-jump apparatus to collect refolding kinetic data for up to 5 ms.[19] We used this new home-built instrument to jump the pressure of the sample from 1.2 kbar to 1 bar (see Materials and Methods). The tryptophan in position 22 (Trp22) was used as a fluorescence probe to study folding of λ*YG after a microsecond pressure drop. Tyrosine in position 33 (Tyr33) was introduced to quench Trp22 fluorescence in the folded state, enhancing the change in fluorescence lifetime upon unfolding.

Titration of λ*YG with GuHCl shows an unfolding midpoint concentration of 1.3 M (Figure S3). Equilibrium pressure denaturation of λ*YG was monitored by fluorescence spectroscopy (Figure S4). Based on these results, we expect no pressure denaturation of λ*YG in 0 M GuHCl up to 1.2 kbar, whereas in 1 M GuHCl, the protein is poised for unfolding when the pressure is increased above 1 bar.

For the P-jump experiments, Trp fluorescence excited at 285 nm was sampled every 12.5 ns and digitized with a time resolution of 100 ps. The fluorescence lifetimes for NATA were then normalized from χ = 0 (before the P-jump) to χ = 1 (4.7 ms after the P-jump) through a linear fitting procedure, and the 0 and 1 M protein samples were analyzed on the same scale for comparison. The dead time of the instrument with a starting pressure of 1.2 kbar was determined from the NATA sample, fitting its step-function-like trace to a single exponential rise of ∼3 μs.

In general, λ variants with a Q33Y mutation exhibit an increase in fluorescence lifetime upon unfolding whether the denaturation is accomplished using temperature, pressure, or a chemical denaturant. This response can be rationalized in terms of nonradiative quenching of the tryptophan fluorescence by tyrosine in the folded state. Figure 3 shows the fluorescence-detected kinetics of NATA, λ*YG without GuHCl, and λ*YG in 1 M GuHCl from 1.2 kbar to 1 bar. λ*YG without GuHCl data is a folded control because λ*YG without GuHCl does not undergo pressure denaturation at 1.2 kbar. The initial lifetime increase for λ*YG in 1 M GuHCl at t = 0 is a factor of 1.4 larger than that observed for λ*YG without GuHCl. This observation indicates that there is a fast burst phase, <3 μs, superimposed on the intrinsic response of tryptophan to the P-jump, during which the lifetime increases. In addition, a slower phase, 1.74 ± 0.02 ms, was observed. This slow phase was not present in the control P-jumps of NATA or λ*YG without GuHCl. The absolute changes in the fluorescence decay of NATA and λ*YG without GuHCl in response to a P-jump are shown in Figure S5.

Figure 3

Pressure-jump of NATA (gray) and λ*YG in 1 M (blue) and 0 M (red) GuHCl from 1.2 kbar to 1 bar probed by tryptophan fluorescence decays. Tryptophan lifetime was normalized from χ = 0 (before the P-jump) to χ = 1 (5 ms after the jump) for NATA. P-jump data of λ*YG were analyzed on the same lifetime scale as NATA data to facilitate direct comparison. The panel on the left shows the data from 20 μs before the P-jump to 100 μs after the P-jump. The P-jump occurs at t = 0 μs. The fast phase was fitted to a single exponential function with time constants of τ = 2.6 ± 0.3 (NATA), 2.3 ± 0.4 (λ*YG with GuHCl), and 3.7 ± 0.5 μs (λ*YG without GuHCl). The panel on the right shows the data from 5 μs after the P-jump to ∼4.7 ms after the P-jump. The millisecond kinetic response of λ*YG with GuHCl was fitted to a single exponential function with a time constant τ = 1.74 ± 0.02 ms (black curve).

Effect of Helix 5 on Overall Protein Folding

Helix 5 does not fold in the simulation. But does that mean helix 5 is unimportant for folding? We checked the expression of two C-terminal truncated versions of λ*HG (a Tyr33His mutant with stability very similar to that of λ*YG). We truncated λ-repressor fragment after either amino acid S72 or S78. Both fragments expressed poorly as compared to the full-length construct. The secondary structure content of both fragments was reduced compared to the wild-type protein, as confirmed by circular dichroism. Neither of the truncated proteins showed a distinct thermal melting transition, as probed by circular dichroism spectroscopy or fluorescence spectroscopy (see SI). It appears that at least a portion of helix 5 is required for successful folding of λ-repressor fragment, even if the C-terminal helix is very flexible under conditions favoring the native state.

Discussion

Our simulation shows that refolding of a relatively large (80 residue) protein domain can be observed successfully after a P-jump in silico. The computed folding process consists of two very fast (∼1–2 μs duration) events, separated by a slower (∼16 μs) conformational search. In the experiment, a burst phase consistent with these fast times was observed; there was also a millisecond phase, likely due to the formation of non-native helix in the pressure-denatured state, and discussed previously.[19] We focus our discussion on the fast refolding induced by P-jump and the microsecond molecular phase which is observed in our simulation.

The two fast events both occur near the “speed limit” of folding proposed for λ-repressor based on T-jump experiments: it was previously observed that as λ repressor fragment is successively stabilized by mutation or lowering the temperature, a ∼1 μs fast phase appears and grows in amplitude.[20,31] This “molecular phase” was attributed to a fraction of the protein population poised at the barrier top, visible to an ensemble experiment only when the protein is stabilized so its folding barrier approaches ∼RT. A similar time scale has also been observed for slower folders by single-molecule experiments resolving the passage of a protein across the transition state.[21]

The molecular time scale corresponds to an average over the complex network of fast dynamics observed in hidden Markov models.[54−59] Our Cα-RMSD probe in Figure 2B should be sensitive to events before, during, or after passage through the bottleneck for folding, and indeed, the two fast events we see correspond to different regions of the free energy landscape.

The first fast event observed here corresponds to relaxation in the denatured basin, or downhill formation within 2 μs of a folding intermediate containing native secondary structure in helices 1, 2, and 4, with 1 and 4 properly aligned. Downhill folding was first invoked to rationalize the very fast appearance of folding intermediates during refolding of phosphoglycerate kinase and RNase H that occurred only under stabilizing conditions.[60,61] A truncated version of the λ repressor fragment with only helices 1 and 4 remaining has been shown to be stable and fold on a time scale of a few microseconds,[62] further supporting the alignment of helices 1 and 4 as an important first step during folding. In a simulation where these helices were incorrectly oriented, the protein remained in a trapped state for a full 60 μs of simulation.[25]

The second fast event can be identified with passage through the bottleneck or transition-state ensemble, as it is preceded by a 16 μs period of unproductive conformational search. This passage is not instantaneous, but takes 0.9 μs during which two helices rearrange and a loop forms. This event represents the speed limit of folding in the absence of unproductive conformational search, unless the three motions that comprise the event can be choreographed even more tightly upon further redesign of the λ repressor fragment sequence.

According to ref (20), the barrier height along a one-dimensional reaction coordinate (e.g., χ = 0–0.5 in the experiment or Cα-RMSD = 2–8 Å in the simulation) can be estimated from the molecular rate km and the activated rate of folding ka asAlthough a single trajectory lacks the statistics to determine these rates, we may approximate km ≈ (0.9 μs)−1 and ka ≈ (19 μs)−1, yielding a barrier of about 3RT. The molecular and activated time scales reported here are in good agreement with measurements of λ*YG by T-jumps at its melting temperature (km = (2 μs)−1 and (ka = 22.5 μs)−1).[30]

The step-by-step folding pathway reported here for λ*YG is in agreement with what has been proposed previously for wild-type λ-repressor from theoretical studies[63,64] and from implicit solvent replica-exchange MD simulation on another fast folding λ-repressor mutant.[53] However, the folding pathway observed here is different from what we reported on the λ*HG mutant.[25] For λ*HG, helices 1–3 assumed their native conformation first. Helix 4 of λ*HG formed individually in the early folding stage, but finding the correct orientation relative to other helices was a slow process that involved kinetic traps.[25] The difference between λ*YG and λ*HG could result from the point mutations altering the energy landscape, from different initial structures (pressure-denatured state vs extended state) altering the initial condition for refolding, or simply from a heterogeneous ensemble of folding pathways. It is clear from past work on λ repressor fragment that different mutations and different solvent conditions lead to different folding pathways,[50,65] but the exponential sensitivity of population to free energy (P1/P2 = exp[−ΔG12/RT]) makes it rather unlikely that the same sequence will fold by several equally Boltzmann-weighted pathways. In terms of a funneled rough energy landscape, it is easy to see how perturbations of the funnel could switch the most likely path, but the random free energy variations introduced by landscape roughness make it unlikely that many equally weighted paths exist. Indeed, relatively few cases of parallel pathway folding are known, such as lysozyme, staphylococcal nuclease, or certain repeat proteins.[66−68]

The same question is much more pressing for the nature of the folding bottleneck: How much does the timing and sequence of the three motions we observed during the 0.9 μs passage through the bottleneck change from trajectory to trajectory? How heterogeneous is the transit through the transition state? The answers to these questions are presently unknown at an atomistic level of resolution and will require multiple trajectories to provide proper sampling. It is in our opinion a priority for computational folding studies during the next several years.

For both λ*YG and λ*YA (the latter has G46A/G48A mutations instead of A37G/A49G), the fast P-jump experiment shows a <3 μs burst phase and a >1 ms slow phase. One of our long simulations (λ*YA) got trapped for the duration of the simulation,[19] whereas the other (λ*YG) folded rapidly. We thus propose that prompt folding is roughly as likely as kinetic trapping in both the P-jump experiments and simulations.

If so, this raises the same question raised by Lapidus and co-workers in their microfluidics experiment and by Pande from simulations:[28,69,70] Are there very slowly interconverting denatured states? For example, a kinetic scheme such aswith U and T initially populated, can explain how some proteins (initially in U) reach the native state rapidly whereas others (initially in T) reach it slowly. If U and T bracket the fluorescence lifetime of N, it is possible to observe microsecond and millisecond phases of opposite sign (Figure 3). It is well documented that pressure unfolding can be very slow due to the positive activation volumes,[7,71,72] and that pressure denaturation forms more compact denatured states than temperature denaturation. Compact states that fold more slowly than highly extended denatured states have been observed. For example, trpzip2 has a heterogeneous denatured population, whose subpopulation with blue-shifted tryptophan fluorescence (less solvent exposed tryptophan) folds more slowly to the native state than the subpopulation with red-shifted tryptophan fluorescence (more solvent exposed tryptophan).[73,74] Thus it remains unproven, but fully consistent with our data and simulations, that the pressure-denatured state of λ-repressor contains slowly interconverting compact traps.

Finally, we investigated what hinders the formation of helix 5 in our simulation. In the native state, helix 5 is stabilized by interacting with helix 4 through a small hydrophobic patch.[25] Formation of this patch is prevented in our simulation by several non-native salt bridges (Figure 4). Since helix 5 has weak helix propensity and needs to form in concert with its tertiary contact with helix 4, the non-native salt bridges that keep the C-terminus away from helix 4 hinder the formation in helix 5. One should not conclude that the non-native salt bridges hinder overall folding. They form within the first 2 μs after the P-jump and remain stable for most of the trajectory (Figure 4). They may play a crucial role in bringing the N- and C-terminus together, facilitating the formation of the main hydrophobic core that involves helices 1 to 4. Particularly the bridges between helix 1 and the C-terminus (Glu13-Arg85 and Arg17-Glu83), may accelerate the folding reaction in its early stage, but slow it down toward the end. In that regard it is noteworthy that λblue1, a two-helix bundle containing only helices 1 and 4, folds rapidly (the conformational search of helices 2 and 3 is eliminated),[62] but truncating helix 5 in the full-length protein (Results) produces a highly unstable protein with low expression level and missing secondary structure.

Figure 4

Non-native salt bridges in the refolding trajectory, monitored by the distance between oxygen bridges.

Outlook

The pressure jump simulation reveals key steps in moving the protein across the transition state during folding. However, limitations do exist. First of all, limited computational resources only allow the production of a single or few trajectories, preventing one from drawing a statistical conclusion, as seen in the case of the molecular rate, folding rate, or the structural heterogeneity of the transition-state ensemble. Second, a bias toward a certain kind of secondary structure has been observed in the current generation of force fields.[26,38,39] We employed the CHARMM22 force field with CMAP correction in our simulations.[35−37] Although it has been shown that this force field has a helical bias,[26,38] we used it successfully to fold a fast folding mutant of λ-repressor (λ*HG),[25] and it seems to explain qualitatively both trapping and prompt refolding after a P-jump (this work). The high temperature (325 K) used in the current study has a destabilizing effect on the helical structure, which can help in balancing the helical bias to some extent. However, the development of force fields that reproduce melting temperatures and hence denatured ensembles of proteins better remains an important goal. Third, a chemical denaturant (GuHCl) was used in the experiment to help unfold the protein at just 1.2 kbar. The simulations, on the other hand, utilized high temperature to accelerate the unfolding process. Such differences between experiment and simulation cannot yet be avoided altogether, but should be minimized as progress is made in both areas. Finally, all force fields are parametrized for ambient temperature and pressure. However, unfolding simulations are widely performed at very high temperature (and pressure in our case) to generate the denatured state for refolding simulation. Since the interatomic interaction strengths are much less sensitive to pressure changes than temperature changes in the ranges chosen for our simulations,[4] we expect that high-pressure simulation at less severe temperatures could yield useful initial states for refolding.