Constant pH replica exchange molecular dynamics in biomolecules using a discrete protonation model

Abstract

A constant pH replica exchange molecular dynamics (REMD) method is proposed and implemented to improve coupled protonation and conformational state sampling. By mixing conformational sampling at constant pH (with discrete protonation states) with a temperature ladder, this method avoids conformational trapping. Our method was tested and applied to seven different biological systems. The constant pH REMD not only predicted pKa correctly for small, model compounds but also converged faster than constant pH molecular dynamics (MD). We further tested our constant pH REMD on a heptapeptide from ovomucoid third domain (OMTKY3). Although constant pH REMD and MD produced very close pKa values, the constant pH REMD showed its advantage in the efficiency of conformational and protonation state samplings.

Figure 1

Diagrams displaying exchanging algorithms in constant-pH REMD. (A): Only molecular structures (denoted as q) are attempted to exchange. In this case, protonation states are not toughed at an exchange attempt; (B): Both molecular structures (denoted as q) and protonation states (denoted as n) are attempted to exchange at the same time. Metropolis criterion is applied in both algorithms to evaluate transitions.

Figure 2

Titration curves of blocked aspartate amino acid from 100 ns MD at 300K and REMD runs. Agreement can be seen between MD and REMD simulations.

Figure 3

Cumulative average protonation fraction of a titratable residue vs Monte Carlo (MC) steps. (A) Aspartic acid reference compound at pH=4. (B) Asp2 in model peptide ADFDA at pH=4.

Figure 4

The titration curves of the model peptide ADFDA at 300K from both MD and REMD simulations. MD simulation time was 100 ns and 10 ns were chosen for each replica for REMD runs.

Figure 5

Backbone dihedral angle (φ, ψ) normalized probability density (Ramachandran plots) for Asp2 at pH 4 in ADFDA. Ramachandran plots at other solution pH values are similar. For Asp2, constant-pH MD and REMD sampled the same local backbone conformational space. Phe3 and Asp4 Ramachandran plots also display the same trend.

Figure 6

Cluster populations of ADFDA at 300K (A): MD vs REMD at pH 4, (B): two REMD runs from different starting structures at pH 4. Large correlation shown in Figure 5B suggests that the REMD runs are converged. Large correlations between two independent REMD runs are also observed at other solution pH values. Correlations between MD and REMD simulations can be found in Table. 3.

Figure 7

(A) and (B) are titration curves of Asp3, Lys5 and Tyr7 in the heptapeptide derived from protein OMTKY3. (C) shows the Hill’s plots of Asp3. The pKa values of Asp3 are found through Hill’s plots.

Figure 7

(A) and (B) are titration curves of Asp3, Lys5 and Tyr7 in the heptapeptide derived from protein OMTKY3. (C) shows the Hill’s plots of Asp3. The pKa values of Asp3 are found through Hill’s plots.

Figure 8

Cumulative average protonation fraction of a titratable residue vs MC steps.

Figure 8

Cumulative average protonation fraction of a titratable residue vs MC steps.

Figure 9

Dihedral angle (φ, ψ) probability densities of Asp3 at pH 4. (A): constant-pH MD results; (B): constant-pH REMD results. All others also show very similar trend.

Figure 10

The root-mean-square deviations (RMSD) between the cumulative (φ, ψ) probability density up to current time and the (φ, ψ) probability density produced by entire simulation. (φ, ψ) probability density convergence behaviors at other pH values also show that REMD runs converge to final distribution faster.

Figure 11

Cluster population at 300 K from constant pH MD and REMD simulations at pH=4. Cluster analysis is performed using the entire simulation. The populations in each cluster from the first and second half of the trajectory are compared and plotted. Ideally, a converged trajectory should yield a correlation coefficient to be 1. (A): constant pH MD (B): constant pH REMD. Much higher correlation coefficient can be seen in constant pH REMD simulation, suggesting much better convergence is achieved by the constant pH REMD run.