MM Study of the Reaction Mechanism of the Dehydratase Domain

Research Article Cite This: ACS Catal. 2018, 8, 10267−10278

pubs.acs.org/acscatalysis

QM/MM Study of the Reaction Mechanism of the Dehydratase Domain from Mammalian Fatty Acid Synthase Fabiola E. Medina, Rui P. P. Neves, Maria J. Ramos, and Pedro A. Fernandes* UCIBIO, REQUIMTE, Departamento de Química e Bioquímica, Faculdade de Ciências, Universidade do Porto, Rua do Campo Alegre s/n, 4169-007 Porto, Portugal

ACS Catal. Downloaded from pubs.acs.org by UNIV OF SUNDERLAND on 10/08/18. For personal use only.

S Supporting Information *

ABSTRACT: Dehydratase (DH) is a catalytic domain of the mammalian fatty acid synthase (mFAS), a multidomain enzyme with seven different active sites that work in tandem to carry out the biosynthesis of palmitic acid for de novo lipogenesis. DH catalyzes the dehydration of the β-hydroxyacyl to an α,β-unsaturated acyl intermediate. We have conducted hybrid QM/MM calculations to clarify the catalytic mechanism for the DH domain at the ONIOM(DFT/Amber) level of theory. The results have shown that the dehydration step occurs in two stages: (i) the His878-imidazole acts as a base deprotonating the Cα of the β-hydroxyacyl (HAC) substrate and (ii) the β-elimination of the β-hydroxyl of HAC proceeds with late protonation of the leaving hydroxide by the Asp1033-carboxylic group, forming a water molecule as a byproduct. The α-deprotonation depends on an oxyanion hole mechanism where the HAC’s α-carbonyl is anchored by two strong hydrogen bonds from the neighboring Gly888 and the intramolecular β-hydroxyl, positioning the Cα of HAC for deprotonation by His878. A positively charged His1037 improves the acidic character of Asp1033 and completes the catalytic triad in DH, because when His1037 is neutral the positively charged His878 behaves as the acid in the β-elimination step. We observe that the positively charged His1037 renders the β-elimination step more thermodynamically favorable (ΔrG of −15.9 kcal·mol−1). The β-elimination step exhibits a Gibbs energy barrier of 14.1 kcal·mol−1 and it is the rate-limiting step of the reaction (in agreement with the experimental barrier of ∼17 kcal·mol−1. Nevertheless, the rate-limiting step does not seem to be dependent on the protonation of His1037. Through evaluation of the electrostatic effect per residue on the rate-limiting step, we concluded also that the electrostatic contribution of the enzyme’s body does not seem significant, even though there are many positively and negatively charged residues close to the leaving β-hydroxyl group of HAC. KEYWORDS: fatty acid synthase (FAS), dehydratase (DH), reaction mechanism, QM/MM ONIOM, DFT, α,β-unsaturated

1. INTRODUCTION Mammalian fatty acid synthase (mFAS) is a cytosolic enzyme that catalyzes the synthesis of the palmitic acid.1−3 The overexpression of mFAS is related to diseases such as cancer, obesity, and infectious diseases,4 which makes it a valuable target for drug design campaigns.5,6 Mammalian FAS (type I FAS) is a homodimeric enzyme with nine domains per monomer, six of them catalytic. In bacteria and plants these domains are freestanding monofunctional enzymes and are classified as type II FAS.1 The chemical reactions catalyzed by FAS are conserved in almost all organisms, regardless of the different architectural arrangement of the domains that catalyze each of them. The overall reaction mechanism of the mFAS is divided into three generic phases: initiation, elongation, and termination. The elongation phase involves four consecutive reaction cycles, repeated cyclically seven times until the palmitic acid product is released. Therefore, the four active monomers participating in the elongation phase carry out a total of 28 catalytic cycles. Each elongation cycle (4 reactions) results in an increment of © XXXX American Chemical Society

an ethyl group (see Figure 1). During each turnover, the acyl carrier protein (ACP) domain, which covalently binds to the substrate through a phosphopantetheine prosthetic group (PNS),7 transports the intermediate produced in each catalytic domain to the active site of the next domain to carry out the following reaction. The initiation phase begins in the malonylacetyl transferase (MAT) domain8−10 where the primer acetyl unit from acetyl−CoA is loaded into MAT and transferred to ACP. The latter loads it into the active site of the β-ketoacyl synthase (KS).11 After that the MAT transfers the extender malonyl unit to KS from the malonyl-CoA, to produce the βketoacyl intermediate through a condensation reaction. The latter is converted to a β-hydroxyacyl intermediate (HAC) by the NADPH-dependent β-ketoreductase (KR) domain.12,13 The next step is the formation of an α,β-unsaturated acyl from the dehydration of the HAC intermediate by the dehydratase Received: July 5, 2018 Revised: September 20, 2018

10267

DOI: 10.1021/acscatal.8b02616 ACS Catal. 2018, 8, 10267−10278

Research Article

ACS Catalysis

Figure 1. Schematic representation of the complete catalytic cycle of mFAS. Turnover is divided into three main stages: initiation, initiation/ loading; elongation, condensation and β-carbon processing; termination, chain termination.

domain (DH).14 This intermediate will undergo a second reduction led by the NADPH-dependent enoyl reductase (ER),15 yielding a fully saturated carbon−carbon bond. The elongation will continue with a new cycle, with a malonyl unit loaded into MAT. The cycle continues until a 16-long carbon chain is formed. Finally, at the termination phase, the final product (palmitic acid) is released from ACP by the thioesterase (TE) domain.1,2,16,17 In this study, we address the mammalian DH (mDH) domain of mFAS, which is a pseudodimeric monomer that contains five α-helix and six antiparallel β-sheets. The β-sheets wrap around the α-helices,14,18,19 as shown in Figure S1, in the Supporting Information (SI). This arrangement is similar to the dimeric type-II FAS DH.2 Despite the fact that the sequences are not similar (13.6% sequence similarity), the fold of the type-II form of DH is superimposable with the mDH domain and the catalytic residues are also conserved (Figure 2). The active site of the mDH is located in the putative interface of the pseudodimer, the same as for the Escherichia coli dimeric DH FAS (ecDH).14,20 However, mDH contains a single active site formed by residues His878 from the Nterminal β-sheet fold and Asp1033 from the C-terminal fold, unlike that of the E. coli DH FAS, which has two independent active sites located between the heterodimer but with the same residues involved. While the side chain groups of an Asp35, Thr39, Glu47, and Ala59 were suggested to be involved in substrate recognition for the type-II FAS (PKS),7 the crystallographic structure of the mFAS (PDB code 2VZ9)2 suggests that substrate access should occur through a hydrophobic pocket. However, there is no X-ray structure where the substrate can be found. On the other hand, the crystal structure of the E. coli (PDB code 4KEH)21 has been reported with a mimetic substrate covalently bonded to the acyl carrier protein (ACP). A general mechanism has been put forward for the homologous bacterial FabA enzyme.14,22 It was proposed that the dehydration is initiated by the deprotonation of the αcarbon in the substrate, concomitantly with the protonation of the β-hydroxyl oxygen atom of the HAC substrate (Scheme 1). The catalytic dyad comprehends a His (His70B for FabA) that acts as a base during the α-deprotonation step and an Asp (Asp84A for FabA) that acts as an acid in the protonation of the β-hydroxyl oxygen.23 Additionally, X-ray studies suggest

Figure 2. mDH obtained from the PDB structure 2VZ9. Catalytic residues are shown in ball-and-stick representation (colored by element, green C atoms). Equivalent residues in the E. coli DH (PDB code 4KEH), the substrate, and the prosthetic group (PNS) are colored by element, with purple C atoms.

that a conserved water molecule in the active site of FabA may participate in catalysis or, alternatively, represent the product of the dehydration reaction.14,18,22 Mutagenesis studies also show that a Gln88A is crucial to anchor the Asp84B in an appropriate position in the active site (see Scheme 1), and the catalytic activity of ecDH is reduced if the Gln is replaced by an Ala or Thr.14 On the contrary, the activity remains unchanged when the Gln is replaced by a His.14 This result is consistent with the structural data from the type-I FAS (mammalian FAS), where in the active site of mDH there is a His (His1037) in an equivalent position to that of glutamine. Hence, this tetrad (His70B, Asp84A, Gln88A, and a water molecule) is suggested to be essential in the catalytic core of FabA.14,19,24 These residues are conserved across different 10268


Research Article

ACS Catalysis

The 1R3 ligand resembled the PNS−fatty acid substrate covalent complex (it exhibits an aliphatic hydrocarbon chain, followed by the PNS that anchors the ACP domain) and was positioned between the proposed catalytic residues (His878 and Asp1033). The HAC substrate was built from the aliphatic chain of 1R3 with the GaussView software26 (see Figure S2 in the SI). We assumed that as the elongation phase of the FAS mechanism undergoes 7 consecutive cycles, the study of the dehydration of the smallest substrate of DH (C4 substrate) should resemble and be representative of that of larger substrates. Hence, following our previous work on the KR domain (Medina et al., 2017),13 we assumed that the binding of the C10 substrate should also be representative of the binding of the C4 substrate, which would simplify the extent and bring more reliability to the modeling of the enzyme:substrate complex. We have also modeled PNS because its presence mimics the steric constraints that are imposed on HAC during the reaction taking place in the mDH domain. The parameters for HAC and PNS were determined from optimized structures at the HF/6-31G(d) level of theory in vacuum. The Antechamber tool, incorporated in the AmberTools of AMBER12,27 was used with the GAFF force field28 to derive the intramolecular and Lennard−Jones parameters, while the RESP (restrained electrostatic potential) method29 was used to derive atomic charges at the HF/631G(d) level of theory30 to be consistent with the AMBER FF99SBildn force field used to parametrize the remaining atoms of the MM portion of the enzyme.31,32 The hydrogen atoms were added using the Xleap module of the AMBER12 software.27 The file parameters are included in the SI. Aside from the proposed catalytic residues, all protonation states were assigned by Propka3.0 (see Table S1 in the SI).33,34 For the proposed catalytic residues, protonation states were assigned based on experimental data (even though some differed from the Propka3.0 predictions):22 Asp1033 (predicted pKa of 5.7) was modeled in the neutral form, His878 (predicted pKa of 6.5) was protonated in the δ-nitrogen,23 and His1037 (predicted pKa of 2.6) was tested in all three possible

Scheme 1. Chemical Reaction Proposed for the DH Domain of the Mammalian FASa

a

The residues in curly brackets are the catalytic residues proposed in the literature for FabA and mFAS (in bold).

species.14,20,25 Accordingly, the His878, Asp1033, and His1037 triad in mDH are the proposed catalytic residues.14,24 In order to gather evidence that supports/denies the proposed reaction mechanism and to understand how specific enzyme−substrate interactions contribute to the activation energy of the rate-limiting step, we performed Molecular Dynamics (MD) and hybrid Quantum Mechanics/Molecular Mechanism (QM/MM) calculations on the DH domain of mFAS to determine the thermodynamic/kinetic profile for the proposed chemical steps. These results provide atomistic-level understanding of the catalytic mechanism of the DH domain, a fundamental piece of the mFAS biosynthetic machinery.

2. METHODS We built the QM/MM model from the DH domain of the subunit A of mFAS (PDB ID 2VZ9). As the substrate was not present in the mDH, it was modeled from a sulfonyl 3decynoyl pantetheinamide probe (1R3) complexed with E. coli DH (PDB code 4KEH), which exhibits a similar fold to mDH (rms of 3.1 Å over 257 pairs of heavy atoms).

Figure 3. Overall QM/MM model. (Left) Cartoon representation of the mDH domain with the PNS prosthetic group in stick representation. Water molecules considered at the MM layer are represented in white with black outline (stick representation). (Right) DFT layer of 94 atoms is represented in a shadowed ball-and-stick representation. 10269


Research Article

ACS Catalysis

urations of His1037. The generalized Born model40 was used to describe the solvent, and the Monte Carlo method was utilized to simulate changes in the protonation state of His1037 for the several CpHMD simulations (details in section 2.3 in the SI). All of the QM/MM calculations were performed using the electrostatic embedding scheme of the ONIOM methodology,41,42 available in Gaussian 09.43 For the DFT layer, the truncated bonds in the DFT frontier were treated with the linkatom method.27,44 All geometry optimizations and linear transit scans were performed at the B3LYP/6-31G(d):FF99SB level of theory.45−47 Previous works have shown that B3LYP can accurately reproduce geometries, even though it is also well known that it provides underestimated activation energies in ca. 2−3 kcal·mol−1.48−51 The structure of the stationary points was accurately determined by fully optimizing each transition state and following its intrinsic reaction coordinate (IRC).52 Hirshfeld atomic charges53−55 and natural bond order (NBO)/ Wiberg bond index,56−58 analysis were calculated at the B3LYP/6-311+G(2d,2p):FF99SB level of theory to discuss the electronic density changes and bond orders in the DFT layer throughout the chemical reaction. Frequency calculations were performed in order to confirm the absence of imaginary frequencies in the stationary points of the reactant, intermediate, and product. In the case of transition states, we confirmed the presence of a single imaginary vibrational frequency which corresponded to a maximum only along the reaction coordinate of the PES. We only accounted for vibrational temperatures larger than 100 K for the calculation of thermal and entropy contributions (as validated elsewhere).59 The interpolation between rotational and harmonic vibrational approximations was evaluated by a weighted function combining the free rotor and harmonic oscillator approximations.60 The QM/MM energies were refined subsequently by singlepoint calculations with five hybrid DFs (B3LYP,46,47,61 BMK,62 M06-2X,63 MPWB1K,64 and ωB97X-D65) using the larger 6311+G(2d,2p) basis set; the results for all of them are shown in Figure S6 in the SI. All of the density functionals showed a similar overall agreement regarding the shape of the PES for the reaction mechanism. Despite small differences, a benchmarking66 study in similar reactions indicates a preference by the ωB97X-D. Hence, the energy calculations for all stationary points were performed at the ωB97X-D/6311+G(2d,2p):FF99SB//B3LYP/6-31G(d):FF99SB level.

protonation states. In particular, for His1037, it is well known that the prediction of pKa of His residues from empirical functions is often erroneous.35,36 Hence, we considered that studying the reaction mechanism with all possible protonation states of His1037 would clarify the role of this residue for the reaction mechanism of mDH. TIP3P water molecules37 were added to solvate the system within 12.0 Å of the solute, in a rectangular water box. We added 11 Na+ counterions to neutralize the charge of the system. The geometry of the system was optimized in four minimization steps from the crystallographic X-ray structure in order to avoid eventual tensions and clashes: (1) the water molecules were freely minimized, and the rest of the system was restrained with harmonic potential force constants of 50 kcal·mol−1·Å−2; (2) the hydrogens were freely minimized, with the remaining system restrained with 50 kcal·mol−1·Å−2 harmonic force constants; (3) then the side chains of the protein were energy minimized, while the remaining solute was kept fixed with 50 kcal·mol−1·Å−2 harmonic force constants; (4) finally, the system was minimized freely, exhibiting a rmsd of 0.65 Å relative to the X-ray structure. All bonds involving hydrogen atoms were constrained with the SHAKE algorithm38 to use an integration step of 2 fs. A progressive heating was performed linearly from 0 to 310 K, for 20 ps, followed by 30 ps at 310 K in the NVT ensemble. A conventional molecular dynamics (cMD) production was then performed for 20 ns in order to analyze the conformations of the active site region (details in section 2.1 in the SI). The cMD showed that the substrate remained in the active site (rmsd ≈ 1.87 Å, see Figure S3); however, the active site adopted conformations that were not appropriate for the reaction to take place through the proposed mechanism (see Figure S4). The QM/MM model was defined from the minimized structure obtained at the FF99SBildn level of theory. It included the whole mDH domain with the HAC substrate and PNS and all water molecules within 6 Å of the HAC substrate and the PNS. The final model comprehended 4083 atoms and an overall charge of −11 au. The DFT layer comprised 94 atoms. We included the side chain of the three catalytic residues proposed in the literature (His878, His1037, and Asp1033). We also included parts of nearby interacting residues: Leu885, Phe886, Pro887, Gly888, and Tyr1003 as well as the entire HAC substrate and the ethylthiol moiety of PNS, as shown in Figure 3. The rest of the system (3989 atoms) was treated at the MM level, and waters at the MM level were kept fixed (see Figure 3, and the coordinates are in the SI). Furthermore, we built an additional QM/MM model adding the side chain of the Leu1036, Met1035, and Ala1034 to the DFT layer to account for the effect of the environment around His1037. This model comprehended 120 atoms in the DFT layer (see Figure S5 in the SI). However, the PES explored did not give relevant additional information (different from the DFT layer with 94 atoms) that could contribute to a better atomistic description of the reaction mechanism (details in section 2.2 in the SI). Finally, the His1037 protonation in its three possible protonation states (δ- and ε-protonated imidazole and imidazolium forms) was also modeled. The implicit solvent constant pH molecular dynamics (CpHMD)39 method, as implemented in AMBER12,27 was used to predict the pKa of His1037 in order to estimate the difference in energy between the neutral and positively charged config-

3. RESULTS AND DISCUSSION In the following sections we discuss the reaction mechanism obtained for the mDH domain with the QM/MM calculations. We also discuss the geometry changes and charge distribution in the stationary states and the thermodynamics of the whole process. The His878 and the Asp1033, proposed as catalytic residues, were modeled in the neutral form. Furthermore, we followed the reaction pathway considering the three protonation states of the His1037 (imidazolium form, positively charged; δ- and ε-nitrogen protonated). Our results indicated that the His1037 was not crucial for the α-deprotonation stage; however, the product obtained in the β-elimination stage with the His1037 positively charged was thermodynamically more favorable than the product obtained with the neutral His1037 (δ- and εnitrogen protonation). Scheme 2 describes the α-deprotonation and the β-hydroxyl elimination stages proposed in this 10270


Research Article

ACS Catalysis Scheme 2. General Mechanism Proposed for the Two Steps Involved in the mDH

work (details of neutral protonations, Scheme S1 in the SI). In addition, we also calculated the relative energy of the δ- and εtautomer of His1037 to the positively charged His1037 to compare the thermodynamics of each of the three mechanisms. We calculated a pKa of 8.34 for His1037, with the constant pH MD (CpHMD) methodology, which indicated that the positively charged form of His1037 was the most favored. In order to estimate the energy difference between the δ- and the ε-tautomer of His1037 (ΔΔEδ→ε), we reoptimized the reactant, interconverting δ-His1037 to ε-His1037 and vice versa, and compared the DFT energies for each transformation with single-point calculations by using the mechanical embedding scheme. 3.1. First Reaction Stage: α-Deprotonation of HAC Substrate. In the reactant (React) state, HAC was bonded through the α-carbonyl to the PNS-thioethyl, forming the HAC-PNS covalent complex. The α-carbonyl was confined in an oxyanion hole formed by an internal hydrogen bond with the β-hydroxyl (1.95 Å) and with the Gly888-amide (2.05 Å), as shown in Figure 4. In addition, the β-hydroxyl oxygen of HAC also formed a hydrogen bond with the Asp1033carboxylic group (1.72 Å). In the geometry-optimized React state, the catalytic His878nitrogen (δ-protonated) was adequately positioned for the αdeprotonation to occur: at 2.04 Å from the α-hydrogen and forming an angle of 153.3° with the α-carbon (see Figure 4). This active form of His878 resembles the distance of 2.5 Å from His70-Nε toward α-carbon in E. coli DH analogues obtained experimentally (FabA and FabZ).22,23 The first transition state (TS1) was observed when the acid hydrogen (Hα) was located at 1.48 Å from the α-carbon and at 1.25 Å from the ε-nitrogen. The imaginary frequency was 1170i cm−1, and it was characterized by an antisymmetric stretching of the Cα−Hα and Hα−Nε bonds. We also observed that the hydrogen bonds established by HAC with the oxyanion hole shortened in 0.12 and 0.13 Å, respectively (see Figure 4 and Figure S7 in the SI), which further supported the effect of the oxyanion hole at stabilizing the charge building up at the α-carbonyl. Additionally, Gly888 was a key residue in the formation of the β-hydroxyl intramolecular hydrogen bond and the resulting oxyanion hole arrangement.67 In the absence of the hydrogen

Figure 4. Optimized coordinates for the DFT layer. (Top) React conformation. (Middle) Transition state (TS1); oxyanion hole region is highlighted in green (shadow square). (Bottom) First intermediate (I1). Colored transparency behind the HAC corresponds to the overall charge variations from positive (blue) to negatively charged (red). Distances are indicated in Angstrsoms.

bond with Gly888, we observed that only the His878imidazole formed a hydrogen bond with the β-hydroxyl oxygen of HAC. Since no base other than His878 could perform the α-deprotonation, the activation energy (ΔE‡) rose above 30 kcal·mol−1, rendering the reaction unfavorable (see Figure S8 in the SI). The unusual kind of oxyanion hole observed, in which one of the hydrogen bonds was donated by the substrate itself, resembles a form of self-promotion. This hypothesis was also corroborated through the replacement of the HAC-hydroxyl and the Gly888-amide to a methyl group in React and TS1. The barrier was increased in 7.3 kcal·mol−1 for the hydroxyl replacement and 4 kcal·mol−1 for the amine group mutation. The double replacement increased the barrier by 11.1 kcal· mol−1 in a cooperative effect of the hydroxyl and amine groups toward the HAC-carbonyl. Thus, the oxyanion hole was a key requisite for the α-deprotonation step. In agreement with our results, studies in others enzymes (thiolase,68 cholinesterases,69 ketosteroid isomerase,70 and transferase10) have revealed that 10271


Research Article

ACS Catalysis the oxyanion hole arrangement stabilizes the transition state in ca. 3−6 kcal·mol−1. The TS1 exhibited a Gibbs activation energy (ΔG1‡) of 9.0 kcal·mol−1 and led to an intermediate product (I1) with a Gibbs reaction energy of 7.1 kcal·mol−1. On I1, the His878imidazolium and a negatively charged HAC (HACδ−) were formed (see Tables S2 and S3 in the SI). In the α-deprotonation we expected a carbanion formation in HACδ− as a product of the deprotonation with the atomic charge located mainly in the α-carbon; however, we observed that the relative atomic charge was distributed between the αcarbon and the adjacent carbonyl. The atomic charge increased in ca. 21% out of the total for each atom in Cα−CO (from React to I1 in HAC). In addition, the Cα−CO bond decreased from 1.50 Å in React to 1.39 Å in I1, resembling the CC bond length in aromatic rings, while the CO bond increased from 1.23 Å in React to 1.28 Å in I1, also similar to the aromatic CO bond length of 1.25 Å. Thus, we verified that no carbanion intermediate was formed, but the formed structure was best described as an enolate species. In fact, in line with previous studies in enzymes such as citrate synthase71 and enzymes of the crotonase family,72 the formation of enolate intermediates is often related to the presence of an oxyanion hole. 3.2. Second Reaction Stage: β-Hydroxyl Elimination. The subsequent step of the reaction required the elimination of the β-hydroxyl group to form the α,β-unsaturated acyl product (CRO, from crotonyl−CoA) and a water molecule as a byproduct.73 For the reaction to take place a structural rearrangement of His878 was observed after I1, in which the His878 made a hydrogen bond with the β-hydroxyl oxygen; at this point His878 could facilitate β-elimination either by electrostatic stabilization or through protonation. Consequently, the endergonic I1 led to an exergonic intermediate (I2) with a Gibbs reaction energy of −4.0 kcal·mol−1 relative to the React state. We observed that upon loss of the interaction between His878 and the Cα, the (CO)Cα−Cβ(OH) bond rotated from −27.3° to 19.9° to favor the interaction between the positively charged His878 and the β-hydroxyl, simultaneously shortening the distance between the α-carbonyl and the β-hydroxyl (from 1.73 to 1.55 Å, see Figure S7, Tables S4 and S5). The enolate species formed in I1 remained after the Cα−Cβ bond rotation, which suggested that I1 was in fact a transient intermediate. Hence, in comparison with I1, I2 presents two acidic residues (His878 and Asp1033) oriented to β-OH, which again suggested that an oxyanion hole could favor the β-elimination. However, despite the stabilization by the newly formed oxyanion hole (His878 and Asp1033), the β-hydroxyl was still covalently bonded to Cβ (Wiberg bond index of 0.22 au, close to that characteristic of single bonds ∼0.25 au). Therefore, the protonation of the β-hydroxyl oxygen should be a requisite for the β-elimination to occur. Contrary to the α-deprotonation, we observed that the acid species promoting β-elimination differed depending on whether His1037 was neutral or positively charged: when His1037 was positively charged, the Asp1033 acted as the acid to protonate the leaving β-hydroxyl oxygen, while when it was in the neutral form, the His878 acted as the acid. In an attempt to clarify how the His1037 protonation influenced the acidic role of both Asp1033 and His878 in the β-elimination step, we characterized the three transition states for this reaction (refer to Figure 5 for the TS2 of the His1037protonated model). At TS2, the β-hydroxyl oxygen was 1.93 Å

Figure 5. Stationary points for the β-elimination stage. Only the DFT layer is shown for simplicity. (Top) Second intermediate (I2). (Middle) Second transition state (TS2) with one imaginary frequency; yellow arrows (displacement vector) indicate the directions with largest atom displacement for the reaction coordinate. (Bottom) Products (Prod) of the reaction, showing the formed water molecule and CRO. Distances are indicated in Angstroms.

far from Cβ, and two identical carbon−carbon bonds (1.40 Å for Cα−Cβ and Cα−CO) were formed in HACδ−. Furthermore, NBO calculations suggested a rehybridization of Cβ from sp3 to sp2 (from I2 to TS2), while the same hybridization (sp2) remained for Cα and CO in the transient enolate structure. The results also revealed that in TS2 the leaving hydroxyl increased the charge and acquired a hydroxide character (atomic charge of −0.12 au vs −0.03 au in I2). When we performed a vibrational analysis on the imaginary frequencies of TS2 we observed that for the positively charged His1037 TS2 was characterized by an imaginary frequency of 313i cm−1, corresponding to an antisymmetric stretching of the Cβ···OβH, His878Hα···OβH, and Asp1033Hδ2···OβH distances (displacement ca. 0.92, 0.44, and 0.65 Å, respectively, see Figure 5). In addition, the displacement of the Asp1033···OβH 10272


Research Article

ACS Catalysis

React state (vs −0.7 and −6.5 kcal·mol−1, δ- and ε-nitrogenprotonated His1037, respectively). In addition, entropy contributions for the Gibbs energy throughout the chemical reaction were not substantial to establish the final thermodynamic profile of DH catalysis (ca. 1 kcal·mol−1, see Table S6 in the SI). Upon inclusion of the estimated energy differences between the δ- and the ε-protonated (ΔΔEδ→ε = ∼2.2 kcal·mol−1, calculated using the QM layer only) and the positively charged His1037 (−1.3 kcal·mol−1), it was clear that the β-elimination step was the rate-limiting step for the dehydration of HAC by DH. Including the energy differences between the several protonation states, the corrected Gibbs activation energies for the δ- and ε-protonated states became 16.0 and 14.8 kcal· mol−1. These barriers are slightly higher than the positively charged His1037; nonetheless, the difference was still small. Finally, we could conclude that even though the Gibbs activation barrier of the rate-limiting step was not dependent on the interactions with His1037, the pathway with Asp1033 and positively charged His1037 was the most favorable: the React state with positively charged His1037 was preferred over the neutral form (pKa of 8.34), and the overall Gibbs reaction energy was considerably more exergonic (>10 kcal·mol−1) than when His1037 was neutral and His878 was the acid. As mentioned, there have been proposals that highlighted the possible role of a water molecule in the β-elimination catalyzed by mDH-like enzymes (e.g., FabA).14,22 Our results did not indicate that a water molecule could be important to improve β-elimination. Nevertheless, experimental work such as kinetic isotope effects could be important to adequately evaluate this proposal. 3.3. Rate-Limiting Step. The enzyme turnover for the E. coli DH FAS (bDH) was proposed as 204 min−1 at 37 °C,14 which corresponds to a Gibbs activation energy of 17.4 kcal· mol−1, as calculated with transition state theory.75 Albeit we cannot make a direct comparison between bDH and mDH, the experimental studies evidenced that both forms shared similar structural and catalytic features, except that the active site in bDH had a Gln1036 instead of the His1037. This would mean that its pathway could be comparable to that found with His1037 in the neutral form (∼16 kcal·mol−1), ca. 1 kcal·mol−1 lower than the bDH experimental value. Furthermore, experimental data14 suggested that the kcat for the bDH remained very similar (changed from 204 to 187 min−1)14 upon Q1036H mutation. In agreement with our results, this observation supported that even though the protonation of His1037 did not influence the role of the catalytic Asp1033 in the reaction mechanism, neither His1037 nor Gln1036 was crucial for the kinetic outcome of the reaction mechanism of the HAC dehydration. Nevertheless, we concluded that His1037 was a key residue to promote the most favorable thermodynamics of the catalysis by mDH. Further mutagenesis studies would be required to clarify the role of other active site residues for the dehydration of HAC, namely, Gly888 or Tyr1003. Gly888 in particular was observed to be a critical residue for the feasibility of the α-deprotonation of HAC. In addition, we point out that our results do not evaluate the kinetics and thermodynamics of the formation of the enzyme:substrate complex or product exit. In this regard we cannot affirm that β-elimination is the rate-limiting step of the whole catalytic process. We have observed, nevertheless, that the β-elimination step is the rate-limiting step of the chemical

distance was ca. 0.2 Å larger than the His878···OβH distance. On the other hand, when His1037 was in the neutral state, TS2 was characterized by imaginary frequencies of 299i and 265i cm−1 (His1037-δ and His1037-ε, respectively) and the displacement of the His878···OβH and the Asp1033···OβH distances was similar for both neutral systems (see Figure S9 in the SI). These observations further supported the fact that positively charged His1037 could increase the acidic character of Asp1033. The chemical reaction proceeded to form the CRO product. We noted that the Cα−Cβ unsaturated bond was formed concomitantly with the hydroxide protonation by the acid species in a synchronous process (see the IRC calculations Figure S10 in the SI). Additionally, when the Asp1033 was the acid, a salt bridge was observed between the Asp1033− carboxylate (upon water formation) and the His1037− imidazolium (see Figure 5). This is a characteristic interaction observed in proteins,74 with a distance of 4.0 Å between the Oδ1 of Asp1033 and the geometric center of the His1037imidazolium (CHis1037) and an angle of ca. 128° between the Oδ1 of Asp1033, the CHis1037, and the Cγ of His1037. All three TSs provided similar Gibbs activation barriers (Figure 6). When Asp1033 was the acid, we calculated a Gibbs

Figure 6. (Top) Thermodynamic profile for the reaction mechanism studied: α-deprotonation and β-elimination. Profile represents the Gibbs energy when the Asp1033 acts as the acid in the β-elimination step (blue) or when His878 acts as the acid (δ-protonated His1037 in green; ε-protonated His1037 in yellow). Energies are in kcal·mol−1. Value highlighted in red corresponds to the estimated energetic cost calculated to tautomerize His1037 from δ- to ε-protonated (ΔΔEδ→ε). Difference between the neutral and the positively charged state of His1037 was calculated as −1.3 kcal·mol−1 (from pKa of 8.34 and a physiological pH of 7.4). (Bottom) Scheme of the stationary states throughout of the dehydration reaction. Atoms highlighted in red correspond to the highest displacement vectors in the TS1 and TS2. Intermediates (I1 and I2) are represented in one structure.

activation barrier of 14.1 kcal·mol−1, while when His878 was the acid we calculated barriers of 14.7 and 11.8 kcal·mol−1 for the δ- and ε-nitrogen-protonated His1037. The same was not observed for the calculated Gibbs reaction energies. We observed that the most exergonic reaction was that in which Asp1033 was the acid and His1037 was positively charged with a Gibbs reaction energy of −15.9 kcal·mol−1 relative to the 10273


Research Article

ACS Catalysis

Figure 7. Electrostatic contribution (ΔΔE‡) of individual residues within 15 Å of the DFT layer (20 Å for charged residues) for the rate-limiting step: (left) mDH domain in cartoon representation. DFT layer (white sticks) and charged residuesnegative sticks in red and pinkpositive sticks in blue and gray. Carbonyl carbon and β-hydroxyl are shown in yellow sticks. Dashed line indicates the separation between residues closer to the α-carbonyl than to the β-hydroxyl (regions I and II) and residues closer to the β-hydroxyl than to the α-carbonyl (regions III and IV). Water molecules were removed for simplicity. (Right) Electrostatic contribution of each residue plotted as a function of the “relative proximity coordinate, d”, where d corresponds to the difference between the distance of the residue to the α-carbonyl carbon (dCO) and the distance to the β-hydroxyl oxygen (dOH). dCO and dOH were calculated as the distance from a reference carbon (Cα for neutral residues, Cγ for Arg and Glu, Cδ for Asp, Cε for Lys) to the substrate α-CO carbon and the substrate β-OH oxygen. Activation energy of the rate-limiting step (TS2) was 13.1 kcal·mol−1.

would have the opposite effect. Our calculations seem to agree with the hypothesis (Figure 7) even though there were more charged residues (both positives and negatives) nearby the βhydroxyl comparatively to the α-carbonyl (see regions III and IV in Figure 7). In summary, the results suggested that positive residues closer to β-hydroxyl and negative residues closer to the α-carbonyl lower the barrier for the rate-limiting step (regions IV and II in Figure 7). The opposite was observed for residues with the inverse charge and relative proximity coordinates (regions III and I in Figure 7). Exceptions were observed for residues which were almost equidistant from both the α-carbonyl and the β-hydroxyl. In those cases, an atom-byatom calculation would be needed to discriminate the relative position of each atom of the residue. We observed this effect in a past work regarding HIV-1 protease,76 suggesting that this qualitative approach could be used for future design of rational mutations that can improve the catalytic performance of enzymes where a significant electronic density shift takes place during the rate-limiting step, as far as the mutations do not introduce significant folding changes in the protein. Furthermore, we noted that most charged residues were found at the outer part of mDH (as verified in Figure 7 and Figure S11 in the SI), which suggests that these residues could be more important to promote protein:protein interactions with neighboring domains. In this regard, the calculated contributions from Asp877 and Asp1002 as well as other exposed residues on mDH should be interpreted with caution. We observed also that there were no charged residues within 10 Å (see Figure S11) of the active site, which is consistent with the hydrophobic environment required to allocate the growing fatty acid.

reaction taking place at the mDH domain of FAS. Future studies can address the energetic details of HAC (un)binding to the active site of the domain as well as the exit of the CRO product, since these processes should also be accounted for to determine the rate-limiting step of the whole catalysis by mDH. In this regard we have also verified how feasible it would be to restore the initial protonation state of the active of mDH (deprotonation of His878 by Asp1033). To do that we tested two hypotheses depending on whether the water byproduct was present in the active site at that moment. Our calculations indicated that this process could easily occur (activation energy in range 3−6 kcal·mol−1); however, both the protonation with/without water byproduct resulted in an endergonic process. Although unexpected, these observations may suggest that restoring the active site of mDH may take place when there is an increase of the solvent available volume in the active site of mDH upon product unbinding. 3.4. Role of the Individual Residues on the Activation Energy Barrier Points to Strategies To Improve the Enzymatic Efficiency. The electrostatic energy contribution (ΔΔE‡) from each individual residue, within 15 Å of the reaction coordinate, to the rate-limiting step was estimated (details in section 3.1 in the SI). A positive ΔΔE‡ contribution indicated residues that increased the energy barrier for the ratelimiting step through electrostatic interactions (delayed the reaction), whereas negative values indicated residues that decreased this barrier through electrostatic interactions (fastened the reaction). Overall, we noted that the electrostatic contributions to the activation energy were mostly lower than ∼1.0 kcal·mol−1, resulting in a cumulative ΔΔE‡ of +7.2 kcal·mol−1,which suggests that the medium−long range electrostatic environment does not accelerate the rate-limiting step, and the enzyme electrostatic environment should mostly confer stability to the protein and the active site. We have observed also, previously, a shift of electron density as the reaction progresses from the negatively charged αcarbonyl carbon in I2 toward the β-hydroxyl at TS2. Consequently, we expected that positive residues closer to the α-carbonyl at I2 would render the reaction more unfavorable and negative residues in the same circumstance

4. CONCLUSION In this work we present the reaction mechanism of the dehydration catalyzed by the mDH domain of FAS which proceeds in two stages: (i) the α-carbon deprotonation of Cα and (ii) the β-carbon elimination of Cβ, concomitant with the formation of the α,β-unsaturated product (CRO). The residues involved in these stages are a catalytic dyad composed by the His878 and the Asp1033. Moreover, the literature proposes that a third residue can play an important role in the 10274


Research Article

ACS Catalysis

determined for the stationary point with five different DFT functional and the 6-311+G(2d,2p) basis set; energies are corrected at the B3LYP/6-31G(d):FF99SB level of theory; general mechanism proposed for the two steps involved in the mammalian DH; alternative mechanism is represented also (His1037 neutral); relevant distance variations between the residues involved in the formation of hydrogen bonds in the protonated system (positively charged); superposition of the two initial structures tested; Gly888-amide hydrogen forms/not forms the hydrogen bond with the HACcarbonyl oxygen; relaxed PES for the reaction coordinate followed for the α-deprotonation stage; variations of the absolute atomic charge along of the reaction pathway for the protonated system studied (positively charged); variations of the absolute atomic charge along the reaction pathway for the amino acids present in the DFT layer, for the δ-protonated and ε-protonated forms; relevant distance variations in the active site for the His1037 protonated in δ- and ε-nitrogen; distance variations in the hydrogen bond network for the His1037 protonated in δ- and ε-nitrogen; rate-limiting step (TS2) structures for the His1037 protonated in δand ε-nitrogen; intrinsic reaction coordinate for the second stage (β-elimination); thermodynamic energies for the α-deprotonation and the β-elimination stages for the three mechanisms followed; role of the individual residues on the activation energy barrier (details) (PDF)

catalysis (His1037). To clarify its influence on catalysis, we calculated the most probable protonation of His1037 with the CpHMD method and followed the mechanism considering its three protonation states. Our results indicate that the positively charged His1037 is the preferred form in mDH (pKa of 8.34). The His878 acts as a base by deprotonating the α-carbon of HAC, and the Asp1033 acts as an acid to promote βelimination of the leaving hydroxyl group of HAC to form the CRO and a water molecule as a byproduct. Alternatively, when His1037 is in the neutral form, His878 acts both as a base in αdeprotonation and as an acid during β-elimination, restoring its initial protonation state. Both mechanisms are plausible because they follow similar kinetics; however, the reaction where His1037 is positively charged and Asp1033 acts as an acid is also more favorable thermodynamically (in about 10 kcal·mol−1). The dehydration reaction was an exergonic process (ΔG of −15.9 kcal·mol−1) with the β-elimination stage as the ratelimiting step (ΔG‡ of 14.1 kcal·mol−1). The dehydration of HAC by mDH is strongly enhanced by an oxyanion hole formed with the α-carbonyl of HAC by its βhydroxyl and the backbone of Gly888 and another one formed between the β-hydroxyl of HAC and the side chain of the catalytic His878 and Asp1033. Our results emphasized that the former oxyanion hole is a key requisite to ensure the feasibility of the α-deprotonation, while the latter oxyanion hole is responsible for the stabilization of the resulting intermediate. Furthermore, we confirm that the active site of mDH is mostly composed by hydrophobic and polar residues, which should stabilize the growing fatty acid substrate, and that the electrostatic contribution of the enzyme environment for the kinetics of HAC’s dehydration is on average lower than 1 kcal· mol−1. Nevertheless, our results showed that there are much more positively and negatively charged residues closer to the βhydroxyl than the α-carbonyl of HAC. In conclusion, our results provide important atomistic, kinetic, and thermodynamic details of the reaction mechanism of fatty acid dehydration by mDH in the FAS multienzyme. Synergistically with experimental evidence built in the past years, we provide relevant rationalization and clarification of the chemical phenomena taking place, which can enhance the future development of drug design campaigns for cancer treatment.

■

DH domain highlighted the charged residues within a 20 Å of the β-hydroxyl nuclear coordinate; Gaussian input file for all optimized stationary states (CoordParmDH.zip)

■

Parameters (DH-pns-hac.inpcrd and DH-pns-hac.prmtop) for DH enzyme, HAC, and PNS molecules (OptDH.zip)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Fabiola E. Medina: 0000-0002-0230-8717 Rui P. P. Neves: 0000-0003-2032-9308 Maria J. Ramos: 0000-0002-7554-8324

ASSOCIATED CONTENT

S Supporting Information *

Notes

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acscatal.8b02616.

The authors declare no competing financial interest.

■

β-Sheet and α-helix representation of the two adjacent pseudodimeric domains in the mDH FAS; catalytic residues of the mammalian DH (PDB code 2VZ9) and the E. coli DH (PDB code 4KEH), classical molecular dynamics of the enzyme:substrate complex (details); RMSd per residue (His878, Asp1033, HAC, and PNS) along of the classical molecular dynamics; active site representations; relevant distances and dihedral angle proposed reaction coordinates along the cMD simulations: between the α-hydrogen of the HAC and the εnitrogen of the His878 and the oxygen hydroxyl of the HAC and the hydrogen of Asp1033; DFT layer of 120 atoms (details); representation of the DFT layers used: 94 and 120 atoms; CpHMD details; single-point energy

ACKNOWLEDGMENTS This work received financial support from the following institutions: European Union (FEDER funds POCI/01/0145/ FEDER/007728) and National Funds (FCT/MEC, Fundaçaõ para a Ciência e Tecnologia and Ministério da Educaçaõ e Ciência) under the Partnership Agreement PT2020 UID/ MULTI/04378/2013, NORTE-01-0145-FEDER-000024, supported by Norte Portugal Regional Operational Programme (NORTE 2020), under the PORTUGAL 2020 Partnership Agreement through the European Regional Development Fund (ERDF), and Fundaçaõ para a Ciência e a Tecnologia (FCT) through project PTDC/QUI-QFI/28714/2017. F.E.M acknowledges the Comisión Nacional de Investigación Cientı ́fica y Tecnológica−Chile (CONICYT) for a Ph.D. scholarship. 10275


Research Article

ACS Catalysis

■

(19) Labonte, J. W.; Townsend, C. A. Active Site Comparisons and Catalytic Mechanisms of the Hot Dog Superfamily. Chem. Rev. 2013, 113, 2182−2204. (20) Li, Y.; Dodge, G. J.; Fiers, W. D.; Fecik, R. A.; Smith, J. L.; Aldrich, C. C. Functional Characterization of a Dehydratase Domain from the Pikromycin Polyketide Synthase. J. Am. Chem. Soc. 2015, 137, 7003−7006. (21) Nguyen, C.; Haushalter, R. W.; Lee, D. J.; Markwick, P. R. L.; Bruegger, J.; Caldara-Festin, G.; Finzel, K.; Jackson, D. R.; Ishikawa, F.; O’Dowd, B.; McCammon, J. A.; Opella, S. J.; Tsai, S. C.; Burkart, M. D. Trapping the dynamic acyl carrier protein in fatty acid biosynthesis. Nature 2014, 505, 427−431. (22) Leesong, M.; Henderson, B. S.; Gillig, J. R.; Schwab, J. M.; Smith, J. L. Structure of a dehydratase-isomerase from the bacterial pathway for biosynthesis of unsaturated fatty acids: Two catalytic activities in one active site. Structure 1996, 4, 253−264. (23) Annand, R. R.; Kozlowski, J. F.; Davisson, V. J.; Schwab, J. M. Mechanism-based inactivation of Escherichia coli.beta.-hydroxydecanoyl thiol ester dehydrase: assignment of the imidazole nitrogen-15 NMR resonances and determination of the structure of the alkylated histidine. J. Am. Chem. Soc. 1993, 115, 1088−1094. (24) Joshi, A. K.; Smith, S. Construction, Expression, and Characterization of a Mutated Animal Fatty-Acid Synthase Deficient in the Dehydratase Function. J. Biol. Chem. 1993, 268, 22508−22513. (25) Biswas, R.; Dutta, A.; Dutta, D.; Hazra, D.; Banerjee, D. R.; Basak, A.; Das, A. K. Crystal structure of dehydratase component HadAB complex of mycobacterial FAS-II pathway. Biochem. Biophys. Res. Commun. 2015, 458, 369−374. (26) Dennington, R.; Keith, T.; Millam, J. GaussView, Version 5; Semichem Inc.: Shawnee Mission KS, 2009. (27) Case, D. A.; Darden, T. A.; Cheatham, T. E.; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Walker, R. C.; Zhang, W.; Merz, K. M.; Roberts, B.; Hayik, S.; Roitberg, A. E.; Seabra, G.; Swails, J.; Götz, A. W.; Kolossváry, I.; Wong, K. F.; Paesani, F.; Vanicek, J.; Wolf, R. M.; Liu, J.; Wu, X.; Brozell, S. R.; Steinbrecher, T.; Gohlke, H.; Cai, Q.; Ye, X.; Wang, J.; Hsieh, M. J.; Cui, G.; Roe, D. R.; Mathews, D. H.; Seetin, M. G.; Salomon-Ferrer, R.; Sagui, C.; Babin, V.; Luchko, T.; Gusarov, S.; Kovalenko, A.; Kollman, P. A. AMBER 12; University of California: San Francisco, 2012. (28) Wang, J. M.; 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. (29) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. A WellBehaved Electrostatic Potential Based Method Using Charge Restraints for Deriving Atomic Charges - the RESP model. J. Phys. Chem. 1993, 97, 10269−10280. (30) Roothaan, C. C. J. New Developments in Molecular Orbital Theory. Rev. Mod. Phys. 1951, 23, 69−89. (31) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Comparison of multiple amber force fields and development of improved protein backbone parameters. Proteins: Struct., Funct., Genet. 2006, 65, 712−725. (32) Lindorff-Larsen, K.; Piana, S.; Palmo, K.; Maragakis, P.; Klepeis, J. L.; Dror, R. O.; Shaw, D. E. Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins: Struct., Funct., Genet. 2010, 78, 1950−1958. (33) Olsson, M. H. M.; Søndergaard, C. R.; Rostkowski, M.; Jensen, J. H. PROPKA3: Consistent Treatment of Internal and Surface Residues in Empirical pKa Predictions. J. Chem. Theory Comput. 2011, 7, 525−537. (34) Søndergaard, C. R.; Olsson, M. H. M.; Rostkowski, M.; Jensen, J. H. Improved Treatment of Ligands and Coupling Effects in Empirical Calculation and Rationalization of pKa Values. J. Chem. Theory Comput. 2011, 7, 2284−2295. (35) Kim, M. O.; Nichols, S. E.; Wang, Y.; McCammon, J. A. Effects of histidine protonation and rotameric states on virtual screening of M. tuberculosis RmlC. J. Comput.-Aided Mol. Des. 2013, 27, 235−246.

REFERENCES

(1) Smith, S.; Tsai, S.-C. The type I fatty acid and polyketide synthases: a tale of two megasynthases. Nat. Prod. Rep. 2007, 24, 1041−1072. (2) Maier, T.; Leibundgut, M.; Ban, N. The crystal structure of a mammalian fatty acid synthase. Science 2008, 321, 1315−1322. (3) Viegas, M. F.; Neves, R. P. P.; Ramos, M. J.; Fernandes, P. A. Modeling of Human Fatty Acid Synthase and in Silico Docking of Acyl Carrier Protein Domain and Its Partner Catalytic Domains. J. Phys. Chem. B 2018, 122, 77−85. (4) McGillick, B. E.; Kumaran, D.; Vieni, C.; Swaminathan, S. betaHydroxyacyl-acyl Carrier Protein Dehydratase (FabZ) from Francisella tularensis and Yersinia pestis: Structure Determination, Enzymatic Characterization, and Cross-Inhibition Studies. Biochemistry 2016, 55, 1091−1099. (5) Zhang, L.; Xiao, J. F.; Xu, J. R.; Fu, T. R.; Cao, Z. W.; Zhu, L.; Chen, H. Z.; Shen, X.; Jiang, H. L.; Zhang, L. Crystal structure of FabZ-ACP complex reveals a dynamic seesaw-like catalytic mechanism of dehydratase in fatty acid biosynthesis. Cell Res. 2016, 26, 1330−1344. (6) Sacco, E.; Slama, N.; Backbro, K.; Parish, T.; Laval, F.; Daffe, M.; Eynard, N.; Quemard, A. Revisiting the Assignment of Rv0241c to Fatty Acid Synthase Type II of Mycobacterium tuberculosis. J. Bacteriol. 2010, 192, 4037−4044. (7) Ishikawa, F.; Haushalter, R. W.; Lee, D. J.; Finzel, K.; Burkart, M. D. Sulfonyl 3-Alkynyi Pantetheinamides as Mechanism-Based CrossLinkers of Acyl Carrier Protein Dehydratase. J. Am. Chem. Soc. 2013, 135, 8846−8849. (8) Bunkoczi, G.; Misquitta, S.; Wu, X. Q.; Lee, W. H.; Rojkova, A.; Kochan, G.; Kavanagh, K. L.; Oppermann, U.; Smith, S. Structural Basis for Different Specificities of Acyltransferases Associated with the Human Cytosolic and Mitochondrial Fatty Acid Synthases. Chem. Biol. 2009, 16, 667−675. (9) Bunkoczi, G.; Pasta, S.; Joshi, A.; Wu, X.; Kavanagh, K. L.; Smith, S.; Oppermann, U. Mechanism and substrate recognition of human holo ACP synthase. Chem. Biol. 2007, 14, 1243−1253. (10) Paiva, P.; Sousa, S. F.; Ramos, M. J.; Fernandes, P. A. Understanding the Catalytic Machinery and the Reaction Pathway of the Malonyl-Acetyl Transferase Domain of Human Fatty Acid Synthase. ACS Catal. 2018, 8, 4860−4872. (11) Robbins, T.; Kapilivsky, J.; Cane, D. E.; Khosla, C. Roles of Conserved Active Site Residues in the Ketosynthase Domain of an Assembly Line Polyketide Synthase. Biochemistry 2016, 55, 4476− 4484. (12) Hardwicke, M. A.; Rendina, A. R.; Williams, S. P.; Moore, M. L.; Wang, L.; Krueger, J. A.; Plant, R. N.; Totoritis, R. D.; Zhang, G.; Briand, J.; Burkhart, W. A.; Brown, K. K.; Parrish, C. A. A human fatty acid synthase inhibitor binds beta-ketoacyl reductase in the ketosubstrate site. Nat. Chem. Biol. 2014, 10, 774−779. (13) Medina, F. E.; Neves, R. P. P.; Ramos, M. J.; Fernandes, P. A. A QM/MM study of the reaction mechanism of human beta-ketoacyl reductase. Phys. Chem. Chem. Phys. 2017, 19, 347−355. (14) Pasta, S.; Witkowski, A.; Joshi, A. K.; Smith, S. Catalytic residues are shared between two pseudosubunits of the dehydratase domain of the animal fatty acid synthase. Chem. Biol. 2007, 14, 1377− 1385. (15) Rosenthal, R. G.; Vogeli, B.; Quade, N.; Capitani, G.; Kiefer, P.; Vorholt, J. A.; Ebert, M. O.; Erb, T. J. The use of ene adducts to study and engineer enoyl-thioester reductases. Nat. Chem. Biol. 2015, 11, 398−400. (16) Keatinge-Clay, A. T. The structures of type I polyketide synthases. Nat. Prod. Rep. 2012, 29, 1050−1073. (17) Maier, T.; Leibundgut, M.; Boehringer, D.; Ban, N. Structure and function of eukaryotic fatty acid synthases. Q. Rev. Biophys. 2010, 43, 373−422. (18) Dillon, S. C.; Bateman, A. The Hotdog fold: wrapping up a superfamily of thioesterases and dehydratases. BMC Bioinf. 2004, 5, 109. 10276


Research Article

ACS Catalysis (36) Uranga, J.; Mikulskis, P.; Genheden, S.; Ryde, U. Can the protonation state of histidine residues be determined from molecular dynamics simulations? Comput. Theor. Chem. 2012, 1000, 75−84. (37) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79, 926−935. (38) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. NumericalIntegration of Cartesian Equations of Motion of a Sysstem with Constrains - Molecular - Dynamics of N-alkanes. J. Comput. Phys. 1977, 23, 327−341. (39) Mongan, J.; Case, D. A.; McCammon, J. A. Constant pH molecular dynamics in generalized born implicit solvent. J. Comput. Chem. 2004, 25, 2038−2048. (40) Bashford, D.; Case, D. A. Generalized born models of macromolecular solvation effects. Annu. Rev. Phys. Chem. 2000, 51, 129−152. (41) Chung, L. W.; Sameera, W. M. C.; Ramozzi, R.; Page, A. J.; Hatanaka, M.; Petrova, G. P.; Harris, T. V.; Li, X.; Ke, Z.; Liu, F.; Li, H.-B.; Ding, L.; Morokuma, K. The ONIOM Method and Its Applications. Chem. Rev. 2015, 115, 5678−5796. (42) Chung, L. W.; Hirao, H.; Li, X.; Morokuma, K. The ONIOM method: its foundation and applications to metalloenzymes and photobiology. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 327− 350. (43) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö .; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, Revision E.01; Gaussian Inc.: Wallingford, CT, 2009. (44) Field, M. J.; Bash, P. A.; Karplus, M. A Combined QuantumMechanical and molecular Mechanical Potential for MolecularDynamics Simulations. J. Comput. Chem. 1990, 11, 700−733. (45) Ditchfield, R.; Hehre, W. J.; Pople, J. A. Self-Consistent Molecular-Orbital Methods.9. Extended Gaussian-Type Basis for Molecular-Orbital Studies of Organic Molecules. J. Chem. Phys. 1971, 54, 724. (46) Lee, C. T.; Yang, W. T.; Parr, R. G. Development of the ColleSalvetti Correlation -Energy Formula into a Functional of the Electron-Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (47) Becke, A. D. Density-Functional Thermochemistry.3. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (48) Ribeiro, A. J. M.; Yang, L.; Ramos, M. J.; Fernandes, P. A.; Liang, Z.-X.; Hirao, H. Insight into Enzymatic Nitrile Reduction: QM/MM Study of the Catalytic Mechanism of QueF Nitrile Reductase. ACS Catal. 2015, 5, 3740−3751. (49) Neves, R. P. P.; Fernandes, P. A.; Ramos, M. J. Unveiling the Catalytic Mechanism of NADP(+)-Dependent Isocitrate Dehydrogenase with QM/MM Calculations. ACS Catal. 2016, 6, 357−368. (50) Calixto, A. R.; Ramos, M. J.; Fernandes, P. A. Influence of Frozen Residues on the Exploration of the PES of Enzyme Reaction Mechanisms. J. Chem. Theory Comput. 2017, 13, 5486−5495. (51) Cantú Reinhard, F. G.; Faponle, A. S.; de Visser, S. P. Substrate Sulfoxidation by an Iron(IV)-Oxo Complex: Benchmarking Computationally Calculated Barrier Heights to Experiment. J. Phys. Chem. A 2016, 120, 9805−9814. (52) Fukui, K. The path of chemical reactions - the IRC approach. Acc. Chem. Res. 1981, 14, 363−368.

(53) Hirshfeld, F. L. Bonded-atom fragments for describing molecular charge densities. Theor. Chim. Acta 1977, 44, 129−138. (54) Ritchie, J. P. Electron-Density Distribution Analysis for Nitromethane, Nitromethide, and Nitramide. J. Am. Chem. Soc. 1985, 107, 1829−1837. (55) Ritchie, J. P.; Bachrach, S. M. Some Methods and Applications of Electron-Density Distribution Analysis. J. Comput. Chem. 1987, 8, 499−509. (56) Mayer, I. Bond order and valence indices: A personal account. J. Comput. Chem. 2007, 28, 204−221. (57) Glendening, E. D.; Weinhold, F. Natural Resonance Theory: II. Natural bond order and valency. J. Comput. Chem. 1998, 19, 610− 627. (58) Wiberg, K. B. Application of Pople-Santry-Segal CNDO Method to cyclopropylcarbinyl and cyclobutyl cation and to bicyclobutane. Tetrahedron 1968, 24, 1083−1096. (59) Neves, R. P. P.; Fernandes, P. A.; Ramos, M. J. Mechanistic insights on the reduction of glutathione disulfide by protein disulfide isomerase. Proc. Natl. Acad. Sci. U. S. A. 2017, 114, E4724−E4733. (60) Grimme, S. Supramolecular Binding Thermodynamics by Dispersion-Corrected Density Functional Theory. Chem. - Eur. J. 2012, 18, 9955−9964. (61) Becke, A. D. Density-Functional Exchange-Energy approximation with Correct Asymptotic-Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (62) Boese, A. D.; Martin, J. M. L. Development of density functionals for thermochemical kinetics. J. Chem. Phys. 2004, 121, 3405−3416. (63) Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120, 215−241. (64) Zhao, Y.; Truhlar, D. G. Hybrid meta density functional theory methods for thermochemistry, thermochemical kinetics, and noncovalent interactions: The MPW1B95 and MPWB1K models and comparative assessments for hydrogen bonding and van der Waals interactions. J. Phys. Chem. A 2004, 108, 6908−6918. (65) Chai, J. D.; Head-Gordon, M. Long-range corrected hybrid density functionals with damped atom-atom dispersion corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (66) Mangiatordi, G. F.; Bremond, E.; Adamo, C. DFT and Proton Transfer Reactions: A Benchmark Study on Structure and Kinetics. J. Chem. Theory Comput. 2012, 8, 3082−3088. (67) Kamerlin, S. C. L.; Chu, Z. T.; Warshel, A. On Catalytic Preorganization in Oxyanion Holes: Highlighting the Problems with the Gas-Phase Modeling of Oxyanion Holes and Illustrating the Need for Complete Enzyme Models. J. Org. Chem. 2010, 75, 6391−6401. (68) Meriläinen, G.; Poikela, V.; Kursula, P.; Wierenga, R. K. The Thiolase Reaction Mechanism: The Importance of Asn316 and His348 for Stabilizing the Enolate Intermediate of the Claisen Condensation. Biochemistry 2009, 48, 11011−11025. (69) Gao, D.; Zhan, C.-G. Modeling Effects of Oxyanion Hole on the Ester Hydrolysis Catalyzed by Human Cholinesterases. J. Phys. Chem. B 2005, 109, 23070−23076. (70) Sigala, P. A.; Kraut, D. A.; Caaveiro, J. M. M.; Pybus, B.; Ruben, E. A.; Ringe, D.; Petsko, G. A.; Herschlag, D. Testing Geometrical Discrimination within an Enzyme Active Site: Constrained Hydrogen Bonding in the Ketosteroid Isomerase Oxyanion Hole. J. Am. Chem. Soc. 2008, 130, 13696−13708. (71) Evans, C. T.; Kurz, L. C.; Remington, S. J.; Srere, P. A. Active site mutants of pig citrate synthase: effects of mutations on the enzyme catalytic and structural properties. Biochemistry 1996, 35, 10661−72. (72) Hamed, R. B.; Batchelar, E. T.; Clifton, I. J.; Schofield, C. J. Mechanisms and structures of crotonase superfamily enzymes − How nature controls enolate and oxyanion reactivity. Cell. Mol. Life Sci. 2008, 65, 2507−2527. 10277


Research Article

ACS Catalysis (73) Kim, S.; Thiessen, P. A.; Bolton, E. E.; Chen, J.; Fu, G.; Gindulyte, A.; Han, L. Y.; He, J. E.; He, S. Q.; Shoemaker, B. A.; Wang, J. Y.; Yu, B.; Zhang, J.; Bryant, S. H. PubChem Substance and Compound databases. Nucleic Acids Res. 2016, 44, D1202−D1213. (74) Donald, J. E.; Kulp, D. W.; DeGrado, W. F. Salt bridges: Geometrically specific, designable interactions. Proteins: Struct., Funct., Genet. 2011, 79, 898−915. (75) Laidler, K. J.; King, M. C. Development of transition-state theory. J. Phys. Chem. 1983, 87, 2657−2664. (76) Ribeiro, A. J. M.; Santos-Martins, D.; Russo, N.; Ramos, M. J.; Fernandes, P. A. Enzymatic Flexibility and Reaction Rate: A QM/ MM Study of HIV-1 Protease. ACS Catal. 2015, 5, 5617−5626.

10278


MM Study of the Reaction Mechanism of the Dehydratase Domain

Recommend Documents