Multiscale Quantum-Classical Simulations of Enzymes

Enzymes are remarkably efficient catalysts evolved to perform well-defined and highly specific chemical transformations. Studying the nature of enzymatic rate enhancements is highly important from several aspects, including the rational design of synthetic catalysts and transition-state inhibitors. Herein, we describe recent progress in our group in the development of multiscale simulation methods and their application to several enzyme systems. In particular, we describe the use of combined quantum mechanics/molecular mechanics (QM/MM) methods in classical and quantum simulations. The development of various novel path-integral methods is reviewed. These methods are tailor-made for enzyme systems, where only a few degrees of freedom involved in the chemistry need to be quantized. The application of the hybrid QM/MM quantum-classical simulation approach to three case studies is presented. The first case involves proton transfer in nitroalkane oxidase, where the enzyme employs tunneling as a catalytic fine-tuning tool. The second case presented involves orotidine 5′-monophosphate decarboxylase, where multidimensional free energy simulations together with kinetic isotope effects are combined in the study of the reaction mechanism. Finally, we discuss the monoterpene cyclase bornyl diphosphate synthase, where non-statistical dynamics is a key component in enzyme function.     

Towards more realistic computational modeling of homogenous catalysis by density functional theory: combined QM/MM and ab initio molecular dynamics

The combined quantum mechanics/molecular mechanics (QM/MM) and the ab initio molecular dynamics methods (AIMD) are fast emerging as viable computational molecular modeling tools. Both methods allow for the incorporation of effects that are often ignored in high level calculations, but may be critical to the real chemistry of the simulated system. In the combined QM/MM method part of the system, say the active site, is treated quantum mechanically whereas the remainder of the system is treated with a faster molecular mechanics force field. This allows high level calculations to be performed where the effects of the environment are incorporated in a computationally tractable manner. With the ab initio molecular dynamics methods, the system is simulated at a finite temperature with no empirical force field. Rather, the forces at each time step are determined with a full electronic structure calculation at the density functional level. Thus, simulations of chemical reactions can be performed where finite temperature effects are realistically represented. In this paper a brief introduction to both methods is given. The methods are further demonstrated with specific applications to modeling homogenous catalytic processes at the molecular level. These applications are our latest efforts to build more realistic computational models of catalytic systems at the density functional level.     

Valence Bond-Based Hybrid Quantum Mechanics Molecular Mechanics Approaches and Proper Inclusion of the Effect of the Surroundings

The pioneering development of multiscale models for complex chemical systems by Karplus, Levitt, and Warshel, including the hybrid quantum mechanics molecular mechanics (QM/MM) approach and its application to enzymes, established a new field in chemistry that allows the modeling of reactivity within complex chemical systems. Inspired by the potential of such methods, many groups developed different QM/MM variants. Valence bond (VB) theory, which always was and still is an important conceptual tool for chemists, is well suited to deal with problems of chemical reactivity. Hence, here we review VB-based QM/MM methods, including the early semi-empirical methods that utilize VB concepts and more recent ab initio VB-based QM/MM methods. Special emphasis is given to the different ways to include effects of the surroundings on the solute. It is shown that within the VB framework, simple mechanical embedding for each diabatic state, followed by mixing of the states, accounts for most of these effects.     

On the Influence of the Protonation States of Active Site Residues on AChE Reactivation: A QM/MM Approach

Acetylcholinesterase (AChE), an enzyme of the serine hydrolase superfamily, is a mediator of signal transmission at cholinergic synapses by catalyzing acetylcholine cleavage into acetate and choline. This enzyme is vulnerable to covalent inhibition by organophosphate compounds (like VX). Covalent inhibition of AChE does not revert spontaneously. Known reactivator compounds have limited action in restoring catalytic activity. QM/MM simulations of VX‐inhibited AChE reactivation by pralidoxime (2‐PAM), a classical reactivator, were performed. These afforded a broad view of the effect of protonation states of active‐site residues, and provide evidence for the role of Glu202, which needs to be protonated for reactivation to occur. In situ deprotonation of 2‐PAM for both protonation states of Glu202 showed that His447 is able to deprotonate 2‐PAM with the assistance of Glu202. Because the active site of serine hydrolases is highly conserved, this work provides new insights on the interplay between the catalytic triad residues and this glutamate, newly identified as protonatable.     

Atomistic-Level Description of the Covalent Inhibition of SARS-CoV-2 Papain-like Protease

Inhibition of the papain-like protease (PLpro) of SARS-CoV-2 has been demonstrated to be a successful target to prevent the spreading of the coronavirus in the infected body. In this regard, covalent inhibitors, such as the recently proposed VIR251 ligand, can irreversibly inactivate PLpro by forming a covalent bond with a specific residue of the catalytic site (Cys¹¹¹), through a Michael addition reaction. An inhibition mechanism can therefore be proposed, including four steps: (i) ligand entry into the protease pocket; (ii) Cys¹¹¹ deprotonation of the thiol group by a Brønsted–Lowry base; (iii) Cys¹¹¹-S⁻ addition to the ligand; and (iv) proton transfer from the protonated base to the covalently bound ligand. Evaluating the energetics and PLpro conformational changes at each of these steps could aid the design of more efficient and selective covalent inhibitors. For this aim, we have studied by means of MD simulations and QM/MM calculations the whole mechanism. Regarding the first step, we show that the inhibitor entry in the PLpro pocket is thermodynamically favorable only when considering the neutral Cys¹¹¹, that is, prior to the Cys¹¹¹ deprotonation. For the second step, MD simulations revealed that His²⁷² would deprotonate Cys¹¹¹ after overcoming an energy barrier of ca. 32 kcal/mol (at the QM/MM level), but implying a decrease of the inhibitor stability inside the protease pocket. This information points to a reversible Cys¹¹¹ deprotonation, whose equilibrium is largely shifted toward the neutral Cys¹¹¹ form. Although thermodynamically disfavored, if Cys¹¹¹ is deprotonated in close proximity to the vinylic carbon of the ligand, then covalent binding takes place in an irreversible way (third step) to form the enolate intermediate. Finally, due to Cys¹¹¹-S⁻ negative charge redistribution over the bound ligand, proton transfer from the initially protonated His²⁷² is favored, finally leading to an irreversibly modified Cys¹¹¹ and a restored His²⁷². These results elucidate the selectivity of Cys¹¹¹ to enable formation of a covalent bond, even if a weak proton acceptor is available, as His²⁷².     

Stereoselective hydride transfer by aryl-alcohol oxidase, a member of the GMC superfamily

Primary alcohol oxidation by aryl‐alcohol oxidase (AAO), a flavoenzyme providing H₂O₂ to ligninolytic peroxidases, is produced by concerted proton and hydride transfers, as shown by substrate and solvent kinetic isotope effects (KIEs). Interestingly, when the reaction was investigated with synthesized (R)‐ and (S)‐α‐deuterated p‐methoxybenzyl alcohol, a primary KIE (≈6) was observed only for the R enantiomer, revealing that the hydride transfer is highly stereoselective. Docking of p‐methoxybenzyl alcohol at the buried crystal active site, together with QM/MM calculations, showed that this stereoselectivity is due to the position of the hydride‐ and proton‐receiving atoms (flavin N5 and His502 Nε, respectively) relative to the alcohol Cα‐substituents, and to the concerted nature of transfer (the pro‐S orientation corresponding to a 6 kcal mol^?1 penalty with respect to the pro‐R orientation). The role of His502 is supported by the lower activity (by three orders of magnitude) of the H502A variant. The above stereoselectivity was also observed, although activities were much lower, in AAO reactions with secondary aryl alcohols (over 98 % excess of the R enantiomer after treatment of racemic 1‐(p‐methoxyphenyl)ethanol, as shown by chiral HPLC) and especially with use of the F501A variant. This variant has an enlarged active site that allow better accommodation of the α‐substituents, resulting in higher stereoselectivity (S/R ratios) than is seen with AAO. High enantioselectivity in a member of the GMC oxidoreductase superfamily is reported for the first time, and shows the potential for engineering of AAO for deracemization purposes.     

Solvent-mediated proton transfer in catalysis by carbonic anhydrase

Considerable attention has been focused on proton transfer through intervening water molecules in complex macromolecules of biological interest, such as bacteriorhodopsin, cytochrome c oxidase, and many others. Proton transfer in catalysis by carbonic anhydrase provides a useful model for the study of the properties of such proton translocations. High-resolution X-ray crystallography in combination with measurements of catalysis have revealed new details of this process. A prominent proton shuttle residue His64 shows evidence of structural mobility, which appears to enhance proton transfer between the active site and bulk solvent. Moreover, the properties of the imidazole side chain of His64, including its conformations and pK(a), are finely tuned by surrounding residues of the active-site cavity. The structure of a network of ordered solvent molecules located between His64 and the active site are also sensitive to surrounding residues. These features combine to provide efficient proton-transfer rates as great as 10(6) s(-1) necessary to sustain rapid catalysis.     

Virtual Alanine Scan of the Main Protease Active Site in Severe Acute Respiratory Syndrome Coronavirus 2

Recently, inhibitors of the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) main protease (Mpro) have been proposed as potential therapeutic agents for COVID-19. Studying effects of amino acid mutations in the conformation of drug targets is necessary for anticipating drug resistance. In this study, with the structure of the SARS-CoV-2 Mpro complexed with a non-covalent inhibitor, we performed molecular dynamics (MD) simulations to determine the conformation of the complex when single amino acid residue in the active site is mutated. As a model of amino acid mutation, we constructed mutant proteins with one residue in the active site mutated to alanine. This method is called virtual alanine scan. The results of the MD simulations showed that the conformation and configuration of the ligand was changed for mutants H163A and E166A, although the structure of the whole protein and of the catalytic dyad did not change significantly, suggesting that mutations in His163 and Glu166 may be linked to drug resistance.     

Drug Repurposing of the Unithiol: Inhibition of Metallo-β-Lactamases for the Treatment of Carbapenem-Resistant Gram-Negative Bacterial Infections

The increasing antibiotic resistance is a clinical problem worldwide. Numerous Gram-negative bacteria have already become resistant to the most widely used class of antibacterial drugs, β-lactams. One of the main mechanisms is inactivation of β-lactam antibiotics by bacterial β-lactamases. Appearance and spread of these enzymes represent a continuous challenge for the clinical treatment of infections and for the design of new antibiotics and inhibitors. Drug repurposing is a prospective approach for finding new targets for drugs already approved for use. We describe here the inhibitory potency of known detoxifying antidote 2,3-dimercaptopropane-1-sulfonate (unithiol) against metallo-β-lactamases. Unithiol acts as a competitive inhibitor of meropenem hydrolysis by recombinant metallo-β-lactamase NDM-1 with the K_I of 16.7 µM. It is an order of magnitude lower than the K_I for l-captopril, the inhibitor of angiotensin-converting enzyme approved as a drug for the treatment of hypertension. Phenotypic methods demonstrate that the unithiol inhibits natural metallo-β-lactamases NDM-1 and VIM-2 produced by carbapenem-resistant K. pneumoniae and P. aeruginosa bacterial strains. The 3D full atom structures of unithiol complexes with NDM-1 and VIM-2 are obtained using QM/MM modeling. The thiol group is located between zinc cations of the active site occupying the same place as the catalytic hydroxide anion in the enzyme–substrate complex. The sulfate group forms both a coordination bond with a zinc cation and hydrogen bonds with the positively charged residue, lysine or arginine, responsible for proper orientation of antibiotics upon binding to the active site prior to hydrolysis. Thus, we demonstrate both experimentally and theoretically that the unithiol is a prospective competitive inhibitor of metallo-β-lactamases and it can be utilized in complex therapy together with the known β-lactam antibiotics.     

Combination of computational prescreening and experimental library construction can accelerate enzyme optimization by directed evolution

Chiral compounds can be produced efficiently by using biocatalysts. However, wild-type enzymes often do not meet the requirements of a production process, making optimization by rational design or directed evolution necessary. Here, we studied the lipase-catalyzed hydrolysis of the model substrate 1-(2-naphthyl)ethyl acetate both theoretically and experimentally. We found that a computational equivalent of alanine scanning mutagenesis based on QM/MM methodology can be applied to identify amino acid positions important for the activity of the enzyme. The theoretical results are consistent with concomitant experimental work using complete saturation mutagenesis and high-throughput screening of the target biocatalyst, a lipase from Bacillus subtilis. Both QM/MM-based calculations and molecular biology experiments identify histidine 76 as a residue that strongly affects the catalytic activity. The experiments demonstrate its important influence on enantioselectivity.     

Poly(ADP-ribose)-polymerase-catalyzed hydrolysis of NAD+: QM/MM simulation of the enzyme reaction

Poly(ADP-ribose) polymerase (PARP) is a nuclear enzyme which uses NAD+ as substrate and catalyzes the transfer of multiple units of ADP-ribose to target proteins. PARP is an attractive target for the discovery of novel therapeutic agents and PARP inhibitors are currently evaluated for the treatment of a variety of pathological conditions such as brain ischemia, inflammation, and cancer. Herein, we use the PARP-catalyzed reaction of NAD+ hydrolysis as a model for gaining insight into the molecular details of the catalytic mechanism of PARP. The reaction has been studied in both the gas-phase and in the enzyme environment through a QM/MM approach. Our results indicate that the cleavage reaction of the nicotinamide-ribosyl bond proceeds through an SN2 dissociative mechanism via an oxacarbenium transition structure. These results confirm the importance of the structural water molecule in the active site and may constitute the basis for the design of transition-state-based PARP inhibitors.     

The mechanism of phosphorylation of natural nucleosides and anti-HIV analogues by nucleoside diphosphate kinase is independent of their sugar substituents

The reaction mechanism of the phosphoryl transfer catalyzed by dinucleoside diphosphate kinase from Dictyostelium discoideum is investigated by semiempirical AM1 molecular orbital computation of an active site model system on the basis of various X-ray crystallographic structures. The computational results suggest that the phosphoryl transfer from adenosine triphosphate to the His122 residue is accompanied by the simultaneous shift of a proton from the histidine residue to one of the oxygen atoms of the gamma phosphate group. This involves a doubly protonated His122 residue whilst this residue is neutral in its ternary complex with ADP and the transition state analogue AlF(3). The proposed mechanism is thus analogous to that of phosphoryl transfer by cyclic adenosine monophosphate dependent protein kinase and uridine/cytidine monophosphate kinase as found in our earlier work and clarifies the role of the ribose 3'-OH group. Furthermore, the energetics of phosphoryl transfer onto other nucleoside analogues such as 3'-azido-3'-deoxythymidine-diphosphate and 2',3'-dideoxy-2',3'-didehydro-thymidine-diphosphate are investigated. The calculated reaction barriers for the phosphorylation of the diphosphates by the enzyme are all within a range of 13.1 kJ mol(-1), which suggests that variations in the activation energies alone cannot account for the experimentally observed differences in enzymatic activity. Consequences for the design of new anti-HIV nucleoside analogues are discussed. Supporting information for this article is available on the WWW under http://www.wiley-vch.de/contents/jc_2268/2002/f360_s.pdf or from the author.     

Mechanisms of Proteolytic Enzymes and Their Inhibition in QM/MM Studies

Experimental evidence for enzymatic mechanisms is often scarce, and in many cases inadvertently biased by the employed methods. Thus, apparently contradictory model mechanisms can result in decade long discussions about the correct interpretation of data and the true theory behind it. However, often such opposing views turn out to be special cases of a more comprehensive and superior concept. Molecular dynamics (MD) and the more advanced molecular mechanical and quantum mechanical approach (QM/MM) provide a relatively consistent framework to treat enzymatic mechanisms, in particular, the activity of proteolytic enzymes. In line with this, computational chemistry based on experimental structures came up with studies on all major protease classes in recent years; examples of aspartic, metallo-, cysteine, serine, and threonine protease mechanisms are well founded on corresponding standards. In addition, experimental evidence from enzyme kinetics, structural research, and various other methods supports the described calculated mechanisms. One step beyond is the application of this information to the design of new and powerful inhibitors of disease-related enzymes, such as the HIV protease. In this overview, a few examples demonstrate the high potential of the QM/MM approach for sophisticated pharmaceutical compound design and supporting functions in the analysis of biomolecular structures.     

High‐Throughput Virtual Screening Using Quantum Mechanical Probes: Discovery of Selective Kinase Inhibitors

A procedure based on semi‐empirical quantum mechanical (QM) calculations of interaction energy is proposed for the rapid screening of compound poses generated by high‐throughput docking. Small molecules (consisting of 2–10 atoms and termed “probes”) are overlapped with polar groups in the binding site of the protein target. The interaction energy values between each compound pose and the probes, calculated by a semi‐empirical Hamiltonian, are used as filters. The QM probe method does not require fixed partial charges and takes into account polarization and charge‐transfer effects which are not captured by conventional force fields. The procedure is applied to screen ～100 million poses (of 2.7 million commercially available compounds) obtained by high‐throughput docking in the ATP binding site of the tyrosine kinase erythropoietin‐producing human hepatocellular carcinoma receptor B4 (EphB4). Three QM probes on the hinge region and one at the entrance pocket are employed to select for binding affinity, while a QM probe on the side chain of the so‐called gatekeeper residue (a hypervariable residue in the kinome) is used to enforce selectivity. The poses with favorable interactions with the five QM probes are filtered further for hydrophobic matching and low ligand strain. In this way, a single‐digit micromolar inhibitor of EphB4 with a relatively good selectivity profile is identified in a multimillion‐compound library upon experimental tests of only 23 molecules.     

D-Tyrosine as a chiral precusor to potent inhibitors of human nonpancreatic secretory phospholipase A2 (IIa) with antiinflammatory activity

Few reported inhibitors of secretory phospholipase A(2) enzymes truly inhibit the IIa human isoform (hnpsPLA(2)-IIa) noncovalently at submicromolar concentrations. Herein, the simple chiral precursor D-tyrosine was derivatised to give a series of potent new inhibitors of hnpsPLA(2)-IIa. A 2.2-A crystal structure shows an inhibitor bound in the active site of the enzyme, chelated to a Ca(2+) ion through carboxylate and amide oxygen atoms, H-bonded through an amide NH group to His48, with multiple hydrophobic contacts and a T-shaped aromatic-group-His6 interaction. Antiinflammatory activity is also demonstrated for two compounds administered orally to rats.     

Computational studies of the primary phototransduction event in visual rhodopsin

This Account addresses recent advances in the elucidation of the detailed molecular rearrangements due to the primary photochemical event in rhodopsin, a prototypical G-protein-coupled receptor (GPCR) responsible for the signal transmission cascade in the vertebrate vision process. The reviewed studies provide fundamental insight on long-standing problems regarding the assembly and function of the individual residues and bound water molecules that form the rhodopsin active site, a center that catalyzes the 11-cis/all-trans isomerization of the retinyl chromophore in the primary step of the phototransduction mechanism. Emphasis is placed on the authors' recent computational studies, based on state-of-the-art quantum mechanics/molecular mechanics (QM/MM) hybrid methods, addressing the structural refinement of the retinyl chromophore binding site in high-resolution X-ray structures of bovine visual rhodopsin, the energy storage mechanism, and the molecular origin of spectroscopic changes due to the primary photochemical event.     

The cooperative effect between active site ionized groups and water desolvation controls the alteration of acid/base catalysis in serine proteases

What is the driving force that alters the catalytic function of His57 in serine proteases between general base and general acid in each step along the enzymatic reaction? The stable tetrahedral complexes (TC) of chymotrypsin with trifluoromethyl ketone transition state analogue inhibitors are topologically similar to the catalytic transition state. Therefore, they can serve as a good model to study the enzyme catalytic reaction. We used DFT quantum mechanical calculations to analyze the effect of solvation and of polar factors in the active site of chymotrypsin on the pK_a of the catalytic histidine in FE (the free enzyme), EI (the noncovalent enzyme inhibitor complex), and TC. We demonstrated that the acid/base alteration is controlled by the charged groups in the active site—the catalytic Asp102 carboxylate and the oxyanion. The effect of these groups on the catalytic His is modulated by water solvation of the active site.     

QM/MM through the 1990s: The First Twenty Years of Method Development and Applications

The 2013 Nobel Prize in Chemistry was awarded to the authors of the first two publications utilizing the concept of combined quantum mechanical and molecular mechanical (QM/MM) methods. In celebrating this great event in computational chemistry, we review the early development of combined QM/MM techniques and the associated events that took place through the mid-1990s. We also offer some prospects for the future development of quantum mechanical techniques for macromolecular systems.     

A new family of beta-hairpin mimetics based on a trypsin inhibitor from sunflower seeds

The ability of proteases to regulate many aspects of cell function and defense accounts for the considerable interest in the design of novel protease inhibitors. There are many naturally occurring proteinaceous serine protease inhibitors, one of which is a 14 amino acid cyclic peptide from sunflower seeds that shows both sequence and conformational similarity with the trypsin-reactive loop of the Bowman-Birk family of serine protease inhibitors. This inhibitor adopts a beta-hairpin conformation when bound at the active site of bovine beta-trypsin. We illustrate here an approach to inhibitor design in which the beta hairpin from the naturally occurring peptide is transplanted onto a hairpin-inducing template. Two mimetics with the sequences RC*TKSIPPIC*F (where C*C* is a disulfide) and TKSIPPI are studied, each mounted onto a D-Pro-L-Pro template. NMR studies revealed a well-defined beta-hairpin conformation for each mimetic in aqueous solution; this conformation is closely related to the trypsin-bound conformation of the natural inhibitor and includes a cis-Ile-Pro peptide bond. Both mimetics inhibit trypsin in the mid nanomolar range. An alanine scan revealed the importance for inhibitory activity of the specificity-determining Lys residue and of the first but not the second Pro residue in the IPPI motif. Since these hairpin mimetics can be prepared by parallel combinatorial synthesis, this family of molecules may be a useful starting point for the discovery of other biologically or medicinally useful serine protease inhibitors.