A generalized hybrid orbital (GHO) method has been developed at the semiempirical level in combined quantum mechanical and molecular mechanical (QM/MM) calculations. In this method, a set of hybrid orbitals is placed on the boundary atom between the QM and MM fragments, and one of the hybrid orbitals participates in the SCF calculation for the atoms in the QM region. The GHO method provides a well-defined potential energy surface for a hybrid QM/MM system and is a significant improvement over the “link-atom” approach by saturating the QM valencies with hydrogen atoms. The method has been tested on small molecules and yields reasonable structural, energetic, and electronic results in comparison with the results of the corresponding QM and MM approximations. The GHO method will greatly increase the applicability of combined QM/MM methods to systems comprising large molecules, such as proteins.

We present a theoretical framework for the calculation of rate constants of enzyme-catalyzed reactions that combines variational optimization of the dynamical bottleneck for overbarrier reactive events and multidimensional quantum mechanical tunneling dynamics for through-barrier reactive events, both in the presence of the protein environment. The theory features a two-zone, three-stage procedure called ensemble-averaged variational transition state theory with multidimensional tunneling (EA-VTST/MT) with the transmission coefficient based on the equilibrium secondary-zone (ESZ) approximation for including the effects of the protein on a catalytic reaction center, called the primary zone. The dynamics is calculated by canonical variational theory with optimized multidimensional tunneling contributions, and the formalism allows for Boltzmann averaging over an ensemble of reactant and transition state conformations. In the first stage of the calculations, we assume that the generalized transition states can be well described by a single progress coordinate expressed in primary-zone internal coordinates; in subsequent steps, the transmission coefficient is averaged over a set of primary-zone reaction paths that depend on the protein configuration, and each reaction path has its own reaction coordinate and optimized tunneling path. We also present a simpler approximation to the transmission coefficient that is called the static secondary-zone (SSZ) approximation. We illustrate both versions of this method by carrying out calculations of the reaction rate constants and kinetic isotope effects for oxidation of benzyl alcoholate to benzaldehyde by horse liver alcohol dehydrogenase. The potential energy surface is modeled by a combined generalized hybrid orbital/quantum mechanical/molecular mechanical/semiempirical valence bond (GHO-QM/MM/SEVB) method. The multidimensional tunneling calculations are microcanonically optimized by employing both the small-curvature tunneling approximation and version 4 of the large-curvature tunneling approximation. We find that the variation of the protein mean force as a function of reaction coordinate is quantitatively significant, but it does not change the qualitative conclusions for the present reaction. We obtain good agreement with experiment for both kinetic isotope effects and Swain−Schaad exponents.

One of the strongest experimental indications of hydrogen tunneling in biology has been the elevated Swain-Schaad exponent for the secondary kinetic isotope effect in the hydride-transfer step catalyzed by liver alcohol dehydrogenase. This process has been simulated using canonical variational transition-state theory for overbarrier dynamics and optimized multidimensional paths for tunneling. Semiclassical quantum effects on the dynamics are included on a 21-atom substrate-enzyme-coenzyme primary zone embedded in the potential of a substrate-enzyme-coenzyme-solvent secondary zone. The potential energy surface is calculated by treating 54 atoms by quantum mechanical electronic structure methods and 5506 protein, coenzyme, and solvent atoms by molecular mechanical force fields. We find an elevated Swain-Schaad exponent for the secondary kinetic isotope effect and generally good agreement with other experimental observables. Quantum mechanical tunneling is calculated to account for ∼60% of the reactive flux, confirming the dominance of tunneling that was inferred from the Swain-Schaad exponent. The calculations provide a detailed picture of the origin of the kinetic isotope effect and the nature of the tunneling process.

We present an overview of new procedures for including quantum mechanical effects in enzyme kinetics. Quantum effects are included in three ways: (1) The electronic structure of the atoms in the catalytic center is treated quantum mechanically in order to calculate a realistic potential energy surface for the bond rearrangement process. (2) The discrete nature of quantum mechanical vibrational energies is incorporated in the treatment of nuclear motion for computing the potential of mean force. (3) Multidimensional tunneling contributions are included. These procedures are illustrated by applications to proton abstractions catalyzed by enolase and methylamine dehydrogenase and hydride-transfer reactions by alcohol dehydrogenase and xylose isomerase.

ABSTRACT:This paper provides an overview of a new method developed to include quantum mechanical effects and free energy sampling in calculations of reaction rates in enzymes. The paper includes an overview of variational transition state theory with optimized multidimensional tunneling for simple gas-phase reactions and then shows how this is extended to incorporate free energy effects and to include protein motions in the reaction coordinate by ensemble averaging. Finally we summarize recent comparisons to experiment for primary and secondary kinetic isotope effects for proton and hydride transfer reactions catalyzed by enzymes.

We have calculated the reaction rate and kinetic isotope effects for conversion of 2-phospho-d-glycerate to phosphoenolpyruvate by yeast enolase. The potential energy surface is modeled by a combined quantum mechanical/molecular mechanical method with generalized hybrid orbitals. The dynamics calculations are carried out by semiclassical variational transition state theory with multidimensional tunneling contributions. Quantum effects are included for a 25-atom cluster consisting of the substrate and part of the protein embedded in a rigid framework consisting of the rest of the protein and water. Quantum effects are important for calculating the absolute rate constant, and variational optimization of the dynamical bottleneck location is important for calculating the kinetic isotope effects. This provides the first evidence that transition state geometries are isotope dependent for enzyme reactions.

We have applied molecular dynamics umbrella-sampling simulation and ensemble-averaged variational transition state theory with multidimensional tunneling (EA-VTST/MT) to calculate the reaction rate of xylose-to- xylulose isomerization catalyzed by xylose isomerase in the presence of two Mg2+ ions. The calculations include determination of the free energy of activation profile and ensemble averaging in the transmission coefficient. The potential energy function is approximated by a combined QM/MM/SVB method involving PM3 for the quantum mechanical (QM) subsystem, CHARMM22 and TIP3P for the molecular mechanical (MM) environment, and a simple valence bond (SVB) local function of two bond distances for the hydride transfer reaction. The simulation confirms the essential features of a mechanism postulated on the basis of kinetics and X-ray data by Whitlow et al. (Whitlow, M.; Howard, A. J.; Finzel, B. C.; Poulos, T. L.; Winborne, E.; Gilliland, G. L. Proteins 1991, 9, 153) and Ringe, Petsko, and coworkers (Labie, A.; Allen, K.-N.; Petsko, G. A.; Ringe, D. Biochemistry 1994, 33, 5469). This mechanism involves a rate-determining 1,2-hydride shift with prior and post proton transfers. Inclusion of quantum mechanical vibrational energy is important for computing the free energy of activation, and quantum mechanical tunneling effects are essential for computing kinetic isotope effects (KIEs). It is found that 85% of the reaction proceeds by tunneling and 15% by overbarrier events. The computed KIE for the ratio of hydride to deuteride transfer is in good agreement with the experimental results. The molecular dynamics simulations reveal that proton and hydride transfer reactions are assisted by breathing motions of the mobile Mg2+ ion in the active site, providing evidence for concerted motion of Mg2+ during the hydride transfer step.

Classical molecular dynamics and Monte Carlo simulations typically exclude quantum effects on the vibrations of reactants and transition states, and this may lead to significant errors in the computed potential of mean force. To correct this deficiency, a simple approximate procedure is proposed for the inclusion of quantum-mechanical vibrational energy in the computation of reactive potentials of mean force in condensed phases. The method is illustrated by a hydrogen atom transfer and a proton transfer reaction in water, in particular, the 1,5-sigmatropic shift in malonaldehyde and the intermolecular proton shift between ammonium ion and ammonia in an encounter complex. In both cases, quantum-mechanical vibrational energy makes significant contributions by reducing the free energy of activation by 2 to 3 kcal/mol. This finding has important implications in developing empirical potential functions for the study of enzyme reactions, and it is essential to quantize vibrational energy in the computed potential of mean force and free energy of activation in order to compare simulations quantitatively with experiment.