An algorithm to find minimum free-energy paths using umbrella integration

Minimum Free Energy Pathways of Reactive Processes with Nudged Elastic Bands

The determination of minimum free energy pathways (MFEP) is one of the most widely used strategies to study reactive processes. For chemical reactions in complex environments, the combination of quantum mechanics (QM) with a molecular mechanics (MM) representation is usually necessary in a hybrid QM/MM framework. However, even within the QM/MM approximation, the affordable sampling of the phase space is, in general, quite restricted. To reduce drastically the computational cost of the simulations, several methods such as umbrella sampling require performing a priori a selection of a reaction coordinate. The quality of the computed results, in an affordable computational time, is intimately related to the reaction coordinate election which is, in general, a nontrivial task. In this work, we provide an approach to model reactive processes in complex environments that does not require the a priori selection of a reaction coordinate. The proposed methodology combines QM/MM simulations with an extrapolation of the nudged elastic bands (NEB) method to the free energy surface (FENEB). We present and apply our own FENEB scheme to optimize MFEP in different reactive processes, using QM/MM frameworks at semiempirical and density functional theory levels. Our implementation is based on performing the FENEB optimization by uncoupling the optimization of the band in a perpendicular and tangential direction. In each step, a full optimization with the spring force is performed, which guarantees that the images remain evenly distributed. The robustness of the method and the influence of sampling on the quality of the optimized MFEP and its associated free energy barrier are studied. We show that the FENEB method provides a good estimation of the reaction barrier even with relatively short simulation times, supporting that its combination with QM/MM frameworks provides an adequate tool to study chemical processes in complex environments.

Curvature-weighted nudged elastic band method using the Riemann curvature

The nudged elastic band (NEB) method is utilized to find reaction paths (RPs) using discretized intermediate structures called "images" between a reactant and a product. In fact, NEB calculations do not always converge because of the bent of RPs. Mathematically, more images are needed for complex curves, and here, we focused on the curvature. In this study, we propose a new method for calculating the curvature of the RPs as well as a method for weighting the spring constant of the NEB with the curvature, which we named the curvature weighted NEB (CW-NEB) method. In addition, we will propose the CW-NEB method with the climbing image (CI) method (CW-CI-NEB). To show the efficiency of our method, calculations for the ene-reaction of the CW-CI-NEB method were performed and compared with those of the CI-NEB methods. The CW-CI-NEB methods converged in shorter iteration times than the CI-NEB calculation, which was found by gathering the bent of RPs.

Multidimensional virtual-system coupled canonical molecular dynamics to compute free-energy landscapes of peptide multimer assembly

An enhanced-sampling method termed multidimensional virtual-system coupled canonical molecular dynamics (mD-VcMD) method is developed. In many cases, generalized-ensemble methods realizing enhanced sampling, for example, adaptive umbrella sampling, apply an effective potential, which is derived from temporarily assumed canonical distribution as a function of one or more arbitrarily defined reaction coordinates. However, it is not straightforward to estimate the appropriate canonical distribution, especially for cases applying multiple reaction coordinates. The current method, mD-VcMD, does not rely on the form of the canonical distribution. Therefore, it is practically useful to explore a high-dimensional reaction-coordinate space. In this article, formulation of mD-VcMD and its evaluation with the simple molecular models consisting of three or four alanine peptides are presented. We confirmed that mD-VcMD efficiently searched 2D and 3D reaction-coordinate spaces defined as interpeptide distances. Direct comparisons with results of long-term canonical MD simulations revealed that mD-VcMD produces correct canonical ensembles. © 2019 Wiley Periodicals, Inc.

Multi-dimensional virtual system introduced to enhance canonical sampling

When an important process of a molecular system occurs via a combination of two or more rare events, which occur almost independently to one another, computational sampling for the important process is difficult. Here, to sample such a process effectively, we developed a new method, named the "multi-dimensional Virtual-system coupled Monte Carlo (multi-dimensional-VcMC)" method, where the system interacts with a virtual system expressed by two or more virtual coordinates. Each virtual coordinate controls sampling along a reaction coordinate. By setting multiple reaction coordinates to be related to the corresponding rare events, sampling of the important process can be enhanced. An advantage of multi-dimensional-VcMC is its simplicity: Namely, the conformation moves widely in the multi-dimensional reaction coordinate space without knowledge of canonical distribution functions of the system. To examine the effectiveness of the algorithm, we introduced a toy model where two molecules (receptor and its ligand) bind and unbind to each other. The receptor has a deep binding pocket, to which the ligand enters for binding. Furthermore, a gate is set at the entrance of the pocket, and the gate is usually closed. Thus, the molecular binding takes place via the two events: ligand approach to the pocket and gate opening. In two-dimensional (2D)-VcMC, the two molecules exhibited repeated binding and unbinding, and an equilibrated distribution was obtained as expected. A conventional canonical simulation, which was 200 times longer than 2D-VcMC, failed in sampling the binding/unbinding effectively. The current method is applicable to various biological systems.

Automated Exploration of Free Energy Landscapes Based on Umbrella Integration

We present a new approach for automated exploration of free energy landscapes on the basis of the umbrella integration (UI) method. The method to search points in the landscape relies on the normal distributions and gradients of the potential of mean force (PMF) obtained from UI calculations. We applied this approach to the alanine dipeptide in solution and demonstrated that the equilibrium and the transition states were efficiently found in the ascending order of the PMF values.

Molecular dynamics coupled with a virtual system for effective conformational sampling

An enhanced conformational sampling method is proposed: virtual-system coupled canonical molecular dynamics (VcMD). Although VcMD enhances sampling along a reaction coordinate, this method is free from estimation of a canonical distribution function along the reaction coordinate. This method introduces a virtual system that does not necessarily obey a physical law. To enhance sampling the virtual system couples with a molecular system to be studied. Resultant snapshots produce a canonical ensemble. This method was applied to a system consisting of two short peptides in an explicit solvent. Conventional molecular dynamics simulation, which is ten times longer than VcMD, was performed along with adaptive umbrella sampling. Free-energy landscapes computed from the three simulations mutually converged well. The VcMD provided quicker association/dissociation motions of peptides than the conventional molecular dynamics did. The VcMD method is applicable to various complicated systems because of its methodological simplicity. © 2018 Wiley Periodicals, Inc.

Toward an Automatic Determination of Enzymatic Reaction Mechanisms and Their Activation Free Energies

We present a combination of the string method and a path collective variable for the exploration of the free energy surface associated to a chemical reaction in condensed environments. The on-the-fly string method is employed to find the minimum free energy paths on a multidimensional free energy surface defined in terms of interatomic distances, which is a convenient selection to study bond forming/breaking processes. Once the paths have been determined, a reaction coordinate is defined as a measure of the advance of the system along these paths. This reaction coordinate can be then used to trace the reaction Potential of Mean Force from which the activation free energy can be obtained. This combination of methodologies has been here applied to the study, by means of Quantum Mechanics/Molecular Mechanics simulations, of the reaction catalyzed by guanidinoacetate methyltransferase. This enzyme catalyzes the methylation of guanidinoacetate by S-adenosyl-l-methionine, a reaction that involves a methyl transfer and a proton transfer and for which different reaction mechanisms have been proposed.

A New Maximum Likelihood Approach for Free Energy Profile Construction from Molecular Simulations

A novel variational method for construction of free energy profiles from molecular simulation data is presented. The variational free energy profile (VFEP) method uses the maximum likelihood principle applied to the global free energy profile based on the entire set of simulation data (e.g from multiple biased simulations) that spans the free energy surface. The new method addresses common obstacles in two major problems usually observed in traditional methods for estimating free energy surfaces: the need for overlap in the re-weighting procedure and the problem of data representation. Test cases demonstrate that VFEP outperforms other methods in terms of the amount and sparsity of the data needed to construct the overall free energy profiles. For typical chemical reactions, only ~5 windows and ~20-35 independent data points per window are sufficient to obtain an overall qualitatively correct free energy profile with sampling errors an order of magnitude smaller than the free energy barrier. The proposed approach thus provides a feasible mechanism to quickly construct the global free energy profile and identify free energy barriers and basins in free energy simulations via a robust, variational procedure that determines an analytic representation of the free energy profile without the requirement of numerically unstable histograms or binning procedures. It can serve as a new framework for biased simulations and is suitable to be used together with other methods to tackle with the free energy estimation problem.