Empirical Valence Bond Simulations Suggest a Direct Hydride Transfer Mechanism for Human Diamine Oxidase

Diamine oxidase (DAO) is an enzyme involved in the regulation of cell proliferation and the immune response. This enzyme performs oxidative deamination in the catabolism of biogenic amines, including, among others, histamine, putrescine, spermidine, and spermine. The mechanistic details underlying the reductive half-reaction of the DAO-catalyzed oxidative deamination which leads to the reduced enzyme cofactor and the aldehyde product are, however, still under debate. The catalytic mechanism was proposed to involve a prototropic shift from the substrate–Schiff base to the product–Schiff base, which includes the rate-limiting cleavage of the Cα–H bond by the conserved catalytic aspartate. Our detailed mechanistic study, performed using a combined quantum chemical cluster approach with empirical valence bond simulations, suggests that the rate-limiting cleavage of the Cα–H bond involves direct hydride transfer to the topaquinone cofactor—a mechanism that does not involve the formation of a Schiff base. Additional investigation of the D373E and D373N variants supported the hypothesis that the conserved catalytic aspartate is indeed essential for the reaction; however, it does not appear to serve as the catalytic base, as previously suggested. Rather, the electrostatic contributions of the most significant residues (including D373), together with the proximity of the Cu2+ cation to the reaction site, lower the activation barrier to drive the chemical reaction.


■ INTRODUCTION
Amine oxidases (AOs) are a group of enzymes that catalyze the degradation of amines into the corresponding aldehydes with the concomitant production of hydrogen peroxide and ammonia, as shown in Scheme 1. They are divided into two structurally distinct classes of enzymes, depending on the cofactor involved in the oxidative deamination reaction: 1−3 (i) the copper-containing AOs [e.g., diamine oxidase (DAO)], and (ii) the flavin adenine dinucleotide-containing AOs [e.g., monoamine oxidase B (MAO B)]. Both DAO and MAO B are involved in histamine degradation in the human body 4 and are thus crucial for the regulation of histamine at physiological concentrations.
In the central nervous system, histamine degradation is initiated by its methylation with histamine N-methyltransferase (HNMT), producing N-methylhistamine, which is inactive at the histamine receptor sites and which is further oxidized by MAO B to the corresponding imine. The imine leaves MAO B and undergoes non-enzymatic hydrolysis to the respective aldehyde, 5,6 which is rapidly further metabolized by aldehyde dehydrogenase to t-methyl-imidazoleacetic acid (Scheme 2a). 7 In contrast, in peripheral tissues, DAO acts directly on histamine (Scheme 2b), yielding imidazoleacetic acid as the end product of histamine degradation. 7 Because of its broad substrate specificity, including activity toward histamine-related compounds such as N-methylhistamine, 8 and others like cadaverine (1,5-diaminopentane), 9 putrescine (1,4-diaminobutane), 8 and the polyamine spermidine, 8 DAO is a frontline enzyme for the degradation of exogenous amines. Hence, reduced levels of both DAO activity and expression have been directly linked to histamine intolerance. 4 This, in turn, leads to numerous unpleasant symptoms, including, for example, diarrhea, headaches, nasal congestion, asthmatoid wheezing, hypotension, and arrhythmia, as described in ref 4 and references cited therein.
The DAO-catalyzed oxidative deamination consists of reductive and oxidative half-reactions. While the former leads to the reduced enzyme cofactor and the aldehyde product, 10 the latter involves the reoxidation of the enzyme cofactor by molecular oxygen, which produces ammonia and hydrogen peroxide (Scheme 1). Measurements of primary kinetic isotope effects indicated that the rate-limiting step of these reactions is the cleavage of the Cα−H bond. 11 Oxidation of pdimethylaminomethylbenzylamine and its deuterated analogue by DAO gave k H /k D values up to 5.5. 11 This led the authors to propose a mechanism that involves the formation of a Schiffbase between the amine group of the substrate and the C5carbonyl atom of the 2,4,5-trihydroxy phenylalanine-quinone (TPQ) cofactor (Scheme 3). The prototropic shift from the substrate−Schiff base to the product−Schiff base includes the rate-limiting cleavage of the Cα−H bond, 11−14 as postulated in a study where 2-hydroxy-5-tert-butyl-1,4-benzoquinone was used as a model for the topaquinone (TPQ) cofactor to study the mechanism of benzylamine oxidation in acetonitrile. 12 According to the proposed mechanism (Scheme 3), benzylamine and the model cofactor (TPQ ox ) form a covalent substrate−Schiff base complex (TPQ s−Sb ), followed by the base-catalyzed proton abstraction from the substrate−Schiff base complex, resulting in the product−Schiff base complex (TPQ p−Sb ). The product−Schiff base complex is finally hydrolyzed to an aldehyde and the aminoresorcinol form of the reduced cofactor (TPQ red ). In the enzyme-catalyzed reaction, a conserved aspartate residue was proposed to act as the catalytic base. 15 Extensive mutagenesis studies have indicated that D383 of Escherichia coli AO (ecAO) is critical for this enzyme's catalytic activity. 13,14 That is, the D383E variant showed very low, but measurable, activity. 13 The D383N and D383A variants were inactive, with the former having drastically reduced affinity for the 2-hydrazinopyridine inhibitor, compared to the wild-type (WT) enzyme. 14 It was, therefore, concluded that D383 performs an essential role in catalysis, acting as the catalytic base, allowing the conversion of a substrate−Schiff to a product−Schiff base complex.
The role of the copper ion has also been extensively studied, 17,18 indicating its necessity for the catalytic mechanism of these AOs. It was shown that the copper-free form of the pig-kidney AO completely lost its activity. 18 Upon the addition of copper, however, the enzyme activity was fully restored. 18 It was also demonstrated that Ni 2+ and Zn 2+ do not bind the protein, whereas Co 2+ can be incorporated to the same extent as Cu 2+ , but the Co-reconstituted enzyme shows only very modest activity. 18 The present study aims to elucidate the catalytic mechanism of the oxidative deamination of histamine by human DAO (hDAO) and to investigate the role of the catalytic D373 residue, which corresponds to D383 in the homologous enzyme ecAO. To address these issues, we performed quantum-mechanical (QM) cluster calculations of the enzyme active site, 19 followed by empirical valence bond (EVB) 20−22 calculations of the full system that allowed us to calculate the reaction free energies. Our calculations show that the ratelimiting step of histamine degradation is the cleavage of the Cα−H bond, as previously described. 11 However, the reaction proceeds via direct hydride transfer from histamine to the TPQ cofactor, thus ruling out the previously proposed Schiff-base formation. Moreover, studies of the D373E and D373N variants indeed prove that the conserved catalytic aspartate is essential for the reaction, but it does not serve as the catalytic base, as previously proposed. 23 By describing the catalytic mechanism of hDAO in detail, we provide a fundamental understanding of this enzyme that is necessary for the design of therapeutics for the treatment of histamine intolerance. 4,24 substrate−Schiff base complex. 12,25 Such models have been used extensively to study enzymatic reaction mechanisms, 19,26 as well as being the basis of the "theozymes" used for de novo enzyme design purposes. 27,28 Our initial model systems for the non-enzymatic reaction included (i) histamine and the TPQ cofactor as well as (ii) an extended model including histamine, TPQ, and a Cu 2+ ion. Both models were optimized at the M11L/6-31G(d) level of theory, with the addition of the SDD basis set for the extended model with copper. Using this approach, it was not possible to obtain an optimized structure for the substrate−Schiff base complex in either the gas phase or continuum solvent [using the conductor-like polarizable continuum model (CPCM), 29 with a dielectric constant of 80]. In all attempts to obtain this complex, potential energy surface scanning resulted in quinone-ring opening, which suggests that it is highly unlikely that such a substrate−Schiff base complex exists in the absence of the enzyme.
To examine whether this mechanism is possible in the hDAO active site, we created a larger QM cluster model of the enzyme. This model consisted of the histamine substrate, the Cu 2+ cation, the TPQ cofactor, and the active site residues (i.e., D373, Y371, W376, Y459, H510, H512, and H675), as shown in Figure 1. All amino acid side chains and the TPQ cluster were truncated at the Cα-atom and capped with methyl groups. The histamine was modeled in its neutral unionized form. Although it has been shown that the substrate binds the active site in its protonated form, 30 the substrate needs to be deprotonated to its unionized form to both restore its nucleophilic character and initiate the reaction shown in Scheme 3. 31,32 Because the physiological pK a value of histamine is around 9.4, the latter can easily be achieved by active site water molecules with a cost of only a few kcal·mol −1 . 33 In addition, the pK a of the 4-hydroxyl group of the TPQ cofactor has been suggested to be as low as ∼3, which is in part due to its proximity to the catalytic Cu 2+ ion which will stabilize the charge on the oxyanion and facilitate its deprotonation. 34 Therefore, deprotonating the TPQ cofactor in the active site will be facile, and in our simulations, it was modeled as an oxyanion, because this form was shown to exist in the active site even before the reaction occurs. 35 The resulting cluster model was embedded into continuum solvent, using the CPCM model with a dielectric constant of 4.
Once the optimized reactant complex was obtained, we performed an energy scan along the distance between the nitrogen atom of the chain amino group of histamine and the C5-carbonyl atom of the TPQ cofactor. As can be seen from Figure S1, it was not possible to locate a transition state along this reaction coordinate. We therefore examined other mechanistic possibilities. Our previous computational work on the related enzyme, MAO B, makes a one electron (radical) pathway highly unlikely; 36 we therefore focused on the direct hydride transfer from the substrate Cα−H moiety to the DAO cofactor, which is analogous to the mechanism proposed for the MAO-catalyzed oxidative deamination of amine substrates. 36−41 We note that it has been recently shown, both experimentally and computationally, that C−H hydrogen bond donors directly transfer the hydride anion to the carbonyl oxygen of 2,3dichloro-5,6-dicyano-1,4-benzoquinone, which is often used model quinone in mechanistic studies. 42,43 Moreover, copperpromoted hydride transfer to a carbonyl group was previously reported as well. 44 Therefore, we modeled the reaction where Cα−H hydride cleavage from the histamine substrate involves a direct transfer of the hydride anion (H − ) to the carbonyl O5atom of the TPQ cofactor ( Figure 2). QM calculations were performed using the same cluster model of the enzyme as described above. For the reaction coordinate, we selected the distance between the donor Cα−H atom of the histamine and the carbonyl O5 atom of the TPQ cofactor. All calculations were performed at the (CPCM)/M11L/6-311++G(2df,2pd)-SDD//(CPCM)/M11L/6-31G(d)-SDD level of theory, as described in the Methodology section.  Table S1, and the associated coordinates are also provided in the Supporting Information. Figure 2. Free energy profile for the direct Cα−H hydride transfer in hDAO-catalyzed deamination of histamine, obtained using a QM cluster model (blue). The feasible, but slightly less favorable, mechanistic alternative involving histamine deprotonation by the side chain of the active site D373 in the second step of the reductive half-reaction is shown in red. The experimental activation energy is displayed in pink. All energies are in kcal·mol −1 . A simplified version of the cluster model used in the present work, showing only the main reacting atoms, is also presented for clarity.

ACS Omega
Article Figure 2 shows the relevant QM-optimized stationary points for the DAO catalyzed oxidative degradation of histamine. The relevant interaction distances are shown in Table S1, and the corresponding coordinates are provided in the Supporting Information. In the reactant state, the histamine is oriented toward the TPQ cofactor in a favorable conformation for the reaction to occur, forming hydrogen bonds between its amino group and the side-chain carboxylate of D373, as well between the histamine-ring imino nitrogen and the side-chain oxygen atom of Y459. The distances between the reacting atoms of histamine and TPQ at the reactant state, i.e., Cα−O5 and Hα− O5, are 3.2 and 2.6 Å, respectively. At the transition state, the hydrogen being transferred is asymmetrically located between the Cα-carbon and the O5-acceptor atom on the TPQ cofactor, with Cα−H and O5−H distances of 1.4 and 1.1 Å, respectively. We note here that the calculated activation free energy for this process is 23.7 kcal·mol −1 (v imag = 1012i cm −1 ), which is significantly higher than the experimental value (17.1 kcal· mol −1 , 8 derived from the experimentally obtained k cat of 139 min −1 , using the Eyring−Polanyi equation 45 ), and therefore, further EVB simulations 20−22 were performed to probe this reaction, as described in the next section. We, however, note that the activation free energy of 23.7 kcal·mol −1 is comparable to the value of 24.4 kcal·mol −1 for the analogous reaction in MAO B that we previously obtained using the same QM-only cluster approach. 33 Natural bond orbital (NBO) 46 analysis of the QM optimized geometries demonstrates that this process is indeed associated with the transfer of a hydride anion, as can be seen from the calculated atomic charges (Table 1). In the reactant state, the sum of the partial charges on the atoms of histamine and the TPQ cofactor are 0.02 and −1.02, respectively. However, at the transition state, the partial charges change to 0.25 and −1.25, respectively, indicating a transfer of negative charge to the TPQ cofactor. Following the hydride abstraction, the system proceeds to the corresponding intermediates, characterized by the semireduced TPQ cofactor and the cationic substrate (Figures 1 and 2). The charge distribution further indicates that the formation of the new O5−H bond increases the total negative charge on the TPQ cofactor to −1.40, while simultaneously the total positive charge on the histamine increases to 0.80. The following step in the amine oxidation involves deprotonation of the substrate amino group by the O2 atom of the TPQ cofactor (for a definition of the atom naming, see Table 1). This step is facilitated by the increase in the negative charge on the O2 atom of the TPQ cofactor during the hydride transfer reaction. This process is calculated to have a very low activation energy of 1.7 kcal·mol −1 (v imag = 1056i cm −1 ) and is thus facile, making the initial hydride abstraction the rate-limiting step of the overall reaction. Finally, the overall reaction free energy (relative to the reactant state) is −2.5 kcal· mol −1 , resulting in a neutral trans-imine and fully reduced TPQ cofactor as the final products ( Figure 2).
To explore the putative role of the active site residue D373 in deprotonating the cationic substrate during the second step of the reaction, as an alternative to the O2 atom of the TPQ cofactor, we also examined this pathway ( Figure 2). Our calculations show that such a pathway is plausible, although the reaction proceeds through the transition state with a notably higher activation energy (4 kcal·mol −1 relative to the intermediate and 9.7 kcal·mol −1 relative to the reactant state) and with a slightly less stable product state (ΔG 0 = 4.0 kcal· mol −1 ). However, both pathways are substantially lower in energy than the initial hydride transfer step and are therefore likely to be viable options to complete the reaction.
EVB Simulations. While our calculations using the QM cluster model strongly suggest that DAO utilizes a hydride transfer mechanism similar to that in MAO B, the calculated activation energy of 23.7 kcal·mol −1 is substantially higher than the experimental value of 17.1 kcal·mol −1 . 8 This discrepancy could, in part, be due to either an incomplete electrostatic treatment in the continuum solvent model or due to insufficient conformational sampling. To address these issues, we complemented our QM cluster calculations with EVB 20−22 simulations using the two-state model described in Figure S2.
A fundamental part of the EVB philosophy is the need for a well-defined reference state, which can be either the nonenzymatic reaction in vacuum or in solution or the WT enzyme against a series of mutants (see detailed discussion in, e.g., refs 21, 22). The reference state is then parameterized to reproduce either experimental data, where available, or higher-level QM calculations, where experimental data is lacking. The obtained

ACS Omega
Article parameters are then transferred to describe the same chemical reaction in different environments. Such an approach is feasible because of the demonstrated phase-independence of the EVB off-diagonal elements, as validated in ref 47. In the present case, due to the lack of experimental data to describe non-enzymatic histamine degradation in aqueous solution, it is necessary to fit the non-enzymatic reference state to energies obtained from higher-level QM calculations. For simplicity, we have performed our reference calculations in vacuum to avoid artifacts introduced by using continuum rather than explicit solvation. We then performed EVB simulations in vacuum, as described in the Methodology section, fitting to the quantum chemical calculations, and then using the same parameters to study the reaction in solution and the DAO active site, as in previous work. 38,41 The vacuum energetics of the reaction were calculated at the M11L/6-31G(d) level of theory, with zeropoint energies and entropies calculated from the vibrational frequencies at the same level of theory as the geometry optimizations. The model for the direct hydride transfer reaction consisted of neutral histamine and TPQ cofactor with a deprotonated hydroxyl group. The obtained transition state was verified by frequency analysis (v imag = 745i cm −1 ), as well as by following the intrinsic reaction coordinate in the reactant and product directions, 48 which yielded an activation free energy, ΔG ⧧ , of 51.5 kcal·mol −1 , and a reaction free energy, ΔG 0 , of 50.8 kcal·mol −1 . These extremely large values are consistent with the expectation that the investigated hydride transfer would be extremely unfavorable in the gas phase.
The EVB free energies, obtained using the calibration parameters shown in Table S2, are given in Figure 3a. As our QM cluster calculations indicated that the hydride transfer is the rate-limiting step of the overall process, we have focused only on this step in our EVB calculations. From this data, it can be seen that, as expected, the inclusion of solvent effects already substantially reduces the hydride transfer barrier even for the non-enzymatic reaction. The barrier is further reduced in the hDAO active site, yielding a calculated activation free energy of 14.8 kcal·mol −1 . Previous studies have shown that protonated histamine is likely bound to the active site; however, deprotonated (neutral) histamine is essential for the reaction to take place. Thus, an intraenzymatic deprotonation of histamine must occur within the enzyme active site, the cost of which would be expected to be 2.7 kcal·mol −1 at physiological pH (assuming a pK a of 9.4 for histamine and using the relationship ΔG PT = 2.303RT(pK a (donor) − pH)). The cost of this deprotonation, which approximately equals the difference between the calculated free energy of activation of 14.8 kcal·mol −1 and the experimental value of 17.1 kcal·mol −1 (derived from k cat = 139 min −1 ), 8 needs to be added to the calculated value to account for the initial histamine deprotonation. This yields an overall calculated free energy of 17.5 kcal·mol −1 , in excellent agreement with the experimental value, thus giving credence to our proposed mechanism (Table  2).
Finally, we note that our EVB calculations predict an activation free energy of 26.6 kcal·mol −1 for the non-enzymatic reaction in aqueous solution, leading to a barrier reduction of 11.8 kcal·mol −1 upon moving to the enzyme active site. This would correspond to a tremendous rate acceleration of 10 9 -fold compared to the non-enzymatic reaction. This substantial rate enhancement is comparable in magnitude to the similarly large rate enhancement observed computationally in MAO Bcatalyzed dopamine oxidation, which has been demonstrated to operate through the same direct hydride transfer mechanism. 38 The hydride transfer is also predicted to be thermodynamically favorable, with a calculated reaction free energy, ΔG 0 , of −2.9 kcal·mol −1 .
Roles of D373 and the Metal Cofactor. Following from our modeling of the WT enzyme, we have also investigated the role of D373 in the catalytic mechanism by mutating this residue both to glutamine and the isostructural, but chemically distinct, asparagine (Figure 4). It was experimentally shown that the D373E mutation leads to low but still detectable activity, whereas the corresponding mutation to asparagine completely abolishes the catalytic activity (Table 2). 14 On the

ACS Omega
Article basis of this, it has been proposed that D373 acts as the catalytic base by abstracting proton from the substrate−Schiff base complex to produce the product−Schiff base complex (Scheme 3, :B − ). However, we demonstrate here ( Table 2) that we obtain good agreement with experiment when modeling this reaction as proceeding through a direct hydride transfer mechanism. Specifically, as can be seen from the obtained EVB energetics, the calculated activation free energy for the rate-limiting step of direct hydride transfer in the D373N hDAO variant is 29.1 kcal·mol −1 , which is a loss in activity of approximately 10 orders of magnitude compared to that calculated for the WT enzyme (with ΔG calc ⧧ = 14.8 kcal·mol −1 ). Also, there is an extremely high thermodynamic barrier to this process, as reflected in the reaction free energy of 19.9 kcal· mol −1 . Therefore, our calculations predict a complete loss of activity for this variant, based on the direct hydride transfer mechanism. Following from this, in the case of the D373E variant, we calculate an activation free energy of 19.4 kcal· mol −1 , in reasonable agreement with the experimental value of 21.5 kcal·mol −1 for ecAO. 14 We also note that in the case of ecAO, this mutation increases the activation free energy by 6.1 kcal·mol −1 . While we do not have the corresponding experimental value for the hDAO mutants, our calculated ΔΔG calc ⧧ of 4.6 kcal·mol −1 for mutant versus WT hDAO is in good agreement with the observations from E. coli (where ΔΔG exp ⧧ is 6.1 kcal·mol −1 ). Here, we take into account that, as shown in Table 2, WT hDAO is less active than the corresponding enzyme from E. coli.
From a structural perspective, our EVB simulations show that, in the calculated transition states (Figure 4a), D373 forms hydrogen bonds with the amino group of histamine (average distance 2.7 ± 0.1 Å and average angle 165.4 ± 7.9°) and with the side chain of Y463 (2.7 ± 0.1 Å). In the transition state of the D373E mutant, however, E373 forms hydrogen bonds primarily with W376 (2.8 ± 0.2 Å) and Y463 (2.6 ± 0.1 Å), thus bridging the two residues and not being engaged in hydrogen bonding with amino group of histamine (3.0 ± 0.2 Å and 142.7 ± 13.8°), which implies that this hydrogen bond is lost at the transition state. At the same time, the average distance between the cofactor TPQ:O2 and Cu 2+ cation is 2.10 ± 0.03 Å for WT DAO in contrast to 3.7 ± 0.2 Å for D373E. This distance further increases to 4.4 ± 0.1 Å for the D373N mutant ( Figure 4). The electrostatic contributions of the most significant residues ( Figure 5) show that the protein residue at position 373 (i.e., D in the wild type, E or N in mutants) affects catalysis the most (while some other residues also make contributions, as annotated in this figure, the contributions are similar between the different enzyme variants, and mostly cancel each out). Thus, the difference in the free energy of activation between the WT and D373E or D373N mutant enzyme may be explained by a combination of the loss of stabilizing interactions from the catalytic Cu 2+ cation, together with (in the case of the N373 variant) a slight loss of stabilizing interaction from D373, which in turn leads to the loss of activity.

■ CONCLUSIONS
In this work, we combined QM cluster calculations and EVB simulations to elucidate the catalytic mechanism of the DAO catalyzed oxidative deamination of histamine. Our results indicate that direct hydride transfer from the substrate Cα−H moiety to the O5 atom of TPQ cofactor is the rate-limiting step of the reductive half-reaction, which is in line with the experimentally confirmed rate-limiting step. EVB simulations of the hydride transfer reaction catalyzed by WT DAO predict an almost 9 orders of magnitude rate enhancement in comparison to the corresponding non-enzymatic reaction in aqueous solution. The same rate enhancement was observed for the equivalent reaction catalyzed by the MAO B enzyme on monoamines. 38 Moreover, we have also shown that while D373 is indeed essential for DAO catalysis, its role is not proton abstraction during the catalytic cycle, as previously proposed. 12−14 Namely, the average hydrogen bond distances at the transition state combined with the electrostatic contribu-

ACS Omega
Article tions of the most significant residues (in terms of energetic contributions) to the calculated activation free energies showed that proximity of the Cu 2+ cation, together with the residue D373, affect catalysis the most, with the loss of key interactions with these moieties leading to the observed loss of activities in the D373E and D373N variants. In conjunction with our previous work, 40 the results presented here suggest that the in vivo degradation of histamine to yield an aldehyde always occurs via the same catalytic mechanism, namely direct hydride transfer from histamine and N-methylhistamine to the corresponding cofactor (TPQ or FAD, respectively), regardless of the histamine catabolic pathway. A detailed understanding of the mechanism of histamine-degrading enzymes is important for the design of therapeutics for the treatment of histamine intolerance, bearing in mind that no medication to date is able to increase DAO activity and therefore keep the level of histamine degradation in balance. Insufficient histamine degradation by DAO however may also disrupt the other histamine-degrading pathway operating through HNMT and MAO B, resulting in the collapse of histamine clearance. 50,51 Therefore, the catalytic mechanism of the degradation of this possibly harmful compound is important for the successful development of drugs that could supplement impaired histamine-degradation pathways.
■ METHODOLOGY QM Calculations. The cluster model of the DAO active site used in this work consisted of the substrate histamine, the Cu 2+ cation, the TPQ cofactor, and the active site residues D373, Y371, W376, Y459, H510, H512, and H675, as shown in Figure  1. Geometries for all key stationary points were optimized using the M11L functional 52 at the M11L/6-31G(d) level of theory, with the thermal Gibbs free energy corrections extracted from frequency calculations at the same level of theory without scaling factors applied. A further single point calculation was performed with the larger 6-311++G(2df,2pd) basis set to refine the electronic energies. The results of these calculations are shown in Figure 2 and Table S4.
System Setup for the EVB Simulations. To calibrate our EVB simulations of the enzyme catalyzed reactions, we performed corresponding simulation of the non-enzymatic reaction in vacuum. The EVB gas-phase shift and coupling parameters (α 2 0 and H 12 , Table S2) were obtained by fitting these simulations to ΔG ⧧ and ΔG 0 obtained from QM calculations in vacuum at the M11L/6-31G(d) level of theory and then using the same parameter set for the reactions in aqueous solution and in the enzyme active site. All QM calculations are performed using the Gaussian 09 package. 53 The starting points for our EVB simulations of DAO were the coordinates of DAO in complex with the inhibitor aminoguanidine (PDB ID: 3MPH). 54 The D373E and D373N mutants of DAO were generated using the Dunbrack backbone-dependent rotamer library, 55 as implemented in Chimera. 56 The position of the inhibitor in this structure served as a reference point for initial manual positioning of the substrate (i.e., histamine) into the hDAO active site. Short MD simulations of 5 ns were then performed to obtain equilibrated and optimized structures, which served as a starting point for the subsequent protein parameterization. The protein model included one subunit of hDAO enclosed in a simulation sphere, with a 30 Å radius, centered at the reactive Cα of the substrate. All protein atoms outside this sphere were kept restrained to their starting positions by applying a 200 kcal·mol −1 ·Å −2 harmonic restraint. All ionized residues in simulations are listed in Table S5. Other ionizable residues outside this region were kept in their neutral state. All simulations were performed using the OPLS-AA force field, 57,58 and the copper ion in this enzyme was reparametrized to reproduce the correct tetrahedral coordination, based on the previously published copper dummy model by Liao et al. 59 This model has been successful in reproducing the Jahn−Teller effect seen in divalent copper and was deemed an excellent starting point for the reparametrization.
All molecular dynamics and EVB simulations in this work were performed using the Q (v5.10) simulation package. 60 The resulting structure of DAO in complex with histamine was equilibrated and gradually heated to 300 K over a short 25 ps simulation. The system was further equilibrated at 300 K, performing 10 ns of MD simulation starting at the transition state, for each reaction. An initial 200 kcal·mol −1 ·Å −2 harmonic restraint was applied on all heavy atoms in the simulation sphere to remove bad contacts. This was gradually reduced during the initial heating phase, while only a 0.5 kcal·mol −1 ·Å −2 restraint was left on the reacting atoms of the substrate histamine and TPQ cofactor. All solvent hydrogen atoms were constrained by the SHAKE algorithm, 61 and system temperatures were regulated using the Berendsen thermostat 62 with a 100 fs bath relaxation time. For the calculation of all nonbonded interactions, a cut-off value of 10 Å was used, with the exception of those involving reacting atoms for which a 99 Å cut-off (i.e., no cut-off) was used. The endpoint structure of the equilibration run was used as a starting point for the subsequent EVB simulations. The equilibrated structure was subjected to 10 parallel FEP runs by randomizing the atomic velocities to perform independent EVB calculations and ensure reliable sampling. The mapping procedure was performed in 51 mapping frames of 100 ps each, yielding a total simulation time of 5.1 ns per replica.
Substrate and Cofactor Parameterization. The program Q, 60 used for the EVB simulations, and GROMACS, 63 which was employed in the original Cu-model paper, use different force field descriptions and the former does not apply particle mesh Ewald summation for the long-range electrostatic interactions. Therefore, the published Cu parameters were further validated and optimized to allow a better fit with the previous experimental data. 64 For the new parameters, both the repulsive and attractive parts of the center particle of the ion model were optimized in respect to the resulting hydration free energies and ion−water radial distribution functions. It was also necessary to reduce the repulsive force on the non-interacting dummy particles from the 0.03 kcal·mol −1 ·Å −6 , used in the previous dummy model paper, 65 to 0.01 kcal·mol −1 ·Å −6 , or the reproduction of the Jahn−Teller distortions was not possible. This setup is in line with the results obtained in the original article with GROMACS, where no force at all was applied to the dummy particles. This approach was not chosen to avoid singularities when calculating the forces in Q.
The parameterization was performed using the TIP3P water model as implemented in Q, with the dummy model being solvated in an 18 Å water sphere subjected to the spherical boundary conditions, according to the surface-constrained allatom solvent model. 66 Hydrogen atoms were subjected to the SHAKE algorithm, with long-range interactions beyond a cutoff of 10 Å being treated using the LRF approach. 67 The metal dummy was kept in the simulation sphere by applying a weak position restraint of 0.5 kcal·mol −1 ·Å −2 . After relaxing the

ACS Omega
Article system by gradually increasing temperature over the course of 60 ps, the final coordinates were used to start five independent calculations of the solvation free energy, using 51 umbrellasampling windows of 100 ps simulation time each (5.1 ns simulation time per trajectory). The analysis of the free energies was performed using the Qfep tool of the Q program suite. The radial distribution function was calculated, from three 1 ns-long independent MD simulations, using VMD (1.9.1), 68 employing a bin size of 0.01 and a cut-off of 10 Å ( Figure S3).
The parameters for nonstandard protein residues (substrate and cofactor) were obtained using the ffld_server as implemented in Schrodinger's Macromodel suite 69 allowing for the construction of the relevant EVB force fields to correctly describe the reactant and product states of the hDAO-catalyzed reaction (Tables S6−S10). For both states, atomic charges were determined by fitting to the electrostatic potentials computed by the QM calculations at the HF/6-31G(d) level of theory according to the RESP scheme, as implemented in AmberTools15. 70 ■ ASSOCIATED CONTENT

* S Supporting Information
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acsomega.8b00346.
Additional computational data and all parameters necessary to reproduce our EVB simulations (PDF)