^{1}and Gerhard Hummer

^{1,a)}

### Abstract

We derive simple analytical expressions for the error and computational efficiency of replica exchange molecular dynamics (REMD) simulations (and by analogy replica exchange Monte Carlo simulations). The theory applies to the important case of systems whose dynamics at long times is dominated by the slow interconversion between two metastable states. As a specific example, we consider the folding and unfolding of a protein. The efficiency is defined as the rate with which the error in an estimated equilibrium property, as measured by the variance of the estimator over repeated simulations, decreases with simulation time. For two-state systems, this rate is in general independent of the particular property. Our main result is that, with comparable computational resources used, the relative efficiency of REMD and molecular dynamics (MD) simulations is given by the ratio of the number of transitions between the two states averaged over all replicas at the different temperatures, and the number of transitions at the single temperature of the MD run. This formula applies if replica exchange is frequent, as compared to the transition times. High efficiency of REMD is thus achieved by including replica temperatures in which the frequency of transitions is higher than that at the temperature of interest. In tests of the expressions for the error in the estimator, computational efficiency, and the rate of equilibration we find quantitative agreement with the results both from kinetic models of REMD and from actual all-atom simulations of the folding of a peptide in water.

We thank Dr. A. Szabo, Dr. N.-V. Buchete, Dr. A. M. Berezhkovskii, Dr. A. B. Adib, and Dr. Y.-C. Kim for many helpful and stimulating discussions. This research used the Biowulf Linux cluster at the NIH and was supported by the Intramural Research Program of the NIDDK, NIH.

I. INTRODUCTION

II. THEORY

A. Rate model of REMD

B. Error and efficiency of MD and REMD simulations

1. MD simulations

2. REMD simulations

C. Limit of fast replica exchange

D. Error of REMD in the limit of fast exchange

E. Optimizing the number and temperature range of replicas

III. RESULTS AND DISCUSSION

A. Convergence of the exact and approximate relaxation times

B. Simulations of penta-alanine

IV. CONCLUDING REMARKS

A. Exchange frequency and number of replicas

B. Temperature range

C. Equilibration

### Key Topics

- Molecular dynamics
- 23.0
- Protein folding
- 15.0
- Proteins
- 10.0
- Statistical properties
- 7.0
- Correlation functions
- 6.0

## Figures

Temperature dependence of folding (, red) and unfolding (, blue) rates and of the relaxation rates (, black) of the -repressor fragment (Ref. 20). The temperature dependent rate constants were calculated based on Table II of Ref. 20 and the temperature dependence of the water viscosity.

Temperature dependence of folding (, red) and unfolding (, blue) rates and of the relaxation rates (, black) of the -repressor fragment (Ref. 20). The temperature dependent rate constants were calculated based on Table II of Ref. 20 and the temperature dependence of the water viscosity.

Kinetic model of replica exchange. Connectivity of all states is illustrated for (a) and (b) replicas. Black arrows represent folding and unfolding processes, while red arrows show replica exchange processes. The corresponding reduced linear kinetic schemes, obtained in the limit of fast replica exchange, are shown below the full kinetic diagrams.

Kinetic model of replica exchange. Connectivity of all states is illustrated for (a) and (b) replicas. Black arrows represent folding and unfolding processes, while red arrows show replica exchange processes. The corresponding reduced linear kinetic schemes, obtained in the limit of fast replica exchange, are shown below the full kinetic diagrams.

Efficiency gain of REMD over MD as a function of the highest replica temperature . The reference temperature is . Plotted is the efficiency in the limit of a large number of replicas, which is independent of , as obtained from Eq. (26). Note that for temperatures above 350 K, we used extrapolated rates.

Efficiency gain of REMD over MD as a function of the highest replica temperature . The reference temperature is . Plotted is the efficiency in the limit of a large number of replicas, which is independent of , as obtained from Eq. (26). Note that for temperatures above 350 K, we used extrapolated rates.

Convergence of the relaxation rate and the efficiency of REMD for the -repressor fragment (Ref. 20). (a) The exchange rate dependence of the rate limiting eigenvalue is shown for the full kinetic problem (blue line), for the exact fast replica exchange limit (black line), and for the approximate Eq. (18) (red line). The eigenvalues corresponding to MD (no replica exchange coupling), and are marked with arrows. (b) REMD efficiency [Eq. (7)] vs the target temperature. The fast-exchange limit is shown as black line with circles, the results using Eq. (25) are shown as red dashed line, the MD simulations correspond to efficiency of 1, shown in blue. The arrow indicates the variation in the efficiency for increasing replica exchange rate. The inset shows the REMD variance for different target temperatures.

Convergence of the relaxation rate and the efficiency of REMD for the -repressor fragment (Ref. 20). (a) The exchange rate dependence of the rate limiting eigenvalue is shown for the full kinetic problem (blue line), for the exact fast replica exchange limit (black line), and for the approximate Eq. (18) (red line). The eigenvalues corresponding to MD (no replica exchange coupling), and are marked with arrows. (b) REMD efficiency [Eq. (7)] vs the target temperature. The fast-exchange limit is shown as black line with circles, the results using Eq. (25) are shown as red dashed line, the MD simulations correspond to efficiency of 1, shown in blue. The arrow indicates the variation in the efficiency for increasing replica exchange rate. The inset shows the REMD variance for different target temperatures.

Convergence of the relaxation rate and the efficiency of REMD for the Arrhenius folding/unfolding rates of . (a) The exchange rate dependence of the rate limiting eigenvalue is shown for the full kinetic problem (blue line), for the exact fast replica exchange limit (black line), and for the approximate Eq. (18) (red line). The eigenvalues corresponding to MD (no replica exchange coupling), and are marked with arrows. (b) REMD efficiency [Eq. (7)] vs the target temperature. The fast-exchange limit is shown as black line with circles, the results using Eq. (25) are shown as red dashed line, the MD simulations correspond to an efficiency of 1, shown in blue. The arrow indicates the variation in the efficiency for increasing replica exchange rate. The inset shows the variance of the fraction folded as a function of the target temperature, scaled for reference and readability to a simulation length of .

Convergence of the relaxation rate and the efficiency of REMD for the Arrhenius folding/unfolding rates of . (a) The exchange rate dependence of the rate limiting eigenvalue is shown for the full kinetic problem (blue line), for the exact fast replica exchange limit (black line), and for the approximate Eq. (18) (red line). The eigenvalues corresponding to MD (no replica exchange coupling), and are marked with arrows. (b) REMD efficiency [Eq. (7)] vs the target temperature. The fast-exchange limit is shown as black line with circles, the results using Eq. (25) are shown as red dashed line, the MD simulations correspond to an efficiency of 1, shown in blue. The arrow indicates the variation in the efficiency for increasing replica exchange rate. The inset shows the variance of the fraction folded as a function of the target temperature, scaled for reference and readability to a simulation length of .

Autocorrelation function of the folding state at different temperatures obtained from REMD simulations of . Also shown are the corresponding for regular MD simulations at (blue) and (green), and the autocorrelation function of the number of folded states (red). The dashed black lines show exponential decays with the slowest relaxation rate obtained from Eq. (18).

Autocorrelation function of the folding state at different temperatures obtained from REMD simulations of . Also shown are the corresponding for regular MD simulations at (blue) and (green), and the autocorrelation function of the number of folded states (red). The dashed black lines show exponential decays with the slowest relaxation rate obtained from Eq. (18).

Variance in the fraction folded calculated from REMD simulations of for different target temperatures. For reference and readability, all results have been scaled to a simulation time of . The red symbols show data obtained from 11 independent 150-ns-long REMD simulations. The solid lines are obtained for the reduced kinetic model (blue) in the fast-exchange limit from the exact solution of the coarse grained kinetic model and in the continuum limit (black symbols) using the formula of Eq. (24). The green line corresponds to the variance for times longer MD simulations.

Variance in the fraction folded calculated from REMD simulations of for different target temperatures. For reference and readability, all results have been scaled to a simulation time of . The red symbols show data obtained from 11 independent 150-ns-long REMD simulations. The solid lines are obtained for the reduced kinetic model (blue) in the fast-exchange limit from the exact solution of the coarse grained kinetic model and in the continuum limit (black symbols) using the formula of Eq. (24). The green line corresponds to the variance for times longer MD simulations.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content