Chemical, physical, and theoretical kinetics of an ultrafast folding protein

Tutorial on how to build non-Markovian dynamic models from molecular dynamics simulations for studying protein conformational changes

Protein conformational changes play crucial roles in their biological functions. In recent years, the Markov State Model (MSM) constructed from extensive Molecular Dynamics (MD) simulations has emerged as a powerful tool for modeling complex protein conformational changes. In MSMs, dynamics are modeled as a sequence of Markovian transitions among metastable conformational states at discrete time intervals (called lag time). A major challenge for MSMs is that the lag time must be long enough to allow transitions among states to become memoryless (or Markovian). However, this lag time is constrained by the length of individual MD simulations available to track these transitions. To address this challenge, we have recently developed Generalized Master Equation (GME)-based approaches, encoding non-Markovian dynamics using a time-dependent memory kernel. In this Tutorial, we introduce the theory behind two recently developed GME-based non-Markovian dynamic models: the quasi-Markov State Model (qMSM) and the Integrative Generalized Master Equation (IGME). We subsequently outline the procedures for constructing these models and provide a step-by-step tutorial on applying qMSM and IGME to study two peptide systems: alanine dipeptide and villin headpiece. This Tutorial is available at https://github.com/xuhuihuang/GME_tutorials. The protocols detailed in this Tutorial aim to be accessible for non-experts interested in studying the biomolecular dynamics using these non-Markovian dynamic models.

Brownian Diffusion of Hexagonal Boron Nitride Nanosheets and Graphene in Two Dimensions

Two-dimensional (2D) nanomaterials have numerous interesting chemical and physical properties that make them desirable building blocks for the manufacture of macroscopic materials. Liquid-phase processing is a common method for forming macroscopic materials from these building blocks including wet-spinning and vacuum filtration. As such, assembling 2D nanomaterials into ordered functional materials requires an understanding of their solution dynamics. Yet, there are few experimental studies investigating the hydrodynamics of disk-like materials. Herein, we report the lateral diffusion of hexagonal boron nitride nanosheets (h-BN and graphene) in aqueous solution when confined in 2-dimensions. This was done by imaging fluorescent surfactant-tagged nanosheets and visualizing them by using fluorescence microscopy. Spectroscopic studies were conducted to characterize the interactions between h-BN and the fluorescent surfactant, and atomic force microscopy (AFM) was conducted to characterize the quality of the dispersion. The diffusion data under different gap sizes and viscosities displayed a good correlation with Kramers' theory. We propose that the yielded activation energies by Kramers' equation express the magnitude of the interaction between fluorescent surfactant tagged h-BN and glass because the energies remain constant with changing viscosity and decrease with increasing confinement size. The diffusion of graphene presented a similar trend with similar activation energy as the h-BN. This relationship suggests that Kramers' theory can also be applied to simulate the diffusion of other 2D nanomaterials.

Toward a Benchmark for Markov State Models: The Folding of HP35

Adopting a 300 μs long MD trajectory of the folding of villin headpiece (HP35) by D. E. Shaw Research, we recently constructed a Markov state model (MSM) based on inter-residue contacts. The model reproduces the folding time and predicts that the native basin and unfolded region consist of metastable substates that are structurally well-characterized. Recognizing the need to establish well-defined benchmark problems, we study to what extent and in what sense this MSM can be employed as a reference model. Hence, we test the robustness of the MSM by comparing it to models that use alternative combinations of features, dimensionality reduction methods, and clustering schemes. The study suggests some main characteristics of the folding of HP35 that should be reproduced by other competitive models. Moreover, the discussion reveals which parts of the MSM workflow matter most for the considered problem and illustrates the promises and pitfalls of state-based models for the interpretation of biomolecular simulations.

Selecting Features for Markov Modeling: A Case Study on HP35

Markov state models represent a popular means to interpret molecular dynamics trajectories in terms of memoryless transitions between metastable conformational states. To provide a mechanistic understanding of the considered biomolecular process, these states should reflect structurally distinct conformations and ensure a time scale separation between fast intrastate and slow interstate dynamics. Adopting the folding of villin headpiece (HP35) as a well-established model problem, here we discuss the selection of suitable input coordinates or "features", such as backbone dihedral angles and interresidue distances. We show that dihedral angles account accurately for the structure of the native energy basin of HP35, while the unfolded region of the free energy landscape and the folding process are best described by tertiary contacts of the protein. To construct a contact-based model, we consider various ways to define and select contact distances and introduce a low-pass filtering of the feature trajectory as well as a correlation-based characterization of states. Relying on input data that faithfully account for the mechanistic origin of the studied process, the states of the resulting Markov model are clearly discriminated by the features, describe consistently the hierarchical structure of the free energy landscape, and─as a consequence─correctly reproduce the slow time scales of the process.

The Wako-Saitô-Muñoz-Eaton Model for Predicting Protein Folding and Dynamics

Despite the recent advances in the prediction of protein structures by deep neutral networks, the elucidation of protein-folding mechanisms remains challenging. A promising theory for describing protein folding is a coarse-grained statistical mechanical model called the Wako-Saitô-Muñoz-Eaton (WSME) model. The model can calculate the free-energy landscapes of proteins based on a three-dimensional structure with low computational complexity, thereby providing a comprehensive understanding of the folding pathways and the structure and stability of the intermediates and transition states involved in the folding reaction. In this review, we summarize previous and recent studies on protein folding and dynamics performed using the WSME model and discuss future challenges and prospects. The WSME model successfully predicted the folding mechanisms of small single-domain proteins and the effects of amino-acid substitutions on protein stability and folding in a manner that was consistent with experimental results. Furthermore, extended versions of the WSME model were applied to predict the folding mechanisms of multi-domain proteins and the conformational changes associated with protein function. Thus, the WSME model may contribute significantly to solving the protein-folding problem and is expected to be useful for predicting protein folding, stability, and dynamics in basic research and in industrial and medical applications.

Correlation-Based Feature Selection to Identify Functional Dynamics in Proteins

To interpret molecular dynamics simulations of biomolecular systems, systematic dimensionality reduction methods are commonly employed. Among others, this includes principal component analysis (PCA) and time-lagged independent component analysis (TICA), which aim to maximize the variance and the time scale of the first components, respectively. A crucial first step of such an analysis is the identification of suitable and relevant input coordinates (the so-called features), such as backbone dihedral angles and interresidue distances. As typically only a small subset of those coordinates is involved in a specific biomolecular process, it is important to discard the remaining uncorrelated motions or weakly correlated noise coordinates. This is because they may exhibit large amplitudes or long time scales and therefore will be erroneously considered important by PCA and TICA, respectively. To discriminate collective motions underlying functional dynamics from uncorrelated motions, the correlation matrix of the input coordinates is block-diagonalized by a clustering method. This strategy avoids possible bias due to presumed functional observables and conformational states or variation principles that maximize variance or time scales. Considering several linear and nonlinear correlation measures and various clustering algorithms, it is shown that the combination of linear correlation and the Leiden community detection algorithm yields excellent results for all considered model systems. These include the functional motion of T4 lysozyme to demonstrate the successful identification of collective motion, as well as the folding of the villin headpiece to highlight the physical interpretation of the correlated motions in terms of a functional mechanism.

Predicting and Simulating Mutational Effects on Protein Folding Kinetics

Mutational perturbations of protein structures, i.e., phi-value analysis, are commonly employed to probe the extent of involvement of a particular residue in the rate-determining step(s) of folding. This generally involves the measurement of folding thermodynamic parameters and kinetic rate constants for the wild-type and mutant proteins. While computational approaches have been reasonably successful in understanding and predicting the effect of mutations on folding thermodynamics, it has been challenging to explore the same on kinetics due to confounding structural, energetic, and dynamic factors. Accordingly, the frequent observation of fractional phi-values (mean of ~0.3) has resisted a precise and consistent interpretation. Here, we describe how to construct, parameterize, and employ a simple one-dimensional free energy surface model that is grounded in the basic tenets of the energy landscape theory to predict and simulate the effect of mutations on folding kinetics. As a proof of principle, we simulate one-dimensional free energy profiles of 806 mutations from 24 different proteins employing just the experimental destabilization as input, reproduce the relative unfolding activation free energies with a correlation of 0.91, and show that the mean phi-value of 0.3 essentially corresponds to the extent of stabilization energy gained at the barrier top while folding.

Folding Intermediates, Heterogeneous Native Ensembles and Protein Function

Single domain proteins fold via diverse mechanisms emphasizing the intricate relationship between energetics and structure, which is a direct consequence of functional constraints and demands imposed at the level of sequence. On the other hand, elucidating the interplay between folding mechanisms and function is challenging in large proteins, given the inherent shortcomings in identifying metastable states experimentally and the sampling limitations associated with computational methods. Here, we show that free energy profiles and surfaces of large systems (>150 residues), as predicted by a statistical mechanical model, display a wide array of folding mechanisms with ubiquitous folding intermediates and heterogeneous native ensembles. Importantly, residues around the ligand binding or enzyme active site display a larger tendency to partially unfold and this manifests as intermediates or excited states along the folding coordinate in ligand binding domains, transcription repressors, and representative enzymes from all the six classes, including the SARS-CoV-2 receptor binding domain (RBD) of the spike protein and the protease M^pro. It thus appears that it is relatively easier to distill the imprints of function on the folding landscape of larger proteins as opposed to smaller systems. We discuss how an understanding of energetic-entropic features in ordered proteins can pinpoint specific avenues through which folding mechanisms, populations of partially structured states and function can be engineered.

Combining Ancestral Reconstruction with Folding-Landscape Simulations to Engineer Heterologous Protein Expression

Obligate symbionts typically exhibit high evolutionary rates. Consequently, their proteins may differ considerably from their modern and ancestral homologs in terms of both sequence and properties, thus providing excellent models to study protein evolution. Also, obligate symbionts are challenging to culture in the lab and proteins from uncultured organisms must be produced in heterologous hosts using recombinant DNA technology. Obligate symbionts thus replicate a fundamental scenario of metagenomics studies aimed at the functional characterization and biotechnological exploitation of proteins from the bacteria in soil. Here, we use the thioredoxin from Candidatus Photodesmus katoptron, an uncultured symbiont of flashlight fish, to explore evolutionary and engineering aspects of protein folding in heterologous hosts. The symbiont protein is a standard thioredoxin in terms of 3D-structure, stability and redox activity. However, its folding outside the original host is severely impaired, as shown by a very slow refolding in vitro and an inefficient expression in E. coli that leads mostly to insoluble protein. By contrast, resurrected Precambrian thioredoxins express efficiently in E. coli, plausibly reflecting an ancient adaptation to unassisted folding. We have used a statistical-mechanical model of the folding landscape to guide back-to-ancestor engineering of the symbiont protein. Remarkably, we find that the efficiency of heterologous expression correlates with the in vitro (i.e., unassisted) folding rate and that the ancestral expression efficiency can be achieved with only 1-2 back-to-ancestor replacements. These results demonstrate a minimal-perturbation, sequence-engineering approach to rescue inefficient heterologous expression which may potentially be useful in metagenomics efforts targeting recent adaptations.

Modern Kinetics and Mechanism of Protein Folding: A Retrospective

Modern experimental kinetics of protein folding began in the early 1990s with the introduction of nanosecond laser pulses to trigger the folding reaction, providing an almost 10⁶-fold improvement in time resolution over the stopped-flow method being employed at the time. These experiments marked the beginning of the "fast-folding" subfield that enabled investigation of the kinetics of formation of secondary structural elements and disordered loops for the first time, as well as the fastest folding proteins. When I started to work on this subject, a fast folding protein was one that folded in milliseconds. There were, moreover, no analytical theoretical models and no atomistic or coarse-grained molecular dynamics simulations to describe the mechanism. Two of the most important discoveries from my lab since then are a protein that folds in hundreds of nanoseconds, as determined from nanosecond laser temperature experiments, and the discovery that the theoretically predicted barrier crossing time is about the same for proteins that differ in folding rates by 10⁴-fold, as determined from single molecule fluorescence measurements. We also developed what has been called the "Hückel model" of protein folding, which quantitatively explains a wide range of equilibrium and kinetic measurements. This retrospective traces the history of contributions to the "fast folding" subfield from my lab until about 3 years ago, when I left protein folding to spend the rest of my research career trying to discover an inexpensive drug for treating sickle cell disease.

MSMPathfinder: Identification of Pathways in Markov State Models

Markov state models represent a popular means to interpret biomolecular processes in terms of memoryless transitions between metastable conformational states. To gain insight into the underlying mechanism, it is instructive to determine all relevant pathways between initial and final states of the process. Currently available methods, such as Markov chain Monte Carlo and transition path theory, are convenient for identifying the most frequented pathways. They are less suited to account for the typically huge amount of pathways with low probability which, though, may dominate the cumulative flux of the reaction. On the basis of a systematic construction of all possible pathways, the here proposed method MSMPathfinder is able to characterize the multitude of unique pathways (say, up to 10¹⁰) in a complex system and to quantitatively calculate their correct weights and associated waiting times with predefined accuracy. Adopting the chiral transitions of a peptide helix and the folding of the villin headpiece as model problems, mechanisms and associated waiting times of these processes are discussed using a kinetic network representation. The analysis reveals that the waiting time distribution may yield only little insight into the diversity of pathways, because the measured folding times do typically not reflect the most probable path lengths but rather the cumulative effect of many different pathways.

Sampling of the conformational landscape of small proteins with Monte Carlo methods

Computer simulation provides an increasingly realistic picture of large-scale conformational change of proteins, but investigations remain fundamentally constrained by the femtosecond timestep of molecular dynamics simulations. For this reason, many biologically interesting questions cannot be addressed using accessible state-of-the-art computational resources. Here, we report the development of an all-atom Monte Carlo approach that permits the modelling of the large-scale conformational change of proteins using standard off-the-shelf computational hardware and standard all-atom force fields. We demonstrate extensive thermodynamic characterization of the folding process of the α-helical Trp-cage, the Villin headpiece and the β-sheet WW-domain. We fully characterize the free energy landscape, transition states, energy barriers between different states, and the per-residue stability of individual amino acids over a wide temperature range. We demonstrate that a state-of-the-art intramolecular force field can be combined with an implicit solvent model to obtain a high quality of the folded structures and also discuss limitations that still remain.

A modified Ising model of Barabási-Albert network with gene-type spins

The central question of systems biology is to understand how individual components of a biological system such as genes or proteins cooperate in emerging phenotypes resulting in the evolution of diseases. As living cells are open systems in quasi-steady state type equilibrium in continuous exchange with their environment, computational techniques that have been successfully applied in statistical thermodynamics to describe phase transitions may provide new insights to the emerging behavior of biological systems. Here we systematically evaluate the translation of computational techniques from solid-state physics to network models that closely resemble biological networks and develop specific translational rules to tackle problems unique to living systems. We focus on logic models exhibiting only two states in each network node. Motivated by the apparent asymmetry between biological states where an entity exhibits boolean states i.e. is active or inactive, we present an adaptation of symmetric Ising model towards an asymmetric one fitting to living systems here referred to as the modified Ising model with gene-type spins. We analyze phase transitions by Monte Carlo simulations and propose a mean-field solution of a modified Ising model of a network type that closely resembles a real-world network, the Barabási–Albert model of scale-free networks. We show that asymmetric Ising models show similarities to symmetric Ising models with the external field and undergoes a discontinuous phase transition of the first-order and exhibits hysteresis. The simulation setup presented herein can be directly used for any biological network connectivity dataset and is also applicable for other networks that exhibit similar states of activity. The method proposed here is a general statistical method to deal with non-linear large scale models arising in the context of biological systems and is scalable to any network size.

A Disordered Loop Mediates Heterogeneous Unfolding of an Ordered Protein by Altering the Native Ensemble

The high flexibility of long disordered or partially structured loops in folded proteins allows for entropic stabilization of native ensembles. Destabilization of such loops could alter the native ensemble or promote alternate conformations within the native ensemble if the ordered regions themselves are held together weakly. This is particularly true of downhill folding systems that exhibit weak unfolding cooperativity. Here, we combine experimental and computational methods to probe the response of the native ensemble of a helical, downhill folding domain PDD, which harbors an 11-residue partially structured loop, to perturbations. Statistical mechanical modeling points to continuous structural changes on both temperature and mutational perturbations driven by entropic stabilization of partially structured conformations within the native ensemble. Long time-scale simulations of the wild-type protein and two mutants showcase a remarkable conformational switching behavior wherein the parallel helices in the wild-type protein sample an antiparallel orientation in the mutants, with the C-terminal helix and the loop connecting the helices displaying high flexibility, disorder, and non-native interactions. We validate these computational predictions via the anomalous fluorescence of a native tyrosine located at the interface of the helices. Our observations highlight the role of long loops in determining the unfolding mechanisms, sensitivity of the native ensembles to mutational perturbations and provide experimentally testable predictions that can be explored in even two-state folding systems.

New Insights into Folding, Misfolding, and Nonfolding Dynamics of a WW Domain

Intermediate states in protein folding are associated with formation of amyloid fibrils, which are responsible for a number of neurodegenerative diseases. Therefore, prevention of the aggregation of folding intermediates is one of the most important problems to overcome. Recently, we studied the origins and prevention of formation of intermediate states with the example of the Formin binding protein 28 (FBP28) WW domain. We demonstrated that the replacement of Leu26 by Asp26 or Trp26 (in ∼15% of the folding trajectories) can alter the folding scenario from three-state folding, a major folding scenario for the FBP28 WW domain (WT) and its mutants, toward two-state or downhill folding at temperatures below the melting point. Here, for a better understanding of the physics of the formation/elimination of intermediates, (i) the dynamics and energetics of formation of β-strands in folding, misfolding, and nonfolding trajectories of these mutants (L26D and L26W) is investigated; (ii) the experimental structures of WT, L26D, and L26W are analyzed in terms of a kink (heteroclinic standing wave solution) of a generalized discrete nonlinear Schrödinger equation. We show that the formation of each β-strand in folding trajectories is accompanied by the emergence of kinks in internal coordinate space as well as a decrease in local free energy. In particular, the decrease in downhill folding trajectory is ∼7 kcal/mol, while it varies between 31 and 48 kcal/mol for the three-state folding trajectory. The kink analyses of the experimental structures give new insights into formation of intermediates, which may become a useful tool for preventing aggregation.

Studying biomolecular folding and binding using temperature-jump mass spectrometry

Characterizing folding and complex formation of biomolecules provides a view into their thermodynamics, kinetics and folding pathways. Deciphering kinetic intermediates is particularly important because they can often be targeted by drugs. The key advantage of native mass spectrometry over conventional methods that monitor a single observable is its ability to identify and quantify coexisting species. Here, we show the design of a temperature-jump electrospray source for mass spectrometry that allows one to perform fast kinetics experiments (0.16-32 s) at different temperatures (10-90 °C). The setup allows recording of both folding and unfolding kinetics by using temperature jumps from high to low, and low to high, temperatures. Six biological systems, ranging from peptides to proteins to DNA complexes, exemplify the use of this device. Using temperature-dependent experiments, the folding and unfolding of a DNA triplex are studied, providing detailed information on its thermodynamics and kinetics.

Identifying metastable states of biomolecules by trajectory mapping and density peak clustering

Efficiently and accurately analyzing high-dimensional time series, such as the molecular dynamics (MD) trajectory of biomolecules, is a long-standing and intriguing task. Two different but related techniques, i.e., dimension reduction methods and clustering algorithms, have been developed and applied widely in this field. Here we show that the combination of these techniques enables further improvement of the analyses, especially with very complicated data. Specifically, we present an approach that combines the trajectory mapping (TM) method, which constructs slow collective variables of a time series, with density peak clustering (DPC) [A. Rodriguez and A. Laio, Science 344, 1492 (2014)SCIEAS0036-807510.1126/science.1242072], which identifies similar data points to form clusters in a static data set. We illustrate the application of the TMDPC approach with hundreds of microseconds of all-atomic MD trajectories of two proteins, the villin headpiece and protein G. The results show that TMDPC is a powerful tool for achieving the metastable states and slow dynamics of these high-dimensional time series due to the efficient consideration of the time successiveness and the geometric distances between data points.

Thermodynamics and folding landscapes of large proteins from a statistical mechanical model

Statistical mechanical models that afford an intermediate resolution between macroscopic chemical models and all-atom simulations have been successful in capturing folding behaviors of many small single-domain proteins. However, the applicability of one such successful approach, the Wako-Saitô-Muñoz-Eaton (WSME) model, is limited by the size of the protein as the number of conformations grows exponentially with protein length. In this work, we surmount this size limitation by introducing a novel approximation that treats stretches of 3 or 4 residues as blocks, thus reducing the phase space by nearly three orders of magnitude. The performance of the 'bWSME' model is validated by comparing the predictions for a globular enzyme (RNase H) and a repeat protein (IκBα), against experimental observables and the model without block approximation. Finally, as a proof of concept, we predict the free-energy surface of the 370-residue, multi-domain maltose binding protein and identify an intermediate in good agreement with single-molecule force-spectroscopy measurements. The bWSME model can thus be employed as a quantitative predictive tool to explore the conformational landscapes of large proteins, extract the structural features of putative intermediates, identify parallel folding paths, and thus aid in the interpretation of both ensemble and single-molecule experiments.

PyFolding: Open-Source Graphing, Simulation, and Analysis of the Biophysical Properties of Proteins

For many years, curve-fitting software has been heavily utilized to fit simple models to various types of biophysical data. Although such software packages are easy to use for simple functions, they are often expensive and present substantial impediments to applying more complex models or for the analysis of large data sets. One field that is reliant on such data analysis is the thermodynamics and kinetics of protein folding. Over the past decade, increasingly sophisticated analytical models have been generated, but without simple tools to enable routine analysis. Consequently, users have needed to generate their own tools or otherwise find willing collaborators. Here we present PyFolding, a free, open-source, and extensible Python framework for graphing, analysis, and simulation of the biophysical properties of proteins. To demonstrate the utility of PyFolding, we have used it to analyze and model experimental protein folding and thermodynamic data. Examples include: 1) multiphase kinetic folding fitted to linked equations, 2) global fitting of multiple data sets, and 3) analysis of repeat protein thermodynamics with Ising model variants. Moreover, we demonstrate how PyFolding is easily extensible to novel functionality beyond applications in protein folding via the addition of new models. Example scripts to perform these and other operations are supplied with the software, and we encourage users to contribute notebooks and models to create a community resource. Finally, we show that PyFolding can be used in conjunction with Jupyter notebooks as an easy way to share methods and analysis for publication and among research teams.

Effect of Mutations on the Global and Site-Specific Stability and Folding of an Elementary Protein Structural Motif

Understanding the folding mechanism of proteins requires detailed knowledge of the roles of individual amino acid residues in stabilization of specific elements and local segments of the native structure. Recently, we have utilized the combination of circular dichroism (CD) and site-specific ¹³C isotopically edited infrared spectroscopy (IR) coupled with the Ising-like model for protein folding to map the thermal unfolding at the residue level of a de novo designed helix-turn-helix motif αtα. Here we use the same methodology to study how the sequence of local thermal unfolding is affected by selected mutations introduced into the most and least stable parts of the motif. Seven different mutants of αtα are screened to find substitutions with the most pronounced effects on the overall stability. Subsequently, thermal unfolding of two mutated αtα sequences is studied with site-specific resolution, using four distinct ¹³C isotopologues of each. The data are analyzed with the Ising-like model, which builds on a previous parametrization for the original αtα sequence and tests different ways of incorporating the amino acid substitution. We show that for both more and less stable mutants only the adjustment of all interaction parameters of the model can yield a satisfactory fit to the experimental data. The stabilizing and destabilizing mutations result, respectively, in a similar increase and decrease of the stability of all probed local segments, irrespective of their position with respect to the mutation site. Consequently, the relative order of their unfolding remains essentially unchanged. These results underline the importance of the interconnectivity of the stabilizing interaction network and cooperativity of the protein structure, which is evident even in a small motif with apparently noncooperative, heterogeneous unfolding. Overall, our findings are consistent with the native structure being the dominant factor in determining the folding mechanism, regardless of the details of its overall or local thermodynamic stabilization.

Perspective: Identification of collective variables and metastable states of protein dynamics

The statistical analysis of molecular dynamics simulations requires dimensionality reduction techniques, which yield a low-dimensional set of collective variables (CVs) {x _i } = x that in some sense describe the essential dynamics of the system. Considering the distribution P( x ) of the CVs, the primal goal of a statistical analysis is to detect the characteristic features of P( x ), in particular, its maxima and their connection paths. This is because these features characterize the low-energy regions and the energy barriers of the corresponding free energy landscape ΔG( x ) = -k _B T ln P( x ), and therefore amount to the metastable states and transition regions of the system. In this perspective, we outline a systematic strategy to identify CVs and metastable states, which subsequently can be employed to construct a Langevin or a Markov state model of the dynamics. In particular, we account for the still limited sampling typically achieved by molecular dynamics simulations, which in practice seriously limits the applicability of theories (e.g., assuming ergodicity) and black-box software tools (e.g., using redundant input coordinates). We show that it is essential to use internal (rather than Cartesian) input coordinates, employ dimensionality reduction methods that avoid rescaling errors (such as principal component analysis), and perform density based (rather than k-means-type) clustering. Finally, we briefly discuss a machine learning approach to dimensionality reduction, which highlights the essential internal coordinates of a system and may reveal hidden reaction mechanisms.

Atomic-Level Description of Protein Folding inside the GroEL Cavity

Chaperonins (ubiquitous facilitators of protein folding) sequester misfolded proteins within an internal cavity, thus preventing protein aggregation during the process of refolding. GroEL, a tetradecameric bacterial chaperonin, is one of the most studied chaperonins, but the role of the internal cavity in the refolding process is still unclear. It has been suggested that rather than simply isolating proteins while they refold, the GroEL cavity actively promotes protein folding. A detailed characterization of the folding dynamics and thermodynamics of protein substrates encapsulated within the cavity, however, has been difficult to obtain by experimental means, due to the system's complexity and the many steps in the folding cycle. Here, we examine the influence of the GroEL cavity on protein folding based on the results of unbiased, atomistic molecular dynamics simulations. We first verified that the computational setup, which uses a recently developed state-of-the-art force field that more accurately reproduces the aggregation propensity of unfolded states, could recapitulate the essential structural dynamics of GroEL. In these simulations, the GroEL tetradecamer was highly dynamic, transitioning among states corresponding to most of the structures that have been observed experimentally. We then simulated a small, unfolded protein both in the GroEL cavity and in bulk solution and compared the protein's folding process within these two environments. Inside the GroEL cavity, the unfolded protein interacted strongly with the disordered residues in GroEL's C-terminal tails. These interactions stabilized the protein's unfolded states relative to its compact states and increased the roughness of its folding free-energy surface, resulting in slower folding compared to the rate in solution. For larger proteins, which are more typical GroEL substrates, we speculate that these interactions may allow substrates to more quickly escape kinetic traps associated with compact, misfolded states, thereby actively promoting folding.

Heterogeneity in the Folding of Villin Headpiece Subdomain HP36

Small single domain proteins that fold on the microsecond time scale have been the subject of intense interest as models for probing the complexity of folding energy landscapes. The villin headpiece subdomain (HP36) has been extensively studied because of its simple three helix structure, ultrafast folding lifetime of a few microseconds, and stable native fold. We have previously shown that folding as measured by a single ¹³C═¹⁸O isotopic label on residue A57 in helix 2 occurs at a different rate than that measured by global probes of folding, indicating noncooperative complexity in the folding of HP36. In order to determine whether this complexity reflects intermediates or parallel pathways over a small activation barrier, ¹³C═¹⁸O labels were individually incorporated at six different positions in HP36, including into all 3 helices. The equilibrium thermal unfolding transitions and the folding/unfolding dynamics were monitored using the unique IR signature of the ¹³C═¹⁸O label by temperature dependent FTIR and temperature jump IR spectroscopy, respectively. Equilibrium experiments reveal that the ¹³C═¹⁸O labels at different positions in HP36 show drastic differences in the midpoint of their transitions ( T_m), ranging from 45 to 67 °C. Heterogeneity is also observed in the relaxation kinetics; there are differences in the microsecond phase when different labeled positions are probed. At a final temperature of 45 °C, the relaxation rate for ¹³C═¹⁸O A57 is 2.4e + 05 s^-1 whereas for ¹³C═¹⁸O L69 HP36 the relaxation rate is 5.1e + 05 s^-1, two times faster. The observation of site-dependent midpoints for the equilibrium unfolding transitions and differences in the relaxation rates of the labeled positions enables us to probe the progressive accumulation of the folded structure, providing insight into the microscopic details of the folding mechanism.

Accurate Protein-Folding Transition-Path Statistics from a Simple Free-Energy Landscape

A central goal of protein-folding theory is to predict the stochastic dynamics of transition paths-the rare trajectories that transit between the folded and unfolded ensembles-using only thermodynamic information, such as a low-dimensional equilibrium free-energy landscape. However, commonly used one-dimensional landscapes typically fall short of this aim, because an empirical coordinate-dependent diffusion coefficient has to be fit to transition-path trajectory data in order to reproduce the transition-path dynamics. We show that an alternative, first-principles free-energy landscape predicts transition-path statistics that agree well with simulations and single-molecule experiments without requiring dynamical data as an input. This "topological configuration" model assumes that distinct, native-like substructures assemble on a time scale that is slower than native-contact formation but faster than the folding of the entire protein. Using only equilibrium simulation data to determine the free energies of these coarse-grained intermediate states, we predict a broad distribution of transition-path transit times that agrees well with the transition-path durations observed in simulations. We further show that both the distribution of finite-time displacements on a one-dimensional order parameter and the ensemble of transition-path trajectories generated by the model are consistent with the simulated transition paths. These results indicate that a landscape based on transient folding intermediates, which are often hidden by one-dimensional projections, can form the basis of a predictive model of protein-folding transition-path dynamics.

Switching Protein Conformational Substates by Protonation and Mutation

Protein modules that regulate the availability and conformational status of transcription factors determine the rapidity, duration, and magnitude of cellular response to changing conditions. One such system is the single-gene product Cnu, a four-helix bundle transcription co-repressor, which acts as a molecular thermosensor regulating the expression of virulence genes in enterobacteriaceae through modulation of its native conformational ensemble. Cnu and related genes have also been implicated in pH-dependent expression of virulence genes. We hypothesize that protonation of a conserved buried histidine (H45) in Cnu promotes large electrostatic frustration, thus disturbing the H-NS, a transcription factor, binding face. Spectroscopic and calorimetric methods reveal that H45 exhibits a suppressed p K_a of ∼5.1, the protonation of which switches the conformation to an alternate native ensemble in which the fourth helix is disordered. The population redistribution can also be achieved through a mutation H45V, which does not display any switching behavior at pH values greater than 4. The Wako-Saitô-Muñoz-Eaton (WSME) statistical mechanical model predicts specific differences in the conformations and fluctuations of the fourth and first helices of Cnu determining the observed pH response. We validate these predictions through fluorescence lifetime measurements of a sole tryptophan, highlighting the presence of both native and non-native interactions in the regions adjoining the binding face of Cnu. Our combined experimental-computational study thus shows that Cnu acts both as a thermo- and pH-sensor orchestrated via a subtle but quantifiable balance between the weak packing of a structural element and protonation of a buried histidine that promotes electrostatic frustration.

Protein folding transition path times from single molecule FRET

The transition path is the tiny segment of a single molecule trajectory when the free energy barrier between states is crossed and for protein folding contains all of the information about the self-assembly mechanism. As a first step toward obtaining structural information during the transition path from experiments, single molecule FRET spectroscopy has been used to determine average transition path times from a photon-by-photon analysis of fluorescence trajectories. These results, obtained for several different proteins, have already provided new and demanding tests that support both the accuracy of all-atom molecular dynamics simulations and the basic postulates of energy landscape theory of protein folding.

How Many Protein Sequences Fold to a Given Structure? A Coevolutionary Analysis

Quantifying the relationship between protein sequence and structure is key to understanding the protein universe. A fundamental measure of this relationship is the total number of amino acid sequences that can fold to a target protein structure, known as the "sequence capacity," which has been suggested as a proxy for how designable a given protein fold is. Although sequence capacity has been extensively studied using lattice models and theory, numerical estimates for real protein structures are currently lacking. In this work, we have quantitatively estimated the sequence capacity of 10 proteins with a variety of different structures using a statistical model based on residue-residue co-evolution to capture the variation of sequences from the same protein family. Remarkably, we find that even for the smallest protein folds, such as the WW domain, the number of foldable sequences is extremely large, exceeding the Avogadro constant. In agreement with earlier theoretical work, the calculated sequence capacity is positively correlated with the size of the protein, or better, the density of contacts. This allows the absolute sequence capacity of a given protein to be approximately predicted from its structure. On the other hand, the relative sequence capacity, i.e., normalized by the total number of possible sequences, is an extremely tiny number and is strongly anti-correlated with the protein length. Thus, although there may be more foldable sequences for larger proteins, it will be much harder to find them. Lastly, we have correlated the evolutionary age of proteins in the CATH database with their sequence capacity as predicted by our model. The results suggest a trade-off between the opposing requirements of high designability and the likelihood of a novel fold emerging by chance.

Structure-Based Prediction of Protein-Folding Transition Paths

We propose a general theory to describe the distribution of protein-folding transition paths. We show that transition paths follow a predictable sequence of high-free-energy transient states that are separated by free-energy barriers. Each transient state corresponds to the assembly of one or more discrete, cooperative units, which are determined directly from the native structure. We show that the transition state on a folding pathway is reached when a small number of critical contacts are formed between a specific set of substructures, after which folding proceeds downhill in free energy. This approach suggests a natural resolution for distinguishing parallel folding pathways and provides a simple means to predict the rate-limiting step in a folding reaction. Our theory identifies a common folding mechanism for proteins with diverse native structures and establishes general principles for the self-assembly of polymers with specific interactions.

When fast is better: protein folding fundamentals and mechanisms from ultrafast approaches

Protein folding research stalled for decades because conventional experiments indicated that proteins fold slowly and in single strokes, whereas theory predicted a complex interplay between dynamics and energetics resulting in myriad microscopic pathways. Ultrafast kinetic methods turned the field upside down by providing the means to probe fundamental aspects of folding, test theoretical predictions and benchmark simulations. Accordingly, experimentalists could measure the timescales for all relevant folding motions, determine the folding speed limit and confirm that folding barriers are entropic bottlenecks. Moreover, a catalogue of proteins that fold extremely fast (microseconds) could be identified. Such fast-folding proteins cross shallow free energy barriers or fold downhill, and thus unfold with minimal co-operativity (gradually). A new generation of thermodynamic methods has exploited this property to map folding landscapes, interaction networks and mechanisms at nearly atomic resolution. In parallel, modern molecular dynamics simulations have finally reached the timescales required to watch fast-folding proteins fold and unfold in silico. All of these findings have buttressed the fundamentals of protein folding predicted by theory, and are now offering the first glimpses at the underlying mechanisms. Fast folding appears to also have functional implications as recent results connect downhill folding with intrinsically disordered proteins, their complex binding modes and ability to moonlight. These connections suggest that the coupling between downhill (un)folding and binding enables such protein domains to operate analogically as conformational rheostats.

Imaging Metastable States and Transitions in Proteins by Trajectory Map

It has been a long-standing and intriguing issue to develop robust methods to identify metastable states and interstate transitions from simulations or experimental data to understand the functional conformational changes of proteins. It is usually hard to define the complicated boundaries of the states in the conformational space using most of the existing methods, and they often lead to parameter-sensitive results. Here, we present a new approach, visualized Trajectory Map (vTM), to identify the metastable states and the rare interstate transitions, by considering both the conformational similarity and the temporal successiveness of conformations. The vTM is able to give a nonambiguous description of slow dynamics. The case study of a β-hairpin peptide shows that the vTM can reveal the states and transitions from all-atom MD trajectory data even when a single observable (i.e, one-dimensional reaction coordinate) is used. We also use the vTM to refine the folding/unfolding mechanism of HP35 in explicit water by analyzing a 125 μs all-atom MD trajectory and obtain folding/unfolding rates of about 1/μs, which are in good agreement with the experimental values.

Toward a quantitative description of microscopic pathway heterogeneity in protein folding

Experimentally consistent statistical modeling of protein folding thermodynamics reveals unprecedented complexity with numerous parallel folding routes in five different proteins.

Cooperativity and modularity in protein folding

A simple statistical mechanical model proposed by Wako and Saitô has explained the aspects of protein folding surprisingly well. This model was systematically applied to multiple proteins by Muñoz and Eaton and has since been referred to as the Wako-Saitô-Muñoz-Eaton (WSME) model. The success of the WSME model in explaining the folding of many proteins has verified the hypothesis that the folding is dominated by native interactions, which makes the energy landscape globally biased toward native conformation. Using the WSME and other related models, Saitô emphasized the importance of the hierarchical pathway in protein folding; folding starts with the creation of contiguous segments having a native-like configuration and proceeds as growth and coalescence of these segments. The Φ-values calculated for barnase with the WSME model suggested that segments contributing to the folding nucleus are similar to the structural modules defined by the pattern of native atomic contacts. The WSME model was extended to explain folding of multi-domain proteins having a complex topology, which opened the way to comprehensively understanding the folding process of multi-domain proteins. The WSME model was also extended to describe allosteric transitions, indicating that the allosteric structural movement does not occur as a deterministic sequential change between two conformations but as a stochastic diffusive motion over the dynamically changing energy landscape. Statistical mechanical viewpoint on folding, as highlighted by the WSME model, has been renovated in the context of modern methods and ideas, and will continue to provide insights on equilibrium and dynamical features of proteins.

Characterization of protein folding by a Φ-value calculation with a statistical-mechanical model

The Φ-value analysis approach provides information about transition-state structures along the folding pathway of a protein by measuring the effects of an amino acid mutation on folding kinetics. Here we compared the theoretically calculated Φ values of 27 proteins with their experimentally observed Φ values; the theoretical values were calculated using a simple statistical-mechanical model of protein folding. The theoretically calculated Φ values reflected the corresponding experimentally observed Φ values with reasonable accuracy for many of the proteins, but not for all. The correlation between the theoretically calculated and experimentally observed Φ values strongly depends on whether the protein-folding mechanism assumed in the model holds true in real proteins. In other words, the correlation coefficient can be expected to illuminate the folding mechanisms of proteins, providing the answer to the question of which model more accurately describes protein folding: the framework model or the nucleation-condensation model. In addition, we tried to characterize protein folding with respect to various properties of each protein apart from the size and fold class, such as the free-energy profile, contact-order profile, and sensitivity to the parameters used in the Φ-value calculation. The results showed that any one of these properties alone was not enough to explain protein folding, although each one played a significant role in it. We have confirmed the importance of characterizing protein folding from various perspectives. Our findings have also highlighted that protein folding is highly variable and unique across different proteins, and this should be considered while pursuing a unified theory of protein folding.

Refined Parameterization of Nonbonded Interactions Improves Conformational Sampling and Kinetics of Protein Folding Simulations

Recent advances in computational technology have enabled brute-force molecular dynamics (MD) simulations of protein folding using physics-based molecular force fields. The extensive sampling of protein conformations afforded by such simulations revealed, however, considerable compaction of the protein conformations in the unfolded state, which is inconsistent with experiment. Here, we show that a set of surgical corrections to nonbonded interactions between amine nitrogen-carboxylate oxygen and aliphatic carbon-carbon atom pairs can considerably improve the realism of protein folding simulations. Specifically, we show that employing our corrections in ∼500 μs all-atom replica-exchange MD simulations of the WW domain and villin head piece proteins increases the size of the denatured proteins' conformations and does not destabilize the native conformations of the proteins. In addition to making the folded conformations a global minimum of the respective free energy landscapes at room temperature, our corrections also make the free energy landscape smoother, considerably accelerating the folding kinetics and, hence, reducing the computational expense of a protein folding simulation.

Energetic and topological determinants of a phosphorylation-induced disorder-to-order protein conformational switch

We show that the phosphorylation of 4E-BP2 acts as a triggering event to shape its folding-function landscape that is delicately balanced between conflicting favorable energetics and intrinsically unfavorable topological connectivity. We further provide first evidence that the fitness landscapes of proteins at the threshold of disorder can differ considerably from ordered domains.

Study on the Application of the Combination of TMD Simulation and Umbrella Sampling in PMF Calculation for Molecular Conformational Transitions

Free energy calculations of the potential of mean force (PMF) based on the combination of targeted molecular dynamics (TMD) simulations and umbrella samplings as a function of physical coordinates have been applied to explore the detailed pathways and the corresponding free energy profiles for the conformational transition processes of the butane molecule and the 35-residue villin headpiece subdomain (HP35). The accurate PMF profiles for describing the dihedral rotation of butane under both coordinates of dihedral rotation and root mean square deviation (RMSD) variation were obtained based on the different umbrella samplings from the same TMD simulations. The initial structures for the umbrella samplings can be conveniently selected from the TMD trajectories. For the application of this computational method in the unfolding process of the HP35 protein, the PMF calculation along with the coordinate of the radius of gyration (R_g) presents the gradual increase of free energies by about 1 kcal/mol with the energy fluctuations. The feature of conformational transition for the unfolding process of the HP35 protein shows that the spherical structure extends and the middle α-helix unfolds firstly, followed by the unfolding of other α-helices. The computational method for the PMF calculations based on the combination of TMD simulations and umbrella samplings provided a valuable strategy in investigating detailed conformational transition pathways for other allosteric processes.

Imprints of function on the folding landscape: functional role for an intermediate in a conserved eukaryotic binding protein

In the computational characterization of single domain protein folding, the effective free energies of numerous microstates are projected onto few collective degrees of freedom that in turn serve as well-defined reaction coordinates. In this regard, one-dimensional (1D) free energy profiles are widely used mainly for their simplicity. Since folding and functional landscapes are interlinked, how well can these reduced representations capture the structural and dynamic features of functional states while being simultaneously consistent with experimental observables? We investigate this issue by characterizing the folding of the four-helix bundle bovine acyl-CoA binding protein (bACBP), which exhibits complex equilibrium and kinetic behaviours, employing an Ising-like statistical mechanical model and molecular simulations. We show that the features of the 1D free energy profile are sufficient to quantitatively reproduce multiple experimental observations including millisecond chevron-like kinetics and temperature dependence, a microsecond fast phase, barrier heights, unfolded state movements, the intermediate structure and average ϕ-values. Importantly, we find that the structural features of the native-like intermediate (partial disorder in helix 1) are intricately linked to a unique interplay between packing and electrostatics in this domain. By comparison with available experimental data, we propose that this intermediate determines the promiscuous functional behaviour of bACBP that exhibits broad substrate specificity. Our results present evidence to the possibility of employing the statistical mechanical model and the resulting 1D free energy profile to not just understand folding mechanisms but to even extract features of functionally relevant states and their energetic origins.

Twist-open mechanism of DNA damage recognition by the Rad4/XPC nucleotide excision repair complex

DNA damage repair starts with the recognition of damaged sites from predominantly normal DNA. In eukaryotes, diverse DNA lesions from environmental sources are recognized by the xeroderma pigmentosum C (XPC) nucleotide excision repair complex. Studies of Rad4 (radiation-sensitive 4; yeast XPC ortholog) showed that Rad4 "opens" up damaged DNA by inserting a β-hairpin into the duplex and flipping out two damage-containing nucleotide pairs. However, this DNA lesion "opening" is slow (˜5-10 ms) compared with typical submillisecond residence times per base pair site reported for various DNA-binding proteins during 1D diffusion on DNA. To address the mystery as to how Rad4 pauses to recognize lesions during diffusional search, we examine conformational dynamics along the lesion recognition trajectory using temperature-jump spectroscopy. Besides identifying the ˜10-ms step as the rate-limiting bottleneck towards opening specific DNA site, we uncover an earlier ˜100- to 500-μs step that we assign to nonspecific deformation (unwinding/"twisting") of DNA by Rad4. The β-hairpin is not required to unwind or to overcome the bottleneck but is essential for full nucleotide-flipping. We propose that Rad4 recognizes lesions in a step-wise "twist-open" mechanism, in which preliminary twisting represents Rad4 interconverting between search and interrogation modes. Through such conformational switches compatible with rapid diffusion on DNA, Rad4 may stall preferentially at a lesion site, offering time to open DNA. This study represents the first direct observation, to our knowledge, of dynamical DNA distortions during search/interrogation beyond base pair breathing. Submillisecond interrogation with preferential stalling at cognate sites may be common to various DNA-binding proteins.

Microscopic interpretation of folding ϕ-values using the transition path ensemble

All-atom molecular dynamics simulations now allow us to create movies of proteins folding and unfolding. However, it is difficult to assess the accuracy of the folding mechanisms observed because experiments cannot yet directly resolve events occurring along the transition paths between unfolded and folded states. Protein folding ϕ-values provide residue-resolved information about folding mechanisms by comparing the effects of mutations on folding rates and stability, but determining ϕ-values by separately simulating mutant proteins would be computationally demanding and prone to large statistical errors. Here we use transition path theory to develop a method for computing ϕ-values directly from the transition path ensemble, without the need for additional simulations. This path-based approach uses the full transition path information available from equilibrium folding and unfolding trajectories, or from transition path sampling, and does not require identification of folding transition states. Applying our approach to a set of simulations of 10 small proteins by Shaw and coworkers [Lindorff-Larsen K, Piana S, Dror RO, Shaw DE (2011) Science 334(6055):517-520; Piana S, Lindorff-Larsen K, Shaw DE (2011) Biophys J100(9):L47-L49; and Piana S, Lindorff-Larsen K, Shaw DE (2013) Proc Natl Acad Sci USA 110(15):5915-5920], we find good agreement with experiments in most cases where data are available. We can further resolve the contributions to fractional ϕ-values coming from partial contact formation versus transition path heterogeneity. Although in some cases, there is substantial heterogeneity of folding mechanism, in others, such as Ubiquitin, the mechanism is strongly conserved.

Ising Model Reprogramming of a Repeat Protein's Equilibrium Unfolding Pathway

Repeat proteins are formed from units of 20-40 aa that stack together into quasi one-dimensional non-globular structures. This modular repetitive construction means that, unlike globular proteins, a repeat protein's equilibrium folding and thus thermodynamic stability can be analysed using linear Ising models. Typically, homozipper Ising models have been used. These treat the repeat protein as a series of identical interacting subunits (the repeated motifs) that couple together to form the folded protein. However, they cannot describe subunits of differing stabilities. Here we show that a more sophisticated heteropolymer Ising model can be constructed and fitted to two new helix deletion series of consensus tetratricopeptide repeat proteins (CTPRs). This analysis, showing an asymmetric spread of stability between helices within CTPR ensembles, coupled with the Ising model's predictive qualities was then used to guide reprogramming of the unfolding pathway of a variant CTPR protein. The designed behaviour was engineered by introducing destabilising mutations that increased the thermodynamic asymmetry within a CTPR ensemble. The asymmetry caused the terminal α-helix to thermodynamically uncouple from the rest of the protein and preferentially unfold. This produced a specific, highly populated stable intermediate with a putative dimerisation interface. As such it is the first step in designing repeat proteins with function regulated by a conformational switch.

Scaling Rules for Vibrational Energy Transport in Globular Proteins

Computational studies of vibrational energy flow in biomolecules have to date mapped out transport pathways on a case-by-case basis. To provide a more general approach, we derive scaling rules for vibrational energy transport in a globular protein, which are identified from extensive nonequilibrium molecular dynamics simulations of vibrational energy flow in the villin headpiece subdomain HP36. We parametrize a master equation based on inter-residue, residue-solvent, and heater-residue energy-transfer rates, which closely reproduces the results of the all-atom simulations. From that fit, two scaling rules emerge, one for energy transport along the protein backbone which relies on a diffusion model and another for energy transport between tertiary contacts, which is based on a harmonic model. Requiring only the calculation of mean and variance of relatively few atomic distances, the approach holds the potential to predict the pathways and time scales of vibrational energy flow in large proteins.

Reduction of All-Atom Protein Folding Dynamics to One-Dimensional Diffusion

Theoretical models have often modeled protein folding dynamics as diffusion on a low-dimensional free energy surface, a remarkable simplification. However, the accuracy of such an approximation and the number of dimensions required were not clear. For all-atom folding simulations of ten small proteins in explicit solvent we show that the folding dynamics can indeed be accurately described as diffusion on just a single coordinate, the fraction of native contacts (Q). The diffusion models reproduce both folding rates, and finer details such as transition-path durations and diffusive propagators. The Q-averaged diffusion coefficients decrease with chain length, as anticipated from energy landscape theory. Although the Q-diffusion model does not capture transition-path durations for the protein NuG2, we show that this can be accomplished by designing an improved coordinate Qopt. Overall, one-dimensional diffusion on a suitable coordinate turns out to be a remarkably faithful model for the dynamics of the proteins considered.

Structural origin of slow diffusion in protein folding

Experimental, theoretical, and computational studies of small proteins suggest that interresidue contacts not present in the folded structure play little or no role in the self-assembly mechanism. Non-native contacts can, however, influence folding kinetics by introducing additional local minima that slow diffusion over the global free-energy barrier between folded and unfolded states. Here, we combine single-molecule fluorescence with all-atom molecular dynamics simulations to discover the structural origin for the slow diffusion that markedly decreases the folding rate for a designed α-helical protein. Our experimental determination of transition path times and our analysis of the simulations point to non-native salt bridges between helices as the source, which provides a quantitative glimpse of how specific intramolecular interactions influence protein folding rates by altering dynamics and not activation free energies.

Sequence, structure, and cooperativity in folding of elementary protein structural motifs

Residue-level unfolding of two helix-turn-helix proteins--one naturally occurring and one de novo designed--is reconstructed from multiple sets of site-specific (13)C isotopically edited infrared (IR) and circular dichroism (CD) data using Ising-like statistical-mechanical models. Several model variants are parameterized to test the importance of sequence-specific interactions (approximated by Miyazawa-Jernigan statistical potentials), local structural flexibility (derived from the ensemble of NMR structures), interhelical hydrogen bonds, and native contacts separated by intervening disordered regions (through the Wako-Saitô-Muñoz-Eaton scheme, which disallows such configurations). The models are optimized by directly simulating experimental observables: CD ellipticity at 222 nm for model proteins and their fragments and (13)C-amide I' bands for multiple isotopologues of each protein. We find that data can be quantitatively reproduced by the model that allows two interacting segments flanking a disordered loop (double sequence approximation) and incorporates flexibility in the native contact maps, but neither sequence-specific interactions nor hydrogen bonds are required. The near-identical free energy profiles as a function of the global order parameter are consistent with expected similar folding kinetics for nearly identical structures. However, the predicted folding mechanism for the two motifs is different, reflecting the order of local stability. We introduce free energy profiles for "experimental" reaction coordinates--namely, the degree of local folding as sensed by site-specific (13)C-edited IR, which highlight folding heterogeneity and contrast its overall, average description with the detailed, local picture.