Catalytic Reaction Mechanism in Native and Mutant Catechol-O

Subscriber access provided by University of South Dakota

B: Biophysics; Physical Chemistry of Biological Systems and Biomolecules

Catalytic Reaction Mechanism in Native and Mutant COMT from the Adaptive String Method and Mean Reaction Force Analysis David Adrian Saez, Kirill Zinovjev, Iñaki Tuñón, and Esteban Vöhringer-Martinez J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b07339 • Publication Date (Web): 30 Aug 2018 Downloaded from http://pubs.acs.org on September 4, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Catalytic Reaction Mechanism in Native and Mutant COMT from the Adaptive String Method and Mean Reaction Force Analysis David Adrian Saez,† Kirill Zinovjev,‡ Iñaki Tuñoń,‡ and Esteban Vöhringer-Martinez∗,† †Departamento de F´ısico-Qu´ımica, Facultad de Ciencias Qu´ımicas, Universidad de Concepción, Chile ‡Departament de Qu´ımica F´ısica, Universitat de València, 46100 Burjassot, Spain E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Catechol-O-Methyltransferase is an enzyme that catalyzes the methylation reaction of dopamine by S-Adenosylmethionine increasing the reaction rate by almost 16 orders of magnitude compared to the reaction in aqueous solution. Here, we combine the recently introduced adaptive string method and the Mean Reaction Force method in combination with structural and electronic descriptors to characterize the reaction mechanism. The catalytic effect of the enzyme is addressed by comparison of the reaction in the human wild-type enzyme, in the less effective Y68A mutant and in aqueous solution. The influence of these different environments at different stages of the chemical process and the significance of key collective variables describing the reaction were quantified. Our results show that the native enzyme limits the access of water molecules to the active site, enhancing the interaction between the reactants and providing a more favorable electrostatic environment to assist the SN 2 methyl transfer reaction.

Introduction Enzymes are biological macromolecules able to catalyze a wide range of chemical reactions. 1 Due to their complex structure and flexibility it is often not trivial to define a reaction coordinate which represents correctly the reaction flux between reactants and products. Even when the reaction coordinate is known, the detailed mechanism by which the enzyme reduces the activation barrier might still be unclear. The main factors contributing to catalysis involve a combination of pre-organization of the reactants, electrostatic stabilization 1–5 or in some special cases even conformational selection. 6 Assuming a correct description of the interactions in the catalytic site and a complete sampling of the configurational space, a method for rational separation of contributions to the activation barrier along the reaction coordinate is needed. Ideally, these contributions should be explained in terms of structural modifications and electronic changes taking place inside the reactive center due to interactions between the 2


Page 2 of 35

Page 3 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


substrates and the enzyme. Recently, we have shown that the Mean Reaction Force method is able to distinguish between structural and electronic contributions to the activation barrier of reactions in aqueous phase and enzymes. 7–10 Distinguishing these factors would contribute to the understanding of the reactants preorganization or the electrostatic stabilization. In this contribution we combine the Mean Reaction Force method with the recently introduced adaptive string method, 11 which identifies the minimum free energy path, to analyze the Catechol-O-Methyltransferase (COMT) enzyme. Comparing the reaction in the wild-type enzyme and in aqueous solution, it is possible to obtain new insights in the reaction and at the same time explain the decreased efficiency of the Y68A mutant. COMT belongs to the family of the Methyltransferases that catalyze the methylation of diverse substrates and are involved in various physiological processes. 12–14 COMT employs S-Adenosylmethionine (SAM) to methylate different catecholic substrates, like neurotransmitters and hormones. This enzyme enhances the reaction rate by up to 16 orders of magnitude with respect to the corresponding reaction occurring in solution 15 without covalently binding SAM or the reactants. 16,17 The early experimental studies provided the first insights into the influence of the chemical environment in the methyl transfer to catecholic substrates. 16–18 These studies reported larger values for the inversion of the secondary kinetic isotope effect (2◦ -KIE) in the enzyme than in water, which was interpreted as a compression of the transition state inside the enzyme and the major factor in the rate enhancement. 16 Computational chemistry approaches provided a more detailed microscopic picture of the catalytic cavity. 19–21 The first studies suggested that the main effect was to select reactive conformations called near-attack-conformations (NAC) that were poorly sampled in solution. These reactive conformations had a smaller distance between the reactants, making a link with the compression hypothesis. 20 Simultaneously, Roca et al. 5,22,23 suggested that the apparent role of NACs could be understood as a side effect of the evolution of the enzyme’s active site to stabilize the enzymatic transition state with respect to the transition state (TS) occurring in solution: the dipoles inside the enzyme

3



are pre-organized to stabilize the TS. In contrast, in aqueous solution the transformation from polar reactants to a less polar transition state requires an energetic cost of solvent reorganization. In the last years, new experimental evidence 24 showed that the mutation of the residue Tyr68 to alanine in the human enzyme decreased the enzyme’s activity, and the inversion of the 2◦ KIE approached the value observed in aqueous solution. Empirical valence bond calculations assessed the activation energy for the wild-type (WT) and Ala68 mutant (Y68A) 25 concluding that the residue in position 68 is not able to compress or push the reactants in any way, since the simulations do not show a sensible variation in the interatomic distances. Further experimental data including Stokes shifts and fluorescence lifetimes addressed the influence of Tyr68. 26 In the wild-type, the observed larger Stokes shift and reduced lifetime for the excited state dipole at Trp143 were interpreted as a more rigid environment surrounding Trp143 after the binding of SAM. The mutation of the residue in position 68, approximately 7˚ A away, resulted in a Stokes shift comparable to the apoenzyme, which would suggest a less structured active site in the mutant. These experimental results and recent advances in electronic structure calculations triggered new computational studies, which addressed the dependence of the reaction energetics on the size of the reactive center treated with the electronic structure methods in a QM/MM approach. 4,27–29 The main conclusion was that the activation barrier is not altered significantly once the magnesium ion and its first coordination shell was added to the QM region together with the substrates and enough sampling of the relevant conformations was achieved. However, larger variations in the free energy difference between reactants and products were observed depending on the number of residues in this region. Most of these calculations had been done with semi-empirical methods such as PM6, which might induce additional errors to the calculation, especially when the coordination of the magnesium cation is considered. Recently, Kulik et al. reported free energy calculations employing a DFT method (PBE/6-31G(d)) and confirmed the conclusions showing additional charge transfer effects in

4


Page 4 of 35

Page 5 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the catalytic cavity. 30 All these studies conclude that these factors are significantly smaller than the catalytic effect and that they are important only for the layer of residues in close contact with the substrate and cofactor, since the free energy barrier is not greatly altered with additional residues included in the QM region. Most of the recent research has been focused on reproducing the experimental value of the activation barrier. Here, we address the reaction in the wild-type (WT) enzyme and the Y68A mutant, and compare the results with the aqueous environment to rationalize the catalytic effect of the enzyme. Before addressing the factors contributing to catalysis, a proper conformational sampling of the enzymes with the bound substrates employing previously parameterized force field parameters for cofactor and substrate 31 was performed to provide representative reactant conformations in the enzyme. Interestingly, these conformations differ between the mutant and WT due to the proximity of water in the active site in the latter. With the reactant state of each system we first used the adaptive string method to identify the Minimum Free Energy Path (MFEP) for each environment in a space of three collective variables (CV) whose changes are linked to the progress of the reaction. 11,32,33 The mechanistic details leading to the catalytic activity of the enzyme were analyzed and quantified from restrained molecular dynamics simulations along the MFEP-based path CV by the Mean Reaction Force. 7 Structural and electronic contributions to the free energy barrier and the decomposition of the free energy barrier in terms of the CVs provided new insights for the reaction in the enzyme the mutant and in aqueous solution, which will contribute to the general understanding of the enzyme catalysis.

Methods System setup The enzyme (WT and Y68A) was centered in a rectangular box of water molecules (TIP3P) maintaining a minimum distance of 10 ˚ A from any protein atom to its edge. The coordinates of the protein atoms, the substrate and SAM were taken from

5



previous work based on the murine enzyme 22,34 . The structure was based on the crystal structure published in 2008 35 (PDB code 3BWM), and modified according to the calculations performed on a similar system by Roca et al. (see Figure 1). 5 All residues were in their standard protonation states. Histidine residues were singly protonated on Nδ1 (57, 193) and N2 (12, 16, 142, 182), respectively. The box contained 7978 solvent molecules and five sodium ions to achieve neutrality. Periodic boundary conditions were used in combination with the Particle Mesh Ewald to account for the long-range electrostatic interactions, using a non-bonded cut-off value of 8˚ A . Parameters for SAM and dopamine were taken from our previous work 31 and GAFF 36 , respectively and the rest of the system was described by the ff14SB force field. 37 All molecular dynamics calculations were performed using tools from the Amber16 and AmberTools17 packages. 38 The mutant and wild-type enzyme were minimized first restraining the position of all non-solvent atoms and then the whole system. Both systems were equilibrated for 200 ps under NVT conditions using a thermostat 39 and for 5 ns applying the Monte Carlo barostat and Langevin dynamics (T=298K, τ =1 ps−1 ) using the pmemd.cuda.MPI program in Amber16. 38 Starting from the output structure of the NPT equilibration, five independent simulations of 100 ns with randomly assigned velocities were ˚−2 ) of 2.5˚ A between the performed maintaining a distance restraint (k = 100 kcal mol−1 A oxygen atom of the DOP hydroxyl group and the magnesium cation. During the dynamics all bonds involving hydrogen were constrained with the SHAKE algorithm and an integration timestep of 2 fs was used. In all simulations a converged RMSD of the backbone atoms and the equilibration of the number of water molecules in the catalytic cavity in the mutant structure was obtained. The number of hydration sites in the catalytic cavity of the WT and Y68A mutant was calculated using SSTMap. 40 For the analysis two representative structures for each system were obtained from a cluster analysis (average-linkage method on RMSD of the carbon α atoms of all protein residues) of the combined trajectories and simulated with position restraints on all protein atoms for 5 ns. In both representative structures the obtained

6


Page 6 of 35

Page 7 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


number of hydration sites were the same for both systems. QM/MM simulations were performed with the sander.MPI module of the AmberTools17 package 38 , minmizing the system first and the reactant structure was equilibrated for 50 ps in the NVT ensemble with the Langevin thermostat at 300 K applying an integration timestep of 1 fs. The QM region comprising DOP, SAM and the magnesium ion was described with the DFTB3 Hamiltonian and the 3ob parameters. 41–44 We have shown previously that this method results in electronic energies of similar quality as high level extrapolated abinitio methods for a model system of this reaction and that it also reproduces experimental activation and reaction free energies in aqueous solution. 10 Remaining atoms were treated classically using the ff14SB force field. 37 The product structure was obtained by steered molecular dynamics transferring the methyl group to the negatively charged oxygen atom of DOP in the equilibrated reactant structure. This structure was then equilibrated further for 50 ps. The resulting reactant and product structures in the WT enzyme were then used to find the MFEP in the enzymatic environment . The same procedure was applied to the mutant Y68A to obtain reactant and product equilibrated structures. The aqueous system consisted of a cubic box with 2271 solvent molecules and one chloride ion to achieve neutrality. The QM/MM partitioning of the substrates was the same as in the enzymatic system, not considering the magnesium ion. The simulation protocol was identical to that used for the enzymatic systems. Minimum Free Energy Path The MFEP was obtained with the string method as implemented in the AmberTools17 package 11,38 following a common scheme for the three systems under study. The progress of the reaction was followed through three collective variables directly involved in the process of the methyl transfer: the distance between the sulfur atom of SAM and the carbon atom of the transferred methyl group (S-C distance or CV1); the distance between the aforementioned carbon atom and the negatively charged oxygen of DOP (C-O distance or CV2) and the hybridization of the carbon atom in the methyl group, expressed as the point-plane distance defined by the three hydrogens and the carbon atom

7



itself (CV3). The semi-empirical DFTB3 method with the 3ob parameters 41–44 was used to describe the atoms in the QM region. Although it has been shown that the size of the QM region can have an effect of up to 5 kcal/mol in the activation barrier, 4,30 here, we are interested in the reaction mechanism and the differences between wild-type, mutant and aqueous systems. The employed QM method in this study, and also the electronic structure methods in the two previous studies (PM6 and PBE/6-31G(d)), can present errors in the range of 5 kcal/mol compared to high level ab-initio methods, 10 which have to be combined with errors from improper sampling in QM/MM simulations. Therefore, the reported deviations of the activation barriers from the experimental value with the QM region might also have other origins than the size itself. The implementation of the string method used to obtain the MFEP is explained in detail elsewhere. 11 Here, the string consisted of 64 nodes, with the initial and final ones defined by the collective variables obtained from the QM/MM minimization and 50 ps equilibration of reactants and products. In the calculation of the MFEP in the WT no nodes were fixed, allowing a complete relaxation of the reactants and products. To keep the systems comparable to the wild-type, the initial (reactant) and final (product) distances between atoms S, C and O were restrained to the values in the converged wild-type enzymatic string, for both the Y68A and AQ systems. The initial points for the string nodes were obtained by linearly interpolating the values of the CVs in reactants and products. The MD replicas were given 5 ps of relaxation, gradually increasing the force constant of the string bias. To keep the system close to the path, an orthogonal potential of 400 kcal mol−1 (amu1/2 ˚ A)−2 was applied during the evolution of the string. Then, every node was allowed to drift, following the direction of the mean force in the selected CV space. To improve the sampling of the configurational space, replica exchange was attempted between neighboring nodes every 20 steps. The convergence of the string was monitored by measuring the mean distance between the past and current positions of the string nodes. Convergence of the MFEP was obtained after 30 ps and the interval between 40 and 50 ps was used to obtain converged

8


Page 8 of 35

Page 9 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


parameters of the string. With these parameters we continued the simulation in the potential of mean force (PMF) stage during 42 ps. This PMF is obtained along a path-CV (the socalled s coordinate) that measures the advance of the system along the MFEP. 33 Once the string in the WT enzyme was converged, the values of the collective variables for each node (time-averaged over the 40-50 ps) were used as an initial guess to define the string in the aqueous and mutant environments. Also here, once the string reached convergence the simulation of the string in the PMF phase was continued during 54 ps (Y68A) and 70 ps (AQ). PMFs were obtained using Umbrella Integration 45 which allowed to calculate the 95% confidence intervals analytically, assuming normal distribution of the sampling along the reaction coordinate in each US window. 11,33 The different lengths of the PMF phase aimed to reach the same value of the 95% confidence interval for the free energy barrier (± 1 kcal/mol). The analysis of the final, converged string was performed over a set of 64 trajectories, one for each node. Each trajectory was processed using the software GROMACS 4.5.4. Different features were evaluated, as average distances and angles between atoms involved in the reaction. For each frame of each trajectory a single point calculation was performed at B3LYP/def2-SV(P) level of theory, using the software ORCA 46 in an QM/MM approach. The obtained electron density was used in the Horton 2.0 package 47 to derive the HirshfeldI charges on every QM atom, for every frame in the short trajectory. Finally, average Hirshfeld-I charges were calculated for every QM-atom in each node. Hirshfeld-I atomic charges Since the reaction involves annihilation of charge, Hirshfeld-I atomic charges were used to monitor the charge transfer process and to assess the environmental effect. 8,10 Hirshfeld-I atomic charges provide a quick and intuitive method to evaluate the electron density distribution in a molecule 48 . Using a modified version of the original scheme proposed by Hirshfeld 49 , it allows the determination of atomic populations based only on the polarized electron density of the system, employing a database of proatomic P densities of charged and neutral atoms: the promolecular density 0B=1 ρ0B (r) is obtained

9



Page 10 of 35

superimposing the calculated proatomic densities of the spherically averaged proatoms ρ0A (r). This promolecular density allows the determination of the contribution of each atom to the molecular density ρ(r). In this way, a density can be assigned to an atom A in the molecule, according to the formula: ρ0 (r) ρA (r) = ρ(r) PN A 0 B=1 ρB (r)

(1)

The charge of this atom is defined through its atomic number ZA as: Z qA = ZA −

ρA (r)dr

(2)

In the Hirshfeld-I method the charge of every proatom is updated until consistence with its calculated charge in the molecule. Proatoms with fractional charges are obtained through linear interpolation between the densities of two proatoms with integer charge. The advantage of these charges is that they have been shown to reproduce the ab initio molecular electrostatic potential 50 and to be conformationally independent. Non-covalent interactions analysis To visualize the interactions between the reactants and the enzyme in the transition state we have used the program NCIPLOT 51,52 on a representative structure extracted from the molecular dynamics simulations. This program calculates the NCI index based on the electron density (ρ(r)) and the reduced density gradient (RDG), defined as:

RDG =

|∇ρ(r)|

1 2(3π 2 )

1 3

4

ρ(r) 3

(3)

Regions of low electron density, like tails located far from the nuclei, show a significant value of RDG. On the other hand, regions of space where a weak interaction is established between two molecular fragments should show a low electron density coupled to a vanishing RDG, due to the presence of a bond critical point. The noncovalent contacts are then identified with the regions of vanishing RDG at low densities. They are located by generating gradient isosurfaces enclosing the corresponding regions in real space. The attractive or repulsive 10


Page 11 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


nature, and also the magnitude of these interactions is identified based on the sign and magnitude of the second eigenvalue of the electron density Hessian matrix. 51 Mean Reaction Force In order to isolate factors involved in the lowering of the activation barrier in the studied reaction, the mean reaction force (MRF) 7,9,53 was used, which divides the free energy profile along a reaction coordinate into four regions. Associated changes in the free energy of each reaction region are mostly related to processes of different nature. 6,7,9,10,53 The first region comprising from reactant to the minimum of the MRF describes mostly structural changes associated with the free energy change W1 (modification of equilibrium length of bonds, changes in the angle between reactants, etc.). W2 is more related to electronic changes as bond rupture and formation occurring close to the transition state without major structural modification. Finally, W3 and W4 describe electronic and structural relaxation from the transition state to the products. Comparison of the free energy involved in each of these regions for different environments yields a more detailed analysis of the enzymatic effect on the reaction mechanism.

Results and Discussion The molecular mechanism by which COMT is able to reduce the activation barrier of the transmethylation reaction with respect to the reaction in aqueous solution was addressed by identifying the Minimum Free Energy Path (MFEP). The reaction along this path was analyzed by the Mean Reaction Force comparing the reaction in aqueous solution with the one in the wild-type enzyme and the less active Y68A mutant. Molecular dynamics simulations of the enzyme were based on the crystal structure of human COMT published in 2008 by Rutherford et al. 35 (see Figure 1). In particular, the dinitrocatechol ligand in the crystal structure was changed to dopamine, leaving only the paraoxygen inside the coordination sphere of magnesium. To keep the coordination number of magnesium constant, one of the carboxylic oxygens of Glu199 interacted with the divalent

11



Figure 1: Representation of the equilibrated reactive complex inside the native (left) and mutated (right) human COMT. The hydration sites within 5˚ A from the carbon and oxygen atoms are shown as water molecules in their most probable orientations.

cation during the initial equilibration process (see Figure 1 left side). This configuration of the active site has been suggested previously 5,22 and is similar to one of the reactant configurations reported recently by Patra et al. using catechol as a substrate. 29 The other reactant conformation of catechol was a bidentate coordination to the magnesium cation but in our study with dopamine this configuration resulted in a disorganization of the WT active site. Therefore, the conformation reported by Roca et al. 5,22 was used to study the reaction of dopamine that resulted in an active site shown in Figure 1 (left). In this conformation the meta-oxygen atom of dopamine is in nucleophilic attack position to the methyl group of SAM in accordance with the reported preferred product. This active site reactant configuration for the WT was stable for 50 ns after which the coordination to the magnesium cation was lost, in agreement with the reactant conformations reported recently for the catechol substrate. 29 To maintain the coordination a distance restraint between the oxygen atom of the DOP hydroxyl group and the magnesium cation was introduced, and five independent 100 ns simulations were performed. Cluster analysis of the five trajectories resulted in rep-

12


Page 12 of 35

Page 13 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


resentative structures which all present the same active site structure shown on the left of Figure 1. For the Y68A mutant the same protocol was used but after the initial equilibration of temperature and pressure water molecules had entered the catalytic cavity (see Figure 1 right side). Additionally, the clustering of the five trajectories provided representative structures, which differed in the conformation of the Lys144 residue in the active site. The most abundant one was similar to the WT as shown on the right of Figure 1 and in the other one Lys144 left the active site and was replaced with additional water molecules (see Figure S1 and S5 in Supporting Information). The representative reactant structures of the WT and the two of the mutant were equilibrated without restraints in DFTB3 QM/MM simulations where the magnesium cation coordination remained constant (see Figure S2-S3).

Free Energy Profile and Mean Reaction Force To study the reaction we traced the minimum free energy path (MFEP) for the transmethylation reaction in the wild-type enzyme (WT), in the Y68A mutant and in aqueous solution (AQ) starting from the structures described above using the adaptive string method implemented in Amber17. 11,38 The activation barriers for the most abundant reactant structures at the DFTB3/ff14SB level are 20.1 kcal mol−1 for the wild type enzyme (WT), 22.5 kcal mol−1 for the Y68A mutant and 26.2 kcal mol−1 for the aqueous system (AQ). One advantage of the Adaptive String Method implementation is the possibility to calculate an uncertainty in the obtained energy profile (shadow area in the Figure 2), which amounts to ±1.0 kcal mol−1 in the activation energy according to a 95% confidence interval. The alternative, less abundant conformation of the mutant without Lys144 in the active site resulted in a slightly smaller activation barrier almost within the error margin of the other conformation (see Figure S5 in Supporting Information). The observed small effect of Lys144 on the MFEP

13



is in agreement with its non-essential role for the methyl-transfer reaction as revealed from experiments showing catalytic activity for the L144A mutant. 54–56 The obtained results have to be compared to experimental activation barriers based on measured rate constants. For the human and murine wild type enzyme with dopamine or dihydroxybenzoic acid the reported values range 18.1-19.1 kcal mol−1 . 25,29,57 Considering our error in the activation energy and the experimental uncertainty the two values lay in the error margin of each other. Recent computational studies have also been able to reach similar values using QM regions that do not include protein residues 25 , and the differences from different QM region size appear to arise when PM6, PM3 or AM1 semi-empirical methods are used. 29,30 With the experimental rate constant of the Y68A mutant an activation energy of 20.0 kcal mol−1 can be estimated. 24 Our calculated activation barrier for the most abundant active site conformation is 22.5 kcal mol−1 , a value between the experimental one and the one calculated by Lameira and Warshel (24.4 kcal mol−1 ). 25 For the experimental activation free energy in aqueous solution often a value of 30 kcal mol−1 is taken as reference, although this value corresponds to the reaction of a model system: phenolate and trimethylsulfonium in methanol. 10 We have shown previously that our methodology yields activation barriers in excellent agreement with experiment for the model system. Because the reactants in this study are different and the initial distance between the reactants is the one found in the enzyme, a comparison to experiment for the uncatalyzed reaction is not straightforward. But, accounting for these differences (model system possesses three methyl groups vs only one in SAM, which requires an entropic correction for activation barrier) the activation barrier found here is in good agreement to our previous result 10 and to the experimental reference. Besides the agreement between the calculated and experimental activation barriers, it is interesting to address different processes that take place during the reaction. The Mean Reaction Force (MRF) naturally divides the whole reaction coordinate in specific regions where different kind of processes take place during the progress of the reaction. 7 The regions

14


Page 14 of 35

Page 15 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


are defined by the derivative of the free energy profile with respect to the reaction coordinate. The minimum, zero, and maximum of the MRF generate the four regions. The activation barrier can be thus separated in an early contribution W1 , which includes the free energy required mainly for the structural changes of the system and a late contribution W2 , which involves changes associated to electronic rearrangement in the proximity of the transition state.

Figure 2: Left: Free Energy Profile, Center: Mean Reaction Force and Right: behavior of collective variables CV1 (continuous line) and CV2 (dashed line) as a function of the path CV for the reaction of methyl transfer in WT (red), Y68A (green) and AQ (blue) environments. Vertical dashed lines represent the minimum, zero and maximum values of the mean reaction force for the WT system. The shaded area shows the 95% confidence interval in the free energy.

Figure 2 (left and center), shows the free energy profile and mean reaction force obtained for the WT (red) and the Y68A (green) mutant enzyme starting from the conformation with the negatively charged, nucleophilic oxygen atom pointing towards SAM (see Figure 1), and in aqueous solution (AQ, blue) starting from the distance found in the native enzyme. Interestingly, the transition state and the position of the MRF extrema are at smaller values of the reaction coordinate in aqueous solution than in the wild type or mutant enzyme. One possible explanation for this difference could be that the reaction in the three systems follows a different path reflected in different distances between the atoms involved in the methyl transfer. However, as can be seen in the right side of Figure 2, the S-C and C-O distances (which were defined as collective variables previously) follow the same trend in the three environments along the MFEP, with a very small variation in the transition state region between aqueous solution and the WT or mutant enzyme. During the first steps, the C-O distance decreases, while the S-C distance remains constant. In the transition state 15



region both distances change and at the end of the reaction the methylated product moves away from SAM. Therefore, the observed difference in the transition state and the MRF extrema between the reaction in water and in the enzyme has to be related to the interactions of the reactants with the surrounding molecules in each environment and the resulting changes in the free energy. These interactions may alter the regions defined through the MRF extrema and their contribution to the activation barrier differently as is shown in Table 1. The energetic difference in W1 and W2 between the WT and AQ systems amounts to approximately 3 kcal mol−1 respectively. This difference for W1 is related to more favorable structural changes in the enzyme due to e.g. preorganization of the WT enzyme, while for W2 the presence of residues stabilizing the less polar (compared to the reactants) transition state might decrease the energy required for the electronic rearrangements. From this first analysis, we can conclude that the wild type enzyme reduces both contributions equally. Table 1 indicates that for the mutant system the contribution to the activation free energy from W1 amounts to 12.1 kcal mol−1 , 1.7 kcal mol−1 larger than the same process in the wild type enzyme, whereas the difference in W2 to the wild type enzyme is relatively small (0.7 kcal mol−1 ). Therefore, one can conclude that the mutant presents similar favorable interactions in the transition state, which reduces the contribution of W2 with respect to aqueous solution. But, the contribution W1 from structural changes is similar to the AQ system and differs considerably from the wild type enzyme. The mutation, therefore, alters mostly the early steps of the reaction, making them energetically similar to the reaction occurring in an aqueous environment. Table 1: Free energy associated to each region defined by the MRF and the activation free energy for the reaction in the three different environments. Units are kcal mol−1 W3 W4 Ea System W1 W2 WT 10.4 9.7 -10.9 -15.1 20.1 Y68A 12.1 10.4 -14.7 -12.4 22.5 AQ 13.3 12.9 -12.1 -19.1 26.2 . 16


Page 16 of 35

Page 17 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


In global terms, the Mean Reaction Force has shown that most of the changes between the mutant and the WT enzyme are located in the first region characterized by the contribution W1 and the structural changes. The catalytic effect of the enzyme is present in the first region (W1 ) and also the transition state region (W2 ), where the AQ system presents a larger contribution to W2 and slightly shifted positions in MRF extrema and transition state. An interesting aspect of the string method is that the contribution of each collective coordinate to the total free energy barrier can be quantified along the MFEP. 11,58 Since the MFEP was defined in a three dimensional space formed by the S-C distance (CV1), C-O distance (CV2) and the hybridization of the methyl group (CV3), the free energy can be decomposed as the sum of individual contributions of free energy calculated for each collective coordinate. The contributions of each CV are shown in Figure 3 for each environment. The largest differences between the enzyme and the aqueous solution are observed in CV1 (left) and CV2 (middle) defined as the S-C distance and the C-O distance respectively. However, CV1 presents the same profile in the wild type enzyme and the mutant with a pronounced difference to the aqueous solution (at smaller scale this also holds for CV3). CV2 in contrast differ throughout the three systems, and the mutant contribution lies in between the one of the enzyme and in water.

Figure 3: Free energy profile for each CV as a function of the path CV. Left: CV1 (S-C distance), center: CV2 (C-O distance), right: CV3 (Methyl hybridization). WT (red), Y68A (green) and AQ (blue) systems.

In the region of structural changes up to the minimum of the MRF (dashed vertical line) the contribution of CV1 to the activation free energy in the enzyme (mutant and wild-type) and in aqueous solution is absent till close before the minimum. But, CV2’s contribution starts 17



earlier and differs between the three environments markedly with the largest contribution in aqueous solution followed by the mutant and the wild-type enzyme. Also CV3 present some changes but they are very small energetically compared to the others. For the transition state region a similar picture emerges where the main differences between the three environments are concentrated in CV2 and the difference between aqueous environment and enzyme is reflected in CV1 and CV3. This observation is in agreement with the result of W1 and W2 from the MRF, as CV2 describes the approximation of the reactants which is more structural in nature whereas CV1 involves electronic changes of bond breaking and bond formation. Following the scheme of the Mean Reaction Force, structural and electronic contributions of each collective variable were used to assign their contribution to the total activation free energy for each environment (see Table S1 in the Supplementary Information). The contributions to W1 are mostly concentrated in CV2, associated with the approximation of the nucleophilic oxygen to the methyl group. The wild type enzyme facilitates this initial step reducing the value by 2.1 kcal mol−1 with respect to the aqueous solution. The mutant, however, present a value that differs only by 0.7 kcal mol−1 from the reaction in water. This confirms that the mutant and the aqueous solution are energetically similar in the structural changes in the beginning of the reaction. Differences in CV1 and CV3 in the structural contributions to the activation barrier are largest between WT and AQ systems in CV1 (0.9 kcal mol−1 ) but overall significantly smaller for the other environments and also smaller compared to the contribution of CV2. This changes when the region between the minimum of the MRF and the transition state is analyzed. In this region associated with electronic rearrangements, breaking and forming of chemical bonds, CV1 (S-C distance) presents the largest contribution to W2 with a difference of the WT system of 1.8 kcal mol−1 to the aqueous solution. For CV1 and CV2 the contribution of the mutant and the wild-type are very similar to each other and differ considerably from the values in aqueous solution (see Table S2 in the Supplementary Information). From this first analysis we can conclude that the enzyme (mutant or wild-type) reduces the

18


Page 18 of 35

Page 19 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


free energy required for the electronic rearrangements in the bond forming and breaking processes considerably. The effect on the structural changes is also evident, but only for the more active wild type enzyme, whereas the mutant presents free energy contributions more similar to the aqueous solution, possibly due to the water molecules entering the catalytic cavity. It is in the first region of the MRF represented with the changes of CV2 (C–O distance) where the mutant differs most from the wild type enzyme. To see if the electronic nature of the oxygen atom differs between the three environments, we studied the Hirshfeld-I atomic charges on sulfur and oxygen atoms involved in the reaction, as done in our previous study 10 . Hirshfeld-I charges are based on the polarized electron density 48 (QM/MM methodology) assigning different portions of the molecular electron density to different atoms in an Atoms in Molecules approach. The advantage over other atomic charges is that they don’t depend on the employed basis set of the electronic structure method and that they reproduce the molecular electrostatic potential accurately. In Figure 4 on the left the charge of the sulfur atom in each environment is shown. In all cases the charge remains constant in the beginning of the reaction until the reactants are properly aligned for the electronic changes to occur. The constant charge is related to the invariance of the distance S-C. During this first approach of the two molecules the charge on the nucleophilic oxygen (Figure 4 right) was subject to a steady increase. These results confirm that little electron transfer and therefore electronic changes occur in this region before the minimum of the MRF. Comparing the absolute values of each charge in this region in the different environments one observes that the charge on the sulfur atom is similar in both enzymes, but differs from the charge in aqueous solution. This suggests that in the beginning of the reaction, the electrostatic environment of SAM is very similar in the wild-type and the mutant enzyme but different to the one found in aqueous solution. The similarity of the enzymatic environments is reinforced by the equal binding isotope effects observed previously for the WT and Y68A enzymes. 26

19



0.2

Sulfur

0.8 0.6

Charge[e]

0.5 0.4 0.3 0.2

0.4 0.5 0.6

0.1 0.0

Oxygen

0.3

0.7 Charge[e]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 35

0

2

4

6 s[amu1/2 Å]

8

10

0.7

12

0

2

4

6 s[amu1/2 Å]

8

10

12

Figure 4: Atomic charge for atoms involved in the methyl transfer as a function of the path CV. Left: sulfur atom, right: oxygen atom. WT (red), Y68A (green) and AQ (blue).

But for the nucleophilic oxygen atom in dopamine the similarity in atomic charge is found between the mutant and the aqueous system, evidencing a different electrostatic environment for this atom in the wild-type enzyme. The similar charge of the oxygen atom in the mutant and in aqueous solution is related to water molecules entering the catalytic cavity of the mutant enzyme. If an electronegative atom is surrounded by water molecules one expects a more negative atomic charge. In the wild type enzyme a reduced number of water molecules would lead to a less negative oxygen atom. To quantify the number and positions of water molecules in the active site of the two enzymes, hydration sites for the two most abundant representative reactant structures from the MM simulation obtained by the SSTMap software are shown in Figure 5. In the WT enzyme only four hydration sites are found whereas in the mutant up to eight sites are identified. The effect of more water molecules in active sites has been discussed recently and was associated with better solvation of charged Michaelis complexes and an increase in the reorganization energy to reach an apolar transition state. Both factors would lead to larger activation barriers, assuming a negligible effect of the water molecules on the transition state. 59–61 In fact in the recent study by Lameira et al. it was found that the reactant complex (for the catechol substrate) is better solvated in the mutant than in the wild type increasing the activation barrier, which is in agreement with our results. 20


Page 21 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


This different molecular environment between the mutant and wild type enzymes also explains the difference in the contribution of CV2 to the structural part of the activation barrier discussed above and the larger value of W1 for the mutant and AQ system. For the sulfur atomic charge, which starts changing close to the transition state region, the mutant and the wild type enzyme present the same evolution, which differs from the AQ system.

Figure 5: Position of water molecules within 5 ˚ A of the reactive oxygen and carbon atoms observed in the hydration site analysis of the most representative reactant conformation. Left: wild-type enzyme, Right: Y68A mutant.

Considering the similarity of the Hirshfeld-I charges on sulfur in the enzyme and the differences observed on the oxygen atomic charge between the three environments together with the analysis of the hydration sites, we conclude that the approach of the reactants requires more free energy in the mutant and in AQ due to the increased desolvation cost of the polar reactant Michaelis complex to reach the less polar transition state. The methyl group transfer reaction of SAM to DOP is directly associated with the transfer of one positive charge. It is in the transition state region where this transfer is mostly concentrated. However, the charges on the oxygen and sulfur atom do not vary between the wild-type enzyme, the mutant and the uncatalyzed reaction. After the maximum of 21



the MRF the change in the atomic charges is small with minimal differences between the molecular environments evidencing a structural relaxation of the products. Although the electron transfer reflected in the charges seems to progress in the same manner, the free energy cost related to this process and especially the contribution of the collective variables to the total free energy is different as discussed above. The decomposition of the free energy in the CVs showed that CV1 has the largest contribution to the electronic rearrangements characterized through W2 and that in this region the aqueous solution (AQ) requires more free energy than in the wild type or mutant enzyme. The differences between the enzyme and aqueous solution are also reflected in the contribution of CV3, although by a much smaller absolute value. This might indicate that the inversion of the methyl group is easier inside the enzymatic environment than in water. But, specific interactions favoring this process in the enzyme should be identified to confirm the catalysis of the inversion process considering the small differences. Recent studies proposed an additional contribution to catalysis of this type of reactions stemming from interactions between methylic hydrogen atoms and protein or water oxygen atoms. 62,63 We examined the active site for residues able to form this type of interaction at this stage of the reaction and assist the inversion of the methyl group. The most favorable interaction is observed between the oxygen atom in the backbone carbonyl group of Asp141 and the hydrogen atoms in the methyl group. This interaction is then compared to the average distance of the oxygen atoms in the water molecules of the aqueous system in Figure 6. In the transition state region of the the wild-type enzyme the distance between oxygen and methylic hydrogen is clearly shorter than the distance found in water (between the dashed lines). For the Y68A mutant this distance starts at levels comparable to the aqueous system and then shortens in the W2 region, probably explaining, along with other electrostatic interactions, the progressive reduction of the MRF in the mutant system and its overlap with the wild-type one (see Figure 2). This type of interactions can assist the hybridization of the methyl group in the transition

22


Page 22 of 35

Page 23 of 35

state region within the enzymes and could also provide an alternative explanation for the experimentally observed inverse 2◦ -KIE difference between wild type, mutant and aqueous solution for this reaction according to Williams et al. 63 To confirm the nature of this interaction in Figure 6 we show the noncovalent interactions (NCI) in one representative structure at the transition state of the wild-type enzyme. Despite the characteristic rings around the two bonds that are broken and formed, stabilizing interactions with backbone oxygen atom of Asp141 and the sulfur atom of Met40 are observed.

3.0 2.9 2.8 Distance[Å]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.7 2.6 2.5 2.4 2.3

0

2

4

6 s[amu1/2 Å]

8

10

12

Figure 6: Left: Average minimal distance along the reaction coordinate between hydrogens in the methyl group and backbone oxygen atoms in Asp141, for WT (red) and Y68A (green), and water oxygen atoms for AQ (blue). Right: Illustration of non-covalent interactions present at the transition state inside the active site of the enzymatic WT system. Red: repulsive/steric interactions, blue: charge-charge interactions, green: hydrogen bonds. The discontinuous blue circles highlight two stabilizing interactions. Some hydrogen atoms were omitted for clarity.

From this analysis one can conclude that the enzyme can assist the hybridization of the methyl group through specific interactions with the oxygen backbone atom of Asp141. These interactions accompanied with an appropriate electrostatic environment reduce the free energy for the electronic rearrangements characterized by W2 and the contribution of CV1 associated with the S-C distance.

23



Conclusions Through a combination of methods we studied the factors which contribute to the catalysis of catechol O-methyltransferase in the methyl transfer reaction from SAM to dopamine, its natural substrate. The adaptive string method identified the minimum free energy path from the reactants to the products in the three environments: wild-type (WT) enzyme, less active Y68A mutant, and in aqueous solution (AQ). The obtained free energy barriers reproduce the relative reactivity among environments and are in good agreement with the experimental values and previous computational studies. Analysis of the inter-atomic distances does not reveal a different MFEP in the three environments and suggests the interactions between the substrate and the surrounding molecular environment as the key factor in the catalysis. The Mean Reaction Force (MRF) together with a combination of structural and electronic descriptors were used to identify the main factors contributing to catalysis. From this analysis one can conclude that the wild-type enzyme reduces the free energy required to fulfill the structural changes, particularly facilitating the approach of the two reactants. This does not hold for the less effective Y68A mutant, which presents more water molecules in the catalytic cavity, which increase the free energy cost of the corresponding structural changes due to the increased desolvation free energy cost of the polar Michaelis complex. This effect on the activation barrier was quantified by the free energy contribution of each collective variable to the structural contribution of the free energy barrier together with the Hirshfeld-I atomic charges on the oxygen atom. The evolution of the atomic charges also allowed to study the electron transfer process in the transition state region between the extrema of the MRF. For these electronic changes the mutant and wild type enzyme provide a favorable electrostatic environment, which reduces the free energy required for the electron transfer. Additionally, specific interactions between the oxygen backbone atom of Asp141 and hydrogen atoms of the transferred methyl group assist the inversion of this methyl group (see Figure 6). In summary, the WT enzyme facilitates the structural changes to reach the transition state region through the elimination of water molecules and assists the electron transfer and the 24


Page 24 of 35

Page 25 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


inversion process of the methyl group through favorable interactions at the transition state. The decreased activity of the mutant is related to water molecules entering the catalytic cavity, which make the approach of the two charged substrates more difficult. This might also explain the observed trend in 2◦ -KIE and the less structured cavity reported by Zhang et al. 26 Furthermore, the less structured cavity from the experimental study is in agreement with our simulations which showed increased fluctuations of the residues in the active site for the Y68A mutant enzyme (see change in root mean square fluctuations (RMSF) per residue between the two enzymes shown in Figure S6 of the Supporting Information). From this first combination of the adaptive string method and the MRF new insights in the catalytic mechanism of the wild type and mutant COMT has been provided. We hope that this approach will be helpful in the future for better understanding of enzymatic catalysis and rational enzyme engineering.

Supporting Information Figure S1-3 displays important distances in the active site in the equilibration of the enzyme, Figure S4 the charge of the methyl group in the reaction, Figure S5 the less abundant reactant conformation of the mutant and Figure S6 the change in RMSF per residue for the WT and Y68A mutant. Table S1-2 summarizes the structural and electronic contributions of the collective variables CV1 and CV2.

Acknowledgement DAS acknowledges a PhD scholarship 21130517 from CONICYT. EVM is thankful for funding support provided by FONDECYT 1160197 and Max-Planck Society for funding the Max-Planck Partner group. Fernanda Duarte is acknowledged for fruitful discussion on the NCI calculations. KZ and IT gratefully acknowledge financial support from Ministerio de Econom´ıa y Competitividad (project CTQ2015-66223-C2-2-P. KZ acknowledges a FPU fel-

25



lowship from Ministerio de Educación. The authors acknowledge computational facilities of the Servei d Informática de la Universitat de València in the “Tirant” supercomputer and Maite Roca for providing the initial dopamine-SAM-enzyme structure.

26


Page 26 of 35

Page 27 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


References (1) Hammes-Schiffer, S. Catalytic Efficiency of Enzymes: A Theoretical Analysis. Biochemistry 2013, 52, 2012–2020. (2) Warshel, A.; Levitt, M. Theoretical Studies of Enzymic Reactions: Dielectric, Electrostatic and Steric Stabilization of the Carbonium Ion in the Reaction of Lysozyme. J. Mol. Biol. 1976, 103, 227–249. (3) Warshel, A. Electrostatic Origin of the Catalytic Power of Enzymes and the Role of Preorganized Active Sites. J. Biol. Chem. 1998, 273, 27035–27038. (4) Jindal, G.; Warshel, A. Exploring the Dependence of QM-MM Calculations of Enzyme Catalysis on the Size of the QM Region. J. Phys. Chem. B 2016, 120, 9913–9921. (5) Roca, M.; Andrés, J.; Moliner, V.; Tu˜ nón, I.; Bertrán, J. On the Nature of the Transition State in Catechol O-Methyltransferase. A Complementary Study Based on Molecular Dynamics and Potential Energy Surface Explorations. J. Am. Chem. Soc. 2005, 127, 10648–10655. (6) Vöhringer-Martinez, E.; Dörner, C. Conformational Substrate Selection Contributes to the Enzymatic Catalytic Reaction Mechanism of Pin1. J. Phys. Chem. B 2016, 120, 12444–12453. (7) Vöhringer-Martinez, E.; Toro-Labbé, A. The Mean Reaction Force: A Method to Study the Influence of the Environment on Reaction Mechanisms. J. Chem. Phys. 2011, 135, 064505. (8) Vöhringer-Martinez, E.; Verstraelen, T.; Ayers, P. W. The Influence of Ser-154, Cys113, and the Phosphorylated Threonine Residue on the Catalytic Reaction Mechanism of Pin1. J. Phys. Chem. B 2014, 118, 9871–9880.

27



(9) Vöhringer-Martinez, E.; Duarte, F.; Toro-Labbé, A. How does Pin1 Catalyze the CisTrans Prolyl Peptide Bond Isomerization? A QM/MM and Mean Reaction Force Study. J. Phys. Chem. B. 2012, 116, 12972–12979. (10) Saez, D. A.; Vogt-Geisse, S.; Inostroza-Rivera, R.; Kubar, T.; Elstner, M.; ToroLabbé, A.; Vöhringer-Martinez, E. The Effect of the Environment on the Methyl Transfer Reaction Mechanism Between Trimethylsulfonium and Phenolate. Phys. Chem. Chem. Phys. 2016, 18, 24033–24042. (11) Zinovjev, K.; Tu˜ nón, I. Adaptive Finite Temperature String Method in Collective Variables. J. Phys. Chem. A 2017, 121, 9764–9772. (12) Kozbial, P. Z.; Mushegian, A. R. Natural History of S-Adenosylmethionine-Binding Proteins. BMC Struct. Biol. 2005, 5, 1–26. (13) Boriack-Sjodin, P. A.; Swinger, K. K. Protein Methyltransferases: A Distinct, Diverse, and Dynamic Family of Enzymes. Biochemistry 2016, 55, 1557–1569. (14) Kiss, L. E.; da Silva, P. S. Medicinal Chemistry of Catechol O-Methyltransferase (COMT) Inhibitors and Their Therapeutic Utility. J. Med. Chem. 2014, 57, 8692– 8717. (15) Mihel, I.; Knipe, J. O.; Coward, J. K.; Schowen, R. L. α-Deuterium Isotope Effects and Transition-State Structure in an Intramolecular Model System for Methyl-Transfer Enzymes. J. Am. Chem. Soc. 1979, 101, 4349–4351. (16) Hegazi, M. F.; Borchardt, R. T.; Schowen, R. L. α-Deuterium and Carbon-13 Isotope Effects for Methyl Transfer Catalyzed by Catechol O-Methyltransferase. SN 2-like Transition State. J. Am. Chem. Soc. 1979, 101, 4359–4365. (17) Gray, C. H.; Coward, J. K.; Schowen, K. B.; Schowen, R. L. α-Deuterium and Carbon13 Isotope Effects for a Simple, Intermolecular Sulfur-to-Oxygen Methyl-Transfer 28


Page 28 of 35

Page 29 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Reaction. Transition-State Structures and Isotope Effects in Transmethylation and Transalkylation. J. Am. Chem. Soc. 1979, 101, 4351–4358. (18) Swain, C. G.; Taylor, L. J. Catalysis by Secondary Valence Forces. J. Am. Chem. Soc. 1962, 84, 2456–2457. (19) Bruice, T. C.; Zheng, Y.-J. A Theoretical Examination of the Factors Controlling the Catalytic Eficiency of a Transmethylation Enzyme: Catechol O-Methyltransferase. J. Am. Chem. Soc. 1997, 119, 8137–8145. (20) Bruice, T. C.; Lau, E. Y. Importance of Correlated Motions in Forming Highly Reactive Near Attack Conformations in Catechol O-Methyltransferase. J. Am. Chem. Soc. 1998, 120, 12387–12394. (21) Bruice, T. C.; Lau, E. Y. Comparison of the Dynamics for Ground-State and TransitionState Structures in the Active Site of Catechol O-Methyltransferase. J. Am. Chem. Soc. 2000, 122, 7165–7171. (22) Roca, M.; Mart´ı, S.; Andrés, J.; Moliner, V.; Tu˜ noń, I.; Bertrán, J.; Williams, I. H. Theoretical Modeling of Enzyme Catalytic Power: Analysis of Cratic and Electrostatic Factors in Catechol O-Methyltransferase. J. Am. Chem. Soc. 2003, 125, 7726–7737. (23) Garc´ıa-Meseguer, R.; Zinovjev, K.; Roca, M.; Ruiz-Pern´ıa, J. J.; Tu˜ nón, I. Linking Electrostatic Effects and Protein Motions in Enzymatic Catalysis. A Theoretical Analysis of Catechol O-Methyltransferase. J. Phys. Chem. B 2015, 119, 873–882. (24) Zhang, J.; Klinman, J. P. Enzymathic Methyl Transfer: Role of an Active Site Residue in Generating Active Site Compaction that Correlates with Catalytic Efficiency. J. Am. Chem. Soc. 2011, 133, 17134–17137.

29



(25) Lameira, J.; Bora, R. P.; Chu, Z. T.; Warshel, A. Methyltransferases do not Work by Compression, Cratic, or Desolvation Effects, but by Electrostatic Preorganization. Proteins: Struct., Funct., Bioinf. 2015, 83, 318–330. (26) Zhang, J.; Kulik, H. J.; Martinez, T. J.; Klinman, J. P. Mediation of DonorAcceptor Distance in an Enzymatic Methyl Transfer Reaction. Proc. Natl. Acad. Sci. U. S. A. 2015, 112, 7954–7959. (27) Kulik, H. J.; Zhang, J.; Klinman, J. P.; Martinez, T. J. How Large Should the QM Region Be in QM/MM Calculations? The Case of Catechol O-Methyltransferase. J. Phys. Chem. B 2016, 120, 11381–11393. (28) Karelina, M.; Kulik, H. J. Systematic Quantum Mechanical Region Determination in QM/MM Simulation. J. Chem. Theory Comput. 2017, 13, 563–576. (29) Patra, N.; Ioannidis, E. I.; Kulik, H. J. Computational Investigation of the Interplay of Substrate Positioning and Reactivity in Catechol O-Methyltransferase. PLoS One 2016, 11, 1–23. (30) Kulik, H. J. Large-Scale QM/MM Free Energy Simulations of Enzyme Catalysis Reveal the Influence of Charge Transfer. Phys. Chem. Chem. Phys. 2018, 20, 20650–20660. (31) Saez, D. A.; Vöhringer-Martinez, E. A Consistent S-Adenosylmethionine Force Field Improved by Dynamic Hirshfeld-I Atomic Charges for Biomolecular Simulation. J. Comput.-Aided Mol. Des. 2015, 29, 951–961. (32) Maragliano, L.; Fischer, A.; Vanden-Eijnden, E.; Ciccotti, G. String Method in Collective Variables: Minimum Free Energy Paths and Isocommittor Surfaces. J. Chem. Phys. 2006, 125, 024106.

30


Page 30 of 35

Page 31 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(33) Zinovjev, K.; Ruiz-Pern´ıa, J. J.; Tu˜ nón, I. Toward an Automatic Determination of Enzymatic Reaction Mechanisms and Their Activation Free Energies. J. Chem. Theory Comput. 2013, 9, 3740–3749. (34) Vidgren, J.; Svensson, L.; Liljas, A. Crystal Structure of Catechol O-Methyltransferase. Nature 1994, 368, 354–358. (35) Rutherford, K.; Le Torng, I.; Stenkamp, R. E.; Parson, W. W. Crystal Structures of Human 108V and 108M Catechol O-Methyltransferase. J. Mol. Biol. 2008, 380, 120– 130. (36) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General Amber Force Field. J. Comput. Chem. 2004, 25, 1157–1174. (37) Maier, J. A.; Martinez, C.; Kasavajhala, K.; Wickstrom, L.; Hauser, K. E.; Simmerling, C. ff14SB: Improving the Accuracy of Protein Side Chain and Backbone Parameters from ff99SB. J. Chem. Theory Comput. 2015, 11, 3696–3713. (38) Case, D.; Cerutti, D.; Cheatham III, T.; Darden, T.; Duke, R.; Giese, T.; Gohlke, H.; Goetz, A.; Greene, D.; Homeyer, N. et al. Amber 2017. University of California, San Francisco, 2017. (39) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684– 3690. (40) Haider, K.; Cruz, A.; Ramsey, S.; Gilson, M. K.; Kurtzman, T. Solvation Structure and Thermodynamic Mapping (SSTMap): An Open-Source, Flexible Package for the Analysis of Water in Molecular Dynamics Trajectories. J. Chem. Theory Comput. 2018, 14, 418–425.

31



(41) Gaus, M.; Cui, Q.; Elstner, M. DFTB3: Extension of the Self-Consistent-Charge Density-Functional Tight-Binding Method (SCC-DFTB). J. Chem. Theory Comput. 2011, 7, 931–948. (42) Gaus, M.; Goez, A.; Elstner, M. Parametrization and Benchmark of DFTB3 for Organic Molecules. J. Chem. Theory Comput. 2013, 9, 338–354. (43) Gaus, M.; Lu, X.; Elstner, M.; Cui, Q. Parameterization of DFTB3/3OB for Sulfur and Phosphorus for Chemical and Biological Applications. J. Chem. Theory Comput. 2014, 10, 1518–1537. (44) Lu, X.; Gaus, M.; Elstner, M.; Cui, Q. Parametrization of DFTB3/3OB for Magnesium and Zinc for Chemical and Biological Applications. J. Phys. Chem. B 2015, 119, 1062– 1082. (45) Kastner, J. Umbrella Integration in Two or More Reaction Coordinates. J. Chem. Phys. 2009, 131, 034109. (46) Neese, F. The ORCA Program System. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 73–78. (47) Verstraelen, T.; Tecmer, P.; Heidar-Zadeh, F.; Boguslawski, K.; Chan, M.; Zhao, Y.; Kim, T. D.; Vandenbrande, S.; Yang, D.; González-Espinoza, C. E. et al. HORTON 2.0.0. 2015; http://theochem.github.com/horton/, (accessed March 7, 2018). (48) Bultinck, P.; Van Alsenoy, C.; Ayers, P. W.; Carbó-Dorca, R. Critical Analysis and Extension of the Hirshfeld Atoms in Molecules. J. Chem. Phys. 2007, 126, 144111. (49) Hirshfeld, F. L. Bonded-Atom Fragments for Describing Molecular Charge Densities. Theor. Chim. Acta 1977, 44, 129–138. (50) Van Damme, S.; Bultinck, P.; Fias, S. Electrostatic Potentials from Self-Consistent Hirshfeld Atomic Charges. J. Chem. Theory Comput. 2009, 5, 334–340. 32


Page 32 of 35

Page 33 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(51) Johnson, E. R.; Keinan, S.; Mori-Snchez, P.; Contreras-Garca, J.; Cohen, A. J.; Yang, W. Revealing Noncovalent Interactions. J. Am. Chem. Soc. 2010, 132, 6498– 6506. (52) Contreras-Garc´ıa, J.; Johnson, E. R.; Keinan, S.; Chaudret, R.; Piquemal, J.-P.; Beratan, D. N.; Yang, W. NCIPLOT: A Program for Plotting Noncovalent Interaction Regions. J. Chem. Theory Comput. 2011, 7, 625–632. (53) Toro-Labbé, A. Characterization of Chemical Reactions from the Profiles of Energy, Chemical Potential and Hardness. J. Phys. Chem. A. 1999, 103, 4398–4403. (54) Law, B. J. C.; Bennett, M. R.; Thompson, M. L.; Levy, C.; Shepherd, S. A.; Leys, D.; Micklefield, J. Effects of Active-Site Modification and Quaternary Structure on the Regioselectivity of Catechol-O-Methyltransferase. Angew. Chem., Int. Ed. 2016, 55, 2683–2687. (55) Ehler, A.; Benz, J.; Schlatter, D.; Rudolph, M. G. Mapping the Conformational Space Accessible to Catechol-O-Methyltransferase. Acta Crystallogr., Sect. D: Struct. Biol. 2014, 70, 2163–2174. (56) Männistö, P. T.; Kaakkola, S. Catechol-O-methyltransferase (COMT): Biochemistry, Molecular Biology, Pharmacology, and Clinical Efficacy of the New Selective COMT Inhibitors. Pharmacol. Rev. 1999, 51, 593–628. (57) Schultz, E.; Nissinen, E. Inhibition of Rat Liver and Duodenum Soluble Catechol-OMethyltransferase by a Tight-Binding Inhibitor OR-462. Biochem. Pharmacol. 1989, 38, 3953–3956. (58) Sanchez-Martinez, M.; Field, M.; Crehuet, R. Enzymatic Minimum Free Energy Path Calculations Using Swarms of Trajectories. J. Phys. Chem. B 2015, 119, 1103–1113.

33



(59) Richard, J. P.; Amyes, T. L.; Goryanova, B.; Zhai, X. Enzyme Architecture: On the Importance of Being in a Protein Cage. Curr. Opin. Chem. Biol. 2014, 21, 1 – 10. (60) Kulkarni, Y. S.; Liao, Q.; Petrovi, D.; Krger, D. M.; Strodel, B.; Amyes, T. L.; Richard, J. P.; Kamerlin, S. C. L. Enzyme Architecture: Modeling the Operation of a Hydrophobic Clamp in Catalysis by Triosephosphate Isomerase. J. Am. Chem. Soc. 2017, 139, 10514–10525. (61) Blaha-Nelson, D.; Kr¨ uger, D. M.; Szeler, K.; Ben-David, M.; Kamerlin, S. C. L. Active Site Hydrophobicity and the Convergent Evolution of Paraoxonase Activity in Structurally Divergent Enzymes: The Case of Serum Paraoxonase 1. J. Am. Chem. Soc. 2017, 139, 1155–1167. (62) Horowitz, S.; Dirk, L. M.; Yesselman, J. D.; Nimtz, J. S.; Adhikari, U.; Mehl, R. A.; Scheiner, S.; Houtz, R. L.; Al-Hashimi, H. M.; Trievel, R. C. Conservation and Functional Importance of Carbon-Oxygen Hydrogen Bonding in AdoMet-Dependent Methyltransferases. J. Am. Chem. Soc. 2013, 135, 15536–15548. (63) Wilson, P. B.; Williams, I. H. Influence of Equatorial CHO Interactions on Secondary Kinetic Isotope Effects for Methyl Transfer. Angew. Chem., Int. Ed. 2016, 55, 3192– 3195.

34


Page 34 of 35

Page 35 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Graphical TOC Entry

35


Catalytic Reaction Mechanism in Native and Mutant Catechol-O

Recommend Documents