Conventional torsion angle potentials used in molecular dynamics (MD) have a singularity problem when three bonded particles are collinearly aligned. This problem is often encountered in ...coarse-grained (CG) simulations. Here, we propose a new form of the torsion angle potential, which introduces an angle-dependent modulating function. By carefully tuning the parameters for this modulating function, our method can eliminate the problematic angle-dependent singularity while being combined with existing models. As an example, we optimized the modulating function of the torsion angle potential for popular CG models of biomolecules based on the statistics over experimental structures deposited in the Protein Data Bank. By applying our method to designed and natural biomolecules, we show that the new torsion angle potential is able to eliminate the singularity problem while maintaining the structural features in the original models. Furthermore, by comparing our design with previous methods, we found that our new potential has advantages in computational efficiency and numerical stability. We strongly recommend the usage of our new potential in the CG simulations of flexible molecules.
Quantum mechanics/molecular mechanics (QM/MM) calculations are applied for anharmonic vibrational analyses of biomolecules and solvated molecules. The QM/MM method is implemented into a molecular ...dynamics (MD) program, GENESIS, by interfacing with external electronic structure programs. Following the geometry optimization and the harmonic normal-mode analysis based on a partial Hessian, the anharmonic potential energy surface (PES) is generated from QM/MM energies and gradients calculated at grid points. The PES is used for vibrational self-consistent field (VSCF) and post-VSCF calculations to compute the vibrational spectrum. The method is first applied to a phosphate ion in solution. With both the ion and neighboring water molecules taken as a QM region, IR spectra of representative hydration structures are calculated by the second-order vibrational quasi-degenerate perturbation theory (VQDPT2) at the level of B3LYP/cc-pVTZ and TIP3P force field. A weight-average of IR spectra over the structures reproduces the experimental spectrum with a mean absolute deviation of 16 cm
. Then, the method is applied to an enzyme, P450 nitric oxide reductase (P450nor), with the NO molecule bound to a ferric (Fe
) heme. Starting from snapshot structures obtained from MD simulations of P450nor in solution, QM/MM calculations have been carried out at the level of B3LYP-D3/def2-SVP(D). The spin state of Fe
(NO) is likely a closed-shell singlet state based on a ratio of N-O and Fe-NO stretching frequencies (ν
and ν
) calculated for closed- and open-shell singlet states. The calculated ν
and ν
overestimate the experimental ones by 120 and 75 cm
, respectively. The electronic structure and solvation of Fe
(NO) affect the structure around the heme of P450nor leading to an increase in ν
and ν
.
Sarcoplasmic reticulum (SR) Ca2+-ATPase transports two Ca2+ ions from the cytoplasm to the SR lumen against a large concentration gradient. X-ray crystallography has revealed the atomic structures of ...the protein before and after the dissociation of Ca2+, while biochemical studies have suggested the existence of intermediate states in the transition between E1P⋅ADP⋅2Ca2+ and E2P. Here, we explore the pathway and free energy profile of the transition using atomistic molecular dynamics simulations with the mean-force string method and umbrella sampling. The simulations suggest that a series of structural changes accompany the ordered dissociation of ADP, the A-domain rotation, and the rearrangement of the transmembrane (TM) helices. The luminal gate then opens to release Ca2+ ions toward the SR lumen. Intermediate structures on the pathway are stabilized by transient sidechain interactions between the A- and P-domains. Lipid molecules between TM helices play a key role in the stabilization. Free energy profiles of the transition assuming different protonation states suggest rapid exchanges between Ca2+ ions and protons when the Ca2+ ions are released toward the SR lumen.
Residue-level coarse-grained (CG) models have become one of the most popular tools in biomolecular simulations in the trade-off between modeling accuracy and computational efficiency. To investigate ...large-scale biological phenomena in molecular dynamics (MD) simulations with CG models, unified treatments of proteins and nucleic acids, as well as efficient parallel computations, are indispensable. In the GENESIS MD software, we implement several residue-level CG models, covering structure-based and context-based potentials for both well-folded biomolecules and intrinsically disordered regions. An amino acid residue in protein is represented as a single CG particle centered at the Cα atom position, while a nucleotide in RNA or DNA is modeled with three beads. Then, a single CG particle represents around ten heavy atoms in both proteins and nucleic acids. The input data in CG MD simulations are treated as GROMACS-style input files generated from a newly developed toolbox, GENESIS-CG-tool. To optimize the performance in CG MD simulations, we utilize multiple neighbor lists, each of which is attached to a different nonbonded interaction potential in the cell-linked list method. We found that random number generations for Gaussian distributions in the Langevin thermostat are one of the bottlenecks in CG MD simulations. Therefore, we parallelize the computations with message-passing-interface (MPI) to improve the performance on PC clusters or supercomputers. We simulate Herpes simplex virus (HSV) type 2 B-capsid and chromatin models containing more than 1,000 nucleosomes in GENESIS as examples of large-scale biomolecular simulations with residue-level CG models. This framework extends accessible spatial and temporal scales by multi-scale simulations to study biologically relevant phenomena, such as genome-scale chromatin folding or phase-separated membrane-less condensations.
Spike (S) protein is the primary antigenic target for neutralization and vaccine development for the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). It decorates the virus surface and ...undergoes large motions of its receptor binding domains (RBDs) to enter the host cell. Here, we observe Down, one-Up, one-Open, and two-Up-like structures in enhanced molecular dynamics simulations, and characterize the transition pathways via inter-domain interactions. Transient salt-bridges between RBD
and RBD
and the interaction with glycan at N343
support RBD
motions from Down to one-Up. Reduced interactions between RBD
and RBD
in one-Up induce RBD
motions toward two-Up. The simulations overall agree with cryo-electron microscopy structure distributions and FRET experiments and provide hidden functional structures, namely, intermediates along Down-to-one-Up transition with druggable cryptic pockets as well as one-Open with a maximum exposed RBD. The inherent flexibility of S-protein thus provides essential information for antiviral drug rational design or vaccine development.
The growing interest in the complexity of biological interactions is continuously driving the need to increase system size in biophysical simulations, requiring not only powerful and advanced ...hardware but adaptable software that can accommodate a large number of atoms interacting through complex forcefields. To address this, we developed and implemented strategies in the GENESIS molecular dynamics package designed for large numbers of processors. Long‐range electrostatic interactions were parallelized by minimizing the number of processes involved in communication. A novel algorithm was implemented for nonbonded interactions to increase single instruction multiple data (SIMD) performance, reducing memory usage for ultra large systems. Memory usage for neighbor searches in real‐space nonbonded interactions was reduced by approximately 80%, leading to significant speedup. Using experimental data describing physical 3D chromatin interactions, we constructed the first atomistic model of an entire gene locus (GATA4). Taken together, these developments enabled the first billion‐atom simulation of an intact biomolecular complex, achieving scaling to 65,000 processes (130,000 processor cores) with 1 ns/day performance. Published 2019. This article is a U.S. Government work and is in the public domain in the USA.
The first atomistic view of an intact gene (GATA4), based on experimental HiC data, consisting of 83 kilobases of double‐stranded helical DNA wrapped around 427 nucleosomes. Including water molecules and ions, the simulation box contains over 1 billion atoms.
The graphics processing unit (GPU) has become a popular computational platform for molecular dynamics (MD) simulations of biomolecules. A significant speedup in the simulations of small- or ...medium-size systems using only a few computer nodes with a single or multiple GPUs has been reported. Because of GPU memory limitation and slow communication between GPUs on different computer nodes, it is not straightforward to accelerate MD simulations of large biological systems that contain a few million or more atoms on massively parallel supercomputers with GPUs. In this study, we develop a new scheme in our MD software, GENESIS, to reduce the total computational time on such computers. Computationally intensive real-space nonbonded interactions are computed mainly on GPUs in the scheme, while less intensive bonded interactions and communication-intensive reciprocal-space interactions are performed on CPUs. On the basis of the midpoint cell method as a domain decomposition scheme, we invent the single particle interaction list for reducing the GPU memory usage. Since total computational time is limited by the reciprocal-space computation, we utilize the RESPA multiple time-step integration and reduce the CPU resting time by assigning a subset of nonbonded interactions on CPUs as well as on GPUs when the reciprocal-space computation is skipped. We validated our GPU implementations in GENESIS on BPTI and a membrane protein, porin, by MD simulations and an alanine-tripeptide by REMD simulations. Benchmark calculations on TSUBAME supercomputer showed that an MD simulation of a million atoms system was scalable up to 256 computer nodes with GPUs.