Peptide folding simulations

Rethinking the protein folding problem from a new perspective

One of the main concerns of Anfinsen was to reveal the connection between the amino-acid sequence and their biologically active conformation. This search gave rise to two crucial questions in structural biology, namely, why the proteins fold and how a sequence encodes its folding. As to the why, he proposes a plausible answer, namely, the thermodynamic hypothesis. As to the how, this remains an unsolved challenge. Consequently, the protein folding problem is examined here from a new perspective, namely, as an 'analytic whole'. Conceiving the protein folding in this way enabled us to (i) examine in detail why the force-field-based approaches have failed, among other purposes, in their ability to predict the three-dimensional structure of a protein accurately; (ii) propose how to redefine them to prevent these shortcomings, and (iii) conjecture on the origin of the state-of-the-art numerical-methods success to predict the tridimensional structure of proteins accurately.

Computing long time scale biomolecular dynamics using quasi-stationary distribution kinetic Monte Carlo (QSD-KMC)

It is a challenge to obtain an accurate model of the state-to-state dynamics of a complex biological system from molecular dynamics (MD) simulations. In recent years, Markov state models have gained immense popularity for computing state-to-state dynamics from a pool of short MD simulations. However, the assumption that the underlying dynamics on the reduced space is Markovian induces a systematic bias in the model, especially in biomolecular systems with complicated energy landscapes. To address this problem, we have devised a new approach we call quasistationary distribution kinetic Monte Carlo (QSD-KMC) that gives accurate long time state-to-state evolution while retaining the entire time resolution even when the dynamics is highly non-Markovian. The proposed method is a kinetic Monte Carlo approach that takes advantage of two concepts: (i) the quasistationary distribution, the distribution that results when a trajectory remains in one state for a long time (the dephasing time), such that the next escape is Markovian, and (ii) dynamical corrections theory, which properly accounts for the correlated events that occur as a trajectory passes from state to state before it settles again. In practice, this is achieved by specifying, for each escape, the intermediate states and the final state that has resulted from the escape. Implementation of QSD-KMC imposes stricter requirements on the lengths of the trajectories than in a Markov state model approach as the trajectories must be long enough to dephase. However, the QSD-KMC model produces state-to-state trajectories that are statistically indistinguishable from an MD trajectory mapped onto the discrete set of states for an arbitrary choice of state decomposition. Furthermore, the aforementioned concepts can be used to construct a Monte Carlo approach to optimize the state boundaries regardless of the initial choice of states. We demonstrate the QSD-KMC method on two one-dimensional model systems, one of which is a driven nonequilibrium system, and on two well-characterized biomolecular systems.

Atomistic Peptide Folding Simulations Reveal Interplay of Entropy and Long-Range Interactions in Folding Cooperativity

Understanding how proteins fold has remained a problem of great interest in biophysical research. Atomistic computer simulations using physics-based force fields can provide important insights on the interplay of different interactions and energetics and their roles in governing the folding thermodynamics and mechanism. In particular, generalized Born (GB)-based implicit solvent force fields can be optimized to provide an appropriate balance between solvation and intramolecular interactions and successfully recapitulate experimental conformational equilibria for a set of helical and β-hairpin peptides. Here, we further demonstrate that key thermodynamic properties and their temperature dependence obtained from replica exchange molecular dynamics simulations of these peptides are in quantitative agreement with experimental results. Useful lessons can be learned on how the interplay of entropy and sequentially long-range interactions governs the mechanism and cooperativity of folding. These results highlight the great potential of high-quality implicit solvent force fields for studying protein folding and large-scale conformational transitions.

Revisiting Chameleon Sequences in the Protein Data Bank

The steady growth of the Protein Data Bank (PDB) suggests the periodic repetition of searches for sequences that form different secondary structures in different protein structures; these are called chameleon sequences. This paper presents a fast (nlog(n)) algorithm for such searches and presents the results on all protein structures in the PDB. The longest such sequence found consists of 20 residues.

A new class of enhanced kinetic sampling methods for building Markov state models

Markov state models (MSMs) and other related kinetic network models are frequently used to study the long-timescale dynamical behavior of biomolecular and materials systems. MSMs are often constructed bottom-up using brute-force molecular dynamics (MD) simulations when the model contains a large number of states and kinetic pathways that are not known a priori. However, the resulting network generally encompasses only parts of the configurational space, and regardless of any additional MD performed, several states and pathways will still remain missing. This implies that the duration for which the MSM can faithfully capture the true dynamics, which we term as the validity time for the MSM, is always finite and unfortunately much shorter than the MD time invested to construct the model. A general framework that relates the kinetic uncertainty in the model to the validity time, missing states and pathways, network topology, and statistical sampling is presented. Performing additional calculations for frequently-sampled states/pathways may not alter the MSM validity time. A new class of enhanced kinetic sampling techniques is introduced that aims at targeting rare states/pathways that contribute most to the uncertainty so that the validity time is boosted in an effective manner. Examples including straightforward 1D energy landscapes, lattice models, and biomolecular systems are provided to illustrate the application of the method. Developments presented here will be of interest to the kinetic Monte Carlo community as well.

Molecular dynamics simulations and applications in computational toxicology and nanotoxicology

Nanotoxicology studies toxicity of nanomaterials and has been widely applied in biomedical researches to explore toxicity of various biological systems. Investigating biological systems through in vivo and in vitro methods is expensive and time taking. Therefore, computational toxicology, a multi-discipline field that utilizes computational power and algorithms to examine toxicology of biological systems, has gained attractions to scientists. Molecular dynamics (MD) simulations of biomolecules such as proteins and DNA are popular for understanding of interactions between biological systems and chemicals in computational toxicology. In this paper, we review MD simulation methods, protocol for running MD simulations and their applications in studies of toxicity and nanotechnology. We also briefly summarize some popular software tools for execution of MD simulations.

Capturing Transition States for tRNA Hybrid-State Formation in the Ribosome

In order to quantitatively describe the energetics of biomolecular rearrangements, it is necessary to identify reaction coordinates that accurately capture the relevant transition events. Here, we perform simulations of A-site tRNA movement (∼20 Å) during hybrid-state formation in the ribosome and quantify the ability of interatomic distances to capture the transition state ensemble. Numerous coordinates are found to be accurate indicators of the transition state, allowing tRNA rearrangements to be described as diffusion across a one-dimensional free-energy surface. In addition to providing insights into the physical-chemical relationship between biomolecular structure and dynamics, these results can help enable single-molecule techniques to probe the free-energy landscape of the ribosome.

Capturing transition paths and transition states for conformational rearrangements in the ribosome

To reveal the molecular determinants of biological function, one seeks to characterize the interactions that are formed in conformational and chemical transition states. In other words, what interactions govern the molecule's energy landscape? To accomplish this, it is necessary to determine which degrees of freedom can unambiguously identify each transition state. Here, we perform simulations of large-scale aminoacyl-transfer RNA (aa-tRNA) rearrangements during accommodation on the ribosome and project the dynamics along experimentally accessible atomic distances. From this analysis, we obtain evidence for which coordinates capture the correct number of barrier-crossing events and accurately indicate when the aa-tRNA is on a transition path. Although a commonly used coordinate in single-molecule experiments performs poorly, this study implicates alternative coordinates along which rearrangements are accurately described as diffusive movements across a one-dimensional free-energy profile. From this, we provide the theoretical foundation required for single-molecule techniques to uncover the energy landscape governing aa-tRNA selection by the ribosome.

Mapping the intrinsically disordered properties of the flexible loop domain of Bcl-2: a molecular dynamics simulation study

Most of the B-cell lymphoma-2 (Bcl-2) protein structure has been elucidated; however, the conformation of its flexible loop domain (FLD) has not yet been experimentally predicted. Its high flexibility under physiological conditions is the reason. FLD behaves as an intrinsically disordered region (IDR) and can adopt regular structures in particular conditions associated with the control of Bcl-2's anti-apoptotic functions. In a previous contribution, we analyzed an engineered Bcl-2 construct (Bcl-2-Δ22Σ3) submitted to 25-ns MD and reported a disordered-to-helix transitions in a region of FLD (rFLD, residues 60-77). However, the conformational preferences in solution of rFLD in the nanosecond to microsecond scale were not analyzed. Herein, an average model was obtained for the native Bcl-2 protein by homology modeling and MD simulation techniques. From this, only the atomic coordinates corresponding to the rFLD were simulated for 1 μs by MD at 310 K. In concordance with previous studies, a disordered-to-helix transitions were exhibited, implying that this "interconversion of folding" in the rFLD suggest a possible set of conformations encoded in its sequence. Principal component analysis (PCA) showed that most of the conformational fluctuation of Bcl-2 is provided by rFLD. Dihedral PCA (dPCA) offered information about all the conformations of rFLD in the μs of the simulation, characterizing a dPCA-based free energy landscape of rFLD, and a conformational ensemble of fast interconverting conformations as other IDRs. Furthermore, despite the conformational heterogeneity of rFLD, the analysis of the dihedral angles (Φ, Ψ) showed that this region does not randomly explore the conformational space in solution.

The "Hot-Solvent/Cold-Solute" Problem Revisited

The temperature steers the equilibrium and nonequilibrium conformational dynamics of macromolecules in solution. Therefore, corresponding molecular dynamics simulations require a strategy for temperature control which should guarantee that the experimental statistical ensemble is also sampled in silico. Several algorithms for temperature control have been proposed. All these thermostats interfere with the macromolecule's "natural" dynamics as given by the Newtonian mechanics. Furthermore, using a single thermostat for an inhomogeneous solute-solvent system can lead to stationary temperature gradients. To avoid this "hot solvent/cold solute" problem, two separate thermostats are frequently applied, one to the solute and one to the solvent. However, such a separate temperature control will perturb the dynamics of the macromolecule much more strongly than a global one and, therefore, can introduce large artifacts into its conformational dynamics. Based on the concept that an explicit solvent environment represents an ideal thermostat concerning the magnitude and time correlation of temperature fluctuations of the solute, we propose a temperature control strategy that, on the one hand, provides a homogeneous temperature distribution throughout the system together with the correct statistical ensemble for the solute molecule while, on the other hand, minimally perturbing its dynamics.

Resolution Exchange Simulation with Incremental Coarsening

We previously developed an algorithm, called "resolution exchange", which improves canonical sampling of atomic resolution models by swapping conformations between high- and low-resolution simulations. Here, we demonstrate a generally applicable incremental coarsening procedure and apply the algorithm to a larger peptide, met-enkephalin. In addition, we demonstrate a combination of resolution and temperature exchange, in which the coarser simulations are also at elevated temperatures. Both simulations are implemented in a "top-down" mode, to allow efficient allocation of CPU time among the different replicas.

Toward a Coarse-Grained Protein Model Coupled with a Coarse-Grained Solvent Model: Solvation Free Energies of Amino Acid Side Chains

Recently, we reported that molecular dynamics (MD) simulations using a coarse-grained (CG) peptide model coupled with a CG water model are able to reproduce many of the structural and thermodynamic features of short peptides with nonpolar side chains at 10(3) times the normal speed (JCTC, 2007, 3, 2146-2161). To further develop a CG protein model for MD simulations, we systematically parametrized the side chains of all 20 naturally occurring amino acids. We developed the parameters by fitting the dihedral potentials of 13 small molecules, the densities and self-solvation free energies of liquids of eight organic molecules, and the hydration free energies of 35 small organic molecules. In a set of 11 classes of compounds (105 in total) including alkanes, alcohols, ethers, ketones/aldehydes, amines, amides, aromatics, carboxylic acids, sulfides/thiols, alkyl ammoniums, and carboxylate ions, the average error in the calculated hydration free energies compared with experimental results is about 1.4 kJ/mol. The average error in the calculated transfer free energies of the 19 side-chain analogues of amino acids from cyclohexane to water is about 2.2 kJ/mol. These results are comparable to the results of all-atom models.

How Efficient Is Replica Exchange Molecular Dynamics? An Analytic Approach

Replica exchange molecular dynamics (REMD) has become a standard technique for accelerating relaxation in biosimulations. Despite its widespread use, questions remain about its efficiency compared with conventional, constant temperature molecular dynamics (MD). An analytic approach is taken to describe the relative efficiency of REMD with respect to MD. This is applied to several simple two-state models and to several real proteins-protein L and the B domain of protein A-to predict the relative efficiency of REMD with respect to MD in actual applications. In agreement with others, we find the following: as long as there is a positive activation energy for folding, REMD is more efficient than MD; the effectiveness of REMD is strongly dependent on the activation enthalpy; and the efficiency of REMD for actual proteins is a strong function of the maximum temperature. Choosing the maximum temperature too high can result in REMD becoming significantly less efficient than conventional MD. A good rule of thumb appears to be to choose the maximum temperature of the REMD simulation slightly above the temperature at which the enthalpy for folding vanishes. Additionally, we find that the number of replicas in REMD, while important for simulations shorter than one or two relaxation times, has a minimal effect on the asymptotic efficiency of the method.

Molecular modeling of peptides

This article presents a review of the field of molecular modeling of peptides. The main focus is on atomistic modeling with molecular mechanics potentials. The description of peptide conformations and solvation through potentials is discussed. Several important computer simulation methods are briefly introduced, including molecular dynamics, accelerated sampling approaches such as replica-exchange and metadynamics, free energy simulations and kinetic network models like Milestoning. Examples of recent applications for predictions of structure, kinetics, and interactions of peptides with complex environments are described. The reliability of current simulation methods is analyzed by comparison of computational predictions obtained using different models with each other and with experimental data. A brief discussion of coarse-grained modeling and future directions is also presented.

Molecular dynamics simulations of glycoproteins using CHARMM

Molecular dynamics simulations are an effective tool to study the structure, dynamics, and thermodynamics of carbohydrates and proteins. However, the simulations of heterogeneous glycoprotein systems have been limited due to the lack of appropriate molecular force field parameters describing the linkage between the carbohydrate and the protein regions as well as the tools to prepare these systems for modeling studies. In this work we outline the recent developments in the CHARMM carbohydrate force field to treat glycoproteins and describe in detail the step-by-step procedures involved in building glycoprotein geometries using CHARMM-GUI Glycan Reader.

Enhanced sampling techniques in molecular dynamics simulations of biological systems

BACKGROUND

Molecular dynamics has emerged as an important research methodology covering systems to the level of millions of atoms. However, insufficient sampling often limits its application. The limitation is due to rough energy landscapes, with many local minima separated by high-energy barriers, which govern the biomolecular motion.

SCOPE OF REVIEW

In the past few decades methods have been developed that address the sampling problem, such as replica-exchange molecular dynamics, metadynamics and simulated annealing. Here we present an overview over theses sampling methods in an attempt to shed light on which should be selected depending on the type of system property studied.

MAJOR CONCLUSIONS

Enhanced sampling methods have been employed for a broad range of biological systems and the choice of a suitable method is connected to biological and physical characteristics of the system, in particular system size. While metadynamics and replica-exchange molecular dynamics are the most adopted sampling methods to study biomolecular dynamics, simulated annealing is well suited to characterize very flexible systems. The use of annealing methods for a long time was restricted to simulation of small proteins; however, a variant of the method, generalized simulated annealing, can be employed at a relatively low computational cost to large macromolecular complexes.

GENERAL SIGNIFICANCE

Molecular dynamics trajectories frequently do not reach all relevant conformational substates, for example those connected with biological function, a problem that can be addressed by employing enhanced sampling algorithms. This article is part of a Special Issue entitled Recent developments of molecular dynamics.

Detailed microscopic unfolding pathways of an α-helix and a β-hairpin: direct observation and molecular dynamics

We present a combined experimental and computational study of unfolding pathways of a model 21-residue α-helical heteropeptide (W1H5-21) and a 16-residue β-hairpin (GB41-56). Experimentally, we measured fluorescence energy transfer efficiency as a function of temperature, employing natural tryptophans as donors and dansylated lysines as acceptors. Secondary structural analysis was performed with circular dichroism and Fourier transform infrared spectroscopy. Our studies present markedly different unfolding pathways of the two elementary secondary structural elements. During thermal denaturation, the helical peptide exhibits an initial decrease in length, followed by an increase, while the hairpin undergoes a systematic increase in length. In the complementary computational part of the project, we performed microsecond length replica-exchange molecular dynamics simulations of the peptides in explicit solvent, yielding a detailed microscopic picture of the unfolding processes. For the α-helical peptide, we found a large heterogeneous population of intermediates that are primarily frayed single helices or helix-turn-helix motifs. Unfolding starts at the termini and proceeds through a stable helical region in the interior of the peptide but shifted off-center toward the C-terminus. The simulations explain the experimentally observed non-monotonic variation of helix length with temperature as due primarily to the presence of frayed-end single-helix intermediate structures. For the β-hairpin peptide, our simulations indicate that folding is initiated at the turn, followed by formation of the hairpin in zipper-like fashion, with Cα···Cα contacts propagating from the turn to termini and hairpin hydrogen bonds forming in parallel with these contacts. In the early stages of hairpin formation, the hydrophobic side-chain contacts are only partly populated. Intermediate structures with low numbers of β-hairpin hydrogen bonds have very low populations. This is in accord with the "broken zipper" model of Scheraga. The monotonic increase in length with temperature may be explained by the zipper-like breaking of the hairpin hydrogen bonds and backbone contacts.

Non-exponential kinetics and a complete folding pathway of an α-helical heteropeptide: direct observation and comprehensive molecular dynamics

We have performed a combined experimental and computational study of the folding of a 21-residue α-helix-forming heteropeptide (WH21). Temperature jump kinetics with improved dynamic range at several temperatures revealed non-exponential relaxation that could be well described with two time constants of 20 and 300 ns at 298 K. In the computational part, we performed multi-microsecond molecular dynamics simulations of WH21 in explicit water, using the AMBER03 and OPLS/AA potentials. The simulations were in good agreement with available experimental data on helix content and relaxation times. On the basis of 70 individual transitions, we identified the main pathways of helix unfolding. Three paths were found in both force fields, with unfolding progressing through (1) N-terminus, C-terminus, and center; (2) C-terminus, N-terminus, and center; and (3) C-terminus, center, and N-terminus. An additional fourth path starting in the central region and expanding to the termini was detected only with AMBER03. Intermediate structures sampled along the pathway included a central helix with frayed termini, an off-center helix, and a helical hairpin. The simulations suggest that the short relaxation should be assigned to partly cooperative fluctuations of several neighboring hydrogen bonds. Overall, by a combination of ultrafast kinetic measurements and detailed microscopic description through comprehensive molecular dynamics, we have obtained important new insights into the helix folding process.

Rational optimization of conformational effects induced by hydrocarbon staples in peptides and their binding interfaces

eIF4E is frequently over-expressed in different cancers and causes increased translation of oncogenic proteins via deregulated cap-dependent translation. Inhibitors of the eIF4E:eIF4G interactions represent an approach that would normalize cap-dependent translation. Stapled peptides represent an emerging class of therapeutics that can target protein: protein interactions. We present here molecular dynamics simulations for a set of rationally designed stapled peptides in solution and in complex with eIF4E, supported with biophysical and crystallographic data. Clustering of the simulated structures revealed the favoured conformational states of the stapled peptides in their bound or free forms in solution. Identifying these populations has allowed us to design peptides with improved affinities by introducing mutations into the peptide sequence to alter their conformational distributions. These studies emphasise the effects that engineered mutations have on the conformations of free and bound peptides, and illustrate that both states must be considered in efforts to attain high affinity binding.

A membrane-proximal, C-terminal α-helix is required for plasma membrane localization and function of the G Protein-coupled receptor (GPCR) TGR5

The C terminus of G protein-coupled receptors (GPCRs) is important for G protein-coupling and activation; in addition, sorting motifs have been identified in the C termini of several GPCRs that facilitate correct trafficking from the endoplasmic reticulum to the plasma membrane. The C terminus of the GPCR TGR5 lacks any known sorting motif such that other factors must determine its trafficking. Here, we investigate deletion and substitution variants of the membrane-proximal C terminus of TGR5 with respect to plasma membrane localization and function using immunofluorescence staining, flow cytometry, and luciferase assays. Peptides of the membrane-proximal C-terminal variants are subjected to molecular dynamics simulations and analyzed with respect to their secondary structure. Our results reveal that TGR5 plasma membrane localization and responsiveness to extracellular ligands is fostered by a long (≥ 9 residues) α-helical stretch at the C terminus, whereas the presence of β-strands or only a short α-helical stretch leads to retention in the endoplasmic reticulum and a loss of function. As a proof-of-principle, chimeras of TGR5 containing the membrane-proximal amino acids of the β2 adrenergic receptor (β2AR), the sphingosine 1-phosphate receptor-1 (S1P1), or the κ-type opioid receptor (κOR) were generated. These TGR5β2AR, TGR5S1P1, or TGR5κOR chimeras were correctly sorted to the plasma membrane. As the exchanged amino acids of the β2AR, the S1P1, or the κOR form α-helices in crystal structures but lack significant sequence identity to the respective TGR5 sequence, we conclude that the secondary structure of the TGR5 membrane-proximal C terminus is the determining factor for plasma membrane localization and responsiveness towards extracellular ligands.

Deducing conformational variability of intrinsically disordered proteins from infrared spectroscopy with Bayesian statistics

As it remains practically impossible to generate ergodic ensembles for large intrinsically disordered proteins (IDP) with molecular dynamics (MD) simulations, it becomes critical to compare spectroscopic characteristics of the theoretically generated ensembles to corresponding measurements. We develop a Bayesian framework to infer the ensemble properties of an IDP using a combination of conformations generated by MD simulations and its measured infrared spectrum. We performed 100 different MD simulations totaling more than 10 µs to characterize the conformational ensemble of αsynuclein, a prototypical IDP, in water. These conformations are clustered based on solvent accessibility and helical content. We compute the amide-I band for these clusters and predict the thermodynamic weights of each cluster given the measured amide-I band. Bayesian analysis produces a reproducible and non-redundant set of thermodynamic weights for each cluster, which can then be used to calculate the ensemble properties. In a rigorous validation, these weights reproduce measured chemical shifts.

On the foldability of tryptophan-containing tetra- and pentapeptides: an exhaustive molecular dynamics study

Short peptides serve as minimal model systems to decipher the determinants of foldability due to their simplicity arising from their smaller size, their ability to echo protein-like structural characteristics, and their direct implication in force field validation. Here, we describe an effort to identify small peptides that can still form stable structures in aqueous solutions. We followed the in silico folding of a selected set of 8640 tryptophan-containing tetra- and pentapeptides through 15 210 molecular dynamics simulations amounting to a total of 272.46 μs using explicit representation of the solute and full treatment of the electrostatics. The evaluation and sorting of peptides is achieved through scoring functions, which include terms based on interatomic vector distances, atomic fluctuations, and rmsd matrices between successive frames of a trajectory. Highly scored peptides are studied further via successive simulation rounds of increasing simulation length and using different empirical force fields. Our method suggested only a handful of peptides with strong foldability prognosis. The discrepancies between the predictions of the various force fields for such short sequences are also extensively discussed. We conclude that the vast majority of such short peptides do not adopt significantly stable structures in water solutions, at least based on our computational predictions. The present work can be utilized in the rational design and engineering of bioactive peptides with desired molecular properties.

Identification of minimally interacting modules in an intrinsically disordered protein

The conformational characterization of intrinsically disordered proteins (IDPs) is complicated by their conformational heterogeneity and flexibility. If an IDP could somehow be divided into smaller fragments and reconstructed later, theoretical and spectroscopic studies could probe its conformational variability in detail. Here, we used replica molecular-dynamics simulations and network theory to explore whether such a divide-and-conquer strategy is feasible for α-synuclein, a prototypical IDP. We characterized the conformational variability of α-synuclein by conducting >100 unbiased all-atom molecular-dynamics simulations, for a total of >10 μs of trajectories. In these simulations, α-synuclein formed a heterogeneous ensemble of collapsed coil states in an aqueous environment. These states were stabilized by heterogeneous contacts between sequentially distant regions. We find that α-synuclein contains residual secondary structures in the collapsed states, and the heterogeneity in the collapsed state makes it feasible to split α-synuclein into sequentially contiguous minimally interacting fragments. This study reveals previously unknown characteristics of α-synuclein and provides a new (to our knowledge) approach for studying other IDPs.

Biomolecular dynamics: order-disorder transitions and energy landscapes

While the energy landscape theory of protein folding is now a widely accepted view for understanding how relatively weak molecular interactions lead to rapid and cooperative protein folding, such a framework must be extended to describe the large-scale functional motions observed in molecular machines. In this review, we discuss (1) the development of the energy landscape theory of biomolecular folding, (2) recent advances toward establishing a consistent understanding of folding and function and (3) emerging themes in the functional motions of enzymes, biomolecular motors and other biomolecular machines. Recent theoretical, computational and experimental lines of investigation have provided a very dynamic picture of biomolecular motion. In contrast to earlier ideas, where molecular machines were thought to function similarly to macroscopic machines, with rigid components that move along a few degrees of freedom in a deterministic fashion, biomolecular complexes are only marginally stable. Since the stabilizing contribution of each atomic interaction is on the order of the thermal fluctuations in solution, the rigid body description of molecular function must be revisited. An emerging theme is that functional motions encompass order-disorder transitions and structural flexibility provides significant contributions to the free energy. In this review, we describe the biological importance of order-disorder transitions and discuss the statistical-mechanical foundation of theoretical approaches that can characterize such transitions.

Structure of a peptide adsorbed on graphene and graphite

Noncovalent functionalization of graphene using peptides is a promising method for producing novel sensors with high sensitivity and selectivity. Here we perform atomic force microscopy, Raman spectroscopy, infrared spectroscopy, and molecular dynamics simulations to investigate peptide-binding behavior to graphene and graphite. We studied a dodecamer peptide identified with phage display to possess affinity for graphite. Optical spectroscopy reveals that the peptide forms secondary structures both in powder form and in an aqueous medium. The dominant structure in the powder form is α-helix, which undergoes a transition to a distorted helical structure in aqueous solution. The peptide forms a complex reticular structure upon adsorption on graphene and graphite, having a helical conformation different from α-helix due to its interaction with the surface. Our observation is consistent with our molecular dynamics calculations, and our study paves the way for rational functionalization of graphene using biomolecules with defined structures and, therefore, functionalities.

MOLECULAR DYNAMIC SIMULATIONS ON THE FOLDING AND CONFORMATIONAL INSIGHTS OF THE TRUNCATED PEPTIDES

A total of 120 ns molecular dynamics simulations was used to study the folding and conformational aspects of six peptides with different lengths (Pep19–25, Pep15–25, Pep1–25, Pep15–39, Pep1–40, and Pep1–50) truncated from the αβ-tubulin dimer. These truncated peptides were found to undergo distinct structural transitions, with Pep1–25 and Pep1–50 folding into their respective stable conformations whereas on the contrary for the others. All the six truncated peptides are more or less compact than the corresponding segments in the αβ-tubulin dimer. The most striking contraction was observed in Pep1–25, which folds in a similar manner of β-hairpin. Pep1–50 has the least contraction and its folded conformation is the closest to that in the αβ-tubulin dimer. Moreover, the same conversions of β12–β23 from helices to hydrogen-bonded turns were witnessed in both Pep1–50 and the αβ-tubulin dimer. The structural instabilities of Pep19–25, Pep15–25, Pep15–39, and Pep1–40 were caused by the lack of long-distance interactions or/and the absence of key residues, with the details given in the discussions. The folding and conformational divergences of six truncated peptides were also observed in their active peptide segments ( Ap 15–25). Ap 15–25 in Pep1–50 achieves the best agreements with the αβ-tubulin dimer, implying that the local structure of Ap 15–25 in the αβ-tubulin dimer can be well reserved in Pep1–50 rather than in the other truncated peptides. The long-distance interactions, especially the key residues (e.g. β48-Arg), play a crucial role in the correct folding of Ap 15–25. The correct folding into the stable conformations is a prerequisite for the peptides to implement their catalytic actions, and therefore the present results are helpful to the future designs of active peptides.

Influence of solvent and intramolecular hydrogen bonding on the conformational properties of o-linked glycopeptides

A detailed investigation of the conformational properties of all the biologically relevant O-glycosidic linkages using the Hamiltonian replica exchange (HREX) simulation methodology and the recently developed CHARMM carbohydrate force field parameters is presented. Fourteen biologically relevant O-linkages between the five sugars N-acetylgalactosamine (GalNAc), N-acetylglucosamine (GlcNAc), D-glucose (Glc), D-mannose (Man), and L-fucose (Fuc) and the amino acids serine and threonine were studied. The force field was tested by comparing the simulation results of the model glycopeptides to various NMR (3)J coupling constants, NOE distances, and data from molecular dynamics with time-averaged restraints (tar-MD). The results show the force field to be in overall agreement with experimental and previous tar-MD simulations, although some small limitations are identified. An in-depth hydrogen bond and bridging water analysis revealed an interplay of hydrogen bonding and bridge water interactions influencing the geometry of the underlying peptide backbone, with the O-linkages favoring extended β-sheet and polyproline type II (PPII) conformations over the compact α(R)-helical conformation. The newly developed parameters were also able to identify hydrogen bonding and water mediated interactions between O-linked sugars and proteins. These results indicate that the newly developed parameters in tandem with HREX conformational sampling provide the means to study glycoproteins in the absence of targeted NMR restraint data.

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.

Using the Wimley-White Hydrophobicity Scale as a Direct Quantitative Test of Force Fields: The MARTINI Coarse-Grained Model

The partitioning of proteins and peptides at the membrane/water interface is a key step in many processes, including the action of antimicrobial peptides, cell-penetrating peptides, and toxins, as well as signaling. To develop a computational model that can be used to accurately represent such systems, the underlying model must be able to quantitatively represent the partitioning preferences of amino acids in the lipid membrane. The MARTINI model provides a consistent set of parameters for building coarse-grained models of systems involving lipids and proteins. Even though MARTINI is parametrized to reproduce the partitioning behavior of small molecules, its ability to reproduce partitioning preferences of amino acids at lipid/water interfaces has never been tested. In this study, we measured the partitioning free energies of side chains of amino acids using alchemical simulations and umbrella sampling. The pentapeptides of sequence Ac-WLXLL were simulated at the POPC/water and cyclohexane/water interfaces using MARTINI, and the computed free energies were compared with the Wimley-White hydrophobicity scale. The free energy values obtained using the free energy perturbation, thermodynamic integration, and umbrella sampling methods were compared to gain insight into the most efficient method and the degree of sampling required to obtain statistically accurate free energies for use with atomistic force fields in future work. With the standard MARTINI water model, the amino acids D, E, K, and R were found to be significantly too favorable in hydrophobic environments, whereas with the polarizable water model, the amino acids D, E, K, and R were found to give correct free energies of partitioning. The amino acids P and F showed significant deviations from the experimental values. This model system will be used in future improvements to the MARTINI model.

Protein folding in a reverse micelle environment: the role of confinement and dehydration

Characterization of the molecular interactions that stabilize the folded state of proteins including hydrogen bond formation, solvation, molecular crowding, and interaction with membrane environments is a fundamental goal of theoretical biophysics. Inspired by recent experimental studies by Gai and co-workers, we have used molecular dynamics simulations to explore the structure and dynamics of the alanine-rich AKA(2) peptide in bulk solution and in a reverse micelle environment. The simulated structure of the reverse micelle shows substantial deviations from a spherical geometry. The AKA(2) peptide is observed to (1) remain in a helical conformation within a spherically constrained reverse micelle and (2) partially unfold when simulated in an unconstrained reverse micelle environment, in agreement with experiment. While aqueous solvation is found to stabilize the N- and C-termini random coil portions of the peptide, the helical core region is stabilized by significant interaction between the nonpolar surface of the helix and the aliphatic chains of the AOT surfactant. The results suggest an important role for nonpolar peptide-surfactant and peptide-lipid interactions in stabilizing helical geometries of peptides in reverse micelle environments.

Order through disorder: hyper-mobile C-terminal residues stabilize the folded state of a helical peptide. a molecular dynamics study

Conventional wisdom has it that the presence of disordered regions in the three-dimensional structures of polypeptides not only does not contribute significantly to the thermodynamic stability of their folded state, but, on the contrary, that the presence of disorder leads to a decrease of the corresponding proteins' stability. We have performed extensive 3.4 µs long folding simulations (in explicit solvent and with full electrostatics) of an undecamer peptide of experimentally known helical structure, both with and without its disordered (four residue long) C-terminal tail. Our simulations clearly indicate that the presence of the apparently disordered (in structural terms) C-terminal tail, increases the thermodynamic stability of the peptide's folded (helical) state. These results show that at least for the case of relatively short peptides, the interplay between thermodynamic stability and the apparent structural stability can be rather subtle, with even disordered regions contributing significantly to the stability of the folded state. Our results have clear implications for the understanding of peptide energetics and the design of foldable peptides.