Predicting extreme pKa shifts in staphylococcal nuclease mutants with constant pH molecular dynamics

Accurate computational methods of determining protein and nucleic acid pK(a) values are vital to understanding pH-dependent processes in biological systems. In this article, we use the recently developed method constant pH molecular dynamics (CPHMD) to explore the calculation of highly perturbed pK(a) values in variants of staphylococcal nuclease (SNase). Simulations were performed using the replica exchange (REX) protocol for improved conformational sampling with eight temperature windows, and yielded converged proton populations in a total sampling time of 4 ns. Our REX-CPHMD simulations resulted in calculated pK(a) values with an average unsigned error (AUE) of 0.75 pK units for the acidic residues in Δ + PHS, a hyperstable variant of SNase. For highly pK(a)-perturbed SNase mutants with known crystal structures, our calculations yielded an AUE of 1.5 pK units and for those mutants based on modeled structures an AUE of 1.4 pK units was found. Although a systematic underestimate of pK shifts was observed in most of the cases for the highly perturbed pK mutants, correlations between conformational rearrangement and plasticity associated with the mutation and error in pK(a) prediction was not evident in the data. This study further extends the scope of electrostatic environments explored using the REX-CPHMD methodology and suggests that it is a reliable tool for rapidly characterizing ionizable amino acids within proteins even when modeled structures are employed.     

Enhanced sampling of peptide and protein conformations using replica exchange simulations with a peptide backbone biasing-potential

During replica exchange molecular dynamics (RexMD) simulations, several replicas of a system are simulated at different temperatures in parallel allowing for exchange between replicas at frequent intervals. This technique allows significantly improved sampling of conformational space and is increasingly being used for structure prediction of peptides and proteins. A drawback of the standard temperature RexMD is the rapid increase of the replica number with increasing system size to cover a desired temperature range. In an effort to limit the number of replicas, a new Hamiltonian-RexMD method has been developed that is specifically designed to enhance the sampling of peptide and protein conformations by applying various levels of a backbone biasing potential for each replica run. The biasing potential lowers the barrier for backbone dihedral transitions and promotes enhanced peptide backbone transitions along the replica coordinate. The application on several peptide cases including in all cases explicit solvent indicates significantly improved conformational sampling when compared with standard MD simulations. This was achieved with a very modest number of 5-7 replicas for each simulation system making it ideally suited for peptide and protein folding simulations as well as refinement of protein model structures in the presence of explicit solvent.     

Adaptive conformational sampling based on replicas

Computer simulations of biomolecules such as molecular dynamics simulations are limited by the time scale of conformational rearrangements. Several sampling techniques are available to search the multi-minima free energy landscape but most efficient, time-dependent methods do generally not produce a canonical ensemble. A sampling algorithm based on a self-regulating ladder of searching copies in the dihedral subspace is developped in this paper. The learning process using short- and long-term memory functions allows an efficient search in phase space while combining a deterministic dynamics and stochastic swaps with the searching copies conserves a canonical limit. The sampling efficiency and accuracy are indicated by comparing the ansatz with conventional molecular dynamics and replica exchange simulations.     

Structure and energy landscape of a photoswitchable peptide: a replica exchange molecular dynamics study

A replica exchange molecular dynamics (REMD) simulation of a bicyclic azobenzene peptide in explicit dimethyl sulfoxide solution is presented in order to characterize the conformational structures and energy landscape of a photoswitchable peptide. It is shown that an enhanced-sampling technique such as the REMD method is essential to obtain a converged conformational sampling of the peptide at room temperature. This is because conventional MD simulations of less than approximately 100-ns length are either trapped in local minima (at 295 K) or-if run at high temperature-do not resemble the room-temperature REMD results. Calculating various nuclear Overhauser effects (NOEs) and (3)J-couplings, a good overall agreement between the REMD simulations and the NMR experiments of Renner et al. (Biopolymers 2000;54:501-514) is found. In particular, the REMD study confirms the general picture drawn by Renner et al. that the trans-isomer of the azobenzene peptide exhibits a well-defined structure, while the cis-isomer is a conformational heterogeneous system; that is, the trans-isomer occurs in 2 well-defined conformers, while the cis-isomer represents an energetically frustrated system that leads to an ensemble of conformational structures. Employing a principal component analysis of the REMD data, the free energy landscape of the systems is studied at various temperatures. The implications for the folding and unfolding pathways of the system are discussed.     

Application of biasing‐potential replica‐exchange simulations for loop modeling and refinement of proteins in explicit solvent

Comparative or homology modeling of a target protein based on sequence similarity to a protein with known structure is widely used to provide structural models of proteins. Depending on the target‐template similarity these model structures may contain regions of limited structural accuracy. In principle, molecular dynamics (MD) simulations can be used to refine protein model structures and also to model loop regions that connect structurally conserved regions but it is limited by the currently accessible simulation time scales. A recently developed biasing potential replica exchange (BP‐REMD) method was used to refine loops and complete decoy protein structures at atomic resolution including explicit solvent. In standard REMD simulations several replicas of a system are run in parallel at different temperatures allowing exchanges at preset time intervals. In a BP‐REMD simulation replicas are controlled by various levels of a biasing potential to reduce the energy barriers associated with peptide backbone dihedral transitions. The method requires much fewer replicas for efficient sampling compared with T‐REMD. Application of the approach to several protein loops indicated improved conformational sampling of backbone dihedral angle of loop residues compared to conventional MD simulations. BP‐REMD refinement simulations on several test cases starting from decoy structures deviating significantly from the native structure resulted in final structures in much closer agreement with experiment compared to conventional MD simulations. Proteins 2010. © 2010 Wiley‐Liss, Inc.     

Ab initio folding of proteins with all-atom discrete molecular dynamics

Discrete molecular dynamics (DMD) is a rapid sampling method used in protein folding and aggregation studies. Until now, DMD was used to perform simulations of simplified protein models in conjunction with structure-based force fields. Here, we develop an all-atom protein model and a transferable force field featuring packing, solvation, and environment-dependent hydrogen bond interactions. We performed folding simulations of six small proteins (20-60 residues) with distinct native structures by the replica exchange method. In all cases, native or near-native states were reached in simulations. For three small proteins, multiple folding transitions are observed, and the computationally characterized thermodynamics are in qualitative agreement with experiments. The predictive power of all-atom DMD highlights the importance of environment-dependent hydrogen bond interactions in modeling protein folding. The developed approach can be used for accurate and rapid sampling of conformational spaces of proteins and protein-protein complexes and applied to protein engineering and design of protein-protein interactions.     

Conformational regions of Boc-Ala-Aib-Ala-OMe. Sampling with molecular dynamics simulations using time averaging of distance restraints.

The method of time averaging of distance restraints in molecular dynamics simulations is applied to Boc-Ala-Aib-Ala-OMe in order to demonstrate the improved sampling properties of this method compared to conventional distance restraining. Two conformational regions, beta-turn type II and gamma-turn, are seen during MD runs at a simulation temperature of 500 K, while in simulations with conventional distance restraining, no conformational transitions could be observed for temperatures up to 1000 K.     

Conformational sampling using high-temperature molecular dynamics

High-temperature molecular dynamics as a method for conformational search was explored on the antigen combining site of McPC 603, a phosphorylcholine binding immunoglobulin. Simulations at temperatures of 500, 800, and 1500 K were run for 111.5, 101.7, and 76.3 ps, respectively. The effectiveness of the search was assessed using a variety of methods. For the shorter hypervariable loops, molecular dynamics explored an appreciable fraction of the conformational space as evidenced by a comparison to a simple theoretical model of the size of the conformational space. However, for the longer loops and the antigen combining site as a whole, the simulation times were too short for a complete search. The simulations at 500 and 800 K both generated conformations that minimized to energies 200 kcal/mole lower than the crystal structure. However, the 1500 K simulation produced higher energy structures, even after minimization; in addition, this highest temperature run had many cis-trans peptide isomerizations. This suggests that 1500 K is too high a temperature for unconstrained conformational sampling. Comparison of the results of high temperature molecular dynamics with a direct conformational search method, [R. E. Bruccoleri & M. Karplus (1987) Biopolymers 26, 137-168]. showed that the two methods did not overlap much in conformational space. Simple geometric measures of the conformational space indicated that the direct method covered more space than molecular dynamics at the lower temperature, but not at 1500 K. The results suggest that high-temperature molecular dynamics can aid in conformational searches.     

Molecular dynamics simulations of peptides and proteins with amplified collective motions

We present a novel method that uses the collective modes obtained with a coarse-grained model/anisotropic network model to guide the atomic-level simulations. Based on this model, local collective modes can be calculated according to a single configuration in the conformational space of the protein. In the molecular dynamics simulations, the motions along the slowest few modes are coupled to a higher temperature by the weak coupling method to amplify the collective motions. This amplified-collective-motion (ACM) method is applied to two test systems. One is an S-peptide analog. We realized the refolding of the denatured peptide in eight simulations out of 10 using the method. The other system is bacteriophage T4 lysozyme. Much more extensive domain motions between the N-terminal and C-terminal domain of T4 lysozyme are observed in the ACM simulation compared to a conventional simulation. The ACM method allows for extensive sampling in conformational space while still restricting the sampled configurations within low energy areas. The method can be applied in both explicit and implicit solvent simulations, and may be further applied to important biological problems, such as long timescale functional motions, protein folding/unfolding, and structure prediction.     

Enhanced sampling of molecular dynamics simulation of peptides and proteins by double coupling to thermal bath

Low sampling efficiency in conformational space is the well-known problem for conventional molecular dynamics. It greatly increases the difficulty for molecules to find the transition path to native state, and costs amount of CPU time. To accelerate the sampling, in this paper, we re-couple the critical degrees of freedom in the molecule to environment temperature, like dihedrals in generalized coordinates or nonhydrogen atoms in Cartesian coordinate. After applying to ALA dipeptide model, we find that this modified molecular dynamics greatly enhances the sampling behavior in the conformational space and provides more information about the state-to-state transition, while conventional molecular dynamics fails to do so. Moreover, from the results of 16 independent 100?ns simulations by the new method, it shows that trpzip2 has one-half chances to reach the naive state in all the trajectories, which is greatly higher than conventional molecular dynamics. Such an improvement would provide a potential way for searching the conformational space or predicting the most stable states of peptides and proteins.     

163 Enhanced sampling of peptides and proteins with a new biasing replica exchange method

Classical MD simulations (cMD) are limited by the sampling of relevant states of the peptides. Replica exchange (REMD) methods aim to search the conformational space of proteins more efficiently (reviewed in Ostermeir & Zacharias, 2013). We have developed a Hamiltonian REMD method that takes advantage of an intrinsic property of proteins, the specific Φ ? dihedral angle combinations along the polymer backbone. By employing a coupled two-dimensional biasing potential the energy barriers along the polymer backbone are reduced more effectively than by a previous approach based on a one-D biasing potential (Kannan & Zacharias, 2007). Thus, adjacent amino acids along the polymers backbone can easily switch between favourable regions in the Ramachandran plot. Additionally, energy barriers of rotameric states of amino acid side chains of proteins are also biased in the replica runs. The method improves the sampling of conformational substates of proteins at a modest number of replicas (nine replicas in the standard set-up with one replica running without biasing potential) compared to much larger numbers necessary in the case of standard temperature (T)-REMD simulations. A further improvement is achieved by a dynamical adjustment of the penalty potential levels in the replicas such that high exchange rates and improved mixing of conformations between different replicas are guaranteed. The biasing potential (BP)-REMD method turns out to be suitable to speed up both the folding of spaghetti-like test peptides and the refinement of loop decoy structures. Starting from extended structures, an α-helical oligo-alanine and β-hairpin chignolin and the Trp-cage protein fold more rapidly in near-native structures than in cMD simulations. The BP-REMD simulations not only accelerate the folding process of test proteins but also enlarge the variety of sampled configurations in conformational space. Since flexible parts of the protein can be penalized selectively, this method provides a precise tool to investigate regions of interest of the protein.     

Simulating force-induced conformational transitions in polysaccharides with the SMD replica exchange method

Conventional steered molecular dynamics (SMD) simulations do not readily reproduce equilibrium conditions of atomic force microscopy (AFM) stretch and release measurements of polysaccharides undergoing force-induced conformational transitions because of the gap between the timescales of computer simulations ( approximately 1 mus) and AFM measurements ( approximately 1 s). To circumvent this limitation, we propose using the replica exchange method (REM) to enhance sampling during SMD simulations. By applying REM SMD to a small polysaccharide system and comparing the results with those from AFM stretching experiments, we demonstrate that REM SMD reproduces the experimental results not only qualitatively but quantitatively, approaching near equilibrium conditions of AFM measurements. As tested in this work, hysteresis and computational time of REM SMD simulations of short polysaccharide chains are significantly reduced as compared to regular SMD simulations, making REM SMD an attractive tool for studying forced-induced conformational transitions of small biopolymer systems.     

A directed essential dynamics simulation of peptide folding

We present a directed essential dynamics (DED) method for peptide and protein folding. DED is a molecular dynamics method based on the essential dynamics sampling and the principal component analysis. The main idea of DED is to use principal component analysis to determine the direction of the most active collective motion of peptides at short intervals of time (20 fs) during the folding process and then add an additional force along it to adjust the folding direction. This method can make the peptides avoid being trapped in the local minima for a long time and enhance the sampling efficiency in conformational space during the simulation. An S-peptide with 15 amino acids is used to demonstrate the DED method. The results show that DED can lead the S-peptide to fold quickly into the native state, whereas traditional molecular dynamics needs more time to do this.     

The consistency of large concerted motions in proteins in molecular dynamics simulations.

A detailed investigation is presented into the effect of limited sampling time and small changes in the force field on molecular dynamics simulations of a protein. Thirteen independent simulations of the B1 IgG-binding domain of streptococcal protein G were performed, with small changes in the simulation parameters in each simulation. Parameters studied included temperature, bond constraints, cut-off radius for electrostatic interactions, and initial placement of hydrogen atoms. The essential dynamics technique was used to reveal dynamic differences between the simulations. Similar essential dynamics properties were found for all simulations, indicating that the large concerted motions found in the simulations are not particularly sensitive to small changes in the force field. A thorough investigation into the stability of the essential dynamics properties as derived from a molecular dynamics simulation of a few hundred picoseconds is provided. Although the definition of the essential modes of motion has not fully converged in these short simulations, the subspace in which these modes are confined is found to be reproducible.     

Structure and dynamics of an RNA tetraloop: a joint molecular dynamics and NMR study

Molecular dynamics simulations of the RNA tetraloop 5'-CGCUUUUGCG-3' with high melting temperature and significant conformational heterogeneity in explicit water solvent are presented and compared to NMR studies. The NMR data allow for a detailed test of the theoretical model, including the quality of the force field and the conformational sampling. Due to the conformational heterogeneity of the tetraloop, high temperature (350 K) and locally enhanced sampling simulations need to be invoked. The Amber98 force field leads to a good overall agreement with experimental data. Based on NMR data and a principal component analysis of the 350 K trajectory, the dynamic structure of the tetraloop is revealed. The principal component free energy surface exhibits four minima, which correspond to well-defined conformational structures that differ mainly by their base stacking in the loop region. No correlation between the motion of the sugar rings and the stacking dynamics of the loop bases is found.     

Salt effects on hydrophobic‐core formation in folding of a helical miniprotein studied by molecular dynamics simulations

We have investigated effects of salt ions on folding events of a helical miniprotein chicken villin headpiece subdomain HP36. Low concentrations of ions alter electrostatic interactions between charged groups of a protein and can change the populations of conformers. Here, we compare two data sets of folding simulations of HP36 in explicit water solvent with or without ions. For efficient sampling of the conformational space of HP36, the multicanonical replica‐exchange molecular dynamics method was employed. Our analyses suggest that salt alters salt‐bridging nature of the protein at later stages of folding at room temperature. Especially, more nonnative, nonlocal salt bridges are formed at near‐native conformations in pure water. Our analyses also show that such salt‐bridge formation hinders the fully native hydrophobic‐core packing at the final stages of folding. Proteins 2014; 82:933–943. © 2013 Wiley Periodicals, Inc.     

Protein flexibility: multiple molecular dynamics simulations of insulin chain B

Multiple molecular dynamics simulations totaling more than 100 ns were performed on chain B of insulin in explicit solvent at 300 K and 400 K. Despite some individual variations, a comparison of the protein dynamics of each simulation showed similar trends and most structures were consistent with NMR experimental values, even at the elevated temperature. The importance of packing interactions in determining the conformational transitions of the protein was observed, sometimes resulting in conformations induced by localized hydrophobic interactions. The high temperature simulation generated a more diverse range of structures with similar elements of secondary structure and populated conformations to the simulations at room temperature. A broad sampling of the conformational space of insulin chain B illustrated a wide range of conformational states with many transitions at room temperature in addition to the conformational states observed experimentally. The T-state conformation associated with insulin activity was consistently present and a possible mechanism of behavior was suggested.     

Headgroup Conformations of Phospholipids from Molecular Dynamics Simulation: Sampling Challenges and Comparison to Experiment

The preferred conformations of the glycerol region of a phospholipid have been explored using replica exchange molecular dynamics (MD) simulations and compared with the results of standard MD approaches and with experiment. We found that due to isomerization rates in key torsions that are slow on the timescale of atomistic MD simulations, standard MD is not able to produce accurate equilibrium conformer distributions from reasonable trajectory lengths (e.g., on the 100 ns) timescale. Replica exchange MD, however, results in quite efficient sampling due to the rapid increase in isomerization rate with temperature. The equilibrium distributions obtained from replica exchange MD have been compared with the results of experimental nuclear magnetic resonance observations. This comparison suggests that the sampling approach demonstrated here is a valuable tool that can be used in evaluating force fields for molecular simulation of lipids.