Machine learning for biomolecular simulations
Authors:   , Jung Hendrik 
Published: January 2025

Many biomolecular processes are characterized by transitions between metastable states separated by large energy barriers. Examples of such metastable states are the inward- and outward-facing configuration of a transporter, the bound and unbound configurations of a protein and small drug-like molecule, or the folded and unfolded states of a protein. A good model for the complex dynamics on the rough high-dimensional energy surfaces of biomolecules is the thermally activated barrier crossing, as shown by the agreement between theory and single-molecule experiments of protein folding. In this model, the transitions between states are very short compared to the average waiting time to observe one transition 𝜏 = 𝑘−1, with 𝑘 the rate constant for the transition. The reason is that 𝜏 is proportional to to the exponential of the barrier height, 𝜏 ∝ exp(𝛽Δ𝐺‡), where Δ𝐺‡ is the barrier height and 𝛽 = (𝑘𝐵𝑇 )−1 is the inverse temperature. On the other hand, the time between leaving one state and entering the other, the transition path time 𝑡TP, is approximately proportional to ln(𝛽Δ𝐺‡) if Δ𝐺‡ > 2/𝛽. This disparity of timescales makes most functional transitions in biomolecules rare events, because they happen very rarely compared to the time they take, i.e. 𝑡TP ≪ 𝜏.
Due to the high temporal and spatial resolution inherent to solving Newtons equation of motion for all 𝑁 particles in the system in small time increments, molecular dynamics (MD) simulations are uniquely positioned to resolve the minute details of the molecular rearrangement during the short TPs connecting the metastable states. At the same time, due to sequential nature of MD and the small time increments needed for an accurate integration, it can be computationally challenging to impossible to reach the timescale of 𝜏, which is needed to observe a single transition event. This fact is known as the "sampling challenge".
A possible solution to the disparity of timescales between the waiting time 𝜏 and the transition path time 𝑡TP is therefore to simulate only the transition events. Unfortunately, due to the stochastic nature of the transitions, it is impossible to know in advance when a transition will occur. Importantly, although it is expected that only a few of the 3𝑁 Cartesian degrees of freedom (DoFs) are relevant to describe the transition, these DoFs are not known a priori either. It is therefore not easily possible to solve the sampling challenge by driving the system along these key DoFs to force the occurrence of TPs. The few relevant DoFs together form the reaction coordinate, which-if known-is sufficient to quantitatively describe the reaction progress along a given transition. This fact is known as the "interpretation challenge".
Solutions to the sampling and the interpretation challenges are mutually dependent. If the sampling challenge was solved, one could statistically infer the key DoFs by observing many TPs. Conversely, if the interpretation challenge was solved, one could easily simulate many TPs by driving the system along the key DoFs that control the transition. The aim of this thesis is therefore to develop and describe a machine learning (ML) algorithm that addresses both challenges simultaneously and autonomously.
The proposed ML algorithm combines transition path sampling (TPS), concepts from statistical physics, and deep learning to learn the succinct features of the transition in an iterative and self-consistent manner. TPS is a Markov chain Monte Carlo method in trajectory space to sample the rare transitions between states by focusing the computational effort on the transition region. The ML algorithm autonomously initializes the short MD simulations used in TPS according to its current best guess of the location of the transition region. It then iteratively refines its understanding of the transition dynamics by analyzing the outcomes of the MD simulations and learning to predict the committor with a neural network (NN). The committor towards a given state is the probability of reaching that state before any other. The committor provides the optimal dimensionality reduction for any transition, because it provides a progress coordinate that monotonically increases from one to the other state and because it preserves the dynamics of the system. By learning to predict the location of the transition state ensemble, the algorithm is able to drastically increase the efficiency of the rare event sampling by almost exclusively initializing MD simulations at the transition state ensemble. Furthermore, by comparing its committor prediction with the later observed true outcomes of the simulations, the algorithm learns to predict the transition dynamics in a self-consistent manner. Finally, at convergence, detailed knowledge about the transition dynamics and the key DoFs determining the transition are encoded in the trained NN. This knowledge is extracted by an input relevance analysis of the NN, and can subsequently be used to build low-dimensional mathematical models of the transition using symbolic regression.
The proposed ML algorithm was applied to molecular systems ranging from ion dissociation in water over reorganizations of membrane proteins to the crystallization of polymers and methane hydrates to show that it is applicable to and relevant for a wide range of rare molecular processes. In addition, the proposed ML algorithm was implemented in two open-source python libraries.