A Computational Study of the Promiscuity of the SAM-Dependent Methyltransferase AtHTMT1

A two-pronged computational approach was taken to study the promiscuity of the SAM+-dependent methyl transferase AtHTMT1 from thale cress with several nucleophiles (Cl–, Br–, I–, NCO–, NCS–). First, enzyme-free methyl transfer reactions were studied with M05/6-311+G(2d,p) DFT calculations and electrostatic continuum models (PCM/SMD) for various chemical environments. Second, QM/MM MD simulations with semiempirical Hamiltonians (PM7, PM6-D3, AM1, PM6-D3H4) and the AMBER 14SB force field were used to study the enzyme catalyzed reaction in silico. The combination of the DFT and MD results shows that reactant desolvation generally accelerates the reaction, but it cannot explain the selectivity of the enzyme. The critical position of H2O molecules at the reactive site favors the reaction of NCS– over Cl– and Br– in agreement with experiments, but not observed in the quantum calculations for the cytosol. The addition of selected H2O molecules to the N terminus of NCS– greatly increases its reactivity, while H2O molecules attached to Cl– slow the reaction. The partial solvation of the nucleophiles in the reactive pouch holds the key to understanding the reactivity of AtHTMT1.


INTRODUCTION
S-adenosyl-L-methionine (SAM)-dependent halide/thiocyanate methyl transferases (HTMTs) catalyze the transfer of a methyl group from the sulfonium cation in SAM + to a nucleophile X − yielding CH 3 X and S-adenosyl-L-homocysteine (SAH) as summarized in eq 1. 1−3 The promiscuity of HTMTs, which accept a wide range of nucleophiles, has been successfully exploited to replace the CH 3 group at the sulfonium cation with other residues R to produce SAM analogues R−SAH + , (2) which then can be used to selectively mark targeted biomolecules with other methyl transferases (MTs) as possible medical diagnostics. 4−7 This versatile promiscuity is not only a useful lab-property of HTMTs. It can also be exploited for the green production of methyl halides from nonfood biomass as valuable industrial and agricultural chemicals. 8 The production of CH 3 Cl can be observed on a large scale in vivo. Methyl chloride is a significant source of atmospheric chlorine (17%), but its emission is not controlled, as about 55% of it is produced by plants. 9−11 The HTMTs linked to the plantal emission of CH 3 Cl typically accept a wide range of nucleophiles, such as SH − , CN − , SCN − , Cl − , Br − , and I − . 12−17 Methyl thiocyanate has been found to be part of the plantal defense against bacterial infections, 18 but no such biological benefit has been found for CH 3 Cl so far.
The analysis of the promiscuity of HTMTs and their catalytic mechanism hold the key to understanding plantal CH 3 Cl production. It has been speculated that the emission of CH 3 Cl can be used to remove excess Cl − from the plant, but SAM is an energetically expensive commodity for the cell. 13,19 The quest of Cl − removal consuming SAM appears to be energetically prohibitively expensive and can be suppressed by feeding NaSCN to the plant. 20 The study of the enzymatic mechanism provides the answer to the question of whether the production of CH 3 Cl is a simple, though expensive, byproduct of the NCS − metabolism.
The kinetic and structural data for the HTMT from thale cress (Arabidopsis thalina, AtHTMT1) lay the foundation for the computational study of HTMT promiscuity with kinetic and structural data from a single source. 14 It has been reported that AtHTMT1 works most efficiently with the NCS − anion in the following order: NCS − > I − > Br − > Cl − . (The minus sign "−" in ambidentate ligands marks the ligand atom later binding directly to the CH 3 group from the SAM molecule.) The reproduction of this order will be used to judge the quality of the computational model. Simulations agreeing with the experiment then can be used to deduce the physical principles controlling the promiscuity of AtHTMT1.
There are theoretical studies on the chemistry of SAM, 21−32 but none of them focuses on AtHTMT1 and its reactivity with small nucleophiles. In this work, starting from the PDB structure, the free energy changes in the rate-determining step (rds) of the formation of CH 3 Cl, CH 3 Br, CH 3 I, CH 3 OCN, CH 3 NCO, CH 3 NCS, and CH 3 SCN are studied with QM/MM MD simulations to explore the details of the underlying mechanism. 33 To quantify which nucleophile benefits the most from the catalytic environment, the MD simulations including AtHTMT1 are compared with enzyme-free QM calculations of the methyl transfer with an electrostatic continuum model using water as solvent to mimic the uncatalyzed reaction in the cytosol.
Comparing results obtained from QM/MM MD simulations with those from DFT calculations with an SMD model for solvent appears to be counterintuitive as the two methods differ dramatically in the number of atoms involved. Nevertheless, both methods strive to reproduce chemical reality starting from different assumptions of the solvation process. The heterogeneous environment is better described by the atomistically detailed MD simulations, while the homogeneous environment of a dilute aqueous solution can be more efficiently handled with a continuum model. The large size of the atomistic model limits the range of computationally affordable electronic structure methods. The principally lower computational costs of a continuum model allow us to use high level quantum models and thereby to analyze individual bonds at a much more sophisticated level. As both methods result in a macroscopically correct picture of the process, their direct comparison yields complementing information for each other of the chemical process studied.
The dielectric constant (ϵ) of the interior of enzymes has been estimated to lie between 4 and 20, 34 and model studies with (CH 3 ) 3 S + mimicking SAM show that the methyl transfer is greatly accelerated in environments with low values for ϵ. 35 A change of the value of ϵ in the electrostatic continuum model from 78 (water) to 4.24 (diethyl ether) in the QM calculations provides an estimate of the magnitude of the electrostatic effect on the reaction rate caused by the desolvation of the reactants in the enzyme. The X-ray structure of AtHTMT1 shows the presence of additional H 2 O molecules at the reactive site, indicating chemical interactions of the nucleophile with these molecules. The interaction between a small solute ion and its environment is addressed in continuum models by the addition of explicit solvent molecules to the quantum calculation. 36 The number and position of these solvent molecules are an integral part of the continuum model and its parametrization. The addition of an explicit H 2 O molecule to the small nucleophile in the continuum solvation models can then be used to identify the effect of an individual hydrogen bond at the reactive site on the catalytic process.
The QM/MM MD and the ab initio QM calculations with an electrostatic continuum handle the effect of the chemical environment differently. Therefore, it is necessary to be prudent comparing their results. The main focus of this work lies on the relative reactivities of various nucleophiles X − and chemical trends presuming that both approaches describe the macroscopic reality correctly. Details of the computational and technical setup will be discussed in the Methods section. The Results section is split into subsections for the results from the computational blocks as laid-out above, and all observed effects will be condensed into one model in the final subsection, which leads to the conclusion about the catalytic mechanism and the promiscuity of AtHTMT1.

METHODS
SAM has an isoelectric point (pH I ) of 7.24 37 and the pK a values of the COOH and NH 3 + groups in the methionine branch are 1.80 and 7.80, respectively. Hence, three possible protonation states of the methione chain in the SAM molecule (anionic a SAM + , bipolar b SAM + , neutral n SAM + ; Scheme 1) were considered for the continuum QM calculations assuming temperate pH values for the cytosol. The acronyms SAM and SAH are extended by a superscript on the left ('a', 'b' or 'n') to indicate the protonation state of the methionine chain. A superscript on the right indicates the positive charge on the sulfur atom in SAM.
The electric charges on the methionine branch of SAM are involved with the methyl transfer reaction. To study this effect, the QM calculations were carried out with all three protonation states. The QM/MM MD simulations, on the other hand, were done exclusively with b SAM + , as a preliminary analysis indicates its dominance in an enzymic environment, vide infra.
2.1. Standard QM Calculations without the Enzyme. The comparatively small number of atoms involved in the reaction between SAM + and a nucleophile X − without the enzyme and an electrostatic continuum model for the chemical environment facilitates the application of midlevel quantum models to the problem. The properties of the sulfonium cation in SAM + (marked red in eq 1) dominate the S N 2 reaction of the methyl transfer. Extensive method testing on the analogue methyl transfer in trimethylsulfonium chloride [(CH 3 ) 3 SCl] shows that the M05/6-311+G(2d,p) computational level can describe the internal S N 2 reaction adequately in various solvents. 35,38 The initial geometry of the SAM + molecule for the geometry optimization is the equilibrated structure of SAM + in aqueous solution taken from our previous QM/MM MD study of the hydrolysis of SAM + . 22 All quantum calculations using continuum models were done with the Gaussian 09 software packages. 39 All results from the geometry optimizations were verified with additional frequency calculations. Minima showed no imaginary frequencies, and the results for transition states showed only one aligned with the S N 2 reaction path. Changes in Gibbs free energy (ΔG°) were calculated for 298.15 K, 1 atm and 1 mol/L. The correction to changes in Gibbs free energy caused by the transfer of the quantum results from data of the gas phase [(1 mol)/(24.5 L)] to solution [1 mol/L] is denoted as δG°, where Δξ is the difference in the sums of stoichiometric coefficients of the reaction. 35,40 The X-ray structure of AtHTMT1 indicates that H 2 O molecules in the reactive cavity near the nucleophile are likely to influence the enzymatic reaction. 14 Hence, it is necessary to add one explicit H 2 O molecule into the continuum calculations. 36 The calculations for Cl − and Br − in water and diethyl ether with and without explicit H 2 O molecules were implemented with the PCM and SMD models. 36 43−46 The H 2 O molecules in the reactive cavity of the enzyme as revealed by X-ray crystallographic analysis 14 can interact with the nucleophile; hence, the correct description of the interactions among the nucleophile and surrounding H 2 O molecules is critical for the computational model. The data in Table 1 indicate that only the SMD model describes this situation adequately. Therefore, we use the SMD model for the continuum calculations.
The calculations for Cl − and Br − also highlight the effect of the SAM + protonation state on the methyl transfer. The absolute barrier heights (Tables S.1 to S.3) show that n SAM + and b SAM + have similar values of Δ TS G env°w ith n SAM + always being the smaller one. The Coulomb repulsion between the halide ion and the overall negative charge on the methionine branch of a SAM + can be seen in the markedly higher barriers for a SAM + . The SMD calculations for the noncatalyzed methyl transfer reaction in the cytosol with all nucleophiles were done with n SAM + to compare the results with the fastest possible reaction in aqueous solution.
2.2. QM/MM MD Simulations with the Enzyme. The structure of AtHTMT1 with SAH at the reactive site was obtained from the Protein Data Bank. 14,47 The protonation state of the enzyme was determined by GROMACS and that of the SAH molecule with CHIMERA. 48,49 The b SAH + molecule was turned into b SAM + by replacing a H 2 O molecule in the cavity with the new methyl group. The influence of the environment on the methyl transfer reaction can be deduced from the comparison of the QM/MM and the QM calculations, as n SAM + and b SAM + yield similar barriers for the methyl transfer to Cl − and Br − .
The AMBER 14SB force field, the general amber force field (GAFF), and the TIP3P force field were used to describe AtHTMT1, b SAM + , and the 94 560 water molecules, respectively. 50−52 Periodic boundary conditions were used on a cubic water-box of (10 nm) 3 . The Coulomb and van der Waals cutoff radii were 0.9 nm, and long-range electrostatic interactions were calculated with the particle mesh Ewald (PME) method. 53 To preserve electrostatic neutrality and to keep the salt concentration close to 0.150 mol/L, Na + and Cl − ions were added to the system. Figure 1 shows a snapshot of the methyl transfer from b SAM + to Cl − in AtHTMT1, where the MM environment and the quantum H 2 O molecules in the solvation shell of b SAM +  and Cl − are omitted for clarity sake. Quantum atoms were described using the semiempirical methods PM7, 54 PM6-D3, 55,56 AM1, 57 and PM6-D3H4 58 as defined in the Gaussian 09/16 39,59 and MOPAC20116 60 program suites. The collective variable λ (PLUMED), 61 defined as the difference between the C−S (d 1 ) and C−Cl (d 2 ) bond lengths was followed to calculate free energy changes using the umbrella sampling/weighted-histogram analysis method (US/ WHAM). 62,63 The MD simulations started with a steepest descent geometry optimization of the whole simulation cell followed by a set of MM MD simulations beginning with a 2 ns NVT simulation [T VR = 300 K (new velocity rescaling), 64 τ = 2 fs (time steps)] with a restraint on the movement of the heavy backbone atoms and d 2 (40 000 kJ/mol/nm 2 ). The restraints on the backbone atoms were then lifted stepwise (1000, 100, 10, 0 kJ/mol/nm 2 ) in a series of 2, 2, 2, 4, and 60 ns NpT simulations [T VR = 300 K, p Ber = 1 bar (Berendsen barostat), τ = 2 fs], which allows the density of the system to relax to its natural state. The last step of the equilibration phase was a 180 ns NpT simulation [T VR = 300 K, p Ber = 1 bar, τ = 2 fs] without any constraints, which yielded a constant rmd (root mean displacement) value for the heavy atoms in the backbone of the enzyme, while the running average of the total energy and its standard deviation have become constant. The equilibrated structure was then used to construct the QM/MM model for the methyl transfer reaction and to determine the number of QM H 2 O molecules for the solvation model by promoting all H 2 O molecules within 500 pm of the nucleophile to the quantum level.
The collective variable λ and the methyl transfer were driven in steps of 0.015 nm (40 ps each) from −0.18 to 0.225 nm or −0.195 to 0.225 nm depending on the nucleophile with a length constraint of 40 000 kJ/mol/nm 2 . The angle of attack of the nucleophile (4000 kJ/mol/rad 2 ) and position of the QM H 2 O molecules in their solvation shells (2000 kJ/mol) was constrained with functions provided by PLUMED to avoid their diffusion away from the nucleophile. The transition region was resampled with smaller λ steps and/or larger constraints, if individual windows in this region sampled poorly. The overall ΔG°curve was then constructed from the individual segments.

RESULTS AND DISCUSSION
3.1. Standard QM Calculations without the Enzyme. Table 2 lists the results of Δ TS G env°( eq 4) from the SMD calculations with the neutral n SAM + molecule. The reactive complex was used as the energetic reference due to its conceptual similarity to the enzyme−substrate complex in the Michaelis−Menten model of enzymatic reactivity. It should be noted that the reactions for Cl − and Br − without explicit H 2 O molecules in diethyl ether and that for Br − with an explicit H 2 O molecule appear to be barrier-less if the solvent separated reactants are taken as the energetic reference, because the formation of the contact ion pair in the reactive complex in ether is associated with a huge gain in Gibbs free energy (Table  S.1), if the solvent separated reactants are chosen as reference.
The results for Δ TS G wat°f rom the SMD calculations with water as solvent and no explicit H 2 O molecule show that the X − with the more nucleophilic atom attacking the methyl group has lower reaction barriers and reacts faster with n SAM + (eq 7) where the ambidentate ions SCN − and OCN − engage in faster reactions with the softer end of them attacking the methyl group favoring S over N and N over O. 45,46 Kinetic experiments with various trimethylsulfonium salts [(CH 3 ) 3 S + · X − ] in methanol (ϵ = 32.6) and ethanol (ϵ = 24.9) show that (CH 3 ) 3 S + ·Br − decomposes at a similar rate to (CH 3 ) 3 S + ·NCS − . Hence, the position of the NCS − ion in eq 7 suggests that more types of interaction between the n SAM + and X − than just the nucleophilic attack have to be considered regarding the reactivity of SAM + . The desolvation of the reactants has been found to be important for the catalysis of the S N 2 methyl transfer. 33,35 The difference Δ DS G°between the Δ TS G env°v alues (eq 4) for water and ether Table 2. Δ TS G°(eq 4) and Δ +1 G°(eq 10) for the Methyl Transfer from the Neutral n SAM + Molecule to Various Nucleophiles Obtained from SMD Calculations with and without an Explicit  (8) is used to quantify the desolvation effect (Table 2) with ether marking the lower end of the effective permittivity range for the interior of an enzyme. 34 The data for Δ DS G°from SMD calculations with no explicit H 2 O molecule show that the methyl transfer from the n SAM + donor to the nucleophile is generally enhanced with the exception of the NCS − ion. Sulfur atoms typically do not form strong hydrogen-bond networks; hence, a small positive value for Δ DS G°is not surprising. The interaction strength of each ion with the environment leads to the following rearrangement of the reactivity order by barrier height The relative softness of the nucleophile appears to be less dominant in an apolar environment, but the halides seem to be generally more reactive. The effect of H 2 O molecules on the methyl transfer can be quantified similarly. The difference between the values for Δ TS G°in the same continuum model solvent (sol) with (sol +1) and without (sol) an explicit H 2 O molecule is defined as The fastest reaction, which has the largest values for v max and k cat , is observed with AtHTMT1 and NCS − ions having the lowest barrier of 82.8 kJ mol −1 . The highest AtHTMT1 barrier of 90.0 kJ mol −1 is observed in the reaction with Cl − . By contrast, the quantum calculation with a continuum solvation model (Table 2) shows that NCS − is the slowest of the trio and Br − is the fastest. Although the quantum calculation with a continuum solvation model indicates that the desolvation of the reactant will increase the rate of the reaction (typically Δ DS G°< 0 and Δ +1 G°> 0), it fails to reproduce the relative order of reactivity of the nucleophiles observed in the experiments with enzymes. Experiments show that the enzymatic environment not only enhances the methyl transfer, but also favors NCS − over other nucleophiles. Therefore, it is necessary to include the enzyme explicitly in the calculations.
3.2. QM/MM MD Simulations with the Enzyme. The CHIMERA analysis of SAH in AtHTM1 indicates that b SAH + is the dominating protonation state, 48,49 while the lowest barriers for the methyl transfer in an aqueous environment are observed with n SAH. The initial analysis with Cl − and Br − as nucleophiles indicates that both protonation states lead to similar barrier height ( Table 1). The docking of SAM + to AtHTMT1 locks its protonation state to b SAM + while the fast, unhindered proton movement in the aqueous environment enables the methyl transfer reaction to pass through the energetically more favorable transition state with n SAM + . The difference in protonation state depending on the chemical

ACS Omega
http://pubs.acs.org/journal/acsodf Article environment ensures that the methyl transfer reaction always passes through its lowest transition state. The current study of methyl transfer in the solution and in the enzyme follows different mechanisms as summarized in eq 11,where E is the enzyme, S is the solution, RC is the reactive complex, TS is the transition state, PC is the product complex, E·SAM + represents SAM + docked to the enzyme, and E·SAH + represents SAH docked to the enzyme. The rate determining step (rds) in which the reaction passes from RC to PC is unimolecular in both cases, and is governed only by the concentration of RC. Hence, the rate determining steps can be compared directly.
The change in free energy in the enzymatic environment (G enz°) related to the methyl transfer as it proceeds along the reaction coordinate λ (eq 6) obtained from QM/MM MD simulations with the PM6-D3H4 Hamiltonian for Cl − , Br − , and NCS − is shown in Figure 2. The first minimum in the G enz°c urves marks the [X − ···(E· b SAM + )] complex with the b SAM + molecule being bound to the enzyme. The structure associated with this minimum in free energy can be compared directly with the reactive complex (RC) in the simple S N 2 reactions with an electrostatic continuum model and the second minimum with the product complex (PC). By focusing solely on the reaction series RC → TS → PC, it is possible to compare the DFT results directly with those of the MD simulations.
The changes in Gibbs free energy G enz°a re plotted relative to the energy of the first minimum G min°, which has been identified as the reactive complex ( Figure 2).
where G min°i s the energetic reference, and the maximum in Δ RC G enz°c urve becomes equivalent to Δ TS G env°( eq 4). Δ TS G enz°i s also the rate determining barrier for the methyl transfer reaction in the enzyme as shown in Figure 2, and the barrier of the second step in the Michaelis-Menten model of enzyme kinetics. Table 4 summarizes the methyl transfer barriers Δ TS G enz°f or all nucleophiles X − obtained with four semiempirical Hamiltonians. No barrier was observed for the methyl transfer with the AM1 Hamiltonian for I − and NCS − . However, experiments on the enzyme kinetics of AtHTMT1 (Table 3) indicate a barrier of 82.8 kJ mol −1 for NCS − . The barrier in the reaction of I − with the other enzymes varies between 81.1 and 88.1 kJ mol −1 . Hence, the results obtained with the AM1 Hamiltonian are dismissed from further discussions. The experiments with AtHTMT1 indicate that I − reacts more rapidly than Br − and Cl − . 14 This observation is not reproduced by the series of MD simulations using the PM7 Hamiltonian, and the PM7 calculations are dismissed from further discussions, too.
The AM1 Hamiltonian has been used successfully in other MD simulations 30−32 of the enzyme catalyzed methyl transfer, but semiempirical Hamiltonians cannot universally transfer due to their empirical nature. The data in Table 4 indicate the necessity to validate their application when changes in the chemical problem studied cannot be ignored.
The MD simulations with the PM6-D3 and PM6-D3H4 Hamiltonians yield the same order of barrier heights Δ TS G enz°o f the experimental trio (Cl − > Br − > NCS − ). The terminal reaction speed v max observed in experiment 14 and the barrier heights deduced from them are given in Table 3. The barriers obtained with the PM6-D3 Hamiltonian are generally closer to barrier heights calculated from the experimental data and hence are more credible. The disagreement between experiment and MD simulation increases with the reaction speed of the nucleophile. It is the largest (23.0 kJ mol −1 ) for NCS − and the smallest (0.5 kJ mol −1 ) for Cl − . The data in Table 3 show   (13) has been studied in detail. 5 The rate constant for the catalyzed step k cat was found to be 364.4 ± 12.7 min −1 . The barrier of the inverse reaction Δ TS inv G enz°i s estimated to be 69.0 kJ mol −1 by transition state theory. 65 Δ TS inv G enz°c an be calculated from MD results by taking the second minimum (products) as the energetic reference (Tables S6 to S9). The PM6-D3H4 and PM6-D3 Hamiltonians yield similar barrier heights of 72.9 and 75.3 kJ mol −1 . Both values are in reasonably good agreement with the experimental values close to the threshold of chemical accuracy of 1 kcal mol −1 . The PM6 Hamiltonians describe the reverse reaction from the uncharged products CH 3 X and SAH to the transition state well, but the calculated barriers for the forward reaction with NCS − and I − appear to be too low compared to the experiment (Table 3), which suggests that the ionic reactive complex RC with these nucleophiles is not well described by these semiempirical Hamiltonians. The simulations with the PM6 Hamiltonians describe the acceleration of the methyl group transfer to NCS − and I − in agreement with experiment, but are likely to overestimate the effect.
The MD simulations for Cl − with the PM6-D3H4 Hamiltonian were done six times. The six values for Δ TS G enz°a re found in the range of 71.6 ± 1.8 kJ mol −1 with an average value of 71.1 kJ mol −1 and a standard deviation of 1.2 kJ mol −1 . The barriers for the studied nucleophiles are separated on average by 7.5 kJ mol −1 so that the relative reactivity orders obtained from a single simulation can be trusted provided that the errors in the barriers for all nucleophiles are similar. The exceptions are Cl − and OCN − , the barriers of which are separated by merely 1.6 kJ mol −1 , and can be regarded as approximately equal. The reactivity order of the seven nucleophiles in the AtHTMT1 catalyzed reaction based on Δ TS G enz°b y QM/MM MD simulations with the PM6 Hamiltonians are as follows   (Table 4) is defined as the number of H 2 O molecules within 500 pm from X − in the enzymatic environment. The nucleophile in the enzymatic environment is always surrounded by H 2 O molecules. The enzymatic environment does not desolvate the nucleophile completely, but the X − ions keep parts of their solvation shell from the aqueous environment. Hence, the nucleophiles can be regarded as partially solvated. The most H 2 O molecules can be found around the halides (Cl − , Br − ) and the least around NCS − and SCN − . This fact is independent of the semiempirical Hamiltonians, and suggests that NCS − and SCN − ions fit the best into the reactive cavity. They are the ions with the smallest number of solvent molecules and hence the most activated nucleophiles, which partially explains why NCS − becomes the most reactive (highest value for v max , Table 3) of the experimental trio.
Both Cl − and Br − have K M values significantly larger than 1 (Table 3), which indicates that the equilibrium for the dissociation of the enzyme−substrate complex lies on the side of the unbound enzyme following the definition of K M by Michaelis and Menten 66 or equivalently that K M is dominated by the dissociation of the enzyme−substrate complex in the steady state interpretation by Briggs and Haldane. 67 The dissociation can be seen in MD simulations for Cl − without constraints on the position of the nucleophile where the Cl − ion quickly escapes from the reactive pouch of AtHTMT1.
Schmidberger et al. reported a site-directed mutagenesis (Tyr172/Phe), which changes Tyrosine 172 to phenylalanine, ACS Omega http://pubs.acs.org/journal/acsodf Article yielded an enzyme showing a reduced "efficiency with chloride ions..., but otherwise a similar efficiency for bromide and thiocyanate". 14 The MD simulations provide a possible explanation for this observation. Figure 3 panels a and b display snapshots from the trajectories for the methyl transfer with Cl − and NCS − in the reactive complex, which is the first minimum in Figure 2. The orientation of the tyrosine OH group changes over time. When it points into the general direction of the SAM methyl group, it can form a chain of hydrogen bonds to the attacking Cl − ion (Figure 3a). These transitory hydrogen-bond chains can stabilize the position of the wandering Cl − ion in the S N 2 transition state and then facilitate the reaction. The directing property of Tyr172 is visible in the decrease of v max from 2.43 nmol min −1 (mg protein) −1 for the native enzyme to 0.91 for the mutant. 14 In the case of the NCS − ion (Figure 3b), this chain points to the S atom of the NCS − ion. Sulfur atoms are generally weak hydrogen-bond acceptors so that a possible hydrogen-bond chain is weak and has only a small influence on the transition state. A chain to the N atom of the NCS − ion perpendicular to the S···HO axis would be longer and less effective. Therefore, the changes in v max caused by the mutation are within the margin of error of the reported experiment. 14 In summary, Δ TS G enz°o btained from MD simulation correctly describes the relative reactivities of the nucleophiles observed in experiment by the value of v max , the regioselectivity of AtHTMT1, and the simulations provide clues for the understanding of the value of K M . The simulations correctly describe the essentials of methyl transfer in thale cress.
3.3. The Role of the Environment in the Catalytic Process. The combination of the results from the standard quantum calculations and the MD simulations provides insight into the mechanism of the enzymatic reaction with AtHTMT1. The number n H 2 O of H 2 O molecules in the vicinity of the nucleophile X − differs greatly in the two approaches. It is very small and fixed in the standard quantum calculation, while it is large and fluctuating in the MD simulations. Hence, the observed solvation shell of X − in the MD simulations is a simulation result in contrast to that in the standard quantum calculations, where it is a rigid, user-defined part of the calculation. The comparison of both can be used to determine the role of solute−solvent interactions in the enzymatic environment.
The MD barriers (Table 4) calculated with the PM6 Hamiltonians are smaller than their SMD counterparts ( Table  2)  The enzyme's preference for the pseudohalides attacking the methyl group with their softer atoms is observed in the change of relative nucleophilic activities. The SMD water calculations with no explicit H 2 O molecules (reactivity order in eq 7) for the cytosol indicate that NCS − would be the least reactive nucleophile of the experimental trio (Cl − , Br − , NCS − ), while the MD simulations for the same reaction catalyzed by AtHTMT1 (reactivity order in eq 14) suggest that NCS − is the most reactive. This change in reactivity order cannot be attributed to a simple desolvation effect because the same order can be observed with ether as the SMD solvent (eq 9). The addition of an explicit H 2 O molecule to the nucleophile does not change the reactivity order either; therefore, the change in the reactivity order has to be credited to the enzymatic environment. Figure 3c shows  Table 2 show that moving a H 2 O molecule from the S to the N atom in NCS − lowers the barrier of the methyl transfer by 10.7 kJ mol −1 . Therefore, AtHTMT1 not only catalyzes the methyl transfer by desolvating the nucleophile to support the underlying S N 2 reaction, but also facilitates the reaction by placing H 2 O molecules near the N atom of the NCS − ion. Such an advantageous arrangement of H 2 O molecules follows directly from the ability of the atoms in the nucleophile to accept hydrogen bonds. SMD calculations with ether as solvent show that the formation of the S···H hydrogen bond (244 pm) in NCS − ···H 2 O is endergonic (ΔG°= +4.9 kJ mol −1 ) and the correction for particle densities (eq 3) makes it slightly exergonic in solution (−3.0 kJ mol −1 ). The corresponding N··· H hydrogen bond in OH 2 ···NCS − (189 pm) is slightly exergonic in the gas phase (−1.4 kJ mol −1 ) and strongly in solution (−9.3 kJ mol −1 ). Hence, not much chemical work is required to desolvate the attacking S atom of NCS − in the enzymatic environment while keeping the N atom solvated to accelerate the reaction. The NCS − ion holds parts of his solvation shell in the catalytic process.
The ability of AtHTMT1 to selectively desolvate the attacking atom in the nucleophile is not very strong for Cl − . The longer Cl···H bond (221 pm) is stronger in solution (−11.5 kJ mol −1 ), and the Cl − ion keeps parts of its solvation shell on entering the enzyme ( Table 4). The H 2 O molecules in the vicinity of Cl − slow the methyl transfer as the quantum calculations indicate (Δ +1 G°, Table 2).
The data for Cl − and NCS − reveal the remarkable finetuning of the desolvation of the nucleophile in AtHTMT1. Nucleophile desolvation generally accelerates the S N 2 reaction ACS Omega http://pubs.acs.org/journal/acsodf Article of the methyl transfer. The partial preservation of the nucleophiles solvation shell boosts the reaction of NCS − , meanwhile it reduces the catalytic effect for the reaction with Cl − . It propels NCS − from the least to the most reactive nucleophile of the experimental trio and hence limits the loss of SAM + to the accidental production of CH 3 Cl.

CONCLUSIONS
A set of quantum SMD calculations and QM/MM MD simulations was done to understand the promiscuity of AtHTMT1 and the production of methyl chloride by thale cress and other plants. Published structural data were used as starting point for the computational analysis, and the companion kinetic data 14 for Cl − , Br − , and NCS − (experimental trio) to validate the results. The combination of the results from both computational methods shows that neither technique alone is sufficient to explain the mechanism of the observed promiscuity. The uncatalyzed reaction in the cytosol was studied with M05/6-311+G(2d,p) DFT calculations and water as SMD solvent, while the QM/MM (PM6-D3H4, AMBER 14SB, TIP3P) MD simulations were used to study the enzyme catalyzed reaction. The direct comparison of the absolute barrier heights obtained from the DFT calculations with those from the MD simulations indicate that NCS − has the largest gain in reaction speed for a nucleophile of the experimental trio, which supports the hypothesis of AtHTMT1 being part of the plantal defense against bacterial infections. 18 The preference for NCS − is apparent in the order of reaction speeds for the experimental trio. NCS − is the slowest reacting nucleophile in water (Table 2) in contrast to the simulations for the enzyme catalyzed one, which predicts NCS − to be the fastest of the trio (Table 4) in agreement with experiments (Table 3).
Reactant desolvation has been suggested to be the acting principle in enzyme catalysis. 33 Water-free calculations with ether as SMD solvent support this idea as barrier heights generally decrease in ether solution, but the relative barrier heights for the experimental trio and the order of reaction rates do not change. Only the MD simulations with AtHTMT1 reproduce the reactivity order correctly. Furthermore, the calculated absolute barrier heights for the rate determining step in the enzymatic reaction are found to be in much better agreement with those obtained from the experimental data ( Table 3). The average difference in barrier height between the MD simulations with the PM6-D3 Hamiltonian and the experimental value is 8.9 kJ mol −1 (PM6-D3H4, 25.1 kJ mol −1 ) with a standard deviation of 12.3 kJ mol −1 (PM6-D3H4, 6.1 kJ mol −1 ) for the experimental trio.
The MD simulations provide the atomistic details necessary to understand the methyl transfer in vivo. The NCS − ion fits best into the reactive site of enzyme since it is surrounded by the least H 2 O molecules. The NCS − and SCN − ions are the most desolvated and most activated ions. The regio-selectivity appears to be controlled by the abilities of the S and N atoms to be engaged in S N 2 reactions.
The H 2 O molecules hover around the N atom in NCS − in the enzymatic process, which reflects the ability of an individual atom in NCS − to engage in hydrogen bonding. The concept of partial desolvation not only focuses on the number of H 2 O molecules close to the nucleophile, but also comprises the hydrogen bonding partners of the H 2 O molecules. No H 2 O molecules can be found close to the S atom of NCS − , but the H 2 O molecules of the partial solvation shell of NCS − can be found close to the N atom. Such a configuration of H 2 O molecules improves the reactivity of the nucleophile as shown by the quantum calculations. The partial solvation shell of the nucleophile in the reactive pouch of the enzyme favors the reactivity of NCS − and consequently the formation of CH 3 SCN. Direct interactions between AtHTMT1 and reactants controlling the methyl transfer are not observed in the MD simulation.
The Cl − and Br − ions amass the most H 2 O molecules around them. The enzyme cannot effectively desolvate the Cl − ion, because the Cl···H hydrogen bond is much stronger than its S···H counterpart. The partial preservation of the solvation shell of the Cl − ion effectively deactivates the nucleophile. It seems that the fine-tuning of the solvation shell around the nucleophile by the enzyme disfavors the small halides.
The quantum chemical analysis 35 of the methyl transfer between a trimethyl sulfonium ion [(CH 3 ) 3 S + ] as stand-in for n SAM + and Cl − shows that the formation of intermediate [(CH 3 ) 3 S + ·Cl − ] ion pairs depends strongly on the permeability of the solvent. Therefore, the partial solvation of the nucleophile in the enzymatic environment also prevents Cl − from binding to the reactive site as observed in MD simulations without spatial constraints. The Cl − ion prefers to move out of the reactive cavity. The tendency of Cl − to leave the reactive pouch is also observed in the experimental K M value being much larger than 1, while it can still engage in the methyl transfer reaction. Hence, the formation of CH 3 Cl by AtHTMT1 may be regarded as an unfortunate byproduct of the general acceleration of S N 2 reactions in the enzymatic environment.
The model of a partial solvation shell for the nucleophile in AtHTMT1 can be verified in experiments with cyanate ions. The model concludes that the softer N atom attacks the methyl group while the H 2 O molecules of its partial solvation shell gather around the O atom. The simulations with the PM6-D3 Hamiltonian yield barrier heights in better agreement with those calculated from experimental data than the simulations with the PM6-D3H4 Hamiltonian. The calculations with the PM6-D3 Hamiltonian show that CH 3 NCO will be formed in such an experiment with a v max value larger than that for the formation of CH 3 Br (eq 14) while the PM6-D3H4 Hamiltonian indicates a lower value. In contrast, the DFT calculations with a H 2 O molecule attached to the nucleophile indicate that the reaction between n SAM + and the halides should be faster than that of OCN − in the uncatalyzed reaction with water as solvent. The proposed model of partial solvation to explain the basic selectivity of AtHTMT1 can assist the rational design of new, more complex biotechnological techniques.
Additional tables; evaluation of published Michaelis− Menten data; geometries, energies, frequencies from QM calculations; geometry data for the QM region in MD simulations (PDF) ACS Omega http://pubs.acs.org/journal/acsodf Article