Search

The ability to make sense of the massive amounts of high-dimensional data generated from molecular dynamics simulations is heavily dependent on the knowledge of a low-dimensional manifold (parameterized by a reaction coordinate or RC) that typically distinguishes between relevant metastable states, and which captures the relevant slow dynamics of interest. Methods based on machine learning and artificial intelligence have been proposed over the years to deal with learning such low-dimensional manifolds, but they are often criticized for a disconnect from more traditional and physically interpretable approaches. To deal with such concerns, in this work we propose a deep learning based state predictive information bottleneck approach to learn the RC from high-dimensional molecular simulation trajectories. We demonstrate analytically and numerically how the RC learnt in this approach is connected to the committor in chemical physics and can be used to accurately identify transition states. A crucial hyperparameter in this approach is the time delay or how far into the future the algorithm should make predictions about. Through careful comparisons for benchmark systems, we demonstrate that this hyperparameter choice gives useful control over how coarse-grained we want the metastable state classification of the system to be. We thus believe that this work represents a step forward in systematic application of deep learning based ideas to molecular simulations. [ABSTRACT FROM AUTHOR]

Artificial intelligence (AI)-based approaches have had indubitable impact across the sciences through the ability to extract relevant information from raw data. Recently, AI has also found use in enhancing the efficiency of molecular simulations, wherein AI derived slow modes are used to accelerate the simulation in targeted ways. However, while typical fields where AI is used are characterized by a plethora of data, molecular simulations, per construction, suffer from limited sampling and thus limited data. As such, the use of AI in molecular simulations can suffer from a dangerous situation where the AI-optimization could get stuck in spurious regimes, leading to incorrect characterization of the reaction coordinate (RC) for the problem at hand. When such an incorrect RC is then used to perform additional simulations, one could start to deviate progressively from the ground truth. To deal with this problem of spurious AI-solutions, here, we report a novel and automated algorithm using ideas from statistical mechanics. It is based on the notion that a more reliable AI-solution will be one that maximizes the timescale separation between slow and fast processes. To learn this timescale separation even from limited data, we use a maximum caliber-based framework. We show the applicability of this automatic protocol for three classic benchmark problems, namely, the conformational dynamics of a model peptide, ligand-unbinding from a protein, and folding/unfolding energy landscape of the C-terminal domain of protein G. We believe that our work will lead to increased and robust use of trustworthy AI in molecular simulations of complex systems. [ABSTRACT FROM AUTHOR]

In this work, we revisit our recent iterative machine learning (ML)–molecular dynamics (MD) technique "Reweighted autoencoded variational Bayes for enhanced sampling" [J. M. L. Ribeiro et al., J. Chem. Phys. 149, 072301 (2018) and Y. Wang, J. M. L. Ribeiro, and P. Tiwary, Nat. Commun. 10, 3573 (2019)] and analyze and formalize some of its approximations. These include (a) the choice of a predictive time-delay, or how far into the future should the ML try to predict the state of a given system output from MD, and (b) that for short time-delays, how much of an error is made in approximating the biased propagator for the dynamics as the unbiased propagator. We demonstrate through a master equation framework as to why the exact choice of time-delay is irrelevant as long as a small non-zero value is adopted. We also derive a correction to reweight the biased propagator, and somewhat to our dissatisfaction but also to our reassurance, we find that it barely makes a difference to the intuitive picture we had previously derived and used. [ABSTRACT FROM AUTHOR]

In this work, we revisit the classic problem of homogeneous nucleation of a liquid droplet in a supersaturated vapor phase. We consider this at different extents of the driving force, or equivalently the supersaturation, and calculate a reaction coordinate (RC) for nucleation as the driving force is varied. The RC is constructed as a linear combination of three order parameters, where one accounts for the number of liquidlike atoms and the other two for local density fluctuations. The RC is calculated from biased and unbiased molecular dynamics (MD) simulations using the spectral gap optimization approach "SGOOP" [P. Tiwary and B. J. Berne, Proc. Natl. Acad. Sci. U. S. A. 113, 2839 (2016)]. Our key finding is that as the supersaturation decreases, the RC ceases to simply be the number of liquidlike atoms, and instead, it becomes important to explicitly consider local density fluctuations that correlate with shape and density variations in the nucleus. All three order parameters are found to have similar barriers in their respective potentials of mean force; however, as the supersaturation decreases, the density fluctuations decorrelate slower and thus carry longer memory. Thus, at lower supersaturations, density fluctuations are non-Markovian and cannot be simply ignored from the RC by virtue of being noise. Finally, we use this optimized RC to calculate nucleation rates in the infrequent metadynamics framework and show that it leads to a more accurate estimate of the nucleation rate with four orders of magnitude acceleration relative to unbiased MD. [ABSTRACT FROM AUTHOR]

Spectral gap optimization of order parameters (SGOOP) [P. Tiwary and B. J. Berne, Proc. Natl. Acad. Sci. U. S. A. 113, 2839 (2016)] is a method for constructing the reaction coordinate (RC) in molecular systems, especially when they are plagued with hard to sample rare events, given a larger dictionary of order parameters or basis functions and limited static and dynamic information about the system. In its original formulation, SGOOP is designed to construct a 1-dimensional RC. Here we extend its scope by introducing a simple but powerful extension based on the notion of conditional probability factorization where known features are effectively washed out to learn additional and possibly hidden features of the energy landscape. We show how SGOOP can be used to proceed in a sequential and bottom-up manner to (i) systematically probe the need for extending the dimensionality of the RC and (ii) if such a need is identified, learn additional coordinates of the RC in a computationally efficient manner. We formulate the method and demonstrate its utility through three illustrative examples, including the challenging and important problem of calculating the kinetics of benzene unbinding from the protein T4L99A lysozyme, where we obtain excellent agreement in terms of dissociation pathway and kinetics with other sampling methods and experiments. In this last case, starting from a larger dictionary of 11 order parameters that are generic for ligand unbinding processes, we demonstrate how to automatically learn a 2-dimensional RC, which we then use in the infrequent metadynamics protocol to obtain 16 independent unbinding trajectories. We believe our method will be a big step in increasing the utility of SGOOP in performing intuition-free sampling of complex systems. Finally, we believe that the utility of our protocol is amplified by its applicability to not just SGOOP but also other generic methods for constructing the RC. [ABSTRACT FROM AUTHOR]

Here we propose the reweighted autoencoded variational Bayes for enhanced sampling (RAVE) method, a new iterative scheme that uses the deep learning framework of variational autoencoders to enhance sampling in molecular simulations. RAVE involves iterations between molecular simulations and deep learning in order to produce an increasingly accurate probability distribution along a low-dimensional latent space that captures the key features of the molecular simulation trajectory. Using the Kullback-Leibler divergence between this latent space distribution and the distribution of various trial reaction coordinates sampled from the molecular simulation, RAVE determines an optimum, yet nonetheless physically interpretable, reaction coordinate and optimum probability distribution. Both then directly serve as the biasing protocol for a new biased simulation, which is once again fed into the deep learning module with appropriate weights accounting for the bias, the procedure continuing until estimates of desirable thermodynamic observables are converged. Unlike recent methods using deep learning for enhanced sampling purposes, RAVE stands out in that (a) it naturally produces a physically interpretable reaction coordinate, (b) is independent of existing enhanced sampling protocols to enhance the fluctuations along the latent space identified via deep learning, and (c) it provides the ability to easily filter out spurious solutions learned by the deep learning procedure. The usefulness and reliability of RAVE is demonstrated by applying it to model potentials of increasing complexity, including computation of the binding free energy profile for a hydrophobic ligand–substrate system in explicit water with dissociation time of more than 3 min, in computer time at least twenty times less than that needed for umbrella sampling or metadynamics. [ABSTRACT FROM AUTHOR]

The ability to predict accurate thermodynamic and kinetic properties in biomolecular systems is of both scientific and practical utility. While both remain very difficult, predictions of kinetics are particularly difficult because rates, in contrast to free energies, depend on the route taken. For this reason, specific enhanced sampling methods are needed to calculate long-time scale kinetics. It has recently been demonstrated that it is possible to recover kinetics through the so-called “infrequent metadynamics” simulations, where the simulations are biased in a way that minimally corrupts the dynamics of moving between metastable states. This method, however, requires the bias to be added slowly, thus hampering applications to processes with only modest separations of time scales. Here we present a frequency-adaptive strategy which bridges normal and infrequent metadynamics. We show that this strategy can improve the precision and accuracy of rate calculations at fixed computational cost and should be able to extend rate calculations for much slower kinetic processes. [ABSTRACT FROM AUTHOR]

We consider a range of model potentials with metastable states undergoing molecular dynamics coupled to a thermal bath in the high friction regime and consider how the optimal reaction coordinate depends on the diffusion anisotropy. For this we use our recently proposed method "spectral gap optimization of order parameters (SGOOP)" [P. Tiwary and B. J. Berne, Proc. Natl. Acad. Sci. U. S. A. 113, 2839 (2016)]. We show how available information about dynamical observables in addition to static information can be incorporated intoSGOOP, which can then be used to accurately determine the "best" reaction coordinate for arbitrary anisotropies.We compare our results with transmission coefficient calculations and published benchmarks wherever applicable or available, respectively. [ABSTRACT FROM AUTHOR]

Condensation of a liquid droplet from a supersaturated vapour phase is initiated by a prototypical nucleation event. As such it is challenging to compute its rate from atomistic molecular dynamics simulations. In fact at realistic supersaturation conditions condensation occurs on time scales that far exceed what can be reached with conventional molecular dynamics methods. Another known problem in this context is the distortion of the free energy profile associated to nucleation due to the small, finite size of typical simulation boxes. In this work the problem of time scale is addressed with a recently developed enhanced sampling method while contextually correcting for finite size effects. We demonstrate our approach by studying the condensation of argon, and showing that characteristic nucleation times of the order of magnitude of hours can be reliably calculated. Nucleation rates spanning a range of 10 orders of magnitude are computed at moderate supersaturation levels, thus bridging the gap between what standard molecular dynamics simulations can do and real physical systems. [ABSTRACT FROM AUTHOR]

We use a recently proposed method called Spectral Gap Optimization of Order Parameters (SGOOP) [P. Tiwary and B. J. Berne, Proc. Natl. Acad. Sci. U. S. A. 113, 2839 (2016)], to determine an optimal 1-dimensional reaction coordinate (RC) for the unbinding of a bucky-ball from a pocket in explicit water. This RC is estimated as a linear combination of the multiple available order parameters that collectively can be used to distinguish the various stable states relevant for unbinding. We pay special attention to determining and quantifying the degree to which water molecules should be included in the RC. Using SGOOP with under-sampled biased simulations, we predict that water plays a distinct role in the reaction coordinate for unbinding in the case when the ligand is sterically constrained to move along an axis of symmetry. This prediction is validated through extensive calculations of the unbinding times through metadynamics and by comparison through detailed balance with unbiased molecular dynamics estimate of the binding time. However when the steric constraint is removed, we find that the role of water in the reaction coordinate diminishes. Here instead SGOOP identifies a good one-dimensional RC involving various motional degrees of freedom. [ABSTRACT FROM AUTHOR]

Metadynamics is a popular enhanced sampling scheme wherein by periodic application of a repulsive bias, one can surmount high free energy barriers and explore complex landscapes. Recently, metadynamics was shown to be mathematically well founded, in the sense that the biasing procedure is guaranteed to converge to the true free energy surface in the long time limit irrespective of the precise choice of biasing parameters. A differential equation governing the post-transient convergence behavior of metadynamics was also derived. In this short communication, we revisit this differential equation, expressing it in a convenient and elegant Riccati-like form. A perturbative solution scheme is then developed for solving this differential equation, which is valid for any generic biasing kernel. The solution clearly demonstrates the robustness of metadynamics to choice of biasing parameters and gives further confidence in the widely used method. [ABSTRACT FROM AUTHOR]

We consider the rate of transition for a particle between two metastable states coupled to a thermal environment for various magnitudes of the coupling strength using the recently proposed infrequent metadynamics approach [P. Tiwary and M. Parrinello, Phys. Rev. Lett. 111, 230602 (2013)]. We are interested in understanding how this approach for obtaining rate constants performs as the dynamics regime changes from energy diffusion to spatial diffusion. Reassuringly, we find that the approach works remarkably well for various coupling strengths in the strong coupling regime, and to some extent even in the weak coupling regime. [ABSTRACT FROM AUTHOR]