Parallel Tempering Simulations Research Articles

Comparison of On-the-Fly Probability Enhanced Sampling and Parallel Tempering Combined with Metadynamics for Atomistic Simulations of RNA Tetraloop Folding.

Atomistic simulations with reliable models are extremely useful in providing exquisitely detailed pictures of biomolecular phenomena that are not always accessible to experiments. One such biomolecular phenomenon is RNA folding, which often requires exhaustive simulations with combined advanced sampling techniques. In this work, we employed the multithermal-multiumbrella on-the-fly probability enhanced sampling (MM-OPES) technique and compared it against combined parallel tempering and metadynamics simulations. We found that MM-OPES simulations were successful in reproducing the free energy surfaces from combined parallel tempering and metadynamics simulations. Importantly, we also investigated a wide range of temperature sets (minimum and maximum temperatures) for MM-OPES simulations in order to identify some guidelines for deciding the temperature limits for an accurate and efficient exploration of the free energy landscapes. We found that most temperature sets yielded almost the same accuracy in reproducing the free energy surface at the ambient conditions as long as (i) the maximum temperature is reasonably high, (ii) the temperature at which we run the simulation is reasonably high (in our simulations, running temperature is defined as [minimum temperature + maximum temperature]/2), and (iii) the effective sample size at the temperature of interest is statistically reasonable. In terms of the computational cost, all MM-OPES simulations were nearly 4 times less costly than the combined parallel tempering and metadynamics simulations. We concluded that the demanding combined parallel tempering and metadynamics simulations can safely be replaced with approximately 4 times less costly MM-OPES simulations (with carefully selected temperature limits) to obtain the same information.

Structures of MERS-CoV macro domain in aqueous solution with dynamics: Impacts of parallel tempering simulation techniques and CHARMM36m and AMBER99SB force field parameters.

A novel virus, severe acute respiratory syndrome Coronavirus 2 (SARS‐CoV‐2), causing coronavirus disease 2019 (COVID‐19) worldwide appeared in 2019. Detailed scientific knowledge of the members of the Coronaviridae family, including the Middle East Respiratory Syndrome Coronavirus (MERS‐CoV) is currently lacking. Structural studies of the MERS‐CoV proteins in the current literature are extremely limited. We present here detailed characterization of the structural properties of MERS‐CoV macro domain in aqueous solution. Additionally, we studied the impacts of chosen force field parameters and parallel tempering simulation techniques on the predicted structural properties of MERS‐CoV macro domain in aqueous solution. For this purpose, we conducted extensive Hamiltonian‐replica exchange molecular dynamics simulations and Temperature‐replica exchange molecular dynamics simulations using the CHARMM36m and AMBER99SB parameters for the macro domain. This study shows that the predicted secondary structure properties including their propensities depend on the chosen simulation technique and force field parameter. We perform structural clustering based on the radius of gyration and end‐to‐end distance of MERS‐CoV macro domain in aqueous solution. We also report and analyze the residue‐level intrinsic disorder features, flexibility and secondary structure. Furthermore, we study the propensities of this macro domain for protein‐protein interactions and for the RNA and DNA binding. Overall, results are in agreement with available nuclear magnetic resonance spectroscopy findings and present more detailed insights into the structural properties of MERS CoV macro domain in aqueous solution. All in all, we present the structural properties of the aqueous MERS‐CoV macro domain using different parallel tempering simulation techniques, force field parameters and bioinformatics tools.

Critical properties of the frustrated Ising model on a honeycomb lattice: A Monte Carlo study

Critical and in the highly frustrated regime also dynamical properties of the J1−J2 Ising model with competing nearest-neighbor J1 and second-nearest-neighbor J2 interactions on a honeycomb lattice are investigated by standard Monte Carlo and parallel tempering simulations. The phase boundary is determined as a function of the coupling ratio for the phase transition between the paramagnetic and ferromagnetic states within R≡J2/|J1|∈[−1/4,0]. It is confirmed that at least for R≥−0.2 the transition remains second-order and complies with the standard Ising universality class. In the highly frustrated regime of R<−0.2 and low temperatures the system tends to freeze to metastable domain states, separated by large energy barriers, which show extremely sluggish dynamics. The resulting huge equilibration and autocorrelation times hinder the analysis of critical properties and thus the character of the transition in this region remains to be determined.

Understanding population annealing Monte Carlo simulations

Population annealing is a recent addition to the arsenal of the practitioner in computer simulations in statistical physics and it proves to deal well with systems with complex free-energy landscapes. Above all else, it promises to deliver unrivaled parallel scaling qualities, being suitable for parallel machines of the biggest caliber. Here we study population annealing using as the main example the two-dimensional Ising model, which allows for particularly clean comparisons due to the available exact results and the wealth of published simulational studies employing other approaches. We analyze in depth the accuracy and precision of the method, highlighting its relation to older techniques such as simulated annealing and thermodynamic integration. We introduce intrinsic approaches for the analysis of statistical and systematic errors and provide a detailed picture of the dependence of such errors on the simulation parameters. The results are benchmarked against canonical and parallel tempering simulations.

Improving Temperature Generator in Parallel Tempering Simulation in the NPT Condition.

Parallel tempering (PT) simulations provide a powerful computational tool to accelerate conformational searches of complex molecular systems. In this method, multiple replica systems assigned to different temperatures run in parallel and undergo tandem replica exchange events to enhance conformational mixing in temperature space. Efficient PT simulations require a uniform acceptance ratio across all the replicas. In the context of a previous energy distribution-based temperature generation (TG) scheme, we propose an improved TG protocol to maintain such a uniform exchange probability in general PT simulations.

Time resolved alkali silicate decondensation by sodium hydroxide solution

Silica is by far the chemical compound the most widespread and used around the world: as a raw product in the buildings and roads industry, as concrete, or as a processed product in the manufacture of glass, ceramics or zeolites. In alkali silicate solutions—often used to synthesize those materials—a complex interplay of decondensation and condensation processes leads to the restructuring of silicate clusters at the atomic scale on a short time-scale. We were able to deconvolute these effects by combining time resolved small angle x-ray scattering, nuclear magnetic resonance, and parallel tempering simulations. We investigated the impact of a dilution by pure water or by a sodium hydroxide solution on the speciation and size of the dissolved silicates in solution. Herein, we show that the silicate clusters are not affected by dilution, suggesting that sodium cations protect the silicate clusters from hydrolysis. Decondensation is triggered by hydroxide ions that weaken and break Si–O bonds. Alongside the decondensation, the evolution of the computed protonation state of the silica species indicates a change in the interaction potential. Our results pave the way towards the investigation at the atomic scale of more complex systems implying alkali silicate solutions in condensation process by the addition of calcium or aluminum to synthesize aluminosilicate binders, hydrogels or zeolites.

Correlation effects in parallel tempering and the role of the swapping frequency.

In this work, I investigate the effect of the swapping time frequency in parallel tempering (or replica exchange molecular dynamics, REMD) on the sampled equilibrium distributions at the different temperatures. By means of a simple deterministic Ising model originally introduced by M. Creutz (Phys. Rev. Lett., 1983, 50, 1411) I show that high frequency swaps can induce a systematic bias on the sampled REMD equilibrium distributions. The origin of this failure is ascribed to correlation effects among close temperature configurations used to evaluate REMD swapping probabilities. These results should serve as a monitor for using too frequent swapping attempts in parallel tempering simulations of generic Hamiltonians, including the ones used in atomistic simulations.

Numerical study of barriers and valleys in the free-energy landscape of spin glasses

We study the problem of glassy relaxations in the presence of an external field in the highly controlled context of a spin-glass simulation. We consider a small spin glass in three dimensions (specifically, a lattice of size L = 8, small enough to be equilibrated through a parallel tempering simulations at low temperatures, deep in the spin glass phase). After equilibrating the sample, an external field is switched on, and the subsequent dynamics is studied. The field turns out to reduce the relaxation time, but huge statistical fluctuations are found when different samples are compared. After taking care of these fluctuations we find that the expected linear regime is very narrow. Nevertheless, when regarded as a purely numerical method, we find that the external field is extremely effective in reducing the relaxation times.

Optimization of replica exchange temperature ladder under the well-tempered ensemble

Selection of temperature sequence and bias factor is a critical step for the preparation of Parallel Tempering simulations in a Well-Tempered Ensemble (PT-WTE). We provide a scheme in this study to generate the temperature sequence and the corresponding bias factor optimally. The number of replicas for a complete coverage of a specific temperature range is adjustable in this scheme, while keeping the average acceptance probability (AAP) between neighboring replica-pairs unchanged. Two series of PT-WTE simulations were tested with the number of replicas as 16 and 8 respectively. This work leads to a better understanding and application of PT-WTE simulations.

Dynamic and thermodynamic crossover scenarios in the Kob-Andersen mixture: Insights from multi-CPU and multi-GPU simulations.

The physical behavior of glass-forming liquids presents complex features of both dynamic and thermodynamic nature. Some studies indicate the presence of thermodynamic anomalies and of crossovers in the dynamic properties, but their origin and degree of universality is difficult to assess. Moreover, conventional simulations are barely able to cover the range of temperatures at which these crossovers usually occur. To address these issues, we simulate the Kob-Andersen Lennard-Jones mixture using efficient protocols based on multi-CPU and multi-GPU parallel tempering. Our setup enables us to probe the thermodynamics and dynamics of the liquid at equilibrium well below the critical temperature of the mode-coupling theory, [Formula: see text]. We find that below [Formula: see text] the analysis is hampered by partial crystallization of the metastable liquid, which nucleates extended regions populated by large particles arranged in an fcc structure. By filtering out crystalline samples, we reveal that the specific heat grows in a regular manner down to [Formula: see text] . Possible thermodynamic anomalies suggested by previous studies can thus occur only in a region of the phase diagram where the system is highly metastable. Using the equilibrium configurations obtained from the parallel tempering simulations, we perform molecular dynamics and Monte Carlo simulations to probe the equilibrium dynamics down to [Formula: see text]. A temperature-derivative analysis of the relaxation time and diffusion data allows us to assess different dynamic scenarios around [Formula: see text]. Hints of a dynamic crossover come from analysis of the four-point dynamic susceptibility. Finally, we discuss possible future numerical strategies to clarify the nature of crossover phenomena in glass-forming liquids.

Water isotope effect on the thermostability of a polio viral RNA hairpin: A metadynamics study.

Oral polio vaccine is considered to be the most thermolabile of all the common childhood vaccines. Despite heavy water (D2O) having been known for a long time to stabilise attenuated viral RNA against thermodegradation, the molecular underpinnings of its mechanism of action are still lacking. Whereas, understanding the basis of D2O action is an important step that might reform the way other thermolabile drugs are stored and could possibly minimize the cold chain problem. Here using a combination of parallel tempering and well-tempered metadynamics simulation in light water (H2O) and in D2O, we have fully described the free energy surface associated with the folding/unfolding of a RNA hairpin containing a non-canonical basepair motif, which is conserved within the 3'-untranslated region of poliovirus-like enteroviruses. Simulations reveal that in heavy water (D2O) there is a considerable increase of the stability of the folded basin as monitored through an intramolecular hydrogen bond (HB), size, shape, and flexibility of RNA structures. This translates into a higher melting temperature in D2O by 41 K when compared with light water (H2O). We have explored the hydration dynamics of the RNA, hydration shell around the RNA surface, and spatial dependence of RNA-solvent collective HB dynamics in the two water systems. Simulation in heavy water clearly showed that D2O strengthens the HB network in the solvent, lengthens inter-residue water-bridge lifetime, and weakens dynamical coupling of the hairpin to its solvation environment, which enhances the rigidity of solvent exposed sites of the native configurations. The results might suggest that like other added osmoprotectants, D2O can act as a thermostabilizer when used as a solvent.

Simulated-quantum-annealing comparison between all-to-all connectivity schemes

Quantum annealing aims to exploit quantum mechanics to speed up the search for the solution to optimization problems. Most problems exhibit complete connectivity between the logical spin variables after they are mapped to the Ising spin Hamiltonian of quantum annealing. To account for hardware constraints of current and future physical quantum annealers, methods enabling the embedding of fully connected graphs of logical spins into a constant-degree graph of physical spins are therefore essential. Here, we compare the recently proposed embedding scheme for quantum annealing with all-to-all connectivity due to Lechner, Hauke and Zoller (LHZ) [Science Advances 1 (2015)] to the commonly used minor embedding (ME) scheme. Using both simulated quantum annealing and parallel tempering simulations, we find that for a set of instances randomly chosen from a class of fully connected, random Ising problems, the ME scheme outperforms the LHZ scheme when using identical simulation parameters, despite the fault tolerance of the latter to weakly correlated spin-flip noise. This result persists even after we introduce several decoding strategies for the LHZ scheme, including a minimum-weight decoding algorithm that results in substantially improved performance over the original LHZ scheme. We explain the better performance of the ME scheme in terms of more efficient spin updates, which allows it to better tolerate the correlated spin-flip errors that arise in our model of quantum annealing. Our results leave open the question of whether the performance of the two embedding schemes can be improved using scheme-specific parameters and new error correction approaches.

Improved estimation of density of states for Monte Carlo sampling via MBAR.

We present a new method to calculate the density of states using the multistate Bennett acceptance ratio (MBAR) estimator. We use a combination of parallel tempering (PT) and multicanonical simulation to demonstrate the efficiency of our method in a statistical model of sampling from a two-dimensional normal mixture and also in a physical model of aggregation of lattice polymers. While MBAR has been commonly used for final estimation of thermodynamic properties, our numerical results show that the efficiency of estimation with our new approach, which uses MBAR as an intermediate step, often improves upon conventional use of MBAR. We also demonstrate that it can be beneficial in our method to use full PT samples for MBAR calculations in cases where simulation data exhibit long correlation.

Superposition-Enhanced Estimation of Optimal Temperature Spacings for Parallel Tempering Simulations.

Effective parallel tempering simulations rely crucially on a properly chosen sequence of temperatures. While it is desirable to achieve a uniform exchange acceptance rate across neighboring replicas, finding a set of temperatures that achieves this end is often a difficult task, in particular for systems undergoing phase transitions. Here we present a method for determination of optimal replica spacings, which is based upon knowledge of local minima in the potential energy landscape. Working within the harmonic superposition approximation, we derive an analytic expression for the parallel tempering acceptance rate as a function of the replica temperatures. For a particular system and a given database of minima, we show how this expression can be used to determine optimal temperatures that achieve a desired uniform acceptance rate. We test our strategy for two atomic clusters that exhibit broken ergodicity, demonstrating that our method achieves uniform acceptance as well as significant efficiency gains.

Direct Mixing of Atomistic Solutes and Coarse-Grained Water.

We present a new dual-resolution approach for coupling atomistic and coarse-grained models in molecular dynamics simulations of hydrated systems. In particular, a coarse-grained point dipolar water model is used to solvate molecules represented with standard all-atom force fields. A unique characteristic of our methodology is that the mixing of resolutions is direct, meaning that no additional or ad hoc scaling factors, intermediate regions, or extra sites are required. To validate the methodology, we compute the hydration free energy of 14 atomistic small molecules (analogs of amino acid side chains) solvated by the coarse-grained water. Remarkably, our predictions reproduce the experimental data as accurately as the predictions from state-of-the-art fully atomistic simulations. We also show that the hydration free energy of the coarse-grained water itself is in comparable or better agreement with the experimental value than the predictions from all but one of the most common multisite atomistic models. The coarse-grained water is then applied to solvate a typical atomistic protein containing both α-helix and β-strand elements. Moreover, parallel tempering simulations are performed to investigate the folding free energy landscape of a representative α helical and a β hairpin structure. For the simulations considered in this work, our dual-resolution method is found to be 3 to 6 times more computationally efficient than corresponding fully atomistic approaches.

Parallel tempering simulation of the three-dimensional Edwards–Anderson model with compact asynchronous multispin coding on GPU

Monte Carlo simulations of the Ising model play an important role in the field of computational statistical physics, and they have revealed many properties of the model over the past few decades. However, the effect of frustration due to random disorder, in particular the possible spin glass phase, remains a crucial but poorly understood problem. One of the obstacles in the Monte Carlo simulation of random frustrated systems is their long relaxation time making an efficient parallel implementation on state-of-the-art computation platforms highly desirable. The Graphics Processing Unit (GPU) is such a platform that provides an opportunity to significantly enhance the computational performance and thus gain new insight into this problem. In this paper, we present optimization and tuning approaches for the CUDA implementation of the spin glass simulation on GPUs. We discuss the integration of various design alternatives, such as GPU kernel construction with minimal communication, memory tiling, and look-up tables. We present a binary data format, Compact Asynchronous Multispin Coding (CAMSC), which provides an additional 28.4% speedup compared with the traditionally used Asynchronous Multispin Coding (AMSC). Our overall design sustains a performance of 33.5 ps per spin flip attempt for simulating the three-dimensional Edwards–Anderson model with parallel tempering, which significantly improves the performance over existing GPU implementations.

MCMC2 (version 1.1): A Monte Carlo code for multiply-charged clusters

MCMC2 (version 1.1): A Monte Carlo code for multiply-charged clusters

Calculation of the “absolute” free energy of a β-hairpin in an all-atom force field

We propose a new approach to calculate the conformational free energy of a macromolecule in a compact stable state in implicit solvent. The free energy is evaluated with respect to an artificial reference system without internal interactions that is confined to a small well-defined multidimensional volume of a regular shape occupying approximately the same part of the conformational space as the macromolecule of interest. We present a practical implementation of our method, successfully apply it to a β-hairpin in all-atom representation, verify the results with direct parallel tempering simulations, and discuss the possibilities of further improvements.

Copper, BDNF and Its N‐terminal Domain: Inorganic Features and Biological Perspectives

Brain-derived neurotrophic factor (BDNF) is a neurotrophin that influences development, maintenance, survival, and synaptic plasticity of central and peripheral nervous systems. Altered BDNF signaling is involved in several neurodegenerative disorders including Alzheimer's disease. Metal ions may influence the BDNF activity and it is well known that the alteration of Cu(2+) homeostasis is a prominent factor in the development of neurological pathologies. The N-terminal domain of BDNF represents the recognition site of its specific receptor TrkB, and metal ions interaction with this protein domain may influence the protein/receptor interaction. In spite of this, no data inherent the interaction of BDNF with Cu(2+) ions has been reported up to now. Cu(2+) complexes of the peptide fragment BDNF(1-12) encompassing the sequence 1-12 of N-terminal domain of human BDNF protein were characterized by means of potentiometry, spectroscopic methods (UV/Vis, CD, EPR), parallel tempering simulations and DFT-geometry optimizations. Coordination features of the acetylated form, Ac-BDNF(1-12), were also characterized to understand the involvement of the terminal amino group. Whereas, an analogous peptide, BDNF(1-12)D3N, in which the aspartate residue was substituted by an asparagine, was synthesized to provide evidence on the possible role of carboxylate group in Cu(2+) coordination. The results demonstrated that the amino group is involved in metal binding and the metal coordination environment of the predominant complex species at physiological pH consisted of one amino group, two amide nitrogen atoms, and one carboxylate group. Noteworthy, a strong decrease of the proliferative activity of both BDNF(1-12) and the whole protein on a SHSY5Y neuroblastoma cell line was found after treatment in the presence of Cu(2+). The effect of metal addition is opposite to that observed for the analogous fragment of nerve growth factor (NGF) protein, highlighting the role of specific domains, and suggesting that Cu(2+) may drive different pathways for the BDNF and NGF in physiological as well as pathological conditions.

Parallel Tempering Simulation on Generalized Canonical Ensemble

AbstractParallel tempering simulation is widely used in enhanced sampling of systems with complex energy surfaces. We hereby introduce generalized canonical ensemble (GCE) instead of the usual canonical ensemble into the parallel tempering to further improve abilities of the simulation technique. GCE utilizes an adapted weight function to obtain a unimodal energy distribution even in phase-coexisting region and then the parallel tempering on GCE yields the steady swap acceptance rates (SARs) instead of the fluctuated SARs in that on canonical ensemble. With the steady SARs, we can facilitate assign the parameters of the parallel tempering simulation to more efficiently reach equilibrium among different phases. We illustrate the parallel tempering simulation on GCE in the phase-coexisting region of 2-dimensional Potts model, a benchmark system for new simulation method developing. The result indicates that the new parallel tempering method is more efficient to estimate statistical quantities (i.e., to sample the conformational space) than the normal parallel tempering, specially in phase-coexisting regions of larger systems.