Multicanonical Method Research Articles

Restricted multicanonical sampling for machine learning potential construction

The restricted multicanonical (rMUCA) ensemble method is developed and combined with the on-the-fly machine learning potential (MLP) generation scheme. The rMUCA simulation performs a random walk in the potential-energy subspace restricted by the selected collective variables and allows us to sample physically relevant configurations without being trapped in local energy minima. No preliminary simulation runs are required to construct a bias potential. The sample structures for training are collected dynamically from the simulations using the MLP itself where the simultaneous error estimation is utilized to judge whether an updated structure should be added to the sample data set or not. The rMUCA formula can be also used for the saddle-point search with minor modification. The method is applied to the oxidation of carbon monoxide on platinum surfaces. The results show that the rMUCA simulation provides an efficient and accurate way to sample rare events for the MLP construction.

Calculation of the residual entropy of Ice Ih by Monte Carlo simulation with the combination of the replica-exchange Wang-Landau algorithm and multicanonical replica-exchange method.

We estimated the residual entropy of Ice Ih by the recently developed simulation protocol, namely, the combination of the replica-exchange Wang-Landau algorithm and multicanonical replica-exchange method. We employed a model with the nearest neighbor interactions on the three-dimensional hexagonal lattice, which satisfied the ice rules in the ground state. The results showed that our estimate of the residual entropy is in accordance with various previous results. In this article, we not only give our latest estimate of the residual entropy of Ice Ih but also discuss the importance of the uniformity of a random number generator in Monte Carlo simulations.

Nonflat histogram techniques for spin glasses.

We study the bimodal Edwards-Anderson spin glass comparing established methods, namely the multicanonical method, the 1/k ensemble, and parallel tempering, to an approach where the ensemble is modified by simulating power-law-shaped histograms in energy instead of flat histograms as in the standard multicanonical case. We show that by this modification a significant speed-up in terms of mean round-trip times can be achieved for all lattice sizes taken into consideration.

Efficient simulation protocol for determining the density of states: Combination of replica-exchange Wang-Landau method and multicanonical replica-exchange method.

By combining two generalized-ensemble algorithms, the replica-exchange Wang-Landau (REWL) method and the multicanonical replica-exchange method (MUCAREM), we propose an effective simulation protocol to determine the density of states with high accuracy. The new protocol is referred to as REWL-MUCAREM, and REWL is first performed and then MUCAREM is performed. In order to verify the effectiveness of our protocol, we performed simulations of a square-lattice Ising model using the three methods, namely REWL, MUCAREM, and REWL-MUCAREM. The results showed that the density of states obtained by REWL-MUCAREM is more accurate than that is estimated by the two methods separately.

Histogram-free multicanonical Monte Carlo sampling to calculate the density of states

We report a new multicanonical Monte Carlo algorithm to obtain the density of states for physical systems with continuous state variables in statistical mechanics. Our algorithm is able to obtain a closed-form expression for the density of states expressed in a chosen basis set, instead of a numerical array of finite resolution as in previous variants of this class of MC methods such as the multicanonical sampling and Wang–Landau sampling. This is enabled by storing the visited states directly and avoiding the explicit collection of a histogram. This practice also has the advantage of avoiding undesirable artificial errors caused by the discretization and binning of continuous state variables. Our results show that this scheme is capable of obtaining converged results with a much reduced number of Monte Carlo steps, leading to a significant speedup over existing algorithms. Program summaryProgram Title: HistogramFreeMUCAProgram Files doi:http://dx.doi.org/10.17632/mtgyvsfzyg.1Licensing provisions: BSD 3-clauseProgramming language: PythonNature of problem: This program implements a novel algorithm to obtain the density of states of a physical system, expanded in a chosen basis set. Unlike existing algorithms that return the density of states as a numerical array, this algorithm avoids binning of a continuous variable and is able to express the density of states as a closed-form expression. It is thus suitable for the study of the statistical mechanics and thermodynamic properties of physical systems where the density of a continuous state variable is of interest.Solution method: The new algorithm presented here is a special reweighting method in classical Monte Carlo approaches. In particular, it can be regarded as a descendant and a hybrid method closely related to the multicanonical method and Wang–Landau sampling.Additional comments: Most updated source code can be found at: https://github.com/yingwaili/HistogramFreeMUCA

Development of a generalized hybrid Monte Carlo algorithm to generate the multicanonical ensemble with applications to molecular systems.

In the present paper, a generalized hybrid Monte Carlo method to generate the multicanonical ensemble has been developed, which is a generalization of the multicanonical hybrid Monte Carlo (HMC) method by Hansmann and co-workers [Chem. Phys. Lett. 259, 321 (1996)]. The generalized hybrid Monte Carlo (GHMC) method is an equations-of-motion guided Monte Carlo combined with partial momentum refreshment. We successfully applied our multicanonical GHMC to dense Lennard-Jones fluids and a coarse grained protein model. It is found that good computational efficiency can be gained in the case of the acceptance ratio around 60% for the models examined. While a large number of molecular dynamics (MD) steps in a single GHMC cycle is needed to yield good computational efficiency at a large mixing ratio of momenta with thermal noise vectors, corresponding to the original multicanonical HMC method, a small number of MD steps are enough to achieve good efficiency at a small mixing ratio. This property is useful to develop a composite algorithm combining the present GHMC method with other Monte Carlo moves.

Multicanonical Monte Carlo evaluation of drag probability distribution functions

The probability distribution of the drag generated by a two-dimensional square/rectangular obstacle is calculated both for quasi-random input flow patterns and for random surface roughness by employing the multicanonical Monte Carlo procedure in conjunction with the lattice Boltzmann method. The results demonstrate that the multicanonical method can estimate the probability distribution function in low-probability regions with far less computational effort than standard techniques.

Massively parallel multicanonical simulations

Generalized-ensemble Monte Carlo simulations such as the multicanonical method and similar techniques are among the most efficient approaches for simulations of systems undergoing discontinuous phase transitions or with rugged free-energy landscapes. As Markov chain methods, they are inherently serial computationally. It was demonstrated recently, however, that a combination of independent simulations that communicate weight updates at variable intervals allows for the efficient utilization of parallel computational resources for multicanonical simulations. Implementing this approach for the many-thread architecture provided by current generations of graphics processing units (GPUs), we show how it can be efficiently employed with of the order of 104 parallel walkers and beyond, thus constituting a versatile tool for Monte Carlo simulations in the era of massively parallel computing. We provide the fully documented source code for the approach applied to the paradigmatic example of the two-dimensional Ising model as starting point and reference for practitioners in the field. Program summaryProgram Title: cudamucaProgram Files doi:http://dx.doi.org/10.17632/tzhfpdymv9.1Licensing provisions: Creative Commons Attribution license (CC BY 4.0)Programming language: C++, CUDAExternal routines/libraries: NVIDIA CUDA Toolkit 6.5 or newerNature of problem: The program determines weights for a multicanonical simulation of the 2D Ising model to result in a flat energy histogram. A final production run with these weights provides an estimate of the density of states of the model.Solution method: The code uses a parallel variant of the multicanonical method employing many parallel walkers that accumulate a common histogram. The resulting histogram is used to determine the weight function for the next iteration. Once the iteration has converged, simulations visit all possible energies with the same probability.Additional comments including restrictions and unusual features: The system size and size of the population of replicas are limited depending on the memory of the GPU device used. Code repository at https://github.com/CQT-Leipzig/cudamuca.

Multicanonical sampling of the space of states of ℋ(2, n)-vector models

Problems of temperature behavior of specific heat are solved by the entropy simulation method for Ising models on a simple square lattice and a square spin ice (SSI) lattice with nearest neighbor interaction, models of hexagonal lattices with short-range (SR) dipole interaction, as well as with long-range (LR) dipole interaction and free boundary conditions, and models of spin quasilattices with finite interaction radius. It is established that systems of a finite number of Ising spins with LR dipole interaction can have unusual thermodynamic properties characterized by several specific-heat peaks in the absence of an external magnetic field. For a parallel multicanonical sampling method, optimal schemes are found empirically for partitioning the space of states into energy bands for Ising and SSI models, methods of concatenation and renormalization of histograms are discussed, and a flatness criterion of histograms is proposed. It is established that there is no phase transition in a model with nearest neighbor interaction on a hexagonal lattice, while the temperature behavior of specific heat exhibits singularity in the same model, in case of LR interaction. A spin quasilattice is found that exhibits a nonzero value of residual entropy.

Population-dynamics method with a multicanonical feedback control.

We discuss the Giardinà-Kurchan-Peliti population dynamics method for evaluating large deviations of time-averaged quantities in Markov processes [Phys. Rev. Lett. 96, 120603 (2006)PRLTAO0031-900710.1103/PhysRevLett.96.120603]. This method exhibits systematic errors which can be large in some circumstances, particularly for systems with weak noise, with many degrees of freedom, or close to dynamical phase transitions. We show how these errors can be mitigated by introducing control forces within the algorithm. These forces are determined by an iteration-and-feedback scheme, inspired by multicanonical methods in equilibrium sampling. We demonstrate substantially improved results in a simple model, and we discuss potential applications to more complex systems.

Knots as a Topological Order Parameter for Semiflexible Polymers.

Using a combination of the multicanonical MonteCarlo algorithm and the replica-exchange method, we investigate the influence of bending stiffness on the conformational phases of a bead-stick homopolymer model and present the pseudophase diagram for the complete range of semiflexible polymers, from flexible to stiff. Although it is a simple model, we observe a rich variety of conformational phases, reminiscent of conformations observed for synthetic polymers or biopolymers. Depending on the bending stiffness, the model exhibits different pseudophases like bent, hairpin, or toroidal. In particular, we find thermodynamically stable knots and unusual transitions into these knotted phases with a clear phase coexistence, but almost constant mean total energy, and hence almost no latent heat.

Investigating Kinetics of Conformational Change using Molecular Dynamics and Milestoning

To understand the functional roles of biomolecules, it is necessary to study both stable and metastable states in conformational space and how they are connected. Because of the advance of hardware and software, molecular dynamics (MD) simulation methods have become a feasible tool for this purpose, resulting in rather quantitative conclusions, which can be compared with experiment. MD simulation methods, however, still have a bottleneck due to “rare events” connecting stable states, and such a rare event can occur with millisecond timescales in a large biomolecule, which hampers a direct application of MD simulations to the rare event of biomolecules. This is why we need some novel ideas and algorithms to overcome this problem of rare events, and many researchers have been developing promising methods, mainly focusing on conformational sampling methods. These include very powerful replica exchange and multicanonical methods, but time information is lost in these methods and dynamic characterization of a rare event is difficult. In this presentation, we apply the milestoning method devised by Elber, Vanden-Eijnden and others to conformational change of a small peptide, chignolin, at the folding temperature (∼420K). The milestoning method is regarded as an alternative of Markov state modeling, and by defining “milestones (collections of small dividing surfaces)” in some order parameter space and counting the crossing of milestones by a trajectory, mean first passage times (MFPTs) between two milestones can be obtained and a kinetic characterization of conformational change can be achieved. Here we will examine how the choice of order parameter space affects the results of MFPTs, and also analyze the effects of non-Markovianity for this process using Zuckerman's method.

Transition Path Times in Protein Folding Studied by Structure-Based Simulation

In trajectories of protein folding and unfolding processes, the transition path that connects denatured and native states is the most interesting part. Recently, single-molecule FRET experiments together with a statistical inference theory successfully identified the transition path times for folding of some small proteins. Yet, its underlying physics is poorly understood. Here we conducted a comprehensive survey of transition path times for 29 small-to-medium two-state-folding proteins using structure-based coarse-grained molecular dynamics simulations. Using the multi-canonical ensemble method, we first identified folding transition temperature accurately. At the transition temperature, we then performed relatively long simulations observing reversible folding and unfolding events, from which we identified and analyzed the transition path times. The distribution of transition path time for each protein can be explained as free diffusion of particle in a reaction coordinate. We sought what reaction coordinates correlate with the transition path time. Among tested coordinates, we found that the average transition path time is most strongly correlated to the difference in the numbers of native contacts, i.e., contact energies, between native and denatured states. These results imply that the transition path time series can be approximated as the nearly free diffusion in the reaction coordinate of native contacts.

Application of the parallel multicanonical method to lattice gas condensation

We present the speedup from a novel parallel implementation of the multicanonical method on the example of a lattice gas in two and three dimensions. In this approach, all cores perform independent equilibrium runs with identical weights, collecting their sampled histograms after each iteration in order to estimate consecutive weights. The weights are then redistributed to all cores. These steps are repeated until the weights are converged. This procedure benefits from a minimum of communication while distributing the necessary amount of statistics efficiently.Using this method allows us to study a broad temperature range for a variety of large and complex systems. Here, a gas is modeled as particles on the lattice, which interact only with their nearest neighbors. For a fixed density this model is equivalent to the Ising model with fixed magnetization. We compare our results to an analytic prediction for equilibrium droplet formation, confirming that a single macroscopic droplet forms only above a critical density.

Multicanonical Determination of the Symbol Error Ratio of WDM Polarization Multiplexed QPSK Systems

AbstractThe multicanonical method is applied to the calculation of the symbol error ratio (SER) as a function of the optical signal to noise ratio (OSNR) at the receiver for a polarization multiplexed quadrature phase shift keying (PM-QPSK) wavelength division multiplexed (WDM) system. We improve upon previous calculations by including polarization mode dispersion (PMD) and subsequently verifying the numerical accuracy of our calculations. Our numerical studies demonstrate that acceptable accuracy can be achieved even when advancing the polarization through the fiber with relatively large propagation step lengths.

Salt Effects on Folding of a Helical Mini Protein Villin Headpiece Subdomain HP36 Studied by Generalized-Ensemble Simulations

Additives dissolved in solvent are important factors that affect proteins' stability and/or folding. In this study we investigated effects of salt ions in solvent on folding events of a helical mini protein HP36. Addition of low concentrations of ions should alter electrostatic interactions among charged groups, so that populations for conformational substates of proteins should be changed. Here we compared two data sets of folding simulations of HP36 with explicit water solvent. For efficient sampling of conformational space of the protein, multicanonical replica-exchange method was adopted.Results of the present analyses suggest that addition of ions reduces the number of nonnative, nonlocal salt bridges in the protein molecule at later stages of folding at room temperature. Especially, nonnative salt bridges between Glu5 and Arg15 and/or another between Asp4 and Lys30 have been kept in the near-native conformations in pure water. Because dehydration of the hydrophobic core of HP36 is completed only at the latest stage of folding where correct hydrophobic-core packing becomes formed, these salt bridges can prevent folding into the fully native structure of HP36 at room temperature.

Salt effects on hydrophobic-core formation in folding of a helical miniprotein studied by molecular dynamics simulations

We have investigated effects of salt ions on folding events of a helical miniprotein chicken villin headpiece subdomain HP36. Low concentrations of ions alter electrostatic interactions between charged groups of a protein and can change the populations of conformers. Here, we compare two data sets of folding simulations of HP36 in explicit water solvent with or without ions. For efficient sampling of the conformational space of HP36, the multicanonical replica-exchange molecular dynamics method was employed. Our analyses suggest that salt alters salt-bridging nature of the protein at later stages of folding at room temperature. Especially, more nonnative, nonlocal salt bridges are formed at near-native conformations in pure water. Our analyses also show that such salt-bridge formation hinders the fully native hydrophobic-core packing at the final stages of folding.

Scaling properties of a parallel implementation of the multicanonical algorithm

The multicanonical method has been proven powerful for statistical investigations of lattice and off-lattice systems throughout the last two decades. We discuss an intuitive but very efficient parallel implementation of this algorithm and analyze its scaling properties for discrete energy systems, namely the Ising model and the 8-state Potts model. The parallelization relies on independent equilibrium simulations in each iteration with identical weights, merging their statistics in order to obtain estimates for the successive weights. With good care, this allows faster investigations of large systems, because it distributes the time-consuming weight-iteration procedure and allows parallel production runs. We show that the parallel implementation scales very well for the simple Ising model, while the performance of the 8-state Potts model, which exhibits a first-order phase transition, is limited due to emerging barriers and the resulting large integrated autocorrelation times. The quality of estimates in parallel production runs remains of the same order at the same statistical cost.

Molecular simulations in generalised ensemble

In this article, we review the generalised-ensemble algorithms. Three well-known methods, namely multicanonical algorithm, simulated tempering and replica-exchange method, are described first. Both Monte Carlo and molecular dynamics versions of the algorithms are presented. We then present further extensions of the above three methods. Finally, we discuss the relations among multicanonical algorithm, Wang–Landau method and metadynamics.

Recent Applications of Replica-Exchange Molecular Dynamics Simulations of Biomolecules

Replica-exchange molecular dynamics (REMD) method is one of the enhanced conformational sampling techniques in MD simulations of proteins or other systems with rugged-energy landscapes. In REMD method, copies of original simulation system at different temperatures are simulated separately and simultaneously. Every few steps, temperatures between neighboring replicas are exchanged if the Metropolis criteria for their instantaneous potential energies are satisfied. Due to its simplicity and high efficiency in parallel computers, the method has been applied to many biological problems including protein folding, aggregation, receptor-ligand binding, and so on. In the last ten years, continuous effort to improve sampling efficiency of REMD simulations for larger biological systems has been carried out by us and other theoretical scientists. In this review article, we introduce two different approaches in REMD simulations to reduce the computational cost. One is the multicanonical replica-exchange method (MUCAREM) for reducing the number of replicas. In this method, each replica has a different multicanonical weight factor and takes a flat energy distribution to cover a wider potential energy space. Another approach is to employ implicit solvent/membrane models for representing surrounding environments of target proteins in REMD simulations. We show two applications of proteinfolding simulations in explicit solvent using the former approach and a structural prediction of a transmembrane protein dimer using the latter. Finally, we discuss possibilities of REMD method to simulate a large-scale conformational change of protein systems using massively parallel supercomputers.