energy landscape theory predicts different 
 scenarios for protein 
 folding . in the case of two - state folding , 

the equilibrium distribution is dominated by only two species ( the 
 folded and unfolded states ) separated by large free - energy barriers .
single - molecule experiments are only now starting 
 to shed some light on the transition paths between these states .
however , in bulk the low population of any species 
 intermediate between the native and unfolded forms in the case of 
 two - state folding results in all of the information about the mechanism 
 being effectively hidden .
conversely , in the downhill folding scenario , 
 where free - energy barriers are comparable to the thermal energy ( kbt ) , a myriad of possible structures 
 between the denatured and native states can be populated as the reaction 
 progresses . hence studies on downhill 
 folding can reveal the details of the conformational motions of the 
 polypeptide chain en route to the native state .
one of the experimental approaches used to reach 
 the downhill regime , 
 first applied to the 80-residue amino - terminal -repressor ( fragment 
 685 ) and later to the pin ww domain , consists 
 of engineering the protein by introducing rate - enhancing mutations . in the seminal fluorescence temperature - jump experiments by yang 
 and
gruebele on -repressor , this resulted in the emergence 
 of a faster kinetic phase .
this new phase 
 was interpreted as the downhill relaxation from a vanishing barrier 
 top , related to the so - called folding  speed limit  . a large range of other mutants
have been studied 
 that exhibit differences in both the stability and the rates and amplitudes 
 of this molecular phase .
mutations have been introduced on the basis 
 of needing a fluorescence probe for folding , altering helix propensity in the helices , removing specific interactions , or 
 reducing backbone flexibility .
still , 
 the experimentally intensive work of designing and expressing protein 
 mutants is to some extent based on chemical intuition .
computational 
 screening methods can be of great help in designing mutations in a 
 systematic way . here
we present one such approach that we apply 
 to a simple simulation 
 model for -repressor .
contrary to most existing methods , the 
 focus here is in engineering the folding cooperativity and dynamics 
 instead of modulating protein stability .
we run dynamics simulations that sample both the folded and unfolded 
 states . assuming that the distribution of states on a reaction 
 coordinate is approximately bimodal at the midpoint ( at least , initially ) , 
 we use the proximity of the two peaks , relative to their broadness , 
 as an indication of cooperativity .
we then estimate the effects on 
 this metric arising from mutation of each individual amino acid residue . 

this description is in the spirit of the  cooperative response  
 defined by freire and coworkers .
finally , 
 the results are confirmed by running simulations on a selected mutant 
 and independently assessing the changes in the dynamics using a markov 
 state model ( msm ) , which does not rely on the reaction coordinate 
 used in the engineering step . for computational 
 tractability
, we use as our reference a coarse - grained , topology - based 
 simulation model of -repressor . using our approach , we engineer 
 the folding from being barrier - limited under midpoint conditions to 
 downhill .
this allows us to identify a new hot spot that , in the context 
 of the simple model , will maximally disrupt the cooperativity of folding . 

the approach that we present is general and can be used to modulate 
 barriers in either direction .
in figure 1a , we present the general workflow 
 for our method . 
 an experimental pdb structure file with the coordinates of the protein 
 atoms is the only input and is used to generate a coarse - grained , 
 structure - based ( go ) model for the protein .
the energy landscape of this reference model is sampled 
 via molecular dynamics ( or monte carlo ) simulation , and the resulting 
 trajectories projected onto a suitable reaction coordinate x. although , for our simulation model , we use the fraction 
 of native contacts ( q ) as a progress variable , for 
 more complicated models the reaction coordinate can be variationally 
 optimized .
using the distribution of values 
 of x for both the native ( n ) and unfolded ( u ) states 
 at the midpoint temperature ( tm ) , we define 
 a cooperativity metric1where xn and varn(x ) are , respectively , 
 the mean and variance of x computed over the native 
 state and similarly for the unfolded state .
intuitively it is easy 
 to see that upon a decrease in the free - energy barriers the mean values 
 of x for n and u will get closer to each other and 
 their distributions will become broader ( figure 1b ) .
hence , the lower the barrier the lower the value of . 
 in our method , we predict the changes in this metric upon point mutations , 
 which are systematically introduced for every amino acid residue independently . 
 on the basis of these results ,
interesting mutants are selected , in 
 this case those with low values of . the folding dynamics of 
 these mutants are sampled by running new simulations , and the effect 
 of the mutations on the dynamics is assessed by using an msm .
the 
 advantage of the msm is that it does not rely on the projection on 
 an order parameter , and hence it allows 
 an independent assessment of the effect of the mutation on the dynamics . 

( b ) theoretical free energy surfaces ( top ) and probability 
 distributions ( bottom ) as a function of an order parameter x under midpoint conditions .
we show illustrative examples 
 corresponding to the shift from bimodal ( two - state , red ) to downhill ( green ) folding .
u and n are the unfolded and native states 
 respectively , and the dashed lines mark the mean values xu / n of the order parameter for the unfolded 
 and native states .
we use a 
 protein model coarse - grained to a single bead per amino acid residue 
 located at the c atom .
we construct the model using the experimental x - ray structure of 
 the -repressor mutant ya ( 3kz3 , see figure 2a ) .
the potential - energy 
 function is a sum of harmonic terms for bonds and angles , a statistical 
 potential for the pseudodihedrals , and terms for nonbonded interactions .
favorable nonbonded interactions are limited 
 to residue pairs that are in contact in the reference structure ( shown 
 in the contact map of figure 2b ) .
for one such 
 pair of residues ( i , j ) the interaction 
 energy is defined as2where rij is the distance between c atoms in the instantaneous conformation , 
 ij is the same distance in the 
 reference structure , and ij is 
 the residue - pair specific interaction 
 energy , taken from the miyazawa 
we 
 run simulations of the go model using a modified version of 
 the gromacs simulation package ( version 4.0.5 ) . to propagate the dynamics based on the langevin equation , we use 
 a leapfrog stochastic integrator with a time step of 10 fs .
the simulation 
 data were projected onto suitable order parameters , mainly the root - mean - square 
 deviation ( rmsd ) and the fraction of native contacts ( q ) .
the data from simulations at different temperatures were then 
 combined using the weighted histogram analysis method ( wham ) .
( a ) cartoon 
 representation of the experimental structure of the 
 ya mutant of the 685 protein fragment ( 3kz3 ) .
we show the heavy 
 atoms of the experimental fluorescent probe w22 ( cyan ) , y33 ( blue ) , 
 and the l18 residue ( yellow ) .
( b ) contact map for the ya mutant 
 including all native interaction pairs , as defined in the structure - based 
 model . to determine the 
 contribution that different residues make to the folding cooperativity , 
 for each mutated residue we reduce the strength of its native interactions 
 ( those defined in the contact map , figure 2b ) by a given amount
( see the results and discussion ) . for each of the mutants we calculate the values of qn and qu that appear in eq 1 by reweighting3where the average
is computed over all of 
 the saved configurations from the reference simulations that are in 
 state s  { u , n}. in eq 3 , e is the difference in the energy between the reference 
 simulation and the mutant , emut = eref + e .
the values 
 of the variance vars(q ) are calculated 
 analogously because vars(q ) = qs  qs2 . to restrict these averages to u or n , we set a dividing 
 line between these states to q = 0.55 , which approximately 
 corresponds to the maximum in the free - energy barrier for the go 
 models ( see the results and discussion ) .
we check the accuracy of this calculation by estimating q and var(q ) from the first and 
 second moments of the probability distribution p(q ) obtained from wham ( see supporting 
 information ( si ) , figure s1 ) . to assess the effects of the dynamics 
 in an independent way , we use an msm methodology .
we assume that the dynamics of the system can be described as a 
 discrete - time markov process using the transition matrix t(t ) .
the elements tji of this matrix are the probabilities that 
 being in microstate i the system will be found in 
 microstate j after a lag time t .
the dynamics of the system can then be expressed using the discrete - time 
 analog of the continuous - time master equation4where 
 the column - vector p(kt ) contains the populations of 
 each microstate at time kt .
we calculate the elements of t using the maximum - likelihood 
 estimator5where n(t ) = ( n)ji is the transition 
 count matrix that we obtain directly from the discretized simulation 
 trajectories ( see later ) . we compute equilibrium populations from 
 the right eigenvector of the stationary mode of the transition matrix 
 ( 0r ) 
 and relaxation times i from its 
 eigenvalues i as i = t / ln i .
we define 
 the microstate space using information from the 
 native contacts alone by first discretizing the conformations into 
 strings and then clustering the strings into microstates .
the native 
 contact map discretization is adequate for a go model , as stable 
 non - native interactions can not be formed ; however , the process can 
 be generalized to other systems by considering the full contact matrix , 
 including non - native interactions .
a recent study has found contact - map - based 
 markov models to be more robust than rmsd - based ones . in the same spirit as previous ising - like 
 models for protein folding , we discretize the simulation 
 data assuming that each of the contacts has two possibilities : either 
 being formed ( 1 ) or not ( 0 ) .
hence , for a protein with n native contacts , there are a total of 2 possible strings of zeros and ones , each corresponding to a contact 
 map . in principle
, this discretization would rapidly become intractable 
 for even a small number of contacts .
however , for -repressor 
 ( with a total of 115 contacts in the contact map , see figure 2b ) , we find that in practice very few states are 
 populated due to the limited length of the simulation runs .
( out of 
 more than 4  10 possibilities , 34  788 
 different strings were visited in the 40  000 frames saved during 
 a 4 s simulation time of ya at the midpoint . ) the discretized trajectory 
 is then clustered using a k - medoids algorithm . instead of randomly initializing the cluster 
 centers , we take advantage of the dynamics trajectory , where contiguous 
 snapshots will usually belong in the same energy basin .
we sequentially 
 assign the time series of strings to an existing cluster when the 
 hamming distance between the instantaneous 
 string conformation and the cluster center is shorter than a certain 
 cutoff .
if the hamming distances to the centers of several existing 
 clusters are below the cutoff , the lowest is chosen .
after all the 
 conformations have been assigned sequentially and the initial clustering 
 has been generated , we optimize it by reading the strings from the 
 trajectory in a randomized order and checking the assignment of each 
 individual string to the clusters . the cluster centers are updated 
 in the event that a string is added , which is more central within 
 the cluster than the existing cluster center .
the procedure is repeated 
 until no strings are reassigned in a randomized parsing of the trajectory . 

we find that this procedure is very efficient in generating structurally 
 meaningful clusters as the initial assignment already produces a very 
 good first guess .
the final number of clusters depends on the 
 value of the hamming distance cutoff , with lower values resulting 
 in a higher number of clusters .
a low cutoff is , in principle , desirable 
 to produce finely resolved clusters .
however , we find that the fraction 
 of the simulation data accounted for by the frequently visited clusters , 
 here defined as those visited for an aggregate simulation time of 
 at least 10 ns , is reduced with decreasing cutoffs .
we therefore choose 
 a hamming distance cutoff as a compromise between resolving multiple 
 clusters in both the unfolded and native states and maximizing the 
 number of trajectory frames included in the frequently visited clusters .
transitions between microstates for the 
 construction of the master equation model are calculated directly . 
 a transition between microstate
i and j is assigned every time that the simulation jumps from one frequently 
 visited cluster to another after a certain lag time t .
when the trajectory reaches an infrequently visited cluster j , we consider that it remains in the initial microstate i. this procedure produces less accurate kinetics than transition - based 
 assignment , but in this case it allows 
 us to capture the qualitative differences between the different models .
we start by analyzing the simulations of 
 the coarse - grained topology - based 
 model of the -repressor .
we use the experimental structure 
 of the ya mutant ( see methods ) that 
 differs from the wild type ( wt ) by only four internal amino acid residues 
 ( 92.5% sequence identity ) and that is also very similar structurally 
 ( rmsd = 0.7  ) .
we choose this structure as ya has been 
 proposed to fold downhill under native conditions and to have a small 
 ( > 3 kbt ) barrier at 
 its tm .
also , it was crystallized as a 
 shorter sequence 
 than the wt and in the absence of dna , which makes it a closer match 
 to the experimental construct . from the projection 
 of our go model simulation trajectories on the fraction of 
 native contacts ( q )
, we find primarily two interconverting 
 states under midpoint conditions ( figure 3a ) . 

the potentials of mean force for the ya mutant ( figure 3c ) indicate that the protein folds downhill under 
 native conditions and is barrier - limited at the simulation midpoint 
 ( tm  310 k ) , consistent with the 
 experimental results . on the basis of 
 the projection on q , the barrier is 2.3 kcal / mol 
 at the midpoint ( i.e. , gu = 3.7kbt ) , in agreement with results by gruebele and coworkers .
we also observe a sparsely populated intermediate 
 state on the folded side of the dominant barrier ( q  0.7 ) , with helix 5 slightly detached from the protein core . 

additional support for its presence comes from 
 explicit solvent simulations , where many contacts that involve helix 
 5 form late in transition paths , the 
 low helical propensity predicted by agadir , and 
 the high b factors .
we note that we have 
 also found this state to appear when the wt structure is used to construct 
 the go model , making this a robust prediction .
( a ) time series for the 
 projection of the simulation trajectory 
 for ya on the fraction of native contacts ( q ) at 310 k. ( b ) heat capacity thermogram calculated from wham .
( c ) 
 potentials of mean force at native ( top ) and midpoint ( bottom ) conditions . in the context of
the simple go 
 model a natural way of simulating 
 mutations is just scaling the contact energies ij in eq 2 .
for this scaling , 
 we choose a value of 0.5 ( i.e. , scaling by half the contact strength 
 of every contact made by the mutated amino acid residue ) to calculate 
 the cooperativity metric . this type of conservative mutation , 
 corresponding to a small decrease in the strength of the interactions 
 formed by the selected residue , is likely to destabilize the folded 
 state ; the effect is similar to the small 
 perturbations used in experimental -value analysis rather than the disruptive effects that might 
 arise from introducing a more bulky or charged group . in figure 4 , we show the resulting values of  of each 
 of the  single - point mutants  normalized by the value 
 for ya , which are generally very close to that for the reference 
 simulation , suggesting very small changes in the cooperativity .
the 
 robustness to mutations is in agreement with a number of studies that 
 suggest that cooperativity and free - energy barriers are carefully 
 selected features of protein - energy landscapes
. however , there are a number of cooperativity hot spots that , according 
 to this calculation , can reduce the folding cooperativity ( shown in 
 green in figure 4a ) .
to calculate the 
 cooperativity metric  , we are relying on the adequacy of q as a reaction coordinate , but the results could be different 
 for alternative progress parameters for the folding reaction . in figure 4b , we compare 
 the estimates of the cooperativity metric from the fraction of native 
 contacts ( q ) and the rmsd ( rmsd ) .
we find that the agreement between the two estimates 
 is extremely good , with a pearson correlation coefficient of r = 0.88 , that increases to 0.97 if we remove the outlier 
 f76 .
the high q and low rmsd for this mutation point to a deviation from the expected 
 behavior of the mutational approach ( i.e. , the model probability distributions 
 from figure 1b ) .
for f76 , the mutation results 
 in a population shift from the native to the intermediate , which differently 
 affects the estimates of q and var(q ) from different reaction coordinates ( see si , figure s2 ) .
it is therefore advisable to 
 examine the distributions that are used to calculate  , as this 
 will allow one to spot possible deviations from the expected behavior . 

although alternative descriptions of  , for example , based on 
 the cooperativity parameter c or the calorimetric criterion , may be able to overcome this limitation , we do not pursue them here .
( a ) normalized 
 value of the cooperativity metric  from the 
 native contacts q upon a single - point mutation for 
 every amino acid residue in the ya sequence .
( b ) correlation between  
 calculated from q and  calculated from the 
 rmsd from native.(c ) cartoon of the structure of -repressor 
 with color code indicating the value of . the color scale goes 
 from high ( red ) to low ( green ) . spheres are shown for the 
 residues with > 10% decrease in . the transparent surface 
 envelops 
 the residues identified by gruebele and coworkers to be the folding 
 core of the protein ( see text ) .
we also analyze the distribution of  values in the 
 3-d structure 
 of the -repressor ( figure 4c ) . a number 
 of hot spots cluster around a central core formed by helices 1 and 
 4 and the turns between helices 34 and 45 .
this indicates 
 that the network of contacts formed by these residues is the most 
 sensitive part of the protein for folding cooperativity .
the cluster 
 of hot spots that we identify is in remarkable agreement with the 
 region that gruebele and coworkers have expressed in the construct 
 blue1 , consisting of the two - helix bundle formed 
 by helices 1 and 4 connected by linkers .
they found that blue1 folds with similar tm and 
 rates as the 685 fragment , indicating that 
 this region comprises the minimal folding core of the protein .
looking at the predicted values of our cooperativity 
 metric , we find that out of the 23 residues with more than a 10% decrease 
 in  , 18 are within the cooperative core proposed by gruebele 
 and coworkers ( figure 4c ) .
although our results 
 also agree with some of the experimental mutations ( e.g. , asp14 ) , 
 in other cases our simple method for changing the energetics fails 
 to predict changes in the barriers derived from analysis of t - jump 
 experiments , indicating that a more detailed energy function will 
 be needed for quantitatively reproducing changes in tm values and free - energy barriers .
our method , however , 
 is a powerful tool for making approximate predictions that can direct 
 mutational analysis from experimentalists .
we focus on the leu18 
 single - point mutant , which according to our 
 calculation would produce the greatest decrease in . ( see figure 4a , b . ) to test this prediction , we run simulations 
 of the go model where we scale the ij of the l18 contacts by 0.5 ( l180.5 ) and where we remove 
 the contact altogether ( l180 ) , leaving only an excluded volume 
 term . in the projection on
q we see that transitions 
 between the folded and unfolded states are much faster for the l180.5 mutant and become entirely diffusive for l180 . 

( see figure 5a . ) from the wham analysis of 
 the simulations , we confirm the results of our prediction : upon mutation , 
 the heat - capacity curves become considerably broader , a characteristic 
 signature of the reduction of the cooperativity ( figure 5b ) .
also , potentials of 
 mean force on q reveal a decrease in the midpoint 
 free - energy barrier from 3.7 kbt ( ya ) to 1.4 kbt ( l180 ) , making the full mutant an apparent downhill 
 folder , at least according to this projection ( figure 5c ) .
the differences in the free - energy barriers are due to 
 small changes in both the enthalpic and entropic contributions to 
 the free energy ( see si , figure s3 ) .
time series 
 for the projection of the simulation trajectory for 
 the ya mutants l180.5 ( blue ) and l180 ( green ) 
 on the fraction of native contacts ( q ) at their midpoint 
 temperatures .
( c ) potentials of mean force under native ( top ) and midpoint 
 ( bottom ) conditions .
the reduction of the free - energy 
 barriers we observe could be due 
 to the choice of a suboptimal progress coordinate for folding .
to 
 assess the effect on the dynamics independently of the projection 
 on q , we construct a transition network from the 
 discretized simulation trajectory ( see the methods ) .
the k - medoids algorithm results in 21 to 23 clusters , 
 accounting for 8993% of the total simulation time . in figure 6 .
we show the clusters represented by the mean and 
 standard deviation of the values of q and rmsd of 
 the corresponding conformations .
when overlaid on the potential of 
 mean force , the clusters appear reasonably well - separated , particularly 
 for the intermediate to high values of q. it is important 
 to note that while for ya there is a gap between clusters on 
 the folded and unfolded sides of the barrier ( q  
 0.55 ) , for the l180 mutant such a gap does not exist as a result 
 of the non - negligible population in the barrier region .
( see figure 6c . ) in figure 6d , we show 
 snapshots corresponding to the seven most populated clusters for the 
 l180 mutant , including folded ( n ) and unfolded ( u ) clusters 
 and multiple clusters in the transition region ( t )
however , this is natural because 
 we have used the native contact - map for the clustering .
the average 
 contact maps clearly indicate the structure formation events that 
 are taking place en route to folding .
potentials of mean force ( in kcal / mol ) 
 for q and 
 the rmsd to native for the different models for -repressor : 
 ya ( a ) , l180.5 ( b ) , and l180 ( c ) . in each case , 
 we overlay ellipsoids centered at the average rmsd and q values for the cluster with principal axes of one standard deviation . 

( d ) snapshots and average contact maps for the seven most populated 
 clusters for the l180 model , with the same color code and 
 names as those in panel c. we construct the msm for the three different models with 
 a lag 
 time t = 1 ns , for which eigenvalues are approximately 
 converged , as required for markovian dynamics ( not shown ) .
we analyze 
 the relaxation times , i , that we 
 calculate from the eigenvalues of the transition matrix , t. ( see figure 7a . )
the slowest mode is about 
 one order of magnitude slower for the ya than it is for the 
 l180 mutant , with the l180.5 mutant being somewhere 
 in between . according to the sign structure of the right eigenvectors 
 1r
, this 
 mode in fact corresponds to exchange between high q and low q microstates ( i.e. , folding , see figure 7b ) .
hence the speedup in 1 is 
 consistent with the reduction in the free - energy barrier in the projection 
 on q. also , for the two mutants we find a reduction 
 of the gap between the first and second slowest modes , a characteristic 
 feature of downhill folding .
interestingly , 
 the next few modes appear at very similar values of the relaxation 
 time for all three models , although they correspond to substantially 
 different processes .
for example , in the case of the ya , the 
 second mode 2r corresponds to the exchange of the native cluster with the 
 helix - detached intermediate , while for the mutants a more diverse 
 set of configurations is involved , including clusters that in the 
 two - state case would be at the top of the barrier ( q  0.55 , see figure 7b ) .
( a ) spectrum 
 of relaxation times for the msm of the -repressor 
 model and the l18 mutants .
( b ) top row shows potentials of mean force 
 ( in kcal / mol ) with overlaid circles indicating the mean values of q for the different microstates of the msm .
the next three rows show the values of the 
 slowest right eigenvectors ir projected on q. shaded areas are 
 shown to illustrate the sign structure of the eigenvectors . to gain further insight into how 
 the dynamical signatures in our 
 model mutants may result in different signals in kinetic experiments , 
 we calculate relaxations using the derived transition matrices .
we 
 use eq 4 to propagate the dynamics from a theoretical 
 initial distribution that we set to be the fully folded microstate 
 ( i.e. , p(0 ) is a vector of zeros but for the fully folded 
 state n ) . to calculate a proxy for the signal of the most usually 
 studied experimental probe , the fluorescence of trp22
, we use the 
 fraction of native contacts of this residue ( qw22 , figure 8a ) .
the relaxation is faster 
 as we go from ya ( figure 8b ) to the 
 engineered models ( figure 8c , d ) .
however , the 
 two - state approximation also starts to become worse as we approach 
 the downhill mutant , which requires three exponentials to fit the 
 data .
the fastest phase in the downhill mutant appears on a time scale 
 similar to the fast eigenmodes of the msm .
this kinetic analysis hence 
 validates the results from the projection approach , confirming the 
 decrease in free - energy barriers , speed - up in rates , and emergence 
 of  strange kinetics  .
residues that form contacts 
 with the w22 side chain are shown in atomic detail for heavy atoms 
 and as transparent spheres centered in the c. ( b ) 
 decay of the normalized fraction of native w22 contacts , qw22 , for ya ( circles ) , with lines corresponding 
 to single ( violet ) , double ( cyan ) , and triple ( magenta ) exponential 
 fitting expressions .
we present a new approach 
 to computational protein engineering 
 that can be very useful as a tool to guide the search for relevant 
 mutants to be studied experimentally .
in particular , we engineer the 
 transition from two - state to downhill folding by tuning the cooperativity 
 of a protein model that folds in a two - state fashion at the midpoint 
 and reaches downhill folding upon single - point modifications in the 
 simplified energy function .
the method is novel in that it focuses 
 on modulating the probability distribution for a given order parameter 
 and therefore accounts for the ensemble nature of protein folding . 
 in this respect ,
our approach is similar to smarteps , an ising - model - based 
 method recently developed to predict changes in stability in protein 
 mutants . in this case
, we assume that q is a good order parameter for folding and modulate the 
 probability distribution for this parameter .
this is a reasonable 
 assumption in the context of a go model . in more complicated 
 scenarios , a step involving the optimization of the reaction
coordinate can be incorporated as part of step 2 in the algorithm ( figure 1a ) .
also , we approximate the probability distribution 
 on q as being bimodal ( see eq 1 ) .
this may break down when the two - state approximation does not 
 hold , particularly if we use more detailed ( i.e. , atomistic ) protein 
 models .
however , we note that for a large subset of single - domain 
 proteins studied with microsecond atomistic md simulations , the distributions 
 along a coordinate are still largely two - state , and many small domains have been shown to fold 
 in a two - state fashion experimentally . while most 
 methods focus on recovering the correct tm , here emphasis
we validate the results 
 of the engineering by comparing the msm of the original protein model 
 with that of the mutants . in the context of a more detailed model 
 it will be possible to make direct comparison with the exact kinetic 
 signatures observed for different mutants .
the work that we 
 present here is based on a simple topology - based 
 model that considers only the interactions present in the native conformation 
 of the protein .
non - native interactions could modify the emerging 
 picture , although these have been found to have only relatively small 
 effects on protein folding cooperativity .
our results do , however , stand by themselves , as there are a number 
 of experimental references that we are able to reproduce .
first , the 
 general features of folding of ya such as two - state folding 
 near the midpoint are captured by the 
 simple go model .
in addition , we are able to locate a cooperative 
 folding core that overlaps with that identified by gruebele and coworkers .
although it is not possible to answer that 
 question in advance , we speculate that they will , as our estimate 
 of the melting temperature for l180 ( tm  290 k ) involves a decrease to 92% 
 of the tm of the original model ( 315 k ) , 
 while the experimental tm of ya 
 is 344 k. our approach is particularly promising in the context 
 of a more 
 detailed description of the effects of mutations in the model for 
 the interactions . for example , we show one possibility here , 
 in which we replace the ij of the 
 residue pairs that change upon mutation with the corresponding values 
 from the miyazawa  jernigan interaction matrix and substitute the knowledge - based torsion terms according 
 to the mutation . using this method for the computational protein engineering , 

the agreement with the experimental tm for a database of 17 different mutants is very good ( see figure 9a ) , although this does not guarantee that the barrier 
 heights will be reproduced .
jernigan 
 contact energies required to produce the mutations from our reference 
 sequence ( figure 9b ) support the scaling factors 
 ( particularly the 0.5 scaling ) used in this study . taken together
, 
 these results suggest the future directions for refining the current 
 approach with a more accurate description of the energetics
. this 
 will be possible by carefully calibrating the results of the model 
 against extensive data sets of mutation effects in the thermodynamics and kinetics of 
 protein folding .
( a ) correlation between experimentally determined tm values and those calculated by substituting 
 the ij involving the mutated residue 
 in the simulations . 

 one 
 strategy for reaching the downhill folding regime , primarily 
 exploited for the 685 protein fragment , 
 consists of cumulatively introducing mutations that speed up folding . 
 
 this is an experimentally demanding process where chemical intuition 
 usually serves as a guide for the choice of amino acid residues to 
 mutate . 
 such an approach can be aided by computational methods that 
 screen for protein engineering hot spots . here 
 we present one such 
 method that involves sampling the energy landscape of the pseudo - wild - type 
 protein and investigating the effect of point mutations on this landscape . 
 using a novel metric for the cooperativity , we identify those residues 
 leading to the least cooperative folding . 
 the folding dynamics of 
 the selected mutants are then directly characterized and the differences 
 in the kinetics are analyzed within a markov - state model framework . 
 
 although the method is general , here we present results for a coarse - grained 
 topology - based simulation model of -repressor , whose barrier 
 is reduced from an initial value of 4kbt at the midpoint to 1kbt , thereby reaching the downhill folding 
 regime . 
Introduction
Methods
Results and Discussion
Conclusions