Hierarchical Nyström methods for constructing Markov state models for conformational dynamics

TLDR: In this paper, the Nystrom method is used to deal with poorly sampled micro-states, which can help spectral clustering identify metastable aggregates with highly populated micro states rather than being distracted by lowly populated states.

Abstract: Markov state models (MSMs) have become a popular approach for investigating the conformational dynamics of proteins and other biomolecules. MSMs are typically built from numerous molecular dynamics simulations by dividing the sampled configurations into a large number of microstates based on geometric criteria. The resulting microstate model can then be coarse-grained into a more understandable macrostate model by lumping together rapidly mixing microstates into larger, metastable aggregates. However, finite sampling often results in the creation of many poorly sampled microstates. During coarse-graining, these states are mistakenly identified as being kinetically important because transitions to/from them appear to be slow. In this paper, we propose a formalism based on an algebraic principle for matrix approximation, i.e., the Nystrom method, to deal with such poorly sampled microstates. Our scheme builds a hierarchy of microstates from high to low populations and progressively applies spectral clustering on sets of microstates within each level of the hierarchy. It helps spectral clustering identify metastable aggregates with highly populated microstates rather than being distracted by lowly populated states. We demonstrate the ability of this algorithm to discover the major metastable states on two model systems, the alanine dipeptide and trpzip2 peptide.

Figures

FIG. 1. Schematic figure for the multilevel analysis of a free energy landscape. (a) a 1D free energy landscape divided into four levels; (b) a cluster tree representing this free energy landscape. At level one, two nodes are formed that correspond to the two deepest free energy minima. At level two, four nodes are identified for four free energy minima. At level three and four, the number of nodes are reduced since some free energy minima are connected.

FIG. 9. Illustration of the overlapping of microstates (or protein conformations) assigned to large macrostates in different lumpings. The optimal lumping (A) is compared with a representative lumping (B) with mutual information at around 2.1. Both of these lumpings contain 9 large macrostates states with population >2%. The joint probability matrix PA, B(i, j) indicates the overlapping of microstates assigned to macrostate i in the optimal lumping A and macrostate j in the representative lumping B. PA, B has large diagonal elements but small off-diagonal elements after permutation. These results indicate that the large macrostates in the two lumpings share a relatively large fraction of identical microstates.

FIG. 8. Histogram of pairwise mutual information (with a mean value of 2.104) between the optimal lumping and all other lumpings (259 of them) obtained from the Nyström method by varying the level sets (red). For comparison, the mutual information between a lumping generated by PCCA (with 13 macrostates) and the optimal lumping is only 0.785 (green). The entropy of the optimal lumping (upper limit of the mutual information) is 3.313 (black), while the averaged mutual information between random lumpings (we produced 259 random lumpings) and the optimal lumping is 0.051.

FIG. 6. (a) An illustration of block diagonal structures of the microstate transition probability matrices. There are 2000 microstates in total, and the matrices are permuted to group microstates that belong to the same macrostate together. The results show that PCCA (left panel) tends to separate nearly disconnected small blocks first, while HNEG (right panel) focuses on identifying the well populated macrostates. The number of macrostates for both lumpings is 11. (b) HNEG successfully identify the large macrostates (around 9) with a much smaller total number of macrostates (<20) compared to PCCA (>80).

FIG. 7. The Hierarchical Nyström method can robustly identify the large metastable macrostates with population greater than 1% (red), 2% (green), and 3% (blue), when varying the fraction of data that are included in the submatrix A (see Sec. II for details). The percentage of data we include(Pm) in the submatrix A is varied from 41% to 99%. The number of large macrostates keeps the same after we include 50% of the data or more.

FIG. 3. Bayesian comparison of MSMs constructed by the Hierarchical Nyström Extension Graph (HNEG) and the PCCA method for the alanine dipeptide. The y axis displays the logarithm of the posterior probability (ln (P(L1|D))) for models generated by HNEG (red) and PCCA (black). The logarithmic Bayes factor lnB= ln (P(L1|D))HNEG − ln (P(L2|D)PCCA) 100, indicating that HNEG consistently provides a better lumping than PCCA. The lag time is 9 ps.