Free energy and hidden barriers of the β-sheet structure of prion protein

Unbiasing Enhanced Sampling on a High-Dimensional Free Energy Surface with a Deep Generative Model

Biased enhanced sampling methods that utilize collective variables (CVs) are powerful tools for sampling conformational ensembles. Due to their large intrinsic dimensions, efficiently generating conformational ensembles for complex systems requires enhanced sampling on high-dimensional free energy surfaces. While temperature-accelerated molecular dynamics (TAMD) can trivially adopt many CVs in a simulation, unbiasing the simulation to generate unbiased conformational ensembles requires accurate modeling of a high-dimensional CV probability distribution, which is challenging for traditional density estimation techniques. Here we propose an unbiasing method based on the score-based diffusion model, a deep generative learning method that excels in density estimation across complex data landscapes. We demonstrate that this unbiasing approach, tested on multiple TAMD simulations, significantly outperforms traditional unbiasing methods and can generate accurate unbiased conformational ensembles. With the proposed approach, TAMD can adopt CVs that focus on improving sampling efficiency and the proposed unbiasing method enables accurate evaluation of ensemble averages of important chemical features.

Molecular Dynamics Simulations of Ion Permeation in Human Voltage-Gated Sodium Channels

The recent determination of cryo-EM structures of voltage-gated sodium (Na_v) channels has revealed many details of these proteins. However, knowledge of ionic permeation through the Na_v pore remains limited. In this work, we performed atomistic molecular dynamics (MD) simulations to study the structural features of various neuronal Na_v channels based on homology modeling of the cryo-EM structure of the human Na_v1.4 channel and, in addition, on the recently resolved configuration for Na_v1.2. In particular, single Na⁺ permeation events during standard MD runs suggest that the ion resides in the inner part of the Na_v selectivity filter (SF). On-the-fly free energy parametrization (OTFP) temperature-accelerated molecular dynamics (TAMD) was also used to calculate two-dimensional free energy surfaces (FESs) related to single/double Na⁺ translocation through the SF of the homology-based Na_v1.2 model and the cryo-EM Na_v1.2 structure, with different realizations of the DEKA filter domain. These additional simulations revealed distinct mechanisms for single and double Na⁺ permeation through the wild-type SF, which has a charged lysine in the DEKA ring. Moreover, the configurations of the ions in the SF corresponding to the metastable states of the FESs are specific for each SF motif. Overall, the description of these mechanisms gives us new insights into ion conduction in human Na_v cryo-EM-based and cryo-EM configurations that could advance understanding of these systems and how they differ from potassium and bacterial Na_v channels.

Identifying hydrophobic protein patches to inform protein interaction interfaces

Interactions between proteins lie at the heart of numerous biological processes and are essential for the proper functioning of the cell. Although the importance of hydrophobic residues in driving protein interactions is universally accepted, a characterization of protein hydrophobicity, which informs its interactions, has remained elusive. The challenge lies in capturing the collective response of the protein hydration waters to the nanoscale chemical and topographical protein patterns, which determine protein hydrophobicity. To address this challenge, here, we employ specialized molecular simulations wherein water molecules are systematically displaced from the protein hydration shell; by identifying protein regions that relinquish their waters more readily than others, we are then able to uncover the most hydrophobic protein patches. Surprisingly, such patches contain a large fraction of polar/charged atoms and have chemical compositions that are similar to the more hydrophilic protein patches. Importantly, we also find a striking correspondence between the most hydrophobic protein patches and regions that mediate protein interactions. Our work thus establishes a computational framework for characterizing the emergent hydrophobicity of amphiphilic solutes, such as proteins, which display nanoscale heterogeneity, and for uncovering their interaction interfaces.

Elucidating Solution Structures of Cyclic Peptides Using Molecular Dynamics Simulations

Protein-protein interactions are vital to biological processes, but the shape and size of their interfaces make them hard to target using small molecules. Cyclic peptides have shown promise as protein-protein interaction modulators, as they can bind protein surfaces with high affinity and specificity. Dozens of cyclic peptides are already FDA approved, and many more are in various stages of development as immunosuppressants, antibiotics, antivirals, or anticancer drugs. However, most cyclic peptide drugs so far have been natural products or derivatives thereof, with de novo design having proven challenging. A key obstacle is structural characterization: cyclic peptides frequently adopt multiple conformations in solution, which are difficult to resolve using techniques like NMR spectroscopy. The lack of solution structural information prevents a thorough understanding of cyclic peptides' sequence-structure-function relationship. Here we review recent development and application of molecular dynamics simulations with enhanced sampling to studying the solution structures of cyclic peptides. We describe novel computational methods capable of sampling cyclic peptides' conformational space and provide examples of computational studies that relate peptides' sequence and structure to biological activity. We demonstrate that molecular dynamics simulations have grown from an explanatory technique to a full-fledged tool for systematic studies at the forefront of cyclic peptide therapeutic design.

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.

Roles of variable linker length in dual acting virucidal entry inhibitors on HIV-1 potency via on-the-fly free energy molecular simulations

The Dual-Acting Virolytic Entry Inhibitors, or DAVEI's, are a class of recombinant chimera fusion proteins consisting of a lectin, a flexible polypeptide linker, and a fragment of the membrane-proximal external region (MPER) of HIV-1 gp41. DAVEIs trigger virolysis of HIV-1 virions through interactions with the trimeric envelope glycoprotein complex (Env), though the details of these interactions are not fully determined as yet. The purpose of this work was to use structural modeling to rationalize a dependence of DAVEI potency on the molecular length of the linker connecting the two components. We used temperature accelerated molecular dynamics and on-the-fly parameterization to compute free energy versus end-to-end distance for two different linker lengths, DAVEI L0 (His6 ) and DAVEI L2 ([Gly4 Ser]2 His6 ). Additionally, an envelope model was created based on a cryo-electron microscopy-derived structure of a cleaved, soluble Env construct, with high-mannose glycans added which served as putative docking locations for the lectin, along with MPER added that served as a putative docking location for the MPER region of DAVEI (MPERDAVEI ). Using MD simulation, distances between the lectin C-terminus and Env gp41 MPER were measured. We determined that none of the glycans were close enough to gp41 MPER to allow DAVEI L0 to function, while one, N448, will allow DAVEI L2 to function. These findings are consistent with the previously determined dependence of lytic function on DAVEI linker lengths. This supports the hypothesis that DAVEI's engage Env at both glycans and the Env MPER in causing membrane poration and lysis.

Unraveling the Molecular Mechanism of Prion H2 C-Terminus Misfolding by Metadynamics Simulations

Conformational transition from the normal cellular form of prion protein (PrP^C) to the pathogenic "scrapie" form (PrP^Sc) is considered to be a key event in the occurrence of prion disease. Additionally, the H2 C-terminus is widely considered to be a vital site for PrP conformational transition, which can be used as an important region to explore the potential mechanism of PrP misfolding. Therefore, to study the misfolding mechanism of PrP, 500 ns well-tempered metadynamics simulations were performed by focusing on the H2 C-terminus of PrP. For comparison, three systems were designed in total, including PrP in neutral and acidic conditions, as well as H187R mutant. The resulting free energy surfaces (FESs) obtained from metadynamics simulations reveal that acidic conditions and H187R mutation can facilitate PrP misfolding by decreasing free energy barriers for conformational transition and forming energy stable conformational states. Further analyses aimed at H2 C-terminus show that due to the increase of positive charge on residue 187 in both acidic and H187R systems, the electrostatic repulsion of residue 187 and R136/R156 increases greatly, which disrupts the electrostatic interaction network around H2 C-terminus and exposes the hydrophobic core to the solvent. Taken together, acidic conditions and H187R mutation can accelerate PrP misfolding mainly by forming more energetically stable metastable conformations with lower free energy barriers, and electrostatic network disruption involving residue 187 drives the initial misfolding of H2 C-terminus. This study provides quantitative insight into the related function of the H2 C-terminus in the PrP misfolding process, which may guide H2 C-terminus mediated drug design in the future.

Combining the biased and unbiased sampling strategy into one convenient free energy calculation method

Constructing a free energy landscape for a large molecule is difficult. One has to use either a high temperature or a strong driving force to enhance the sampling on the free energy barriers. In this work, we propose a mixed method that combines these two kinds of acceleration strategies into one simulation. First, it applies an adaptive biasing potential to some replicas of the molecule. These replicas are particularly accelerated in a collective variable space. Second, it places some unbiased and exchangeable replicas at various temperature levels. These replicas generate unbiased sampling data in the canonical ensemble. To improve the sampling efficiency, biased replicas transfer their state variables to the unbiased replicas after equilibrium by Monte Carlo trial moves. In comparison to previous integrated methods, it is more convenient for users. It does not need an initial reference biasing potential to guide the sampling of the molecule. And it is also unnecessary to insert many replicas for the requirement of passing the free energy barriers. The free energy calculation is accomplished in a single stage. It samples the data as fast as a biased simulation and it processes the data as simple as an unbiased simulation. The method provides a minimalist approach to the construction of the free energy landscape. © 2019 Wiley Periodicals, Inc.

Polymorphism at 129 dictates metastable conformations of the human prion protein N-terminal β-sheet

We study the thermodynamic stability of the native state of the human prion protein using a new free-energy method, replica-exchange on-the-fly parameterization. This method is designed to overcome hidden-variable sampling limitations to yield nearly error-free free-energy profiles along a conformational coordinate. We confirm that all four (M129V, D178N) polymorphs have a ground-state conformation with three intact β-sheet hydrogen bonds. Additionally, they are observed to have distinct metastabilities determined by the side-chain at position 129. We rationalize these findings with reference to the prion "strain" hypothesis, which links the variety of transmissible spongiform encephalopathy phenotypes to conformationally distinct infectious prion forms and classifies distinct phenotypes of sporadic Creutzfeldt-Jakob disease based solely on the 129 polymorphism. Because such metastable structures are not easily observed in structural experiments, our approach could potentially provide new insights into the conformational origins of prion diseases and other pathologies arising from protein misfolding and aggregation.

Mechanism of Unfolding of Human Prion Protein

Misfolding and aggregation of prion proteins are associated with several neurodegenerative diseases. Therefore, understanding the mechanism of the misfolding process is of enormous interest in the scientific community. It has been speculated and widely discussed that the native cellular prion protein (PrP^C) form needs to undergo substantial unfolding to a more stable PrP^C* state, which may further oligomerize into the toxic scrapie (PrP^Sc) form. Here, we have studied the mechanism of the unfolding of the human prion protein (huPrP) using a set of extensive well-tempered metadynamics simulations. Through multiple microsecond-long metadynamics simulations, we find several possible unfolding pathways. We show that each pathway leads to an unfolded state of lower free energy than the native state. Thus, our study may point to the signature of a PrP^C* form that corresponds to a global minimum on the conformational free-energy landscape. Moreover, we find that these global minima states do not involve an increased β-sheet content, as was assumed to be a signature of PrP^Sc formation in previous simulation studies. We have further analyzed the origin of metastability of the PrP^C form through free-energy surfaces of the chopped helical segments to show that the helices, particularly H2 and H3 of the prion protein, have the tendency to form either a random coil or a β-structure. Therefore, the secondary structural elements of the prion protein are only weakly stabilized by tertiary contacts and solvation forces so that relatively weak perturbations induced by temperature, pressure, pH, and so forth can lead to substantial unfolding with characteristics of intrinsically disordered proteins.

Calculation of the Local Free Energy Landscape in the Restricted Region by the Modified Tomographic Method

The free energy landscape is the most important information in the study of the reaction mechanisms of the molecules. However, it is difficult to calculate. In a large collective variable space, a molecule must take a long time to obtain the sufficient sampling during the simulation. To save the calculation quantity, decreasing the sampling region and constructing the local free energy landscape is required in practice. However, the restricted region in the collective variable space may have an irregular shape. Simply restricting one or more collective variables of the molecule cannot satisfy the requirement. In this paper, we propose a modified tomographic method to perform the simulation. First, it divides the restricted region by some hyperplanes and connects the centers of hyperplanes together by a curve. Second, it forces the molecule to sample on the curve and the hyperplanes in the simulation and calculates the free energy data on them. Finally, all the free energy data are combined together to form the local free energy landscape. Without consideration of the area outside the restricted region, this free energy calculation can be more efficient. By this method, one can further optimize the path quickly in the collective variable space.

Sampling errors in free energy simulations of small molecules in lipid bilayers

Free energy simulations are a powerful tool for evaluating the interactions of molecular solutes with lipid bilayers as mimetics of cellular membranes. However, these simulations are frequently hindered by systematic sampling errors. This review highlights recent progress in computing free energy profiles for inserting molecular solutes into lipid bilayers. Particular emphasis is placed on a systematic analysis of the free energy profiles, identifying the sources of sampling errors that reduce computational efficiency, and highlighting methodological advances that may alleviate sampling deficiencies. This article is part of a Special Issue entitled: Biosimulations edited by Ilpo Vattulainen and Tomasz Róg.