Hisashi Okumura

Temperature and pressure denaturation of chignolin: Folding and unfolding simulation by multibaric-multithermal molecular dynamics method

A multibaric‐multithermal molecular dynamics (MD) simulation of a 10‐residue protein, chignolin, was performed. All‐atom model with the Amber parm99SB force field was used for the protein and the TIP3P model was used for the explicit water molecules. This MD simulation covered wide ranges of temperature between 260 and 560 K and pressure between 0.1 and 600 MPa and sampled many conformations without getting trapped in local‐minimum free‐energy states. Folding events to the native β‐hairpin structure occurred five times and unfolding events were observed four times. As the temperature and/or pressure increases, fraction of folded chignolin decreases. The partial molar enthalpy change ΔH and partial molar volume change ΔV of unfolding were calculated as ΔH = 24.1 ± 4.9 kJ/mol and ΔV = −5.6 ± 1.5 cm3/mol, respectively. These values agree well with recent experimental results. Illustrating typical local‐minimum free‐energy conformations, folding and unfolding pathways were revealed. When chignolin unfolds from the β‐hairpin structure, only the C terminus or both C and N termini open first. It may undergo an α‐helix or 310‐helix structure and finally unfolds to the extended structure. Difference of the mechanism between temperature denaturation and pressure denaturation is also discussed. Temperature denaturation is caused by making the protein transferred to a higher entropy state and making it move around more with larger space. The reason for pressure denaturation is that water molecules approach the hydrophobic residues, which are not well hydrated at the folded state, and some hydrophobic contacts are broken. Proteins 2012;.

Replica-Permutation Method with the Suwa–Todo Algorithm beyond the Replica-Exchange Method

We propose a new method for molecular dynamics and Monte Carlo simulations, which is referred to as the replica-permutation method (RPM), to realize more efficient sampling than the replica-exchange method (REM). In RPM, not only exchanges between two replicas but also permutations among more than two replicas are performed. Furthermore, instead of the Metropolis algorithm, the Suwa-Todo algorithm is employed for replica-permutation trials to minimize its rejection ratio. We applied RPM to particles in a double-well potential energy, Met-enkephalin in a vacuum, and a C-peptide analog of ribonuclease A in explicit water. For comparison purposes, replica-exchange molecular dynamics simulations were also performed. As a result, RPM sampled not only the temperature space but also the conformational space more efficiently than REM for all systems. From our simulations of C-peptide, we obtained the α-helix structure with salt bridges between Gly2 and Arg10, which is known in experiments. Calculating its free-energy landscape, the folding pathway was revealed from an extended structure to the α-helix structure with the salt bridges. We found that the folding pathway consists of the two steps: The first step is the salt-bridge formation step, and the second step is the α-helix formation step.

Amyloid Fibril Disruption by Ultrasonic Cavitation: Nonequilibrium Molecular Dynamics Simulations

We describe the disruption of amyloid fibrils of Alzheimers amyloid-β peptides by ultrasonic cavitation. For this purpose, we performed nonequilibrium all-atom molecular dynamics simulations with sinusoidal pressure and visualized the process with movies. When the pressure is negative, a bubble is formed, usually at hydrophobic residues in the transmembrane region. Most β-strands maintain their secondary structures in the bubble. When the pressure becomes positive, the bubble collapses, and water molecules crash against the hydrophilic residues in the nontransmembrane region to disrupt the amyloid. Shorter amyloids require longer sonication times for disruption because they do not have enough hydrophobic residues to serve as a nucleus to form a bubble. These results agree with experiments in which monodispersed amyloid fibrils were obtained by ultrasonication.

Pressure-Induced Helical Structure of a Peptide Studied by Simulated Tempering Molecular Dynamics Simulations.

It is known experimentally that an AK16 peptide forms more α-helix structures with increasing pressure while proteins unfold in general. In order to understand this abnormality, molecular dynamics (MD) simulations with the simulated tempering method for the isobaric-isothermal ensemble were performed in a wide pressure range from 1.0 × 10(-4) GPa to 1.4 GPa. From the results of the simulations, it is found that the fraction of the folded state decreases once and increases after that with increasing pressure. The partial molar volume change from the folded state to unfolded state increases monotonically from a negative value to a positive value with pressure. The behavior under high pressure conditions is consistent with the experimental results. The radius of gyration of highly helical structures decreases with increasing pressure, which indicates that the helix structure shrinks with pressure. This is the reason why the fraction of the folded state increases as pressure increases.

Dimerization process of amyloid-β(29-42) studied by the Hamiltonian replica-permutation molecular dynamics simulations.

The amyloid-β peptides form amyloid fibrils which are associated with Alzheimers disease. Amyloid-β(29-42) is its C-terminal fragment and a critical determinant of the amyloid formation rate. This fragment forms the amyloid fibril by itself. However, the fragment conformation in the fibril has yet to be determined. The oligomerization process including the dimerization process is also still unknown. The dimerization process corresponds to an early process of the amyloidogenesis. In order to investigate the dimerization process and conformations, we applied the Hamiltonian replica-permutation method, which is a better alternative to the Hamiltonian replica-exchange method, to two amyloid-β(29-42) molecules in explicit water solvent. At the first step of the dimerization process, two amyloid-β(29-42) molecules came close to each other and had intermolecular side chain contacts. When two molecules had the intermolecular side chain contacts, the amyloid-β(29-42) tended to have intramolecular secondary structures, especially β-hairpin structures. The two molecules had intermolecular β-bridge structures by coming much closer at the second step of the dimerization process. Formation of these intermolecular β-bridge structures was induced by the β-hairpin structures. The intermolecular β-sheet structures elongated at the final step. Structures of the amyloid-β(29-42) in the monomer and dimer states are also shown with the free-energy landscapes, which were obtained by performing efficient sampling in the conformational space in our simulations.

Coulomb replica-exchange method: Handling electrostatic attractive and repulsive forces for biomolecules

We propose a new type of the Hamiltonian replica‐exchange method (REM) for molecular dynamics (MD) and Monte Carlo simulations, which we refer to as the Coulomb REM (CREM). In this method, electrostatic charge parameters in the Coulomb interactions are exchanged among replicas while temperatures are exchanged in the usual REM. By varying the atom charges, the CREM overcomes free‐energy barriers and realizes more efficient sampling in the conformational space than the REM. Furthermore, this method requires only a smaller number of replicas because only the atom charges of solute molecules are used as exchanged parameters. We performed Coulomb replica‐exchange MD simulations of an alanine dipeptide in explicit water solvent and compared the results with those of the conventional canonical, replica exchange, and van der Waals REMs. Two force fields of AMBER parm99 and AMBER parm99SB were used. As a result, the CREM sampled all local‐minimum free‐energy states more frequently than the other methods for both force fields. Moreover, the Coulomb, van der Waals, and usual REMs were applied to a fragment of an amyloid‐β peptide (Aβ) in explicit water solvent to compare the sampling efficiency of these methods for a larger system. The CREM sampled structures of the Aβ fragment more efficiently than the other methods. We obtained β‐helix, α‐helix, 310‐helix, β‐hairpin, and β‐sheet structures as stable structures and deduced pathways of conformational transitions among these structures from a free‐energy landscape.

Decomposition-order effects of time integrator on ensemble averages for the Nosé-Hoover thermostat

Decomposition-order dependence of time development integrator on ensemble averages for the Nosé-Hoover dynamics is discussed. Six integrators were employed for comparison, which were extensions of the velocity-Verlet or position-Verlet algorithm. Molecular dynamics simulations by these integrators were performed for liquid-argon systems with several different time steps and system sizes. The obtained ensemble averages of temperature and potential energy were shifted from correct values depending on the integrators. These shifts increased in proportion to the square of the time step. Furthermore, the shifts could not be removed by increasing the number of argon atoms. We show the origin of these ensemble-average shifts analytically. Our discussion can be applied not only to the liquid-argon system but also to all MD simulations with the Nosé-Hoover thermostat. Our recommended integrators among the six integrators are presented to obtain correct ensemble averages.

Length Dependence of Polyglycine Conformations in Vacuum

We performed replica-exchange molecular dynamics simulations of polyglycines in vacuum to investigate their conformational difference due to different numbers of residues. We employed the polyglycines of which the numbers of residues are 1 (PG1), 5 (PG5), 10 (PG10), and 15 (PG15). We show and discuss the conformations of the polyglycine molecules, which have the lowest potential energies in local minimum states. The polyglycines PG5 and PG10 often have helical structures. The helical structures of the polyglycine PG10 are β-helix structures. The PG15 have complicated tertiary structures. The tertiary structures have two β-hairpins in the N- and C-terminal regions. A parallel β-sheet structure is also formed between the N-terminal side of the N-terminal β-hairpin and the C-terminal side of the C-terminal β-hairpin.

Manifold Correction Method for the Nose-Hoover and Nose-Poincare Molecular Dynamics Simulations

We introduce the manifold correction method to molecular dynamics (MD) simulations with the Nose–Hoover and Nose–Poincare thermostats. The manifold correction method was originally developed in astronomy, as an accurate numerical method for many body systems. Because the Nose–Hoover thermostat is not a symplectic algorithm, the quantity which is conserved analytically is not conserved but increases in actual MD simulations. Using the manifold correction method, this quantity is completely conserved, and it makes the MD simulation stable. Because the conservation of this quantity is required in the proof that the Nose–Hoover thermostat gives the canonical ensemble, the manifold correction method guarantees to provide the correct statistical ensemble. Although the time development of the Nose–Poincare thermostat is described as a symplectic algorithm, if the interatomic potential energy is truncated, the Nose–Poincare thermostat is no longer symplectic. In this case, the Hamiltonian increases, and temperatu...

Cutoff Effect in the Nosé–Poincaré and Nosé–Hoover Thermostats

We performed molecular dynamics (MD) simulations of a Lennard-Jones system and investigated the effect of potential cutoff in the Nose–Poincare and Nose–Hoover thermostats. The Nose–Poincare thermostat is the symplectic algorithm of the Nose thermostat, while the Nose–Hoover thermostat is not a symplectic algorithm. If the potential energy is twice or more differentiable, the Hamiltonian was conserved well in the Nose–Poincare thermostat. If the potential energy is once or less differentiable, however, the Hamiltonian was not conserved, but increased because the continuity of potential energy is required in a symplectic MD simulation. The increase in the Hamiltonian caused the increase in instantaneous temperature, and physical quantities cannot be obtained correctly. It is because the difference in the Hamiltonian effectively increases the set temperature in the equations of motion. On the other hand, the Hamiltonian was not conserved for any cutoff method in the Nose–Hoover thermostat because it is not ...