Molecular dynamics simulations of protein folding from the transition state

Repurposing of SARS-CoV-2 compounds against Marburg Virus using MD simulation, mm/GBSA, PCA analysis, and free energy landscape

The significant mortality rate associated with Marburg virus infection made it the greatest hazard among infectious diseases. Drug repurposing using in silico methods has been crucial in identifying potential compounds that could prevent viral replication by targeting the virus's primary proteins. This study aimed at repurposing the drugs of SARS-CoV-2 for identifying potential candidates against the matrix protein VP40 of the Marburg virus. Virtual screening was performed where the control compound, Nilotinib, showed a binding score of -9.99 kcal/mol. Based on binding scores, hit compounds 9549298, 11960895, 44545852, 51039094, and 89670174 were selected that had a lower binding score than the control. Subsequent molecular dynamics (MD) simulation revealed that compound 9549298 consistently formed a hydrogen bond with the residue Gln290. This was observed both in molecular docking and MD simulation poses, indicating a strong and significant interaction with the protein. 11960895 had the most stable and consistent RMSD pattern exhibited in 100 ns simulation, while 9549298 had the most identical RMSD plot compared to the control molecule. MM/PBSA analysis showed that the binding free energy (ΔG) of 9549298 and 11960895 was lower than the control, with -30.84 and -38.86 kcal/mol, respectively. It was observed by the PCA (principal component analysis) and FEL (free energy landscape) analysis that compounds 9549298 and 11960895 had lesser conformational variation. Overall, this study proposed 9549298 and 11960895 as potential binders of VP40 MARV that can cause its inhibition, however it inherently lacks experimental validation. Furthermore, the study proposes in-vitro experiments as the next step to validate these computational findings, offering a practical approach to further explore these compounds' potential as antiviral agents.Communicated by Ramaswamy H. Sarma.

Calculation of Protein Folding Thermodynamics Using Molecular Dynamics Simulations

Despite advances in artificial intelligence methods, protein folding remains in many ways an enigma to be solved. Accurate computation of protein folding energetics could help drive fields such as protein and drug design and genetic interpretation. However, the challenge of calculating the state functions governing protein folding from first-principles remains unaddressed. We present here a simple approach that allows us to accurately calculate the energetics of protein folding. It is based on computing the energy of the folded and unfolded states at different temperatures using molecular dynamics simulations. From this, two essential quantities (ΔH and ΔCp) are obtained and used to calculate the conformational stability of the protein (ΔG). With this approach, we have successfully calculated the energetics of two- and three-state proteins, representatives of the major structural classes, as well as small stability differences (ΔΔG) due to changes in solution conditions or variations in an amino acid residue.

Molecular Dynamics and Docking Simulations of Homologous RsmE Methyltransferases Hints at a General Mechanism for Substrate Release upon Uridine Methylation on 16S rRNA

In this study, molecular dynamics (MD) and docking simulations were carried out on the crystal structure of Neisseria Gonorrhoeae RsmE aiming at free energy of binding estimation (ΔGbinding) of the methyl transfer substrate S-adenosylmethionine (SAM), as well as its homocysteine precursor S-adenosylhomocysteine (SAH). The mechanistic insight gained was generalized in view of existing homology to two other crystal structures of RsmE from Escherichia coli and Aquifex aeolicus. As a proof of concept, the crystal poses of SAM and SAH were reproduced reflecting a more general pattern of molecular interaction for bacterial RsmEs. Our results suggest that a distinct set of conserved residues on loop segments between β12, α6, and Met169 are interacting with SAM and SAH across these bacterial methyltransferases. Comparing molecular movements over time (MD trajectories) between Neisseria gonorrhoeae RsmE alone or in the presence of SAH revealed a hitherto unknown gatekeeper mechanism by two isoleucine residues, Ile171 and Ile219. The proposed gating allows switching from an open to a closed state, mimicking a double latch lock. Additionally, two key residues, Arg221 and Thr222, were identified to assist the exit mechanism of SAH, which could not be observed in the crystal structures. To the best of our knowledge, this study describes for the first time a general catalytic mechanism of bacterial RsmE on theoretical ground.

Effect of Charge State on the Equilibrium and Kinetic Properties of Mechanically Interlocked [5]Rotaxane: A Molecular Dynamics Study

Rotaxanes can exhibit stimuli-responsive behavior by allowing positional fluctuations of their rota groups in response to physiochemical conditions such as the changes in solution pH. However, ionic strength of the solution also affects the molecular conformation by altering the charge state of the entire molecule, coupling the stimuli-responsiveness of rotaxanes with their conformation. A molecular-scale investigation on a model system can allow the decoupling and identification of various effects and can greatly benefit applications of such molecular switches. By using atomistic molecular dynamics simulations, we study equilibrium and kinetics properties of various charge states of the [5]rotaxane, which is a supramolecular moiety with four rotaxanes bonded to a porphyrin core. We model various physiochemical charge states, each of which can be realized at various solution pH levels as well as several exotic charge distributions. By analyzing molecular configurations, hydrogen bonding, and energetics of single molecules in salt-free water and its polyrotaxanated network at the interface of water and chloroform, we demonstrate that charge-neutral and negatively charged molecules often tend to collapse in a way that they can expose their porphyrin core. Contrarily, positively charged moieties tend to take more extended molecular configurations blocking the core. Further, sudden changes in the charge states emulating the pH alterations in solution conditions lead to rapid, sub-10 ns level, changes in the molecular conformation of [5]rotaxane via shuttling motion of CB6 rings along axles. Finally, simulations of 2D [5]rotaxane network structures support our previous findings on a few nanometer-thick film formation at oil-water interfaces. Overall, our results suggest that rotaxane-based structures can exhibit a rich spectrum of molecular configurations and kinetics depending on the ionic strength of the solution.

FSATOOL 2.0: An integrated molecular dynamics simulation and trajectory data analysis program

Molecular dynamics simulation is important in the computational study of the biomolecules. In this paper, we upgrade our previous FSATOOL to version 2.0. It is no longer a plugin as before. Besides the existed enhanced sampling and Markov state model analysis module, FSATOOL 2.0 has three new features now. First, it contains a molecular dynamics simulation engine on both CPU and GPU device. The engine works with an embedded enhanced sampling module. Second, it can do the free energy calculation by various practical methods, including the weighted histogram analysis method and Gaussian mixture model. Third, it has many subroutines to process the trajectory data, such as principal component analysis, time-structure based independent component analysis, contact analysis, and Φ-value analysis. Most importantly, all these calculations are integrated into one package. The trajectory data format is compatible with all the modules. With a proper input parameter file, users can do the molecular dynamics simulation and data analysis work by only a few simplified commands. The capabilities and theoretical backgrounds of FSATOOL 2.0 are introduced in the paper.

Investigating the folding mechanism of the N-terminal domain of ribosomal protein L9

Protein folding is a popular topic in the life science. However, due to the limited sampling ability of experiments and simulations, the general folding mechanism is not yet clear to us. In this work, we study the folding of the N-terminal domain of ribosomal protein L9 (NTL9) in detail by a mixing replica exchange molecular dynamics method. The simulation results are close to previous experimental observations. According to the Markov state model, the folding of the protein follows a nucleation-condensation path. Moreover, after the comparison to its 39-residue β-α-β motif, we find that the helix at the C-terminal has a great influence on the folding process of the intact protein, including the nucleation of the key residues in the transition state ensemble and the packing of the hydrophobic residues in the native state.

A comparative study for the intermediate states of myelin oligodendrocyte glycoprotein in the absence and presence of glycan - A computational approach

Myelin Oligodendrocyte glycoprotein (MOG) is found to play an important role in providing structural integrity to myelin sheath at the same time it acts as an auto-antigen which might lead to Multiple Sclerosis (MS). What causes this specific property of being an auto-antigen is still not known. Here we present molecular dynamics simulation studies of unfolding and folding of the protein MOG in both the absence and presence of N-glycan in order to understand the role of glycosylation in the stability and flexibility of the protein. The main results from these studies show that the glycosylation increases the stability of the protein MOG and inhibits the complete unfolding of MOG in the SMD. From the folding studies using TMD, it was observed that the glycan helps the protein to attain the near-native folded conformation. However, it was also observed from the direct TMD studies that the pathway of protein folding was enhanced by the trace-back of intermediate states in the presence of glycan.

Deep clustering of protein folding simulations

BACKGROUND

We examine the problem of clustering biomolecular simulations using deep learning techniques. Since biomolecular simulation datasets are inherently high dimensional, it is often necessary to build low dimensional representations that can be used to extract quantitative insights into the atomistic mechanisms that underlie complex biological processes.

RESULTS

We use a convolutional variational autoencoder (CVAE) to learn low dimensional, biophysically relevant latent features from long time-scale protein folding simulations in an unsupervised manner. We demonstrate our approach on three model protein folding systems, namely Fs-peptide (14 μs aggregate sampling), villin head piece (single trajectory of 125 μs) and β- β- α (BBA) protein (223 + 102 μs sampling across two independent trajectories). In these systems, we show that the CVAE latent features learned correspond to distinct conformational substates along the protein folding pathways. The CVAE model predicts, on average, nearly 89% of all contacts within the folding trajectories correctly, while being able to extract folded, unfolded and potentially misfolded states in an unsupervised manner. Further, the CVAE model can be used to learn latent features of protein folding that can be applied to other independent trajectories, making it particularly attractive for identifying intrinsic features that correspond to conformational substates that share similar structural features.

CONCLUSIONS

Together, we show that the CVAE model can quantitatively describe complex biophysical processes such as protein folding.

VP40 of the Ebola Virus as a Target for EboV Therapy: Comprehensive Conformational and Inhibitor Binding Landscape from Accelerated Molecular Dynamics

The first account of the dynamic features of the loop region of VP40 of the Ebola virus was studied using accelerated molecular dynamics simulations and reported herein. Among the proteins of the Ebola virus, the matrix protein (VP40) plays a significant role in the virus lifecycle thereby making it a promising therapeutic target. Of interest is the newly elucidated N-terminal domain loop region of VP40 comprising residues K127, T129, and N130 which when mutated to alanine have demonstrated an unrecognized role for N-terminal domain-plasma membrane interaction for efficient VP40-plasma membrane localization, oligomerization, matrix assembly, and egress. The molecular understanding of the conformational features of VP40 in complex with a known inhibitor still remains elusive. Using accelerated molecular dynamics approaches, we conducted a comparative study on VP40 apo and bound systems to understand the conformational features of VP40 at the molecular level and to determine the effect of inhibitor binding with the aid of a number of post-dynamic analytical tools. Significant features were seen in the presence of an inhibitor as per molecular mechanics/generalized born surface area binding free energy calculations. Results revealed that inhibitor binding to VP40 reduces the flexibility and mobility of the protein as supported by root mean square fluctuation and root mean square deviation calculations. The study revealed a characteristic "twisting" motion and coiling of the loop region of VP40 accompanied by conformational changes in the dimer interface upon inhibitor binding. We believe that results presented in this study will ultimately provide useful insight into the binding landscape of VP40 which could assist researchers in the discovery of potent Ebola virus inhibitors for anti-Ebola therapies.

Identification of simple reaction coordinates from complex dynamics

Reaction coordinates are widely used throughout chemical physics to model and understand complex chemical transformations. We introduce a definition of the natural reaction coordinate, suitable for condensed phase and biomolecular systems, as a maximally predictive one-dimensional projection. We then show that this criterion is uniquely satisfied by a dominant eigenfunction of an integral operator associated with the ensemble dynamics. We present a new sparse estimator for these eigenfunctions which can search through a large candidate pool of structural order parameters and build simple, interpretable approximations that employ only a small number of these order parameters. Example applications with a small molecule's rotational dynamics and simulations of protein conformational change and folding show that this approach can filter through statistical noise to identify simple reaction coordinates from complex dynamics.

Mechanism of Inhibition of Hsp90 Dimerization by Gyrase B Inhibitor Coumermycin A1 (C-A1) Revealed by Molecular Dynamics Simulations and Thermodynamic Calculations

Heat shock protein (Hsp) 90 an emerging and attracting target in the anti-HIV drug discovery process due to the key role it plays in the pathogenicity of HIV-1 virus. In this research study, long-range all-atom molecular dynamics simulations were engaged for the bound and the unbound proteins to enhance the understanding of the molecular mechanisms of the Hsp90 dimerization and inhibition. Results evidently showed that coumermycin A1 (C-A1), a recently discovered Hsp90 inhibitor, binds at the dimer's active site of the Hsp90 protein and leads to a substantial parting between dimeric opposed residues, which include Arg591.B, Lys594.A, Ser663.A, Thr653.B, Ala665.A, Thr649.B, Leu646.B and Asn669.A. Significant differences in magnitudes were observed in radius of gyration, root-mean-square deviation and root-mean-square fluctuation, which confirms a reasonably more flexible state in the apo conformation associated with it dimerization. In contrast, the bound conformer of Hsp90 showed less flexibility. This visibly highpoints the inhibition process resulting from the binding of the ligand. These findings were further validated by principal component analysis. We believe that the detailed dynamic analyses of Hsp90 presented in this study, would give an imperative insight and better understanding to the function and mechanisms of inhibition. Furthermore, information obtained from the binding mode of the inhibitor would be of great assistance in the design of more potent inhibitors against the HIV target Hsp90.

Conformational selection and induced fit for RNA polymerase and RNA/DNA hybrid backtracked recognition

RNA polymerase catalyzes transcription with a high fidelity. If DNA/RNA mismatch or DNA damage occurs downstream, a backtracked RNA polymerase can proofread this situation. However, the backtracked mechanism is still poorly understood. Here we have performed multiple explicit-solvent molecular dynamics (MD) simulations on bound and apo DNA/RNA hybrid to study backtracked recognition. MD simulations at room temperature suggest that specific electrostatic interactions play key roles in the backtracked recognition between the polymerase and DNA/RNA hybrid. Kinetics analysis at high temperature shows that bound and apo DNA/RNA hybrid unfold via a two-state process. Both kinetics and free energy landscape analyses indicate that bound DNA/RNA hybrid folds in the order of DNA/RNA contracting, the tertiary folding and polymerase binding. The predicted Φ-values suggest that C7, G9, dC12, dC15, and dT16 are key bases for the backtracked recognition of DNA/RNA hybrid. The average RMSD values between the bound structures and the corresponding apo ones and Kolmogorov-Smirnov (KS) P-test analyses indicate that the recognition between DNA/RNA hybrid and polymerase might follow an induced fit mechanism for DNA/RNA hybrid and conformation selection for polymerase. Furthermore, this method could be used to relative studies of specific recognition between nucleic acid and protein.

Simulations and experiments in protein folding

The interplay between simulations and experiments of protein folding has largely contributed to the elucidation of many important aspects of the phenomenon. In this chapter, I briefly describe the experiments which provide information on the kinetics of the protein folding process, and help to characterize the folding transition state. Then, I show how to probe the kinetics of protein folding using molecular dynamics simulations, how to compare the simulations with the experiments and how to help and rationalize the latter, ultimately offering a molecular picture of the process. After the production of suitable molecular dynamics simulation data in the form of trajectories, the procedure involves sequentially the identification of the stable states of the protein, the identification of the transition pathways connecting the stable states, the identification of the transition state conformations, comparison with experimental results, and finally, the identification of the molecular determinants or reaction coordinates of the folding process, that is, the features that clearly help distinguishing the transition state from the stable states.

Recent developments in methods for identifying reaction coordinates

In the study of rare events in complex systems with many degrees of freedom, a key element is to identify the reaction coordinates of a given process. Over recent years, a number of methods and protocols have been developed to extract the reaction coordinates based on limited information from molecular dynamics simulations. In this review, we provide a brief survey over a number of major methods developed in the past decade, some of which are discussed in greater detail, to provide an overview of the problems that are partially solved and challenges that still remain. A particular emphasis has been placed on methods for identifying reaction coordinates that are related to the committor.

Characterization of the free-energy landscapes of proteins by NMR-guided metadynamics

The use of free-energy landscapes rationalizes a wide range of aspects of protein behavior by providing a clear illustration of the different states accessible to these molecules, as well as of their populations and pathways of interconversion. The determination of the free-energy landscapes of proteins by computational methods is, however, very challenging as it requires an extensive sampling of their conformational spaces. We describe here a technique to achieve this goal with relatively limited computational resources by incorporating nuclear magnetic resonance (NMR) chemical shifts as collective variables in metadynamics simulations. As in this approach the chemical shifts are not used as structural restraints, the resulting free-energy landscapes correspond to the force fields used in the simulations. We illustrate this approach in the case of the third Ig-binding domain of protein G from streptococcal bacteria (GB3). Our calculations reveal the existence of a folding intermediate of GB3 with nonnative structural elements. Furthermore, the availability of the free-energy landscape enables the folding mechanism of GB3 to be elucidated by analyzing the conformational ensembles corresponding to the native, intermediate, and unfolded states, as well as the transition states between them. Taken together, these results show that, by incorporating experimental data as collective variables in metadynamics simulations, it is possible to enhance the sampling efficiency by two or more orders of magnitude with respect to standard molecular dynamics simulations, and thus to estimate free-energy differences among the different states of a protein with a k(B)T accuracy by generating trajectories of just a few microseconds.

Global conformational selection and local induced fit for the recognition between intrinsic disordered p53 and CBP

The transactivation domain (TAD) of tumor suppressor p53 can bind with the nuclear coactivator binding domain (NCBD) of cyclic-AMP response element binding protein (CBP) and activate transcription. NMR experiments demonstrate that both apo-NCBD and TAD are intrinsic disordered and bound NCBD/TAD undergoes a transition to well folded. The recognition mechanism between intrinsic disordered proteins is still hotly debated. Molecular dynamics (MD) simulations in explicit solvent are used to study the recognition mechanism between intrinsic disordered TAD and NCBD. The average RMSD values between bound and corresponding apo states and Kolmogorov-Smirnov P test analysis indicate that TAD and NCBD may follow an induced fit mechanism. Quantitative analysis indicates there is also a global conformational selection. In summary, the recognition of TAD and NCBD might obey a local induced fit and global conformational selection. These conclusions are further supported by high-temperature unbinding kinetics and room temperature landscape analysis. These methods can be used to study the recognition mechanism of other intrinsic disordered proteins.

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.

The redundancy of NMR restraints can be used to accelerate the unfolding behavior of an SH3 domain during molecular dynamics simulations

1 Abstract

Improvement of structure-based potentials for protein folding by native and nonnative hydrogen bonds

Pure Gō models (where every native interaction equally stabilizes the folded state) have widely proved their convenience in the computational investigation of protein folding. However, a chemistry-based description of the real interactions also provides a desirable tune in the analysis of the folding process, and thus some hybrid Gō potentials that combine both aspects have been proposed. Among all the noncovalent interactions that contribute to protein folding, hydrogen bonds are the only ones with a partial covalent character. This feature makes them directional and, thus, more difficult to model as part of the coarse-grained descriptions that are typically employed in Gō models. Thanks to a simplified but rigorous representation of backbone hydrogen bonds that we have recently proposed, we present in this article a combined potential (Gō + backbone hydrogen bond) to study the thermodynamics of protein folding in the frame of very simple simulation models. We show that the explicit inclusion of hydrogen bonds leads to a systematic improvement in the description of protein folding. We discuss a representative set of examples (from two-state folders to downhill proteins, with different types of native structures) that reveal a relevant agreement with experimental data.

Structured pathway across the transition state for peptide folding revealed by molecular dynamics simulations

Small globular proteins and peptides commonly exhibit two-state folding kinetics in which the rate limiting step of folding is the surmounting of a single free energy barrier at the transition state (TS) separating the folded and the unfolded states. An intriguing question is whether the polypeptide chain reaches, and leaves, the TS by completely random fluctuations, or whether there is a directed, stepwise process. Here, the folding TS of a 15-residue β-hairpin peptide, Peptide 1, is characterized using independent 2.5 μs-long unbiased atomistic molecular dynamics (MD) simulations (a total of 15 μs). The trajectories were started from fully unfolded structures. Multiple (spontaneous) folding events to the NMR-derived conformation are observed, allowing both structural and dynamical characterization of the folding TS. A common loop-like topology is observed in all the TS structures with native end-to-end and turn contacts, while the central segments of the strands are not in contact. Non-native sidechain contacts are present in the TS between the only tryptophan (W11) and the turn region (P7-G9). Prior to the TS the turn is found to be already locked by the W11 sidechain, while the ends are apart. Once the ends have also come into contact, the TS is reached. Finally, along the reactive folding paths the cooperative loss of the W11 non-native contacts and the formation of the central inter-strand native contacts lead to the peptide rapidly proceeding from the TS to the native state. The present results indicate a directed stepwise process to folding the peptide.

The folding dynamics of many small protein/peptides investigated recently are in terms of simple two-state model in which only two populations exist (folded and unfolded), separated by a single free energy barrier with only one kinetically important transition state (TS). However, dynamical characterization of the folding TS is challenging. We have used independent unbiased atomistic molecular dynamics simulations with clear folding-unfolding transitions to characterize structural and dynamical features of transition state ensemble of Peptide 1. A common loop-like topology is observed in all TS structures extracted from multiple simulations. The trajectories were used to examine the mechanism by which the TS is reached and subsequent events in folding pathways. The folding TS is reached and crossed in a directed stagewise process rather than through random fluctuations. Specific structures are formed before, during, and after the transition state, indicating a clear structured folding pathway.

Computational model for protein unfolding simulation

The protein folding problem is one of the fundamental and important questions in molecular biology. However, the all-atom molecular dynamics studies of protein folding and unfolding are still computationally expensive and severely limited by the time scale of simulation. In this paper, a simple and fast protein unfolding method is proposed based on the conformational stability analyses and structure modeling. In this method, two structure-based conditions are considered to identify the unstable regions of proteins during the unfolding processes. The protein unfolding trajectories are mimicked through iterative structure modeling according to conformational stability analyses. Two proteins, chymotrypsin inhibitor 2 (CI2) and α -spectrin SH3 domain (SH3) were simulated by this method. Their unfolding pathways are consistent with the previous molecular dynamics simulations. Furthermore, the transition states of the two proteins were identified in unfolding processes and the theoretical Φ values of these transition states showed significant correlations with the experimental data (the correlation coefficients are >0.8). The results indicate that this method is effective in studying protein unfolding. Moreover, we analyzed and discussed the influence of parameters on the unfolding simulation. This simple coarse-grained model may provide a general and fast approach for the mechanism studies of protein folding.

Turbulent phenomena in protein folding

Protein folding and hydrodynamic turbulence are two long-standing challenges, in molecular biophysics and fluid dynamics, respectively. The theories of these phenomena have been developed independently and used different formalisms. Here we show that the protein folding flows can be surprisingly similar to turbulent fluid flows. Studying a benchmark model protein (an SH3 domain), we have found that the flows for the slow folding trajectories of the protein, in which a partly formed N- and C-terminal β sheet hinders the RT loop from attaching to the protein core, have many properties of turbulent flows of a fluid. The flows are analyzed in a three-dimensional (3D) space of collective variables, which are the numbers of native contacts between the terminal β strands, between the RT loop and the protein core, and the rest of the native contacts. We have found that the flows have fractal nature and are filled with 3D eddies; the latter contain strange attractors, at which the tracer flow paths behave as saddle trajectories. Two regions of the space increment have been observed, in which the flux variations are self-similar with the scaling exponent h=1/3, in surprising agreement with the Kolmogorov inertial range theory of turbulence. In one region, the cascade of protein rearrangements is directed from larger to smaller scales (net folding), and in the other, it is oppositely directed (net unfolding). Folding flows for the fast trajectories are essentially "laminar" and do not have the property of self-similarity. Based on the results of our study, we infer, and support this inference by simulations, that the origin of the similarity between the protein folding and turbulent motion of a fluid is in a cascade mechanism of structural transformations in the systems that underlies these phenomena.

Induced fit or conformational selection for RNA/U1A folding

The hairpin II of U1 snRNA can bind U1A protein with high affinity and specificity. NMR spectra suggest that the loop region of apo-RNA is largely unstructured and undergoes a transition from unstructured to well-folded upon U1Abinding. However, the mechanism that RNA folding coupled protein binding is poorly understood. To get an insight into the mechanism, we have performed explicit-solvent molecular dynamics (MD) to study the folding kinetics of bound RNA and apo-RNA. Room-temperature MD simulations suggest that the conformation of bound RNA has significant adjustment and becomes more stable upon U1A binding. Kinetic analysis of high-temperature MD simulations shows that bound RNA and apo-RNA unfold via a two-state process, respectively. Both kinetics and free energy landscape analyses indicate that bound RNA folds in the order of RNA contracting, U1A binding, and tertiary folding. The predicted Phi-values suggest that A8, C10, A11, and G16 are key bases for bound RNA folding. Mutant Arg52Gln analysis shows that electrostatic interaction and hydrogen bonds between RNA and U1A (Arg52Gln) decrease. These results are in qualitative agreement with experiments. Furthermore, this method could be used in other studies about biomolecule folding upon receptor binding.

Temperature-induced unfolding of epidermal growth factor (EGF): insight from molecular dynamics simulation

Thermal disruption of protein structure and function is a potentially powerful therapeutic vehicle. With the emerging nanoparticle-targeting and femtosecond laser technology, it is possible to deliver heating locally to specific molecules. It is therefore important to understand how fast a protein can unfold or lose its function at high temperatures, such as near the water boiling point. In this study, the thermal damage of EGF was investigated by combining the replica exchange (136 replicas) and conventional molecular dynamics simulations. The REMD simulation was employed to rigorously explore the free-energy landscape of EGF unfolding. Interestingly, besides the native and unfolded states, we also observed a distinct molten globule (MG) state that retained substantial amount of native contacts. Based on the understanding that which the unfolding of EGF is a three-state process, we have examined the unfolding kinetics of EGF (N-->MG and MG-->D) with multiple 20-ns conventional MD simulations. The Arrhenius prefactors and activation energy barriers determined from the simulation are within the range of previously studied proteins. In contrast to the thermal damage of cells and tissues which take place on the time scale of seconds to hours at relatively low temperatures, the denaturation of proteins occur in nanoseconds when the temperature of heat bath approaches the boiling point.

Dynameomics: a consensus view of the protein unfolding/folding transition state ensemble across a diverse set of protein folds

The Dynameomics project aims to simulate a representative sample of all globular protein metafolds under both native and unfolding conditions. We have identified protein unfolding transition state (TS) ensembles from multiple molecular dynamics simulations of high-temperature unfolding in 183 structurally distinct proteins. These data can be used to study individual proteins and individual protein metafolds and to mine for TS structural features common across all proteins. Separating the TS structures into four different fold classes (all proteins, all-alpha, all-beta, and mixed alpha/beta and alpha +beta) resulted in no significant difference in the overall protein properties. The residues with the most contacts in the native state lost the most contacts in the TS ensemble. On average, residues beginning in an alpha-helix maintained more structure in the TS ensemble than did residues starting in beta-strands or any other conformation. The metafolds studied here represent 67% of all known protein structures, and this is, to our knowledge, the largest, most comprehensive study of the protein folding/unfolding TS ensemble to date. One might have expected broad distributions in the average global properties of the TS relative to the native state, indicating variability in the amount of structure present in the TS. Instead, the average global properties converged with low standard deviations across metafolds, suggesting that there are general rules governing the structure and properties of the TS.

Insights into protein aggregation by NMR characterization of insoluble SH3 mutants solubilized in salt-free water

Protein aggregation in vivo has been extensively associated with a large spectrum of human diseases. On the other hand, mechanistic insights into protein aggregation in vitro were incomplete due to the inability in solubilizing insoluble proteins for high-resolution biophysical investigations. However, a new avenue may be opened up by our recent discovery that previously-thought insoluble proteins can in fact be solubilized in salt-free water. Here we use this approach to study the NMR structural and dynamic properties of an insoluble SH3 mutant with a naturally-occurring insertion of Val22 at the tip of the diverging turn. The obtained results reveal: 1) regardless of whether the residue is Val, Ala, Asp or Arg, the insertion will render the first hNck2 SH3 domain to be insoluble in buffers. Nevertheless, all four mutants could be solubilized in salt-free water and appear to be largely unfolded as evident from their CD and NMR HSQC spectra. 2) Comparison of the chemical shift deviations reveals that while in V22-SH3 the second helical region is similarly populated as in the wild-type SH3 at pH 2.0, the first helical region is largely unformed. 3) In V22-SH3, many non-native medium-range NOEs manifest to define non-native helical conformations. In the meanwhile a small group of native-like long-range NOEs still persists, indicating the existence of a rudimentary native-like tertiary topology. 4) Although overall, V22-SH3 has significantly increased backbone motions on the ps-ns time scale, some regions still own restricted backbone motions as revealed by analyzing ¹⁵N relaxation data. Our study not only leads to the establishment of the first high-resolution structural and dynamic picture for an insoluble protein, but also shed more light on the molecular events for the nonhierarchical folding mechanism. Furthermore, a general mechanism is also proposed for in vivo protein aggregation triggered by the genetic mutation and posttranslational modification.

The nucleation mechanism of protein folding: a survey of computer simulation studies

The nucleation mechanism of protein folding, originally proposed by Baldwin in the early 1970s, was firstly observed by Shakhnovich and co-workers two decades later in the context of Monte Carlo simulations of a simple lattice model. At about the same time the extensive use of φ-value analysis provided the first experimental evidence that the folding of Chymotrypsin-inhibitor 2, a small single-domain protein, which folds with two-state kinetics, is also driven by a nucleation mechanism. Since then, the nucleation mechanism is generally considered the most common form of folding mechanism amongst two-state proteins. However, recent experimental data has put forward the idea that this may not necessarily be so, since the accuracy of the experimentally determined φ values, which are used to identify the critical (i.e. nucleating) residues, is typically poor. Here, we provide a survey of in silico results on the nucleation mechanism, ranging from simple lattice Monte Carlo to more sophisticated off-lattice molecular dynamics simulations, and discuss them in light of experimental data.

Molecular dynamics simulation of phosphorylated KID post-translational modification

Background

Kinase-inducible domain (KID) as transcriptional activator can stimulate target gene expression in signal transduction by associating with KID interacting domain (KIX). NMR spectra suggest that apo-KID is an unstructured protein. After post-translational modification by phosphorylation, KID undergoes a transition from disordered to well folded protein upon binding to KIX. However, the mechanism of folding coupled to binding is poorly understood.

Methodology

To get an insight into the mechanism, we have performed ten trajectories of explicit-solvent molecular dynamics (MD) for both bound and apo phosphorylated KID (pKID). Ten MD simulations are sufficient to capture the average properties in the protein folding and unfolding.

Conclusions

Room-temperature MD simulations suggest that pKID becomes more rigid and stable upon the KIX-binding. Kinetic analysis of high-temperature MD simulations shows that bound pKID and apo-pKID unfold via a three-state and a two-state process, respectively. Both kinetics and free energy landscape analyses indicate that bound pKID folds in the order of KIX access, initiation of pKID tertiary folding, folding of helix α_B, folding of helix α_A, completion of pKID tertiary folding, and finalization of pKID-KIX binding. Our data show that the folding pathways of apo-pKID are different from the bound state: the foldings of helices α_A and α_B are swapped. Here we also show that Asn139, Asp140 and Leu141 with large Φ-values are key residues in the folding of bound pKID. Our results are in good agreement with NMR experimental observations and provide significant insight into the general mechanisms of binding induced protein folding and other conformational adjustment in post-translational modification.

CHARMM: the biomolecular simulation program

CHARMM (Chemistry at HARvard Molecular Mechanics) is a highly versatile and widely used molecular simulation program. It has been developed over the last three decades with a primary focus on molecules of biological interest, including proteins, peptides, lipids, nucleic acids, carbohydrates, and small molecule ligands, as they occur in solution, crystals, and membrane environments. For the study of such systems, the program provides a large suite of computational tools that include numerous conformational and path sampling methods, free energy estimators, molecular minimization, dynamics, and analysis techniques, and model-building capabilities. The CHARMM program is applicable to problems involving a much broader class of many-particle systems. Calculations with CHARMM can be performed using a number of different energy functions and models, from mixed quantum mechanical-molecular mechanical force fields, to all-atom classical potential energy functions with explicit solvent and various boundary conditions, to implicit solvent and membrane models. The program has been ported to numerous platforms in both serial and parallel architectures. This article provides an overview of the program as it exists today with an emphasis on developments since the publication of the original CHARMM article in 1983.

Characterisation of transition state structures for protein folding using 'high', 'medium' and 'low' {Phi}-values

It has been suggested that Phi-values, which allow structural information about transition states (TSs) for protein folding to be obtained, are most reliably interpreted when divided into three classes (high, medium and low). High Phi-values indicate almost completely folded regions in the TS, intermediate Phi-values regions with a detectable amount of structure and low Phi-values indicate mostly unstructured regions. To explore the extent to which this classification can be used to characterise in detail the structure of TSs for protein folding, we used Phi-values divided into these classes as restraints in molecular dynamics simulations. This type of procedure is related to that used in NMR spectroscopy to define the structure of native proteins from the measurement of inter-proton distances derived from nuclear Overhauser effects. We illustrate this approach by determining the TS ensembles of five proteins and by showing that the results are similar to those obtained by using as restraints the actual numerical Phi-values measured experimentally. Our results indicate that the simultaneous consideration of a set of low-resolution Phi-values can provide sufficient information for characterising the architecture of a TS for folding of a protein.

Abeta-globulomers are formed independently of the fibril pathway

Soluble A beta-oligomers are currently discussed as the major causative species for the development of Alzheimer's disease (AD). Consequently, the beta-amyloid cascade hypothesis was extended by A beta-oligomers and their central neuropathogenic role in AD. However, the molecular structure of A beta-oligomers and their relation to amyloid fibril formation remains elusive. Previously we demonstrated that incubation of A beta(1-42) with SDS or fatty acids induces the formation of a homogeneous globular A beta-oligomer termed A beta-globulomer. In this study we investigated the role of A beta-globulomers in the aggregation pathway of A beta-peptide. We used in vitro assays such as thioflavin-T binding and aggregation inhibitors like Congo red to reveal that A beta-peptide in its A beta-globulomer conformation is a structural entity which is independent from amyloid fibril formation. In addition, cellular Alzheimer's-like plaque forming assays show the resistance of A beta-globulomers to deposition as amyloid plaques. We hypothesize that a conformational switch of A beta is decisive for either fibril formation or alternatively and independently A beta-globulomer formation.

Molecular dynamics simulations from putative transition states of alpha-spectrin SH3 domain

A series of molecular dynamics simulations in explicit solvent were started from nine structural models of the transition state of the SH3 domain of alpha-spectrin, which were generated by Lindorff-Larsen et al. (Nat Struct Mol Biol 2004;11:443-449) using molecular dynamics simulations in which experimental Phi - values were incorporated as restraints. Two of the nine models were simulated 10 times for 200 ns and the remaining models simulated two times for 200 ns. Complete folding was observed in one case, while in the other simulations partial folding or unfolding events were observed, which were characterized by a regularization of elements of secondary structure. These results are consistent with recent experimental evidence that the folding of SH3 domains involves low populated intermediate states.

A one-dimensional reaction coordinate for identification of transition states from explicit solvent P(fold)-like calculations

A properly identified transition state ensemble (TSE) in a molecular dynamics (MD) simulation can reveal a tremendous amount about how a protein folds and offer a point of comparison to experimentally derived Phi(F) values, which reflect the degree of structure in these transient states. In one such method of TSE identification, dubbed P(fold), MD simulations of individual protein structures taken from an unfolding trajectory are used to directly assess an input structure's probability of folding before unfolding, and P(fold) is, by definition, 0.5 for the TSE. Other, less computationally intensive methods, such as multidimensional scaling (MDS) of the pairwise root mean-squared deviation (RMSD) matrix of the conformations sampled in a thermal unfolding trajectory, have also been used to identify the TSE. Identification of the TSE is made from the original MD simulation without the need to run further simulations. Here we present a P(fold)-like study and describe methods for identification of the TSE through the derivation of a high fidelity, bounded, one-dimensional reaction coordinate for protein folding. These methods are applied to the engrailed homeodomain. The TSE identified by this approach is essentially identical to the TSE identified previously by MDS of the pairwise RMSD matrix. However, the cost of performing P(fold), or even our reduced P(fold)-like calculations, is at least 36,000 times greater than the MDS method.

Effect of Different Salt Ions on the Propensity of Aggregation and on the Structure of Alzheimer’s Aβ(1-40) Amyloid Fibrils

The formation of amyloid fibrils and other polypeptide aggregates depends strongly on the physico-chemical environment. One such factor affecting aggregation is the presence and concentration of salt ions. We have examined the effects of salt ions on the aggregation propensity of Alzheimer's Abeta(1-40) peptide and on the structure of the dissolved and of the fibrillar peptide. All salts examined promote aggregation strongly. The most pronounced effect is seen within the cationic series, i.e. for MgCl2. Evaluation of different possible explanations suggests that Abeta(1-40) aggregation depends on direct interaction between ions and Abeta(1-40) peptide, and correlates with ion-induced changes of the surface tension. Salts have profound effects on the fibril structure. In the presence of salts, fibrils are associated with smaller diameters, narrower crossover distances and lower amide I maxima. Since Abeta(1-40) aggregation responds to salts in a manner unlike that for other polypeptides, such as glucagon, beta2-microglobulin or alpha-synuclein; these data argue that there is no fully uniform way in which salts affect aggregation of different polypeptide chains. These observations are important for understanding and predicting aggregation on the basis of simple physico-chemical properties.

Parallel folding pathways in the SH3 domain protein

The transition-state ensemble (TSE) is the set of protein conformations with an equal probability to fold or unfold. Its characterization is crucial for an understanding of the folding process. We determined the TSE of the src-SH3 domain protein by using extensive molecular dynamics simulations of the Go model and computing the folding probability of a generated set of TSE candidate conformations. We found that the TSE possesses a well-defined hydrophobic core with variable enveloping structures resulting from the superposition of three parallel folding pathways. The most preferred pathway agrees with the experimentally determined TSE, while the two least preferred pathways differ significantly. The knowledge of the different pathways allows us to design the interactions between amino acids that guide the protein to fold through the least preferred pathway. This particular design is akin to a circular permutation of the protein. The finding motivates the hypothesis that the different experimentally observed TSEs in homologous proteins and circular permutants may represent potentially available pathways to the wild-type protein.

From coarse-grain to all-atom: Toward multiscale analysis of protein landscapes

Multiscale methods are becoming increasingly promising as a way to characterize the dynamics of large protein systems on biologically relevant time-scales. The underlying assumption in multiscale simulations is that it is possible to move reliably between different resolutions. We present a method that efficiently generates realistic all-atom protein structures starting from the C(alpha) atom positions, as obtained for instance from extensive coarse-grain simulations. The method, a reconstruction algorithm for coarse-grain structures (RACOGS), is validated by reconstructing ensembles of coarse-grain structures obtained during folding simulations of the proteins src-SH3 and S6. The results show that RACOGS consistently produces low energy, all-atom structures. A comparison of the free energy landscapes calculated using the coarse-grain structures versus the all-atom structures shows good correspondence and little distortion in the protein folding landscape.