J. Phys. Chem. B 1999, 103, 2535-2542
2535
A Study on Local-Global Cooperativity in Protein Collapse James J. Chou and Eugene I. Shakhnovich* Department of Chemistry and Chemical Biology, HarVard UniVersity, 12 Oxford Street, Cambridge, Massachusetts 02138 ReceiVed: September 30, 1998; In Final Form: January 26, 1999
We study the effect of secondary structure formation on cooperativity of protein collapse using a simple homopolymer model on a “210” lattice that allows a geometrically consistent representation of R-helical conformation. Monte carlo simulations were carried out in the range of temperatures and energetic parameters that characterize relative strength of global (i.e., between monomers that are far apart in sequence) and local (helix-stabilizing) contacts R ) ηg/ηl. The complete phase diagram presented here exhibits a narrow region of temperatures and R-values in which compact conformations with considerable helical content dominate the equilibrium. It was shown that at certain values of R, the collapse transition occurs in a narrow temperature range and is accompanied by a pronounced increase of helical content. However, the simulations do not support a two-state transition scenario at which transition occurs between two (metha) stable states corresponding to free energy minima. Rather, we observe at all temperatures and all values of R a single stable state that evolves, as temperature gets lower, toward more compact conformations. We argue that additional factors such as sequence specificity and/or side-chain packing should be taken into account to explain the two-state character of folding transition observed in real proteins.
1. Introduction The thermodynamic and kinetics aspects of the protein folding problem remain to be a vibrant area of research. Considerable progress has been achieved in understanding the basic physics of folding (reviewed in refs 1-3). Many theoretical studies in protein folding have been focused on the analysis of simple models for which complete folding simulations and analysis are possible. These studies proved useful in elucidating the origin of protein cooperativity in thermodynamics and its implications for kinetics. Simplicity of these models and their ability to capture essential physics of folding (Levinthal paradox, cooperativity) makes them very appealing. Consequently, many important results were derived from their study. However, simplifications inherent to manysbut not alls lattice models have been criticized recently by Honig and Cohen.4 These authors argued that simplified (side chain only in their terminology) models do not take into account the protein backbone structure and thus have a potential to miss certain important features of protein folding, in particular the role of secondary structure in thermodynamics and kinetics of folding. In the response note, Shakhnovich5 argued that lattice models are useful to focus our thinking about protein folding on fundamental issues, since they reproduce several universal, invariant features of this phenomenon. It is clear that better understanding of protein folding and the origin of protein stability must be gained by studying models that incorporate more details of protein structure. In particular, the origin of protein cooperativity should be addressed. It was argued that a likely reason for folding cooperativity in heteropolymer protein-like models is due to special selection of sequences that provide pronounced energy minimum for the * Author to whom correspondence should be addressed. Dept. of Chemistry Harvard, University, 12 Oxford Street, Cambridge, MA 02138. Phone: (617) 495-4130. Fax: (617) 496-5948. E-mail: eugene@diamond. harvard.edu.
native conformation.6-10 However, there may be other reasons of cooperativity in macromolecules which may be essential for proteins too. An important example is coil-to-globule transition in stiff homopolymers which was shown both in theory and simulations11-16 to be a cooperative one. Cooperative coil-to-globule transitions in homopolymers should be distinguished from cooperative folding transitions in protein-like heteropolymers. In the latter case, transition occurs between the multitude of unfolded conformations and unique native structure, while in the former case of stiff homopolymers a possible cooperative transition may occur between the coil state and compact globular state which may not have a unique structure but rather represents a thermodynamic state where a multitude of compact conformations dominates. Apparently both mechanisms may contribute to cooperativity of experimentally observed folding transition.17 In this study we address the issue of how formation of secondary structure (neglected in some earlier lattice simulation studies) can influence the character of collapse transition in polypeptides. This issue is important since secondary structures such as R-helices and β-sheets are present in most folded polypeptides. They are important members of the protein hierarchic structure, which is likely to be essential for most proteins to be stable and function properly. From the physical viewpoint, secondary structure is simply a short-range order partly stabilized by favorable interactions of nearby polymer links. Formation of stable secondary structures must also satisfy local geometrical constraints of biopolymers. In the case of polypeptides, the stereochemistry permits the existence of only few stable short-range orders such as helices, sheets, and turns.18 Therefore, in the absence of effects from the rest of the chain, secondary structure can be determined by local geometrical constraints and interactions among nearby elements on the chain. However, secondary structure formation during the folding of the entire polypeptide chain becomes much more complicated.
10.1021/jp9839192 CCC: $18.00 © 1999 American Chemical Society Published on Web 03/12/1999
2536 J. Phys. Chem. B, Vol. 103, No. 13, 1999 In this case, tertiary interactions, or the interactions between elements located far from one another along the chain, play important roles in determining the folding outcome. Their effects on secondary structure formation can be either positive or negative. Hence local geometrical and energetic properties of a polypeptide do not determine the formation of secondary structure; they only introduce local conformational propensity to the polypeptide. Since the folding outcome is a compromise between volume interactions and local conformational propensity, it is important to explore the relation between these two factors and thereby obtain a better understanding of the globallocal cooperativity in protein folding. A number of studies addressed the issue of the relative role of local and tertiary interactions in protein folding. These studies mostly centered around the discussion of the relative importance of local and tertiary interactions both in the folding process and in stabilization of the native state. However, the results from these studies are somewhat contradictory. Using the results from 2D-lattice simulation, Dill and co-workers19,20 argued that helical propensities may be weak determinants of globular protein structures. Govindarajan and Goldstein21 also reported that the optimal conditions for folding are achieved when the contributions from local interactions to the stability of the native state are small. Experimental studies suggest that there must be an efficient balance between local and nonlocal interactions to provide efficient folding.22 In a related development, Zwanzig and co-workers23 found, on the basis of mathematical analysis of a simple model, that pronounced local bias toward the correct conformation can greatly reduce the folding time compared to the Levinthal time of exhaustive search. A similar point was made by Luthey-Schulten and co-workers,24 who studied analytically a phenomenological statistical-mechanical model of a heteropolymer where formation of secondary structure was taken into account within a Bragg and Zimm approximation along with simplified treatment of the effect of confinement in the globular state, excluded volume and geometrical packing of the helices. Using this model Luthey-Schulten and co-workers showed that formation of helices may substantially decrease the total number of conformations making exhaustive search feasible for small proteins. The differences between the views of different authors may have resulted from significant differences among the models used in these studies. For instance, the definition of local propensity and local order would have profound effects on folding. Since, on average, secondary structures form a large portion of a protein, it is clear that both local and tertiary interactions are essential in determining the native state. To estimate the contribution from local interactions in protein folding, a simulation model with a more detailed implementation of protein backbone stereochemistry may be necessary. Hence, instead of focusing on the relative importance of local and tertiary interactions, it is perhaps more informative to explore the intrinsic properties associated with local and tertiary interactions in polymer folding. A simple and relevant model for studying the aspect of cooperativity related to the chain collapse and the role of local and nonlocal interactions is polyamino acids. Several of the ordered structures found in proteins, such as the R-helix and the pleated sheet or the β-structure, have been detected in polyamino acids.18 In fact, polyamino acid such as polyalanine, though having only one type of side-chain, can exhibit some very interesting folding phenomenon. The stereochemistry and thermodynamics of polyalanine favor R-helix. Moreover, since the alanine side-chains are nonpolar, the R-helices can be
Chou and Shakhnovich stabilized by β1-R4 hydrophobic bonds in straight helical segments. However, nonpolar homopolymers in water tend to collapse to a compact form, which disrupts the formation of any extended helices. Clearly, compromise must be made between the two forces if any stable structures should form. Alternatively, it is possible that tertiary interactions between nonpolar elements enhance helix formation. Indeed, a Monte Carlo simulation study on homopolypeptides has reported that the secondary structure seen in the folded state is the result of an interplay between the short- and long-range interactions.13 Experimental investigation on the conformation of a homopolypeptides, carried out by Scheraga and co-workers,25 further suggested the notion of “cooperativity yield stability”. In their work, Scheraga and co-workers used optical rotatory dispersion (ORD) measurement to observe R-helix formation in poly-L-alanine, a polyamino acid known to exist in the righthanded R-helical form in the solid state, in solution in nonaqueous solvents, and in water as a block copolymer. Since nonpolar polyamino acids are insoluble in water, the authors incorporated the nonpolar poly-L-alanine chain as a central block polymer between two blocks of water-soluble poly-D,L-lysine to increase the solubility. Note that the two lysine blocks repel against each other in water, causing the alanine portion to be extended. According to their optical rotation measurements, the helix-coil transition temperature of the polyalanine portion is more than 25 °C higher in 0.2 M NaCl than it is in pure water. In addition, the same effect can be obtained by raising the pH from 6 to 12. Since the electrostatic repulsion between the lysine end blocks is greatly reduced under these conditions (0.2 M NaCl and pH 12), the finding strongly indicates that in the absence of “stretching” by the lysine blocks, poly- L-alanine acquires more stable conformations, resulting from hydrophobic interactions among helical segments. The lattice model has been a computationally convenient platform for studying conformational statistics of polymers and protein folding. While it is difficult to model the secondary structure of a protein using a simple cubic lattice model, more detailed lattice models were suggested, most notable of them is the so-called 210 lattice model.13 In fact, the 210 lattice model had been used with much success in both simulating secondary structure formation13 and studying the helix-coil transition of polypeptides.13 In this paper we use the 210 lattice to study the effect of formation of secondary structure on cooperativity of collapse transition. To avoid confusion with heteropolymer cooperativity, we study a homopolymer model with the understanding that collapse transition in this case is not a true folding transition to the unique ground state but rather to the ensemble of compact conformations. In the current study, polyamino acids are modeled using a 40-mer homopolymer with local conformational propensity for helix. In this paper, we refer the terms helix to R -helix and contact to a tertiary interaction for conciseness, unless specified otherwise. The implementation of local geometrical and energetic details of a real polypeptide, which favors helix formation, was accomplished by introducing a pseudodihedral potential throughout the polymer chain. Dynamic Monte Carlo simulations of the 40-mer represented in the 210 lattice were performed at various conditions, and the average values of order parameters such as helicity and number of contacts were recorded. By varying the relative contribution of global and local energy to the total energy of the system, a sharp transition in folding characteristics was identified. At this transition, simulation yielded compact conformations in which helices are stabilized
Local-Global Cooperativity in Protein Collapse by contacts. While the above transition region, in which local order and tertiary contacts coexist, seems to suggest cooperativity, we have shown on the basis of melting curves and distribution of order parameters that for homopolymer, folding transition in this region is noncooperative. Our results are summarized in the phase diagram that outlines the conditions at which different forms of the polymer chain exist. 2. Model 210 Lattice. The 210 lattice was designed to account for the most important conformational features of polypeptides, but at the same time, keeping the computational simplicity of a lattice model.13 In this model, polymer chains are still built with simple cubic lattice points. However the adjacent beads of a polymer are connected by vectors of the type ((2, (1, 0), rather than the type ((1, 0, 0) of cubic lattice model. Since the ((2, (1, 0) vector has 24 orientations, compared to 6 of the ((1, 0, 0) vector, polymers built in the 210 lattice have much more conformational flexibility. It has been shown that this model is able to fit with reasonable accuracy the local and global geometry of polypeptides.13 In the 210 lattice, the lattice distance between two consecutive beads, or the bond length, is equal to x5 (in units of simple cubic lattice). If the unit length of the simple cubic lattice is equal to 1.7 Å, the bond length becomes 3.785 Å, the distance between two consecutive R-carbons in real protein. To account for the excluded volume of the polypeptide backbone, each bead is taken to be a sphere of radius 1. This way, the six lattice points around the bead are no longer accessible and thus a backbone envelope is conveniently introduced. In the 210 lattice, there are 10 possible values of the bond angle θ: 36.9°, 53.1°, 66.4°, 78.5°, 90°, 101.5°, 113.6°, 126.8°, 141.1°, and 180°. However, θ values of 36.9°, 53.1°, 66.4°, and 180° are forbidden since they are not possible in real protein. Therefore, other than the 210 lattice itself, the two other major limitations to conformational flexibility of the polymer are excluded volume and the bond angle constraints. Pseudodihedrals. The pseudodihedral was introduced as a simplified representation to the traditionally used dihedral angle ψ. This representation is most efficiently used with the simple polymer model such as a chain of R-carbons. The precise definition of pseudodihedral angle is as follows. Consider four consecutive R-carbon atoms, CRi , connected by three vectors Vi ) Ri+1 - Ri. The pseudodihedral angle Φi+1,i+2 is the angle between the normal vectors N1 ) ||V1 × V2|| and N2 ) ||V2 × V3||. The sign of the pseudodihedral is defined to be the opposite of the projection of V3 on N1. Though pseudodihedral alone is still not sufficient to define the conformation of the chain, its simple vector representation using four R-carbons can be conveniently handled in the 210 lattice. Energy. For the study of local-global cooperativity, a hydrophobic homopolymer with local conformational propensity toward R-helix is used. This homopolymer is implemented in the 210 lattice as a chain of beads, in which the beads can be considered as the R-carbons of a polypeptide. To mimic the effect of hydrophobicity, any two beads more than two bonds away are subjected to hard-core attractive potential. In folding simulation on the lattice, hard-core potential is the standard potential used to mimic the effect of electrostatic and hydrophobic interactions among particles. This potential is crucial for the formation and stabilization of a tertiary fold. Although R-carbon trace of a polypeptide may be conveniently used to describe the general backbone conformation, it cannot give local structural details such as the orientations of the rigid planes
J. Phys. Chem. B, Vol. 103, No. 13, 1999 2537 formed by atoms CR, C′, N, and CR. More importantly the information on interactions between side-chains of adjacent residues is also missing. These details are indispensable in characterizing the local backbone potential of the polypeptide. A simple and efficient way to introduce local conformational propensity to a chain of beads is to apply pseudodihedral potential to all pairs of adjacent beads. A typical value for the pseudodihedral angle of an R-helix is 45°.27 In the 210 lattice, the value of Φ closest to 45° is 48°. It must be mentioned that pseudodihedral potential alone is not enough to specifically favor the formation of R-helix, which must also be stabilized by (βi, βi+3) contacts. Hence, local conformational propensity is implemented by applying both the pseudodihedral potential which favors R-helix formation to all successive pairs of beads (βi+1, βi+2), and hard-core attractive potential between βi and βi+3. In summary, the total energy of this polymer may be separated into two parts. The first part is the sum of all tertiary contact energies between beads separated by at least four bonds. It is named as the global energy, Eg, and defined as N
Eg ) where
{
∞
Ecij
) -ηg 0
∑ Ecij i)j-3
(1)
r2ij e 2 2 < r2ij < 6, |i - j| g 4 otherwise
(2)
In the above equation, the negative sign in front of ηg indicates hydrophobic attraction. The second part consists of pseudodihedral energy attributed to R -helix conformation and contact energy between βi and βi+3; it is named as the local energy, El, and defined as N-4
El )
N
Epi +1,i+2(Φ) + ∑ Ehij ∑ i)0 i)j-3
(3)
where
{
l 45°< Φ < 50° Epi +1,i+2(Φ) ) -η /2 0 otherwise
{
∞
Ehij ) -ηl/2 0
r2ij e 2 2 < r2ij < 6, |i - j| ) 3 otherwise
(4)
(5)
Note the above local energy function gives a strong preference for R -helix formation. The (βi, βi+3) contact energy is included in the local energy since it is a crucial factor in local R-helix formation. Finally the total energy of a chain of N beads, labeled β0, ..., βN-1, is given by
Etot ) Eg + El
(6)
To study the relative importance of local and nonlocal interactions for the features of collapse transition, the simulations were performed with various ratios of ηg to ηl. The ratio ηg/ηl was varied in such a way that ηg + ηl ) 10, the condition that keeps total interaction energy approximately constant. Simulation Parameters. For evaluating volume interaction and compactness of the homopolymer, the number of pairwise contacts, NC, is used as a simulation parameter. In this study, a pairwise contact is defined much in consistency with the
2538 J. Phys. Chem. B, Vol. 103, No. 13, 1999
Chou and Shakhnovich
definition of contact energy in eq 2. A pair of monomers are considered to be in contact if 2 < rij < 6 and |i - j| g 4, where rij is the lattice distance between the beads βi and βj. Here the second condition, |i - j| g 4, is applied to discount local contacts that are not informative about the volume interactions, and exclude particularly the (βi, βi+3) contacts that can appear along with the formation of helices. To measure the degree of helix formation, the parameter named helicity is introduced. As mentioned previously, in the 210 lattice a segment of helix is formed only if all pseudodihedral angles within the segment take on the value of 48°. The helix length m is defined to be the total number of beads included in the helical segment. The helicity, Ξ, of a conformation is thus defined as
Ξ)
∑i m2i /N2
(7)
where mi is the length of the i th helical segment, and N is the length of the entire chain. The advantage of such definition is that a complete helical conformation, or the entire chain is a rigid helix with no interrupts, would yield a helicity of 1. Dynamics. Folding simulation involves a set of dynamic moves performed on a polymer of N beads, labeled successively as β0, ..., βN-1. For all movements, the bead on which the move is to be performed is chosen at random and the move cannot violate excluded volume or bond angle constraints. While the moves are restricted by backbone constraints, they must also not violate the topological constraint that a segment of chain cannot cut through another segment (or bond cutting). This topological constraint would greatly reduce the chance of unrealistic knot formation. Three types of moves are used in all simulations. (1) Single bead flips, or “spike” motions.28 For a detailed description of this type of move in the 210 lattice, see the paper by Skolnick and Kolinski.29 This move does not violate the topological constraint because the thickness of the backbone hard-core envelope is large enough to prevent any bond cutting. (2) Two bead end flips, in which two end bonds are randomly transformed to a new set of vectors.29 It has been shown that moves (1) and (2) are sufficient to satisfy the correct dynamics for the athermal random coil state in the absence of hydrodynamics interactions,30 i.e., they obey Rouse-like dynamics with excluded volume.13 Moreover, they even give the correct character of local conformational relaxation as probed by orientational correlation functions. However, due to the nature of the 210 lattice and backbone constraints, some simple move such as the flipping of a square loop cannot be accomplished by any simple combination of single bead flips unless the loop is stretched first. Hence, in order to get around such local minima faster, chain rotations are also applied. This considerably shortens the equilibration time for a simulation. (3) Chain rotation. The first bead of a segment (βi f βi+L+1) is randomly chosen for rotation, where 0 e L e N - i is a random integer and L is defined to be the rotation segment length. Note the L + 2 beads on this segment are connected by vectors Vi, ..., Vi+L. For computational simplicity, only Vi and Vi+L are used and combined to perform single bead flips. The rest of the vectors in the segment remain unchanged. Hence, in general, this is not a symmetrical rotation with the exception of L ) 2 (the three-bond loop). Unlike the single bead flip, chain rotation involves moves of much larger scale and is thus capable of bond cutting. Hence during simulation, any rotation that leads to bond cutting is rejected. To statistically favor the rotation of shorter segment, an exponential distribution function e-ξL is introduced. With this function, both move (1) and (3)
Figure 1. The average values of NC and helicity are plotted for simulations at various ηg/ηl, the ratio of tertiary interaction energy to local (helical) interaction energy. In the plot, open circles represent 〈helicity〉 and filled circle represent 〈NC〉. The three regions I, II, and III are separated by two dashed lines. All simulations are performed at T ) 1.2.
may be considered as chain rotation with the special case of L)1, which then becomes the single bead flip. Clearly the parameter ξ determines the mobility of the chain. Smaller value of ξ leads to larger scale chain rotations. In this study, ξ ) 0.03 since this value gives reasonably efficient collapsing of a 40-mer hydrophobic chain. A single simulation move consists of the following steps: • A bead βi is randomly selected for movement. • If i ) 0, or N - 1, perform move (2). • If i * 0, or N - 1, randomly select a segment length L, obeying the condition 0 e L e N - i, for chain rotation. If L ) 1, perform move (1). Otherwise, generate a random real number 0 < R < 1. If R < eξL, perform move (3), else regenerate R and repeat this process until the condition is satisfied and move (3) is performed. After each simulation move, the energy of the new conformation is calculated and compared to the energy of the previous conformation, and a standard asymmetric Metropolis sampling scheme is employed. In this study, the collapse of a 40-mer is simulated at various temperatures and ηg/ηl ratios. Since equilibration time required by systems of different parameters may vary significantly, all simulations are performed far beyond the initial equilibrium point to ensure validity of the results. For each simulation, the expectation value of NC and helicity were calculated. 3. Results To examine the effect of relative contribution of local and global energy on the formation of helices and contacts, Monte Carlo simulations were performed at T ) 1.2 (a relatively low temperature) with various ηg/ηl ratios ranging from 0.2/9.8 to 3.0/7.0. The plot of 〈Ξ〉 and 〈NC〉 vs ηg/ηl ratio is shown in Figure 1. For visual interpretation, characteristic conformation traces for each simulation were also generated and shown in Figure 2. The simulation starts at ηg/ηl ) 0.2/9.8, a ratio that gives a dominant local energy term and almost no global energy. As expected, the folding outcome is a rigid, extended helix (Figure 2a). Interestingly enough, 〈Ξ〉 remains to be large as ηg increases from 0.2 to 1.4, and then drops rapidly in a narrow transition region around ηg/ηl ) 2.0/8.0. In this region (region II in Figure 1), a high degree of local-global cooperativity is implied by the coexistence of significant amounts of helices and contacts. A characteristic conformation at ηg/ηl ) 2.0/8.0
Local-Global Cooperativity in Protein Collapse
J. Phys. Chem. B, Vol. 103, No. 13, 1999 2539
Figure 2. Characteristic conformations at simulation equilibrium. (a) The characteristic conformation at three different stages of melting for the case ηg/ηl ) 0.2/9.8, in which local interactions dominate. (b) Cooperativity is observed for the case ηg/ηl ) 2.0/8.0, in which three helix bundle is formed at T ) 1.1. (c) In the case of ηg/ηl ) 2.6/7.4, contact formation dominates, and the characteristic conformation at T ) 1.2 is typical for the collapse of homopolymer without local conformational propensities.
Figure 3. Phase diagram of ηg/ηl ratio vs T showing four different phases. The HC phase is helical and compact, and is defined such that 〈NC〉 g 14 and 〈Ξ〉 g 0.14. The letter C represents the compact phase in which 〈NC〉 g 45 and 〈Ξ〉 < 0.08. The letter H represents the helical phase in which 〈Ξ〉 g 0.2 and 〈NC〉 < 9. NN represents the nonhelical and noncompact phase in which none of the above three criteria is satisfied.
(Figure 2b) shows a three-helix bundle, in which each helix is stabilized by contacts with the other two helices. As ηg/ηl increases further to 2.4/7.6, helices vanish completely and contacts become dominant. The characteristic conformation in Figure 2c shows that the 40-mer is collapsed into a compact form with about 65 contacts. In addition to Figure 1, a more complete 2D phase diagram (shown in Figure 3) was constructed by varying both the ηg/ηl ratio and temperature. The phase diagram exhibits a very narrow region (HC) corresponding to compact helical globules. As can be seen in Figure 2, at other conditions stable conformations
Figure 4. The melting curves for the three cases: ηg/ηl ) 0.2/9.8 (represented by open diamond), ηg/ηl ) 2.0/8.0 (represented by open circles), and ηg/ηl ) 2.6/7.4 (represented by filled circles). (a) The melting curves in terms of total energy. (b) The melting curves in terms of NC. (c) The melting curves in terms of helicity.
correspond to either compact states with little or no helicity or to helical but noncompact conformations or even random coils. From Figure 1, three ηg/ηl ratios that result in three cases of distinct folding characteristics were selected for further simulation studies. These three cases are (I) ηg/ηl ) 0.2/9.8, dominant helix formation; (II) ηg/ηl ) 2.0/8.0, coexistence of helix and contacts; and (III) ηg/ηl ) 2.6/7.4, dominant contacts. Melting curves in terms of 〈Etot〉, 〈Ξ〉, and 〈NC〉 for above three cases were generated and shown in Figure 4. The 〈Etot〉 vs temperature curves for all three cases are similar and have similar slope of transition. Figure 4, parts b and c, provides a more detailed look at the melting curves in terms of 〈Ξ〉 and 〈NC〉. From these two
2540 J. Phys. Chem. B, Vol. 103, No. 13, 1999 plots, it is clear that the melting curve in case I corresponds to the melting of helices, and similarly, the melting curve in case III corresponds to the disruption of contacts. For case II, the melting curves in terms of 〈Ξ〉 and 〈NC〉 illustrate a more complicated melting transition. In Figure 3, parts b and c, initially there is a significant drop in 〈NC〉 at T ) 1.3, while 〈Ξ〉 remains much the same or even increased by a little. Then 〈NC〉 levels off at T ) 1.5, whereas 〈Ξ〉 decreases rapidly. This suggests that the initial stage of melting involves the disruption of local-global cooperativity. The characteristic conformations on this melting curve (Figure 2b) show also that the initial melting of the three-helix bundle involves the separation of one helix from the other two by disrupting the contacts between the helices. Then the extended helices are gradually melted as temperature increases. In terms of total energy, the slopes of melting transition in the three cases are similar. To further explore the nature of the melting transitions in the above three cases, frequency distribution of 〈Ξ〉, 〈NC〉, and 〈Etot〉 at various temperatures were generated (Figure 5). For case I, the distribution of Ξ (Figure 5a) at low temperature (T ) 1.2) shows a narrow unimodal peak at high helicity (Ξ ) 0.98), indicating a rigid, extended helical conformation. As the polymer is denatured by increasing temperature, the distribution of Ξ becomes wider and the maximum moves to a lower value of Ξ. At transition (T ) 2.0), the 〈Etot〉 distribution is symmetric and unimodal, indicating that the transition is noncooperative and diffusive. In case II (Figure 5b), distributions of 〈NC〉 and 〈Ξ〉 at low temperature (T ) 1.1) are clearly bimodal. However, the distribution of 〈Etot〉 is unimodal. This suggests the presence of multiple energetically similar conformations at low temperature, which is consistence with the expectation that homopolymers lack an unique native state. At transition (T ) 1.6), 〈Etot〉 distribution is also symmetrical and unimodal, indicating that the transition is noncooperative. For case III, 〈Etot〉 distribution at transition (T ) 1.9) is slightly bimodal (Figure 5c), suggesting weak cooperativity. 4. Discussion As the local-global contribution ratio changes, three regions with distinct collapse characteristics are identified (see Figure 1). Among them, the most intriguing case is region II, a narrow transition region in which folding simulations yielded conformations with significant amount of both helices and contacts. Around this region, increasing the helical propensity resulted in rapid decrease in global contacts and increase in helicity. Similarly, increasing the magnitude of nonlocal interactions lead to rapid increase in the number of contacts and disappearance of local order. This finding suggests that the collapse transition in the homopolymer chain with secondary structure features a certain degree of cooperativity as manifested in the sigmoidal shape of some of the transition curves shown in Figure 1. Our simulations clearly suggest that the collapse transition occurs in a narrow temperature interval only for a certain optimal ratio of energetic contributions from local and nonlocal interactions. Further, our simulations show that the range of this ratio in which a significant amount of cooperativity exists is very narrow, and cooperativity reduces rapidly as small changes (in either way) are made to the relative contribution from local and tertiary interactions. This is in line with earlier theoretical predictions by Grosberg.31 However, a deeper look at the distributions of different order parameters, such as total number of contacts, degree of helix formation and, finally energy, (Figure 5) shows that in all cases the collapse transition is not really two-state (first-order-like).
Chou and Shakhnovich Indeed the distribution of energies is always monomodal for all studied temperatures and ratios of energies of local to nonlocal contacts. In this sense the collapse transition that is observed for this model is very different from the cooperative folding transition observed in proteins. Clearly other factors such as sequence heterogeneity, evolutionary selection, and specific side-chain packing may play an important role in determining cooperative first-order-like transition in proteins.6,7,32-35 It is interesting to note that for a certain ratio of local vs nonlocal contact energies the distribution of the number of contacts and degree of helix formation are bimodal at lowenough temperature (see Figure 5b). However distribution of energy at this temperature is clearly monomodal and is shifted to the low energy domain. This clearly shows that the reason for the bimodality in the distribution of the order parameter is that several compact conformations having somewhat different degrees of helicity are equally probable at low temperature. This observation is consistent with the fact that homopolymers may not have unique structure. To connect our results more directly to real proteins we estimate what ratio of local to nonlocal contact energy corresponds to any of the real polypeptides. Since the details of interactions in lattice model varies significantly from that in real polypeptides, it is difficult to directly compare its energy with the experimental values for real proteins. However, a rather general comparison can be made on the basis of relative conformational propensities. Our goal is to estimate the parameter ηg/ηl for polypeptides. It turns out that for polyalanine this parameter is close to 0.2. Indeed ηg can be evaluated from MJ parameters,36 while ηl can be evaluated from the R -helical propensity of polyalanine ηl ) kT log sR where sR is a Zimm-Bragg parameter for helix-coil transition.37,18 The result is close to 0.2. However, an essential caveat of this estimate is that absolute values of interactions cannot be derived unambiguously from statistics in quasichemical approximation because the energetic factor (“effective temperature”) is not known,38 therefore this estimate should be taken with caution. There has been some controversy in the literature that stemmed from the argument by Chan and Dill39-41 that compactness alone is a significant driving force for the formation of secondary structure in globular proteins. On the basis of exhaustive searches of conformational states on 2D and some 3D cubic lattice, Dill and co-workers reported significant increase of secondary structure when the polymers are closely packed. In testing this conjecture, Gregoret and Cohen42 carried out simulations using a nonlattice model and found that upon packing, in a density range from as low as of a polymer from unconstrained random walk to the density of typical proteins, secondary structure does not increase significantly (
