Smooth solvation method for d-orbital semiempirical calculations of biological reactions. 2. Application to transphosphorylation thio effects in solution

Density-functional and semiempirical quantum methods and continuum dielectric and explicit solvation models are applied to study the role of solvation on the stabilization of native and thio-substituted transphosphorylation reactions. Extensive comparison is made between results obtained from the different methods. For the semiempirical methods, explicit solvation was treated using a hybrid quantum mechanical/molecular mechanical (QM/MM) approach and the implicit solvation was treated using a recently developed smooth solvation model implemented into a d-orbital semiempirical framework (MNDO/d-SCOSMO) within CHARMM. The different quantum and solvation methods were applied to the transesterification of a 3'-ribose,5'-methyl phosphodiester that serves as a nonenzymatic model for the self-cleavage reaction catalyzed by the hammerhead and hairpin ribozymes. Thio effects were studied for a double sulfur substitution at the nonbridging phosphoryl oxygen positions. The reaction profiles of both the native and double sulfur-substituted reactions from the MNDO/d-SCOSMO calculations were similar to the QM/MM results and consistent with the experimentally observed trends. These results underscore the need for a d-orbital semiempirical representation for phosphorus and sulfur for the study of experimentally observed thio effects in enzymatic and nonenzymatic phosphoryl transfer reactions. One of the major advantages of the present approach is that it can be applied to model chemical reactions at a significantly lower computational cost than either the density-functional calculations with implicit solvation or the semiempirical QM/MM simulations with explicit solvent.     

Hybrid QM/MM study of thio effects in transphosphorylation reactions: the role of solvation

Transphosphorylation thio effects in solution are studied using hybrid QM/MM calculations with a d-orbital semiempirical Hamiltonian. Activated dynamics simulations were performed for a 3' ribose-phosphate model in an explicit 20 A sphere of TIP3P water surrounded by a solvent boundary potential, and free energy analysis was performed using the weighted histogram analysis method. Single thio-substitutions at all of the phosphoryl oxygen positions and a double thio-substitution at the nonbridging positions were considered. The reaction free energy profiles are compared with available experimental data, and the role of solvation on the barrier heights and reaction coordinate is discussed. These results provide an important step in the characterization of thio effects in reactions of biological phosphates that may aid in the interpretation of kinetic data and ultimately help to unravel the catalytic mechanisms of ribozymes.     

Density functional study of the in-line mechanism of methanolysis of cyclic phosphate and thiophosphate esters in solution: insight into thio effects in RNA transesterification

Density functional calculations of thio effects on the in-line mechanism of methanolysis of ethylene phosphate (a reverse reaction model for RNA phosphate transesterification) are presented. A total of 12 reaction mechanisms are examined using the B3LYP functional with large basis sets, and the effects of solvation were treated using the PCM, CPCM, and SM5 solvation models. Single thio substitutions at all of the distinct phosphoryl oxygen positions (2', 3', 5', pro-R) and a double thio substitution at the nonbridging (pro-R/pro-S) positions were considered. Profiles for each reaction were calculated in the dianionic and monoanionic/monoprotic states, corresponding to reaction models under alkaline and nonalkaline conditions, respectively. These models provide insight into the mechanisms of RNA transesterification thio effects and serve as a set of high-level quantum data that can be used in the design of new semiempirical quantum models for hybrid quantum mechanical/molecular mechanical simulations and linear-scaling electronic structure calculations.     

Thio effects on the departure of the 3'-linked ribonucleoside from diribonucleoside 3',3'-phosphorodithioate diesters and triribonucleoside 3',3',5'-phosphoromonothioate triesters: implications for ribozyme catalysis

To provide a solid chemical basis for the mechanistic interpretations of the thio effects observed for large ribozymes, the cleavage of triribonucleoside 3',3',5'-phosphoromonothioate triesters and diribonucleoside 3',3'-phosphorodithioate diesters has been studied. To elucidate the role of the neighboring hydroxy group of the departing 3'-linked nucleoside, hydrolysis of 2',3'-O-methyleneadenosin-5'-yl bis[5'-O-methyluridin-3'-yl] phosphoromonothioate (1 a) has been compared to the hydrolysis of 2',3'-O-methyleneadenosin-5'-yl 5'-O-methyluridin-3'-yl 2',5'-di-O-methyluridin-3'-yl phosphoromonothioate (1 b) and the hydrolysis of bis[uridin-3'-yl] phosphorodithioate (2 a) to the hydrolysis of uridin-3'-yl 2',5'-di-O-methyluridin-3'-yl phosphorodithioate (2 b). The reactions have been followed by RP HPLC over a wide pH range. The phosphoromonothioate triesters 1 a,b undergo two competing reactions: the starting material is cleaved to a mixture of 3',3'- and 3',5'-diesters, and isomerized to 2',3',5'- and 2',2',5'-triesters. With phosphorodithioate diesters 2 a,b, hydroxide-ion-catalyzed cleavage of the P--O3' bond is the only reaction detected at pH >6, but under more acidic conditions desulfurization starts to compete with the cleavage. The 3',3'-diesters do not undergo isomerization. The hydroxide-ion-catalyzed cleavage reaction with both 1 a and 2 a is 27 times as fast as that compared with their 2'-O-methylated counterparts 1 b and 2 b. The hydroxide-ion-catalyzed isomerization of the 3',3',5'-triester to 2',3',5'- and 2',2',5'-triesters with 1 a is 11 times as fast as that compared with 1 b. These accelerations have been accounted for by stabilization of the anionic phosphorane intermediate by hydrogen bonding with the 2'-hydroxy function. Thio substitution of the nonbridging oxygens has an almost negligible influence on the cleavage of 3',3'-diesters 2 a,b, but the hydrolysis of phosphoromonothioate triesters 1 a,b exhibits a sizable thio effect, k(PO)/k(PS)=19. The effects of metal ions on the rate of the cleavage of diesters and triesters have been studied and discussed in terms of the suggested hydrogen-bond stabilization of the thiophosphorane intermediates derived from 1 a and 2 a.     

A quantum mechanical/molecular mechanical study on the catalysis of the pyridoxal 5'-phosphate-dependent enzyme L-serine dehydratase

The catalytic mechanism of a pyridoxal 5'-phosphate-dependent enzyme, l-serine dehydratase, has been investigated using ab initio quantum mechanical/molecular mechanical (QM/MM) methods. New insights into the chemical steps have been obtained, including the chemical role of the substrate carboxyl group in the Schiff base formation step and a proton-relaying mechanism involving the phosphate of the cofactor in the beta-hydroxyl-leaving step. The latter step is of no barrier and follows sequentially after the elimination of the alpha-proton, leading to a single but sequential alpha, beta-elimination step. The rate-limiting transition state is specifically stabilized by the enzyme environment. At this transition state, charges are localized on the substrate carboxyl group, as well as on the amino group of Lys41. Specific interactions of the enzyme environment with these groups are able to lower the activation barrier significantly. One major difficulty associated with studies of complicated enzymatic reactions using ab initio QM/MM models is the appropriate choices of reaction coordinates. In this study, we have made use of efficient semiempirical models and pathway optimization techniques to overcome this difficulty.     

Peptide hydrolysis catalyzed by matrix metalloproteinase 2: a computational study

The MMP-2 reaction mechanism is investigated by using different computational methodologies. First, quantum mechanical (QM) calculations are carried out on a cluster model of the active site bound to an Ace-Gly approximately Ile-Nme peptide. Along the QM reaction path, a Zn-bound water molecule attacks the Gly carbonyl group to give a tetrahedral intermediate. The breaking of the C-N bond is completed thanks to the Glu 404 residue that shuttles a proton from the water molecule to Ile-N atom. The gas-phase QM energy barrier is quite low ( approximately 14 kcal/mol), thus suggesting that the essential catalytic machinery is included in the cluster model. A similar reaction path occurs in the MMP-2 catalytic domain bound to an octapeptide substrate according to hybrid QM and molecular mechanical (QM/MM) geometry optimizations. However, the rupture of the Gly( P 1) approximately Ile( P 1') amide bond is destabilized in the static QM/MM calculations, owing to the positioning of the Ile( P 1') side chain inside the MMP-2 S 1' pocket and to the inability of simple energy miminization methodologies to properly relax complex systems. Molecular dynamics simulations show that these steric limitations are overcome easily through structural fluctuations. The energetic effect of structural fluctuations is taken into account by combining QM energies with average MM Poisson-Boltzmann free energies, resulting in a total free energy barrier of 14.8 kcal/mol in good agreement with experimental data. The rate-determining event in the MMP-2 mechanism corresponds to a H-bond rearrangement involving the Glu 404 residue and/or the Glu 404-COOH --> N-Ile( P 1') proton transfer. Overall, the present computational results and previous experimental data complement each other well in order to provide a detailed view of the MMPs catalytic mechanism.     

Mechanistic Insights into a Stibene Cleavage Oxygenase NOV1 from Quantum Mechanical/Molecular Mechanical Calculations

NOV1, a stilbene cleavage oxygenase, catalyzes the cleavage of the central double bond of stilbenes to two phenolic aldehydes, using a 4-His Fe(II) center and dioxygen. Herein, we use in-protein quantum mechanical/molecular mechanical (QM/MM) calculations to elucidate the reaction mechanism of the central double bond cleavage of phytoalexin resveratrol by NOV1. Our results showed that the oxygen molecule prefers to bind to the iron center in a side-on fashion, as suggested from the experiment. The quintet Fe−O₂ complex with the side-on superoxo antiferromagnetic coupled to the resveratrol radical is identified as the reactive oxygen species. The QM/MM results support the dioxygenase mechanism involving a dioxetane intermediate with a rate-limiting barrier of 10.0 kcal mol⁻¹. The alternative pathway through an epoxide intermediate is ruled out due to a larger rate-limiting barrier (26.8 kcal mol⁻¹). These findings provide important insight into the catalytic mechanism of carotenoid cleavage oxygenases and also the dioxygen activation of non-heme enzymes.     

A quantum mechanical investigation of possible mechanisms for the nucleotidyl transfer reaction catalyzed by DNA polymerase beta

Several quantum mechanical (QM) and hybrid quantum/molecular mechanical (QM/MM) studies have been employed recently to analyze the nucleotidyl transfer reaction in DNA polymerase beta (pol beta). Our examination reveals strong dependence of the reported mechanism on the initial molecular model. Thus, we explore here several model systems by QM methods to investigate pol beta's possible pathway variations. Although our most favorable pathway involves a direct proton transfer from O3'(primer) to O2alpha(Palpha), we also discuss other initial proton-transfer steps--to an adjacent water, to triphosphate, or to aspartic units--and the stabilizing effect of crystallographic water molecules in the active site. Our favored reaction route has an energetically undemanding initial step of less than 1.0 kcal/mol (at the B3LYP/6-31G(d,p) level), and involves a slight rearrangement in the geometry of the active site. This is followed by two major steps: (1) direct proton transfer from O3'(primer) to O2alpha(Palpha) leading to the formation of a pentavalent, trigonal bipyramidal Palpha center, via an associative mechanism, at a cost of about 28 kcal/mol, and (2) breakage of the triphosphate unit (exothermic process, approximately 22 kcal/mol) that results in the full transfer of the nucleotide to the DNA and the formation of pyrophosphate. These energy values are expected to be lower in the physical system when full protein effects are incorporated. We also discuss variations from this dominant pathway, and their impact on the overall repair process. Our calculated barrier for the chemical reaction clearly indicates that chemistry is rate-limiting overall for correct nucleotide insertion in pol beta, in accord with other studies. Protonation studies on relevant intermediates suggest that, although protonation at a single aspartic residue may occur, the addition of a second proton to the system significantly disturbs the active site. We conclude that the active site rearrangement step necessary to attain a reaction-competent geometry is essential and closely related to the "pre-chemistry" avenue described recently as a key step in the overall kinetic cycle of DNA polymerases. Thus, our work emphasizes the many possible ways for DNA polymerase beta's chemical reaction to occur, determined by the active site environment and initial models.     

Hydrolysis of chemically distinct sites of human serum albumin by polyoxometalate: A hybrid QM/MM (ONIOM) study

In this study, mechanisms of hydrolysis of all four chemically diverse cleavage sites of human serum albumin (HSA) by [Zr(OH)(PW₁₁O₃₉)]⁴⁻ (ZrK) have been investigated using the hybrid two-layer QM/MM (ONIOM) method. These reactions have been proposed to occur through the following two mechanisms: internal attack (IA) and water assisted (WA). In both mechanisms, the cleavage of the peptide bond in the Cys392-Glu393 site of HSA is predicted to occur in the rate-limiting step of the mechanism. With the barrier of 27.5 kcal/mol for the hydrolysis of this site, the IA mechanism is found to be energetically more favorable than the WA mechanism (barrier = 31.6 kcal/mol). The energetics for the IA mechanism are in line with the experimentally measured values for the cleavage of a wide range of dipeptides. These calculations also suggest an energetic preference (Cys392-Glu393, Ala257-Asp258, Lys313-Asp314, and Arg114-Leu115) for the hydrolysis of all four sites of HSA. © 2018 Wiley Periodicals, Inc.     

Variational electrostatic projection (VEP) methods for efficient modeling of the macromolecular electrostatic and solvation environment in activated dynamics simulations

New methods for the calculation of electrostatic interactions between the active dynamical region and surrounding external solvated macromolecular environment in hybrid quantum mechanical/molecular mechanical (QM/MM) simulations are presented. The variational electrostatic projection (VEP) method, and related variational reverse-mapping procedure (VEP-RVM) utilize an expansion in Gaussian surface elements for representation of electrostatic interactions. The use of a discretized surface that surrounds the active dynamical region greatly reduces the number of interactions with the particles of the external environment. The methods are tested on two catalytic RNA systems: the hammerhead and the hairpin ribozymes. It is shown that with surface elements numbering from 302 to 1202 points the direct VEP and VEP-RVM methods are able to obtain relative force errors in the range of 0.5-0.05% and 0.09-0.0001%, respectively, using a 4.0 A projection buffer. These results are encouraging and provide an essential step in the development of new variational macromolecular solvent-boundary methods for QM/MM calculations of enzyme reactions.     

Ground-state structures and vertical excitations for the kindling fluorescent protein asFP595

Geometry configurations of a large fraction of the kindling fluorescent protein asFP595 around the chromophore region were optimized by using the effective fragment potential quantum mechanical-molecular mechanical (QM/MM) method. The initial coordinates of heavy atoms were taken from the structure from the Protein Data Bank archive corresponding to the dark-adapted state of the Ala143 --> Gly mutant of asFP595. Optimization of geometry parameters was performed for all internal coordinates in the QM part composed of the chromophore unit and the side chains of His197, Glu215, and Arg92 as well as for positions of effective fragments constituting the MMpart. The structures corresponding to the anion trans, anion cis, and zwitterion trans moieties were considered among various alternatives for the chromophore unit inside the protein matrix. The QM/MM simulations show that the protein environment provides stabilization for the trans-zwitterion isomer compared to the gas-phase conditions. By using the multiconfigurational CASSCF and the time-dependent density functional theory calculations, we estimated positions of spectral bands corresponding to vertical S(0)-S(1) transitions. The results of simulations support the assumption that the dark state of asFP595 corresponds to the anionic or zwitterionic trans-conformation, while the kindled state corresponds to the anionic cis-conformation.     

Performance assessment of semiempirical molecular orbital methods in describing halogen bonding: quantum mechanical and quantum mechanical/molecular mechanical-molecular dynamics study

The performance of semiempirical molecular-orbital methods--MNDO, MNDO-d, AM1, RM1, PM3 and PM6--in describing halogen bonding was evaluated, and the results were compared with molecular mechanical (MM) and quantum mechanical (QM) data. Three types of performance were assessed: (1) geometrical optimizations and binding energy calculations for 27 halogen-containing molecules complexed with various Lewis bases (Two of the tested methods, AM1 and RM1, gave results that agree with the QM data.); (2) charge distribution calculations for halobenzene molecules, determined by calculating the solvation free energies of the molecules relative to benzene in explicit and implicit generalized Born (GB) solvents (None of the methods gave results that agree with the experimental data.); and (3) appropriateness of the semiempirical methods in the hybrid quantum-mechanical/molecular-mechanical (QM/MM) scheme, investigated by studying the molecular inhibition of CK2 protein by eight halobenzimidazole and -benzotriazole derivatives using hybrid QM/MM molecular-dynamics (MD) simulations with the inhibitor described at the QM level by the AM1 method and the rest of the system described at the MM level. The pure MM approach with inclusion of an extra point of positive charge on the halogen atom approach gave better results than the hybrid QM/MM approach involving the AM1 method. Also, in comparison with the pure MM-GBSA (generalized Born surface area) binding energies and experimental data, the calculated QM/MM-GBSA binding energies of the inhibitors were improved by replacing the G(GB,QM/MM) solvation term with the corresponding G(GB,MM) term.     

Quantum mechanics/molecular mechanics minimum free-energy path for accurate reaction energetics in solution and enzymes: sequential sampling and optimization on the potential of mean force surface

To accurately determine the reaction path and its energetics for enzymatic and solution-phase reactions, we present a sequential sampling and optimization approach that greatly enhances the efficiency of the ab initio quantum mechanics/molecular mechanics minimum free-energy path (QM/MM-MFEP) method. In the QM/MM-MFEP method, the thermodynamics of a complex reaction system is described by the potential of mean force (PMF) surface of the quantum mechanical (QM) subsystem with a small number of degrees of freedom, somewhat like describing a reaction process in the gas phase. The main computational cost of the QM/MM-MFEP method comes from the statistical sampling of conformations of the molecular mechanical (MM) subsystem required for the calculation of the QM PMF and its gradient. In our new sequential sampling and optimization approach, we aim to reduce the amount of MM sampling while still retaining the accuracy of the results by first carrying out MM phase-space sampling and then optimizing the QM subsystem in the fixed-size ensemble of MM conformations. The resulting QM optimized structures are then used to obtain more accurate sampling of the MM subsystem. This process of sequential MM sampling and QM optimization is iterated until convergence. The use of a fixed-size, finite MM conformational ensemble enables the precise evaluation of the QM potential of mean force and its gradient within the ensemble, thus circumventing the challenges associated with statistical averaging and significantly speeding up the convergence of the optimization process. To further improve the accuracy of the QM/MM-MFEP method, the reaction path potential method developed by Lu and Yang [Z. Lu and W. Yang, J. Chem. Phys. 121, 89 (2004)] is employed to describe the QM/MM electrostatic interactions in an approximate yet accurate way with a computational cost that is comparable to classical MM simulations. The new method was successfully applied to two example reaction processes, the classical SN2 reaction of Cl-+CH3Cl in solution and the second proton transfer step of the reaction catalyzed by the enzyme 4-oxalocrotonate tautomerase. The activation free energies calculated with this new sequential sampling and optimization approach to the QM/MM-MFEP method agree well with results from other simulation approaches such as the umbrella sampling technique with direct QM/MM dynamics sampling, demonstrating the accuracy of the iterative QM/MM-MFEP method.     

Mechanism of proteolysis in matrix metalloproteinase‐2 revealed by QM/MM modeling

The mechanism of enzymatic peptide hydrolysis in matrix metalloproteinase‐2 (MMP‐2) was studied at atomic resolution through quantum mechanics/molecular mechanics (QM/MM) simulations. An all‐atom three‐dimensional molecular model was constructed on the basis of a crystal structure from the Protein Data Bank (ID: 1QIB), and the oligopeptide Ace‐Gln‐Gly～Ile‐Ala‐Gly‐Nme was considered as the substrate. Two QM/MM software packages and several computational protocols were employed to calculate QM/MM energy profiles for a four‐step mechanism involving an initial nucleophilic attack followed by hydrogen bond rearrangement, proton transfer, and C? N bond cleavage. These QM/MM calculations consistently yield rather low overall barriers for the chemical steps, in the range of 5–10 kcal/mol, for diverse QM treatments (PBE0, B3LYP, and BB1K density functionals as well as local coupled cluster treatments) and two MM force fields (CHARMM and AMBER). It, thus, seems likely that product release is the rate‐limiting step in MMP‐2 catalysis. This is supported by an exploration of various release channels through QM/MM reaction path calculations and steered molecular dynamics simulations. © 2015 Wiley Periodicals, Inc.     

Antibiotic binding to dizinc beta-lactamase L1 from Stenotrophomonas maltophilia: SCC-DFTB/CHARMM and DFT studies

A dizinc beta-lactamase (L1 from Stenotrophomonas maltophilia) complexed with an antibiotic compound (moxalactam) has been studied using a hybrid quantum mechanical/molecular mechanical (QM/MM) approach. The QM region is described by the self-consistent charge-density functional tight binding (SCC-DFTB) model while the MM by CHARMM. The Michaelis complex, which is constructed from a recent X-ray structure of the L1 enzyme with the hydrolyzed moxalactam, is simulated by molecular dynamics. The simulation yields valuable insights into substrate-enzyme interaction, whose implications in the enzyme catalysis are discussed. Finally, the QM/MM results are compared with a high-level density functional theory study of a truncated active-site model and the agreement provides strong support for the SCC-DFTB treatment of the QM region.     

Reaction path potential for complex systems derived from combined ab initio quantum mechanical and molecular mechanical calculations

Combined ab initio quantum mechanical and molecular mechanical calculations have been widely used for modeling chemical reactions in complex systems such as enzymes, with most applications being based on the determination of a minimum energy path connecting the reactant through the transition state to the product in the enzyme environment. However, statistical mechanics sampling and reaction dynamics calculations with a combined ab initio quantum mechanical (QM) and molecular mechanical (MM) potential are still not feasible because of the computational costs associated mainly with the ab initio quantum mechanical calculations for the QM subsystem. To address this issue, a reaction path potential energy surface is developed here for statistical mechanics and dynamics simulation of chemical reactions in enzymes and other complex systems. The reaction path potential follows the ideas from the reaction path Hamiltonian of Miller, Handy and Adams for gas phase chemical reactions but is designed specifically for large systems that are described with combined ab initio quantum mechanical and molecular mechanical methods. The reaction path potential is an analytical energy expression of the combined quantum mechanical and molecular mechanical potential energy along the minimum energy path. An expansion around the minimum energy path is made in both the nuclear and the electronic degrees of freedom for the QM subsystem internal energy, while the energy of the subsystem described with MM remains unchanged from that in the combined quantum mechanical and molecular mechanical expression and the electrostatic interaction between the QM and MM subsystems is described as the interaction of the MM charges with the QM charges. The QM charges are polarizable in response to the changes in both the MM and the QM degrees of freedom through a new response kernel developed in the present work. The input data for constructing the reaction path potential are energies, vibrational frequencies, and electron density response properties of the QM subsystem along the minimum energy path, all of which can be obtained from the combined quantum mechanical and molecular mechanical calculations. Once constructed, it costs much less for its evaluation. Thus, the reaction path potential provides a potential energy surface for rigorous statistical mechanics and reaction dynamics calculations of complex systems. As an example, the method is applied to the statistical mechanical calculations for the potential of mean force of the chemical reaction in triosephosphate isomerase.     

Prediction of the mechanisms of enzyme-catalysed reactions using hybrid quantum mechanical/molecular mechanical methods

The use of hybrid methods, involving both quantum mechanics and molecular mechanics, to model the mechanism of enzyme-catalysed reactions, is discussed. Two alternative approaches to treating the electrostatic interactions between the quantum mechanical and molecular mechanical regions are studied, involving either the inclusion of this term in the electronic Hamiltonian (QM/MM), or evaluating it purely classically (MO + MM). In the latter scheme, possible problems of using force fields that are standard for macromolecular modelling are identified. The use of QM/MM schemes to investigate the mechanism of the enzymes thymidine phosphorylase (ThdPase) and protein tyrosine phosphatase (PTP) is described. For both systems, transition states have been identified using a PM3 Hamiltonian. For ThdPase, concerted motion of the enzyme during the course of the reaction is suggested and, for PTP, a two-step dephosphorylation reaction is indicated, both with quite low barriers.