Molecular dynamics simulations of beta-hairpin folding

Molecular dynamics simulations of beta-hairpin folding have been carried out with a solvent-referenced potential at 274 K. The model peptide V4DPGV4 formed stable beta-hairpin conformations and the beta-hairpin ratio calculated by the DSSP algorithm was about 56% in the 50-ns simulation. Folding into beta-hairpin conformations is independent of the initial conformations. The simulations provided insights into the folding mechanism. The hydrogen bond often formed in a beta-turn first, and then propagated by forming more hydrogen bonds along the strands. Unfolding and refolding occurred repeatedly during the simulations. Both the hydrogen bonding and the hydrophobic interaction played important roles in forming the ordered structure. Without the hydrophobic effect, stable beta-hairpin conformations did not form in the simulations. With the same energy functions, the alanine-based peptide (AAQAA)3Y folded into helical conformations, in agreement with experiments. Folding into an alpha-helix or a beta-hairpin is amino acid sequence-dependent.     

NMR study and molecular dynamics simulations of optimized beta-hairpin fragments of protein G

The stability and structure of several beta-hairpin peptide variants derived from the C-terminus of the B1 domain of protein G were investigated by a number of experimental and computational techniques. Our analysis shows that the structure and stability of this hairpin can be greatly affected by one or a few simple mutations. For example, removing an unfavorable charge near the N-terminus of the peptide (Glu42 to Gln or Thr) or optimization of the N-terminal charge-charge interactions (Gly41 to Lys) both stabilize the peptide, even in water. Furthermore, a simple replacement of a charged residue in the turn (Asp47 to Ala) changes the beta-turn conformation. Finally, we show that the effects of combining these single mutations are additive, suggesting that independent stabilizing interactions can be isolated and evaluated in a simple model system. Our results indicate that the structure and stability of this beta-hairpin peptide can be modulated in numerous ways and thus contributes toward a more complete understanding of this important model beta-hairpin as well as to the folding and stability of larger peptides and proteins.     

Molecular dynamics simulations of folding processes of a beta-hairpin in an implicit solvent

Computer simulations of beta-hairpin folding are relatively difficult, especially those based on the explicit water model. This greatly limits the complete analysis and understanding of their folding mechanisms. In this paper, we use the generalized Born/solvent accessible implicit solvent model to simulate the folding processes of a nine-residue beta-hairpin. We find that the beta-hairpin can fold into its native structure very easily, even using the traditional molecular dynamics method. This allows us to extract 21 complete folding events and investigate the folding process sufficiently. Our results show that there exist four most stable states on the free energy landscape of the short peptide, one native state and three intermediates. We find that two of the non-native stable states have almost the same potential energy as the native state but with lower entropy. This suggests that the native state can be stabilized entropically. Furthermore, we find that the folding processes of this peptide have common features: to fold into its native state, the peptide undergoes a continuous collapsing-extending-recollapsing process to adjust the positions of the side chains in order to form the native middle inter-strand hydrogen bonds. The formations of these bonds are the key step of the folding process. Once these bonds are formed, the peptide can fold into the native state quickly.     

Understanding beta-hairpin formation by molecular dynamics simulations of unfolding

We have studied the mechanism of formation of a 16-residue beta-hairpin from the protein GB1 using molecular dynamics simulations in an aqueous environment. The analysis of unfolding trajectories at high temperatures suggests a refolding pathway consisting of several transient intermediates. The changes in the interaction energies of residues are related with the structural changes during the unfolding of the hairpin. The electrostatic energies of the residues in the turn region are found to be responsible for the transition between the folded state and the hydrophobic core state. The van der Waals interaction energies of the residues in the hydrophobic core reflect the behavior of the radius of gyration of the core region. We have examined the opposing influences of the protein-protein (PP) energy, which favors the native state, and the protein-solvent (PS) energy, which favors unfolding, in the formation of the beta-hairpin structure. It is found that the behavior of the electrostatic components of PP and PS energies reflects the structural changes associated with the loss of backbone hydrogen bonding. Relative changes in the PP and PS van der Waals interactions are related with the disruption of the hydrophobic core of a protein. The results of the simulations support the hydrophobic collapse mechanism of beta-hairpin folding.     

Molecular dynamics simulations and membrane protein structure quality

Despite a growing repertoire of membrane protein structures (currently ∼120 unique structures), considerations of low resolution and crystallization in the absence of a lipid bilayer require the development of techniques to assess the global quality of membrane protein folds. This is also the case for assessment of, e.g. homology models of human membrane proteins based on structures of (distant) bacterial homologues. Molecular dynamics (MD) simulations may be used to help evaluate the quality of a membrane protein structure or model. We have used a structure of the bacterial ABC transporter MsbA which has the correct transmembrane helices but an incorrect handedness and topology of their packing to test simulation methods of quality assessment. An MD simulation of the MsbA model in a lipid bilayer is compared to a simulation of another bacterial ABC transporter, BtuCD. The latter structure has demonstrated good conformational stability in the same bilayer environment and over the same timescale (20 ns) as for the MsbA model simulation. A number of comparative analyses of the two simulations were performed to assess changes in the structural integrity of each protein. The results show a significant difference between the two simulations, chiefly due to the dramatic structural deformations of MsbA. We therefore propose that MD could become a useful quality control tool for membrane protein structural biology. In particular, it provides a way in which to explore the global conformational stability of a model membrane protein fold.     

Molecular dynamics simulations of Trp side-chain conformational flexibility in the gramicidin A channel

Gramicidin A (gA) is prototypical peptide antibiotic and a model ion channel former. Configured in the solid-state NMR beta(6.5)-helix channel conformation, gA was subjected to 1-ns molecular dynamics (MD) gas phase simulations using the all-atom charmm22 force field to ascertain the conformational stability of the Trp side chains as governed by backbone and neighboring side-chain contacts. Three microcanonical trajectories were computed using different initial atomic velocities for each of twenty different initial structures. For each set, one of the four Trp side chains in each monomer was initially positioned in one of the five non-native conformations (A. E. Dorigo et al., Biophysical Journal, 1999, Vol. 76, 1897-1908), the other Trps being positioned in the native state, o1. In three additional control simulations, all Trps were initiated in the native conformation. After equilibration, constraints were removed and subsequent conformational changes of the initially constrained Trp were measured. The chi(1) was more flexible than chi(2.1). The energetically optimal orientation, o1 (Dorigo et al., 1999), was the most stable in all four Trp positions (9, 11, 13, 15) and remained unchanged for the entire 1 ns simulation in 19 of 24 trials. Changes in chi(1) from each of the 5 suboptimal states occur readily. Two of the non-native conformations reverted readily to o1, whereas the other three converted to an intermediate state, i2. There were frequent interconversions between i2 and o1. We speculate that experimentally observed Trp stability is caused by interactions with the lipid-water interface, and that stabilization of one of the suboptimal conformations in gA, such as i2, by lipid headgroups could produce a secondary, metastable conformational state. This could explain recent experimental studies of differences in the channel conductance dispersity between gA and a Trp-to-Phe gA analog, gramicidin M (gM, J. C. Markham et al., Biochimica et Biophysica Acta, 2001, Vol. 1513, 185-192).     

Role of hydrophilic and hydrophobic contacts in folding of the second beta-hairpin fragment of protein G: molecular dynamics simulation studies of an all-atom model

Predicting the folding mechanism of the second beta-hairpin fragment of the Ig-binding domain B of streptococcal protein G is unexpectedly challenging for simplified reduced models because the models developed so far indicated a different folding mechanism from what was suggested from high-temperature unfolding and equilibrium free-energy surface analysis based on established all-atom empirical force fields in explicit or implicit solvent. This happened despite the use of empirical residue-based interactions, multibody hydrophobic interactions, and inclusions of hydrogen bonding effects in the simplified models. This article employs a recently developed all-atom (except nonpolar hydrogens) model interacting with simple square-well potentials to fold the peptide fragment by molecular dynamics simulation methods. In this study, 193 out of 200 trajectories are folded at two reduced temperatures (3.5 and 3.7) close to the transition temperature T* approximately 4.0. Each simulation takes <7 h of CPU time on a Pentium 800-MHz PC. Folding of the new all-atom model is found to be initiated by collapse before the formation of main-chain hydrogen bonds. This verifies the mechanism proposed from previous all-atom unfolding and equilibrium simulations. The new model further predicts that the collapse is initiated by two nucleation contacts (a hydrophilic contact between D46 and T49 and a hydrophobic contact between Y45 and F52), in agreement with recent NMR measurements. The results suggest that atomic packing and native contact interactions play a dominant role in folding mechanism.     

Molecular dynamics simulations as a tool for improving protein stability

Haloalkane dehalogenase (DhlA) was used as a model protein to explore the possibility to use molecular dynamics (MD) simulations as a tool to identify flexible regions in proteins that can serve as a target for stability enhancement by introduction of a disulfide bond. DhlA consists of two domains: an alpha/beta-hydrolase fold main domain and a cap domain composed of five alpha-helices. MD simulations of DhlA showed high mobility in a helix-loop-helix region in the cap domain, involving residues 184-211. A disulfide cross-link was engineered between residue 201 of this flexible region and residue 16 of the main domain. The mutant enzyme showed substantial changes in both thermal and urea denaturation. The oxidized form of the mutant enzyme showed an increase of the apparent transition temperature from 47.5 to 52.5 degrees C, whereas the T(m,app) of the reduced mutant decreased by more than 8 degrees C compared to the wild-type enzyme. Urea denaturation results showed a similar trend. Measurement of the kinetic stability showed that the introduction of the disulfide bond caused a decrease in activation free energy of unfolding of 0.43 kcal mol(-1) compared to the wild-type enzyme and also indicated that the helix-loop-helix region was involved early in the unfolding process. The results show that MD simulations are capable of identifying mobile protein domains that can successfully be used as a target for stability enhancement by the introduction of a disulfide cross-link.     

Molecular dynamics simulations of the TrkH membrane protein

TrkH is a transmembrane protein that mediates uptake of K(+) through the cell membrane. Despite the recent determination of its crystallographic structure, the nature of the permeation mechanism is still unknown, that is, whether K(+) ions move across TrkH by active transport or passive diffusion. Here, molecular dynamics simulations and the umbrella sampling technique have been employed to shed light on this question. The existence of binding site S3 and two alternative binding sites have been characterized. Analysis of the coordination number renders values that are almost constant, with a full contribution from the carbonyls of the protein only at S3. This observation contrasts with observations of K(+) channels, where the contribution of the protein to the coordination number is roughly constant in all four binding sites. An intramembrane loop is found immediately after the selectivity filter at the intracellular side of the protein, which obstructs the permeation pathway, and this is reflected in the magnitude of the energy barriers.     

Molecular dynamics simulations on parallel and antiparallel C.G*G triplexes

Molecular dynamics (MD) studies have been carried out on the Hoogsteen hydrogen bonded parallel and the reverse Hoogsteen hydrogen bonded antiparallel C.G*G triplexes. Earlier, the molecular mechanics studies had shown that the parallel structure was energetically more favourable than the antiparallel structure. To characterize the structural stability of the two triplexes and to investigate whether the antiparallel structure can transit to an energetically more favourable structure, due to the local fluctuations in the structure during the MD simulation, the two structures were subjected to 200ps of constant temperature vacuum MD simulations at 300K. Initially no constraints were applied to the structures and it was observed that for the antiparallel triplex, the structure showed a large root mean square deviation from the starting structure within the first 12ps and the N4-H41--O6 hydrogen bond in the WC duplex got distorted due to a high propeller twist and a moderate increase in the opening angle in the basepairs. Starting from an initial value of 30 degrees , helical twist of the average structure from this simulation had a value of 36 degrees , while the parallel structure stabilized at a twist of 33 degrees. In spite of the hydrogen bond distortions in the antiparallel triplex, it was energetically comparable to the parallel triplex. To examine the structural characteristics of an undistorted structure, another MD simulation was performed on the antiparallel triplex by constraining all the hydrogen bonds. This structure stabilized at an average twist of 33 degrees. In the course of the dynamics though the energy of the molecule - compared to the initial structure - improved, it did not become comparable to the parallel structure. Energy minimization studies performed in the presence of explicit water and counterions also showed the two structures to be equally favourable energetically. Together these results indicate that the parallel C.G*G triplex with Hoogsteen hydrogen bonds also represents a stereochemically and energetically favourable structure for this class of triplexes.     

Molecular dynamics simulations of a protein on hydrophobic and hydrophilic surfaces.

Molecular dynamics simulations have been used to investigate the behavior of the peripheral membrane protein, cytochrome c, covalently tethered to hydrophobic (methyl-terminated) and hydrophilic (thiol-terminated) self-assembled monolayers (SAMs). The simulations predict that the protein will undergo minor structural changes when it is tethered to either surface, and the structures differ qualitatively on the two surfaces: the protein is less spherical on the hydrophilic SAM where the polar surface residues reach out to interact with the SAM surface. The protein is completely excluded from the hydrophobic SAM but partially dissolves in the hydrophilic SAM. Consequently, the surface of the thiol-terminated SAM is considerably less ordered than that of the methyl-terminated SAM, although a comparable, high degree of order is maintained in the bulk of both SAMs: the chains exhibit collective tilts in the nearest-neighbor direction at angles of 20 degrees and 17 degrees with respect to the surface normal in the hydrophobic and the hydrophilic SAMs, respectively. On the hydrophobic SAM the protein is oriented so that the heme plane is more nearly parallel to the surface, whereas on the hydrophilic surface it is more nearly perpendicular. The secondary structure of the protein, dominated by alpha helices, is not significantly affected, but the structure of the loops as well as the helix packing is slightly modified by the surfaces.     

Molecular dynamics simulations of a protein model in uniform and elongational flows

We present a molecular dynamics study of the conformational deformation of a minimalist beta-barrel protein model in two different types of hydrodynamic flows: uniform and elongational. We investigate the characteristics of protein stretching, paying special attention to the unfolding intermediate states and their relationship to the protein folding/unfolding problem. In the uniform flow simulations, one end of the modeled protein was tethered to a fixed point in space and the forced unfolding process was observed. The unfolding takes place via a few stages involving one or two intermediate states, depending on which end is tethered. The calculated force-extension curves show plateau regimes and hysteresis as the protein is stretched and refolded, in qualitative agreement with the experimental measurements. The physical behavior observed in our numerical simulations of the forced unfolding in an elongational flow is very different from that in uniform flow. The protein unfolds abruptly from the globular state to a fully stretched state without going through any observable intermediate states. From these observation, we stress that protein unfolding pathways under the influence of an external force are highly dependent on the mechanism of the exerted force.     

Molecular dynamics simulations of N-terminal peptides from a nucleotide binding protein

Molecular dynamics (MD) simulations of N-terminal peptides from lactate dehydrogenase (LDH) with increasing length and individual secondary structure elements were used to study their stability in relation to folding. Ten simulations of 1–2 ns of different peptides in water starting from the coordinates of the crystal structure were performed. The stability of the peptides was compared qualitatively by analyzing the root mean square deviation (RMSD) from the crystal structure, radius of gyration, secondary and tertiary structure, and solvent accessible surface area. In agreement with earlier MD studies, relatively short (< 15 amino acids) peptides containing individual secondary structure elements were generally found to be unstable; the hydrophobic α₁-helix of the nucleotide binding fold displayed a significantly higher stability, however. Our simulations further showed that the first βαβ supersecondary unit of the characteristic dinucleotide binding fold (Rossmann fold) of LDH is somewhat more stable than other units of similar length and that the α₂-helix, which unfolds by itself, is stabilized by binding to this unit. This finding suggests that the first βαβ unit could function as an N-terminal folding nucleus, upon which the remainder of the polypeptide chain can be assembled. Indeed, simulations with longer units (βαβα and βαβαββ) showed that all structural elements of these units are rather stable. The outcome of our studies is in line with suggestions that folding of the N-terminal portion of LDH in vivo can be a cotranslational process that takes place during the ribosomal peptide synthesis.     

Molecular dynamics simulations of a helicase

Helicases are ubiquitous enzymes involved in nucleic acid metabolism. The PcrA DNA helicase is an essential bacterial protein involved in rolling circle plasmid replication and DNA repair. Recent crystal structures of PcrA bound to DNA indicate that a flexible loop mediates a functionally important rigid-body-domain rotation. In this study, we report stochastic boundary molecular dynamics simulations focused on this region for wild-type and mutants designed to increase the rigidity of the region. Residues in the region that were helix-disfavoring, such as glycine, threonine, and others, were mutated to alanine. The simulated dynamics, analyzed with a variety of measures of structure and mobility, indicate that a few point mutations will substantially increase helix formation in this region. Subnanosecond stochastic boundary molecular dynamics simulations at several temperatures offer a rapid protocol for assessing large numbers of mutants and provides a novel strategy for the design of experiments to test the role of this flexible loop region in the function of PcrA.     

Prediction and evaluation of side-chain conformations for protein backbone structures

A common approach to protein modeling is to propose a backbone structure based on homology or threading and then to attempt to build side chains onto this backbone. A fast algorithm using the simple criteria of atomic overlap and overall rotamer probability is proposed for this purpose. The method was first tested in the context of exhaustive searches of side chain configuration space in protein cores and was then applied to all side chains in 49 proteins of known structure, using simulated annealing to sample space. The latter procedure obtains the correct rotamer for 57% and the correct χ₁ value for 74% of the 6751 residues in the sample. When low-temperature Monte-Carlo simulations are initiated from the results of the simulated-annealing processes, consensus configurations are obtained which exhibit slightly more accurate predictions. The Monte-Carlo procedure also allows converged side chain entropies to be calculated for all residues. These prove to be accurate indicators of prediction reliability. For example, the correct rotamer is obtained for 79% and the correct χ₁ value is obtained for 84% of the half of the sample residues exhibiting the lowest entropies. Side chain entropy and predictability are nearly completely uncorrelated with solvent-accessible area. Some precedents for and implications of this observation are discussed. © 1996 Wiley-Liss, Inc.     

Comparison between the backbone dynamics of an 11-amino acid peptide sequence in alpha-helical and beta-hairpin structural contexts

To investigate the relationship between backbone motions and the structural environment of a peptide sequence, we have used (15)N NMR relaxation data to characterize the backbone motions of the "chameleon-alpha" (Chm-alpha) and "chameleon-beta" (Chm-beta) proteins designed previously by Minor and Kim [Minor, D. L., Jr., and Kim, P. S. (1996) Nature 380, 730-734]. These two proteins contain an identical 11-amino acid sequence (dubbed the "chameleon" peptide sequence) in alpha-helix and beta-hairpin conformations, respectively, within the B1 domain of protein G. When placed in an alpha-helical context, the chameleon peptide shows very limited backbone motions, but some remote regions of the protein are induced to undergo conformational exchange motions, apparently due to modification of packing interactions with the chameleon peptide. In contrast, within a beta-hairpin context, the chameleon peptide displays substantial motions on both picosecond and microsecond-to-millisecond time scales, suggesting that it cannot be readily accommodated within the native reverse turn structure. These observations are consistent with the relatively low stability of the Chm-beta protein and can be rationalized in terms of native turn-stabilizing interactions that may be disrupted in the Chm-beta protein.     

Molecular dynamics simulations of the photoactive protein nitrile hydratase

Nitrile hydratase (NHase) is an enzyme used in the industrial biotechnological production of acrylamide. The active site, which contains nonheme iron or noncorrin cobalt, is buried in the protein core at the interface of two domains, α and β. Hydrogen bonds between βArg-56 and αCys-114 sulfenic acid (αCEA114) are important to maintain the enzymatic activity. The enzyme may be inactivated by endogenous nitric oxide (NO) and activated by absorption of photons of wavelength λ < 630 nm. To explain the photosensitivity and to propose structural determinants of catalytic activity, differences in the dynamics of light-active and dark-inactive forms of NHase were investigated using molecular dynamics (MD) modeling. To this end, a new set of force field parameters for nonstandard NHase active sites have been developed. The dynamics of the photodissociated NO ligand in the enzyme channel was analyzed using the locally enhanced sampling method, as implemented in the MOIL MD package. A series of 1 ns trajectories of NHases shows that the protonation state of the active site affects the dynamics of the catalytic water and NO ligand close to the metal center. MD simulations support the catalytic mechanism in which a water molecule bound to the metal ion directly attacks the nitrile carbon.     

Molecular dynamics simulations of the whey protein beta-lactoglobulin.

Molecular dynamics simulations have been used to model the motions and conformational behavior of the whey protein bovine beta-lactoglobulin. Simulations were performed for the protein by itself and complexed to a single retinol ligand located in a putative interior binding pocket. In the absence of the retinol ligand, the backbone loops around the opening of this interior pocket shifted inward to partially close off this cavity, similar to the shifts observed in previously reported molecular dynamics simulations of the uncomplexed form of the homologous retinol binding protein. The protein complexed with retinol does not exhibit the same conformational shifts. Conformational changes of this type could serve as a recognition signal allowing in vivo discrimination between the free and retinol complexed forms of the beta-lactoglobulin molecule. The unusual bending of the single alpha-helix observed in the simulations of retinol binding protein were not observed in the present calculations.