J. Phys. Chem. B 2008, 112, 13079–13090
13079
Peptide Folding Using Multiscale Coarse-Grained Models Ian F. Thorpe, Jian Zhou,† and Gregory A. Voth* Center for Biophysical Modeling and Simulation, Department of Chemistry, UniVersity of Utah, Salt Lake City, Utah 84112-0850 ReceiVed: February 22, 2008; ReVised Manuscript ReceiVed: July 4, 2008
The multiscale coarse-graining (MS-CG) method has been previously used to describe the equilibrium properties of peptides. The present study reveals that MS-CG models of R-helical polyalanine and the β-hairpin V5PGV5 possess the capacity to efficiently refold in simulations initiated from unfolded configurations. The MS-CG peptides exhibit free energy landscapes that are funneled toward folded configurations and two-state folding behavior, consistent with the known characteristics of small, rapidly folding peptides. Moreover, the models demonstrate enhanced sampling capabilities when compared to systems with full atomic detail. The significance of these observations with respect to the theoretical basis of the MS-CG approach is discussed. The MS-CG peptides were used to reconstruct atomically detailed configurations in order to evaluate the extent to which MS-CG ensembles embody all-atom peptide free energy landscapes. Ensembles obtained from these reconstructed configurations display good agreement with the all-atom simulation data used to generate the MS-CG models and also corroborate the presence of features observed in the MS-CG peptide free energy landscapes. These findings suggest that MS-CG models may be of significant utility in the study of peptide folding. 1. Introduction Coarse grained (CG) modeling in the field of molecular simulation has recently undergone a significant expansion.1-3 CG models typically omit the explicit description of certain components of a molecular system, reducing the dimensionality of the underlying phase space. What remains is an essential subspace that allows one to focus attention on the molecular properties that are of greatest interest; other components need only be described implicitly. Diverse approaches to develop and characterize CG models have been explored. Representative studies include those for polymers,4 membranes and membrane proteins,5-13 protein folding1,14-18 and the dynamics of protein assemblies.19-22 These efforts generally involve the use of “effective” interactions to incorporate the effects of those parts of the system that are described implicitly. The basic concepts involved have been appreciated for some time: the use of united atoms for molecular dynamics is one well-known example of a CG approach. One could also consider implicit solvent models as CG models in which explicit solvent has been coarse-grained into effective solute-solute interactions.23,24 The impetus behind such efforts has primarily been to save time and computational resources by focusing these resources on the components of a system that are of greatest interest. Coarse-graining is therefore one approach by which the length and time regimes accessible to molecular simulations can be expanded. A particularly promising new class of CG methodology is inherently multiscale in character,2 and seeks to cohesively and consistently incorporate information from different length and time regimes. Recently our research group introduced a novel approach to obtaining effective interactions for CG molecular systems that is referred to as the multiscale coarse-graining (MS* Corresponding author. E-mail:
[email protected]. Telephone: (801) 581-5419. Fax: (801) 581-4353. † Current address: School of Chemical Engineering, South China University of Technology, GuangZhou, GuangDong Province, 510640, P. R. China
CG) method.3,9,10 The motivation for this approach is to employ information obtained from molecular simulations carried out with full atomic detail in order to systematically derive CG models. An important tool used to obtain these effective interactions is the so-called force matching methodology.9 In force matching, forces computed for CG sites using detailed all-atom simulations are used to derive effective CG forces. A noteworthy benefit of this approach is that it offers a systematic way to obtain effective interactions without making empirical assumptions about their underlying functional form. In principle this offers a rigorous and systematic multiscale route toward CG interactions. In previous studies MS-CG representations of R-helical polyalanine and the β-hairpin V5PGV5 were generated.25 These two peptides represent basic building blocks of protein structure: it is anticipated that information obtained for these systems will also find application in the studies of larger peptides and proteins. The models were shown to capably reproduce the equilibrium structural properties of the two peptides. In contrast, the present studies investigate the folding properties of these peptide models. Understanding the process of protein folding is key to comprehending the molecular origin of protein structure and function; this field remains a highly active research area.26 As reported herein, the MS-CG peptides have been found to possess the capacity to efficiently refold unfolded configurations. The peptide models also exhibit two-state folding behavior, consistent with the known characteristics of small, rapidly folding peptides. This finding is unanticipated given that the MS-CG models were developed using simulations containing only folded peptides. Moreover, the MS-CG peptides possess free energy landscapes that are funneled toward the folded state and demonstrate enhanced sampling capabilities when compared to systems with full atomic detail. The theoretical basis of these observations with respect to the unique features of the MS-CG approach is discussed below. Benefits of using the MSCG approach to study folding will also be outlined. The MS-
10.1021/jp8015968 CCC: $40.75 2008 American Chemical Society Published on Web 09/20/2008
13080 J. Phys. Chem. B, Vol. 112, No. 41, 2008
Figure 1. Coarse grained models for small peptides representing distinct secondary structural motifs. All-atom (left-most images in each panel) and CG models of (a) the A15 helix and (b) the V5PGV5 hairpin. For all-atom representations the red, blue and cyan beads represent the oxygen, nitrogen and carbon atoms in the CO, NH and CH backbone groups respectively; white beads indicate the carbon atom nearest to the side-chain center of mass. For CG models red, blue, and cyan beads represent the OBB, NBB and CBB backbone groups respectively while side-chain groups are displayed as white beads.
CG peptide models were employed to rebuild all-atom detail into the models. This allowed us to determine the extent to which the MS-CG ensembles correspond to all-atom peptide free energy landscapes. The reconstructed configurations allowed us to generate ensembles that agree well with the original allatom data from which the MS-CG models were generated and lend support to the presence of features identified in the MSCG free energy landscapes. These observations suggest that MSCG models may be of significant utility in the study of peptide folding, especially given the enhanced sampling capabilities of these models. 2. Methods 2.1. Coarse-Graining Scheme. Explicit TIP3P27 water is represented by a single CG site (CGW) located at the geometrical center of each molecule. Each amino acid residue is represented by four CG sites; three for the peptide backbone and one for each side-chain. Similar “intermediate resolution” peptide coarse-graining has been previously employed by other authors.28 Separate CG sites are used to represent the CO, NH, and CH groups of the peptide backbone and are referred to as OBB, NBB, and CBB respectively. The OBB site was placed at the carbonyl oxygen position, the NBB site was placed at the amino hydrogen position and CBB was placed at the CR position (see Supporting Information). The backbone N atom of proline was considered a separate group while the CR of glycine as well as the two hydrogen atoms bonded to it were represented by a single site located at the CR position. A single CG site was used to represent the side chain groups of each residue; these were placed at the center of mass of the constituent atoms (see Figure 1). Different types of CG side chain sites were employed to represent the side chain groups of alanine, valine and proline. 2.2. Coarse-Grained Force Field. The CG effective force field interactions are divided into two categories: “bonded” and “non-bonded” energy terms. Bonded terms involve sites separated by 1, 2, or 3 CG bonds. CG bonds are defined based on chemical connectivity, with bonds connecting different CG sites
Thorpe et al. that contain chemically bonded atoms. Bonded energy terms employ standard functional forms used in conventional molecular mechanics force fields. These include harmonic terms for intersite separations or angles and cosine series for dihedral angles. Energy terms Ux defining fluctuations along a given coordinate x were obtained by fitting the probability distribution along this coordinate, P(x), occurring during the original allatom simulations via: P(x) ) Cx exp[-βUx], where Cx is an arbitrary constant. Each Ux is defined by a standard force field energy expression using parameters (such as the average bond length or angle) fit directly from the all-atom trajectories. Nonbonded interactions were evaluated using the MS-CG force matching methodology developed in our research group. The core of this approach lies in minimizing the difference between computed atomistic forces projected unto CG sites and predicted CG forces. The resulting effective forces are radial, pairwise decomposable and depend only on intersite separations. Detailed descriptions of these methods have been presented in previous studies.9,10,25 A list of the parameters defining the CG energy functions in each of the peptide models has also been presented previously.25 Two different sets of simulations were employed to generate the CG force fields for A15 (polyalanine) and V5PGV5. The simulation system used to generate each CG peptide model consisted of five copies of each folded peptide molecule, all immersed in a single box of water. The choice to use five peptides rather than a single folded peptide to generate the effective interactions was made in order to facilitate adequate sampling of nonbonded interactions. Due to the constraints imposed by the peptides’ well-defined folded structure, intersite separations exist only within very narrow ranges for the folded state. However, in order to obtain a model that can be used for folding, nonbonded interactions that occur throughout a much wider distance range are necessary. Including five peptides in the simulation system allows nonbonded interactions to be well sampled at a wide range of intersite separations while maintaining probability distributions along internal degrees of freedom that are appropriate for the peptides’ folded state.25 The assumption made in this treatment is that nonbonded interactions for CG groups within a peptide chain are the same as those for the same CG groups in different peptide chains. All CG sites of a given type were treated equivalently and interactions averaged over sites at different locations in the simulations to construct the final force field. It is possible that CG sites of a given type may experience slightly different environments at various locations in the system. This procedure averages over such differences to ensure that the set of interactions obtained is the best fit for diverse environments. The goal of this scheme is to generate effective interactions that are as transferable as possible. A detailed description of the development of these models has been presented previously.25 2.3. Simulation Details. All-Atom Simulations. After generating the peptide models as described above, a number of additional all-atom simulations were carried out in order to validate the models. Configurations of the peptides A15 (polyalanine) and V5PGV5 were generated in full atomic detail using the CHARMM molecular simulation program.29 For folded simulations the initial configuration of A15 was set to be an idealized R-helix while V5PGV5 was generated as a β-hairpin with the proline and glycine residues forming the midpoint of the β-turn. Note that the “D” isomer of proline was employed. This peptide is a variant of the hairpin originally studied by Sung, who demonstrated that this molecule intrinsically favors
Coarse-Grained Peptide Folding
J. Phys. Chem. B, Vol. 112, No. 41, 2008 13081
Figure 2. Free energy landscape in the vicinity of the folded basin for each of the two peptides. Representative snapshots of the folded state are inserted into the plots. The color scale employed is shown below the free energy contours. Contours are drawn at intervals of kBT/2 (kB ) Boltzmann’s constant, T ) temperature). Data from (left) all-atom and (center) CG trajectories. Both sets of simulations were initiated from folded configurations. Apart from the slightly smaller Rg for CG V5PGV5, the folded basin is well preserved by the CG models. (Right) Regions of CG refolding landscapes proximal to the folded basin, magnified to show detail. One can observe that the CG folded basin is well-reproduced when starting from either folded or unfolded configurations.
the β-hairpin configuration.30 Either the DLPOLY31 or CHARMM molecular dynamics (MD) engines were employed in concert with the CHARMM 27 all-atom force field to simulate the peptides at 310 K in cuboid unit cells filled with TIP3P water. Simulations were performed in the constant NVT ensemble using a Nose´-Hoover thermostat with a relaxation time of 0.5 ps. Bonds containing hydrogen were held rigid using the SHAKE32 method with a tolerance of 10-6. Periodic boundary conditions were applied and the particle mesh Ewald method used for electrostatic interactions. A 10 Å cutoff for nonbonded interactions was employed. A 2 fs time step was used and configurations written every 2 ps. Trajectories of 4 ns were propagated for each folded peptide. In the case of refolding trajectories initial configurations were generated by carrying out high temperature (1000 K) unfolding simulations of each of the peptides in vacuum. Nine conformations for the helix and ten conformations for the hairpin were randomly selected to be solvated in TIP3P solvent and used to initiate refolding trajectories at a temperature of 310 K. One additional configuration was also employed in which each peptide was generated in a fully extended configuration and immersed in a large solvent box. Simulations were continued until peptide free energy landscapes computed from the trajectory data (see section 2.4) showed no further change. Simulation lengths for all trajectories are provided as Supporting Information. CG Simulations. Simulations of the CG peptides were started from CG versions of the same initial configurations as the allatom MD simulations. Except where noted below, simulation conditions were identical to those present for the all-atom
simulations. As a consequence of the coarse-graining procedure, all nonbonded (electrostatics and van der Waals) interactions were subsumed into a short-range effective force field. These effective interactions were precomputed and tabulated for use in the CG simulations. Due to the short-range nature of these interactions it was possible to use a longer cutoff than the allatom MD simulations (12 Å). 2.4. Data Analysis. Free Energy Landscapes. Relative free energy plots were computed by inverting the probability distributions obtained from molecular dynamics simulations using standard techniques, with the free energy F given by F ) -kBT[ln P(x,y)], where kB is Boltzmann’s constant, T is temperature, and P(x,y) is the normalized probability distribution occurring along variables x and y. To facilitate comparisons, the landscapes were projected unto the radius of gyration (Rg) and the root-mean-squared deviation (rmsd) from the folded configuration (Figures 2 and 3). The reference state for the rmsd calculations was a representative snapshot from equilibrated simulations of the folded state. rmsd was computed using the peptide backbone sites only. These two variables provide a concise description of the overall shape of the peptide molecules as well as the degree of similarity to the folded state. Note that the panels labeled “All-atom” in each of the plots denote the results from all-atom trajectories that were coarse-grained, allowing the requisite quantities to be computed for the corresponding CG sites from the all-atom data. Backbone Dihedral Angles. Due to the CG scheme employed, it is possible to compute angles analogous to the Φ and Ψ dihedrals conventionally used to characterize peptide
13082 J. Phys. Chem. B, Vol. 112, No. 41, 2008
Thorpe et al.
Figure 3. Peptide refolding landscapes. Representative peptide snapshots (not to scale) are inserted into the plots. (Left) All-atom simulations initiated from random initial configurations, (Center) MS-CG simulations started from coarse-grained versions of the same initial configuration employed for the left panel. (Right) All-atom simulations initiated from MS-CG configurations inset on the right-most side of the center panels. These MS-CG configurations were used to reconstruct atomic detail as described in the text before initiating the all-atom simulations. Note that the all-atom refolding landscapes do not span the same range of rmsd as the CG panels because none of the all-atom refolding trajectories revisit the folded state. In contrast, each of the CG trajectories are able to recover the folded basin (dark blue area in lower left corner of CG panels, also see Figure 2). The MS-CG models are able to span configuration space (including the folded basin) that the all-atom systems do not sample: these regions are circled in red in the CG refolding landscapes. To aid in visualization, constant offsets of 2.5 and 5.5 kBT have been added to the A15 and V5PGV5 all-atom refolding landscapes, respectively, in order to make the color scale more consistent with that present in the MS-CG landscapes.
configuration space by evaluating dihedrals for CG groups that contain the relevant atoms. For example, the OBB-NBB-CBBOBB dihedral corresponds to Φ while the NBB-CBB-OBBNBB angle corresponds to Ψ. Because the positions of the relevant all-atom and CG groups are quite similar, distributions of the CG backbone dihedrals are comparable to those obtained for conventional Φ and Ψ dihedrals. Thus, it is observed that the CG R-helical peptide exhibits dihedrals corresponding to the R-R region of dihedral map while the β-hairpin displays dihedrals appropriate for the β region. Because of this observation the CG backbone dihedrals are also referred to using Φ and Ψ terminology. However, these are denoted “pseudo” to highlight that they are computed for the CG representation (see Figure 4). As noted above, the panels labeled “All-atom” in Figure 4 are results from all-atom trajectories that were coarsegrained, allowing the pseudo dihedrals for the appropriate CG groups to be computed from the all-atom data. Thus, the differences in the left and right panels represent intrinsic differences between the all-atom and MS-CG simulations rather than differences in the dihedral angle definitions. Hydrogen Bonding. The MS-CG models no longer contain groups that can engage in conventional hydrogen bonds (HBs), as these structures at the all-atom level have been coarse grained away. However, the pairwise MS-CG interaction between NBB and OBB peptide groups implicitly incorporates hydrogen bonding effects in a spherically averaged manner, as do the interactions between these groups and CGW. In this work, two criteria are employed for peptide backbone groups to assess hydrogen bonding that conforms to that found in the folded
configuration (i.e., “native” HBs). The first criterion is the distance between NBB and OBB peptide groups. The distribution of interparticle distances that occur in CG versions of the original all-atom simulations revealed that distances less than or equal to 3 Å are consistent with native HBs. However, this simple distance criterion on its own is insufficient to uniquely identify native HB configurations. Consequently, dihedral angles involving nearby peptide backbone groups were employed in tandem with the distance requirement (see Figure 5). This ensures that the local configuration of the peptide backbone is consistent with a native HB. Allowed dihedral intervals for individual HBs were determined by computing the average value assumed by each dihedral during all-atom simulations of the folded peptides. A given dihedral angle was defined to be consistent with a HB configuration if it differed by no more than one standard deviation from the mean value computed from the all-atom trajectories. Hydrogen bonding probability distributions were employed to generate free energy landscapes in a manner similar to that described above (see Figure 6). 2.5. Revisiting the All-Atom Ensemble. All-Atom Reconstruction. Using the MS-CG peptide models it is possible to reconstruct the missing details of the all-atom systems, allowing CG configurations to be used to explore corresponding regions of the all-atom free energy landscape. The approach taken will be described briefly below. For a given snapshot from the MSCG peptide trajectories CG solvent molecules were first removed. Due to the method of coarse-graining employed (see section 2.1), the MS-CG models immediately provide the locations of the carbonyl oxygen, amino hydrogen and CR atoms
Coarse-Grained Peptide Folding
J. Phys. Chem. B, Vol. 112, No. 41, 2008 13083
Figure 5. Schematic depicting criteria used to identify hydrogen bonds (HBs) consistent with the folded state (i.e., “native” HBs). The distance between NBB and OBB groups (thick dashed lines) in concert with dihedral angles involving nearby peptide backbone groups (thin dotted lines) were employed. Figure 4. Peptide free energy landscape projected onto backbone dihedral space. The color scheme employed is the same as that used for Figures 2 and 3. Due to the CG scheme employed (see Methods), it is possible to compute angles analogous to conventional Φ and Ψ dihedrals for the CG peptides: these are denoted “pseudo” in the plot. The panels labeled “All-atom” on the left denote the results of allatom trajectories that have been coarse-grained; the pseudo dihedrals were then computed for the CG configurations obtained from these trajectories. Thus, discrepancies between the left and right-hand panels reflect intrinsic differences in the properties of the MS-CG and allatom peptides rather than differences that result from employing the CG representation. Trajectories for polyalanine are displayed in the upper panels while those for the hairpin are shown in lower panels. Note that there is considerable agreement between the all-atom and CG dihedral distributions for folded peptides, indicating the degree to which the native structure is faithfully represented by the models. Also note that the all-atom refolding simulations explore considerably more configuration space than the corresponding CG models.
for the peptide backbone. The positions of the remaining atoms in the peptide backbone were estimated using simple geometric relationships involving the locations of known atoms. These estimates were based on standard features of peptide bonding geometry. Estimates for the locations of a single side-chain site for each residue were obtained by assigning the atom closest to the side-chain center of mass to the position of the MS-CG peptide side-chain site. The CHARMM module INTCOR was then employed to generate positions for the remaining atoms by using the known and estimated atomic positions to produce internal coordinate tables based on CHARMM 27 topology files. The generated structure was subjected to 1000 steps of steepest descent minimization to correct for coordinates with poor initial estimates. During this initial minimization, positions of the carbonyl oxygen, amino hydrogen and CR atoms for the peptide backbone were kept fixed since the locations of these sites are known with certainty. All constraints were then removed and another 500 steps of steepest descent minimization carried out. When tested on previously minimized all-atom peptide configurations that had been coarse-grained, this reconstruction
procedure lead to a deviation in backbone rmsd of approximately 1 × 10-3 Å. When the contribution of side-chains was included the rmsd was found to be less than 0.9 Å for A15 and within 2 Å for V5PGV5. Further details are provided as Supporting Information. Simulations with Reconstructed Configurations. After minimization the reconstructed peptide was placed in a new allatom MD simulation cell of the same size and geometry as the original MS-CG system. The simulation cell was filled with a number of atomically detailed solvent molecules equal to the number of MS-CG solvent molecules in the original MS-CG trajectory. The protein was held fixed and the solvent allowed to equilibrate around the structure during 50 ps of MD. This configuration was then used to initiate additional all-atom trajectories with the same simulation conditions as the all-atom simulations from which the original MS-CG models were obtained. 3. Results and Discussion 3.1. MS-CG Effective Interactions. Much insight is gained into the unique features of multiscale methodology by invoking principles of statistical mechanics.3 A general theoretical and algorithmic framework for this approach has been described in recent publications.33,34 In the discussions that follow, CG variables are displayed in upper-case and vector quantities are shown in bold. For an atomic system composed of n particles with coordinates rn, the process of coarse-graining the system into a set of N variables with coordinates RN involves determining the probability p(RN) of observing all-atom structures that map unto a given CG structure, where: N
p(RN) )
∫ drnp(rn)∏ δ[MR (rn) - RI] I)1
I
(1)
In this expression p(rn) ) e-βu(rn)/zn is the probability distribution function for the all-atom system, u(rn) is the original potential energy function, zn is the canonical configuration
13084 J. Phys. Chem. B, Vol. 112, No. 41, 2008
Thorpe et al. TABLE 1: Hydrogen Bonding Distributions from Peptide Simulationsc A15
V5PGV5
all-atom folded CG folded all-atom refolding CG refolding all-atom folded CG folded all-atom refolding CG refolding
〈HB〉a
〈contig〉b
0.76 0.90 0.00 0.51 0.29 0.45 0.00 0.38
5.69 7.93 1.00 4.75 1.36 1.59 1.00 1.49
a Average fraction of native hydrogen bonds b Average number of residues in longest contiguous hydrogen bonding segment c Note that the all-atom refolding trajectories never reach the folded state, which leads to an average native hydrogen bonding content of zero; no contiguous hydrogen bonding segments longer than two consecutive residues are observed for these trajectories.
integral for the all-atom system, and β ) 1/(kBT). MRI(rn)is a mapping operator that generates the coordinates for CG particle I that is appropriate for a given all-atom configuration: MRI(rn) ) ∑in ) 1 cIiri, where the coefficients cIi depend on the nature of the CG scheme employed. The coarse-graining procedure can then be framed as the search for a CG probability distribution function P(RN) that is consistent with the all-atom distribution: P(RN) ) e-βUCG(RN)/ZN ) p(RN), where ZN is the canonical configuration integral for the CG model and the original potential energy function u(rn) is replaced by an (effective) CG potential UCG(RN). It is thus apparent that UCG(RN) is given by:
UCG(RN) ) -kBT ln[p(RN)] + C
Figure 6. Peptide free energy landscapes projected unto hydrogen bonding coordinates. The color scale is identical to that used in Figures 2-4. The y axis denotes average fraction of hydrogen bonds (HBs) in each peptide while the x axis displays the average number of residues in the longest contiguous HB segment. These plots illustrate the ease of HB formation in the different simulation systems. Note that MS-CG models have higher propensity to form HBs than all-atom peptides. The contours indicate that the average number of hydrogen bonds increases in concert with the average length of HB segments and implies that HB formation is cooperative. Otherwise, increased HB content would tend to generate a large number of (isolated) short hydrogen bonding fragments. Similar plots are not shown for all-atom refolding trajectories because the low number of native HBs observed during those simulations leads to relatively featureless free energy landscapes (see also Table 1).
(2)
In the expression above, C is a constant that depends on the CG configuration integral. When formulated in this way it is apparent that the MS-CG effective potential is in principle a many-body free energy function that depends only on the CG variables. In other words, it is a many-body potential of mean force (PMF) that differs from the conventional two-body PMF. This was clearly demonstrated in our previous MS-CG peptide studies.25 An additional (but not obligatory) assumption made in the MS-CG approach is that UCG(RN) can be represented by a sum of “bonded” and pairwise “non-bonded” energy terms. We found it convenient to compute the energy terms for bonds, angles and dihedrals directly from the probability distributions obtained from all-atom simulations (section 2.2). This choice was made because the distributions of these variables in CG representations of the all-atom trajectories were found to be well described by standard functional forms commonly used in molecular simulations.25 In contrast, there is no easy way to predict a priori what form the effective nonbonded interactions should be. However, by using the force matching procedure it is not necessary to assume a given functional form for these interactions. Subject only to the simple assumptions that the effective interactions are pairwise additive and spherically symmetric, the force matching approach allows these interactions to take whatever functional form is dictated by the properties of the system. The MS-CG force matching procedure used to evaluate nonbonded interactions essentially computes NB the gradients of the nonbonded components UCG (RN) with NB respect to the positions of each CG particle I: -dUCG(RN)/dRI ) FINB, where the term on the right is the CG nonbonded force evaluated for the Ith particle in the system. As seen from eqs 1 and 2, obtaining UCG(RN) in principle involves evaluating a complete integral over the configuration space domain determined by rn. In practice the complete
Coarse-Grained Peptide Folding configuration space is unlikely to be sampled. Thus, the domain of configuration space actually explored implicitly plays a role in determining the values of the CG interactions obtained. As noted in section 2.2, the MS-CG models were parametrized using MD simulations containing five peptides together in a single simulation cell and many of the interactions were averaged over CG amino acid sites. Consequently, the domain of integration described in eq 1 should be understood to implicitly include these conditions, which differs from the situation in which a single folded peptide is employed and each CG site treated as having distinct interactions. This approach was utilized under the assumption that the presence of multiple peptides only serves to allow the exploration of intersite separations which are inaccessible to a single folded peptide and does not affect the properties that govern peptide folding. There is indeed evidence to support this assumption. For example, the all-atom peptides maintained folded configurations and did not form unstructured aggregates during the simulations. Additionally, our previous studies have shown that probability distributions along internal degrees of freedom for the MS-CG peptides agree well with those observed for single, folded, allatom peptides.25 This is also demonstrated in Figure 4, where the pseudo dihedral space explored by CG representations of single peptide, all-atom trajectories and the corresponding MSCG simulations is seen to be similar. Finally, the MS-CG peptides are able to correctly reproduce folded configurations,25 indicating that the effective interactions perform as intended. These observations suggest that the folded state which defines the desired domain of integration in eq 1 is not fundamentally altered by the presence of multiple all-atom peptides in the original simulation system, at least in the cases studied. Despite the evidence to the contrary, it does remain possible that the presence of multiple peptides in the MS-CG force-matching procedure has subtle, unanticipated effects on the physical properties of the resulting MS-CG peptides. Future studies will focus on this issue. In the context of the above discussion, it is useful to employ the theoretical framework described in eqs 1 and 2 to understand the relationship between observed features of the MS-CG models and characteristics of the all-atom systems. One aspect of the models that will become apparent in subsequent sections is that the MS-CG peptides exhibit a preference for folded regions of configuration space. This observation is a consequence of the models being developed using information from folded peptides. As a result, effective interactions in the CG systems are to some degree incompatible with regions of configuration space which differ from those typically explored by folded peptides. This effect biases sampled CG configurations toward the original domain of all-atom integration by facilitating fluctuations that convey the system down the relevant free energy gradients. As a result unfolded configurations are disfavored and progress toward the folded basin accelerated. As will be shown below, this characteristic allows the MS-CG peptide simulations to move more efficiently toward the folded basin when initiated in unfolded configurations (see Figures 2 and 3). In contrast, the all-atom MD simulations require significantly more simulation time to traverse the underlying configuration space. In fact, they are not observed to fold on the time-scale of the present all-atom simulations. Apart from the sampling bias inherent to the CG effective interactions, roughness of the underlying landscape and increased friction of solvent also contribute to limit configurational sampling in the all-atom systems. This observation is consistent with the results of our previous study, which suggested the potential of
J. Phys. Chem. B, Vol. 112, No. 41, 2008 13085 the MS-CG models to enhance sampling.25 The funneled folding landscapes obtained are similar to those produced by Goj models in that they are smoothly directed toward the native state.35 In contrast to Goj models, however, the MS-CG representations were not explicitly constructed to be smoothly funneled toward the native state; this characteristic arose as a consequence of the CG methodology. Moreover, the MS-CG models were obtained in a more systematic, rigorous way that allows for transparent connections between the properties of the model and the molecular interactions underlying the original all-atom system. This feature of the method may facilitate a more complete understanding of the physical properties of the system. 3.2. Realization of the Folded Basin. In our previous study,25 two different peptide CG schemes were described. The first, applied only to polyalanine, employed two CG sites for each amino acid: one each for the side-chain and peptide backbone (2-site model). The other CG scheme employed was that used in the present study (4-site model). As the folding properties of these models were studied, it was found that the 2-site model was able to refold into the right-handed helix present in the all-atom MD simulations from which the MSCG model was derived. However, approximately half of the refolding trajectories were observed to fold into left-handed helices. While it is technically possible to populate this region of the dihedral map, such configurations are found quite rarely in the Protein Data Bank and conventional MD simulations do not often report this phenomenon. Moreover, such helices were never observed for refolding trajectories of the 4-site model. The relatively high percentage of left-handed helices observed for the 2-site model likely stems from the fact that, when all interactions involving a given CG type are set to be the same, the resolution of this model is too low to differentiate between stereoisomers (such as helices of opposite handedness). This emphasizes the care that must be taken in choosing a particular CG scheme, even given the many inherent strengths of the MSCG approach. Due to this observation, further characterization of folding was only carried out for the 4-site model; only these results are discussed here. To test the folding efficacy of the models, different initial configurations were used to start trajectories for each of the two peptides as described in Methods; henceforth these will be referred to as refolding trajectories. In each of these cases the all-atom peptides were not observed to refold on the time scale of the MD simulations. In contrast, the MS-CG models were observed to successfully reconstitute folded peptide configurations in each of the CG refolding MD trajectories for both R and β peptides. Average refolding times were 5.9 ns for A15 and 1.3 ns for V5PGV5. The refolding time for each trajectory was defined as the time required for the peptide to achieve an rmsd of e1 Å from the folded configuration. Note that these times are not directly comparable to all-atom time scales due to the inherent sampling advantages of the MS-CG models. However, these numbers provide a qualitative illustration of the overall roughness of the folding free energy landscapes in each peptide model. Refolding times for individual trajectories are provided as Supporting Information. To facilitate comparison of all-atom and MS-CG data, the all-atom trajectories were first coarse-grained to generate CG configurations. The resulting free energy landscapes were projected unto the radius of gyration (Rg) and the root-mean-squared deviation (rmsd) from the folded state (Figures 2 and 3). Other variables were also investigated including the fraction of native hydrogen bonds and the number of contiguous hydrogen bonding segments
13086 J. Phys. Chem. B, Vol. 112, No. 41, 2008 (e.g., see Figure 6). The general observations described below are unaffected by employing alternative combinations of coordinates. Figure 2 compares all-atom and MS-CG simulations of the folded basin while Figure 3 compares refolding trajectories. One can observe that all-atom and CG landscapes for the folded basin are quite similar. The plots only differ appreciably for V5PGV5, which is slightly more compact in its CG representation than its all-atom counterpart. This is apparent from the Rg at the free energy minimum. However, inspection of representative snapshots (drawn to scale) in the inset of the contour plots shows that the overall peptide structure is well preserved. It is also worth noting the degree to which the MS-CG models are able to faithfully represent the native state, as refolded structures are typically found to be within 1 Å rmsd of a representative folded configuration. The MS-CG description of the folded basin is the same irrespective of the specific initial conditions employed; this can be seen by comparing CG simulations initiated from the folded basin with magnified views of the folded basin in CG refolding trajectories (see right-most panels in Figure 2). This result demonstrates the robust nature of the MS-CG models in representing the folded state. 3.3. MS-CG Refolding Landscapes Exhibit Two-State Behavior. A more complete picture of the process of folding may be obtained by examining the refolding landscape, as is done in the Figure 3. Recall that the all-atom peptides never reach the folded state on the time scale of the MD refolding trajectories. Also note that overall shapes of the all-atom and CG refolding landscapes appear qualitatively similar in their regions of overlap. The MS-CG refolding landscapes reflect the distinctive two-state folding signature of small, rapidly folding peptides, with a single dominant non-native minimum present in addition to the native basin. This behavior has been demonstrated for a number of small helices36-38 and hairpins39-41 and appears to be a general feature of small peptides with a relatively stable structure. This characteristic can be attributed to the collective nature of the transition that arises when folding depends on the formation of contacts distributed throughout the peptide molecule. It is somewhat remarkable that the MS-CG models are able to accurately represent this feature of folding given that they were not explicitly designed to do so. As noted in section 3.2, the all-atom trajectories employed in this study do not refold on the time scale of the simulations. However, comparisons of the folding landscapes presented in the studies cited above with the MS-CG landscapes presented in Figure 3 reveal similar general features including a broad, relatively shallow non-native minimum and a concerted transition to the folded state (which is the global minimum). The concerted transition toward the folded state is more pronounced in the refolding of polyalanine (Figure 6). The characteristic non-native minimum observed stems from two main sources. The first is nonspecific collapse of the peptide due to the hydrophobic effect. The second involves favorable peptidesolvent interactions that stabilize the resulting compact nonnative state. For hydrophobic peptides such as those probed in this study, these interactions primarily include HBs between the peptide backbone and solvent. The fact that the MS-CG peptides eventually fold during the limited MD simulation time suggests that peptide-peptide HBs are more favorable than HBs between peptide groups and water in the MS-CG models. This is consistent with results presented in our previous work25 demonstrating that the attractive well for MS-CG pair interactions is more favorable for peptide-peptide HBs than for the corresponding peptide-solvent interactions (see Supporting
Thorpe et al. Information of that study). This observation provides another example of the ability of the MS-CG methodology to incorporate the fundamental physical properties of the underlying molecular system without ad hoc parametrization. The results presented in this study suggest that one can realistically describe the processes that govern peptide folding with the MS-CG models. Note that the observed two-state folding behavior might not have been predicted a priori if one considers the discussion presented in section 3.1. Because the non-native minimum observed is outside the original domain of integration (i.e., the set of folded peptide MD configurations used for the development of the MS-CG models), the ability of these models to represent this feature of peptide folding is unanticipated. A possible explanation for this observation is that the effective interactions are, to some degree, transferable between configuration spaces. One would expect transferability to increase in regions of configuration space that are more similar to the folded basin. For both peptides it is apparent that the native and nonnative minima are not widely separated (Figure 3). A relatively high degree of similarity between the two configuration spaces may be one reason non-native minima are reasonably welldescribed by the MS-CG approach, particularly for polyalanine. While farther apart with respect to the variables employed in Figure 3, the native and non-native minima of the β-hairpin peptide exhibit shared characteristics which may allow a single set of interactions to be representative of both configuration spaces. This issue will be examined further in section 3.5. Transferability of the effective CG interactions may also be enhanced because the MS-CG approach effectively incorporates crucial atomistic-scale three-body correlations into the CG interactions.42 This feature of the methodology may provide the flexibility needed to more correctly represent different peptide configurations and allow for greater transferability. It is often difficult for CG systems employing a single set of pairwise effective interactions to represent different thermodynamic state points.43 However, such difficulties may arise from a lack of the higher order correlations which are included in the MS-CG approach.42 The extent and origin of transferability for the MSCG interactions is being actively investigated and will be further examined in future studies. 3.4. Hydrogen Bonding Is Enhanced in the MS-CG Peptides. The MS-CG models are able to effectively describe hydrogen bonding via NBB-OBB interactions within the peptides as well as via the relevant peptide-solvent and solventsolvent interactions. Hydrogen bonding between peptide groups is, however, more favored in the MS-CG models than in the all-atom systems (Figure 6). Note that that the average HB percentage (which, for polyalanine, roughly corresponds to average helicity) is significantly increased in the MS-CG refolding simulations relative to the corresponding all-atom simulations (Table 1). Enhanced HB formation is also reflected in the increased average length of contiguous HB segments for the MS-CG peptides. The larger propensity to form HBs is manifest in the MSCG refolding landscape of polyalanine as a sharp all-or-nothing folding transition. Hydrogen bonding in CG polyalanine is seen to be more cooperative than that observed for the CG hairpin (see Figure 6). This is likely due to topological features of the two peptides. The structure of the helix is determined by a larger number of interactions distributed throughout the peptide. The presence of helical configurations along the peptide chain favors HBs in adjoining regions of the chain by preorganizing the peptide backbone into conformations amenable to HB formation. Moreover, these backbone configurations strengthen HBs by
Coarse-Grained Peptide Folding excluding solvent from the local environment.44 While these factors also play a role in determining the collectivity of β-hairpin HB formation, their effects are diminished in such peptides. One reason for this observation is that the β-hairpin topology is inherently more flexible compared to the relatively rigid R-helix. As a result, formation of a given HB in a β-hairpin does not reduce the entropy of nearby residues as much as it does in an R-helix.36 In addition, solvent is not excluded from the local environment of a hairpin HB to the same extent as it is for the helix. This situation may be why only a small fraction of HBs is observed even in the folded all-atom hairpin (Figure 6 and Table 1). Despite the small number of HBs present in any given snapshot, overall hairpin structure is preserved during the course of the trajectory because different HBs participate in stabilizing the structure at various points in time. The tendency for enhanced HB formation in CG refolding trajectories relative to the all-atom simulations may be due to the bias inherent to the MS-CG models. This bias may induce the peptides to stabilize HBs in order to more closely emulate the folded configurations from which the effective MS-CG interactions were obtained (see section 3.1). It is interesting to note that Han and Wu have made similar observations with respect to increased HB propensity in a recent CG peptide study.45 The CG peptide models employed by these authors were composed primarily of alanine residues and had a resolution similar to that used in the present investigation. Potential energy terms for the CG models developed in that study were fit to engender agreement between CG simulations and a combination of observables from all-atom simulations and experiment. The functional form of bond and angle terms employed in their potential energy function is similar to that described in this study. These were parametrized using databases of X-ray crystallographic data. The dihedral energy terms used by these authors also possess the same functional form described in our work. Han and Wu optimized these terms to refold a polyalanine helix as well as to sample distinct regions of the dihedral map. Nonbonded interactions were parametrized by attempting to match physical properties (such as density, compressibility, diffusion constants, partition coefficients and solvation free energies) of pure liquids and liquid mixtures. The manner in which the latter interactions were obtained is similar to methods developed by Marrink and co-workers.6 The presence of enhanced HB formation observed in both sets of CG studies suggests there are fundamental elements which underlie both models. To different degrees, both sets of CG models incorporate information from folded peptide configurations. This may be the common thread that links the two sets of observations. Despite this shared observation, the approaches employed by Han and Wu and in the present study have quite different philosophies. In effect Han and Wu have parametrized their model to provide the characteristics they desire by employing a phenomenological approach. While the desired phenomenon (reproducible folding) is achieved, the parameters for their CG potential energy function are obtained from a variety of sources. Thus, it is not clear whether the balance between individual terms represents that found in realistic peptide systems. In contrast, each of the energy terms employed in the present study is systematically obtained in a consistent way from the same set of all-atom data. Thus, the results noted in this work arise naturally from the MS-CG procedure without explicitly constructing the models to achieve these effects.
J. Phys. Chem. B, Vol. 112, No. 41, 2008 13087 3.5. Configurational Sampling Is Facilitated in the MSCG Simulations. The ability of the MS-CG models to sample regions of configuration space not explored by the all-atom refolding trajectories is notable. The degree to which the CG systems reproduce the probability distributions found for the all-atom peptides can be assessed by comparing the all-atom and CG landscapes in their regions of overlap. These areas appear very similar for the β-peptide, where the unfolded peptide distribution is observed to exhibit roughly the same shape in both all-atom and CG landscapes. Note that the hairpin exhibits a small region of relatively low free energy adjacent to the main folded basin. This region is due to disorder in the termini of the β-strands that increases backbone rmsd without dramatically affecting the peptide’s overall topology. For the purposes of this discussion this region will be considered part of the folded basin. For the helical peptide all-atom and MS-CG refolding landscapes appear more dissimilar. The non-native minimum observed for CG polyalanine occupies a relatively small region of the landscape while that exhibited by the CG β-hairpin is extended and is distributed across a range of configuration space (Figure 3). This observation is consistent with a sharp twostate folding transition for the helical peptide and a less abrupt transition for the hairpin. Apart from the topological features of the two peptides discussed in section 3.4, another source of these observations could be that the driving force for refolding in the hairpin is less than that for the helical peptide. Unfolded peptide configurations are likely to be more extended (i.e., possess a larger Rg) than folded configurations in general. Thus, they tend to be more similar to configurations displayed by peptides exploring the β region of the dihedral map than those exploring R-helical dihedral space. This can be confirmed by examining the dihedral angle distributions presented in Figure 4. The unfolded all-atom hairpin explores regions of dihedral space that are already well sampled in the hairpin’s folded basin. In contrast the dihedral space explored by the folded all-atom helical peptide is only a small subset of the space sampled by the unfolded peptide. Thus, a corollary of the arguments presented in section 3.1 is that unfolded configurations are not as strongly disfavored in the MS-CG hairpin due to the similarities in the local environment of the peptide backbone for folded and unfolded configurations. In contrast, the MSCG helical peptide greatly disfavors extended configurations because they are quite dissimilar from the peptide’s folded state. Consequently, when this peptide is unfolded there is a strong bias that favors the folded basin, causing the model to spend very little time in extended configurations. This is likely why the non-native minimum present in the refolding landscape of the CG helix is much more compact than that observed for the corresponding all-atom landscape. This is demonstrated by comparing representative structures from the two sets of simulations as shown in Figure 3. In contrast, configurations from the non-native minima for refolding trajectories of the hairpin are similar for both MS-CG and all-atom systems. Another source of the bias toward the native basin for the MS-CG helical peptide may be that the helical native basin possesses a free energy minimum that is more stable than the β-hairpin. The β-hairpin possesses only five native HBs compared to ten for polyalanine and thus will be much less stabilized by HB interactions. If one examines the contour plots for MS-CG simulations presented in Figure 2 it is apparent that the free energy increases much more quickly for the helix than for the hairpin as distance from the folded minimum is increased. Thus, the helical peptide is subject to much steeper free energy
13088 J. Phys. Chem. B, Vol. 112, No. 41, 2008 gradients as it moves away from the folded state. Evaluating effective interactions using folded configurations and then employing these interactions to represent unfolded configuration essentially extrapolates the MS-CG interactions beyond the region of configuration space from which they were originally developed. If one linearly extrapolates the free energy change for MS-CG simulations of the folded state with respect to either of the variables employed for the contour plots displayed in Figure 2, it is apparent that a much larger free energy is obtained at Rg or rmsd values consistent with unfolded configurations for the helix than would be obtained for the hairpin. For example, using the contours presented in Figure 2, and assuming a small displacement from the folded basin (e.g., 0.2) in either rmsd or Rg dimensions, one estimates the free energy change of unfolding per unit of Rg or rmsd to always be ∼2 kBT more for the helical peptide than for the β-hairpin. This corresponds to ∼1.2 kcal/mol at 310 K. These simple considerations suggest a larger driving force for refolding of the CG helix than the CG hairpin at any given region in the landscape. As pointed out previously, the all-atom refolding simulations are never observed to visit the folded basin. In the case of polyalanine, the all-atom simulations are not even able to explore the compact non-native minimum observed for the CG models. Estimates for the folding time of small peptides such as these range from hundreds of nanoseconds to microseconds.46 This is significantly longer than the average of 5.9 or 1.3 ns needed to observe folding for MS-CG polyalanine or V5PGV5 respectively. One reason the MS-CG models fold more readily is likely that the driving force propelling the MS-CG peptides toward the folded state is larger than that for the all-atom systems. Unlike the MS-CG peptides, a bias toward folded configurations has not been expressly incorporated into the allatom models. This makes the configuration space accessible to the all-atom peptides much greater than that available to the MS-CG models. Reduction in accessible conformation space for the MS-CG systems serves to enhance sampling along preferred directions (i.e., toward the folded state). A separate but related issue is the reduction in dimensionality for CG systems compared to the original all-atom systems. The accompanying decrease in entropy facilitates the conformational search for the folded state. However, previous analyses indicate that this reduced dimensionality would only account for a speed up of approximately a factor of 3.25 This modest acceleration is not sufficient to account for the sampling differences demonstrated in Figure 2. If the average refolding times are multiplied by three to estimate an effective time scale for MSCG refolding, the resulting effective folding times of ∼18 ns and ∼8 ns for A15 and V5PGV5 respectively are still well below experimental estimates. An additional source of sampling enhancement is that the CG simulations possess smoother free energy landscapes due to the averaging inherent to the MS-CG approach. One component of this smoothing is the reduced roughness of the free energy landscape associated with the CG solvent representation. Because peptide fluctuations are determined to a large extent by solvent fluctuations,47 it is likely that this feature is a significant determinant of the sampling enhancement observed. The enhanced sampling capabilities of MS-CG models suggest a number of potential applications such as examining large scale conformational changes, protein folding studies and protein structure prediction. A possible use for the models described in this work is to study processes coupled to peptide folding such as folding mediated protein association or ligand binding. Due to their inherent bias toward sampling the folded
Thorpe et al. state, models such as these could be employed to drive such processes in order to better study them. This can be done efficiently due to the low computational cost of the models. It is quite noteworthy that none of the observations described above were explicitly parametrized into the models. Instead, these properties arise naturally from the MS-CG methodology and are a consequence of its theoretical foundation and the physical considerations that motivate the approach. There are a number of ways to combine information from both all-atom and CG systems in a multiscale manner in order to exploit the benefits of each type of model. The approach presented in this work is to first use the MS-CG energy function to extensively sample configurations and then select specific configurations of interest in order to reconstruct all-atom representations. This method takes advantage of the speed and efficiency of the CG models for sampling while incorporating the detailed description afforded by all-atom mdoels. Similar approaches have been explored previously.48-51 The MS-CG models described in this work are particularly well-suited to this task. At this intermediate level of resolution it remains a relatively straightforward undertaking to reconstruct the missing degrees of freedom necessary to reconstitute all-atom configurations. 3.6. Ensembles Generated from Reconstructed All-Atom Configurations Are Consistent with the MS-CG Landscapes. In principle, the ability to reconstruct missing all-atom detail provides the ultimate validation for any CG model. This process allows one to determine whether the features observed in the MS-CG simulations are compatible with the properties of the original all-atom system. In this way a concrete link between different levels of resolution can be maintained at all times. It is not always a simple matter to relate minima identified in the MS-CG refolding simulations with features observed in the allatom refolding landscapes because all-atom trajectories may be unable to reach these regions in a reasonable time frame. In this study, this situation was apparent particularly for the polyalanine simulations. Since the all-atom trajectories never fully explore the non-native minimum observed in the MS-CG simulations of this system, it is hard to determine whether this feature is a realistic representation of the peptide landscape or is merely an artifact due to the MS-CG model. To resolve this issue, a CG peptide structure from the non-native minimum observed for each MS-CG peptide was used to reconstruct atomically detailed configurations which were then employed to generate a number of different trajectories using the original all-atom energy function (see Supporting Information). It was observed that free energy landscapes generated in this way exhibit reasonable stability in the vicinity of the minima originally observed in the MS-CG simulations (Figure 3). This result demonstrates that the MS-CG models do realistically represent peptide refolding landscapes and provides further validation for the models. However, unlike the MS-CG simulations, the reconstructed all-atom simulations are still unable to surmount the barrier to folding and do not realize the folded basin. Instead these trajectories avoid the barrier region observed in the MS-CG landscapes and spread out throughout the non-native basin. This result suggests that a barrier does exist between the native and non-native basins for the all-atom peptides that the MS-CG peptides were able to surmount. Thus the barriers present in MS-CG simulations appear to be at approximately correct locations in the refolding landscapes, though presumably are of lower magnitude than those present in the all-atom simulations. Barriers observed for the MS-CG simulations have a height on the order of 6kBT for A15 and only 2kBT for V5PGV5
Coarse-Grained Peptide Folding with respect to the non-native minimum observed for each peptide. The barrier to folding in the all-atom peptides is likely to be significantly larger. Moreover, it is anticipated that the MS-CG peptides will be able to traverse the folding barrier in a more facile manner than all-atom models due to the reduced energetic frustration that results from a smoother landscape. This frustration may make it very difficult for the all-atom trajectories to locate the narrow channel between the folded basin and the non-native minimum present in each of the MS-CG refolding landscapes (Figure 3). In addition to validating the MS-CG simulations, the ability to connect MS-CG configurations with all-atom ensembles enhances the sampling capabilities of approach. In this context MS-CG models can be employed to explore configuration space efficiently. The CG peptides could then be used to reconstruct all-atom configurations when regions of interest are identified. These regions could then be studied exhaustively with the detailed all-atom model to fully characterize their properties. As a result, less time and computational effort is wasted exploring uninteresting regions of configuration space with potentially computationally intensive and inefficient (from a sampling standpoint) all-atom simulations. In principle it should be possible to reconstruct the entire landscape obtained from the MS-CG simulations, generating information about the barriers to refolding along the pathways revealed by the MSCG models. One could further extend this concept by augmenting the sampling capabilities of the MS-CG model with other enhanced sampling techniques such as replica exchange simulations. This allows for the most effective use of resources by focusing computational effort on models that are already more efficient from a sampling perspective. A final advantage of the reconstruction procedure is that it can correct for deficiencies in the original CG model. Thus, even a CG model that has limited overlap with the configuration space of the all-atom system can be employed. As long as the CG model traverses regions of configuration space that are at least approximately similar to those explored by the all-atom system, reconstructed all-atom configurations can provide ensembles corresponding to the original energy function. The all-atom energy function provides the “exact” result in the context of the MS-CG models since these models were generated from and designed to replicate all-atom ensembles. The results observed in this study as well as the possible applications described in the previous paragraphs indicate the powerful potential of this methodology. Converting between allatom and CG representations as necessary ensures that a concrete link between descriptions of the system at each scale or resolution level can be maintained. We are currently engaged in exploring potential applications of this approach.
J. Phys. Chem. B, Vol. 112, No. 41, 2008 13089 from the three factors. First, the MS-CG models possess fewer degrees of freedom, reducing the search space for peptide folding. Second, the roughness of the free energy landscape is reduced in the MS-CG models. Finally, due to the physical basis of the models, the resulting effective interactions incorporate a bias toward the regions of configuration space from which the interactions were originally obtained. These combined effects are enough to reduce the barrier to peptide folding such that it is readily observable on simulation time scales. The MS-CG peptide models were used to reconstruct corresponding all-atom configurations; all-atom simulations initiated from these configurations reiterate the observations made for the MS-CG refolding simulations. Not only does the reconstruction approach provide validation for the MS-CG simulations, it further enhances the sampling capabilities of the models by providing access to the “correct” underlying all-atom distributions. In this way, the best attributes of both MS-CG and atomistic models can be utilized: MS-CG simulations for enhanced sampling and all-atom detail for more accurate thermodynamic and kinetic properties. This approach provides a general strategy by which MS-CG methods and atomically detailed models can be combined in a truly multiscale manner. The properties of the MS-CG method provide a concrete link between characteristics of the models and properties of fully detailed trajectories, facilitating interpretation of the observations. Moreover, the systematic nature of the MS-CG approach facilitates fine-tuning of the method for a range of applications. It is anticipated that future studies will employ models such as these to provide insight into the physical principles that govern protein structure and function. Acknowledgment. This research was supported by the National Science Foundation Collaborative Research in Chemistry program (CHE-0628257). Computations were carried out using an allocation of advanced computing resources supported by the National Science Foundation (No. MCA94P017) via the TeraGrid and provided by the National Center for Supercomputing Applications and Indiana University. The authors thank Dr. C. Mark Maupin for assistance with preparation of Figure 3. Supporting Information Available: Text giving a more detailed description of the all-atom reconstruction method for the MS-CG peptide models as well as figures showing additional free energy landscapes for simulations initiated from reconstructed configurations and a table giving a detailed listing of simulation and refolding times. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes
4. Concluding Remarks In this study, multiscale coarse-grained (MS-CG) peptide models25 were shown to efficiently refold. The peptide models demonstrate realistic refolding landscapes despite being generated using all-atom information from folded configurations. This observation suggests the potential for transferability of the MSCG effective interactions. The basis of transferability may stem from the capacity of the methodology to incorporate higher order structural correlations into these interactions. An added benefit of the MS-CG models is that the approach is based on a systematic theoretical framework. The MS-CG peptides demonstrate the ability to enhance sampling of the underlying configuration space, visiting regions that were not readily attainable using all-atom simulations. These observations stem
(1) Tozzini, V. Curr. Opin. Struct. Biol. 2005, 15, 144. (2) Ayton, G. S.; Noid, W. G.; Voth, G. A. Curr. Opin. Struct. Biol. 2007, 17, 192. (3) Ayton, G. S.; Noid, W. G.; Voth, G. A. MRS Bull. 2007, 32, 929. (4) Muller-Plathe, F. ChemPhysChem 2002, 3, 755. (5) Marrink, S. J.; de Vries, A. H.; Mark, A. E. J. Phys. Chem. B 2004, 108, 750. (6) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. J. Phys. Chem. B 2007, 111, 7812. (7) Shelley, J. C.; Shelley, M. Y.; Reeder, R. C.; Bandyopadhyay, S.; Klein, M. L. J. Phys. Chem. B 2001, 105, 4464. (8) Shelley, J. C.; Shelley, M. Y.; Reeder, R. C.; Bandyopadhyay, S.; Moore, P. B.; Klein, M. L. J. Phys. Chem. B 2001, 105, 9785. (9) Izvekov, S.; Voth, G. A. J. Phys. Chem. B 2005, 109, 2469. (10) Izvekov, S.; Voth, G. A. J. Chem. Phys. 2005, 123, 134105. (11) Stevens, M. J. J. Chem. Phys. 2004, 121, 11942. (12) Shih, A. Y.; Freddolino, P. L.; Arkhipov, A.; Schulten, K. J. Struct. Biol. 2007, 157, 579.
13090 J. Phys. Chem. B, Vol. 112, No. 41, 2008 (13) Bond, P. J.; Holyoake, J.; Ivetac, A.; Khalid, S.; Sansom, M. S. P. J. Struct. Biol. 2007, 157, 593. (14) Buchete, N. V.; Straub, J. E.; Thirumalai, D. Protein Sci. 2004, 13, 862. (15) Das, P.; Moll, M.; Stamati, H.; Kavraki, L. E.; Clementi, C. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 9885. (16) Khalili, M.; Liwo, A.; Scheraga, H. A. J. Mol. Biol. 2006, 355, 536. (17) Murarka, R. K.; Liwo, A.; Scheraga, H. A. J. Chem. Phys. 2007, 127. (18) Marianayagam, N. J.; Fawzi, N. L.; Head-Gordon, T. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 16684. (19) Chu, J. W.; Voth, G. A. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 13111. (20) Chu, J. W.; Voth, G. A. Biophys. J. 2006, 90, 1572. (21) Arkhipov, A.; Freddolino, P. L.; Imada, K.; Namba, K.; Schulten, K. Biophys. J. 2006, 91, 4589. (22) Kim, Y. C.; Hummer, G. J. Mol. Biol. 2008, 375, 1416. (23) Feig, M.; Onufriev, A.; Lee, M. S.; Im, W.; Case, D. A.; Brooks, C. L., III J. Comput. Chem. 2004, 25, 265. (24) Lu, B. Z.; McCammon, J. A. J. Chem. Theory Comput. 2007, 3, 1134. (25) Zhou, J.; Thorpe, I. F.; Izvekov, S.; Voth, G. A. Biophys. J. 2007, 92, 4289. (26) Onuchic, J. N.; Wolynes, P. G. Curr. Opin. Struct. Biol. 2004, 14, 70. (27) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (28) Smith, A. V.; Hall, C. K. Proteins: Struct., Funct., Bioinf. 2001, 44, 376. (29) MacKerell, A. D., Jr.; Brooks, B.; Brooks, C. L., III; Nilsson, L.; Roux, B.; Won, Y.; Karplus, M. CHARMM: The Energy Function and its Parameterization with an Overview of the Program. In The Encyclopedia of Computational Chemistry; Schleyer, P. v. R., Ed.; John Wiley and Sons: Chichester, U.K., 1998; Vol. 1; p 271. (30) Sung, S. S. Biophys. J. 1999, 76, 164.
Thorpe et al. (31) Smith, W.; Forester, T. R.; Todorov, I. T.; Leslie, M. DL _POLY, 2.16 ed.; CCLRC Daresbury Laboratory: Daresbury, Warrington, U.K., 2006. (32) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (33) Noid, W. G.; Chu, J.-W.; Ayton, G. S.; Krishna, V.; Izvekov, S.; Voth, G. A.; Das, A.; Andersen, H. C. J. Chem. Phys. 2008, 128, 244114. (34) Noid, W. G.; Liu, P.; Wang, Y. T.; Chu, J.-W.; Ayton, G. S.; Izvekov, S.; Andersen, H. C.; Voth, G. A. J. Chem. Phys. 2008, 128, 244115. (35) Go, N.; Abe, H. Biopolymers 1981, 20, 991. (36) Onuchic, J. N.; Socci, N. D.; Luthey-Schulten, Z.; Wolynes, P. G. Fold. Des. 1996, 1, 441. (37) Levy, Y.; Jortner, J.; Becker, O. M. Proc. Natl. Acad. Sci. U.S.A. 2001, 98, 2188. (38) Sorin, E. J.; Pande, V. S. Biophys. J. 2005, 88, 2472. (39) Snow, C. D.; Qiu, L.; Du, D.; Gai, F.; Hagen, S. J.; Pande, V. S. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 4077. (40) Streicher, W. W.; Makhatadze, G. I. J. Am. Chem. Soc. 2006, 128, 30. (41) Baumketner, A.; Shea, J. E. Theor. Chem. Acc. 2006, 116, 262. (42) Noid, W. G.; Chu, J.-W.; Ayton, G. S.; Voth, G. A. J. Phys. Chem. B 2007, 111, 4116. (43) Johnson, M. E.; Head-Gordon, T.; Louis, A. A. J. Chem. Phys. 2007, 126. (44) Wang, Z. X.; Duan, Y. J. Theor. Comput. Chem. 2005, 4, 689. (45) Han, W.; Wu, Y.-D. J. Chem. Theory Comput. 2007, 3, 2146. (46) Lednev, I. K.; Karnoup, A. S.; Sparrow, M. C.; Asher, S. A. J. Am. Chem. Soc. 1999, 121, 8074. (47) Frauenfelder, H.; Fenimore, P. W.; Chen, G.; McMahon, B. H. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 15469. (48) De Mori, G. M.; Colombo, G.; Micheletti, C. Proteins 2005, 58, 459. (49) Lyman, E.; Ytreberg, F. M.; Zuckerman, D. M. Phys. ReV. Lett. 2006, 96, 028105. (50) Heath, A. P.; Kavraki, L. E.; Clementi, C. Proteins 2007, 68, 646. (51) Liu, P.; Voth, G. A. J. Chem. Phys. 2007, 126, 045106.
JP8015968
