Intermediate scattering function for polymer molecules: An approach based on relaxation mode analysis.

The theory of polymer dynamics describes the intermediate scattering function for a polymer molecule in terms of relaxation modes defined by normal coordinates for the corresponding coarse-grained model. However, due to the difficulty of defining the normal coordinates for arbitrary polymer molecules, it is generally challenging to express the intermediate scattering function for a polymer molecule in terms of relaxation modes. To overcome this challenge, we propose a general method to calculate the intermediate scattering function for a polymer molecule on the basis of a relaxation mode analysis approach [Takano and Miyashita, J. Phys. Soc. Jpn. 64, 3688 (1995)]. In the proposed method, relaxation modes defined by eigenfunctions in a Markov process are evaluated on the basis of the simulation results for a polymer molecule and used to calculate the intermediate scattering function for that molecule. To demonstrate the effectiveness of the present method, we simulate the dynamics of a linear polymer molecule in a dilute solution and apply it to the calculation of the intermediate scattering function for the polymer molecule. The evaluation results regarding the relaxation modes reasonably describe the intermediate scattering function on the length scale of the radius of gyration of the polymer molecule. Accordingly, we examine the contributions of the pure relaxation and oscillatory relaxation processes to the entire intermediate scattering function.

Structural and Computational Insights into Dynamics and Intermediate States of Orexin 2 Receptor Signaling.

Orexin 2 receptor (OX2R) is a G protein-coupled receptor (GPCR) whose activation is crucial to regulation of the sleep-wake cycle. Recently, inactive and active state structures were determined from X-ray crystallography and cryo-electron microscopy single particle analysis, and the activation mechanisms have been discussed based on these static data. GPCRs have multiscale intermediate states during activation, and insights into these dynamics and intermediate states may aid the precise control of intracellular signaling by ligands in drug discovery. Molecular dynamics (MD) simulations are used to investigate dynamics induced in response to thermal perturbations, such as structural fluctuations of main and side chains. In this study, we proposed collective motions of the TM domain during activation by performing 30 independent microsecond-scale MD simulations for various OX2R systems and applying relaxation mode analysis. The analysis results suggested that TM3 had a vertical structural movement relative to the membrane surface during activation. In addition, we extracted three characteristic amino acid residues on TM3, i.e., Q1343.32, V1423.40, and R1523.50, which exhibited large conformational fluctuations. We quantitatively evaluated the changes in their equilibrium during activation in relation to the movement of TM3. We also discuss the regulation of ligand binding recognition and intracellular signal selectivity by changes in the equilibrium of Q1343.32 and R1523.50, respectively, according to MD simulations and GPCR database. Additionally, the OX2R-Gi signaling complex is stabilized in the conformation resembling a non-canonical (NC) state, which was previously proposed as an intermediate state during activation of neurotensin 1 receptor. Insights into the dynamics and intermediate states during activation gained from this study may be useful for developing biased agonists for OX2R.

Effect of Main and Side Chains on the Folding Mechanism of the Trp-Cage Miniprotein.

Proteins that do not fold into their functional native state have been linked to diseases. In this study, the influence of the main and side chains of individual amino acids on the folding of the tryptophan cage (Trp-cage), a designed 20-residue miniprotein, was analyzed. For this purpose, we calculated the solvation free energy (SFE) contributions of individual atoms by using the 3D-reference interaction site model with the atomic decomposition method. The mechanism by which the Trp-cage is stabilized during the folding process was examined by calculating the total energy, which is the sum of the conformational energy and SFE. The folding process of the Trp-cage resulted in a stable native state, with a total energy that was 62.4 kcal/mol lower than that of the unfolded state. The solvation entropy, which is considered to be responsible for the hydrophobic effect, contributed 31.3 kcal/mol to structural stabilization. In other words, the contribution of the solvation entropy accounted for approximately half of the total contribution to Trp-cage folding. The hydrophobic core centered on Trp6 contributed 15.6 kcal/mol to the total energy, whereas the solvation entropy contribution was 6.3 kcal/mol. The salt bridge formed by the hydrophilic side chains of Asp9 and Arg16 contributed 10.9 and 5.0 kcal/mol, respectively. This indicates that not only the hydrophobic core but also the salt bridge of the hydrophilic side chains gain solvation entropy and contribute to stabilizing the native structure of the Trp-cage.

Analysis of Protein Folding Simulation with Moving Root Mean Square Deviation.

We apply moving root-mean-square deviation (mRMSD), which does not require a reference structure, as a method for analyzing protein dynamics. This method can be used to calculate the root-mean-square deviation (RMSD) of structure between two specified time points and to analyze protein dynamics behavior through time series analysis. We applied this method to the Trp-cage trajectory calculated by the Anton supercomputer and found that it shows regions of stable states as well as the conventional RMSD. In addition, we extracted a characteristic structure in which the side chains of Asp1 and Arg16 form hydrogen bonds near the most stable structure of the Trp-cage. We also determined that ≥20 ns is an appropriate time interval to investigate protein dynamics using mRMSD. Applying this method to NuG2 protein, we found that mRMSD can be used to detect regions of metastable states in addition to the stable state. This method can be applied to molecular dynamics simulations of proteins whose stable structures are unknown.

Non-Markov-Type Analysis and Diffusion Map Analysis for Molecular Dynamics Trajectory of Chignolin at a High Temperature.

We apply the non-Markov-type analysis of state-to-state transitions to nearly microsecond molecular dynamics (MD) simulation data at a folding temperature of a small artificial protein, chignolin, and we found that the time scales obtained are consistent with our previous result using the weighted ensemble simulations, which is a general path-sampling method to extract the kinetic properties of molecules. Previously, we also applied diffusion map (DM) analysis, which is one of a manifold of learning techniques, to the same trajectory of chignolin in order to cluster the conformational states and found that DM and relaxation mode analysis give similar results for the eigenvectors. In this paper, we divide the same trajectory into shorter pieces and further apply DM to such short-length trajectories to investigate how the obtained eigenvectors are useful to characterize the conformational change of chignolin.

Computational Analysis of the SARS-CoV-2 RBD-ACE2-Binding Process Based on MD and the 3D-RISM Theory.

The binding process of angiotensin-converting enzyme 2 (ACE2) to the receptor-binding domain (RBD) of the severe acute respiratory syndrome-like coronavirus 2 spike protein was investigated using molecular dynamics simulation and the three-dimensional reference interaction-site model theory. The results suggested that the protein-binding process consists of a protein-protein approaching step, followed by a local structural rearrangement step. In the approaching step, the interprotein interaction energy decreased as the proteins approached each other, whereas the solvation free energy increased. As the proteins approached, the glycan of ACE2 first established a hydrogen bond with the RBD. Thereafter, the number of interprotein hydrogen bonds increased rapidly. The solvation free energy increased because of the desolvation of the protein as it approached its partner. The spatial distribution function of the solvent revealed the presence of hydrogen bonds bridged by water molecules on the RBD-ACE2 interface. Finally, principal component analysis revealed that ACE2 showed a pronounced conformational change, whereas there was no significant change in RBD.

Characteristic structural difference between inactive and active states of orexin 2 receptor determined using molecular dynamics simulations.

To investigate the dynamics of the orexin 2 receptor, which is a class A G protein-coupled receptor, we recently performed several microsecond-scale molecular dynamics simulations of the wild-type protein, of a mutant that stabilizes the inactive state, and of constitutively active mutants of the class A G protein-coupled receptors. Herein, we review the results of these molecular dynamics simulations of the orexin 2 receptor. In these simulations, characteristic conformational changes were observed in the V3096.40Y mutant. The conformational changes were related to the outward movement of the transmembrane helix 6 and the inward movement of the transmembrane helix 7, which are common structural changes in the activation of G protein-coupled receptors. The index for the quantitative evaluation of the active and inactive states of class A G protein-coupled receptors and the mechanism of the inward movement of the transmembrane helix 7 were examined. In this review, we also discuss the activation mechanism by comparing the structures obtained from the molecular dynamics simulations with the structure of the active state of the orexin 2 receptor clarified by cryo-electron microscopy in the recent years.

Molecular Dynamics Simulations for the Determination of the Characteristic Structural Differences between Inactive and Active States of Wild Type and Mutants of the Orexin2 Receptor.

The orexin2 receptor (OX2R), which is classified as a class A G protein-coupled receptor (GPCR), is the target of our study. We performed over 20 several-microsecond-scale molecular dynamics simulations of the wild type and mutants of OX2R to extract the characteristics of the structural changes taking place in the active state. We introduced mutations that exhibited the stable inactive state and the constitutively active state in class A GPCRs. In these simulations, significant characteristic structural changes were observed in the V3096.40Y mutant, which corresponded to a constitutively active mutant. These conformational changes include the outward movement of the transmembrane helix 6 (TM6) and the inward movement of TM7, which are common structural changes in the activation of GPCRs. In addition, we extracted a suitable index for the quantitative evaluation of the active and inactive states of GPCRs, namely, the inter-atomic distance of Cα atoms between x(3.46) and Y(7.53). The structures of the inactive and active states solved by X-ray crystallography and cryo-electron microscopy can be classified using the inter-atomic distance. Furthermore, we clarified that the inward movement of TM7 requires the swapping of M3056.36 on TM6 and L3677.56 on TM7. Finally, we discussed the structural advantages of TM7 inward movement for GPCR activation.

Universal Relation between Instantaneous Diffusivity and Radius of Gyration of Proteins in Aqueous Solution.

Protein conformational fluctuations are highly complex and exhibit long-term correlations. Here, molecular dynamics simulations of small proteins demonstrate that these conformational fluctuations directly affect the protein's instantaneous diffusivity D_{I}. We find that the radius of gyration R_{g} of the proteins exhibits 1/f fluctuations that are synchronous with the fluctuations of D_{I}. Our analysis demonstrates the validity of the local Stokes-Einstein-type relation D_{I}∝1/(R_{g}+R_{0}), where R_{0}∼0.3 nm is assumed to be a hydration layer around the protein. From the analysis of different protein types with both strong and weak conformational fluctuations, the validity of the Stokes-Einstein-type relation appears to be a general property.

Mutation-induced change in chignolin stability from π-turn to α-turn.

Chignolin, which consists of 10 amino acids, adopts two stable states in simulations at room temperature at 1 atm: the native and misfolded states. The sequence of chignolin is optimized to form a stable π-turn and thus the native state has a π-turn from Asp3 to Thr8. On the other hand, the misfolded state adopts an α-turn from Asp3 to Gly7. We previously investigated the differences in the stability mechanism of the two states using computational techniques. Our previous detailed energy analysis implied that the native state was stabilized by hydrogen bonding between the side chain atoms of Thr6 and Thr8, and Thr8 was not involved in stabilization of the misfolded state. Thus, we predicted that mutation of Thr8 to a neutral amino acid could stabilize the misfolded structure over the native structure. In the present work, we performed 4 µs molecular dynamics simulations for 19 mutants of the 8th residue. Among them, the T8I, T8F, T8P, T8N, and T8Y mutants, in which the 8th residue was changed to a neutral residue, formed only the misfolded structure at room temperature. Even at high temperature, for the T8P mutant, the native structure was not observed, as the T8P mutant cannot form the native structure because of steric hindrance caused by the distinctive cyclic structure of proline. Interestingly, the T8P mutant at high temperature has trans and cis conformations in the Gly7-Pro8 sequence, with the trans conformation corresponding to the misfolded state. NMR analysis of the T8P mutant supported our results.

Folding pathways of NuG2-a designed mutant of protein G-using relaxation mode analysis.

Dynamic analysis methods are important for analyzing long simulations such as folding simulations. Relaxation mode analysis, which approximately extracts slow modes and rates, has been applied in molecular dynamics (MD) simulations of protein systems. Previously, we showed that slow modes are suitable for analyzing simulations in which large conformational changes occur. Here, we applied relaxation mode analysis to folding simulations of a designed mutant of protein G, NuG2, to investigate its folding pathways. The folding simulations of NuG2 were previously performed for this mutant with Anton. In the present study, the free energy surfaces were calculated by projecting the coordinates on the axis of the slow relaxation modes obtained from relaxation mode analysis. We classified various characteristic states such as native, nativelike, intermediate, and random states and clarified two main folding pathways. In the early folding process, the first and second ß strands formed an N-terminal ß-sheet. After the early folding process, the fourth ß strand formed along the first ß strand in the same or opposite direction as the native structure; two characteristic intermediate states were identified. Finally, the intermediate structures folded to the native structure in the folding process. Relaxation mode analysis can be applied to folding simulations of complex proteins to investigate their folding processes.

Identification of slow relaxation modes in a protein trimer via positive definite relaxation mode analysis.

Recently, dynamic analysis methods in signal processing have been applied to the analysis of molecular dynamics (MD) trajectories of biopolymers. In the context of a relaxation mode analysis (RMA) method, based on statistical physics, it is explained why the signal-processing methods work well for the simulation trajectories of biopolymers. A distinctive difference between the RMA method and the signal-processing methods is the introduction of an additional parameter, called an evolution time parameter. This parameter enables us to better estimate the relaxation modes and rates, although it increases computational difficulty. In this paper, we propose a simple and effective extension of the RMA method, which is referred to as the positive definite RMA method, to introduce the evolution time parameter robustly. In this method, an eigenvalue problem for the time correlation matrix of physical quantities relevant to slow relaxation in a system is first solved to find the subspace in which the matrix is numerically positive definite. Then, we implement the RMA method in the subspace. We apply the method to the analysis of a 3-µs MD trajectory of a heterotrimer of an erythropoietin protein and two of its receptor proteins, and we demonstrate the effectiveness of the method.

Analysis of molecular dynamics simulations of 10-residue peptide, chignolin, using statistical mechanics: Relaxation mode analysis and three-dimensional reference interaction site model theory.

Molecular dynamics simulation is a fruitful tool for investigating the structural stability, dynamics, and functions of biopolymers at an atomic level. In recent years, simulations can be performed on time scales of the order of milliseconds using special purpose systems. Since the most stable structure, as well as meta-stable structures and intermediate structures, is included in trajectories in long simulations, it is necessary to develop analysis methods for extracting them from trajectories of simulations. For these structures, methods for evaluating the stabilities, including the solvent effect, are also needed. We have developed relaxation mode analysis to investigate dynamics and kinetics of simulations based on statistical mechanics. We have also applied the three-dimensional reference interaction site model theory to investigate stabilities with solvent effects. In this paper, we review the results for designing amino-acid substitution of the 10-residue peptide, chignolin, to stabilize the misfolded structure using these developed analysis methods.

Conformational change of a biomolecule studied by the weighted ensemble method: Use of the diffusion map method to extract reaction coordinates.

We simulate the nonequilibrium ensemble dynamics of a biomolecule using the weighted ensemble method, which was introduced in molecular dynamics simulations by Huber and Kim and further developed by Zuckerman and co-workers. As the order parameters to characterize its conformational change, we here use the coordinates derived from the diffusion map (DM) method, one of the manifold learning techniques. As a concrete example, we study the kinetic properties of a small peptide, chignolin in explicit water, and calculate the conformational change between the folded and misfolded states in a nonequilibrium way. We find that the transition time scales thus obtained are comparable to those using previously employed hydrogen-bond distances as the order parameters. Since the DM method only uses the 3D Cartesian coordinates of a peptide, this shows that the DM method can extract the important distance information of the peptide without relying on chemical intuition. The time scales are compared well with the previous results using different techniques, non-Markovian analysis and core-set milestoning for a single long trajectory. We also find that the most significant DM coordinate turns out to extract a dihedral angle of glycine, and the previously studied relaxation modes are well correlated with the most significant DM coordinates.

Relaxation mode analysis for molecular dynamics simulations of proteins.

Molecular dynamics simulation is a powerful method for investigating the structural stability, dynamics, and function of biopolymers at the atomic level. In recent years, it has become possible to perform simulations on time scales of the order of milliseconds using special hardware. However, it is necessary to derive the important factors contributing to structural change or function from the complicated movements of biopolymers obtained from long simulations. Although some analysis methods for protein systems have been developed using increasing simulation times, many of these methods are static in nature (i.e., no information on time). In recent years, dynamic analysis methods have been developed, such as the Markov state model and relaxation mode analysis (RMA), which was introduced based on spin and homopolymer systems. The RMA method approximately extracts slow relaxation modes and rates from trajectories and decomposes the structural fluctuations into slow relaxation modes, which characterize the slow relaxation dynamics of the system. Recently, this method has been applied to biomolecular systems. In this article, we review RMA and its improved versions for protein systems.

Analysis of Structural Stability of Chignolin.

We discuss the stability of an entire protein and the influence of main chains and side chains of individual amino acids to investigate the protein-folding mechanism. For this purpose, we calculated the solvation free-energy contribution of individual atoms using the three-dimensional reference interaction site model with the atomic decomposition method. We generated structures of chignolin miniprotein by a molecular dynamics simulation and classified them into six types: native 1, native 2, misfolded 1, misfolded 2, intermediate, and unfolded states. The total energies of the native (-171.1 kcal/mol) and misfolded (-171.2 kcal/mol) states were almost the same and lower than those of the intermediate (-158.5 kcal/mol) and unfolded (-148.1 kcal/mol) states; however, their components were different. In the native state, the side-chain interaction between Thr6 and Thr8 is important for the formation of π-turn. On the other hand, the hydrogen bonds between the atoms of the main chains in the misfolded state become stronger than those in the intermediate state.

Application of reference-modified density functional theory: Temperature and pressure dependences of solvation free energy.

Recently, we proposed a reference-modified density functional theory (RMDFT) to calculate solvation free energy (SFE), in which a hard-sphere fluid was introduced as the reference system instead of an ideal molecular gas. Through the RMDFT, using an optimal diameter for the hard-sphere reference system, the values of the SFE calculated at room temperature and normal pressure were in good agreement with those for more than 500 small organic molecules in water as determined by experiments. In this study, we present an application of the RMDFT for calculating the temperature and pressure dependences of the SFE for solute molecules in water. We demonstrate that the RMDFT has high predictive ability for the temperature and pressure dependences of the SFE for small solute molecules in water when the optimal reference hard-sphere diameter determined for each thermodynamic condition is used. We also apply the RMDFT to investigate the temperature and pressure dependences of the thermodynamic stability of an artificial small protein, chignolin, and discuss the mechanism of high-temperature and high-pressure unfolding of the protein. © 2017 Wiley Periodicals, Inc.

Stability of Unfolded and Folded Protein Structures Using a 3D-RISM with the RMDFT.

Protein stability is determined by the characteristics of the protein itself as well as the surrounding solvent. Herein, we discuss the stability of the folded and unfolded structures of proteins obtained from Anton's long simulations (Lindorff-Larsen, K.; Piana, S.; Dror, R. O.; Shaw, D.E. Science, 2011, 334, 517-520). Specifically, the stabilities of CLN025, the WW domain variant GTT, the triple mutant of the redesigned protein G variant NuG2, and the de novo-designed three-helix bundle protein are investigated. The solvation free energy of the structures is calculated using the three-dimensional reference interaction site model with the reference-modified density functional theory. The total energy is given by the sum of the conformational energy and the solvation free energy, and their balance results in the stabilization of protein structure, as demonstrated by the correspondence between structures with the lowest total energy of all proteins to their native structures. Overall, these findings indicate that the total energy function is appropriate for evaluating the stability of protein folding systems. Moreover, decomposing the energy terms reveals that proteins achieve their stabilities from the balance between the conformational energy and the solvation free energy. In particular, the solvation entropy is the main contributor to the process of folding from more extended structures to compact structures. The native structure is more stable than the compact structure owing to competition between intramolecular and intermolecular interactions.

A reference-modified density functional theory: An application to solvation free-energy calculations for a Lennard-Jones solution.

In the conventional classical density functional theory (DFT) for simple fluids, an ideal gas is usually chosen as the reference system because there is a one-to-one correspondence between the external field and the density distribution function, and the exact intrinsic free-energy functional is available for the ideal gas. In this case, the second-order density functional Taylor series expansion of the excess intrinsic free-energy functional provides the hypernetted-chain (HNC) approximation. Recently, it has been shown that the HNC approximation significantly overestimates the solvation free energy (SFE) for an infinitely dilute Lennard-Jones (LJ) solution, especially when the solute particles are several times larger than the solvent particles [T. Miyata and J. Thapa, Chem. Phys. Lett. 604, 122 (2014)]. In the present study, we propose a reference-modified density functional theory as a systematic approach to improve the SFE functional as well as the pair distribution functions. The second-order density functional Taylor series expansion for the excess part of the intrinsic free-energy functional in which a hard-sphere fluid is introduced as the reference system instead of an ideal gas is applied to the LJ pure and infinitely dilute solution systems and is proved to remarkably improve the drawbacks of the HNC approximation. Furthermore, the third-order density functional expansion approximation in which a factorization approximation is applied to the triplet direct correlation function is examined for the LJ systems. We also show that the third-order contribution can yield further refinements for both the pair distribution function and the excess chemical potential for the pure LJ liquids.

Relaxation mode analysis and Markov state relaxation mode analysis for chignolin in aqueous solution near a transition temperature.

It is important to extract reaction coordinates or order parameters from protein simulations in order to investigate the local minimum-energy states and the transitions between them. The most popular method to obtain such data is principal component analysis, which extracts modes of large conformational fluctuations around an average structure. We recently applied relaxation mode analysis for protein systems, which approximately estimates the slow relaxation modes and times from a simulation and enables investigations of the dynamic properties underlying the structural fluctuations of proteins. In this study, we apply this relaxation mode analysis to extract reaction coordinates for a system in which there are large conformational changes such as those commonly observed in protein folding/unfolding. We performed a 750-ns simulation of chignolin protein near its folding transition temperature and observed many transitions between the most stable, misfolded, intermediate, and unfolded states. We then applied principal component analysis and relaxation mode analysis to the system. In the relaxation mode analysis, we could automatically extract good reaction coordinates. The free-energy surfaces provide a clearer understanding of the transitions not only between local minimum-energy states but also between the folded and unfolded states, even though the simulation involved large conformational changes. Moreover, we propose a new analysis method called Markov state relaxation mode analysis. We applied the new method to states with slow relaxation, which are defined by the free-energy surface obtained in the relaxation mode analysis. Finally, the relaxation times of the states obtained with a simple Markov state model and the proposed Markov state relaxation mode analysis are compared and discussed.