Why Computed Protein Folding Landscapes Are Sensitive to the Water Model

Computational studies on the catalytic potential of the double active site for enzyme engineering

Proteins possessing double active sites have the potential to revolutionise enzyme design strategies. This study extensively explored an enzyme that contains both a natural active site (NAS) and an engineered active site (EAS), focusing on understanding its structural and functional properties. Metadynamics simulations were employed to investigate how substrates interacted with their respective active sites. The results revealed that both the NAS and EAS exhibited similar minimum energy states, indicating comparable binding affinities. However, it became apparent that the EAS had a weaker binding site for the substrate due to its smaller pocket and constrained conformation. Interestingly, the EAS also displayed dynamic behaviour, with the substrate observed to move outside the pocket, suggesting the possibility of substrate translocation. To gain further insights, steered molecular dynamics (SMD) simulations were conducted to study the conformational changes of the substrate and its interactions with catalytic residues. Notably, the substrate adopted distinct conformations, including near-attack conformations, in both the EAS and NAS. Nevertheless, the NAS demonstrated superior binding minima for the substrate compared to the EAS, reinforcing the observation that the engineered active site was less favourable for substrate binding due to its limitations. The QM/MM (Quantum mechanics and molecular mechanics) analyses highlight the energy disparity between NAS and EAS. Specifically, EAS exhibited elevated energy levels due to its engineered active site being located on the surface. This positioning exposes the substrate to solvents and water molecules, adding to the energy challenge. Consequently, the engineered enzyme did not provide a significant advantage in substrate binding over the single active site protein. Further, the investigation of internal channels and tunnels within the protein shed light on the pathways facilitating transport between the two active sites. By unravelling the complex dynamics and functional characteristics of this double-active site protein, this study offers valuable insights into novel strategies of enzyme engineering. These findings establish a solid foundation for future research endeavours aimed at harnessing the potential of double-active site proteins in diverse biotechnological applications.

Transferable Implicit Solvation via Contrastive Learning of Graph Neural Networks

Implicit solvent models are essential for molecular dynamics simulations of biomolecules, striking a balance between computational efficiency and biological realism. Efforts are underway to develop accurate and transferable implicit solvent models and coarse-grained (CG) force fields in general, guided by a bottom-up approach that matches the CG energy function with the potential of mean force (PMF) defined by the finer system. However, practical challenges arise due to the lack of analytical expressions for the PMF and algorithmic limitations in parameterizing CG force fields. To address these challenges, a machine learning-based approach is proposed, utilizing graph neural networks (GNNs) to represent the solvation free energy and potential contrasting for parameter optimization. We demonstrate the effectiveness of the approach by deriving a transferable GNN implicit solvent model using 600,000 atomistic configurations of six proteins obtained from explicit solvent simulations. The GNN model provides solvation free energy estimations much more accurately than state-of-the-art implicit solvent models, reproducing configurational distributions of explicit solvent simulations. We also demonstrate the reasonable transferability of the GNN model outside of the training data. Our study offers valuable insights for deriving systematically improvable implicit solvent models and CG force fields from a bottom-up perspective.

Visualization of Self-Assembly and Hydration of a β-Hairpin through Integrated Small and Wide-Angle Neutron Scattering

Fundamental understanding of the structure and assembly of nanoscale building blocks is crucial for the development of novel biomaterials with defined architectures and function. However, accessing self-consistent structural information across multiple length scales is challenging. This limits opportunities to exploit atomic scale interactions to achieve emergent macroscale properties. In this work we present an integrative small- and wide-angle neutron scattering approach coupled with computational modeling to reveal the multiscale structure of hierarchically self-assembled β hairpins in aqueous solution across 4 orders of magnitude in length scale from 0.1 Å to 300 nm. Our results demonstrate the power of this self-consistent cross-length scale approach and allows us to model both the large-scale self-assembly and small-scale hairpin hydration of the model β hairpin CLN025. Using this combination of techniques, we map the hydrophobic/hydrophilic character of this model self-assembled biomolecular surface with atomic resolution. These results have important implications for the multiscale investigation of aqueous peptides and proteins, for the prediction of ligand binding and molecular associations for drug design, and for understanding the self-assembly of peptides and proteins for functional biomaterials.

Transferable Coarse Graining via Contrastive Learning of Graph Neural Networks

Coarse-grained (CG) force fields are essential for molecular dynamics simulations of biomolecules, striking a balance between computational efficiency and biological realism. These simulations employ simplified models grouping atoms into interaction sites, enabling the study of complex biomolecular systems over biologically relevant timescales. Efforts are underway to develop accurate and transferable CG force fields, guided by a bottom-up approach that matches the CG energy function with the potential of mean force (PMF) defined by the finer system. However, practical challenges arise due to many-body effects, lack of analytical expressions for the PMF, and limitations in parameterizing CG force fields. To address these challenges, a machine learning-based approach is proposed, utilizing graph neural networks (GNNs) to represent CG force fields and potential contrasting for parameterization from atomistic simulation data. We demonstrate the effectiveness of the approach by deriving a transferable GNN implicit solvent model using 600,000 atomistic configurations of six proteins obtained from explicit solvent simulations. The GNN model provides solvation free energy estimations much more accurately than state-of-the-art implicit solvent models, reproducing configurational distributions of explicit solvent simulations. We also demonstrate the reasonable transferability of the GNN model outside the training data. Our study offers valuable insights for building accurate coarse-grained models bottom-up.

Understanding the role of water on temperature-dependent structural modifications of SARS CoV-2 main protease binding sites

Thermally stable and labile proteases are found in microorganisms. Protease mediates the cleavage of polyproteins in the virus replication and transcription process. 6 µs MD simulations were performed for monomer/dimer SARS CoV-2 main protease system in both SPC/E and mTIP3P water model to analyse the temperature-dependent behaviour of the protein. It is found that maximum conformational changes are observed at 348 K which is near the melting temperature. Network distribution of evolved conformations shows an increase in the number of communities with the rise in the temperature. The global conformation of the protein was found to be intact whereas a local conformational space evolved due to thermal fluctuations. The global conformational change in the free energy ΔΔG value for the monomer and the dimer between 278 K and 383 K is found to be 2.51 and 2.10 kJ/mol respectively. A detailed analysis was carried out on the effect of water on the temperature-dependent structural modifications of four binding pockets of SARS CoV-2 main protease namely, catalytic dyad, substrate-binding site, dimerization site and allosteric site. It is found that the water structure around the binding sites is altered with temperature. The water around the dimer sites is more ordered than the monomer sites regardless of the rise in temperature due to structural rigidity. The energy expense of binding the small molecules at substrate binding is less compared to the allosteric site. The water-water hydrogen bond lifetime is found to be more near the cavity of His41. Also, it is observed that mTIP3P water molecules have a similar effect to that of SPC/E water molecules on the main protease.

Water's Variable Role in Protein Stability Uncovered by Liquid-Observed Vapor Exchange NMR

Water is essential to protein structure and stability, yet our understanding of how water shapes proteins is far from thorough. Our incomplete knowledge of protein-water interactions is due in part to a long-standing technological inability to assess experimentally how water removal impacts local protein structure. It is now possible to obtain residue-level information on dehydrated protein structures via liquid-observed vapor exchange (LOVE) NMR, a solution NMR technique that quantifies the extent of hydrogen-deuterium exchange between unprotected amide protons of a dehydrated protein and D₂O vapor. Here, we apply LOVE NMR, Fourier transform infrared spectroscopy, and solution hydrogen-deuterium exchange to globular proteins GB1, CI2, and two variants thereof to link mutation-induced changes in the dehydrated protein structure to changes in solution structure and stability. We find that a mutation that destabilizes GB1 in solution does not affect its dehydrated structure, whereas a mutation that stabilizes CI2 in solution makes several regions of the protein more susceptible to dehydration-induced unfolding, suggesting that water is primarily responsible for the destabilization of the GB1 variant but plays a stabilizing role in the CI2 variant. Our results indicate that changes in dehydrated protein structure cannot be predicted from changes in solution stability alone and demonstrate the ability of LOVE NMR to uncover the variable role of water in protein stability. Further application of LOVE NMR to other proteins and their variants will improve the ability to predict and modulate protein structure and stability in both the hydrated and dehydrated states for applications in medicine and biotechnology.

Systematic Comparison of the Structural and Dynamic Properties of Commonly Used Water Models for Molecular Dynamics Simulations

Water is a unique solvent that is ubiquitous in biology and present in a variety of solutions, mixtures, and materials settings. It therefore forms the basis for all molecular dynamics simulations of biological phenomena, as well as for many chemical, industrial, and materials investigations. Over the years, many water models have been developed, and it remains a challenge to find a single water model that accurately reproduces all experimental properties of water simultaneously. Here, we report a comprehensive comparison of structural and dynamic properties of 30 commonly used 3-point, 4-point, 5-point, and polarizable water models simulated using consistent settings and analysis methods. For the properties of density, coordination number, surface tension, dielectric constant, self-diffusion coefficient, and solvation free energy of methane, models published within the past two decades consistently show better agreement with experimental values compared to models published earlier, albeit with some notable exceptions. However, no single model reproduced all experimental values exactly, highlighting the need to carefully choose a water model for a particular study, depending on the phenomena of interest. Finally, machine learning algorithms quantified the relationship between the water model force field parameters and the resulting bulk properties, providing insight into the parameter-property relationship and illustrating the challenges of developing a water model that can accurately reproduce all properties of water simultaneously.

Machine learning implicit solvation for molecular dynamics

Accurate modeling of the solvent environment for biological molecules is crucial for computational biology and drug design. A popular approach to achieve long simulation time scales for large system sizes is to incorporate the effect of the solvent in a mean-field fashion with implicit solvent models. However, a challenge with existing implicit solvent models is that they often lack accuracy or certain physical properties compared to explicit solvent models as the many-body effects of the neglected solvent molecules are difficult to model as a mean field. Here, we leverage machine learning (ML) and multi-scale coarse graining (CG) in order to learn implicit solvent models that can approximate the energetic and thermodynamic properties of a given explicit solvent model with arbitrary accuracy, given enough training data. Following the previous ML-CG models CGnet and CGSchnet, we introduce ISSNet, a graph neural network, to model the implicit solvent potential of mean force. ISSNet can learn from explicit solvent simulation data and be readily applied to molecular dynamics simulations. We compare the solute conformational distributions under different solvation treatments for two peptide systems. The results indicate that ISSNet models can outperform widely used generalized Born and surface area models in reproducing the thermodynamics of small protein systems with respect to explicit solvent. The success of this novel method demonstrates the potential benefit of applying machine learning methods in accurate modeling of solvent effects for in silico research and biomedical applications.

Polarizable and non-polarizable force fields: Protein folding, unfolding, and misfolding

Molecular dynamics simulations are an invaluable tool to characterize the dynamic motions of proteins in atomistic detail. However, the accuracy of models derived from simulations inevitably relies on the quality of the underlying force field. Here, we present an evaluation of current non-polarizable and polarizable force fields (AMBER ff14SB, CHARMM 36m, GROMOS 54A7, and Drude 2013) based on the long-standing biophysical challenge of protein folding. We quantify the thermodynamics and kinetics of the β-hairpin formation using Markov state models of the fast-folding mini-protein CLN025. Furthermore, we study the (partial) folding dynamics of two more complex systems, a villin headpiece variant and a WW domain. Surprisingly, the polarizable force field in our set, Drude 2013, consistently leads to destabilization of the native state, regardless of the secondary structure element present. All non-polarizable force fields, on the other hand, stably characterize the native state ensembles in most cases even when starting from a partially unfolded conformation. Focusing on CLN025, we find that the conformational space captured with AMBER ff14SB and CHARMM 36m is comparable, but the ensembles from CHARMM 36m simulations are clearly shifted toward disordered conformations. While the AMBER ff14SB ensemble overstabilizes the native fold, CHARMM 36m and GROMOS 54A7 ensembles both agree remarkably well with experimental state populations. In addition, GROMOS 54A7 also reproduces experimental folding times most accurately. Our results further indicate an over-stabilization of helical structures with AMBER ff14SB. Nevertheless, the presented investigations strongly imply that reliable (un)folding dynamics of small proteins can be captured in feasible computational time with current additive force fields.

Implementation of the Freely Jointed Chain Model to Assess Kinetics and Thermodynamics of Thermosensitive Coil-Globule Transition by Markov States

We revived and implemented a method developed by Kuhn in 1934, originally only published in German, that is, the so-called “freely jointed chain” model. This approach turned out to be surprisingly useful for analyzing state-of-the-art computer simulations of the thermosensitive coil–globule transition of N-Isopropylacrylamide 20-mer. Our atomistic computer simulations are orders of magnitude longer than those of previous studies and lead to a reliable description of thermodynamics and kinetics at many different temperatures. The freely jointed chain model provides a coordinate system, which allows us to construct a Markov state model of the conformational transitions. Furthermore, this guarantees a reliable reconstruction of the kinetics in back-and-forth directions. In addition, we obtain a description of the high diversity and variability of both conformational states. Thus, we gain a detailed understanding of the coil–globule transition. Surprisingly, conformational entropy turns out to play only a minor role in the thermodynamic balance of the process. Moreover, we show that the radius of gyration is an unexpectedly unsuitable coordinate to comprehend the transition kinetics because it does not capture the high conformational diversity within the different states. Consequently, the approach presented here allows for an exhaustive description and resolution of the conformational ensembles of arbitrary linear polymer chains.

Effect of Water Models on Transmembrane Self-Assembled Cyclic Peptide Nanotubes

Self-assembling cyclic peptide nanotubes can form nanopores when they are inserted in lipid bilayers, acting as ion and/or water permeable channels. In order to improve the versatility of these systems, it is possible to specifically design cyclic peptides with a combination of natural and non-natural amino acids, enabling the control of the nature of the inner cavity of the channels. Here, the behavior of two types of self-assembling peptide motifs, alternating α-amino acids with γ- or δ-aminocycloalkanecarboxylic acids, is studied via molecular dynamics (MD) simulations. The behavior of water molecules in nanopores is expected to affect the properties of these channels and therefore merits detailed examination. A number of water models commonly used in MD simulations have been validated by how well they reproduce bulk water properties. However, it is less clear how these water models behave in the nanoconfined condition inside a channel. The behavior of four different water models-TIP3P, TIP4P, TIP4P/2005, and OPC-are evaluated in MD simulations of self-assembled cyclic peptide nanotubes of distinct composition and diameter. The dynamic behavior of the water molecules and ions in these designed artificial channels depends subtly on the water model used. TIP3P water molecules move faster than those of TIP4P, TIP4P/2005, and OPC. This demeanor is clearly observed in the filling of the nanotube, in water diffusion within the pore, and in the number and stability of hydrogen bonds of the peptides with water. It was also shown that the water model influences the simulated ion flux through the nanotubes, with TIP3P producing the greatest ion flux. Additionally, the two more recent models, TIP4P/2005 and OPC, which are known to reproduce the experimental self-diffusion coefficient of bulk water quite well, exhibit very similar results under the nanoconfined conditions studied here. Because none of these models have been parametrized specifically for waters confined in peptide nanotubes, this study provides a point of reference for further validation.

Sampling of the conformational landscape of small proteins with Monte Carlo methods

Computer simulation provides an increasingly realistic picture of large-scale conformational change of proteins, but investigations remain fundamentally constrained by the femtosecond timestep of molecular dynamics simulations. For this reason, many biologically interesting questions cannot be addressed using accessible state-of-the-art computational resources. Here, we report the development of an all-atom Monte Carlo approach that permits the modelling of the large-scale conformational change of proteins using standard off-the-shelf computational hardware and standard all-atom force fields. We demonstrate extensive thermodynamic characterization of the folding process of the α-helical Trp-cage, the Villin headpiece and the β-sheet WW-domain. We fully characterize the free energy landscape, transition states, energy barriers between different states, and the per-residue stability of individual amino acids over a wide temperature range. We demonstrate that a state-of-the-art intramolecular force field can be combined with an implicit solvent model to obtain a high quality of the folded structures and also discuss limitations that still remain.

Significant compaction of H4 histone tail upon charge neutralization by acetylation and its mimics, possible effects on chromatin structure

The intrinsically disordered, positively charged H4 histone tail is important for chromatin structure and function. We have explored conformational ensembles of human H4 tail in solution, with varying levels of charge neutralization via acetylation or amino-acid substitutions such as K→Q. We have employed an explicit water model shown recently to be well suited for simulations of intrinsically disordered proteins. Upon progressive neutralization of the H4, its radius of gyration decreases linearly with the tail charge q, the trend is explained using a simple polymer model. While the wild type state (q=+8) is essentially a random coil, hyper-acetylated H4 (q=+3) is virtually as compact and stable as a globular protein of the same number of amino-acids. Conformational ensembles of acetylated H4 match the corresponding K→X substitutions only approximately: based on the ensemble similarity, we propose K→M as a possible alternative to the commonly used K→Q. Possible effects of the H4 tail compaction on chromatin structure are discussed within a qualitative model in which the chromatin is highly heterogeneous, easily inter-converting between various structural forms. We predict that upon progressive charge neutralization of the H4 tail, the least compact sub-states of chromatin de-condense first, followed by de-condensation of more compact structures, e.g. those that harbor a high fraction of stacked di-nucleosomes. The predicted hierarchy of DNA accessibility increase upon progressive acetylation of H4 might be utilized by the cell for selective DNA accessibility control.

Water in Nanopores and Biological Channels: A Molecular Simulation Perspective

This Review explores the dynamic behavior of water within nanopores and biological channels in lipid bilayer membranes. We focus on molecular simulation studies, alongside selected structural and other experimental investigations. Structures of biological nanopores and channels are reviewed, emphasizing those high-resolution crystal structures, which reveal water molecules within the transmembrane pores, which can be used to aid the interpretation of simulation studies. Different levels of molecular simulations of water within nanopores are described, with a focus on molecular dynamics (MD). In particular, models of water for MD simulations are discussed in detail to provide an evaluation of their use in simulations of water in nanopores. Simulation studies of the behavior of water in idealized models of nanopores have revealed aspects of the organization and dynamics of nanoconfined water, including wetting/dewetting in narrow hydrophobic nanopores. A survey of simulation studies in a range of nonbiological nanopores is presented, including carbon nanotubes, synthetic nanopores, model peptide nanopores, track-etched nanopores in polymer membranes, and hydroxylated and functionalized nanoporous silica. These reveal a complex relationship between pore size/geometry, the nature of the pore lining, and rates of water transport. Wider nanopores with hydrophobic linings favor water flow whereas narrower hydrophobic pores may show dewetting. Simulation studies over the past decade of the behavior of water in a range of biological nanopores are described, including porins and β-barrel protein nanopores, aquaporins and related polar solute pores, and a number of different classes of ion channels. Water is shown to play a key role in proton transport in biological channels and in hydrophobic gating of ion channels. An overall picture emerges, whereby the behavior of water in a nanopore may be predicted as a function of its hydrophobicity and radius. This informs our understanding of the functions of diverse channel structures and will aid the design of novel nanopores. Thus, our current level of understanding allows for the design of a nanopore which promotes wetting over dewetting or vice versa. However, to design a novel nanopore, which enables fast, selective, and gated flow of water de novo would remain challenging, suggesting a need for further detailed simulations alongside experimental evaluation of more complex nanopore systems.

Using the generalized Born surface area model to fold proteins yields more effective sampling while qualitatively preserving the folding landscape

Protein folding is a long-standing problem and has been widely investigated using molecular dynamics simulations with both explicit and implicit solvents. However, to what extent the folding mechanisms observed in two water models agree remains an open question. In this study, ab initio folding simulations of ten proteins with different topologies are performed in two combinations of force fields and water models (ff14SB+TIP3P and ff14SBonlysc+GB-Neck2). Interestingly, the latter combination not only folds more proteins but also provides a better balance of different secondary structures than the former in the same number of integration time steps. More importantly, the folding pathways found in the two types of simulations are conserved and they may only differ in their weights. Our results suggest that simulations with an implicit solvent may also be suitable for the investigation of the mechanism of protein folding.

Pressure dependence of vibrational optical activity of model biomolecules. A computational study

Change of molecular properties with pressure is an attracting means to regulate molecular reactivity or biological activity. However, the effect is usually small and so far explored rather scarcely. To obtain a deeper insight and estimate the sensitivity of vibrational optical activity spectra to pressure-induced conformational changes, we investigate small model molecules. The Ala-Ala dipeptide, isomaltose disaccharide and adenine-uracil dinucleotide were chosen to represent three different biomolecular classes. The pressure effects were modeled by molecular dynamics and density functional theory simulations. The dinucleotide was found to be the most sensitive to the pressure, whereas for the disaccharide the smallest changes are predicted. Pressure-induced relative intensity changes in vibrational circular dichroism and Raman optical activity spectra are predicted to be 2-3-times larger than for non-polarized IR and Raman techniques.

Split the Charge Difference in Two! A Rule of Thumb for Adding Proper Amounts of Ions in MD Simulations

Despite the relevance of properly setting ionic concentrations in Molecular Dynamics (MD) simulations, methods or practical rules to set ionic strength are scarce and rarely documented. Based on a recently proposed thermodynamics method we provide an accurate rule of thumb to define the electrolytic content in simulation boxes. Extending the use of good practices in setting up MD systems is promptly needed to ensure reproducibility and consistency in molecular simulations.

Exploring Conformational Change of Adenylate Kinase by Replica Exchange Molecular Dynamic Simulation

Replica exchange molecular dynamics (REMD) simulation is a popular enhanced sampling method that is widely used for exploring the atomic mechanism of protein conformational change. However, the requirement of huge computational resources for REMD, especially with the explicit solvent model, largely limits its application. In this study, the availability and efficiency of a variant of velocity-scaling REMD (vsREMD) was assessed with adenylate kinase as an example. Although vsREMD achieved results consistent with those from conventional REMD and experimental studies, the number of replicas required for vsREMD (30) was much less than that for conventional REMD (80) to achieve a similar acceptance rate (∼0.2), demonstrating high efficiency of vsREMD to characterize the protein conformational change and associated free-energy profile. Thus, vsREMD is a highly efficient approach for studying the large-scale conformational change of protein systems.

Exploring optimization strategies for improving explicit water models: Rigid n-point model and polarizable model based on Drude oscillator

Rigid n-point water models are widely used in atomistic simulations, but have known accuracy drawbacks. Increasing the number of point charges, as well as adding electronic polarizability, are two common strategies for accuracy improvements. Both strategies come at considerable computational cost, which weighs heavily against modest possible accuracy improvements in practical simulations. In an effort to provide guidance for model development, here we have explored the limiting accuracy of "electrostatically globally optimal" n-point water models in terms of their ability to reproduce properties of water dimer-a mimic of the condensed state of water. For a given n, each model is built upon a set of reference multipole moments (e.g. ab initio) and then optimized to reproduce water dimer total dipole moment. The models are then evaluated with respect to the accuracy of reproducing the geometry of the water dimer. We find that global optimization of the charge distribution alone can deliver high accuracy of the water model: for n = 4 or n = 5, the geometry of the resulting water dimer can be almost within 50 of the ab initio reference, which is half that of the experimental error margin. Thus, global optimization of the charge distribution of classical n-point water models can lead to high accuracy models. We also find that while the accuracy improvement in going from n = 3 to n = 4 is substantial, the additional accuracy increase in going from n = 4 to n = 5 is marginal. Next, we have explored accuracy limitations of the standard practice of adding electronic polarizability (via a Drude particle) to a "rigid base"-pre-optimization rigid n-point water model. The resulting model (n = 3) shows a relatively small improvement in accuracy, suggesting that the strategy of merely adding the polarizability to an inferior accuracy water model used as the base cannot fix the defects of the latter. An alternative strategy in which the parameters of the rigid base model are globally optimized along with the polarizability parameter is much more promising: the resulting 3-point polarizable model out-performs even the 5-point optimal rigid model by a large margin. We suggest that future development efforts consider 3- and 4-point polarizable models where global optimization of the "rigid base" is coupled to optimization of the polarizability to deliver globally optimal solutions.

Bayesian-Maximum-Entropy Reweighting of IDP Ensembles Based on NMR Chemical Shifts

Bayesian and Maximum Entropy approaches allow for a statistically sound and systematic fitting of experimental and computational data. Unfortunately, assessing the relative confidence in these two types of data remains difficult as several steps add unknown error. Here we propose the use of a validation-set method to determine the balance, and thus the amount of fitting. We apply the method to synthetic NMR chemical shift data of an intrinsically disordered protein. We show that the method gives consistent results even when other methods to assess the amount of fitting cannot be applied. Finally, we also describe how the errors in the chemical shift predictor can lead to an incorrect fitting and how using secondary chemical shifts could alleviate this problem.

Computational Estimation of Microsecond to Second Atomistic Folding Times

Despite the development of massively parallel computing hardware including inexpensive graphics processing units (GPUs), it has remained infeasible to simulate the folding of atomistic proteins at room temperature using conventional molecular dynamics (MD) beyond the microsecond scale. Here, we report the folding of atomistic, implicitly solvated protein systems with folding times τ ranging from ∼10 μs to ∼100 ms using the weighted ensemble (WE) strategy in combination with GPU computing. Starting from an initial structure or set of structures, WE organizes an ensemble of GPU-accelerated MD trajectory segments via intermittent pruning and replication events to generate statistically unbiased estimates of rate constants for rare events such as folding; no biasing forces are used. Although the variance among atomistic WE folding runs is significant, multiple independent runs are used to reduce and quantify statistical uncertainty. Folding times are estimated directly from WE probability flux and from history-augmented Markov analysis of the WE data. Three systems were examined: NTL9 at low solvent viscosity (yielding τf = 0.8-9 μs), NTL9 at water-like viscosity (τf = 0.2-2 ms), and Protein G at low viscosity (τf = 3-200 ms). In all cases, the folding time, uncertainty, and ensemble properties could be estimated from WE simulation; for Protein G, this characterization required significantly less overall computing than would be required to observe a single folding event with conventional MD simulations. Our results suggest that the use and calibration of force fields and solvent models for precise estimation of kinetic quantities is becoming feasible.

General Purpose Water Model Can Improve Atomistic Simulations of Intrinsically Disordered Proteins

Unconstrained atomistic simulations of intrinsically disordered proteins and peptides (IDP) remain a challenge: widely used, "general purpose" water models tend to favor overly compact structures relative to experiment. Here we have performed a total of 93 μs of unrestrained MD simulations to explore, in the context of IDPs, a recently developed "general-purpose" 4-point rigid water model OPC, which describes liquid state of water close to experiment. We demonstrate that OPC, together with a popular AMBER force field ff99SB, offers a noticeable improvement over TIP3P in producing more realistic structural ensembles of three common IDPs benchmarks: 55-residue apo N-terminal zinc-binding domain of HIV-1 integrase ("protein IN"), amyloid β-peptide (Aβ₄₂) (residues 1-42), and 26-reside H4 histone tail. As a negative control, computed folding profile of a regular globular miniprotein (CLN025) in OPC water is in appreciably better agreement with experiment than that obtained in TIP3P, which tends to overstabilize the compact native state relative to the extended conformations. We employed Aβ₄₂ peptide to investigate the possible influence of the solvent box size on simulation outcomes. We advocate a cautious approach for simulations of IDPs: we suggest that the solvent box size should be at least four times the radius of gyration of the random coil corresponding to the IDP. The computed free energy landscape of protein IN in OPC resembles a shallow "tub" - conformations with substantially different degrees of compactness that are within 2 k_B T of each other. Conformations with very different secondary structure content coexist within 1 k_B T of the global free energy minimum. States with higher free energy tend to have less secondary structure. Computed low helical content of the protein has virtually no correlation with its degree of compactness, which calls into question the possibility of using the helicity as a metric for assessing performance of water models for IDPs, when the helicity is low. Predicted radius of gyration ( R _g) of H4 histone tail in OPC water falls in-between that of a typical globular protein and a fully denatured protein of the same size; the predicted R _g is consistent with two independent predictions. In contrast, H4 tail in TIP3P water is as compact as the corresponding globular protein. The computed free energy landscape of H4 tail in OPC is relatively flat over a significant range of compactness, which, we argue, is consistent with its biological function as facilitator of internucleosome interactions.