Preorganization and reorganization as related factors in enzyme catalysis: the chorismate mutase case

In this paper a deeper insight into the chorismate-to prephenate-rearrangement, catalyzed by Bacillus subtilis chorismate mutase, is provided by means of a combination of statistical quantum mechanics/molecular mechanics simulation methods and hybrid potential energy surface exploration techniques. The main aim of this work is to present an estimation of the preorganization and reorganization terms of the enzyme catalytic rate enhancement. To analyze the first of these, we have studied different conformational equilibria of chorismate in aqueous solution and in the enzyme active site. Our conclusion is that chorismate mutase preferentially binds the reactive conformer of the substrate--that presenting a structure similar to the transition state of the reaction to be catalyzed--with shorter distances between the carbon atoms to be bonded and more diaxial character. With respect to the reorganization effect, an energy decomposition analysis of the potential energies of the reactive reactant and of the reaction transition state in aqueous solution and in the enzyme shows that the enzyme structure is better adapted to the transition structure. This means not only a more negative electrostatic interaction energy with the transition state but also a low enzyme deformation contribution to the energy barrier. Our calculations reveal that the structure of the enzyme is responsible for stabilizing the transition state structure of the reaction, with concomitant selection of the reactive form of the reactants. This is, the same enzymatic pattern that stabilizes the transition structure also promotes those reactant structures closer to the transition structure (i.e., the reactive reactants). In fact, both reorganization and preorganization effects have to be considered as the two faces of the same coin, having a common origin in the effect of the enzyme structure on the energy surface of the substrate.     

A hybrid elastic band string algorithm for studies of enzymatic reactions

A common challenge in theoretical biophysics is the identification of a minimum energy path (MEP) for the rearrangement of a group of atoms from one stable configuration to another. The structure with maximum energy along the MEP approximates the transition state for the process and the energy profile itself permits estimation of the transition rates. In this work we describe a computationally efficient algorithm for the identification of minimum energy paths in complicated biosystems. The algorithm is a hybrid of the nudged elastic band (NEB) and string methods. It has been implemented in the pDynamo simulation program and tested by examining elementary steps in the reaction mechanisms of three enzymes: citrate synthase, RasGAP, and lactate dehydrogenase. Good agreement is found for the energies and geometries of the species along the reaction profiles calculated using the new algorithm and previous versions of the NEB and string techniques, and also those obtained by the common method of adiabatic exploration of the potential energy surface as a function of predefined reaction coordinates. Precisely refined structures of the saddle points along the paths may be subsequently obtained with the climbing image variant of the NEB algorithm. Directions in which the utility of the methods that we have implemented can be further improved are discussed.     

A hybrid potential reaction path and free energy study of the chorismate mutase reaction.

We present a combination of two techniques--QM/MM statistical simulation methods and QM/MM internal energy minimizations--to get a deeper insight into the reaction catalyzed by the enzyme chorismate mutase. Structures, internal energies and free energies, taken from the paths of the reaction in solution and in the enzyme have been analyzed in order to estimate the relative importance of the reorganization and preorganization effects. The results we obtain for this reaction are in good agreement with experiment and show that chorismate mutase achieves its catalytic efficiency in two ways; first, it preferentially binds the active conformer of the substrate and, second, it reduces the free energy of activation for the reaction relative to that in solution by providing an environment which stabilizes the transition state.     

How enzyme dynamics helps catalyze a reaction in atomic detail: a transition path sampling study

We have applied the Transition Path Sampling algorithm to the reaction catalyzed by the enzyme Lactate Dehydrogenase. This study demonstrates the ease of scaling Transition Path Sampling for applications on many degree of freedom systems, whose energy surface is a complex terrain of valleys and saddle points. As a Monte Carlo importance sampling method, transition path sampling is capable of surmounting barriers in path phase space and focuses simulation on the rare event of enzyme catalyzed atom transfers. Generation of the transition path ensemble, for this reaction, resolves a paradox in the literature in which some studies exposed the catalytic mechanism of hydride and proton transfer by lactate dehydrogenase to be concerted and others stepwise. Transition path sampling has confirmed both mechanisms as possible paths from reactants to products. With the objective to identify a generalized, reduced reaction coordinate, time series of both donor-acceptor distances and residue distances from the active site have been examined. During the transition from pyruvate to lactate, residues located behind the transferring hydride collectively compress toward the active site causing residues located behind the hydride acceptor to relax away. It is demonstrated that an incomplete compression/relaxation transition across the donor-acceptor axis compromises the reaction.     

Ab initio quantum mechanical/molecular mechanical molecular dynamics simulation of enzyme catalysis: the case of histone lysine methyltransferase SET7/9

To elucidate enzyme catalysis through computer simulation, a prerequisite is to reliably compute free energy barriers for both enzyme and solution reactions. By employing on-the-fly Born-Oppenheimer molecular dynamics simulations with the ab initio quantum mechanical/molecular mechanical approach and the umbrella sampling method, we have determined free energy profiles for the methyl-transfer reaction catalyzed by the histone lysine methyltransferase SET7/9 and its corresponding uncatalyzed reaction in aqueous solution, respectively. Our calculated activation free energy barrier for the enzyme catalyzed reaction is 22.5 kcal/mol, which agrees very well with the experimental value of 20.9 kcal/mol. The difference in potential of mean force between a corresponding prereaction state and the transition state for the solution reaction is computed to be 30.9 kcal/mol. Thus, our simulations indicate that the enzyme SET7/9 plays an essential catalytic role in significantly lowering the barrier for the methyl-transfer reaction step. For the reaction in solution, it is found that the hydrogen bond network near the reaction center undergoes a significant change, and there is a strong shift in electrostatic field from the prereaction state to the transition state, whereas for the enzyme reaction, such an effect is much smaller and the enzyme SET7/9 is found to provide a preorganized electrostatic environment to facilitate the methyl-transfer reaction. Meanwhile, we find that the transition state in the enzyme reaction is a little more dissociative than that in solution.     

String method in collective variables: minimum free energy paths and isocommittor surfaces

A computational technique is proposed which combines the string method with a sampling technique to determine minimum free energy paths. The technique only requires to compute the mean force and another conditional expectation locally along the string, and therefore can be applied even if the number of collective variables kept in the free energy calculation is large. This is in contrast with other free energy sampling techniques which aim at mapping the full free energy landscape and whose cost increases exponentially with the number of collective variables kept in the free energy. Provided that the number of collective variables is large enough, the new technique captures the mechanism of transition in that it allows to determine the committor function for the reaction and, in particular, the transition state region. The new technique is illustrated on the example of alanine dipeptide, in which we compute the minimum free energy path for the isomerization transition using either two or four dihedral angles as collective variables. It is shown that the mechanism of transition can be captured using the four dihedral angles, but it cannot be captured using only two of them.     

几个P-ylide反应机理的量子拓扑研究

采用MP2(FC)/6-311++G(d,p)对磷叶立德和类磷叶立德自由基反应机理进行了探讨.优化了中间体、过渡态和产物的几何构型,并采用内禀反应坐标法进行追踪.侧重从量子拓扑学的角度,对反应过程中各点进行电子密度拓扑分析,讨论了反应过程中化学键的断裂、生成和化学键的变化规律.上述两个反应都经历三员环过渡态,找到了这类反应的能量过渡态和结构过渡态.     

Exploring potential energy surfaces for chemical reactions: an overview of some practical methods

Potential energy surfaces form a central concept in the application of electronic structure methods to the study of molecular structures, properties, and reactivities. Recent advances in tools for exploring potential energy surfaces are surveyed. Methods for geometry optimization of equilibrium structures, searching for transition states, following reaction paths and ab initio molecular dynamics are discussed. For geometry optimization, topics include methods for large molecules, QM/MM calculations, and simultaneous optimization of the wave function and the geometry. Path optimization methods and dynamics based techniques for transition state searching and reaction path following are outlined. Developments in the calculation of ab initio classical trajectories in the Born-Oppenheimer and Car-Parrinello approaches are described.     

Exploring reaction pathways with transition path and umbrella sampling: application to methyl maltoside

The transition path sampling (TPS) method is a powerful approach to study chemical reactions or transitional properties on complex potential energy landscapes. One of the main advantages of the method over potential of mean force methods is that reaction rates can be directly accessed without knowledge of the exact reaction coordinate. We have investigated the complementary nature of these two differing approaches, comparing transition path sampling with the weighted histogram analysis method to study a conformational change in a small model system. In this case study, the transition paths for a transition between two rotational conformers of a model disaccharide molecule, methyl beta-D-maltoside, were compared with a free energy surface constrained by the two commonly used glycosidic (phi,psi) torsional angles. The TPS method revealed a reaction channel that was not apparent from the potential of mean force method, and the suitability of phi and psi as reaction coordinates to describe the isomerization in vacuo was confirmed by examination of the transition path ensemble. Using both transition state theory and transition path sampling methods, the transition rate was estimated. We have estimated a characteristic time between transitions of approximately 160 ns for this rare isomerization event between the two conformations of the carbohydrate. We conclude that transition path sampling can extract subtle information about the dynamics not apparent from the potential of mean force method. However, in calculating the reaction rate, the transition path sampling method required 27.5 times the computational effort than was needed by the potential of mean force method.     

Comparison of reaction pathways calculated by different algorithms for disilane and water trimer

There still exists some confusion in the literature concerning the definition of a minimum energy pathway and the coordinate system in which it is calculated. Here we compare steepest-descent and eigenvector-following pathways, both with and without a mass-weighted metric. The systems studied are disilane and the water trimer, and we employ various basis sets at the SCF level of theory. We find that paths calculated using eigenvector-following and steepestdescent are practically the same, at least in terms of the reaction mechanism.We find that for the mass-weighted metric the pathways are similar, although in principle they do not have to be identical. Finally, we verify that the geometrical symmetry selection rules hold for a pathway mediated by a recently discovered transition state of the disilane system.     

On the nature of the transition state in catechol O-methyltransferase. A complementary study based on molecular dynamics and potential energy surface explorations

The way in which enzymes influence the rate of chemical processes is still a question of debate. The protein promotes the catalysis of biochemical processes by lowering the free energy barrier in comparison with the reference uncatalyzed reaction in solution. In this article we are reporting static and dynamic aspects of the enzyme catalysis in a bimolecular reaction, namely a methyl transfer from S-adenosylmethionine to the hydroxylate oxygen of a substituted catechol catalyzed by catechol O-methyltransferase. From QM/MM optimizations, we will first analyze the participation of the environment on the transition vector. The study of molecular dynamics trajectories will allow us to estimate the transmission coefficient from a previously localized transition state as the maximum in the potential of mean force profile. The analysis of the reactive and nonreactive trajectories in the enzyme environment and in solution will also allow studying the geometrical and electronic changes, with special attention to the chemical system movements and the coupling with the environment. The main result, coming from both analyses, is the approximation of the magnesium cation to the nucleophilic and the hydroxyl group of the catecholate as a result of a general movement of the protein, stabilizing in this way the transition state. Consequently, the free energy barrier of the enzyme reaction is dramatically decreased with respect to the reaction in solution.     

Theoretical Studies of PCl3/H2 Gas Chromatography Mechanism

Density functional theory B3LYP method was firstly applied to fully optimize the geometrical configuration of each stable point on PCl3/H2 gas chromatography reaction potential energy surface on the 6-311G** basis set, and single point energy was computed at the QCISD(T)/ 6-311G** level, then the transition state was validated by analyzing the unique imaginary vibration modes of each transition state and calculating intrinsic reaction coordinate (IRC), and the major reaction and competing reaction paths of PCl3/H2 gas chromatography reaction were presented through comparing active energy barrier, and phosphor was finally gained from the reaction of PH and PCl.     

Analysis of classical and quantum paths for deprotonation of methylamine by methylamine dehydrogenase.

The hydrogen-transfer reaction catalysed by methylamine dehydrogenase (MADH) with methylamine (MA) as substrate is a good model system for studies of proton tunnelling in enzyme reactions--an area of great current interest--for which atomistic simulations will be vital. Here, we present a detailed analysis of the key deprotonation step of the MADH/MA reaction and compare the results with experimental observations. Moreover, we compare this reaction with the related aromatic amine dehydrogenase (AADH) reaction with tryptamine, recently studied by us, and identify possible causes for the differences observed in the measured kinetic isotope effects (KIEs) of the two systems. We have used combined quantum mechanics/molecular mechanics (QM/MM) techniques in molecular dynamics simulations and variational transition state theory with multidimensional tunnelling calculations averaged over an ensemble of paths. The results reveal important mechanistic complexity. We calculate activation barriers and KIEs for the two possible proton transfers identified-to either of the carboxylate oxygen atoms of the catalytic base (Asp428beta)-and analyse the contributions of quantum effects. The activation barriers and tunnelling contributions for the two possible proton transfers are similar and lead to a phenomenological activation free energy of 16.5+/-0.9 kcal mol(-1) for transfer to either oxygen (PM3-CHARMM calculations applying PM3-SRP specific reaction parameters), in good agreement with the experimental value of 14.4 kcal mol(-1). In contrast, for the AADH system, transfer to the equivalent OD1 was found to be preferred. The structures of the enzyme complexes during reaction are analysed in detail. The hydrogen bond of Thr474beta(MADH)/Thr172beta(AADH) to the catalytic carboxylate group and the nonconserved active site residue Tyr471beta(MADH)/Phe169beta(AADH) are identified as important factors in determining the preferred oxygen acceptor. The protein environment has a significant effect on the reaction energetics and hence on tunnelling contributions and KIEs. These environmental effects, and the related clearly different preferences for the two carboxylate oxygen atoms (with different KIEs) in MADH/MA and AADH/tryptamine, are possible causes of the differences observed in the KIEs between these two important enzyme reactions.     

Simulation studies of the fidelity of biomolecular structure ensemble recreation

We examine the ability of Bayesian methods to recreate structural ensembles for partially folded molecules from averaged data. Specifically we test the ability of various algorithms to recreate different transition state ensembles for folding proteins using a multiple replica simulation algorithm using input from "gold standard" reference ensembles that were first generated with a Go-like Hamiltonian having nonpairwise additive terms. A set of low resolution data, which function as the "experimental" phi values, were first constructed from this reference ensemble. The resulting phi values were then treated as one would treat laboratory experimental data and were used as input in the replica reconstruction algorithm. The resulting ensembles of structures obtained by the replica algorithm were compared to the gold standard reference ensemble, from which those "data" were, in fact, obtained. It is found that for a unimodal transition state ensemble with a low barrier, the multiple replica algorithm does recreate the reference ensemble fairly successfully when no experimental error is assumed. The Kolmogorov-Smirnov test as well as principal component analysis show that the overlap of the recovered and reference ensembles is significantly enhanced when multiple replicas are used. Reduction of the multiple replica ensembles by clustering successfully yields subensembles with close similarity to the reference ensembles. On the other hand, for a high barrier transition state with two distinct transition state ensembles, the single replica algorithm only samples a few structures of one of the reference ensemble basins. This is due to the fact that the phi values are intrinsically ensemble averaged quantities. The replica algorithm with multiple copies does sample both reference ensemble basins. In contrast to the single replica case, the multiple replicas are constrained to reproduce the average phi values, but allow fluctuations in phi for each individual copy. These fluctuations facilitate a more faithful sampling of the reference ensemble basins. Finally, we test how robustly the reconstruction algorithm can function by introducing errors in phi comparable in magnitude to those suggested by some authors. In this circumstance we observe that the chances of ensemble recovery with the replica algorithm are poor using a single replica, but are improved when multiple copies are used. A multimodal transition state ensemble, however, turns out to be more sensitive to large errors in phi (if appropriately gauged) and attempts at successful recreation of the reference ensemble with simple replica algorithms can fall short.     

A theoretical study of the catalytic mechanism of formate dehydrogenase

A theoretical study of the hydride transfer between formate anion and nicotinamide adenine dinucleotide (NAD(+)) catalyzed by the enzyme formate dehydrogenase (FDH) has been carried out by a combination of two hybrid quantum mechanics/molecular mechanics techniques: statistical simulation methods and internal energy minimizations. Free energy profiles, obtained for the reaction in the enzyme active site and in solution, allow obtaining a comparative analysis of the behavior of both condensed media. Moreover, calculations of the reaction in aqueous media can be used to probe the dramatic differences between reactants state in the enzyme active site and in solution. The results suggest that the enzyme compresses the substrate and the cofactor into a conformation close to the transition structure by means of favorable interactions with the amino acid residues of the active site, thus facilitating the relative orientation of donor and acceptor atoms to favor the hydride transfer. Moreover, a permanent field created by the protein reduces the work required to reach the transition state (TS) with a concomitant polarization of the cofactor that would favor the hydride transfer. In contrast, in water the TS is destabilized with respect to the reactant species because the polarity of the solute diminishes as the reaction proceeds, and consequently the reaction field, which is created as a response to the change in the solute polarity, is also decreased. Therefore protein structure is responsible of both effects; substrate preorganization and TS stabilization thus diminishing the activation barrier. Because of the electrostatic features of the catalyzed reaction, both media preferentially stabilize the ground-state, thus explaining the small rate constant enhancement of this enzyme, but FDH does so to a much lower extent than aqueous solution. Finally, a good agreement between experimental and theoretical kinetic isotope effects is found, thus giving some credit to our results.     

Tunneling and classical paths for proton transfer in an enzyme reaction dominated by tunneling: oxidation of tryptamine by aromatic amine dehydrogenase

Proton tunneling dominates the oxidative deamination of tryptamine catalyzed by the enzyme aromatic amine dehydrogenase. For reaction with the fast substrate tryptamine, a H/D kinetic isotope effect (KIE) of 55 +/- 6 has been reported-one of the largest observed in an enzyme reaction. We present here a computational analysis of this proton-transfer reaction, applying combined quantum mechanics/molecular mechanics (QM/MM) methods (PM3-SRP//PM3/CHARMM22). In particular, we extend our previous computational study (Masgrau et al. Science 2006, 312, 237) by using improved energy corrections, high-level QM/MM methods, and an ensemble of paths to estimate the tunneling contributions. We have carried out QM/MM molecular dynamics simulations and variational transition state theory calculations with small-curvature tunneling corrections. The results provide detailed insight into the processes involved in the reaction. Transfer to the O2 oxygen of the catalytic base, Asp128beta, is found to be the favored reaction both thermodynamically and kinetically, even though O1 is closer in the reactant complex. Comparison of quantum and classical models of proton transfer allows estimation of the contribution of hydrogen tunneling in lowering the barrier to reaction in the enzyme. A reduction of the activation free energy due to tunneling of 3.1 kcal mol-1 is found, which represents a rate enhancement due to tunneling by 2 orders of magnitude. The calculated KIE of 30 is significantly elevated over the semiclassical limit, in agreement with the experimental observations; a semiclassical value of 6 is obtained when tunneling is omitted. A polarization of the C-H bond to be broken is observed due to the close proximity of the catalytic aspartate and the (formally) positively charged imine nitrogen. A comparison is also made with the related quinoprotein methylamine dehydrogenase (MADH)-the much lower KIE of 11 that we obtain for the MADH/methylamine system is found to arise from a more endothermic potential energy surface for the MADH reaction.     

Theoretical modeling of enzyme catalytic power: analysis of "cratic" and electrostatic factors in catechol O-methyltransferase

A comparative theoretical study of a bimolecular reaction in aqueous solution and catalyzed by the enzyme catechol O-methyltransferase (COMT) has been carried out by a combination of two hybrid QM/MM techniques: statistical simulation methods and internal energy minimizations. In contrast to previous studies by other workers, we have located and characterized transition structures for the reaction in the enzyme active site, in water and in a vacuum, and our potential of mean force calculations are based upon reaction coordinates obtained from features of the potential energy surfaces in the condensed media, not from the gas phase. The AM1/CHARMM calculated free energy of activation for the reaction of S-adenosyl methionine (SAM) with catecholate catalyzed by COMT is 15 kcal mol(-1) lower the AM1/TIP3P free-energy barrier for the reaction of the trimethylsulfonium cation with the catecholate anion in water at 300 K, in agreement with previous estimates. The thermodynamically preferred form of the reactants in the uncatalyzed model reaction in water is a solvent-separated ion pair (SSIP). Conversion of the SSIP into a contact ion pair, with a structure resembling that of the Michaelis complex (MC) for the reaction in the COMT active site, is unfavorable by 7 kcal mol(-1), largely due to reorganization of the solvent. We have considered alternative ways to estimate the so-called "cratic" free energy for bringing the reactant species together in the correct orientation for reaction but conclude that direct evaluation of the free energy of association by means of molecular dynamics simulation with a simple standard-state correction is probably the best approach. The latter correction allows for the fact that the size of the unit cell employed with the periodic boundary simulations does not correspond to the standard state concentration of 1 M. Consideration of MC-like species allows a helpful decomposition of the catalytic effect into preorganization and reorganization phases. In the preorganization phase, the substrates are brought together into the MC-like species, either in water or in the enzyme active site. In the reorganization phase, the roles of the enzymic and aqueous environments may be compared directly because reorganization of the substrate is about the same in both cases. Analysis of the electric field along the reaction coordinate demonstrates that in water the TS is destabilized with respect to the MC-like species because the polarity of the solute diminishes and consequently the reaction field is also decreased. In the enzyme, the electric field is mainly a permanent field and consequently there is only a small reorganization of the environment. Therefore, destabilization of the TS is lower than in solution, and the activation barrier is smaller.     

Generation of initial trajectories for transition path sampling of chemical reactions with ab initio molecular dynamics

Transition path sampling is an innovative method for focusing a molecular dynamics simulation on a reactive event. Although transition path sampling methods can generate an ensemble of reactive trajectories, an initial reactive trajectory must be generated by some other means. In this paper, the authors have evaluated three methods for generating initial reactive trajectories for transition path sampling with ab initio molecular dynamics. The authors have tested each of these methods on a set of chemical reactions involving the breaking and making of covalent bonds: the 1,2-hydrogen elimination in the borane-ammonia adduct, a tautomerization, and the Claisen rearrangement. The first method is to initiate trajectories from the potential energy transition state, which was effective for all reactions in the test set. Assigning atomic velocities found using normal mode analysis greatly improved the success of this method. The second method uses a high temperature molecular dynamics simulation and then iteratively reduces the total energy of the simulation until a low temperature reactive trajectory is found. This was effective in generating a low temperature trajectory from an initial trajectory run at 3000 K of the tautomerization reaction, although it failed for the other two. The third uses an orbital based bias potential to find a reactive trajectory and uses this trajectory to initiate an unbiased trajectory. The authors found that a highest occupied molecular orbital-lowest unoccupied molecular orbital bias could be used to find a reactive trajectory for the Claisen rearrangement, although it failed for the other two reactions. These techniques will help make it practical to use transition path sampling to study chemical reaction mechanisms that involve bond breaking and forming.     

Gibbs free-energy estimates from direct path-sampling computations

We have implemented a path-sampling scheme enabling a direct estimation of Gibbs free energy. This scheme consists of a Monte Carlo sampling of constant-pressure Langevin paths, followed by an ensemble averaging carried out over the Markov chain of paths. In practice, we sample an umbrella path ensemble, which requires to rigorously define a statistical weight for the paths, equivalent of the Boltzmann weight. This statistical weight is a function of an effective work related to the path. The umbrella ensemble is chosen so that its work histogram overlaps with the histograms corresponding to the ensembles of forward and backward paths. We have finally investigated the relations between numerical efficiency and overlapping properties of the various work histograms. This analysis yields a built-in criterion for diagnosing the convergence during a single-run simulation.