Reproducible protein folding with the stochastic tunneling method

Design the RNA aptamer of PCA3 long non-coding ribonucleic acid by the coarse-grained molecular mechanics

The stochastic tunneling-basin hopping-discrete molecular dynamics (STUN-BH-DMD) method was applied to predict the tertiary structure of the prostate cancer marker PCA3 using two respective secondary structures predicted by the Vienna RNA package and Mathews lab package. The RNA CG force field with the geometrical restraints for maintaining PCA3 secondary structures is used. For each secondary structure, 5000 PCA3 structures were predicted by using 5000 independent initial structures. These structures were then evaluated by a scoring function, considering the contributions from the radius of gyration, contact energy, and surface fraction of complementary nucleotides to ASO683 and ASO735 used in the related experiment. For each secondary structure, the PCA3 structures with the highest three scores were selected for aptamer design and further adsorption simulation. The ASOs complementary to PCA3 surface segments possessing relatively higher RMSF values are selected to be the potential PCA3 aptamers. After the adsorption simulation, the adsorption energies of ASO961, ASO3181, ASO3533, and ASO3595 are higher than or comparable to those of ASO683 and ASO735 used in the experiment. The NEB method was used to obtain MEPs for the adsorption process of all predicted ASOs onto PCA3. The adsorption barriers range between 29 ∼ 39 kcal/mol, while the desorption barriers range between 112 ∼ 352 kcal/mol, indicating these aptamer/PCA3 complexes are very stable. Using PCA3 surface segments with relatively higher RMSF values, longer ASOs can be also obtained and most longer ASOs possess lower binding energy, ranging between -486.1 and -618.2 kcal/mol.Communicated by Ramaswamy H. Sarma.

Exploring the most stable aptamer/target molecule complex by the stochastic tunnelling-basin hopping-discrete molecular dynamics method

The stochastic tunnelling-basin hopping-discrete molecular dynamics (STUN-BH-DMD) method was applied to the search for the most stable biomolecular complexes in water by using the MARTINI coarse-grained (CG) model. The epithelial cell adhesion molecule (EpCAM, PDB code: 4MZV) was used as an EpCAM adaptor for an EpA (Apt_EpA) benchmark target molecule. The effects of two adsorption positions on the EpCAM were analysed, and it is found that the Apt_EpA adsorption configuration located within the EpCAM pocket-like structure is more stable and the energy barrier is lower due to the interaction with water. By the root mean square deviation (RMSD), the configuration of EpCAM in water is more conservative when the Apt_EpA binds to EpCAM by attaching to the pocket space of the EpCAM dimer. For Apt_EpA, the root mean square fluctuation (RMSF) analysis result indicates Nucleobase 1 and Nucleobase 2 display higher flexibility during the CGMD simulation. Finally, from the binding energy contour maps and histogram plots of EpCAM and each Apt_EpA nucleobase, it is clear that the binding energy adsorbed to the pocket-like structure is more continuous than that energy not adsorbed to the pocket-like structure. This study has proposed a new numerical process for applying the STUN-BH-DMD with the CG model, which can reduce computational details and directly find a more stable Apt_EpA/EpCAM complex in water.

Including residual contact information into replica-exchange MD simulations significantly enriches native-like conformations

Proteins are complex biomolecules which perform critical tasks in living organisms. Knowledge of a protein's structure is essential for understanding its physiological function in detail. Despite the incredible progress in experimental techniques, protein structure determination is still expensive, time-consuming, and arduous. That is why computer simulations are often used to complement or interpret experimental data. Here, we explore how in silico protein structure determination based on replica-exchange molecular dynamics (REMD) can benefit from including contact information derived from theoretical and experimental sources, such as direct coupling analysis or NMR spectroscopy. To reflect the influence from erroneous and noisy data we probe how false-positive contacts influence the simulated ensemble. Specifically, we integrate varying numbers of randomly selected native and non-native contacts and explore how such a bias can guide simulations towards the native state. We investigate the number of contacts needed for a significant enrichment of native-like conformations and show the capabilities and limitations of this method. Adhering to a threshold of approximately 75% true-positive contacts within a simulation, we obtain an ensemble with native-like conformations of high quality. We find that contact-guided REMD is capable of delivering physically reasonable models of a protein's structure.

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.

Predicting the Most Stable Aptamer/Target Molecule Complex Configuration Using a Stochastic-Tunnelling Basin-Hopping Discrete Molecular Dynamics Method: A Novel Global Minimum Search Method for a Biomolecule Complex

This study proposed a novel global minimum search method for predicting the most stable biomolecule complex, which combines the strengths of three global minimum search methods (stochastic tunnelling, basin hopping, and discrete molecular dynamics) to efficiently improve the spatial domain search ability of the stochastic tunnelling-basin hopping (STUN-BH) method from our previous study. The epithelial cell adhesion molecule (EpCAM, PDB code: 4MZV) was used as a benchmark target molecule for the EpCAM aptamer EpA (Apt_EpA). For the most stable Apt_EpA/EpCAM complex predicted by our new method, the Apt_EpA was attached to the entangling loop fragments of the two EpCAM molecules with the most Apt_EpA residues. After the Apt_EpA/EpCAM complex had equilibrated with the water environment through a molecular dynamics simulation at 300 K for 10 ns, stable hydrogen bonds formed between the bases of Apt_EpA and EpCAM residues of the secondary structures, which included the alpha helix and beta sheet becoming less stable in the water environment. Those hydrogen bonds formed between the bases of Apt_EpA and EpCAM loop fragment residues remained stable in the water environment.

A Hybrid Hamiltonian for the Accelerated Sampling along Experimental Restraints

In this article, we present an enhanced sampling method based on a hybrid Hamiltonian which combines experimental distance restraints with a bias dependent from multiple path-dependent variables. This simulation method determines the bias-coordinates on the fly and does not require a priori knowledge about reaction coordinates. The hybrid Hamiltonian accelerates the sampling of proteins, and, combined with experimental distance information, the technique considers the restraints adaptively and in dependency of the system's intrinsic dynamics. We validate the methodology on the dipole relaxation of two water models and the conformational landscape of dialanine. Using experimental NMR-restraint data, we explore the folding landscape of the TrpCage mini-protein and in a second example apply distance restraints from chemical crosslinking/mass spectrometry experiments for the sampling of the conformation space of the Killer Cell Lectin-like Receptor Subfamily B Member 1A (NKR-P1A). The new methodology has the potential to adaptively introduce experimental restraints without affecting the conformational space of the system along an ergodic trajectory. Since only a limited number of input- and no-order parameters are required for the setup of the simulation, the method is broadly applicable and has the potential to be combined with coarse-graining methods.

Enriched Conformational Sampling of DNA and Proteins with a Hybrid Hamiltonian Derived from the Protein Data Bank

In this article, we present a method for the enhanced molecular dynamics simulation of protein and DNA systems called potential of mean force (PMF)-enriched sampling. The method uses partitions derived from the potentials of mean force, which we determined from DNA and protein structures in the Protein Data Bank (PDB). We define a partition function from a set of PDB-derived PMFs, which efficiently compensates for the error introduced by the assumption of a homogeneous partition function from the PDB datasets. The bias based on the PDB-derived partitions is added in the form of a hybrid Hamiltonian using a renormalization method, which adds the PMF-enriched gradient to the system depending on a linear weighting factor and the underlying force field. We validated the method using simulations of dialanine, the folding of TrpCage, and the conformational sampling of the Dickerson⁻Drew DNA dodecamer. Our results show the potential for the PMF-enriched simulation technique to enrich the conformational space of biomolecules along their order parameters, while we also observe a considerable speed increase in the sampling by factors ranging from 13.1 to 82. The novel method can effectively be combined with enhanced sampling or coarse-graining methods to enrich conformational sampling with a partition derived from the PDB.

Replica-exchange-with-tunneling for fast exploration of protein landscapes

While the use of replica-exchange molecular dynamics in protein simulations has become ubiquitous, its utility is limited in many practical applications. We propose to overcome some shortcomings that hold back its use in settings such as multi-scale or explicit solvent simulations by integrating ideas of hybrid MC/MD into the replica-exchange protocol. This Replica-Exchange-with-Tunneling method is tested by simulating the Trp-cage protein, a system often used in molecular biophysics for testing sampling techniques.

PACE Force Field for Protein Simulations. 2. Folding Simulations of Peptides

We present the application of our recently developed PACE force field to the folding of peptides. These peptides include α-helical (AK17 and Fs), β-sheet (GB1m2 and Trpzip2), and mixed helical/coil (Trp-cage) peptides. With replica exchange molecular dynamics (REMD), our force field can fold the five peptides into their native structures while maintaining their stabilities reasonably well. Our force field is also able to capture important thermodynamic features of the five peptides that have been observed in previous experimental and computational studies, such as different preferences for a helix-turn-helix topology for AK17 and Fs, the relative contribution of four hydrophobic side chains of GB1p to the stability of β-hairpin, and the distinct role of a hydrogen bond involving Trp-Hε and a D9/R16 salt bridge in stabilizing the Trp-cage native structure. Furthermore, multiple folding and unfolding events are observed in our microsecond-long normal MD simulations of AK17, Trpzip2, and Trp-cage. These simulations provide mechanistic information such as a "zip-out" pathway of the folding mechanism of Trpzip2 and the folding times of AK17 and Trp-cage, which are estimated to be about 51 ± 43 ns and 270 ± 110 ns, respectively. A 600 ns simulation of the peptides can be completed within one day. These features of our force field are potentially applicable to the study of thermodynamics and kinetics of real protein systems.

Coevolutionary Landscape Inference and the Context-Dependence of Mutations in Beta-Lactamase TEM-1

The quantitative characterization of mutational landscapes is a task of outstanding importance in evolutionary and medical biology: It is, for example, of central importance for our understanding of the phenotypic effect of mutations related to disease and antibiotic drug resistance. Here we develop a novel inference scheme for mutational landscapes, which is based on the statistical analysis of large alignments of homologs of the protein of interest. Our method is able to capture epistatic couplings between residues, and therefore to assess the dependence of mutational effects on the sequence context where they appear. Compared with recent large-scale mutagenesis data of the beta-lactamase TEM-1, a protein providing resistance against beta-lactam antibiotics, our method leads to an increase of about 40% in explicative power as compared with approaches neglecting epistasis. We find that the informative sequence context extends to residues at native distances of about 20 Å from the mutated site, reaching thus far beyond residues in direct physical contact.

A fluorescence polarization assay for the experimental validation of an in silico model of the chemokine CXCL8 binding to receptor-derived peptides

Peptide based inhibitors of protein-protein interactions are of great interest in proteomics, structural biology and medicinal chemistry. Optimized inhibitors can be designed by experimental approaches or by computational prediction. Ideally, computational models are adjusted to the peptide-protein complex of interest according to experimental data obtained in specific binding experiments. The chemokine CXCL8 (interleukin-8) is an interesting target for drug discovery due to its role in inflammatory diseases. Given the available structural data and information on its receptor interactions it constitutes a basis for the rational design of inhibitor peptides. Starting from the reported structure of CXCL8 in complex with a peptide derived from its receptor CXCR1 we developed a computational docking procedure to estimate the changes in binding energy as a function of individual amino acid exchanges. This indicates whether the respective amino acid residue must be preserved or can be substituted to maintain or improve affinity, respectively. To validate and improve the assumptions made in this docking simulation we established a fluorescence polarization assay for receptor-derived peptides binding to CXCL8. A peptide library was tested comprising selected mutants characterized by docking simulations. A number of predictions regarding electrostatic interactions were confirmed by these experiments and it was revealed that the model needed to be corrected for backbone flexibility. Therefore, the assay presented here is a promising tool to systematically improve the computational model by iterative cycles of modeling, experimental validation and refinement of the algorithm, leading to a more reliable model and peptides with improved affinity.

Protein structure prediction: assembly of secondary structure elements by basin-hopping

The prediction of protein tertiary structure from primary structure remains a challenging task. One possible approach to this problem is the application of basin-hopping global optimization combined with an all-atom force field. In this work, the efficiency of basin-hopping is improved by introducing an approach that derives tertiary structures from the secondary structure assignments of individual residues. This approach is termed secondary-to-tertiary basin-hopping and benchmarked for three miniproteins: trpzip, trp-cage and ER-10. For each of the three miniproteins, the secondary-to-tertiary basin-hopping approach successfully and reliably predicts their three-dimensional structure. When it is applied to larger proteins, correctly folded structures are obtained. It can be concluded that the assembly of secondary structure elements using basin-hopping is a promising tool for de novo protein structure prediction.

Using D-Amino Acids to Delineate the Mechanism of Protein Folding: Application to Trp-cage

Using the miniprotein Trp-cage as a model, we show that D-amino acids can be used to facilitate the delineation of protein folding mechanism. Specifically, we study the folding-unfolding kinetics of three Trp-cage mutants where the native glycine residue near the C-terminus of the α-helix is replaced by a D-amino acid. A previous study showed that these mutations increase the Trp-cage stability, due to a terminal capping effect. Our results show that the stabilizing effect of D-asparagine and D-glutamine originates almost exclusively from a decrease in the unfolding rate, while the D-alanine mutation results in a similar decrease in the unfolding rate, but it also increases the folding rate. Together, these results support a folding mechanism wherein the α-helix formation in the transition state is nucleated at the N-terminus, whereas those long-range native interactions stabilizing this helix are developed at the downhill side of the folding free energy barrier.

Protein structure prediction using global optimization by basin-hopping with NMR shift restraints

Computational methods that utilize chemical shifts to produce protein structures at atomic resolution have recently been introduced. In the current work, we exploit chemical shifts by combining the basin-hopping approach to global optimization with chemical shift restraints using a penalty function. For three peptides, we demonstrate that this approach allows us to find near-native structures from fully extended structures within 10,000 basin-hopping steps. The effect of adding chemical shift restraints is that the α and β secondary structure elements form within 1000 basin-hopping steps, after which the orientation of the secondary structure elements, which produces the tertiary contacts, is driven by the underlying protein force field. We further show that our chemical shift-restraint BH approach also works for incomplete chemical shift assignments, where the information from only one chemical shift type is considered. For the proper implementation of chemical shift restraints in the basin-hopping approach, we determined the optimal weight of the chemical shift penalty energy with respect to the CHARMM force field in conjunction with the FACTS solvation model employed in this study. In order to speed up the local energy minimization procedure, we developed a function, which continuously decreases the width of the chemical shift penalty function as the minimization progresses. We conclude that the basin-hopping approach with chemical shift restraints is a promising method for protein structure prediction.

Further optimization of a hybrid united-atom and coarse-grained force field for folding simulations: Improved backbone hydration and interactions between charged side chains

PACE, a hybrid force field which couples united-atom protein models with coarse-grained (CG) solvent, has been further optimized, aiming to improve itse ciency for folding simulations. Backbone hydration parameters have been re-optimized based on hydration free energies of polyalanyl peptides through atomistic simulations. Also, atomistic partial charges from all-atom force fields were combined with PACE in order to provide a more realistic description of interactions between charged groups. Using replica exchange molecular dynamics (REMD), ab initio folding using the new PACE has been achieved for seven small proteins (16 - 23 residues) with different structural motifs. Experimental data about folded states, such as their stability at room temperature, melting point and NMR NOE constraints, were also well reproduced. Moreover, a systematic comparison of folding kinetics at room temperature has been made with experiments, through standard MD simulations, showing that the new PACE may speed up the actual folding kinetics 5-10 times. Together with the computational speedup benefited from coarse-graining, the force field provides opportunities to study folding mechanisms. In particular, we used the new PACE to fold a 73-residue protein, 3D, in multiple 10 - 30 μs simulations, to its native states (C(α) RMSD ~ 0.34 nm). Our results suggest the potential applicability of the new PACE for the study of folding and dynamics of proteins.

Mirror images as naturally competing conformations in protein folding

Evolution has selected a protein's sequence to be consistent with the native state geometry, as this configuration must be both thermodynamically stable and kinetically accessible to prevent misfolding and loss of function. In simple protein geometries, such as coiled-coil helical bundles, symmetry produces a competing, globally different, near mirror image with identical secondary structure and similar native contact interactions. Experimental techniques such as circular dichroism, which rely on probing secondary structure content, cannot readily distinguish these folds. Here, we want to clarify whether the native fold and mirror image are energetically competitive by investigating the free energy landscape of three-helix bundles. To prevent a bias from a specific computational approach, the present study employs the structure prediction forcefield PFF01/02, explicit solvent replica exchange molecular dynamics (REMD) with the Amber94 forcefield, and structure-based simulations based on energy landscape theory. We observe that the native fold and its mirror image have a similar enthalpic stability and are thermodynamically competitive. There is evidence that the mirror fold has faster folding kinetics and could function as a kinetic trap. All together, our simulations suggest that mirror images might not just be a computational annoyance but are competing folds that might switch depending on environmental conditions or functional considerations.

Velocity-scaling optimized replica exchange molecular dynamics of proteins in a hybrid explicit/implicit solvent

We propose a scheme for replica exchange molecular dynamics of proteins in explicit solvent that minimizes the number of required replicas using velocity rescaling. Our approach relies on a hybrid method where the protein evolves at each temperature in an explicit solvent, but replica exchange moves utilize an implicit solvent term. The two terms are coupled through the velocity rescaling. We test the efficiency of this approach for a common test case, the trp-cage protein.

A WAVELET APPROACH FOR THE ANALYSIS OF FOLDING TRAJECTORY OF PROTEIN TRP-CAGE

Understanding how protein folds into a functional native structure is arguably one of the most challenging problems remaining in computational biology. Currently, the protein folding mechanism is often characterized by calculating the free energy landscape in terms of various reaction coordinates such as the fraction of native contacts, the radius of gyration, the RMS deviation from the native and so on. In this paper, we present a wavelet approach towards understanding the global state changes during protein folding. The approach is based on the wavelet analysis on the trajectories of various reaction coordinates to identify the significant intermediate states or structural motifs in the folding process. We demonstrate through an example protein Trp-cage that this approach extracts crucial information about protein folding intermediate states as well as the time correlation among these states. Furthermore, the current approach reveals a meaningful structural pattern that had been overlooked in previous works, which provides a better understanding of the folding mechanism as well as the limitation of the current force fields.

Influence of water-protein hydrogen bonding on the stability of Trp-cage miniprotein. A comparison between the TIP3P and TIP4P-Ew water models

We report extensive replica exchange molecular dynamics (REMD) simulations on the folding/unfolding equilibrium of Trp-cage miniprotein using the Amber ff99SB all atom forcefield and TIP3P and TIP4P-Ew explicit water solvent models. REMD simulation-lengths in the 500 ns to the microsecond regime per replica are required to adequately sample the folding/unfolding equilibrium. We observe that this equilibrium is significantly affected by the choice of the water model. Compared with experimental data, simulations using the TIP3P solvent describe the stability of the Trp-cage quite realistically, providing a melting point which is just a few Kelvins above the experimental transition temperature of 317 K. The TIP4P-Ew model shifts the equilibrium towards the unfolded state and lowers the free energy of unfolding by about 3 kJ mol(-1) at 280 K, demonstrating the need to fine-tune the protein-forcefield depending on the chosen water model. We report evidence that the main difference between the two water models is mostly due to the different solvation of polar groups of the peptide. The unfolded state of the Trp-cage is stabilized by an increasing number of hydrogen bonds, destabilizing the α-helical part of the molecule and opening the R-D salt bridge. By reweighting the strength of solvent-peptide hydrogen bonds by adding a hydrogen bond square well potential, we can fully recover the effect of the different water models and estimate the shift in population as due to a difference in hydrogen bond-strength of about 0.4 kJ mol(-1) per hydrogen bond.

Optimal salt bridge for Trp-cage stabilization

Gai and co-workers [Bunagan, M. R., et al. (2006) J. Phys. Chem. B 110, 3759-3763] reported computational design studies suggesting that a D9E mutation would stabilize the Trp-cage. Experimental studies for this mutation were reported in 2008 [Hudaky, P., et al. (2008) Biochemistry 47, 1007-1016]; the authors suggested that [D9E]-TC5b presented a more compact and melting resistant structure because of the "optimal distance between the two sides of the molecule". Nonetheless, the authors reported essentially the same circular dichroism (CD) melting temperature, 38 ± 0.3 °C, for TC5b and its [D9E] mutant. In this study, a more stable Trp-cage, DAYAQ WLKDG GPSSG RPPPS, was examined by nuclear magnetic resonance and CD with the following mutations: [D9E], [D9R,R16E], [R16O], [D9E,R16O], [R16K], and [D9E,R16K]. Of these, the [D9E] mutant displayed the smallest acidification-induced change in the apparent T(m). In analogy to the prior study, the CD melts of TC10b and its [D9E] mutant were, however, very similar; all of the other mutations were significantly fold destabilizing by all measures. A detailed analysis indicates that the original D9-R16 salt bridge is optimal with regard to fold cooperativity and fold stabilization. Evidence of salt bridge formation is also provided for a swapped pair, the [D9R,R16E] mutant. Model systems reveal that an ionized aspartate at the C-terminus of a helix significantly decreases intrinsic helicity, a requirement for Trp-cage fold stability. The CD evidence that was cited as supporting increased fold stability for [D9E]-TC5b at higher temperatures appears to be a reflection of increased helix stability in both the folded and unfolded states rather than a more favorable salt bridge. Our study also provides evidence of other Trp-cage stabilizing roles of the R16 side chain.

Balance between alpha and beta structures in ab initio protein folding

Despite initial successes in folding of proteins by molecular simulation, it is becoming increasingly evident that current energy functions (force fields) tend to favor either alpha or beta secondary structure, such that the choice of force field is governed by the structural class of the protein. Here, we study the folding of peptides with either predominantly alpha (Trp cage) or beta (GB1 hairpin) structure with a modified version of the Amber ff03 force field, optimized to reproduce structural propensity in a helix-forming peptide. Using extensive replica exchange molecular dynamics simulations starting from completely unfolded configurations, we obtain the correct folded structure for each peptide, in close agreement with the experimental native structure (<1.5 A all-atom root-mean-square deviation). We obtain converged equilibrium distributions, with folded populations at standard conditions (approximately 300 K), in remarkable accord with experiment. Further comparison to experimental data from NMR spectroscopy and FRET suggests that although the folded structures are accurately reproduced, the unfolded state remains too structured and compact. Our results suggest that the backbone correction results in a force field that is transferable to the folding of proteins from different structural classes.

Modelling proteins: conformational sampling and reconstruction of folding kinetics

In the last decades biomolecular simulation has made tremendous inroads to help elucidate biomolecular processes in-silico. Despite enormous advances in molecular dynamics techniques and the available computational power, many problems involve long time scales and large-scale molecular rearrangements that are still difficult to sample adequately. In this review we therefore summarise recent efforts to fundamentally improve this situation by decoupling the sampling of the energy landscape from the description of the kinetics of the process. Recent years have seen the emergence of many advanced sampling techniques, which permit efficient characterisation of the relevant family of molecular conformations by dispensing with the details of the short-term kinetics of the process. Because these methods generate thermodynamic information at best, they must be complemented by techniques to reconstruct the kinetics of the process using the ensemble of relevant conformations. Here we review recent advances for both types of methods and discuss their perspectives to permit efficient and accurate modelling of large-scale conformational changes in biomolecules. This article is part of a Special Issue entitled: Protein Dynamics: Experimental and Computational Approaches.

Ab Initio Protein Folding Using a Cooperative Swarm of Molecular Dynamics Trajectories

The use of atomistic simulation techniques to directly resolve the protein tertiary structure from the primary amino acid sequence is hindered by the rough topology of the protein free energy surface and the resulting simulation time scales required. We explore here the use of a molecular dynamics technique based on swarm intelligence to identify the native states of two peptides and a Trp-cage miniprotein. In all cases, the presence of cooperative swarm interactions significantly enhanced the efficiency of molecular dynamics simulations in predicting the native conformation.

Robustness and generalization of structure-based models for protein folding and function

Functional dynamics of native proteins share the energy landscape that guides folding into the native state. Folding simulations of structure-based protein models, using an minimally frustrated energy landscape dominated by native interactions, can describe the geometrical aspects of the folding mechanism. Technical limitations imposed by the fixed shape of conventional contact potentials are a key obstacle toward advanced applications of structure-based models like allostery or ligand binding, which require multiple stable conformations. Generalizations of existing models, commonly using Lennard-Jones-like potentials, lead to inevitable clashes between their repulsive branches. To resolve these challenges, a new contact potential is developed that combines an attractive part based on Gaussians with a separate repulsive term allowing flexibility for adjustments of the potential shape. With this new model multiple minima for studies of functional transitions can be introduced easily and consistently. A sensitivity analysis for five small proteins confirms the robust behavior of structure-based models with our adaptable potential and explores their capacity for quantitative adjustment of the folding thermodynamics. We demonstrate its ability to incorporate alternative contact distances in simulations of structural transitions for the well-studied ROP dimer. Individual contact pairs can switch between distinct states to match the competing syn and anti structures. The flexibility of the new potential facilitates advanced uses of structure-based models. Depending on the application, features can be chosen from physical considerations or to match experiments. Generalized models can be built from multiple structures to study structural transitions or effects of disorder.

High-resolution protein complexes from integrating genomic information with molecular simulation

Bacteria use two-component signal transduction systems (TCS) extensively to sense and react to external stimuli. In these, a membrane-bound sensor histidine kinase (SK) autophosphorylates in response to an environmental stimulus and transfers the phosphoryl group to a transcription factor/response regulator (RR) that mediates the cellular response. The complex between these two proteins is ruled by transient interactions, which provides a challenge to experimental structure determination techniques. The functional and structural homolog of an SK/RR pair Spo0B/Spo0F, however, has been structurally resolved. Here, we describe a method capable of generating structural models of such transient protein complexes. By using existing structures of the individual proteins, our method combines bioinformatically derived contact residue information with molecular dynamics simulations. We find crystal resolution accuracy with existing crystallographic data when reconstituting the known system Spo0B/Spo0F. Using this approach, we introduce a complex structure of TM0853/TM0468 as an exemplary SK/RR TCS, consistent with all experimentally available data.

Microcanonical replica exchange molecular dynamics simulation of proteins

We present microcanonical replica exchange molecular dynamics simulations as an alternative to canonical ones. Its advantage is the easily tunable high acceptance rate for replica exchange. We present the theory, comment on its actual implementation, and demonstrate its application for a common test case, the trp-cage protein.

A kinetic model of trp-cage folding from multiple biased molecular dynamics simulations

Trp-cage is a designed 20-residue polypeptide that, in spite of its size, shares several features with larger globular proteins.Although the system has been intensively investigated experimentally and theoretically, its folding mechanism is not yet fully understood. Indeed, some experiments suggest a two-state behavior, while others point to the presence of intermediates. In this work we show that the results of a bias-exchange metadynamics simulation can be used for constructing a detailed thermodynamic and kinetic model of the system. The model, although constructed from a biased simulation, has a quality similar to those extracted from the analysis of long unbiased molecular dynamics trajectories. This is demonstrated by a careful benchmark of the approach on a smaller system, the solvated Ace-Ala3-Nme peptide. For theTrp-cage folding, the model predicts that the relaxation time of 3100 ns observed experimentally is due to the presence of a compact molten globule-like conformation. This state has an occupancy of only 3% at 300 K, but acts as a kinetic trap.Instead, non-compact structures relax to the folded state on the sub-microsecond timescale. The model also predicts the presence of a state at Calpha-RMSD of 4.4 A from the NMR structure in which the Trp strongly interacts with Pro12. This state can explain the abnormal temperature dependence of the Pro12-delta3 and Gly11-alpha3 chemical shifts. The structures of the two most stable misfolded intermediates are in agreement with NMR experiments on the unfolded protein. Our work shows that, using biased molecular dynamics trajectories, it is possible to construct a model describing in detail the Trp-cage folding kinetics and thermodynamics in agreement with experimental data.

A free-energy approach for all-atom protein simulation

All-atom free-energy methods offer a promising alternative to kinetic molecular mechanics simulations of protein folding and association. Here we report an accurate, transferable all-atom biophysical force field (PFF02) that stabilizes the native conformation of a wide range of proteins as the global optimum of the free-energy landscape. For 32 proteins of the ROSETTA decoy set and six proteins that we have previously folded with PFF01, we find near-native conformations with an average backbone RMSD of 2.14 A to the native conformation and an average Z-score of -3.46 to the corresponding decoy set. We used nonequilibrium sampling techniques starting from completely extended conformations to exhaustively sample the energy surface of three nonhomologous hairpin-peptides, a three-stranded beta-sheet, the all-helical 40 amino-acid HIV accessory protein, and a zinc-finger beta beta alpha motif, and find near-native conformations for the minimal energy for each protein. Using a massively parallel evolutionary algorithm, we also obtain a near-native low-energy conformation for the 54 amino-acid engrailed homeodomain. Our force field thus stabilized near-native conformations for a total of 20 proteins of all structure classes with an average RMSD of only 3.06 A to their respective experimental conformations.

An all-atom structure-based potential for proteins: bridging minimal models with all-atom empirical forcefields

Protein dynamics take place on many time and length scales. Coarse-grained structure-based (Go) models utilize the funneled energy landscape theory of protein folding to provide an understanding of both long time and long length scale dynamics. All-atom empirical forcefields with explicit solvent can elucidate our understanding of short time dynamics with high energetic and structural resolution. Thus, structure-based models with atomic details included can be used to bridge our understanding between these two approaches. We report on the robustness of folding mechanisms in one such all-atom model. Results for the B domain of Protein A, the SH3 domain of C-Src Kinase, and Chymotrypsin Inhibitor 2 are reported. The interplay between side chain packing and backbone folding is explored. We also compare this model to a C(alpha) structure-based model and an all-atom empirical forcefield. Key findings include: (1) backbone collapse is accompanied by partial side chain packing in a cooperative transition and residual side chain packing occurs gradually with decreasing temperature, (2) folding mechanisms are robust to variations of the energetic parameters, (3) protein folding free-energy barriers can be manipulated through parametric modifications, (4) the global folding mechanisms in a C(alpha) model and the all-atom model agree, although differences can be attributed to energetic heterogeneity in the all-atom model, and (5) proline residues have significant effects on folding mechanisms, independent of isomerization effects. Because this structure-based model has atomic resolution, this work lays the foundation for future studies to probe the contributions of specific energetic factors on protein folding and function.

Replica exchange molecular dynamics simulations of coarse-grained proteins in implicit solvent

Current approaches aimed at determining the free energy surface of all-atom medium-size proteins in explicit solvent are slow and are not sufficient to converge to equilibrium properties. To ensure a proper sampling of the configurational space, it is preferable to use reduced representations such as implicit solvent and/or coarse-grained protein models, which are much lighter computationally. Each model must be verified, however, to ensure that it can recover experimental structures and thermodynamics. Here we test the coarse-grained implicit solvent OPEP model with replica exchange molecular dynamics (REMD) on six peptides ranging in length from 10 to 28 residues: two alanine-based peptides, the second beta-hairpin from protein G, the Trp-cage and zinc-finger motif, and a dimer of a coiled coil peptide. We show that REMD-OPEP recovers the proper thermodynamics of the systems studied, with accurate structural description of the beta-hairpin and Trp-cage peptides (within 1-2 A from experiments). The light computational burden of REMD-OPEP, which enables us to generate many hundred nanoseconds at each temperature and fully assess convergence to equilibrium ensemble, opens the door to the determination of the free energy surface of larger proteins and assemblies.

Computing the stability diagram of the Trp-cage miniprotein

We report molecular dynamics simulations of the equilibrium folding/unfolding thermodynamics of an all-atom model of the Trp-cage miniprotein in explicit solvent. Simulations are used to sample the folding/unfolding free energy difference and its derivatives along 2 isochores. We model the DeltaG(u)(P,T) landscape using the simulation data and propose a stability diagram model for Trp-cage. We find the proposed diagram to exhibit features similar to globular proteins with increasing hydrostatic pressure destabilizing the native fold. The observed energy differences DeltaE(u) are roughly linearly temperature-dependent and approach DeltaE(u) = 0 with decreasing temperature, suggesting that the system approached the region of cold denaturation. In the low-temperature denatured state, the native helical secondary structure elements are largely preserved, whereas the protein conformation changes to an "open-clamp" configuration. A tighter packing of water around nonpolar sites, accompanied by an increasing solvent-accessible surface area of the unfolded ensemble, seems to stabilize the unfolded state at elevated pressures.

Conformational landscape of the HIV-V3 hairpin loop from all-atom free-energy simulations

Small beta hairpins have many distinct biological functions, including their involvement in chemokine and viral receptor recognition. The relevance of structural similarities between different hairpin loops with near homologous sequences is not yet understood, calling for the development of methods for de novo hairpin structure prediction and simulation. De novo folding of beta strands is more difficult than that of helical proteins because of nonlocal hydrogen bonding patterns that connect amino acids that are distant in the amino acid sequence and there is a large variety of possible hydrogen bond patterns. Here we use a greedy version of the basin hopping technique with our free-energy forcefield PFF02 to reproducibly and predictively fold the hairpin structure of a HIV-V3 loop. We performed 20 independent basin hopping runs for 500 cycles corresponding to 7.4 x 10(7) energy evaluations each. The lowest energy structure found in the simulation has a backbone root mean square deviation (bRMSD) of only 2.04 A to the native conformation. The lowest 9 out of the 20 simulations converged to conformations deviating less than 2.5 A bRMSD from native.

The Trp-cage: optimizing the stability of a globular miniprotein

The Trp-cage, as the smallest miniprotein, remains the subject of numerous computational and experimental studies of protein folding dynamics and pathways. The original Trp-cage (NLYIQWLKDGGPSSGRPPPS, Tm = 42 degrees C) can be significantly stabilized by mutations; melting points as high as 64 degrees C are reported. In helical portions of the structure, each allowed replacement of Leu, Ile, Lys or Ser residues by Ala results in a 1.5 (+/-0.35) kJ/mol fold stabilization. No changes in structure or fluxionality of the core results upon stabilization. Contrary to the initial hypothesis, specific Pro/Trp interactions are not essential for core formation. The entropic advantage of Pro versus Ala (DeltaDeltaS(U) = 11 +/- 2 J/mol K) was measured at the solvent-exposed P17 site. Pro-Ala mutations at two of the three prolines (P12 and P18) that encage the indole ring result in less fold destabilization (2.3-3.4 kJ/mol). However, a P19A mutation reduces fold stability by 16 kJ/mol reflecting a favorable Y3/P19 interaction as well as Trp burial. The Y3/P19 hydrophobic staple interaction defines the folding motif as an 18-residue unit. Other stabilizing features that have been identified include a solvent-exposed Arg/Asp salt bridge (3.4-6 kJ/mol) and a buried H-bonded Ser side chain ( approximately 10 kJ/mol).

Mutations as trapdoors to two competing native conformations of the Rop-dimer

Conformational transitions play a central role in regulating protein function. Structure-based models with multiple basins have been used to understand the mechanisms governing these transitions. A model able to accommodate multiple folding basins is proposed to explore the mutational effects in the folding of the Rop-dimer (Rop). In experiments, Rop mutants show unusually strong increases in folding rates with marginal effects on stability. We investigate the possibility of two competing conformations representing a parallel (P) and the wild-type antiparallel (AP) arrangement of the monomers as possible native conformations. We observe occupation of both distinct states and characterize the transition pathways. An interesting observation from the simulations is that, for equivalent energetic bias, the transition to the P basin (non-wild-type basin) shows a lower free-energy barrier. Thus, the rapid kinetics observed in experiments appear to be the result of two competing states with different kinetic behavior, triggered upon mutation by the opening of a trapdoor arising from Rop's symmetric structure. The general concept of having competing conformations for the native state goes beyond explaining Rop's mutational behaviors and can be applied to other systems. A switch between competing native structures might be triggered by external factors to allow, for example, allosteric control or signaling.