Minimum Free Energy Path Research Articles

Nucleation and interfacial adsorption in ternary systems.

Nucleation is studied in incompressible ternary fluids by examining the topology of the overall landscape of the energy surface. Minimum free energy paths for nucleation (MFEPs) of a single nucleus in an infinite matrix are computed with the string method in the framework of the continuum theory of nucleation for the regular solution. Properties of the critical nucleus are compared with the predictions of the classical nucleation theory. MFEPs are found to exhibit complex nucleation pathways with non-monotonic variations of compositions in the interfacial region, specifically adsorption of a component. In the symmetric regular solution, the minority component is found to segregate at the interface during nucleation with a concomitant depletion of the nucleus core, resulting in unpredicted partition of the non-selective component. Despite increasing the gradient energy, such inhomogeneity in composition is shown to lower the nucleation barrier.

Regulation of the Molecular Brake Region in FGFR2 Kinase

Unregulated activity of tyrosine kinases (TKS) is responsible for numerous developmental musculoskeletal diseases and cancers. The regulatory mechanisms of tks and how they are affected by point mutations are still in need of elucidation. Fibroblast growth factor receptor (FGFR) kinase is a TK whose regulatory element, the “molecular brake,” is thought to prevent constitutive activation of the receptor, despite not being near the kinase's active site. germline mutations in this “molecular brake” region in FGFR cause musculoskeletal disorders such as craniosynostosis and dwarfism, while somatic mutations can lead to uninhibited cell growth, causing an array of cancers. furthermore, it is thought that the majority of pathogenic activating mutations in protein kinases are in “molecular brake”-like regions. However, the mechanism through which the “molecular brake” of FGFR kinase is disengaged, allowing for signal transduction, is still unknown. using path-based molecular dynamics (MD) simulations, we have shown that the inward motion of the kinase's activation loop upon tyrosine phosphorylation is correlated with disengagement of the “molecular brake.” results of the string method in collective variables, which finds the minimum free energy path connecting the loop-out and loop-in states, indicate found that inward motion of loop occurs within 1 frame of “molecular brake” disengagement. Furthermore, umbrella sampling shows that in loop-out conformations, the disengaged state of the “molecular brake” is less stable, while in loop-in conformations, the disengaged state becomes more stable. Our results suggest that the “molecular brake” is allosterically controlled by the activation loop motion. An understanding of the “molecular brake” disengagement mechanism will provide avenues for targeted therapies to counteract pathologic over-activation of fgfr kinase leading to uncontrolled signal transduction and uninhibited cell growth.

A fast tomographic method for searching the minimum free energy path.

Minimum Free Energy Path (MFEP) provides a lot of important information about the chemical reactions, like the free energy barrier, the location of the transition state, and the relative stability between reactant and product. With MFEP, one can study the mechanisms of the reaction in an efficient way. Due to a large number of degrees of freedom, searching the MFEP is a very time-consuming process. Here, we present a fast tomographic method to perform the search. Our approach first calculates the free energy surfaces in a sequence of hyperplanes perpendicular to a transition path. Based on an objective function and the free energy gradient, the transition path is optimized in the collective variable space iteratively. Applications of the present method to model systems show that our method is practical. It can be an alternative approach for finding the state-to-state MFEP.

Enzymatic minimum free energy path calculations using swarms of trajectories.

The development of approaches for simulating rare events in complex molecular systems is a central concern in chemical physics. In recent work, Roux and co-workers proposed a novel, swarms of trajectories (SoT) method for determining the transition paths of such events. It consists of the dynamical refinement on the system's free energy surface of a putative transition path that is parametrized in terms of a set of collective variables (CVs) that are identified as being important for the transition. In this work, we have implemented the SoT method and used it to investigate the catalytic mechanisms of two enzymatic reactions using hybrid QM/MM potentials. Our aim has been to test the performance of SoT for enzyme systems and to devise robust simulation protocols that can be employed in future studies of this type. We identify the conditions under which converged results can be obtained using inertial and Brownian dynamical evolutions of the CVs, show that the inclusion of several CVs can give significant additional insight into the mechanisms of the reactions, and show that the use of minimum energy paths as starting guesses can greatly accelerate path refinement.

Defect removal in the course of directed self-assembly is facilitated in the vicinity of the order-disorder transition.

The stability of prototypical defect morphologies in thin films of symmetric diblock copolymers on chemically patterned substrates is investigated by self-consistent field theory. The excess free energy of defects and barriers of defect-removal mechanisms are obtained by computing the minimum free-energy path. Distinct defect-removal mechanisms are illustrated demonstrating that (i) defects will become unstable at a characteristic value of incompatibility χN* above the order-disorder transition and (ii) the kinetics is accelerated at weak segregation. Numerical findings are placed in the context of physical mechanisms, and implications for directed self-assembly are discussed.

Validation of the reaction thermodynamics paths associated with LiK(BH4)2 by detection of metastable reaction paths from first-principles calculations

LiK(BH4)2, with a 10.65 wt% H2 capacity, is a promising hydrogen storage candidate material. We describe all the possible metastable paths for three reactant mixtures associated with LiK(BH4)2, using first-principles thermodynamic calculations. However, with the first-principles calculations inherent disadvantage is absent from the vibrational entropy contribution to the prediction of the reaction mechanism, which indicates that the reaction mechanism of the reaction path prediction may not be the most stable. The study focuses on examining the possible metastable paths for the three reactants, using first-principles calculations. The objective of this study is to find the minimum free energy path among all the possible paths of the three reactants mixture. A clear conclusion is that the minimum free energy path predicted from first principles thermodynamic calculations is the most stable reaction paths at an appropriate H2 pressure range for all cases. An additional aim of paper is to assess whether the predicted mechanisms can be unambiguously described under the DFT uncertainty.

Computing the Free Energy Barriers for Less by Sampling with a Coarse Reference Potential while Retaining Accuracy of the Target Fine Model.

Proposed in this contribution is a protocol for calculating fine-physics (e.g., ab initio QM/MM) free-energy surfaces at a high level of accuracy locally (e.g., only at reactants and at the transition state for computing the activation barrier) from targeted fine-physics sampling and extensive exploratory coarse-physics sampling. The full free-energy surface is still computed but at a lower level of accuracy from coarse-physics sampling. The method is analytically derived in terms of the umbrella sampling and the free-energy perturbation methods which are combined with the thermodynamic cycle and the targeted sampling strategy of the paradynamics approach. The algorithm starts by computing low-accuracy fine-physics free-energy surfaces from the coarse-physics sampling in order to identify the reaction path and to select regions for targeted sampling. Thus, the algorithm does not rely on the coarse-physics minimum free-energy reaction path. Next, segments of high-accuracy free-energy surface are computed locally at selected regions from the targeted fine-physics sampling and are positioned relative to the coarse-physics free-energy shifts. The positioning is done by averaging the free-energy perturbations computed with multistep linear response approximation method. This method is analytically shown to provide results of the thermodynamic integration and the free-energy interpolation methods, while being extremely simple in implementation. Incorporating the metadynamics sampling to the algorithm is also briefly outlined. The application is demonstrated by calculating the B3LYP//6-31G*/MM free-energy barrier for an enzymatic reaction using a semiempirical PM6/MM reference potential. These modifications allow computing the activation free energies at a significantly reduced computational cost but at the same level of accuracy compared to computing full potential of mean force.

Comparison between Mean Forces and Swarms-of-Trajectories String Methods.

The original formulation of the string method in collective variable space is compared with a recent variant called string method with swarms-of-trajectories. The assumptions made in the original method are revisited and the significance of the minimum free energy path (MFEP) is discussed in the context of reactive events. These assumptions are compared to those made in the string method with swarms-of-trajectories, and shown to be equivalent in a certain regime: in particular an expression for the path identified by the swarms-of-trajectories method is given and shown to be closely related to the MFEP. Finally, the algorithmic aspects of both methods are compared.

Conformational Transitions in ATP-Driven Calcium Pump SERCA

Sarco/endoplasmic reticulum Ca2+-ATPase (SERCA) is an integral membrane protein that uses ATP hydrolysis as a source of free energy to pump two calcium ions per ATP molecule from calcium poor cytoplasm of the muscle cell to the calcium rich lumen of the sarcoplasmic reticulum, thereby maintaining a ten thousand fold concentration gradient. Detailed structural studies of the pump under different conditions provided analogues of various intermediates in the reaction cycle and revealed important changes in the tertiary structure of the protein both in the cytoplasmic and in the transmembrane parts. Two major outstanding issues are the pathways of the ions to and from the transmembrane binding sites and a detailed understanding of the large scale conformational changes among various functionally relevant states. We have applied all-atom molecular dynamics (MD) and string method with swarms-of-trajectories to study transition pathways among various experimental structures. The size of the system (291,000) and the time-scale involved in the large scale motions prohibit the use of brute-force unbiased simulations to obtain statistically meaningful information. To circumvent this challenge, the string method with swarms-of-trajectories is used to discover the optimal minimum free energy path between the two end-states. The path, which is represented as a chain of states or images in the space of relevant collective variables, is optimized by iterating two steps: moving each image along the drift calculated from an ensemble (swarm) of short unbiased MD trajectories initiated from the image and a reparametrization procedure that keeps all the images equidistant. We hope to understand the molecular details of the complex motions involved in the conformational transitions by studying the pathways produced by the string method and correlating these observations to available structural and biochemical experiments.

Wild-type and molten globular chorismate mutase achieve comparable catalytic rates using very different enthalpy/entropy compensations

The origin of the catalytic power of enzymes with a meta-stable native state, e.g. molten globular state, is an unsolved challenging issue in biochemistry. To help understand the possible differences between this special class of enzymes and the typical ones, we report here computer simulations of the catalysis of both the well-folded wild-type and the molten globular mutant of chorismate mutase. Using the ab initio quantum mechanical/molecular mechanical minimum free-energy path method, we determined the height of reaction barriers that are in good agreement with experimental measurements. Enzyme-substrate interactions were analyzed in detail to identify factors contributing to catalysis. Computed angular order parameters of backbone N-H bonds and side-chain methyl groups suggested site-specific, non-uniform rigidity changes of the enzymes during catalysis. The change of conformational entropy from the ground state to the transition state revealed distinctly contrasting entropy/enthalpy compensations in the dimeric wild-type enzyme and its molten globular monomeric variant. A unique catalytic strategy was suggested for enzymes that are natively molten globules: some may possess large conformational flexibility to provide strong electrostatic interactions to stabilize the transition state of the substrate and compensate for the entropy loss in the transition state. The equilibrium conformational dynamics in the reactant state were analyzed to quantify their contributions to the structural transitions enzymes needed to reach the transition states. The results suggest that large-scale conformational dynamics make important catalytic contributions to sampling conformational regions in favor of binding the transition state of substrate.

String method for calculation of minimum free-energy paths in Cartesian space in freely-tumbling systems.

The string method is a molecular-simulation technique that aims to calculate the minimum free-energy path of a chemical reaction or conformational transition, in the space of a pre-defined set of reaction coordinates that is typically highly dimensional. Any descriptor may be used as a reaction coordinate, but arguably the Cartesian coordinates of the atoms involved are the most unprejudiced and intuitive choice. Cartesian coordinates, however, present a non-trivial problem, in that they are not invariant to rigid-body molecular rotations and translations, which ideally ought to be unrestricted in the simulations. To overcome this difficulty, we reformulate the framework of the string method to integrate an on-the-fly structural-alignment algorithm. This approach, referred to as SOMA (String method with Optimal Molecular Alignment), enables the use of Cartesian reaction coordinates in freely tumbling molecular systems. In addition, this scheme permits the dissection of the free-energy change along the most probable path into individual atomic contributions, thus revealing the dominant mechanism of the simulated process. This detailed analysis also provides a physically-meaningful criterion to coarse-grain the representation of the path. To demonstrate the accuracy of the method we analyze the isomerization of the alanine dipeptide in vacuum and the chair-to-inverted-chair transition of β-D mannose in explicit water. Notwithstanding the simplicity of these systems, the SOMA approach reveals novel insights into the atomic mechanism of these isomerizations. In both cases, we find that the dynamics and the energetics of these processes are controlled by interactions involving only a handful of atoms in each molecule. Consistent with this result, we show that a coarse-grained SOMA calculation defined in terms of these subsets of atoms yields nearidentical minimum free-energy paths and committor distributions to those obtained via a highly-dimensional string.

Erratum: On‐the‐path random walk sampling for efficient optimization of minimum free energy path

Published in Journal of Computational Chemistry 2009, 30, 1649–1653; DOI: 10.1002/jcc.21311 The authors have notified us of errors in the article. The corrected text and equation are shown below. We apologize for any inconvenience this may have caused. Reference 2 should read, “E, W. N.; Ren, W. Q.; Vanden-Eijnden, E. Phys. Rev. B 2002, 66, 052301”; and Reference 8 should read, ”E, W. N.; Ren. W. Q.; Vanden-Eijnden, E. J. Phys. Chem. 2005, 109, 6688–6693”.

Hidden Conformation Events in DNA Base Extrusions: A Generalized Ensemble Path Optimization and Equilibrium Simulation Study.

DNA base extrusion is a crucial component of many biomolecular processes. Elucidating how bases are selectively extruded from the interiors of double-strand DNAs is pivotal to accurately understanding and efficiently sampling this general type of conformational transitions. In this work, the on-the-path random walk (OTPRW) method, which is the first generalized ensemble sampling scheme designed for finite-temperature-string path optimizations, was improved and applied to obtain the minimum free energy path (MFEP) and the free energy profile of a classical B-DNA major-groove base extrusion pathway. Along the MFEP, an intermediate state and the corresponding transition state were located and characterized. The MFEP result suggests that a base-plane-elongation event rather than the commonly focused base-flipping event is dominant in the transition state formation portion of the pathway; and the energetic penalty at the transition state is mainly introduced by the stretching of the Watson-Crick base pair. Moreover to facilitate the essential base-plane-elongation dynamics, the surrounding environment of the flipped base needs to be intimately involved. Further taking the advantage of the extended-dynamics nature of the OTPRW Hamiltonian, an equilibrium generalized ensemble simulation was performed along the optimized path; and based on the collected samples, several base-flipping (opening) angle collective variables were evaluated. In consistence with the MFEP result, the collective variable analysis result reveals that none of these commonly employed flipping (opening) angles alone can adequately represent the base extrusion pathway, especially in the pre-transition-state portion. As further revealed by the collective variable analysis, the base-pairing partner of the extrusion target undergoes a series of in-plane rotations to facilitate the base-plane-elongation dynamics. A base-plane rotation angle is identified to be a possible reaction coordinate to represent these in-plane rotations. Notably, these in-plane rotation motions may play a pivotal role in determining the base extrusion selectivity.

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.

Discontinuous Bubble Nucleation Due to a Metastable Condensation Transition in Polymer–CO2 Mixtures

We combine a newly developed density-functional theory with the string method to calculate the minimum free energy path of bubble nucleation in compressible polymer-CO2 mixtures. Nucleation is initiated by saturating the polymer liquid with high pressure CO2 and subsequently reducing the pressure to ambient condition. Below a critical temperature, we find that there is a discontinuous drop in the nucleation barrier with increased initial CO2 pressure, as a result of an underlying metastable transition from a CO2-rich-vapor phase to a CO2-rich-liquid phase. This phenomenon is different from previously proposed nucleation mechanisms involving metastable transitions.

Efficiently finding the minimum free energy path from steepest descent path

Minimum Free Energy Path (MFEP) is very important in computational biology and chemistry. The barrier in the path is related to the reaction rate, and the start-to-end difference gives the relative stability between reactant and product. All these information is significant to experiment and practical application. But finding MFEP is not an easy job. Lots of degrees of freedom make the computation very complicated and time consuming. In this paper, we use the Steepest Descent Path (SDP) to accelerate the sampling of MFEP. The SHAKE algorithm and the Lagrangian multipliers are used to control the optimization of both SDP and MFEP. These strategies are simple and effective. For the former, it is more interesting. Because as we known, SHAKE algorithm was designed to handle the constraints in molecular dynamics in the past, has never been used in geometry optimization. Final applications on ALA dipeptide and 10-ALA peptide show that this combined optimization method works well. Use the information in SDP, the initial path could reach the more optimal MFEP. So more accurate free energies could be obtained and the amount of computation time could be saved.

Insight into the molecular mechanism of water evaporation via the finite temperature string method

The process of water's evaporation at its liquid/air interface has proven challenging to study experimentally and, because it constitutes a rare event on molecular time scales, presents a challenge for computer simulations as well. In this work, we simulated water's evaporation using the classical extended simple point charge model water model, and identified a minimum free energy path for this process in terms of 10 descriptive order parameters. The measured free energy change was 7.4 kcal/mol at 298 K, in reasonable agreement with the experimental value of 6.3 kcal/mol, and the mean first-passage time was 1375 ns for a single molecule, corresponding to an evaporation coefficient of 0.25. In the observed minimum free energy process, the water molecule diffuses to the surface, and tends to rotate so that its dipole and one O-H bond are oriented outward as it crosses the Gibbs dividing surface. As the water molecule moves further outward through the interfacial region, its local density is higher than the time-averaged density, indicating a local solvation shell that protrudes from the interface. The water molecule loses donor and acceptor hydrogen bonds, and then, with its dipole nearly normal to the interface, stops donating its remaining hydrogen bond. At that point, when the final, accepted hydrogen bond is broken, the water molecule is free. We also analyzed which order parameters are most important in the process and in reactive trajectories, and found that the relative orientation of water molecules near the evaporating molecule, and the number of accepted hydrogen bonds, were important variables in reactive trajectories and in kinetic descriptions of the process.

Activated pathways for the directed insertion of patterned nanoparticles into polymer membranes.

We combine the string method with self-consistent field theory to compute the most probable transition pathway, i.e. the minimum free energy path, for the insertion of Janus and protein-like nanoparticles into a polymer membrane bilayer. The method makes no assumptions in the reaction coordinate and overcomes the long timescales challenge associated with simulating rare events. Our study suggests that one approach to building functional polymer–nanoparticle composite membranes with oriented nanoparticles is through electrostatic interactions. In particular, hydrophobic Janus nanoparticles with an asymmetric charge distribution can be made to directionally insert into charged membranes. This process is kinetically driven, and involves overcoming a thermally surmountable activation barrier, which requires favorable interactions between the nanoparticle and the hydrophilic block of the membrane. In contrast, the insertion of protein-like nanoparticles with alternating hydrophilic–hydrophobic–hydrophilic domains into polymer membranes does not occur as a thermally activated event.

Bubble nucleation in polymer–CO2 mixtures.

We combine density-functional theory with the string method to calculate the minimum free energy path of bubble nucleation in two polymer–CO2 mixture systems, poly(methyl methacrylate) (PMMA)–CO2 and polystyrene (PS)–CO2. Nucleation is initiated by saturating the polymer liquid with high pressure CO2 and subsequently reducing the pressure to ambient condition. Below a critical temperature (Tc), we find that there is a discontinuous drop in the nucleation barrier as a function of increased initial CO2 pressure (P0), as a result of an underlying metastable transition from a CO2-rich-vapor phase to a CO2-rich-liquid phase. The nucleation barrier is generally higher for PS–CO2 than for PMMA–CO2 under the same temperature and pressure conditions, and both higher temperature and higher initial pressure are required to lower the nucleation barrier for PS–CO2 to experimentally relevant ranges. Classical nucleation theory completely fails to capture the structural features of the bubble nucleus and severely underestimates the nucleation barrier.

Ferryl–Oxo Species Produced from Fenton’s Reagent via a Two-Step Pathway: Minimum Free-Energy Path Analysis

A mixture of ferrous ions and hydrogen peroxide, known as Fenton's reagent, is an effective oxidant and has been widely used in various industrial applications; however, there is still controversy about what the oxidizing agents are and how they are produced. In this study, we have determined minimum free-energy paths (MFEPs) from Fenton's reagent to possible oxidizing agents such as hydroxyl radicals and ferryl-oxo species by combining ab initio molecular dynamics simulations and an MFEP search method. Along the MFEPs, representative free-energy profiles of the Fenton reaction were elucidated. On the basis of the free-energy profiles, we revealed that the reaction producing ferryl-oxo species from Fenton's reagent is more energetically favorable than that yielding a free hydroxyl radical, by 24.4 kcal mol(-1), which indicates that the ferryl-oxo species is the primary oxidizing agent in reactions of Fenton's reagent. Moreover, we clarified that the ferryl-oxo species is favorably formed via a two-step reaction pathway, which reaches the product through a dihydroxyiron(IV) intermediate. The energetics charting the free-energy profiles provided valuable information for a comprehensive understanding of Fenton reactions. We concluded that a ferryl-oxo species produced from Fenton's reagent serves as the primary oxidizing agent in the Fenton reaction.