J. Phys. Chem. B 2009, 113, 9669–9680
9669
Searching the Conformational Space of Cyclic β-Amino Acid Peptides Fredy Sussman,* M. Carmen Villaverde, Juan Carlos Este´vez, and Ramo´n J. Este´vez Departamento de Quı´mica Orga´nica, Facultad de Quı´mica, UniVersidad de Santiago de Compostela, 15782 Santiago de Compostela, Spain ReceiVed: December 22, 2008; ReVised Manuscript ReceiVed: May 20, 2009
There is an increasing interest in the secondary structure of β-amino-acid-containing peptides, since these compounds exhibit an intrinsic propensity to form stable folds even for short peptides, a feature that is rarely observed in R-amino-acid-containing peptides. In this work, we use a multiple trajectory molecular dynamics approach to study a panel of cyclic β-amino-acid-containing peptides with a variety of motifs that differ in the ring size, ring substituents, and terminal protecting groups. We find a reasonable agreement between the predicted and the experimentally observed structures, in spite of the simple solvent representation used, indicating that in most cases the folding proceeds energetically downhill and it is driven to a great extent by structural preferences coded in the internal degrees of freedom, a result supported by our energy partition analysis. Our results also show that when the N-terminal end is unprotected, it is likely to be charged in a protic polar solvent. In that case, we find that only a molecular dynamics simulation with an “all atom” solvent representation is capable of reproducing the experimentally observed secondary structure of the peptide. Finally, the time evolution analysis of the hydrogen-bond-induced turns as well as of the root-mean-square deviation from the observed structure indicates that some peptides could have a higher intrinsic flexibility than others, within a given fold, a result that correlates to some degree with our molecular mechanics energy analysis. Introduction A growing body of research on nonnatural amino acids (e.g., β- and γ-amino acids) oligomers has identified specific secondary structures with a strong propensity to adopt compact ordered structures for small peptides, a feature rarely observed in R-amino acids. The intrinsic structuring observed in these molecules (dubbed foldamers on this account) represents an advantage in the design of larger systems with tertiary structure (e.g., receptors and enzymes) built upon these kinds of amino acids.1-11 Moreover, some of the members of these new families of peptides have found applications in the search for peptidomimetic drug leads, since they have been shown to have an advantageous bioavailability profile when compared to their R-amino-acids-based counterparts.3 Other groups have explored the possibility of building supramolecular structures (e.g., nanotubes, pores) with a desired functionality, starting from peptides made of unnatural amino acids with restricted conformational flexibility.9 The ability of these peptides to aggregate into supramolecular entities, given the right complementary sequence, has led to the design, synthesis, and characterization of β-peptide-based helical bundles with a built in functionality.10 The research leading to the identification of a growing number of peptide folds in β-amino-acid-based oligomers has experienced a notable surge.1-8,10,11 In the case of the cyclic β-amino acids, it has been shown that the chemical nature of the cyclic side chains modulates the fold exhibited by the peptides. For instance, the group headed by Gellman has shown that transcycloalkane containing β-amino acids like the trans-2-aminocyclohexanecarboxylic acids (ACHC) form peptides that fold as R-helices with a 14-atom turn while trans-2-aminocyclopentacarboxylic acids (ACPC) containing peptides form a 12-atom * To whom correspondence
[email protected].
should
be
addressed.
E-mail:
turn R-helix.4,5 In some cases (e.g., in ACHC-based peptides), the number of amino acids in a peptide as well as the terminal blocking groups have a strong effect on the turn size of the resulting helix.6 The group headed by Fleet demonstrated that even shorter turns can be produced by four-membered ring (oxetane) containing amino acids.8 It also has been shown that changing the conformation in the additional torsional angle (θ) in β-peptides may profoundly affect the secondary structure fold of these molecules. For instance, the group headed by Fu¨lo¨p has shown that ACPC-based peptides with a cis conformation at the θ torsion angle lead to an extended strand with a zigzag pattern, rather than the helix fold observed when θ is in the trans conformation.7 The computer-assisted prediction of β-peptides in general and of cyclic amino acids containing β-peptides in particular has become a field of increasing interest.4-8,11-14 The design of the first β-amino acid oligomers was based on short and simple molecular mechanics/molecular dynamics protocols without the inclusion of an explicit solvent model, as in the case of the trans-ACPC and ACHC peptides.4 Nevertheless, not all the molecular mechanics based protocols reproduce well the resulting secondary structure of β-amino acid homo-oligomers. For instance, not even a Monte Carlo based full molecular mechanics structural search was able to predict the folding pattern for oxetane-containing oligomers.8 The most sophisticated calculations have been carried out at a quantum mechanical (QM) level.6,11 However, the QM calculations have been restricted to rank some of the possible folding patterns through the comparison of the energetics of these resulting secondary structures instead of simulating the folding of these peptides from scratch. While this approach may be adequate for some homo-oligomers, it is less acceptable for small heterogeneous polypeptides where there may be coexisting structural motifs in the same structure, which might not be known a priori. For instance, Seebach et
10.1021/jp811321n CCC: $40.75 2009 American Chemical Society Published on Web 06/25/2009
9670
J. Phys. Chem. B, Vol. 113, No. 29, 2009
Sussman et al.
TABLE 1: Compounds Studied and Predicted Folds
a SA denotes the simulated annealing protocol. The other protocols are identified by the temperature value in K and by the dielectric constant value (DC). b Most predicted folds agree with the experimentally reported structures. See text for details on compounds 1a and 3. For compounds 5a and 5b there are not experimentally reported structures. c Frayed at the C-terminal end.
al. have shown that hetero-oligomer β-peptides built from acyclic β-amino acids display folds that contain more than one kind of turn.3b Many of the previous molecular mechanics based predictions were carried out with a variety of protocols that differ in many aspects (e.g., force fields and solvent representation, etc). The central goal of this work is to determine whether a multiple copy molecular dynamics protocol is able to “fold” a set of chemically diverse cyclic β-amino-acid-containing peptides (see Table 1), starting from an extended conformation. Our protocols use different starting conditions and/or different starting conformations in order to improve the chances for a more successful sampling of the conformational space of the peptides. We show that the proposed approach is able to reproduce the observed fold of many of the peptides studied, and in some cases to shed light onto the physical interactions that help direct the folding process. For instance, it is known that both R-amino-acid-based peptides and proteins fold on account of the so-called hydrophobic effect, where the solvent plays a central role, while the folding of the cyclic amino-acid-containing β-peptides into a given secondary structure is believed to be directed by a strong internal propensity generated by the additional backbone degree of freedom (torsion angle θ), and only partially modulated by the solvent, a premise behind all of the MD and QM calculations carried out until now. Our energy-partitioning analysis supports this rationale. It shows that for most peptides studied (with one
exception) the internal degrees of freedom energy components favor energetically the helical structure over the extended ones. Our calculations show that not all peptides fold with equal ease and that in some cases the solvent may play an important role in reaching the observed fold, especially when the peptide has an unprotected N-terminal end and it is solvated by a polar protic medium. For that case (compound 1a), our protocol has allowed to predict a charged N-terminal end, and a slightly different turn assignment, features that could not be obtained directly from the original NMR data.6 Finally, we expect that our structural search protocols could be of help as a first tier tool in the prediction of the fold of β-amino acid peptides libraries, as part of our research on the design, synthesis, and structural characterization15 of β-peptides with interesting folding propensities. Methods The cyclic β-amino-acid-containing peptides studied in this work, which include four-, five-, and six-membered rings with various substituents and N- and C-terminal capping groups, are listed in Table 1. We used two protocols for the conformational search aimed at predicting the highest structural fold for each of the peptides displayed in Table 1. The starting point for the first protocol, a simulated annealing (SA) approach, is a 1 ns molecular
Cyclic β-Amino Acid Peptides
J. Phys. Chem. B, Vol. 113, No. 29, 2009 9671
Figure 1. MD generated time evolution of the HB distances in the 14-atom turns for compound 1e. The predicted HB pattern is shown schematically at the top. All distances are in Å.
dynamics (MD) simulation protocol at high temperature (600 K) that started with an energy minimization followed by a 2 ps equilibration stage. The high-temperature MD trajectory rendered several conformations (collected every 200 ps) that were then subjected to 5000 steps of energy minimization, using a steepest descent (SD) protocol, followed by a 1 ns molecular dynamics (MD) simulation at 300 K, using the leapfrog algorithm. In line with previous protocols for the simulation of β-peptides, a distance-dependent dielectric constant is used to simulate implicitly the effect of the solvent.4 The second protocol starts with 5000 steps of energy minimization, using (as above) the SD protocol with a gradient tolerance of 0.01. Each compound is then subjected to various parallel MD simulations that differ in temperature (the temperature values range from 300 to 340 K) as well as dielectric constants (DC), with DC values of vacuum (1) and the most commonly used solvents in NMR structure determination of β-peptides, ranging from chloroform (4.7) to methanol (33.0), and DMSO (48.0). Protocols with MD searches at various temperatures have been previously used for the prediction of the folding of acyclic β-amino-acid-containing peptides, probably to enhance sampling and make it possible to surmount the small energy barriers present in these peptides.14 By the same token, the aim of using various DC values was to reduce the likelihood of structures being trapped at local minima by energy barriers set up by electrostatic interactions.
The CVFF force field16 and the Discover program17 were used throughout this work. We choose this force field because it has wide parameter availability for a large range of functional groups. Both protocols used in this work could be construed as a multiple copy molecular dynamics (MCMD) approach. Taking into account the ergodicity of the system, every peptide listed in Table 1 had its conformational space searched for well over a 10 ns time span, and we expect the MCMD conformational sampling to be more exhaustive than the one that starts with a single conformation and lasts the same amount of time. The set of simulations for every single peptide does not take more than 8 h in a single processor SGI Fuel workstation. Peptides 1b and 3 could have their N-terminal ends charged when solvated in a protic polar medium like methanol (see Results and Discussion). To overcome the spurious effects that a charged group could induce in an implicit solvent simulation, we proceeded to simulate the structure of one of these compounds (1a) in an explicit sphere of 640 methanol molecules, with a density of 0.79 g/cm3. The first step in this protocol was a MD simulation having the solute fixed, lasting 0.6 ns with the purpose of randomizing the initial solvent structure. The next MD stage left all atoms free and lasted for 6 ns. The temperature used in these simulations was 300 K and the integration time step was 1 fs. For these calculations we used the CVFF93, a close relative of the CVFF potential used
9672
J. Phys. Chem. B, Vol. 113, No. 29, 2009
Sussman et al.
Figure 3. MD-generated time evolution of the two HB distances that define the 14-atom turns in compound 1c. The predicted HB pattern is shown at the top of the figure.
found in InsightII. New insights into the folding process leading to an R-helix were obtained by calculating the difference in the energy components between the starting extended structures and the folded structures chosen to have the closest rmsd to an ideal folded structure. For comparison’s sake we also performed the same kind of analysis for an R-Ala octapeptide. The extended and R-helical conformations for this peptide were built by using the tools found in the Biopolymer module in InsightII.17 Results and Discusion
Figure 2. Snapshots of a MD trajectory for ACHC-based oligomers. The upper panel displays the octamer (1e), the middle panel depicts the tetramer with the N-terminal end protected by a t-Boc group (1c), and the lower panel displays the fold of compound 1a, the tetramer with unprotected charged N-terminal end. The black thin lines depict the hydrogen bonds present in these structures.
above, and the simulations were performed with the help of the CHARMm software suite.18 Since the structural studies for compound 3 were carried out in DMSO, an aprotic polar solvent, we felt no need for explicit solvent simulations for this compound. In order to understand the folding process and to gain insight into the forces that lead to a folded structure, we have performed several types of analysis. For instance, the helix turn formation and its variability were gauged by following the time evolution of the distances characteristic for a given hydrogen bond turn, and by evaluating the root-mean-square deviation (rmsd) from an ideal helix structure built a priori with the help of the tools
The main aim of the protocol introduced above is the search for the most structured fold that a given peptide may generate, starting from a variety of initial conditions (e.g., temperature, solvent representation, etc.) and initial conformations. In some cases, protocol variants lead to different alternative structures (for the same peptide), while others show that the predicted fold has resiliency to the MD protocol used. The resulting structures could be classified according to the degree of internal order. Table 1 summarizes the most structured folds predicted by our calculations for the peptides studied in this work, as well as the specific protocol by which it was obtained. As we shall see below, there is a reasonable agreement between the highest ordered structures from our MD simulations and the folds observed from experiment. A. Folding Prediction. In what follows we will discuss the resulting folds predicted by our simulation protocols for each of the structures studied. I. trans-ACHC Oligomers. We studied the resulting folds of trans-ACHC peptides of various lengths (ranging from tetramers to octamers) with assorted capping groups (see compounds 1a-1e). The longest peptide (compound 1e) is an octamer and it folds in a 14-atom helix in most simulations
Cyclic β-Amino Acid Peptides
Figure 4. MD-generated time evolution of the three HB distances that define the 10-14-atom turns observed in compound 1a. The first two panels depict the time profile for the HBs formed by the charged N-terminal nitrogen. The HB pattern is shown at the top.
regardless of the starting conditions. The turn type predicted by our MD simulation can be deduced from the analysis of the time evolution of the distance between the backbone polar atoms that are prone to form any possible hydrogen bond in these molecules. Figure 1 displays the time evolution for the six possible 14-atom turns in compound 1e. As seen from this figure, the peptide takes only a few picoseconds to form the 14-atom hydrogen bonds (HB) and then does not deviate grossly from these HB motifs with the exception of the one closest to the C-terminal end. Reducing the size of the peptide to a hexamer (compound 1d) does not modify the resulting secondary structure fold. The overall folding motif is in agreement with experimental results of the trans-ACHC containing hexamer with a t-Boc and OBn protecting groups at the N- and C-terminal ends, respectively.4,5 Moreover, our simulations on even shorter trans-ACHC-based peptides (i.e., a tetramer capped by a N-terminal protecting group that contains an extra amide bond (compound 1c))4-6 predict the same folding motif as for the longer oligomers (1d,
J. Phys. Chem. B, Vol. 113, No. 29, 2009 9673 1e), again in agreement with experiment.4,5 Snapshots of a MD simulation for compound 1c, a tetramer with its N-terminal protected with t-Boc, as well as for compound 1e, are shown in Figure 2, clearly displaying the helicoidal character of these molecules. Perusal of the conformational data in the literature led to the notion that the main-chain torsional angles for the 10- and the 14-atom turns with the same helicity occupy adjacent positions in the Ramachandran plot,6 allowing for a smooth rearrangement between two conformations underpinned by these HB motifs. Actually, equilibrium between the 14- and 10-atom turns have been proposed for an acyclic β-hexapeptide, based on the results of CD experiments and MD simulations, by the Seebach group.19 We have calculated the HB time evolution for the putative 10atom turns that could be formed by compound 1c. The results are shown in Supporting Information (Figure S1). As seen from the time evolution, there are a number of configurations that adopt the 10-atom turn HB motif. Nevertheless, the frequency of their occurrence is much smaller, and the oscillations out of these conformations are much bigger than the structures with 14-atom HB turns displayed in Figure 3. These results taken as a whole indicate that, although the 10-atom helix configurations are energetically accessible, they are a minority when compared to the population of 14-atom turn helix conformations. NMR experiments suggest that deprotecting the ACHC tetramer on the N-terminal end (see compounds 1a and 1b) stabilizes a helix conformation with 10-atom turns. The results of the HB analysis of our MD simulations on the neutral form of the unprotected tetramer (compound 1b) indicate that the number of 10-atom turns increase substantially upon leaving the N-terminal end unprotected and neutral (see Supporting Information, Figure S1), although they still represent a minority when compared to the 14-atom turn helix in compound 1c. A possible explanation for the difference between the experimental and simulation results could reside in the interplay between the polarity of the solvent and the net charge of the solute, which could play a substantial role in the actual fold of this peptide. The solvent used in all NMR structural studies of ACHC is deuterated methanol, a protic polar molecule. Such environment could induce the N-terminal aliphatic amine to increase its pKa value, and hence the likelihood of being protonated in this solvent. There are some precedents for this rationale. For instance, it has been observed that a heterogeneous β-amino acid hexamer changes its conformation from a mixed 10/12/10 helix to a 14-atom turn helix when the N-terminal end in this peptide is unprotected.3b These results were rationalized by proposing that the N-terminal amine in the uncapped peptide has a positive charge in methanol. Later, MD-based calculations on this and other peptides have shown that an N-terminal charge helps to stabilize the helix through a charge dipole interaction, provided the helix has an N f C polarity.20 In order to investigate the effect of the charge on the peptide fold of the ACHC containing tetramer, we have carried out MD simulations on the unprotected charged Nterminal molecule (compound 1a). The resulting MD trajectory conformations indicate the possible existence of a mixed 14/10 structure with a hitherto unobserved HB motif, in which the charged N-terminal nitrogen switches HBs with the peptide backbone carbonyl groups generating 10- and/or 14-atom turns. The HB distance-time evolution for this folding motif is shown in Figure 4. As seen from this figure, the most stable hydrogen bonds are the ones that form the 10-atom turn generated by the N-terminal nitrogen and two 14-atom turns generated by both the N-terminal and the nitrogen in the second amino acid. A MD snapshot view of the predicted fold for compound 1a,
9674
J. Phys. Chem. B, Vol. 113, No. 29, 2009
Sussman et al.
Figure 5. Time evolution of HB distances that underlie all possible 10- and 14-atom turns in compound 1a, resulting from the last 4 ns of a MD simulation with an “all atom” solvent explicit representation. The HB distance pattern is shown at the top of the figure.
displaying this unusual HB motifs is shown in Figure 2. The NMR spectroscopic data suggest that the actual secondary structure of this peptide is a 10-atom turn helix with a frayed N-terminal end.6 Nevertheless, some of NMR data presented in the original work, specifically the ratio of the NOE cross peaks (NHi-CRHi-1/NHi-CβHi) for the second residue (i.e., i ) 2), indicates that this part of the structure is closer to a 14atom helix than to a 10-atom helix.6 Given the possible charged nature of this peptide, we felt necessary to revisit the folding prediction with a MD protocol with an “all atom” solvent explicit representation, in order to rule on the actual folding motif. It is well-known that in “dry” MD simulations the charged molecules tend to generate additional interactions between the charged fragment and polar groups in the peptide in order to compensate for the shortfall in interactions with the solvent.21 For this sake we have carried out preliminary MD studies of the unprotected ACHC-containing tetramer (1a) in methanol, as described in the Methods section. We have chosen as a starting peptide conformation for this MD simulation one of the structures that results from the implicit solvation simulation described above, in which the N-terminal end is involved in
Figure 6. Snapshot of a MD simulation for the trans-ACPC-based oligomer (compound 2). The dark thin lines indicate the HB pattern.
Cyclic β-Amino Acid Peptides
J. Phys. Chem. B, Vol. 113, No. 29, 2009 9675
Figure 7. MD-generated time evolution of the three HB 12-atom turns predicted for the trans-ACPC-containing peptide (compound 2), obtained from our simulated annealing protocol. The HB pattern is shown at the top of the figure.
two HB interactions (see HB scheme in Figure 4). We have analyzed all the HBs that could be produced during the MD trajectory, and their time evolution is shown in Figure 5. Perusal of this figure indicates that the N-terminal nitrogen loses the HBs that are necessary for 10- or 14-atom turns that were present in the starting MD structure, and hence the N-terminal end becomes frayed (see the two upper panels in Figure 5). The second residue nitrogen seems to be able to form the HBs that are required for 10- and 14-atom turns, although the latter turn seems to be more populated than the former, in agreement with the aforementioned NOE ratio for this residue. Finally, we find that the closest turn to the C-terminal end contains 10 atoms. These results are in closer agreement with experiment than those obtained with an implicit solvent and are supporting evidence that the N-terminal amine group on this uncapped peptide is charged. We are performing more extensive calculations of compound 1a in an explicit solvent in order to get more information on this issue. II. trans-ACPC Oligomer. While various MD simulations presented here may result in folds with a variety of secondary structure arrangements for most peptides, a large number of
the MD simulations for the trans-ACPC-containing octamer (2) only fold into a 12-atom turn R-helix, independently of the MD conditions and starting parameters used in this work. The persistence of the folding pattern in the MD simulations indicates that this peptide together with the trans-ACHC octapeptide are some of the foldamers with the highest intrinsic proclivity to a given fold. The most frequent fold predicted by our MD calculations is a 12-atom turn R-helix (Figure 6) in agreement with NMR and X-ray experiments4,5 The MD-predicted turns are anchored by HBs between the carbonyl oxygens closer to the N-terminal end and the nitrogens three residues toward the C-terminal end, leading to an opposite helix polarity to that of the ACHC containing compounds, again in agreement with experiment, and with previous calculations.4,5 The MD HB time profiles for one of the trajectories is shown in Figure 7. III. cis-ACPC Oligomer. The change of trans-ACPC to cisACPC (compound 3) in a peptide leads to the only oligomer that is not predicted to fold into a R-helix, among the compounds studied in this work. The most frequent fold observed in MD simulations at high dielectric constant is an extended strand with a zigzag shape (see Figure 8, upper panel). When simulations
9676
J. Phys. Chem. B, Vol. 113, No. 29, 2009
Figure 8. Snapshots of the MD trajectory for compound 3, for high dielectric constant simulation (upper panel) and low dielectric constant (lower panel). The simulation temperature was 300 K in both cases.
are performed at low dielectric constant, the peptide bends in a hairpin-like structure whose branches are close to the straight strands observed at the higher dielectric MD simulations (see Figure 8, lower panel). The results are very much in line with structural IR and NMR studies of these compounds, which predict zigzag like extended strands in solvents like MeOD.7 It is known that a highly polar solvent screens electrostatic interactions and breaks the HBs that may be formed in the hairpin-like structures.7 These differences in conformation that result from dielectric constant change can be readily seen in the peptide’s end to end distance (EED) time profile (see Supporting Information, Figure S2). As seen from this figure, the end-to-end distance in the simulation at higher dielectric constants reaches much higher values than those at low dielectric constant, predicting that these structures are more extended in a polar solvent. To quantify the conformational preponderance of the extended structure at high and low dielectric constants values, we have evaluated the amount of conformations that have an end to end distance above certain cutoff values. The results indicate that, at low dielectric constant, the conformations that have EED values higher than 12 Å are clearly a minority, reaching values that are not beyond 1.2%. At high dielectric constant the percentage of extended conformations at the same EED values, can reach up to 40%, clearly indicating that the high DC values favor elongated conformations. IV. Xylofuranoic Acid Based Peptides. The configuration around the CR-Cβ bond (discussed above) is not the only structural variable that has a large effect on the resulting fold. It has been observed that a five-membered ring (like the xylofuranoic acid) based peptide leads to a 14-atom helix, in spite of the presence of a cis CR-Cβ bond.22 Our simulated annealing MD results for the hexamer with an N-terminal t-Boc protected group (compound 4) predict that this peptide folds
Sussman et al. into a 14-atom turn helix with four HBs, in agreement with experiment. These results indicate that the ring conformation space is of the utmost importance in the peptide resulting fold and that it can be modulated by the ring composition and its pattern substitution. The time evolution of the distance for the four HBs in this structure is shown in Figure 9. As seen in this figure, the HB distance profiles are more variable than that of the bare cyclopentanic ring based peptides (compound 2). Moreover, as seen from this figure, the highest variability is observed for the C-terminal end, indicating that this peptide segment is frayed, in agreement with experiment.22 V. Oxetane-Containing Peptides. Ηelix folds with turns smaller than 12 atoms have not been as frequently observed.11 One exception is the fold afforded by peptides containing β-amino acids with a substituted four membered ring (oxetane), which were observed to form a 10 atom turn R-helix.8 We have sampled the conformational space of this kind of peptides with and without substitution motifs (compounds 5). The unsubstituted oxetane peptide (compound 5a) displays several possible conformations. The highest organized structure is a R-helix with 10 atom turns (see HB time evolution in Figure 10). Some other HB motifs present in these MD simulations are generated by 14-atom turns (results not shown). The intermingling of these two turn motifs could be due to their structural proximity in the Ramachandran plot, as discussed above for short transACHC peptides. The overall fold does not change when the rings are substituted with an OBn group (compound 5b), although the actual structural determination was carried out with a longer chain substituent (R ) CH2OBn). We have performed MCMD simulations on the latter compound as well. The results point to various helicoidal structures, but the predicted turns differ from the experimental values.8 It is worth noticing that the group that produced the experimental structure also performed a full torsional angle search powered by a genetic algorithm, and that the resulting predictions were at variance with the experimental results.8 We believe that the difficulties in predicting the observed fold for this peptide are related to the “roughness” of the folding landscape (as we will discuss below). B. Further Characterization of the Folding of Cyclic β-Amino-Acid-Based Peptides. Our molecular simulations offer the possibility of gaining insight into the interactions that produce stable folds for short β-amino acid peptides. One of the open questions in this field relates to the intrinsic ability of β-amino-acid-based peptides in general and cyclic β-amino acid peptides in particular to generate stable folds even for as short peptides as tetramers.6 On the contrary, it is well-known that only in rare occasions R-amino acid peptides fold with less than 10 amino acids in length.23 The folding rates for the formation of a R-helix fold could in principle be roughly gauged from the MD time evolution of the root-mean-square deviation (rmsd) from an experimentally based R-helix, starting from an extended structure. As representative examples of the peptides studied here, we have calculated this rmsd time evolution for the peptide backbone atoms in compounds 1e, 2, 4, and 5a. The results are shown in Figure 11. As seen from this figure, at the beginning, the structure of the molecule deviates by more than 4 Å from an “optimal” R-helix, and then it tends to reach a plateau in less than 100 ps, a value which can be construed as a preliminary estimate for the folding rates of the compounds studied here. These predicted folding rates are much faster than the ones observed for most R-amino acid peptides, which have been determined experimentally to be in the nanosecond range or
Cyclic β-Amino Acid Peptides
J. Phys. Chem. B, Vol. 113, No. 29, 2009 9677
Figure 9. Time evolution for the four HB distances that underlie the 14-atom turns predicted for compound 4. The predicted HB network is schematically drawn at the top. Notice the predicted fraying at the C-terminal end.
slower with few exceptions.23 Only very recently, there has been a report about the folding rate in the picosecond regime.24 No folding rates has been determined experimentally for acyclic or cyclic β-amino acid containing peptides, probably because it falls out of the time regimes attainable with the experimental protocols available until now. We hope that with the very recent advent of novel spectroscopic methods, this hurdle could be overcome.24 On the other hand, our protocols may overestimate the speed of the folding process, because the MD simulations where performed in the absence of a friction force provided by an explicit solvent, or by a Langevin-type molecular dynamics protocol. The plots shown in Figure 11 have corroborated the fold assignation based partly on the HB time evolution plots (Figures 1, 7, 9, and 10) and have provided additional information about the extent of conformational sampling and the departure from the optimal helical structure afforded by our MD protocol. Analysis of these figures indicates that most of the rmsd values fall below 2 Å and many fall below 1 Å away from the ideal structures. For instance, compound 1e has 99.9% of its conformations below 2 Å and 66.3% below 1 Å. The physical basis that underlies the folding process of a β-amino-acid-containing peptide is not as well understood as those of an R-amino acid peptide. To gain some insight into the interactions that favor an R-helix fold, we have calculated the components of the molecular mechanics total energy for the starting extended structure and a folded structure for some of those β-amino-acid-containing peptides that are predicted to fold into an R-helix (1e, 2, 4, and 5a). For the sake of
comparison, we also performed the same calculation for a R-Ala octamer. The folded structures chosen for our calculations were the ones that had the lowest rmsd with respect to an ideal R-helix. The differences in the molecular mechanics energy components between the extended and folded states are listed in Table 2. As seen from this table, the β-amino acid peptide folding into an R-helix structure is favored by the long-range nonbonded and electrostatic interactions (as in the case of R-amino acid peptides), an expected outcome since the folded molecule contain additional HB and van der Waals contacts. What differentiates most of the β-amino acid oligomers from the R-amino acid counterparts is that the localized internal energy components, like bond stretching, valence angle bending, and torsional angle distortions, favor the folded conformation. These results indicate that most β-amino acid peptides studied here have their resulting fold “coded” into some of their degrees of freedom, a property that could be the rationale behind the natural folding proclivity of these molecules. That is not the case for R-amino acid peptides: as seen from the results for the poly-R-Ala peptide shown in Table 2, the bond stretching and valence angle energy differences are positive for this peptide, a result that favors the extended structure. Interestingly enough, the oxetane-containing compound (5a) displays positive values for valence angle and torsion angle energy differences as well, raising the possibility that this compound does not have the strong intrinsic proclivity to folding found in the other β-amino acid peptides, and this feature could help explain the difficulties that molecular mechanics protocols face in the prediction of oxetane-containing amino acids (5) when R )
9678
J. Phys. Chem. B, Vol. 113, No. 29, 2009
Sussman et al.
Figure 10. Distance-time profile for the 10-atom turns in oxetane containing compound 5a. The predicted HB pattern is shown schematically at the top.
-CH2-OBn. In line with this rationale, the larger structural variability that can be observed in the HB time evolution plots for compound 5a (see Figure 10) could be the result of the lack of an intrinsic propensity to a helix fold. Moreover, the peptide that displays the most favorable local interactions (see Table 2) when folded as an R-helix in reference to an extended chain, (compound 1e), also displays the lowest structural variability as observed from its HB time evolution (see Figures 1 and 11). The reduced conformational landscape for this kind of peptides produce a less “frustrated” funnel, where the tendency to increasing entropy is countered by the intrinsic fold coded into the structure. Hence, the energy partition analysis of the kind shown in Table 2 may be of help in the prediction of peptides that produce the most stable folds. Nevertheless, a final gauge of the size of the entropic factors into the folding process would have to wait for simulations with an explicit solvent, since these are vital for estimating the hydrophobic effect, which is a solvent-driven entropic factor. Conclusions In this work, we sampled the conformational space of a group of diverse cyclic β-amino-acid-containing peptides that differ
in the length, ring size, and chemical structure as well as in the N- and C-terminal capping groups, through a set of MD runs that vary in the starting conformations and in the MD initial conditions. We choose this approach in order to enhance the amount of the conformational space sampled. In many cases our protocol was able to reproduce the secondary structure, in spite of a very simple implicit solvent representation. This outcome indicates that the β-peptide folding is not guided by the hydrophobic-driven solvent effect, like in R-peptides. Our energy component analysis (see Table 2) indicates that most cyclic β-amino acid peptides may be driven to a given fold by the propensity coded into the internal degrees of freedom like the valence and torsional angle components. The solvent affects the resulting structure through the screening of the electrostatic interactions, and hence it modulates the breaking or forming of hydrogen bonds like in the case (discussed above) of the cis-ACPC-based peptides (compound 3). Nevertheless, when the resulting peptide is highly structured, the solvent representation and other initial conditions do not seem to have a large effect on the resulting secondary structure. However, we have found out that in the case of charged peptides (e.g., compound 1a), a full atom explicit solvent representation
Cyclic β-Amino Acid Peptides
J. Phys. Chem. B, Vol. 113, No. 29, 2009 9679
Figure 11. Rmsd time evolution (for the backbone atoms) from an optimal R-helix, as obtained from MD simulations for compounds 1e, 2, 4, and 5a.
TABLE 2: Difference in Energy Components between the Extended and r-Helix Conformationsa energy component
(Ala)8
1e
2
4
5a
bond stretching valence angle torsion angle out-of-plane torsion nonbonded Coulombic
4.0 5.3 0.0 0.1 -28.7 -22.2
-4.0 -11.1 -16.0 -0.4 -33.2 -15.7
-1.3 -11.3 -6.1 -1.6 -90.5 -21.9
-0.4 -10.1 -7.2 -0.2 -15.1 -18.8
-0.7 2.1 2.8 0.2 -23.9 -12.4
a
All energy values are in kcal mol-1.
is needed to avoid spurious interactions and to reproduce the experimental peptide fold. Comparisons of the time evolution profiles for the hydrogen bonds that define a turn motif for those compounds that fold into a R-helix indicate that some peptides have more conformational variability than others. For instance, the HB time evolution for some of the ACHC-containing peptides indicate that the distances that define the HB turn vary to a smaller extent than in those peptides that include a heteroatom in the cycle like oxetane-containing peptides (compounds 5). These results are also borne out by the analysis of the rmsd time profile from
an optimal structure, which clearly indicate the conformational variability predicted by our algorithm for every peptide studied. These plots also confirm that our MD-based protocol produce conformations that are structurally close to the observed structures. Our energy component breakdown for these molecules indicates that the internal degree of freedom energy components favor the R-helix over an extended structure. This outcome supports the notion that for many β-amino acid peptides, the folding “frustration” is reduced by this internally coded energy propensity, favoring the folding of even short β-amino acid peptides. For this reason, we believe the analysis based on the MM energy breakdown could have a predictive value. Nevertheless, it is possible that in this case the solvent may play an additional role by either providing a frictional force that could help to keep their fold intact, and our calculations cannot rule out that the solvent could additionally help to compact the structure through a hydrophobic like force. As a whole, the results presented here indicate that the multiple MD approach introduced therein is able to reproduce
9680
J. Phys. Chem. B, Vol. 113, No. 29, 2009
in many cases the folding structure of these peptides and to shed some light into the physical forces that lead to a given fold. The protocols used here could be regarded as a useful and affordable first tier approach for the screening and/or design of cyclic β-amino acid peptide libraries as a step prior to their synthesis and structural characterization. Acknowledgment. This work was supported by a Ministerio de Educacio´n y Ciencia (MEC) grant to R.J.E. and Xunta de Galicia financial aid to F.S. We also thank the Supercomputing Center of Galicia (CESGA) for computer time. Supporting Information Available: This section contains two figures (S1 and S2). The first figure describes the time evolution for the possible 10-atom turns in peptides 1c (upper panel) and 1b (lower panel) from a MD trajectory, while the second figure describes the time evolution for the distance between the N-terminal nitrogen and C-terminal carbon in compound 3, resulting from a MD trajectory at low (black) and high dielectric constants (red). This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Gellman, S. H. Acc. Chem. Res. 1998, 31, 173–180. (2) Hill, D. J.; Mio, M. J.; Prince, R. B.; Hughes, T. S.; Moore, J. S. Chem. ReV. 2001, 101, 3893–4011. (3) (a) Seebach, D.; Beck, A. K.; Bierbaum, D. J. Chem. BiodiVersity 2004, 1, 1111–1239. (b) Seebach, D.; Abele, S.; Gademann, K.; Guichard, G.; Hintermann, T.; Jaun, B.; Matthews, J. L.; Schreiber, J. V. HelV. Chim. Acta 1998, 81, 932–982. (4) Barchi, J. J., Jr.; Huang, X.; Appella, D. H.; Christianson, L. A.; Durell, S. R.; Gellman, S. H. J. Am. Chem. Soc. 2000, 122, 2711–2718. (5) Apella, D. H.; Christianson, L. A.; Karle, I. L.; Powell, D. R.; Gellman, S. H. J. Am. Chem. Soc. 1996, 118, 13071–13072. (6) Hete´nyi, A.; Ma´ndity, I. M.; Martinek, T. A.; To´th, G. K.; Fu¨lo¨p, F. J. Am. Chem. Soc. 2005, 127, 547–553. (7) Martinek, T. A.; To´th, G. K.; Vass, E.; Hollo´si, M.; Fu¨lo¨p, F. Angew. Chem., Int. Ed. 2002, 41, 1718–1721.
Sussman et al. (8) Claridge, T. D. W.; Goodman, J. M.; Moreno, A.; Angus, D.; Barker, S. F.; Taillefumier, C.; Watterson, M. P.; Fleet, G. W. J. Tetrahedron Lett. 2001, 42, 4251–4255. (9) (a) Alema´n, C.; Zanuy, D.; Jime´nez, A. I.; Cativiela, C.; Haspel, N.; Zheng, J.; Casanovas, J.; Wolfson, H.; Nussinov, R. Phys. Biol. 2006, 3, S54–S62. (b) Zanuy, D.; Nussinov, R.; Alema´n, C. Phys. Biol. 2006, 3, S80-S90. (c) Clark, T. D.; Buehler, L. K.; Ghadiri, M. R. J. Am. Chem. Soc. 1998, 120, 651–656. (10) Goodman, J. L.; Petersson, E. J.; Daniels, D. S.; Qiu, J. X.; Schepartz, A. J. Am. Chem. Soc. 2007, 129, 14746–14751. (11) Beke, T.; Somlai, C.; Perczel, A. J. Comput. Chem. 2006, 27, 20– 38. (12) Rathore, N.; Gellman, S. H.; de Pablo, J. J. Biophys. J. 2006, 91, 3425–3435. (13) Zhu, X.; Koenig, P.; Gellman, S. H.; Yethiraj, A.; Cui, Q. J. Phys. Chem. B 2008, 112, 5439–5448. (14) (a) Boned, R.; van Gunsteren, W. F.; Daura, X. Chem.sEur. J. 2008, 14, 5039–5046. (b) Daura, X.; van Gunsteren, W. F.; Mark, A. E. Proteins: Struct. Funct., Genet. 1999, 34, 269–280. (c) Daura, X.; Jaun, B.; Seebach, D.; van Gunsteren, W. F.; Mark, A. E. J. Mol. Biol. 1998, 280, 925–932. (15) Soengas, R. G.; Este´vez, J. C.; Este´vez, R. J. Org. Lett. 2003, 5, 1423–1425. (16) Dauber-Osguthorpe, P.; Roberts, V. A.; Osguthorpe, D. J.; Wolff, J.; Genest, M.; Hagler, A. T. Proteins: Struct., Funct., Genet. 1988, 4, 31– 47. (17) InsightII (version 2000), Biopolymer and Discover are trademarks from Accelrys Inc. San Diego, CA. (18) (a) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187–217. (b) CHARMM 34 was provided by Harvard University, Cambridge, MA. (19) Seebach, D.; Schreiber, J. V.; Abele, S.; Daura, X.; van Gunsteren, F. W. HelV. Chim. Acta 2000, 83, 34–57. (20) Gee, P. J.; van Gunsteren, W. F. Proteins: Struct., Funct., Bioinf. 2006, 63, 136–143. (21) Okur, A.; Wickstrom, L.; Simmerling, C. J. Chem. Theory Comput. 2008, 4, 488–498. (22) Chandrasekhar, S.; Reddy, M. S.; Jagadeesh, B.; Prabhakar, A.; Ramana Rao, M. H. V.; Jagannadh, B. J. Am. Chem. Soc. 2004, 126, 13586– 13587. (23) Thompson, P. A.; Eaton, W. A.; Hofrichter, J. Biochemistry 1997, 36, 9200–9210. (24) Mohammed, O. F.; Jas, G. S.; Lin, M. M.; Zewail, A. H. Angew. Chem., Int. Ed. 2009, 48, DOI: 10.1002/anie.200900513.
JP811321N