Folding a 20 amino acid αβ peptide with the diffusion process-controlled Monte Carlo method

Multi-scale simulations of biological systems using the OPEP coarse-grained model

Biomolecules are complex machines that are optimized by evolution to properly fulfill or contribute to a variety of biochemical tasks in the cellular environment. Computer simulations based on quantum mechanics and atomistic force fields have been proven to be a powerful microscope for obtaining valuable insights into many biological, physical, and chemical processes. Many interesting phenomena involve, however, a time scale and a number of degrees of freedom, notably if crowding is considered, that cannot be explored at an atomistic resolution. To bridge the gap between reality and simulation, many different advanced computational techniques and coarse-grained (CG) models have been developed. Here, we report some applications of the CG OPEP protein model to amyloid fibril formation, the response of catch-bond proteins to two types of fluid flow, and interactive simulations to fold peptides with well-defined 3D structures or with intrinsic disorder.

High-Resolution Structures of the Amyloid-β 1-42 Dimers from the Comparison of Four Atomistic Force Fields

The dimer of the amyloid-β peptide Aβ of 42 residues is the smallest toxic species in Alzheimer's disease, but its equilibrium structures are unknown. Here we determined the equilibrium ensembles generated by the four atomistic OPLS-AA, CHARMM22*, AMBER99sb-ildn, and AMBERsb14 force fields with the TIP3P water model. On the basis of 144 μs replica exchange molecular dynamics simulations (with 750 ns per replica), we find that the four force fields lead to random coil ensembles with calculated cross-collision sections, hydrodynamics properties, and small-angle X-ray scattering profiles independent of the force field. There are, however, marked differences in secondary structure, with the AMBERsb14 and CHARMM22* ensembles overestimating the CD-derived helix content, and the OPLS-AA and AMBER99sb-ildn secondary structure contents in agreement with CD data. Also the intramolecular beta-hairpin content spanning residues 17-21 and 30-36 varies between 1.5% and 13%. Overall, there are significant differences in tertiary and quaternary conformations among all force fields, and the key finding, irrespective of the force field, is that the dimer is stabilized by nonspecific interactions, explaining therefore its possible transient binding to multiple cellular partners and, in part, its toxicity.

Lattice model for amyloid peptides: OPEP force field parametrization and applications to the nucleus size of Alzheimer's peptides

Coarse-grained protein lattice models approximate atomistic details and keep the essential interactions. They are, therefore, suitable for capturing generic features of protein folding and amyloid formation at low computational cost. As our aim is to study the critical nucleus sizes of two experimentally well-characterized peptide fragments Aβ16-22 and Aβ37-42 of the full length Aβ1-42 Alzheimer's peptide, it is important that simulations with the lattice model reproduce all-atom simulations. In this study, we present a comprehensive force field parameterization based on the OPEP (Optimized Potential for Efficient protein structure Prediction) force field for an on-lattice protein model, which incorporates explicitly the formation of hydrogen bonds and directions of side-chains. Our bottom-up approach starts with the determination of the best lattice force parameters for the Aβ16-22 dimer by fitting its equilibrium parallel and anti-parallel β-sheet populations to all-atom simulation results. Surprisingly, the calibrated force field is transferable to the trimer of Aβ16-22 and the dimer and trimer of Aβ37-42. Encouraged by this finding, we characterized the free energy landscapes of the two decamers. The dominant structure of the Aβ16-22 decamer matches the microcrystal structure. Pushing the simulations for aggregates between 4-mer and 12-mer suggests a nucleus size for fibril formation of 10 chains. In contrast, the Aβ37-42 decamer is largely disordered with mixed by parallel and antiparallel chains, suggesting that the nucleus size is >10 peptides. Our refined force field coupled to this on-lattice model should provide useful insights into the critical nucleation number associated with neurodegenerative diseases.

Folding Atomistic Proteins in Explicit Solvent Using Simulated Tempering

Following a previous report on a coarse-grained protein model in implicit solvent, we applied simulated tempering (ST) with on-the-fly Helmholtz free energy (weight factors) determination to the folding or aggregation of seven proteins with the CHARMM, OPLS, and AMBER protein, and the SPC and TIP3P water force fields. For efficiency and reliability, we also performed replica exchange molecular dynamics (REMD) simulations on the alanine di- and deca-peptide, and the dimer of the Aβ16-22 Alzheimer's fragment, and used experimental data and previous simulation results on the chignolin, beta3s, Trp-cage, and WW domain peptides of 10-37 amino acids. The sampling with ST is found to be more efficient than with REMD for a much lower CPU cost. Starting from unfolded or extended conformations, the WW domain and the Trp-cage peptide fold to their NMR structures with a backbone RMSD of 2.0 and 1 Å. Remarkably, the ST simulation explores transient non-native topologies for Trp-cage that have been rarely discussed by other simulations. Our ST simulations also show that the CHARMM22* force field has limitations in describing accurately the beta3s peptide. Taken together, these results open the door to the study of the configurations of single proteins, protein aggregates, and any molecular systems at atomic details in explicit solvent using a single normal CPU. They also demonstrate that our ST scheme can be used with any force field ranging from quantum mechanics to coarse-grain and atomistic.

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.

The OPEP protein model: from single molecules, amyloid formation, crowding and hydrodynamics to DNA/RNA systems

The OPEP coarse-grained protein model has been applied to a wide range of applications since its first release 15 years ago. The model, which combines energetic and structural accuracy and chemical specificity, allows the study of single protein properties, DNA-RNA complexes, amyloid fibril formation and protein suspensions in a crowded environment. Here we first review the current state of the model and the most exciting applications using advanced conformational sampling methods. We then present the current limitations and a perspective on the ongoing developments.

Effect of the Tottori familial disease mutation (D7N) on the monomers and dimers of Aβ40 and Aβ42

Recent experiments have shown that the mutation Tottori (D7N) alters the toxicity, assembly and rate of fibril formation of the wild type (WT) amyloid beta (Aβ) Aβ40 and Aβ42 peptides. We used all-atom molecular dynamics simulations in explicit solvent of the monomer and dimer of both alloforms with their WT and D7N sequences. The monomer simulations starting from a random coil and totaling 3 μs show that the D7N mutation changes the fold and the network of salt bridges in both alloforms. The dimer simulations starting from the amyloid fibrillar states and totaling 4.4 μs also reveal noticeable changes in terms of secondary structure, salt bridge, and topology. Overall, this study provides physical insights into the enhanced rate of fibril formation upon D7N mutation and an atomic picture of the D7N-mediated conformational change on Aβ40 and Aβ42 peptides.

PRIMO: A Transferable Coarse-grained Force Field for Proteins

We describe here the PRIMO (PRotein Intermediate Model) force field, a physics-based fully transferable additive coarse-grained potential energy function that is compatible with an all-atom force field for multi-scale simulations. The energy function consists of standard molecular dynamics energy terms plus a hydrogen-bonding potential term and is mainly parameterized based on the CHARMM22/CMAP force field in a bottom-up fashion. The solvent is treated implicitly via the generalized Born model. The bonded interactions are either harmonic or distance-based spline interpolated potentials. These potentials are defined on the basis of all-atom molecular dynamics (MD) simulations of dipeptides with the CHARMM22/CMAP force field. The non-bonded parameters are tuned by matching conformational free energies of diverse set of conformations with that of CHARMM all-atom results. PRIMO is designed to provide a correct description of conformational distribution of the backbone (ϕ/ψ) and side chains (χ1) for all amino acids with a CMAP correction term. The CMAP potential in PRIMO is optimized based on the new CHARMM C36 CMAP. The resulting optimized force field has been applied in MD simulations of several proteins of 36-155 amino acids and shown that the root-mean-squared-deviation of the average structure from the corresponding crystallographic structure varies between 1.80 and 4.03 Å. PRIMO is shown to fold several small peptides to their native-like structures from extended conformations. These results suggest the applicability of the PRIMO force field in the study of protein structures in aqueous solution, structure predictions as well as ab initio folding of small peptides.

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.

Structural, thermodynamical, and dynamical properties of oligomers formed by the amyloid NNQQ peptide: insights from coarse-grained simulations

Characterizing the early formed oligomeric intermediates of amyloid peptides is of particular interest due to their links with neurodegenerative diseases. Here we study the NNQQ peptide, known to display parallel β-strands in amyloid fibrils by x-ray microcrystallography, and investigate the structural, thermodynamical, and dynamical properties of 20 NNQQ peptides using molecular dynamics and replica exchange molecular dynamics simulations coupled to a coarse-grained force field. All simulations are initiated from randomized and fully dispersed monomeric conformations. Our simulations reveal that the phase transition is characterized by a change in the oligomer and β-sheet size distributions and the percentage of mixed parallel/antiparallel β-strands when the sheets are formed. At all temperatures, however, the fraction of parallel β-strands remains low, though there are many association/fragmentation events. This work and a growing body of computational studies provide strong evidence that the critical nucleus goes beyond 20 chains and reordering of the β-strands occurs in larger oligomers.

Structures of Aβ17-42 trimers in isolation and with five small-molecule drugs using a hierarchical computational procedure

The amyloid-β protein (Aβ) oligomers are believed to be the main culprits in the cytoxicity of Alzheimer's disease (AD) and p3 peptides (Aβ17-42 fragments) are present in AD amyloid plaques. Many small-molecule or peptide-based inhibitors are known to slow down Aβ aggregation and reduce the toxicity in vitro, but their exact modes of action remain to be determined since there has been no atomic level of Aβ(p3)-drug oligomers. In this study, we have determined the structure of Aβ17-42 trimers both in aqueous solution and in the presence of five small-molecule inhibitors using a multiscale computational study. These inhibitors include 2002-H20, curcumin, EGCG, Nqtrp, and resveratrol. First, we used replica exchange molecular dynamics simulations coupled to the coarse-grained (CG) OPEP force field. These CG simulations reveal that the conformational ensemble of Aβ17-42 trimer can be described by 14 clusters with each peptide essentially adopting turn/random coil configurations, although the most populated cluster is characterized by one peptide with a β-hairpin at Phe19-Leu31. Second, these 14 dominant clusters and the less-frequent fibril-like state with parallel register of the peptides were subjected to atomistic Autodock simulations. Our analysis reveals that the drugs have multiple binding modes with different binding affinities for trimeric Aβ17-42 although they interact preferentially with the CHC region (residues 17-21). The compounds 2002-H20 and Nqtrp are found to be the worst and best binders, respectively, suggesting that the drugs may interfere at different stages of Aβ oligomerization. Finally, explicit solvent molecular dynamics of two predicted Nqtrp-Aβ17-42 conformations describe at atomic level some possible modes of action for Nqtrp.

Flexibility and binding affinity in protein-ligand, protein-protein and multi-component protein interactions: limitations of current computational approaches

The recognition process between a protein and a partner represents a significant theoretical challenge. In silico structure-based drug design carried out with nothing more than the three-dimensional structure of the protein has led to the introduction of many compounds into clinical trials and numerous drug approvals. Central to guiding the discovery process is to recognize active among non-active compounds. While large-scale computer simulations of compounds taken from a library (virtual screening) or designed de novo are highly desirable in the post-genomic area, many technical problems remain to be adequately addressed. This article presents an overview and discusses the limits of current computational methods for predicting the correct binding pose and accurate binding affinity. It also presents the performances of the most popular algorithms for exploring binary and multi-body protein interactions.

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.

Structure and folding of disulfide-rich miniproteins: insights from molecular dynamics simulations and MM-PBSA free energy calculations

The fold of small disulfide-rich proteins largely relies on two or more disulfide bridges that are main components of the hydrophobic core. Because of the small size of these proteins and their high cystine content, the cysteine connectivity has been difficult to ascertain in some cases, leading to uncertainties and debates in the literature. Here, we use molecular dynamics simulations and MM-PBSA free energy calculations to compare similar folds with different disulfide pairings in two disulfide-rich miniprotein families, namely the knottins and the short-chain scorpion toxins, for which the connectivity has been discussed. We first show that the MM-PBSA approach is able to discriminate the correct knotted topology of knottins from the laddered one. Interestingly, a comparison of the free energy components for kalata B1 and MCoTI-II suggests that cyclotides and squash inhibitors, although sharing the same scaffold, are stabilized through different interactions. Application to short-chain scorpion toxins suggests that the conventional cysteine pairing found in many homologous toxins is significantly more stable than the unconventional pairing reported for maurotoxin and for spinoxin. This would mean that native maurotoxin and spinoxin are not at the lowest free energy minimum and might result from kinetically rather than thermodynamically driven oxidative folding processes. For both knottins and toxins, the correct or conventional disulfide connectivities provide lower flexibilities and smaller deviations from the initial conformations. Overall, our work suggests that molecular dynamics simulations and the MM-PBSA approach to estimate free energies are useful tools to analyze and compare disulfide bridge connectivities in miniproteins.

Energy landscapes of the monomer and dimer of the Alzheimer's peptide Abeta(1-28)

The cytotoxicity of Alzheimer's disease has been linked to the self-assembly of the 4042 amino acid of the amyloid-beta (Abeta) peptide into oligomers. To understand the assembly process, it is important to characterize the very first steps of aggregation at an atomic level of detail. Here, we focus on the N-terminal fragment 1-28, known to form fibrils in vitro. Circular dichroism and NMR experiments indicate that the monomer of Abeta(1-28) is alpha-helical in a membranelike environment and random coil in aqueous solution. Using the activation-relaxation technique coupled with the OPEP coarse grained force field, we determine the structures of the monomer and of the dimer of Abeta(1-28). In agreement with experiments, we find that the monomer is predominantly random coil in character, but displays a non-negligible beta-strand probability in the N-terminal region. Dimerization impacts the structure of each chain and leads to an ensemble of intertwined conformations with little beta-strand content in the region Leu17-Ala21. All these structural characteristics are inconsistent with the amyloid fibril structure and indicate that the dimer has to undergo significant rearrangement en route to fibril formation.

The complex folding pathways of protein A suggest a multiple-funnelled energy landscape

Folding proteins into their native states requires the formation of both secondary and tertiary structures. Many questions remain, however, as to whether these form into a precise order, and various pictures have been proposed that place the emphasis on the first or the second level of structure in describing folding. One of the favorite test models for studying this question is the B domain of protein A, which has been characterized by numerous experiments and simulations. Using the activation-relaxation technique coupled with a generic energy model (optimized potential for efficient peptide structure prediction), we generate more than 50 folding trajectories for this 60-residue protein. While the folding pathways to the native state are fully consistent with the funnel-like description of the free energy landscape, we find a wide range of mechanisms in which secondary and tertiary structures form in various orders. Our nonbiased simulations also reveal the presence of a significant number of non-native beta and alpha conformations both on and off pathway, including the visit, for a non-negligible fraction of trajectories, of fully ordered structures resembling the native state of nonhomologous proteins.

Potential energy and free energy landscapes

Familiar concepts for small molecules may require reinterpretation for larger systems. For example, rearrangements between geometrical isomers are usually considered in terms of transitions between the corresponding local minima on the underlying potential energy surface, V. However, transitions between bulk phases such as solid and liquid, or between the denatured and native states of a protein, are normally addressed in terms of free energy minima. To reestablish a connection with the potential energy surface we must think in terms of representative samples of local minima of V, from which a free energy surface is projected by averaging over most of the coordinates. The present contribution outlines how this connection can be developed into a tool for quantitative calculations. In particular, stepping between the local minima of V provides powerful methods for locating the global potential energy minimum, and for calculating global thermodynamic properties. When the transition states that link local minima are also sampled we can exploit statistical rate theory to obtain insight into global dynamics and rare events. Visualizing the potential energy landscape helps to explain how the network of local minima and transition states determines properties such as heat capacity features, which signify transitions between free energy minima. The organization of the landscape also reveals how certain systems can reliably locate particular structures on the experimental time scale from among an exponentially large number of local minima. Such directed searches not only enable proteins to overcome Levinthal's paradox but may also underlie the formation of "magic numbers" in molecular beams, the self-assembly of macromolecular structures, and crystallization.

Computational study of the Trp-cage miniprotein based on the ECEPP/3 force field

Using a newly developed Monte Carlo global optimization method called basin paving, we have performed an ab initio computation for the structure of Trp-cage based on the ECEPP/3 force field in vacuo. The lowest energy minimum has been located. Its corresponding configuration is comparable to the native structure of Trp-cage (PDB code 1L2Y) with a backbone root mean square deviation of 2.24 A.

Structural characterization of VGVAPG, an elastin-derived peptide

Elastic fibers are an important component of the extracellular matrix, providing elasticity and resilience to tissues that require the ability to deform repetitively and reversibly. Among the elastin-derived peptides, the Val-Gly-Val-Ala-Pro-Gly (VGVAPG) hexapeptide is known for its chemotactic activity and metalloproteinases upregulation properties. As other elastin-derived peptides, having homologous similar sequences, do not exhibit any biological activity, the following question arises: Does the peptide-receptor interaction need a specific active conformation? Previous experimental studies including NMR and CD spectroscopies did not clearly identify the conformations adopted by the VGVAPG peptide in solution. However, structural predictions made on VGVAPG and related XGXXPG peptides suggested a folded beta-turn conformation. So we undertook a theoretical and experimental study of the VGVAPG peptide. The work presented here, which gives an overall structural description of VGVAPG behavior in water, also provides an additional insight into its structure-activity relationship. Both theoretical and experimental results suggest the existence of an ensemble of rather extended and folded conformations in solution. All the folded structures obtained exhibit a type VIII beta-turn spanning the GVAP sequence. In the lack of any structural information concerning the elastin receptor, these results suggest that such a conformation could be relevant for the peptide-receptor interaction and thus for biological activity.

The Antitumor Properties of the α3(IV)-(185-203) Peptide from the NC1 Domain of Type IV Collagen (Tumstatin) Are Conformation-dependent

Tumor progression may be controlled by various fragments derived from noncollagenous 1 (NC1) C-terminal domains of type IV collagen. We demonstrated previously that a peptide sequence from the NC1 domain of the alpha3(IV) collagen chain inhibits the in vitro expression of matrix metalloproteinases in human melanoma cells through RGD-independent binding to alpha(v)beta(3) integrin. In the present paper, we demonstrate that in a mouse melanoma model, the NC1 alpha3(IV)-(185-203) peptide inhibits in vivo tumor growth in a conformation-dependent manner. The decrease of tumor growth is the result of an inhibition of cell proliferation and a decrease of cell invasive properties by down-regulation of proteolytic cascades, mainly matrix metalloproteinases and the plasminogen activation system. A shorter peptide comprising the seven N-terminal residues 185-191 (CNYYSNS) shares the same inhibitory profile. The three-dimensional structures of the CNYYSNS and NC1 alpha3(IV)-(185-203) peptides show a beta-turn at the YSNS (188-191) sequence level, which is crucial for biological activity. As well, the homologous MNYYSNS heptapeptide keeps the beta-turn and the inhibitory activity. In contrast, the DNYYSNS heptapeptide, which does not form the beta-turn at the YSNS level, is devoid of inhibitory activity. Structural studies indicate a strong structure-function relationship of the peptides and point to the YSNS turn as necessary for biological activity. These peptides could act as potent and specific antitumor antagonists of alpha(v)beta(3) integrin in melanoma progression.

The chicken-egg scenario of protein folding revisited

What is the first step in protein folding - hydrophobic collapse (compaction) or secondary structure formation? It is still not clear if the major driving force in protein folding is hydrogen bonding or hydrophobic interactions or both. We analyzed data on the conformational characteristics of 41 globular proteins in native and partially folded conformational states. Our analysis shows that a good correlation exists between relative decrease in hydrodynamic volume and increase in secondary structure content. No compact equilibrium intermediates lacking secondary structure, or highly ordered non-compact species, were found. This correlation provides experimental support for the hypothesis that hydrophobic collapse occurs simultaneously with formation of secondary structure in the early stages of the protein folding.

Generating ensemble averages for small proteins from extended conformations by Monte Carlo simulations

Since it is not feasible to determine the structure of every protein by experiment, algorithms delivering the folded conformation of a protein solely from its amino acid sequence are desirable. Here the diffusion-process controlled-Monte Carlo approach has been applied to generating ensemble averages for three small proteins with 31, 36, and 46 residues. Starting from extended conformations and using an energy model that was developed on other protein models, the simulations find nativelike structures deviating by 3 A rms from the experimental structures for the main chain atoms. The balance between long-range and short-range interactions is discussed briefly in the context of stability and folding.