Journal of Molecular Biology
Computational Modeling of the Catalytic Reaction in Triosephosphate Isomerase
Introduction
Triosephosphate isomerase (TIM) is one of the most extensively characterized enzymes in the chemical literature. TIM catalyses the reversible isomerization of dihydroxyacetone phosphate (DHAP) to d-glyceraldehyde 3-phosphate (GAP) with exceptionally high efficiency, while suppressing elimination of orthophosphate. To accomplish this isomerization, TIM first extracts a pro-R hydrogen from the C1 of DHAP and then stereospecifically introduces a proton onto the C2 atom (Figure 1). A wealth of thermodynamic, spectroscopic, kinetic, and structural studies have been carried out, providing extensive experimental data with which to compare theoretical calculations. Furthermore, the well-known large-scale catalytic loop motion in TIM (loop 6, residues 166–177), which alternates between open and closed positions, poses a substantial challenge to computational modeling efforts, one that is different from modeling enzyme reactions with relatively small-scale structural changes in the protein over the course of the catalytic cycle. The closed form of the loop is observed whenever the enzyme is ligated, and is responsible for properly aligning active-site residues fundamental to the catalytic process and for the stabilization of reaction intermediates.1
The conversion of DHAP to GAP, a thermodynamically unfavorable but biologically significant reaction, occurs with a catalytic rate and a Michaelis constant The rate-limiting step appears to be a slow conformational change as described by a variety of isotope data.2 In the reverse direction, with GAP as a substrate, kcat=8.7(±0.3)×103 s−1, and the rate-limiting step appears to be a proton transfer.3., 4. Therefore, the opening of the loop on DHAP appears to be faster than 8×103 s−1, whereas on GAP it must be comparable to 7×102 s−1. Recently, these kinetic measurements have been complemented with a variety of spectroscopic kinetic probes of the loop-opening rates in the presence of a variety of ligands.5., 6., 7. Information on the chemical mechanism and the key players at the active-site has been provided from a variety of enzymological, mutagenic and spectroscopic studies. Mutational changes of Lys128 and Glu1659 have important effects in both Km and kcat, whereas the His95Q mutant reduces the catalytic rate by 200-fold.10., 11. A large collection of crystal structures, including the recent TIM-DHAP Michaelis complex at near-atomic resolution, 1NEY.pdb,12 indicate a highly compact active-site. Spectroscopic and X-ray crystallographic data indicate that the substrate is primarily in the form of the ketone DHAP in situ. Short-range interactions to all three key residues, Glu165, His95 and Lys12, are seen in the ground state or Michaelis complex. It is interesting, based on this high-resolution picture of the ground state, to ask how stabilization is provided for the key intermediates. This issue is best addressed through computational studies, and must combine predictions of quantum chemical pathways with predictions of conformational freedom for key residues.
Several theoretical studies, using semiempirical and DFT-based QM/MM methods to study TIM catalysis, have provided an initial atomic-level picture of the various reaction steps.13., 14., 15., 16., 17. However, these studies were not able to utilize the recent TIM-DHAP crystal structure discussed above as a starting point. The resulting structural uncertainties, in conjunction with issues concerning the accuracy of the computational technology employed in the calculations, made it difficult to confidently draw quantitative conclusions. Furthermore, no investigation of the catalytic side-chains or loop motion was carried out. Thus, a complete, robust theoretical picture of TIM catalysis has not been produced.
Here, we use a combination of QM/MM and protein structure modeling techniques to provide a comprehensive, and quantitatively reliable, description of the TIM catalytic cycle. Our QM/MM methods enable routine treatment of quantum regions of the active-site as large as 150–200 atoms, for which full DFT geometry optimizations can be carried out with sufficiently large basis sets to guarantee accurate structures. The QM/MM interface parameterization that we employ has been tested extensively in a wide range of prior calculations on small model systems and proteins;18., 19. we have demonstrated that the errors arising from the interface are less than 1 kcal/mol Furthermore, the overall accuracy in calculating free energies of intermediates and transition states for enzymatic reactions appears to be of the order of 2–3 kcal/mol. This technology represents a significant improvement in accuracy when compared with models and methods used in previous theoretical studies, thus enabling reliable quantitative conclusions to be drawn concerning barrier heights, protonation states, and other details of the reaction. As discussed below, the protonation state of the phosphate group on the ligand is of particular importance, and we have invested considerable effort in considering the alternatives and assessing them according to total energies obtained.
We investigate the catalytic loop motion in TIM using methods for loop and side-chain prediction,20., 21., 22., 23. based on a protein molecular mechanics force-field and a surface Generalized Born (SGB) continuum solvation model that we have developed several years ago.24 The methods have been validated through an extensive series of tests in predicting loop and side-chain structures in the PDB, and have been shown to yield a remarkably high level of accuracy, even for relatively large loops (the catalytic TIM loop is 11 residues, which is within the validation range of our test suite). QM/MM calculations are used to generate charges for the ligand, which is likely to be a substantial improvement over using generic force-field parameters, particularly for the charged phosphate group.
The structural prediction methods are able to successfully generate open and closed forms of the loop, with and without the ligand present, with thermodynamics that are consistent with the entire range of experimental data. These calculations provide an initial atomic level picture of the loop opening with the ligand present, which has been proposed as the rate-determining step in the TIM catalytic cycle.
The specific reaction catalyzed by TIM that we have chosen to study is given below in Equation (1):The thermodynamic stability of the product and reactant in aqueous solution has been measured experimentally, and yields a ΔG of 1.9 kcal/mol.25 Thus, the products are somewhat less stable than the reactants. In vivo, the product is essential for optimal throughput in the glycolytic pathway. Depletion of product continues to drive transformation of the reactants. Without catalysis, the reaction rate in solution is ∼10−7 s−1,26 which would be inadequate to serve the relevant biological function. TIM accelerates the rate by about nine orders of magnitude.
The key residues involved in the catalysis are believed to be Glu165, His95, and Lys12 (Figure 1). It is widely accepted that the first step involves the transfer of a proton from C1 of DHAP to the catalytic base, Glu165 (atom notation can be followed with Figure 1, Figure 2). From this point, three different paths have been proposed and discussed by several experimental and theoretical groups. Paths A and C involve the formation of an enediol intermediate by transfer of a proton from His95 (path A) or Glu165 (path C) to the substrate O2 atom. The enediol intermediate then forms enediolate by the transfer of a proton from O1 to His95 (path A) or Glu165 (path C). Path A, the “classical” mechanism, was proposed by Knowles et al.4 and has been accepted widely. Path C has been supported recently by NMR experiments,27 and by ab initio active-site theoretical models.17 Path B13 forms an enediolate intermediate directly by internal proton transfer from O1 to O2. The final catalytic step, common to all paths, involves a proton transfer from Glu165 to C2 to form the product GAP. The proposed steps for the reaction should account for central experimental evidence: the substantial isotope washout when 2H and 3H are used as substrate.28
In previous theoretical studies, the enzymatic cycle was studied using the dianionic (D state) form of the ligand phosphate group. This is the thermodynamically most stable form of the ligand in neutral solution, as the second pKa of the phosphate group is ∼6.5. However, the free energy difference between the singly (monoanionic, M state) and doubly charged species would be less than 1 kcal/mol; furthermore, partial burial of the phosphate group in the enzyme active-site would tend to favor the protonated species, as the outer shell Born energy of the dianion (substantially larger than that of a single charged species) would be diminished. Thus, a pKa shift from the solution value is expected. In what follows, we study the above reactions using both the singly (M) and doubly (D) charged forms of the phosphate group. We did briefly investigate the neutral species as well, which produce a substantial reduction in activation barrier; however, the pKa of this last proton is extremely low, ∼2–3, making it highly unlikely that this species contributes substantially to the reaction mechanism. Note that even if the effect of the protein on pKa is negligibly small and the dianion substrate is the major bound complex,29 the singly charged form can still dominate the optimum reaction pathway; all that is required is that the overall activation free energy of this pathway is lower than that achieved by going through the doubly charged species exclusively. In formulating this argument, we assume that protonation is diffusion-controlled, so that the kinetic barrier to forming the monoanion from the dianion at pH 7 would be significantly lower than that observed for the overall reaction.
Section snippets
Results
We first introduce the QM/MM study of the catalytic chemistry and discuss the different proposed paths. We then describe our predictions of the opening and closing of the catalytic loop for several different enzymatic structures, in the absence of substrate and in the presence of DHAP and GAP.
QM/MM modeling of the catalytic chemistry
We first discuss the effects of QM/MM geometry optimization of the crystal structure. These calculations reveal a closer interaction between Glu165 and the substrate than that observed in the crystal structure, positioning the carboxylate group of Glu165 in an optimal geometry for the abstraction of H1. The experimental structure indicated a well-defined electron density, with a low Debye Waller factor for both Lys12 and His95, in agreement with the minor adjustment of both residues in the
Conclusions
We have presented a comprehensive picture of TIM-mediated catalysis of the substrates DHAP and GAP at an atomic level of detail. Excellent agreement with the experimental data is obtained at both a qualitative and (when available) quantitative level. Among the key new findings of this work are: (1) the identification of the monoanionic form of the substrate phosphate group as the principal channel for abstraction of the initial proton from DHAP; (2) an energy-based description of the opening
QM/MM protocol
The QM/MM model has been derived from the recent experimental structure, 1NEY.pdb.12 The full protein, a dimer, has been included in the calculations, although only the active-site of the B chain has been treated at the QM level. The total system contains approximately 8000 atoms. All crystallographic water molecules within 20 Å of the substrate are incorporated into the QM/MM system. Both the dianionic and the monoanionic, D and M, states are obtained initially after a 20 ps molecular dynamics
Acknowledgements
The authors thank Dr Sharon Rozovsky for helpful advice in the preparation of the models. This work was supported, in part, by a grant from the NIH to R.A.F. (GM40526 and GM 52018) and to A.M. (GM66388). M.P.J. acknowledges support from NSF grant 0302445.
References (38)
- et al.
The time scale of the catalytic loop motion in triosephosphate isomerase
J. Mol. Biol.
(2001) - et al.
Solution-state NMR investigations of triosephosphate isomerase active site loop motion: ligand release in relation to active site loop dynamics
J. Mol. Biol.
(2001) - et al.
Computer simulation of the triosephosphate isomerase catalyzed reaction
J. Biol. Chem.
(1996) - et al.
On the role of the crystal environment in determining protein side-chain conformations
J. Mol. Biol.
(2002) - et al.
Extending the accuracy limits of prediction for side-chain conformations
J. Mol. Biol.
(2001) Enzyme catalysis—not different, just better
Nature
(1991)- et al.
Energetics of triosephosphate isomerase—appearance of solvent tritium in substrate dihydroxyacetone phosphate and in product
Biochemistry
(1976) - et al.
Triosephosphate isomerase—energetics of the reaction catalyzed by the yeast enzyme expressed in Escherichia coli
Biochemistry
(1988) - et al.
Free-energy profile for reaction catalyzed by triosephosphate isomerase
Biochemistry
(1976) - et al.
Active site loop motion in triosephosphate isomerase: T-jump relaxation spectroscopy of thermal activation
Biochemistry
(2003)
Triosephosphate isomerase requires a positively charged active-site—the role of lysine-12
Biochemistry
Crystal structure of the mutant yeast triosephosphate isomerase in which the catalytic base glutamic acid 165 is changed to aspartic acid
Biochemistry
Electrophilic catalysis in triosephosphate isomerase—the role of histidine 95
Biochemistry
The structural basis for pseudoreversion of the H95N lesion by the secondary S96P mutation in triosephosphate isomerase
Biochemistry
Optimal alignment for enzymatic proton transfer: structure of the Michaelis complex of triosephosphate isomerase at 1.2 Å resolution
Proc. Natl Acad. Sci. USA
The intramolecular mechanism for the second proton transfer in triosephosphate isomerase (TIM): A QM/FE approach
J. Comput. Chem.
Computer-simulation and analysis of the reaction pathway of triosephosphate isomerase
Biochemistry
Quantum mechanical/molecular mechanical studies of the triosephosphate isomerase-catalyzed reaction: verification of methodology and analysis of reaction mechanisms
J. Phys. Chem. B
Ab initio models for receptor-ligand interactions in proteins. 4. Model assembly study of the catalytic mechanism of triosephosphate isomerase
Proteins: Struct. Funct. Genet.
Cited by (74)
Comments on 'Flexibility of enzymatic transitions as a hallmark of optimized enzyme steady-state kinetics and thermodynamics'
2021, Computational Biology and ChemistryChemical Kinetics: From Molecular Structure to Chemical Reactivity
2021, Chemical Kinetics: From Molecular Structure to Chemical ReactivityMapping Conformational Dynamics to Individual Steps in the TEM-1 β-Lactamase Catalytic Mechanism
2018, Journal of Molecular BiologyCitation Excerpt :For instance, CPMG relaxation dispersion NMR (arguably the most powerful approach for monitoring active enzymes) requires that the reaction being monitored have an equilibrium constant close to 1 so that steady state turnover can be maintained indefinitely while measurements are taken [2–6]. Computational models have yielded insights to a much broader set of reactions, but are often limited by short timescales and a lack of experimental parameterization [7–9]. The inability to acquire a detailed picture of catalysis-linked dynamics has real implications for drug development, limiting the ability to design inhibitors targeted at specific, catalytically required dynamic modes.
Porphyrin binding to Gun4 protein, facilitated by a flexible loop, controls metabolite flow through the chlorophyll biosynthetic pathway
2015, Journal of Biological ChemistryCitation Excerpt :This is actually what happens in Gun4, where simulations with crystal mates produce loops with a ∼1–2 Å alpha carbon RMSD with respect to the crystallographically determined structures (hardly distinguishable from the experimental structures). Taken together, these results on both systems point to a strong bias in the loop conformation mediated by crystallographic contacts; thus, the conformation of functionally important loops in Gun4 cannot be assessed purely on the basis of crystallographic data, and study/correction with in silico methods (now a well established practice; Ref. 35) is required. Importantly, this bias, originating from the structural data for crystallographic Gun4, produces a conformation with the loop restricting the access to the expected binding pocket.
Mechanism of hydrogen production in [Fe-Fe]-hydrogenases: A quantum mechanics/molecular mechanics study
2010, International Journal of Hydrogen EnergyCitation Excerpt :QSite [38] couples Jaguar (for the QM region) and the protein-modeling program IMPACT [38] (for the MM region). The theoretical foundation for QSite’s methodology [43–45] enabled numerous applications [46–60]. The QM-MM boundary at a covalent bond is described by frozen orbitals [45].
Quantum chemical descriptors based on semiempirical methods for large biomolecules
2023, Journal of Chemical Physics
- †
Present address: V. Guallar, Department of Biochemistry and Molecular Biophysics, Washington School of Medicine, 660 South Euclid Avenue, St. Louis, MO 63110, USA.