Implied timescales¶

Frank Noe frank.noe@fu-berlin.de
FU Berlin, Arnimallee 6, 14195 Berlin, Germany

Summary: Implied timescales refers to the relaxation timescales of a molecule implied by a Markov model transition matrix estimated at a lag time $\tau$. Since $\tau$ is a model parameter but the relaxation timescales is a physical property of the simulated system, it is expected that the implied timescales are independent of $\tau$. This is not the case for very small values of $\tau$, where the eigenvector approximation error and the spectral error of the Markov model are large; nor it is the case for very large values of $\tau$, where the estimation is dominated by numerical errors. For a good state space discretization, however, an intermediate range of $\tau$ values can be found in which the largest implied timescales should be approximately constant. Computing the implied timescales and plotting them is thus a useful tool to decide for a lag time $\tau$ at which the Markov model will be estimated. Testing whether the implied timescales are constant over a range of $\tau$ within statistical error is also a weak version of making a Chapman-Kolmogorow test.

One of the most interesting kinetic properties of molecular systems are the relaxation timescales of the system. They can be both experimentally accessed via relaxation or correlation functions that are measurable with various spectroscopic techniques [3][1][7][8], but can also be directly calculated from the Markov model eigenvalues as implied timescales: $$t_{i}(k\tau_{0})=-\frac{k\tau_{0}}{\ln\lambda_{i}(k\tau_{0})}. \:\:\:\:(1)$$ Where $\tau_{0}$ is some arbitrary reference lag time and $\tau=k\tau_{0}$ is the lagtime used. [14] has argued: if dynamics are exactly Markovian, then $\lambda_{i}(k\tau)=\lambda_{i}(\tau)^{k}$, and hence: $$t_{i}(k\tau_{0})=t_{i}=-\frac{\tau_{0}}{\ln\lambda_{i}(\tau_{0})}. \:\:\:\:(2)$$ must be independent of the lagtime chosen. Ref. [14] has therefore suggested to test if the slowest implied timescales, $t_{i}(k\tau)$, are approximately constant in the lagtime $k\tau$. Note that this is not a sufficient test of Markovianity, because it does not test whether the eigenvectors are also independent of $k\tau_{0}$. However, it has been empirically observed that constant timescales are a very strong indicator of having an MSM which approximates the underlying dynamics well. Hence, implied timescale tests are frequently used in order to select between different discretizations (clusterings) of state space.

If the discretized dynamics are not exactly Markovian, how do the implied timescales behave in $k\tau_{0}$? In [11], we have derived a tight bound for the slowest relaxation rate $\kappa_{2}=t_{2}^{-1}$. From this bound follows that the relative error of the largest implied timescale is bounded by: $$\begin{aligned} \frac{\hat{t}_{2}-t_{2}}{t_{2}} & \le & \frac{\ln\frac{1}{\alpha}}{\frac{\tau}{t_{2}}+\ln\frac{1}{\alpha}} \:\:\:\:(3)\end{aligned}$$ Where $\alpha=\langle\psi_{2},\hat{\psi}_{2}\rangle_{\mu}$ is the discretization quality with respect to the second propagator eigenfunction. In simple words, if $\alpha=1$, the state space discretization resolves the slowest process perfectly, while if $\alpha=0$, the slowest process is completely concealed by the discretization. For the limit $\frac{\tau}{t_{2}}\gg\ln\frac{1}{\alpha}$, i.e. when either the state space discretization is very good, or the lagtime is very large, the error becomes $$\begin{aligned} \frac{\hat{t}_{2}-t_{2}}{t_{2}} & \lessapprox & \frac{t_{2}}{\tau}\ln\frac{1}{\alpha}. \:\:\:\:(4)\end{aligned}$$ We observe that: (1) the implied timescales are well approximated if the state space discretization is very good ($\alpha=1\rightarrow\ln\alpha^{-1}=0$), (2) the implied timescale estimate converges towards the true implied timescale when the lagtime $\tau$ is increased. Unfortunately this convergence is very slow, with $\tau^{-1}$. Following [6], we can derive that also for all other relaxation processes, the corresponding implied timescale converges to its true value when the discretization quality increases or the lagtime $\tau$ increases:

$$\lim_{\tau\rightarrow\infty}\left|t_{j}(\tau)-\hat{t}_{j}(\tau)\right|=0, \:\:\:\:(5)$$

and also

$$\lim_{\delta_{j}\rightarrow0}\left|t_{j}(\tau)-\hat{t}_{j}(\tau)\right|=0. \:\:\:\:(6)$$

where $\delta_{j}=1-\alpha_{j}$ is the projection error of the state space discretization with respect to the $j$th dynamical process. This fact has been empirically observed in many previous studies [14][13][4][5][2][9][10].

From the mathematical facts above, the following rationale to assess the quality of the state space discretization follows:

For a given state space discretization, estimate a series of transition matrices $\mathbf{T}(\tau_{k}=k\Delta t)$, where $\Delta t$ is the time step between saved trajectory frames and $k$ is a variable integer, using the methods described in lecture_ml_nonrev and lecture_ml_rev.
Compute the $m$ largest eigenvalues of $\mathbf{T}(\tau_{k})$, and from these the $m$ slowest implied timescales $t_{i}(\tau_{k})$ depending on lag time $\tau_{k}$.
When the implied timescales $t_{i}$ reach an approximately constant value for lagtimes $\tau_{k}$, the state space discretization is sufficiently good to resolve the dynamics in these slowest processes. Usually, it is also expected that the lagtimes for which this approximate constant value is reached are significantly smaller than the timescales $t_{i}$ of interest.
Select the minimal lagtime $\tau$ at which $t_{i}$ are approximately constant, and use $\mathbf{T}(\tau)$ as Markov model.

Two notes of caution must be made at this point: (1) The argument above does not include the effect of statistics and is thus strictly only valid for the limit of good sampling. In many practical cases, statistics are insufficient and the implied timescales do not show the expected monotonous behavior that permits the quality of the discretization to be assessed. In this case, additional sampling is needed. (2) Observing convergence of the slowest implied timescales in $\tau$ is not a strict test of Markovianity. While Markovian dynamics implies constancy of implied timescales in $\tau$ [9][13], the reverse is not true and would require the eigenvectors to be constant as well. However, observing the lag time-dependence of the implied timescales is a useful approach to asses the quality of the discretization, and to choose a lag time $\tau$ at which $\mathbf{T}(\tau)$ shall be estimated, but this model needs to be validated subsequently (see lecture_chapman_kolmogorow).

Figure 1: Illustration of the eigenfunction approximation error $\delta_{2}$ on the slow transition in the diffusion in a double well (top, black line). The slowest eigenfunction is shown in the lower four panels (black), along with the step approximations (green) of the partitions (vertical black lines) at $x=50$; $x=40$; $x=10,20,...,80,90$; and $x=40,45,50,55,60$. The eigenfunction approximation error $\delta_{2}$ is shown as red area and its norm is printed.

Figure 2: Illustration of the eigenfunction approximation errors $\delta_{2}$ and $\delta_{3}$ on the two slowest processes in a two-dimensional three-well diffusion model (see Supplementary Information for model details). The columns from left to right show different state space discretizations with white lines as state boundaries: (i) 3 states with maximum metastability, (ii) the metastable states were further subdivided manually to better resolve the transition region, resulting in a partition where no individual state is metastable, (iii)/(iv) Voronoi partition using 25/100 randomly chosen centers, respectively. (a) Potential, (b) The exact eigenfunctions of the slow processes, $\psi_{2}(\mathbf{x})$ and $\psi_{3}(\mathbf{x})$, (c) The approximation of eigenfunctions with discrete states, $Q\psi_{2}(\mathbf{x})$ and $Q\psi_{3}(\mathbf{x})$, (d) Approximation errors $|\psi_{2}-Q\psi_{2}|(\mathbf{x})$ and $|\psi_{3}-Q\psi_{3}|(\mathbf{x})$. The error norms $\delta_{2}$ and $\delta_{3}$ are given. As illustrative examples, consider the one- and two-dimensional diffusion processes and different state space discretization shown in Figs. 3 and 1. Fig. 3 shows the slowest implied timescale $t_{2}$ for the diffusion in a two-well potential (see Fig. 1) with discretizations shown in Fig. 1. The two-state partition at $x=50$ requires a lag time of $\approx2000$ steps in order to reach an error of $<3\%$ with respect to the true implied timescale, which is somewhat slower than $t_{2}$ itself. When the two-state partition is distorted by shifting the discretization border to $x=40$, this quality is not reached before the process itself has relaxed. Thus, in this system two states are not sufficient to build a Markov model that is at the same time precise and has a time resolution good enough to trace the decay of the slowest process. By using more states and particularly a finer discretization of the transition region, the same approximation quality is obtained with only $\tau\approx1500$ (blue) and $\tau\approx500$ steps (green).

Figure 3: Convergence of the slowest implied timescale $t_{2}=-\tau/\ln\lambda_{2}(\tau)$ of the diffusion in a double-well potential depending on the MSM discretization. The metastable partition (black, solid) has greater error than non-metastable partitions (blue, green) with more states that better trace the change of the slow eigenfunction near the transition state. Fig. 4 shows the two slowest implied timescales $t_{2}$, $t_{3}$ for the diffusion in a two-dimensional three-well potential with discretizations shown in Fig. 2a. The metastable 3-state partition requires a lag time of $\approx1000$ steps in order to reach an error of $<3\%$ with respect to the true implied timescale, which is somewhat shorter than the slow but longer than the fast timescale, while refining the discretization near the transition states achieves the same precision with $\tau\approx200$ using only 12 states. A k-means clustering with $k=25$ is worse than the metastable partition, as some clusters cross over the transition region and fail to resolve the slow eigenfunctions. Increasing the number of clusters to $k=100$ improves the result significantly, but is still worse than the 12 states that have been manually chosen so as to well resolve the transition region. This suggests that excellent MSMs could be built with rather few states when an adaptive algorithm that more finely partitions the transition region is employed.

Figure 4: Implied timescales for the two slowest processes in the two-dimensional three-well diffusion model (see 2a for potential and Supplementary Information for details). The colors black, red, yellow, green correspond to the four choices of discrete states shown in columns 1 to 4 of Fig. 2. A fine discretization of the transition region clearly gives the best approximation to the timescales at small lag times.

Acknowledgements:¶

The major part of this article has been published in [12].

Citing implied timescales:¶

The idea of using implied timescales to check for Markovianity was originally introduced in [13].

Bibliography¶

[1]: Bieri, Oliver, Wirz, Jakob, Hellrung, Bruno, Schutkowski, Mike, Drewello, Mario and Kiefhaber, Thomas: The speed limit for protein folding measured by triplet-triplet energy transfer. Proc. Natl. Acad. Sci. USA 96, 9597-9601 (1999).

[2]: N. V. Buchete and G. Hummer: Coarse Master Equations for Peptide Folding Dynamics. J. Phys. Chem. B 112, 6057-6069 (2008).

[3]: Chan, Chi-Kin, Hu, Yi, Takahashi, Satoshi, Rousseau, Denis L., Eaton, William A. and Hofrichter, James: Submillisecond protein folding kinetics studied by ultrarapid mixing. Proc. Natl. Acad. Sci. USA 94, 1779-1784 (1997).

[4]: Chodera, J. D., Dill, K. A., Singhal, N., Pande, V. S., Swope, W. C. and Pitera, J. W.: Automatic discovery of metastable states for the construction of Markov models of macromolecular conformational dynamics. J. Chem. Phys. 126, 155101 (2007).

[5]: J. D. Chodera and F. Noé: Probability distributions of molecular observables computed from Markov models. II: Uncertainties in observables and their time-evolution. J. Chem. Phys. 133, 105102 (2010).

[6]: N. Djurdjevac, M. Sarich and C. Schütte: Estimating the eigenvalue error of Markov State Models. Multiscale Model. Simul. 10, 61-81 (2012).

[7]: Jäger, Marcus, Nguyen, Houbi, Crane, Jason C., Kelly, Jeffery W. and Gruebele, Martin: The folding mechanism of a beta-sheet: the WW domain. J. Mol. Biol. 311, 373-393 (2001).

[8]: Neuweiler, Hannes, Doose, Sören and Sauer, Markus: A microscopic view of miniprotein folding: Enhanced folding efficiency through formation of an intermediate. Proc. Natl. Acad. Sci. USA 102, 16650-16655 (2005).

[9]: Noé, F., Horenko, I., Schütte, C. and Smith, J. C.: Hierarchical Analysis of Conformational Dynamics in Biomolecules: Transition Networks of Metastable States. J. Chem. Phys. 126, 155102 (2007).

[10]: F. Noé, C. Schütte, E. Vanden-Eijnden, L. Reich and T.R. Weikl: Constructing the Full Ensemble of Folding Pathways from Short Off-Equilibrium Simulations. Proc. Natl. Acad. Sci. USA 106, 19011-19016 (2009).

[11]: J.-H. Prinz, J. D. Chodera and F. Noé: Spectral rate theory for two-state kinetics. Phys. Rev. X 4, 011020 (2014).

[12]: J.-H. Prinz, H. Wu, M. Sarich, B. Keller, M. Senne, M. Held, J. D. Chodera, C. Schütte and F. Noé: Markov models of molecular kinetics: Generation and Validation. J. Chem. Phys. 134, 174105 (2011).

[13]: Swope, W. C., Pitera, J. W. and Suits, F.: Describing protein folding kinetics by molecular dynamics simulations: 1. Theory. J. Phys. Chem. B 108, 6571-6581 (2004).

[14]: Swope, W. C., Pitera, J. W., Suits, F., Pitman, M. and Eleftheriou, M.: Describing protein folding kinetics by molecular dynamics simulations: 2. Example applications to alanine dipeptide and beta-hairpin peptide. Journal of Physical Chemistry B 108, 6582-6594 (2004).