Parallel Tempering Simulations Research Articles

Dynamics and optimal number of replicas in parallel tempering simulations

We study the dynamics of parallel tempering simulations, also known as the replica exchange technique, which has become the method of choice for simulation of proteins and other complex systems. Recent results for the optimal choice of the control parameter discretization allow a treatment independent of the system in question. By analyzing mean first passage times across the control parameter space, we find an expression for the optimal number of replicas in simulations covering a given temperature range. Our results suggest a particular protocol to optimize the number of replicas in actual simulations.

Optimized folding simulations of protein A

We describe optimized parallel tempering simulations of the 46-residue B-fragment of protein A. Native-like configurations with a root-mean-square deviation of approximately 3 A to the experimentally determined structure (Protein Data Bank identifier 1BDD) are found. However, at biologically relevant temperatures such conformations appear with only approximately 10 % frequency in our simulations. Possible shortcomings in our energy function are discussed.

Efficient parallel tempering for first-order phase transitions

We present a Monte Carlo algorithm that facilitates efficient parallel tempering simulations of the density of states g(E) . We show that the algorithm eliminates the supercritical slowing down in the case of the Q=20 and Q=256 Potts models in two dimensions, typical examples for systems with extreme first-order phase transitions. As recently predicted, and shown here, the microcanonical heat capacity along the calorimetric curve has negative values for finite systems.

Equilibrium density of states and thermodynamic properties of a model glass former

This paper presents an analysis of the thermodynamics of a model glass former. We have performed equilibrium sampling of a popular binary Lennard-Jones model, employing parallel tempering Monte Carlo to cover the crystalline, amorphous, and liquid regions of configuration space. Disconnectivity graphs are used to visualize the potential energy landscape in the vicinity of a crystalline geometry and in an amorphous region of configuration space. The crystalline global minimum is separated from the bulk of the minima by a large potential energy gap, leading to broken ergodicity in conventional simulations. Our sampling reveals crystalline global minima that are lower in potential energy than some of the previous candidates. We present equilibrium thermodynamic properties based on parallel tempering simulations, including heat capacities and free energy profiles, which depend explicitly on the crystal structure. We also report equilibrium melting temperatures.

Improving an all-atom force field

Experimentally well-characterized proteins that are small enough to be computationally tractable provide useful information for refining existing all-atom force fields. This is used by us for reparametrizing a recently developed all-atom force field. Relying on high statistics parallel tempering simulations of a designed 20 residue beta-sheet peptide, we propose incremental changes that improve the force field's range of applicability.

Coarse-Grained Molecular Dynamics Modeling of Strongly Associating Fluids: Thermodynamics, Liquid Structure, and Dynamics of Symmetric Binary Mixture Fluids

Symmetric binary mixtures capable of strong association via a highly directional and saturable specific interaction between unlike molecules are investigated by canonical molecular dynamics simulations. The specific interaction of the molecules is defined in a new coarse-grained pair potential that can be applied in continuous molecular dynamics as well as in Monte Carlo simulations. The thermodynamic, structural, and dynamic properties of the associating mixture fluids are investigated as a function of density, temperature, and association strength of the specific interaction. Detailed analysis of the simulation data confirms a two-stage mechanism in the formation of specific bonds with increasing interaction strength, including a fast dimerization process and a subsequent stage of perfecting the bonds. A large heat capacity peak is found during the formation or breaking of the bonds, reflecting the large energy fluctuation introduced by the strong association. The fractions of nonbonded molecules obtained from the simulations as a function of density, temperature, and interaction strength are in excellent agreement with the predictions of Wertheim's thermodynamic perturbation theory. The translational and rotational dynamics of the Tmer mixture are effectively retarded with increasing association strength and are analyzed in terms of autocorrelation functions and a non-Gaussian parameter for the translational dynamics. The lifetimes of molecules in bonded and nonbonded states provide detailed information about the transformation of molecules between the bonded state and the nonbonded state. Finally, simulation sampling problems inherent to strongly interacting systems are easily overcome using the parallel tempering simulation technique. This latter result confirms that with the new continuous coarse-grained simulation potential we have a versatile and flexible interaction potential that can be used with many available molecular dynamics and Monte Carlo algorithms under various ensembles.

Rigid quantum Monte Carlo simulations of condensed molecular matter: Water clusters in the n=2→8 range

The numerical advantage of quantum Monte Carlo simulations of rigid bodies relative to the flexible simulations is investigated for some simple systems. The results show that if high frequency modes in molecular condensed matter are predominantly in the ground state, the convergence of path integral simulations becomes nonuniform. Rigid body quantum parallel tempering simulations are necessary to accurately capture thermodynamic phenomena in the temperature range where the dynamics are influenced by intermolecular degrees of freedom; the stereographic projection path integral adapted for quantum simulations of asymmetric tops is a significantly more efficient strategy compared with Cartesian coordinate simulations for molecular condensed matter under these conditions. The reweighted random series approach for stereographic path integral Monte Carlo is refined and implemented for the quantum simulation of water clusters treated as an assembly of rigid asymmetric tops.

ALPHA-HELIX FORMATION IN C-PEPTIDE RNASE-A INVESTIGATED BY PARALLEL TEMPERING SIMULATIONS

We have performed parallel tempering simulations of a 13-residue peptide fragment of ribonuclease-A, c-peptide, in implicit solvent with constant dielectric permittivity. This peptide has a strong tendency to form α-helical conformations in solvent as suggested by circular dichroism (CD) and nuclear magnetic resonance (NMR) experiments. Our results demonstrate that 5th and 8–12 residues are in the α-helical region of the Ramachandran map for global minimum energy state in solvent environment. Effects of salt bridge formation on stability of α-helix structure are discussed.

Modelling the thermal stability of precursor nanoparticles in zeolite synthesis

We have modelled the thermal stability of silica nanoparticles commonly observed as precursors in the synthesis of zeolites. We performed canonical Monte Carlo and parallel tempering simulations on a lattice model that describes the self-assembly of nanoparticles under conditions at which they are observed experimentally. The effect of heating on the relative stability of the phases of the model was analysed by running simulations at various temperatures. At low temperature, the model yields a metastable multi-particle phase with a characteristic size distribution, which is separated by an energy barrier from the true equilibrium phase, a dense silica solid. As temperature increases, the system enters a transition region and eventually reaches the bulk phase. This transition is reminiscent of the experimentally observed transition from nanoparticles to zeolite. The transition temperature scales with the inverse of the system volume, approaching an asymptotic value for large system sizes. This indicates the transition temperature is a reproducible macroscopic property of the system. The transition temperature in the model is within the range of temperatures at which nanoparticles form zeolite crystals in experiments.

Use of the Weighted Histogram Analysis Method for the Analysis of Simulated and Parallel Tempering Simulations.

The growing adoption of generalized-ensemble algorithms for biomolecular simulation has resulted in a resurgence in the use of the weighted histogram analysis method (WHAM) to make use of all data generated by these simulations. Unfortunately, the original presentation of WHAM by Kumar et al. is not directly applicable to data generated by these methods. WHAM was originally formulated to combine data from independent samplings of the canonical ensemble, whereas many generalized-ensemble algorithms sample from mixtures of canonical ensembles at different temperatures. Sorting configurations generated from a parallel tempering simulation by temperature obscures the temporal correlation in the data and results in an improper treatment of the statistical uncertainties used in constructing the estimate of the density of states. Here we present variants of WHAM, STWHAM and PTWHAM, derived with the same set of assumptions, that can be directly applied to several generalized ensemble algorithms, including simulated tempering, parallel tempering (better known as replica-exchange among temperatures), and replica-exchange simulated tempering. We present methods that explicitly capture the considerable temporal correlation in sequentially generated configurations using autocorrelation analysis. This allows estimation of the statistical uncertainty in WHAM estimates of expectations for the canonical ensemble. We test the method with a one-dimensional model system and then apply it to the estimation of potentials of mean force from parallel tempering simulations of the alanine dipeptide in both implicit and explicit solvent.

Folding of Proteins with Diverse Folds

Using parallel tempering simulations with high statistics, we investigate the folding and thermodynamic properties of three small proteins with distinct native folds: the all-helical 1RIJ, the all-sheet beta3s, and BBA5, which has a mixed helix-sheet fold. In all three cases, simulations with our energy function find the native structures as global minima in free energy at experimentally relevant temperatures. However, the folding process strongly differs for the three molecules, indicating that the folding mechanism is correlated with the form of the native structure.

Optimized parallel tempering simulations of proteins

We apply a recently developed adaptive algorithm that systematically improves the efficiency of parallel tempering or replica exchange methods in the numerical simulation of small proteins. Feedback iterations allow us to identify an optimal set of temperatures/replicas which are found to concentrate at the bottlenecks of the simulations. A measure of convergence for the equilibration of the parallel tempering algorithm is discussed. We test our algorithm by simulating the 36-residue villin headpiece subdomain HP-36 where we find a lowest-energy configuration with a root-mean-square deviation of less than 4 A to the experimentally determined structure.

Optimal Estimates of Free Energies from Multistate Nonequilibrium Work Data

We derive the optimal estimates of the free energies of an arbitrary number of thermodynamic states from nonequilibrium work measurements; the work data are collected from forward and reverse switching processes and obey a fluctuation theorem. The maximum likelihood formulation properly reweights all pathways contributing to a free energy difference and is directly applicable to simulations and experiments. We demonstrate dramatic gains in efficiency by combining the analysis with parallel tempering simulations for alchemical mutations of model amino acids.

Improved Efficiency of Replica Exchange Simulations through Use of a Hybrid Explicit/Implicit Solvation Model.

The use of parallel tempering or replica exchange molecular dynamics (REMD) simulations has facilitated the exploration of free energy landscapes for complex molecular systems, but application to large systems is hampered by the scaling of the number of required replicas with increasing system size. Use of continuum solvent models reduces system size and replica requirements, but these have been shown to provide poor results in many cases, including overstabilization of ion pairs and secondary structure bias. Hybrid explicit/continuum solvent models can overcome some of these problems through an explicit representation of water molecules in the first solvation shells, but these methods typically require restraints on the solvent molecules and show artifacts in water properties due to the solvation interface. We propose an REMD variant in which the simulations are performed with a fully explicit solvent, but the calculation of exchange probability is carried out using a hybrid model, with the solvation shells calculated on the fly during the fully solvated simulation. The resulting reduction in the perceived system size in the REMD exchange calculation provides a dramatic decrease in the computational cost of REMD, while maintaining a very good agreement with results obtained from the standard explicit solvent REMD. We applied several standard and hybrid REMD methods with different solvent models to alanine polymers of 1, 3, and 10 residues, obtaining ensembles that were essentially independent of the initial conformation, even with explicit solvation. Use of only a continuum model without a shell of explicit water provided poor results for Ala3 and Ala10, with a significant bias in favor of the α-helix. Likewise, using only the solvation shells and no continuum model resulted in ensembles that differed significantly from the standard explicit solvent data. Ensembles obtained from hybrid REMD are in very close agreement with explicit solvent data, predominantly populating polyproline II conformations. Inclusion of a second shell of explicit solvent was found to be unnecessary for these peptides.

Equilibrium thermodynamics from basin-sampling

We present a "basin-sampling" approach for calculation of the potential energy density of states for classical statistical models. It combines a Wang-Landau-type uniform sampling of local minima and a novel approach for approximating the relative contributions from local minima in terms of the volumes of basins of attraction. We have employed basin-sampling to study phase changes in atomic clusters modeled by the Lennard-Jones potential and for ionic clusters. The approach proves to be efficient for systems involving broken ergodicity and has allowed us to calculate converged heat capacity curves for systems that could previously only be treated using the harmonic superposition approximation. Benchmarks are also provided by comparison with parallel tempering and Wang-Landau simulations, where these proved feasible.

Effect of Single-Point Sequence Alterations on the Aggregation Propensity of a Model Protein

Sequences of contemporary proteins are believed to have evolved through a process that optimized their overall fitness, including their resistance to deleterious aggregation. Biotechnological processing may expose therapeutic proteins to conditions that are much more conducive to aggregation than those encountered in a cellular environment. An important task of protein engineering is to identify alternative sequences that would protect proteins when processed at high concentrations without altering their native structure associated with specific biological function. Our computational studies exploit parallel tempering simulations of coarse-grained model proteins to demonstrate that isolated amino acid residue substitutions can result in significant changes in the aggregation resistance of the protein in a crowded environment while retaining protein structure in isolation. A thermodynamic analysis of protein clusters subject to competing processes of folding and association shows that moderate mutations can produce effects similar to those caused by changes in system conditions, including temperature, concentration, and solvent composition, that affect the aggregation propensity. The range of conditions where a protein can resist aggregation can therefore be tuned by sequence alterations, although the protein generally may retain its generic ability for aggregation.

All-exchanges parallel tempering

An alternative exchange strategy for parallel tempering simulations is introduced. Instead of attempting to swap configurations between two randomly chosen but adjacent replicas, the acceptance probabilities of all possible swap moves are calculated a priori. One specific swap move is then selected according to its probability and enforced. The efficiency of the method is illustrated first on the case of two Lennard-Jones (LJ) clusters containing 13 and 31 atoms, respectively. The convergence of the caloric curve is seen to be at least twice as fast as in conventional parallel tempering simulations, especially for the difficult case of LJ31. Further evidence for an improved efficiency is reported on the ergodic measure introduced by Mountain and Thirumalai [J. Phys. Chem. 93, 6975 (1989)], calculated here for LJ13 close to the melting point. Finally, tests on two simple spin systems indicate that the method should be particularly useful when a limited number of replicas are available.

Parallel tempering: Theory, applications, and new perspectives

We review the history of the parallel tempering simulation method. From its origins in data analysis, the parallel tempering method has become a standard workhorse of physicochemical simulations. We discuss the theory behind the method and its various generalizations. We mention a selected set of the many applications that have become possible with the introduction of parallel tempering, and we suggest several promising avenues for future research.

Optimal allocation of replicas in parallel tempering simulations

We have studied the efficiency of parallel tempering simulations for a variety of systems including a coarse-grained protein, an atomistic model polypeptide, and the Lennard-Jones fluid. A scheme is proposed for the optimal allocation of temperatures in these simulations. The method is compared to the existing empirical approaches used for this purpose. Accuracy associated with the computed thermodynamic quantities such as specific heat is also computed and their dependence on the trial-exchange acceptance rate is reported.

Optimal Allocation of Replicas to Processors in Parallel Tempering Simulations

The optimal allocation of replicas to a homogeneous or heterogeneous set of processors is derived for parallel tempering simulations on multiprocessor machines. In the general case, it is possible without substantially increasing wall clock time to achieve nearly perfect utilization of CPU time. Random fluctuations in the execution time of each replica do not significantly degrade the performance of the scheduler.