INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION NOTE ADDED IN PROOF REFERENCES

Advances in single-molecule manipulation techniques have recently made it possible to mechanically unfold single polypeptides and characterize their force-extension curves (Rief et al., 1997; Kellermayer et al., 1997; Oberhauser et al., 1998). Although studies to date have focused on molecules whose elastic properties are of physiological relevance, it is an exciting possibility that mechanical unfolding experiments on arbitrary proteins may yield information about forces that drive folding, the microscopic distribution of folding intermediates and substrates, and the kinetic barriers in folding pathways. Mechanical unfolding experiments offer a rare opportunity to follow a complex stochastic process on the level of single molecules; ultimately, they may provide the first microscopic test of statistical mechanical models of protein folding. However, the relationships between the observables (force and extension) and microscopic molecular events are not self-evident, nor is the relationship between mechanically induced unfolding pathways and the folding and unfolding pathways of proteins free in solution. These questions have begun to be addressed by theoretical lattice model studies (Socci et al., 1999; Klimov and Thirumalai, 1999) as well as by all-atom molecular dynamics (Lu et al., 1998; Lu and Schulten, 1999) and molecular mechanics (Rohs et al., 1999) simulations.

A recent comparison of mechanical and chemical denaturation of Ig domains (Carrion-Vasquez et al., 1999) yielded encouraging results for the relevance of atomic force microscopy (AFM) unfolding experiments to the folding pathways of untethered proteins: unfolding rates obtained by AFM were extrapolated back to zero force, and found to agree with the unfolding rates calculated by extrapolating bulk-solution results back to zero denaturant concentration. Moreover, the fractional extension of the force-induced transition state was found to agree with the fractional solvent exposure (based on m-value) of the denaturant-induced transition state, suggesting that a similar unfolding pathway is followed in the two pathways. It remains to be seen how general these correlations will be; Socci et al. (1999) have argued that, based on lattice models for folding, subjecting proteins to high forces should shift the transition state closer to the native state, and it may be that the titin domain studied by Carrion-Vasquez et al. (1999), whose transition-state m-value is unusually low, behaves in an exceptional fashion due to optimization for mechanical stability (Lu et al., 1998; Lu and Schulten, 1999). Whether or not mechanical unfolding experiments are found to mimic untethered unfolding, the projection of the high-dimensional protein folding process onto the arbitrary coordinates of force and extension presents a formidable interpretational challenge. It is therefore helpful, at these early stages, to correlate atomically detailed mechanical unfolding simulations with simulated traces of the experimental observables (Lu et al., 1998; Lu and Schulten, 1999; Rohs et al., 1999).

We have performed molecular dynamics (MD) to investigate the mechanical unfolding pathway of the C-terminal -hairpin fragment of protein G. This 16-residue peptide has many of the folding characteristics of larger proteins: it adopts a unique native conformation (Kobayashi et al., 1993), folds cooperatively (Muñoz et al., 1997), and contains both specific secondary structure and a cluster of aromatic sidechains that pack into a structure reminiscent of a hydrophobic core. Its unfolding and refolding trajectories have been extensively characterized by high-temperature molecular dynamics and transition-state analysis (Pande and Rokhsar, 1999), allowing a point of reference for comparison with mechanically induced events. Equilibrium-free energy surfaces have also been calculated by the multicanonical Monte Carlo method using an effective solvent approximation (Dinner et al., 1999). The small size of the peptide allows us to use pulling velocities more than an order of magnitude smaller than previous MD mechanical unfolding simulations of the immunoglobulin domain (Lu et al., 1998; Lu and Schulten, 1999), although they are still several orders of magnitude faster than the pulling velocities of AFM experiments (Rief et al., 1997). Compared to previous molecular dynamics studies (Lu et al., 1998; Lu and Schulten, 1999), we have used relatively soft cantilever stiffnesses (0.4 or 2.0 kcal/mol·Å = 0.28 or 1.4 N/m), approaching the order of those used experimentally in AFM studies (0.05 N/m; Guoliang Yang, personal communication). Although this leads to relatively large fluctuations in distance, it also contributes to smaller perturbing forces. We chose these parameters because we were primarily interested in qualitatively investigating the minimally perturbed unfolding pathway, rather than precisely mapping out a potential along the pulling coordinate.

An experiment in which measurements are made on the mechanical unfolding of a single -hairpin may eventually become feasible. However, our purpose is not to predict the results of such a technical tour-de-force. Rather, by comparing the thermal and mechanical unfolding pathways of the simulated peptide, and by correlating features of the force-extension curve with microscopic events, we hope to gain general insight applicable to the design and interpretation of mechanical protein unfolding experiments.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION NOTE ADDED IN PROOF REFERENCES

Molecular dynamics was carried out using NAMD2 (Nelson et al., 1996, Kale et al., 1999) with the CHARMM19 potential set (Brooks et al., 1993). The application of SHAKEH constraints (Ryckaert et al., 1979; Kale et al., 1999) to all hydrogens allowed the use of 2-fs timesteps. All simulations were carried out at 300K, with temperature rescaling performed every 10 timesteps. Periodic boundary conditions were used corresponding to a box of dimensions 51.5 × 34.75 × 34.75 Å. To generate the solvent, this box was filled with 1936 molecules of TIP3P (Jorgenson et al., 1983) water to a density of 0.93 g/cm³, and equilibrated via 100 ps of MD at 300K. The initial conformation of the C-terminal fragment of protein G (GEWTYDDATKTFTVE) was either the initial conformation used by (Pande and Rokhsar, 1999) (the result of a 100-ps equilibration, with water, starting from the first Protein Data Bank entry for 1GB1), or the PDB entry itself (see Fig. 2, A and C), with waters allowed to relax for 100 ps. The hairpin was added to the box in different orientations to facilitate pulling from terminal or core residues (Fig. 1). For each of these orientations, overlapping waters were removed using X-PLOR (Brünger, 1992), and the solvent/protein system was equilibrated at 300K for a further 100 ps of MD. During the equilibration of the system, the  carbon of F52 or E56 was held fixed, and that of Y45 or G41, respectively, was harmonically constrained to its starting position. The final conformations of these equilibrations were used to begin each unfolding simulation. Each simulation was initiated by assigning a Maxwell velocity distribution at 300K, followed by 5-10 ps of MD using the constraints described above. (In one simulation, labeled "*" in Fig. 4 C, the protein was being pulled during this brief interval, but the results were indistinguishable from other runs using the same spring constant and velocity.) For the remainder of the simulation, the harmonic restraint point was moved at a constant velocity in the direction of the vector connecting the two constrained atoms, while the temperature was maintained at 300K by velocity rescaling. The pulling velocity was chosen to be 2.5 m/s, 2.26 m/s, or 0.9 m/s. Spring constants of 0.4 kcal/mol·Å² or 2.0 kcal/mol·Å² were used.

View larger version (63K): [in this window] [in a new window]   FIGURE 1   Attachment points used for mechanical unfolding experiments. (A) Pulling from the core residues: the  carbon of Phe-56 is held fixed while that of Tyr-45 is harmonically constrained to a moving restraint point. (B) Pulling from the termini: the C-terminal -carbon is held fixed while the N-terminal -carbon is attached to a moving harmonic restraint.

Two measures of unfolding were calculated from the molecular structures along the unfolding pathway. The radius of gyration of the core was defined as the root mean square distance of the sidechain atoms of the core tyrosine, phenylalanine, and tryptophan residues from their collective center of mass. The number of backbone hydrogen bonds was calculated by counting the number of donor-acceptor pairs in which an amide hydrogen was within 2.5 Å of a carbonyl oxygen. Unlike the previous study (Pande and Rokhsar, 1999), only canonical hydrogen-bonding partners (Muñoz et al., 1997) were considered.

	RESULTS

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION NOTE ADDED IN PROOF REFERENCES

We performed nine simulations of mechanical unfolding of the  hairpin (Table 1), varying the pulling speed, spring constant, and attachment points (Fig. 1). In both of the slow pulling speed simulations (soft/core/slow, soft/ends/slow), which were expected to produce the smallest perturbations in the thermally activated behavior of the peptide, unfolding was found to follow a stepwise pathway (Fig. 2) that passes through an intermediate state that lacks any backbone hydrogen bonds but maintains a dense hydrophobic core. This agrees well with the unfolding pathway proposed earlier based on high-temperature simulations (Pande and Rokhsar, 1999), and suggests that, under mild mechanical stresses, the force-induced unfolding trajectory can reflect the thermally accessible energy landscape.

View this table: [in this window] [in a new window]   TABLE 1   -Hairpin mechanical unfolding simulations

View larger version (49K): [in this window] [in a new window]   FIGURE 2   Unfolding pathway of the -hairpin at low pulling speed. Force-time curves for (A) soft/core/slow and (B) soft/ends/slow simulations are overlaid with plots of the number of hydrogen bonds and the radius of gyration of the core sidechains. All traces have been smoothed by boxcar averaging using a 21-ps window, with the raw data sampled once per picosecond (this smoothing accounts for apparent discrepancies with Table 1). Hydrogen bond rupture is completed (blue arrows) before dissociation of the core (green arrows, identifiable as the beginning of the increase in core radius of gyration). Force peaks indicate the resistance of hydrogen bonds (black arrows), the hydrophobic core (red arrows), or the central kink (purple arrow). (C, D) Molecular dynamics snapshots demonstrating the unfolding pathway under tension applied at (C) the core residues or (D) the termini. An intermediate state can be seen in each case [(C) 900 ps; (D) 2041 ps], which largely retains the hydrophobic core but does not contain hydrogen bonds.

The overall unfolding pathway is shared between these two simulations in spite of radically different pulling geometry (Fig. 1), but there are significant differences in the details of the trajectories. In each case, the metastable hydrophobic core involves only three of the four hydrophobic residues. In the soft/core/slow simulation (Fig. 2, A and C), the tyrosine residue, which is under direct tension, flips out early. This causes only a minor rise in the radius of gyration of the core, due to the compaction of the remaining aromatic groups. When pulling from the ends (Fig. 2, B and D), the valine is removed early from the core; the rearrangement of the three aromatic groups into a tight nucleus actually yields an initial decrease in the core radius of gyration (R_G), due to the exclusion of the valine atoms from the core R_G calculation, following Pande et al. (1999). The trajectories also differ in the order and duration of hydrogen bond breakage: in the soft/core/slow simulation (Fig. 2 A), rupture progresses from the core outward (data not shown) and occurs within the first 700 ps, whereas the hydrogen bonds in the soft/ends/slow simulation unzip one by one from the ends upward, over the course of 2000 ps.

In both simulations, one or more peaks in the force (Fig. 2, arrows) mark barriers between the three major folding species (folded, metastable core, and unfolded). A histogram of the time the peptide spends at each extension in the soft/core/slow simulation (Fig. 3) lends support to the interpretation of the force peaks as transition states. In the soft/core/slow simulation, a final barrier can be seen (Fig. 2 A, purple arrow), which involves a transition to an unfolded state that was not observed in high-temperature simulations (Pande and Rokhsar, 1999); in this state (Fig. 2 C, final frame), the central bend of the hairpin is straightened. In the period preceding its unbending, the turn appears to be stabilized by an interaction between the sidechain carboxyl of Asp-46 and the amide hydrogens of Lys-50 and Thr-49 (data not shown). The unbent hyper-unfolded state is unlikely to be significantly sampled in the absence of external forces.

View larger version (40K): [in this window] [in a new window]   FIGURE 3   Force-extension curve of the soft/core/slow simulation. Extension is defined as the distance between the two constrained alpa carbons. The underlying histogram displays undervisited barrier regions between states. In the lower panel, the smoothed number of hydrogen bonds and radius of gyration of the core are shown as a function of extension.

Under both pulling geometries, the highest force peak (Fig. 2, red arrows) occurs in the absence of any hydrogen bonds, and appears to represent the resistance of the hydrophobic core, as can be seen by the rise in the radius of gyration of the core following the force peak. In this respect, the hairpin seems to differ from the Class I -sandwich domains whose mechanical unfolding was simulated by Lu and Schulten (1999): they determined that cooperative hydrogen bonds were mostly or entirely responsible for the large resistance of Class I domains to mechanical unfolding.

In the faster pulling speed simulations (Fig. 4), hydrogen bond rupture and dissociation of the hydrophobic core are less well-separated events. This is particularly true of the core-pulling trajectories, in which an initial steep buildup of force is followed by rapid and simultaneous loss of hydrogen bonds and the compact hydrophobic cluster. The high forces, short peak distances, and cooperative unfolding induced by pulling quickly on a -hairpin from the core residues are qualitatively reminiscent of simulation results obtained by Lu et al. (1998) and Lu and Schulten (1999) for Class Ia (titin-like) -sandwich domains, and of interpretations of pulling data obtained experimentally on molecules of the same fold (Rief et al., 1997; Carrion-Vasquez et al., 1999). In simulations of the hairpin using a fivefold higher spring constant (Fig. 4 C), considerably larger forces are generated, and the trajectories are relatively monotonic, perhaps reflecting high-force swamping of the thermal behavior of the peptide.

View larger version (28K): [in this window] [in a new window]   FIGURE 4   Unfolding pathways and force-extension curves of simulations at higher pulling speed (2.5 m/s). (A) In soft/ends/fast simulations, hydrogen bond rupture can overlap with core dissociation. (B) In soft/core/fast simulations, core dissociation occurs concurrently with cooperative hydrogen bond rupture. (C) Force-extension curves of soft/core/fast simulations show variable peak forces but generally repeatable features, including a steep initial rise in force (the lower dotted gray curve, which differs qualitatively from the others, represents a run in which partial hydrogen bond breakage allowed some relaxation before the main force peak). In the hard/core/fast simulation (upper dotted curve), using a stiffer spring constant (2.0 kcal/mol·Å2 = 1.4 N/m), higher forces are generated and the peaks are less sharp.

	DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION NOTE ADDED IN PROOF REFERENCES

We have simulated the mechanical unfolding of a model protein. The peak forces (2.5-5.2 kcal/mol·Å = 170-360 pN) generated in our simulations using a soft cantilever were of the same order of magnitude as those experimentally measured in AFM unfolding of larger proteins (Rief et al., 1997; Oberhauser et al., 1998; Carrion-Vasquez et al., 1999). This is consistent with kinetic models for mechanical unfolding (Carrion-Vasquez et al., 1999), given that our choice of pulling velocities accelerates hairpin unfolding by a factor of ~10³ (Muñoz et al., 1997), similar to the factor by which unfolding of Ig domains was accelerated in the slowest AFM experiments (Carrion-Vasquez et al., 1999). We found that, under the mildest mechanical stresses, the C-terminal -hairpin fragment of protein G unfolded via a stepwise pathway in which hydrogen bond rupture precedes core dissociation. In addition to the states (unfolded, collapsed core, and folded) observed in high-temperature unfolding simulations (Pande and Rokhsar, 1999), a final barrier was observed to the formation of a distinct hyper-unfolded state in which the central kink of the hairpin is straightened. This state is likely to be very poorly populated in the absence of external forces, and is an example of the ability of mechanical unfolding experiments to probe regions of the folding free energy landscape that are difficult to access by other means.

These simulations add to growing theoretical support (Pande and Rokhsar, 1999; Dinner et al., 1999) for a model of hairpin formation in which an early step is the formation of a hydrophobic cluster, independent of hydrogen bond formation. This picture was noted first by Muñoz et al. (1997), who favored an alternative scenario in which hairpin folding is initiated by specific interactions at the turn. They proposed a simple statistical mechanical zippering model that was in quantitative agreement with their experimental data. Experiments to distinguish the two models have been suggested by Muñoz et al. (1997), Pande and Rokhsar (1999), and Dinner et al. (1999), but the distinction remains unresolved.

Although these simulations provide some optimistic results for the relevance of mechanical unfolding experiments to the untethered folding process, some caveats can also be extracted. First, the use of stiff cantilevers and/or high pulling speeds can hamper the detection of folding intermediates by forcing the protein across high barriers, as seen in Fig. 4: pulling quickly from the core residues led to simultaneous core dissociation and hydrogen bond rupture, bypassing the intermediate seen in simulations of thermal unfolding and generally yielding a single enlarged force peak.

Second, even under milder conditions, the details of the pathway are dependent on the geometry of the experiment: the order of hydrogen bond breakage and the conformational nature of the intermediate collapsed core differed between simulations of pulling from core or end residues. Internal coordinate molecular mechanics simulations have already pointed out that striking differences in rupture forces can exist between different pulling geometries for the same structure (Rohs et al., 1999). Therefore, mechanical unfolding experiments should be repeated using different geometries, to extract the features that are in common between them, before making conclusions about the untethered pathway.

It is unlikely that atomic force microscopes will attain sufficient temporal and spatial resolution in the near future to measure transitions between intermediates in -hairpin formation via the measurement of force-extension curves such as those simulated herein, but such measurements on larger single-domain proteins may soon be realized. By pulling sufficiently slowly, and using alternative methods of attachment to the polymerization strategies currently in vogue, it should be possible to resolve features in the force-extension curves, analogous to those presented here, which correspond to transitions between intermediate states. By analyzing the effects of mutations on these features, or by the simultaneous use of fluorescent structural probes, the conformational details of these intermediates may be elucidated.

	NOTE ADDED IN PROOF

TOP ABSTRACT INTRODUCTION METHODS RESULTS DISCUSSION NOTE ADDED IN PROOF REFERENCES

In recent studies of mechanical unfolding in tandem repeats of titin modules, Marszalek et al. (1999) have observed a previously undetected shoulder in the force-extension curve. Their analysis, including molecular dynamics simulations and mutagenesis, suggests that this feature represents a transition to a specific mechanical unfolding intermediate.