Effective energy function for proteins in solution

A Gaussian solvent-exclusion model for the solvation free energy is developed. It is based on theoretical considerations and parametrized with experimental data. When combined with the CHARMM 19 polar hydrogen energy function, it provides an effective energy function (EEF1) for proteins in solution. The solvation model assumes that the solvation free energy of a protein molecule is a sum of group contributions, which are determined from values for small model compounds. For charged groups, the self-energy contribution is accounted for primarily by the exclusion model. Ionic side-chains are neutralized, and a distance-dependent dielectric constant is used to approximate the charge-charge interactions in solution. The resulting EEF1 is subjected to a number of tests. Molecular dynamics simulations at room temperature of several proteins in their native conformation are performed, and stable trajectories are obtained. The deviations from the experimental structures are similar to those observed in explicit water simulations. The calculated enthalpy of unfolding of a polyalanine helix is found to be in good agreement with experimental data. Results reported elsewhere show that EEF1 clearly distinguishes correctly from incorrectly folded proteins, both in static energy evaluations and in molecular dynamics simulations and that unfolding pathways obtained by high-temperature molecular dynamics simulations agree with those obtained by explicit water simulations. Thus, this energy function appears to provide a realistic first approximation to the effective energy hypersurface of proteins.     

Reversible scaling of dihedral angle barriers during molecular dynamics to improve structure prediction of cyclic peptides.

Peptide cyclization or chemical cross-linking has frequently been used to restrict the conformational freedom of a peptide, for example, to enhance its capacity for selective binding to a target receptor molecule. Structure prediction of cyclic peptides is important to evaluate possible conformations prior to synthesis. Because of the conformational constraints imposed by cyclization low energy conformations of cyclic peptides can be separated by large energy barriers. In order to improve the conformational search properties of molecular dynamics (MD) simulations a potential scaling method has been designed. The approach consists of several consecutive MD simulations with a specific lowering of dihedral energy barriers and reduced nonbonded interactions between atoms separated by three atoms followed by gradually scaling the potential until the original barriers are reached. Application to four cyclic penta- and hexa-peptide test cases and a protein loop of known structure indicates that the potential scaling method is more efficient and faster in locating low energy conformations than standard MD simulations. Combined with a generalized Born implicit solvation model the low energy cyclic peptide conformations and the loop structure are in good agreement with experiment. Applications in the presence of explicit water molecules during the simulations showed also improved convergence to structures close to experiment compared with regular MD.     

Structure refinement of protein model decoys requires accurate side‐chain placement

In this study, the application of temperature‐based replica‐exchange (T‐ReX) simulations for structure refinement of decoys taken from the I‐TASSER dataset was examined. A set of eight nonredundant proteins was investigated using self‐guided Langevin dynamics (SGLD) with a generalized Born implicit solvent model to sample conformational space. For two of the protein test cases, a comparison of the SGLD/T‐ReX method with that of a hybrid explicit/implicit solvent molecular dynamics T‐ReX simulation model is provided. Additionally, the effect of side‐chain placement among the starting decoy structures, using alternative rotamer conformations taken from the SCWRL4 modeling program, was investigated. The simulation results showed that, despite having near‐native backbone conformations among the starting decoys, the determinant of their refinement is side‐chain packing to a level that satisfies a minimum threshold of native contacts to allow efficient excursions toward the downhill refinement regime on the energy landscape. By repacking using SCWRL4 and by applying the RWplus statistical potential for structure identification, the SGLD/T‐ReX simulations achieved refinement to an average of 38% increase in the number of native contacts relative to the original I‐TASSER decoy sets and a 25% reduction in values of C_α root‐mean‐square deviation. The hybrid model succeeded in obtaining a sharper funnel to low‐energy states for a modeled target than the implicit solvent SGLD model; yet, structure identification remained roughly the same. Without meeting a threshold of near‐native packing of side chains, the T‐ReX simulations degrade the accuracy of the decoys, and subsequently, refinement becomes tantamount to the protein folding problem. Proteins 2013. 2012 Published by Wiley Periodicals, Inc.     

The dynamics,structure, and conformational free energy of proline-containing antifreeze glycoprotein

Recent NMR studies of the solution structure of the 14-amino acid antifreeze glycoprotein AFGP-8 have concluded that the molecule lacks long-range order. The implication that an apparently unstructured molecule can still have a very precise function as a freezing inhibitor seems startling at first consideration. To gain insight into the nature of conformations and motions in AFGP-8, we have undertaken molecular dynamics simulations augmented with free energy calculations using a continuum solvation model. Starting from 10 different NMR structures, 20 ns of dynamics of AFGP were explored. The dynamics show that AFGP structure is composed of four segments, joined by very flexible pivots positioned at alanine 5, 8, and 11. The dynamics also show that the presence of prolines in this small AFGP structure facilitates the adoption of the poly-proline II structure as its overall conformation, although AFGP does adopt other conformations during the course of dynamics as well. The free energies calculated using a continuum solvation model show that the lowest free energy conformations, while being energetically equal, are drastically different in conformations. In other words, this AFGP molecule has many structurally distinct and energetically equal minima in its energy landscape. In addition, conformational, energetic, and hydrogen bond analyses suggest that the intramolecular hydrogen bonds between the N-acetyl group and the protein backbone are an important integral part of the overall stability of the AFGP molecule. The relevance of these findings to the mechanism of freezing inhibition is discussed.     

Conformational space exploration of Met- and Leu-enkephalin using the MOLS method, molecular dynamics, and Monte Carlo simulation--a comparative study

We report here a comparative study of the molecular conformational energy landscape generated using the mutually orthogonal Latin squares (MOLS) method, molecular dynamics (MD), and Monte Carlo (MC) simulation. The MOLS method, as described earlier from our laboratory, uses an experimental design technique to rapidly and exhaustively sample the low energy conformations of a molecule. MD and MC simulations have been used to perform similar tasks. In the comparison reported here, the three methods were applied to a pair of neuropeptides, namely Met- and Leu-enkephalin. A set of 1500 conformations of these enkephalins were generated using these methods with CHARMM22 force field, and the resulting samples were analyzed to determine the extent and nature of coverage of the conformational space. The results indicate that the MOLS method samples a larger number of possible conformations and identifies conformations closer to the experimental structures than the MD and MC simulations.     

Brownian dynamics theory for predicting internal and external blockages of tetraethylammonium in the KcsA potassium channel

The theory of Brownian dynamics is used to model permeation and the blocking of KcsA potassium channels by tetraethylammonium (TEA). A novel Brownian dynamics simulation algorithm is implemented that comprises two free energy profiles; one profile is seen by the potassium ions and the other by the TEA molecules whose shape is approximated by a sphere. Our simulations reveal that internally applied TEA blocks the passage of K⁺ ions by physically occluding the pore. A TEA molecule in the external reservoir encounters an attractive energy-well created by four tyrosine residues at position 82, in addition to all other attractive and repulsive forces impinging on it. Using Brownian dynamics, we investigate how deep the energy-well needs to be to reproduce the experimentally determined inhibitory constant k_i for the TEA blockade of KcsA or the mutant Shaker T449Y. The one-dimensional free energy profile obtained from molecular dynamics is first converted into a one-dimensional potential energy profile, and is then transformed into a three-dimensional free energy profile in Brownian dynamics by adding the short-range potential from the channel walls. When converted, the free energy profile calculated from molecular dynamics gives a well-depth of ∼10 kT. We systematically alter the depths of the profiles, and then use Brownian dynamics simulations to numerically determine the current versus TEA-concentration curves. We show that the sequence of binding and unbinding events of the TEA molecule to the binding pocket can be modeled by a first-order Markov process. The Brownian dynamics simulations also reveal that the probability of a TEA molecule binding to the binding pocket in KcsA potassium channels increases exponentially with TEA concentration and depends also on the applied potential and the K⁺ concentration in the simulation assembly.     

Monte Carlo simulation of the conformational behaviour of a polypeptide chain near a charged surface

The conformational behaviour of a short polypeptide chain in the neighbourhood of a charged plane is simulated using a Monte Carlo method. In this approach, the plane is taken as a model of an interface separating a hydrophobic region from a hydrophilic one. It is shown that in the neighbourhood of the plane, folded molecular conformations are prevailing whilst stretched conformations are preferred far from the interface. When the plane is not charged, the molecule adjusts itself parallel to the interface. For a given position of the molecule with respect to the plane, when the charge density of the plane is increased, the molecule tends to turn perpendicular to the plane. The surface may either attract or repulse the molecule depending on the value of the charge density (the plane is always negatively charged).     

Conformation and hydration of aspartame

Conformational free energy calculations using an empirical potential (ECEPP/2) and the hydration shell model were carried out on the neutral, acidic, zwitterionic, and basic forms of aspartame in the hydrated state. The results indicate that as the molecule becomes more charged, the number of low energy conformations becomes smaller and the molecule becomes less flexible. The calculated free energies of hydration of charged aspartames show that hydration has a significant effect on the conformation in solution. Only two feasible conformations were found for the zwitterionic form, and these are consistent with the conformations deduced from NMR and X-ray diffraction experiments. The calculated free energy difference between these two conformations was 1.25 kcal/mol. The less favored of the two solvated conformations can be expected to be stabilized by hydrophobic interaction of the phenyl groups in the crystal.     

Molecular dynamics simulations of a guaiacyl beta-O-4 lignin model compound: examination of intramolecular hydrogen bonding and conformational flexibility

The dynamical conformational behavior of a guaiacyl beta-O-4 lignin model compound has been investigated by molecular simulations. The potential energy surface of the molecule in vacuum has been examined by means of an adiabatic map, showing a large accessible conformational space with multiple energy minima separated by low barriers. Molecular dynamics simulations have been performed in vacuum and with explicit solvent molecules for 10 and 2.1 ns, respectively. Molecular dynamics trajectories recorded in vacuum have shown the molecule to be flexible and to visit a large number of conformations. Many intramolecular H-bonds have been observed, existing for more than 90% of the total simulation time. The presence of explicit solvent molecules induces a significant broadening of some regions of the accessible conformational space and also largely reduces the statistical significance of intramolecular H-bonding. Intramolecular H-bonds observed in vacuum do not persist significantly and are preferentially exchanged with intermolecular H-bonds to the surrounding solvent molecules. The theoretical results are in good agreement with experimental NMR data that do not support the existence of strong and persistent intramolecular H-bonds in solution but instead indicate that H-bonds to solvent predominate. Finally, both molecular modeling and NMR approaches predict the guaiacyl beta-O-4 structure to be flexible and indicate that intramolecular H-bonds are not strong and persistent enough to confer rigidity to the molecule in solution.     

Molecular dynamics simulations of the interaction of phospholipid bilayers with polycaprolactone

The molecular interaction between common polymer chains and the cell membrane is unknown. Molecular dynamics simulations offer an emerging tool to characterise the nature of the interaction between common degradable polymer chains used in biomedical applications, such as polycaprolactone, and model cell membranes. Herein we characterise with all-atomistic and coarse-grained molecular dynamics simulations the interaction between single polycaprolactone chains of varying chain lengths with a phospholipid membrane. We find that the length of the polymer chain greatly affects the nature of interaction with the membrane, as well as the membrane properties. Furthermore, we next utilise advanced sampling techniques in molecular dynamics to characterise the two-dimensional free energy surface for the interaction of varying polymer chain lengths (short, intermediate, and long) with model cell membranes. We find that the free energy minimum shifts from the membrane-water interface to the hydrophobic core of the phospholipid membrane as a function of chain length. Finally, we perform coarse-grained molecular dynamics simulations of slightly larger membranes with polymers of the same length and characterise the results as compared with all-atomistic molecular dynamics simulations. These results can be used to design polymer chain lengths and chemistries to optimise their interaction with cell membranes at the molecular level.     

Binding of Net-Fla,a Netropsin-Flavin Hybrid Molecule,to DNA: Molecular Mechanics and Dynamics Studies in vacuo and in Water Solution

Abstract

We have studied the binding of the hybrid netropsin-flavin (Net-Fla) molecule onto four sequences containing four A.T base pairs. Molecular mechanics minimizations in vacuo show numerous minimal conformations separated by one base pair. 400 ps molecular dynamics simulations in vacuo have been performed using the lowest minima as the starting conformations. During these simulations, the flavin moiety of the drug makes two hydrogen bonds with an amino group of a neighboring guanine. A 200 ps molecular dynamics simulation in explicit water solution suggests that the binding of Net-Fla upon the DNA substrate is enhanced by water bridges. A water molecule bridging the amidinium of Net-Fla to the N3 atom of an adenine seems to be stuck in the dmg-DNA complex during the whole simulation. The fluctuations of the DNA helical parameters and of the torsion angles of the sugar-phosphate backbone are very similar in the simulations in vacuo and in water. The time auto-correlation functions for the DNA helical parameters decrease rapidly in the picosecond range in vacuo. The same functions computed from the water solution molecular dynamics simulations seem to have two modes: the rapid mode is similar to the behavior in vacuo, and is followed by a slower mode in the 10 ps range.     

On Calculation of the Electrostatic Potential of a Phosphatidylinositol Phosphate-Containing Phosphatidylcholine Lipid Membrane Accounting for Membrane Dynamics

Many signaling events require the binding of cytoplasmic proteins to cell membranes by recognition of specific charged lipids, such as phosphoinositol-phosphates. As a model for a protein-membrane binding site, we consider one charged phosphoinositol phosphate (PtdIns(3)P) embedded in a phosphatidylcholine bilayer. As the protein-membrane binding is driven by electrostatic interactions, continuum solvent models require an accurate representation of the electrostatic potential of the phosphoinositol phosphate-containing membrane. We computed and analyzed the electrostatic potentials of snapshots taken at regular intervals from molecular dynamics simulations of the bilayer. We observe considerable variation in the electrostatic potential of the bilayer both along a single simulation and between simulations performed with the GAFF or CHARMM c36 force fields. However, we find that the choice of GAFF or CHARMM c36 parameters has little effect on the electrostatic potential of a given configuration of the bilayer with a PtdIns(3)P embedded in it. From our results, we propose a remedian averaging method for calculating the electrostatic potential of a membrane system that is suitable for simulations of protein-membrane binding with a continuum solvent model.     

Folding Trp-cage to NMR resolution native structure using a coarse-grained protein model

We develop a coarse-grained protein model with a simplified amino acid interaction potential. Using this model, we perform discrete molecular dynamics folding simulations of a small 20-residue protein--Trp-cage--from a fully extended conformation. We demonstrate the ability of the Trp-cage model to consistently reach conformations within 2-angstroms backbone root-mean-square distance from the corresponding NMR structures. The minimum root-mean-square distance of Trp-cage conformations in simulations can be <1 angstroms. Our findings suggest that, at least in the case of Trp-cage, a detailed all-atom protein model with a molecular mechanics force field is not necessary to reach the native state of a protein. Our results also suggest that the success of folding Trp-cage in our simulations and in the reported all-atom molecular mechanics simulation studies may be mainly due to the special stabilizing features specific to this miniprotein.     

A combined 2D-NMR and molecular dynamics analysis of the structure of the actinomycin D: d(ATGCAT)2 complex

We present a comparative analysis of an NMR experiment and molecular and harmonic dynamics simulations of an actinomycin D: d(ATGCAT)2 complex. A comparison of NOE measurements and 1/R6 weighted proton-proton distances confirm the general correctness of the Actinomycin D-DNA model proposed by Sobell. There are, however, some substantial differences between the proton-proton distances inferred from the NOE results and the molecular and harmonic dynamics simulations. The remaining discrepancies could either come from contributions of other conformations to the average properties of the complex or from uncertainties in the NMR distance analysis. An analysis of the molecular dynamics helix properties, sugar puckers, hydrogen bonding, rms fluctuations and torsional properties are qualitatively consistent with those from previous simulations, but the presence of an intercalated drug leads to some new structural and dynamical features.     

Molecular dynamics simulation of amyloid beta dimer formation

Recent experiments with amyloid beta (Abeta) peptide indicate that formation of toxic oligomers may be an important contribution to the onset of Alzheimer's disease. The toxicity of Abeta oligomers depends on their structure, which is governed by assembly dynamics. Due to limitations of current experimental techniques, a detailed knowledge of oligomer structure at the atomic level is missing. We introduce a molecular dynamics approach to study Abeta dimer formation. 1), We use discrete molecular dynamics simulations of a coarse-grained model to identify a variety of dimer conformations; and 2), we employ all-atom molecular mechanics simulations to estimate thermodynamic stability of all dimer conformations. Our simulations of a coarse-grained Abeta peptide model predicts 10 different planar beta-strand dimer conformations. We then estimate the free energies of all dimer conformations in all-atom molecular mechanics simulations with explicit water. We compare the free energies of Abeta(1-42) and Abeta(1-40) dimers. We find that 1), dimer conformations have higher free energies compared to their corresponding monomeric states; and 2), the free-energy difference between the Abeta(1-42) and the corresponding Abeta(1-40) dimer conformation is not significant. Our results suggest that Abeta oligomerization is not accompanied by the formation of thermodynamically stable planar beta-strand dimers.     

Conformational Particularities of the Tachykinin-Like Decapeptide Sialokinin I

This paper reports the study of the three-dimensional structure of the tachykinin-like decapeptide sialokinin I molecule using molecular mechanics and molecular dynamics. Fragment analysis was used to identify the stable spatial structures of sialokinin I, which may be present as a set of conformations that are characterized by the relatively labile N-terminal tripeptide and the conformationally rigid C-terminal heptapeptide. It has been demonstrated that the sialokinin I molecule tends to adopt nearly isoenergetic conformations with different structural types at the N-end of the peptide chain, which change into the alpha-helix turn at its C-end. Molecular dynamics was employed to model the sialokinin I molecule mobility in its stable conformations both in a vacuum and when surrounded by water molecules.

     

DFTMD studies of glucose and epimers: anomeric ratios, rotamer populations, and hydration energies

Results are presented from density functional molecular dynamics (DFTMD) simulations, based on constant energy dynamics, of glucose and its cyclic form of 6-carbon epimers. Both in vacuo and an implicit solvent method (COSMO) were examined, including simulations of low-energy conformations of each molecule. Analysis of the DFTMD results includes the following: energies averaged over the simulation time, calculated anomeric ratios, hydroxyl and hydroxymethyl rotamer populations, and hydration energies. Hydrogen-bonding networks persistence times were examined, and the effects of solvation on rotamer populations were described. Anomeric ratios calculated from energy optimization of an ensemble of low-energy conformers are compared to those obtained from ensemble averages from molecular dynamics, with dynamics simulations giving populations in best agreement with experimental anomeric ratios. Ensemble results in vacuo were not in agreement with experimental anomeric ratios or hydroxymethyl populations, producing in some cases reversal of the α:β ratios. The difficulty in obtaining correct α:β ratios increases with the number of axial groups; the mono-axial epimers being best represented, epimers with two axial groups being more difficult, and the epimers with three axial hydroxyl groups being most difficult to analyze, the result of a large number of very strong hydrogen-bonding networks that form the ensemble of low-energy conformations in the multi-axial structures.     

Simulated Q-annealing: conformational search with an effective potential

We have tested a version of the generalized simulated annealing algorithm based on molecular dynamics simulations with effective potential suggested by Tsallis statistics. The generalized annealing method, termed “simulated Q-annealing (SQ),” is applied to the simulations of a synthetic 11-residue peptide segment (1AQG). In SQ, the energy barriers between local minima change as the parameter q is varied and specific degrees of freedom can be selectively heated up and annealed. Conformational dynamics obtained by ordinary simulated annealing (SA) and SQ simulations are compared in order to illustrate the effectiveness of the SQ approach in conformational searching. We show that SQ can navigate the potential energy surface efficiently with a simple annealing protocol and demonstrate that conformations sampled by SQ can represent the funnel-like free energy surface.     

Conformations of Gly(n)H+ and Ala(n)H+ peptides in the gas phase

High-resolution ion mobility measurements and molecular dynamics simulations have been used to probe the conformations of protonated polyglycine and polyalanine (Gly(n)H and Ala(n)H+, n = 3-20) in the gas phase. The measured collision integrals for both the polyglycine and the polyalanine peptides are consistent with a self-solvated globule conformation, where the peptide chain wraps around and solvates the charge located on the terminal amine. The conformations of the small peptides are governed entirely by self-solvation, whereas the larger ones have additional backbone hydrogen bonds. Helical conformations, which are stable for neutral Alan peptides, were not observed in the experiments. Molecular dynamics simulations for Ala(n)H+ peptides suggest that the charge destabilizes the helix, although several of the low energy conformations found in the simulations for the larger Ala(n)H+ peptides have small helical regions.