Replica exchanging self-guided Langevin dynamics for efficient and accurate conformational sampling

Abstract

This work presents a replica exchanging self-guided Langevin dynamics (RXSGLD) simulation method for efficient conformational searching and sampling. Unlike temperature-based replica exchanging simulations, which use high temperatures to accelerate conformational motion, this method uses self-guided Langevin dynamics (SGLD) to enhance conformational searching without the need to elevate temperatures. A RXSGLD simulation includes a series of SGLD simulations, with simulation conditions differing in the guiding effect and/or temperature. These simulation conditions are called stages and the base stage is one with no guiding effect. Replicas of a simulation system are simulated at the stages and are exchanged according to the replica exchanging probability derived from the SGLD partition function. Because SGLD causes less perturbation on conformational distribution than high temperatures, exchanges between SGLD stages have much higher probabilities than those between different temperatures. Therefore, RXSGLD simulations have higher conformational searching ability than temperature based replica exchange simulations. Through three example systems, we demonstrate that RXSGLD can generate target canonical ensemble distribution at the base stage and achieve accelerated conformational searching. Especially for large systems, RXSGLD has remarkable advantages in terms of replica exchange efficiency, conformational searching ability, and system size extensiveness.

Figure 1

Thermal motions of an alanine dipeptide in Langevin dynamics and in self-guided Langevin dynamics simulations. Atoms are drawn as sticks. Carbon, oxygen, nitrogen, and hydrogen atoms are colored as grey, red, blue, and white, respectively. In LD, all motions have kinetic energies equivalent to a temperature of T, while in SGLD, low frequency motion gain more kinetic energy, becoming hotter, and high frequency motions lose some kinetic energy, becoming cooler. It is the low frequency motion that controls the conformational searching, therefore, SGLD achieves enhanced conformational searching ability without raising temperature.

Figure 2

A basic scheme describing the replica exchanging self-guided Langevin dynamics simulation. There are multiple stages with different simulation conditions, T⁽ⁱ⁾and $T_{\mathrm{SG}}^{\left(i\right)}$. The base stage has the simulation condition of interest, T⁽⁰⁾and $T_{\mathrm{SG}}^{\left(0\right)}=T^{\left(0\right)}$. A simulation system is replicated to many copies, called replicas. On each stage there are one or more replicas. Between stages, a pair of randomly chosen replicas are exchanged according to the exchanging probability. A TRXLD simulation has different stage temperatures, T⁽ⁱ⁾ ⩾ T⁽⁰⁾, but no guiding force, $T_{\mathrm{SG}}^{\left(i\right)}=T^{\left(0\right)}$, and a RXSGLD simulation has difference self-guiding temperatures, $T_{\mathrm{SG}}^{\left(i\right)}\ge T^{\left(0\right)}$.

Figure 3

The skewed double well potentials in the y-dimension. The potential function is shown in Eq. 27. The skew parameter, s, controls the relative depths of the two wells. The energy barrier between the two wells is about 10 kT₀.

Figure 4

Average stage energies in the TRXLD and RXSGLD simulations. The solutions from Eq. 29a are shown as dashed lines. The TRXLD simulations has 8 stages with T = 50/100 K, and the RXSGLD simulations have 8 stages with T_SG = 50/100 K and T = T₀ = 50 K.

Figure 5

The y-coordinate distributions at stage 0, 4, and 7 in the TRXLD (T = 50/100 K), RXSGLD (T_SG = 50/100 K and T = T₀ = 50 K) simulations. The data is for the skewed double well system with s = 2kT₀.

Figure 6

The x-z distributions at stage 0, 4, and 7 in the TRXLD (T = 50/100 K), RXSGLD (T_SG = 50/100Kand T = T₀ = 50K) simulations. The data is for the skewed double well system with s = 2kT₀.

Figure 7

Distribution fractions of well 1, x₁, on the base stages (T₀ = 50 K) from the TRXLD (T = 50/100 K) and RXSGLD (T_SG = 50/100 K and T = T₀ = 50 K) simulations. The solutions of Eq. 30 are shown as dashed lines.

Figure 8

The distribution of the 16 backbone dihedral angles in the 9-reside β-hairpin folding peptide calculated from the replica exchange simulations. One to three regions are defined as labeled for each dihedral angle for the subset index clustering.

Figure 9

Total numbers of clusters searched in the TRXLD (T = 274/500 K) and RXSGLD (T_SG = 274/400 K and T = T₀ = 274 K) simulations. The 9-residue β-hairpin folding peptide is simulated with an implicit salvation model.

Figure 10

Comparison of the cluster populations between the TRXLD (T = 274/400 K) result and the RXSGLD (T_SG = 274/400 K and T = T₀ = 274 K) result. For convenience in plotting in the logarithm scale, the cluster population counts start from 1. The 9-residue β-hairpin folding peptide is simulated with an implicit salvation model. The subset indexing clustering method is described in the text.

Figure 11

Acceptance ratios at each stage in the TRXLD and RXSGLD simulations. The 9-residue β-hairpin folding peptide is simulated with explicit water.

Figure 12

Stage identities visited by replica 0 in the TRXLD and RXSGLD simulations. The 9-residue β-hairpin folding peptide is simulated with explicit water.

Figure 13

Potential energy distributions at each stage in the TRXLD (T = 274/400 K) and RXSGLD (T_SG = 274/400 K and T = T₀ = 274 K) simulations. The 9-residue β-hairpin folding peptide is simulated with explicit water.

Figure 14

The total numbers of clusters searched during the TRXLD and RXSGLD simulations. The 9-residue β-hairpin folding peptide is simulated with explicit water.