High Level QM/MM Modeling of the Formation of the Tetrahedral Intermediate in the Acylation of Wild Type and K73A Mutant TEM-1 Class A β-Lactamase

Exploring Covalent Docking Mechanisms of Boron-Based Inhibitors to Class A, C and D β-Lactamases Using Time-dependent Hybrid QM/MM Simulations

Recently, molecular covalent docking has been extensively developed to design new classes of inhibitors that form chemical bonds with their biological targets. This strategy for the design of such inhibitors, in particular boron-based inhibitors, holds great promise for the vast family of β-lactamases produced, inter alia, by Gram-negative antibiotic-resistant bacteria. However, the description of covalent docking processes requires a quantum-mechanical approach, and so far, only a few studies of this type have been presented. This study accurately describes the covalent docking process between two model inhibitors - representing two large families of inhibitors based on boronic-acid and bicyclic boronate scaffolds, and three β-lactamases which belong to the A, C, and D classes. Molecular fragments containing boron can be converted from a neutral, trigonal, planar state with sp2 hybridization to the anionic, tetrahedral sp3 state in a process sometimes referred to as morphing. This study applies multi-scale modeling methods, in particular, the hybrid QM/MM approach which has predictive power reaching well beyond conventional molecular modeling. Time-dependent QM/MM simulations indicated several structural changes and geometric preferences, ultimately leading to covalent docking processes. With current computing technologies, this approach is not computationally expensive, can be used in standard molecular modeling and molecular design works, and can effectively support experimental research which should allow for a detailed understanding of complex processes important to molecular medicine. In particular, it can support the rational design of covalent boron-based inhibitors for β-lactamases as well as for many other enzyme systems of clinical relevance, including SARS-CoV-2 proteins.

Unraveling the energetic significance of chemical events in enzyme catalysis via machine-learning based regression approach

The bacterial enzyme class of β-lactamases are involved in benzylpenicillin acylation reactions, which are currently being revisited using hybrid quantum mechanical molecular mechanical (QM/MM) chain-of-states pathway optimizations. Minimum energy pathways are sampled by reoptimizing pathway geometry under different representative protein environments obtained through constrained molecular dynamics simulations. Predictive potential energy surface models in the reaction space are trained with machine-learning regression techniques. Herein, using TEM-1/benzylpenicillin acylation reaction as the model system, we introduce two model-independent criteria for delineating the energetic contributions and correlations in the predicted reaction space. Both methods are demonstrated to effectively quantify the energetic contribution of each chemical process and identify the rate limiting step of enzymatic reaction with high degrees of freedom. The consistency of the current workflow is tested under seven levels of quantum chemistry theory and three non-linear machine-learning regression models. The proposed approaches are validated to provide qualitative compliance with experimental mutagenesis studies.

Acylation and deacylation mechanism and kinetics of penicillin G reaction with Streptomyces R61 DD-peptidase

Two quantum mechanical (QM)-cluster models are built for studying the acylation and deacylation mechanism and kinetics of Streptomyces R61 DD-peptidase with the penicillin G at atomic level detail. DD-peptidases are bacterial enzymes involved in the cross-linking of peptidoglycan to form the cell wall, necessary for bacterial survival. The cross-linking can be inhibited by antibiotic beta-lactam derivatives through acylation, preventing the acyl-enzyme complex from undergoing further deacylation. The deacylation step was predicted to be rate-limiting. Transition state and intermediate structures are found using density functional theory in this study, and thermodynamic and kinetic properties of the proposed mechanism are evaluated. The acyl-enzyme complex is found lying in a deep thermodynamic sink, and deacylation is indeed the severely rate-limiting step, leading to suicide inhibition of the peptidoglycan cross-linking. The usage of QM-cluster models is a promising technique to understand, improve, and design antibiotics to disrupt function of the Streptomyces R61 DD-peptidase.

Benchmarking quantum mechanical methods for calculating reaction energies of reactions catalyzed by enzymes

To assess the accuracy of different quantum mechanical methods for biochemical modeling, the reaction energies of 20 small model reactions (chosen to represent chemical steps catalyzed by commonly studied enzymes) were calculated. The methods tested included several popular Density Functional Theory (DFT) functionals, second-order Møller Plesset perturbation theory (MP2) and its spin-component scaled variant (SCS-MP2), and coupled cluster singles and doubles and perturbative triples (CCSD(T)). Different basis sets were tested. CCSD(T)/aug-cc-pVTZ results for all 20 reactions were used to benchmark the other methods. It was found that MP2 and SCS-MP2 reaction energy calculation results are similar in quality to CCSD(T) (mean absolute error (MAE) of 1.2 and 1.3 kcal mol−1, respectively). MP2 calculations gave a large error in one case, and are more subject to basis set effects, so in general SCS-MP2 calculations are a good choice when CCSD(T) calculations are not feasible. Results with different DFT functionals were of reasonably good quality (MAEs of 2.5–5.1 kcal mol−1), whereas popular semi-empirical methods (AM1, PM3, SCC-DFTB) gave much larger errors (MAEs of 11.6–14.6 kcal mol−1). These results should be useful in guiding methodological choices and assessing the accuracy of QM/MM calculations on enzyme-catalyzed reactions.

An Efficient Computational Assay for β-Lactam Antibiotic Breakdown by Class A β-Lactamases

Class A β-lactamases cause clinically relevant resistance to β-lactam antibiotics. Carbapenem degradation is a particular concern. We present an efficient QM/MM molecular simulation protocol that accurately predicts the activity of β-lactamases against carbapenems. Simulations take less than 24 CPU hours, a greater than 99% reduction, and do not require fitting against experimental data or significant parametrization. This computational assay also reveals mechanistic details of β-lactam breakdown and should assist in evaluating emerging β-lactamase variants and developing new antibiotics.

β-Lactamases and β-Lactamase Inhibitors in the 21st Century

The β-lactams retain a central place in the antibacterial armamentarium. In Gram-negative bacteria, β-lactamase enzymes that hydrolyze the amide bond of the four-membered β-lactam ring are the primary resistance mechanism, with multiple enzymes disseminating on mobile genetic elements across opportunistic pathogens such as Enterobacteriaceae (e.g., Escherichia coli) and non-fermenting organisms (e.g., Pseudomonas aeruginosa). β-Lactamases divide into four classes; the active-site serine β-lactamases (classes A, C and D) and the zinc-dependent or metallo-β-lactamases (MBLs; class B). Here we review recent advances in mechanistic understanding of each class, focusing upon how growing numbers of crystal structures, in particular for β-lactam complexes, and methods such as neutron diffraction and molecular simulations, have improved understanding of the biochemistry of β-lactam breakdown. A second focus is β-lactamase interactions with carbapenems, as carbapenem-resistant bacteria are of grave clinical concern and carbapenem-hydrolyzing enzymes such as KPC (class A) NDM (class B) and OXA-48 (class D) are proliferating worldwide. An overview is provided of the changing landscape of β-lactamase inhibitors, exemplified by the introduction to the clinic of combinations of β-lactams with diazabicyclooctanone and cyclic boronate serine β-lactamase inhibitors, and of progress and strategies toward clinically useful MBL inhibitors. Despite the long history of β-lactamase research, we contend that issues including continuing unresolved questions around mechanism; opportunities afforded by new technologies such as serial femtosecond crystallography; the need for new inhibitors, particularly for MBLs; the likely impact of new β-lactam:inhibitor combinations and the continuing clinical importance of β-lactams mean that this remains a rewarding research area.

The reaction mechanism of Zika virus NS2B/NS3 serine protease inhibition by dipeptidyl aldehyde: a QM/MM study

Zika virus (ZIKV) infection has become a global public health problem, associated with microcephaly in newborns and Guillain-Barré syndrome in adults. Currently, there are no commercially available anti-ZIKV drugs. The viral protease NS2B/NS3, which is involved in viral replication and maturation, is a potential drug target. Peptidomimetic aldehyde inhibitors bind covalently to the catalytic S135 of the NS3 protease. Here, we apply hybrid quantum mechanics/molecular mechanics (QM/MM) free-energy simulations at the PDDG-PM3/ff14SB level to investigate the inhibition mechanism of the ZIKV protease by a dipeptidyl aldehyde inhibitor (acyl-KR-aldehyde). The results show that proton transfer from the catalytic S135 to H51 occurs in concert with nucleophilic addition on the aldehyde warhead by S135. The anionic covalent complex between the dipeptidyl aldehyde and the ZIKV protease is analogous to the tetrahedral intermediate for substrate hydrolysis. Spontaneous protonation by H51 forms the hemiacetal. In addition, we use correlated ab initio QM/MM potential energy path calculations at levels up to LCCSD(T)/(aug)-cc-pVTZ to obtain accurate potential energy profiles of the reaction, which also support a concerted mechanism. These results provide detailed insight into the mechanism of ZIKV protease inhibition by a peptidyl aldehyde inhibitor, which will guide in the design of inhibitors.

Semi-Empirical Born-Oppenheimer Molecular Dynamics (SEBOMD) within the Amber Biomolecular Package

Semi-empirical quantum methods from the neglect of differential diatomic overlap (NDDO) family such as MNDO, AM1, or PM3 are fast albeit approximate quantum methods. By combining them with linear scaling methods like the divide & conquer (D&C) method, it is possible to quickly evaluate the energy of systems containing hundreds to thousands of atoms. We here present our implementation in the Amber biomolecular package of a SEBOMD module that provides a way to run semi-empirical Born-Oppenheimer molecular dynamics. At each step of a SEBOMD, a fully converged self-consistent field (SCF) calculation is performed to obtain the semiempirical quantum potential energy of a molecular system encaged or not in periodic boundary conditions. We describe the implementation and the features of our SEBOMD implementation. We show the requirements to conserve the total energy in NVE simulations, and how to accelerate SCF convergence through density matrix extrapolation. Specific ways of handling periodic boundary conditions using mechanical embedding or electrostatic embedding through a tailored quantum Ewald summation is developed. The parallel performance of SEBOMD simulations using the D&C scheme are presented for liquid water systems of various sizes, and a comparison between the traditional full diagonalization scheme and the D&C approach for the reproduction of the structure of liquid water illustrates the potentiality of SEBOMD to simulate molecular systems containing several hundreds of atoms for hundreds of picoseconds with a quantum mechanical potential in a reasonable amount of CPU time.

On-surface synthesis on a bulk insulator surface

On-surface synthesis has rapidly emerged as a most promising approach to prepare functional molecular structures directly on a support surface. Compared to solution synthesis, performing chemical reactions on a surface offers several exciting new options: due to the absence of a solvent, reactions can be envisioned that are otherwise not feasible due to the insolubility of the reaction product. Perhaps even more important, the confinement to a two-dimensional surface might enable reaction pathways that are not accessible otherwise. Consequently, on-surface synthesis has attracted great attention in the last decade, with an impressive number of classical reactions transferred to a surface as well as new reactions demonstrated that have no classical analogue. So far, the majority of the work has been carried out on conducting surfaces. However, when aiming for electronic decoupling of the resulting structures, e.g. for the use in future molecular electronic devices, non-conducting surfaces are highly desired. Here, we review the current status of on-surface reactions demonstrated on the (10.4) surface of the bulk insulator calcite. Besides thermally induced C-C coupling of halogen-substituted aryls, photochemically induced [2  +  2] cycloaddition has been proven possible on this surface. Moreover, experimental evidence exists for coupling of terminal alkynes as well as diacetylene polymerization. While imaging of the resulting structures with dynamic atomic force microscopy provides a direct means of reaction verification, the detailed reaction pathway often remains unclear. Especially in cases where the presence of metal atoms is known to catalyze the corresponding solution chemistry reaction (e.g. in the case of the Ullmann reaction), disclosing the precise reaction pathway is of importance to understand and generalize on-surface reactivity on a bulk insulator surface. To this end, density-functional theory calculations have proven to provide atomic-scale insights that have greatly contributed to unravelling the details of on-surface synthesis on a bulk insulator surface.

Excess positional mutual information predicts both local and allosteric mutations affecting beta lactamase drug resistance

MOTIVATION

Bacterial resistance to antibiotics, particularly plasmid-encoded resistance to beta lactam drugs, poses an increasing threat to human health. Point mutations to beta-lactamase enzymes can greatly alter the level of resistance conferred, but predicting the effects of such mutations has been challenging due to the large combinatorial space involved and the subtle relationships of distant residues to catalytic function. Therefore we desire an information-theoretic metric to sensitively and robustly detect both local and distant residues that affect substrate conformation and catalytic activity.

RESULTS

Here, we report the use of positional mutual information in multiple microsecond-length molecular dynamics (MD) simulations to predict residues linked to catalytic activity of the CTX-M9 beta lactamase. We find that motions of the bound drug are relatively isolated from motions of the protein as a whole, which we interpret in the context of prior theories of catalysis. In order to robustly identify residues that are weakly coupled to drug motions but nonetheless affect catalysis, we utilize an excess mutual information metric. We predict 31 such residues for the cephalosporin antibiotic cefotaxime. Nine of these have previously been tested experimentally, and all decrease both enzyme rate constants and empirical drug resistance. We prospectively validate our method by testing eight high-scoring mutations and eight low-scoring controls in bacteria. Six of eight predicted mutations decrease cefotaxime resistance greater than 2-fold, while only one control shows such an effect. The ability to prospectively predict new variants affecting bacterial drug resistance is of great interest to clinical and epidemiological surveillance.

AVAILABILITY AND IMPLEMENTATION

Excess mutual information code is available at https://github.com/kassonlab/positionalmi CONTACT: kasson@virginia.edu.

Insight into the Mechanism of Hydrolysis of Meropenem by OXA-23 Serine-β-lactamase Gained by Quantum Mechanics/Molecular Mechanics Calculations

The fast and constant development of drug resistant bacteria represents a serious medical emergency. To overcome this problem, the development of drugs with new structures and modes of action is urgently needed. In this work, we investigated, at the atomistic level, the mechanisms of hydrolysis of Meropenem by OXA-23, a class D β-lactamase, combining unbiased classical molecular dynamics and umbrella sampling simulations with classical force field-based and quantum mechanics/molecular mechanics potentials. Our calculations provide a detailed structural and dynamic picture of the molecular steps leading to the formation of the Meropenem-OXA-23 covalent adduct, the subsequent hydrolysis, and the final release of the inactive antibiotic. In this mechanistic framework, the predicted activation energy is in good agreement with experimental kinetic measurements, validating the expected reaction path.

QM/MM simulations as an assay for carbapenemase activity in class A β-lactamases

Carbapenems, 'last resort' antibiotics for many bacterial infections, can now be broken down by several class A β-lactamases (i.e. carbapenemases). Here, carbapenemase activity is predicted through QM/MM dynamics simulations of acyl-enzyme deacylation, requiring only the 3D structure of the apo-enzyme. This may assist in anticipating resistance and future antibiotic design.

In silico screening of 393 mutants facilitates enzyme engineering of amidase activity in CalB

Our previously presented method for high throughput computational screening of mutant activity (Hediger et al., 2012) is benchmarked against experimentally measured amidase activity for 22 mutants of Candida antarctica lipase B (CalB). Using an appropriate cutoff criterion for the computed barriers, the qualitative activity of 15 out of 22 mutants is correctly predicted. The method identifies four of the six most active mutants with ≥3-fold wild type activity and seven out of the eight least active mutants with ≤0.5-fold wild type activity. The method is further used to screen all sterically possible (386) double-, triple- and quadruple-mutants constructed from the most active single mutants. Based on the benchmark test at least 20 new promising mutants are identified.

Combined quantum mechanics/molecular mechanics (QM/MM) methods in computational enzymology

Computational enzymology is a rapidly maturing field that is increasingly integral to understanding mechanisms of enzyme-catalyzed reactions and their practical applications. Combined quantum mechanics/molecular mechanics (QM/MM) methods are important in this field. By treating the reacting species with a quantum mechanical method (i.e., a method that calculates the electronic structure of the active site) and including the enzyme environment with simpler molecular mechanical methods, enzyme reactions can be modeled. Here, we review QM/MM methods and their application to enzyme-catalyzed reactions to investigate fundamental and practical problems in enzymology. A range of QM/MM methods is available, from cheaper and more approximate methods, which can be used for molecular dynamics simulations, to highly accurate electronic structure methods. We discuss how modeling of reactions using such methods can provide detailed insight into enzyme mechanisms and illustrate this by reviewing some recent applications. We outline some practical considerations for such simulations. Further, we highlight applications that show how QM/MM methods can contribute to the practical development and application of enzymology, e.g., in the interpretation and prediction of the effects of mutagenesis and in drug and catalyst design.

QM/MM modelling of ketosteroid isomerase reactivity indicates that active site closure is integral to catalysis

Ketosteroid isomerase (Δ⁵-3-keto steroid isomerase or steroid Δ-isomerase) is a highly efficient enzyme at the centre of current debates on enzyme catalysis. We have modelled the reaction mechanism of the isomerization of 3-oxo-Δ⁵-steroids into their Δ⁴-conjugated isomers using high-level combined quantum mechanics/molecular mechanics (QM/MM) methods, and semi-empirical QM/MM molecular dynamics simulations. Energy profiles were obtained at various levels of QM theory (AM1, B3LYP and SCS-MP2). The high-level QM/MM profile is consistent with experimental data. QM/MM dynamics simulations indicate that active site closure and desolvation of the catalytic Asp38 occur before or during formation of dienolate intermediates. These changes have a significant effect on the reaction barrier. A low barrier to reaction is found only when the active site is closed, poising it for catalysis. This conformational change is thus integral to the whole process. The effects on the barrier are apparently largely due to changes in solvation. The combination of high-level QM/MM energy profiles and QM/MM dynamics simulation shows that the reaction involves active site closure, desolvation of the catalytic base, efficient isomerization and re-opening of the active site. These changes highlight the transition between the ligand binding/releasing form and the catalytic form of the enzyme. The results demonstrate that electrostatic interactions (as a consequence of pre-organization of the active site) are crucial for stabilization during the chemical reaction step, but closure of the active site is essential for efficient catalysis to occur.

A computational methodology to screen activities of enzyme variants

We present a fast computational method to efficiently screen enzyme activity. In the presented method, the effect of mutations on the barrier height of an enzyme-catalysed reaction can be computed within 24 hours on roughly 10 processors. The methodology is based on the PM6 and MOZYME methods as implemented in MOPAC2009, and is tested on the first step of the amide hydrolysis reaction catalyzed by the Candida Antarctica lipase B (CalB) enzyme. The barrier heights are estimated using adiabatic mapping and shown to give barrier heights to within 3 kcal/mol of B3LYP/6-31G(d)//RHF/3-21G results for a small model system. Relatively strict convergence criteria (0.5 kcal/(molÅ)), long NDDO cutoff distances within the MOZYME method (15 Å) and single point evaluations using conventional PM6 are needed for reliable results. The generation of mutant structures and subsequent setup of the semiempirical calculations are automated so that the effect on barrier heights can be estimated for hundreds of mutants in a matter of weeks using high performance computing.

Molecular dynamics of class A β-lactamases-effects of substrate binding

The effects of substrate binding on class A β-lactamase dynamics were studied using molecular dynamics simulations of two model enzymes; 40 100-ns trajectories of the free and substrate-bound forms of TEM-1 (with benzylpenicillin) and PSE-4 (with carbenicillin) were recorded (totaling 4.0 μs). Substrates were parameterized with the CHARMM General Force Field. In both enzymes, the Ω loop exhibits a marked flexibility increase upon substrate binding, supporting the hypothesis of substrate gating. However, specific interactions that are formed or broken in the Ω loop upon binding differ between the two enzymes: dynamics are conserved, but not specific interactions. Substrate binding also has a global structuring effect on TEM-1, but not on PSE-4. Changes in TEM-1's normal modes show long-range effects of substrate binding on enzyme dynamics. Hydrogen bonds observed in the active site are mostly preserved upon substrate binding, and new, transient interactions are also formed. Agreement between NMR relaxation parameters and our theoretical results highlights the dynamic duality of class A β-lactamases: enzymes that are highly structured on the ps-ns timescale, with important flexibility on the μs-ms timescale in regions such as the Ω loop.

Crystal structure of a preacylation complex of the β-lactamase inhibitor sulbactam bound to a sulfenamide bond-containing thiol-β-lactamase

The rise of inhibitor-resistant and other β-lactamase variants is generating an interest in developing new β-lactamase inhibitors to complement currently available antibiotics. To gain insight into the chemistry of inhibitor recognition, we determined the crystal structure of the inhibitor preacylation complex of sulbactam, a clinical β-lactamase inhibitor, bound in the active site of the S70C variant of SHV-1 β-lactamase, a resistance enzyme that is normally present in Klebsiella pneumoniae. The S70C mutation was designed to affect the reactivity of that catalytic residue to allow for capture of the preacylation complex. Unexpectedly, the 1.45 Å resolution inhibitor complex structure revealed that residue C70 is involved in a sulfenamide bond with K73. Such a covalent bond is not present in the wild-type SHV-1 or in an apo S70C structure also determined in this study. This bond likely contributed significantly to obtaining the preacylation complex with sulbactam due to further decreased reactivity toward substrates. The intact sulbactam is positioned in the active site such that its carboxyl moiety interacts with R244, S130, and T235 and its carbonyl moiety is situated in the oxyanion hole. To our knowledge, in addition to being the first preacylation inhibitor β-lactamase complex, this is also the first observation of a sulfenamide bond between a cysteine and lysine in an active site. Not only could our results aid, therefore, structure-based inhibitor design efforts in class A β-lactamases, but the sulfenamide-bond forming approach to yield preacylation complexes could also be applied to other classes of β-lactamases and penicillin-binding proteins with the SXXK motif.

Enumerating pathways of proton abstraction based on a spatial and electrostatic analysis of residues in the catalytic site

The pathways of proton abstraction (PA), a key aspect of most catalytic reactions, is often controversial and highly debated. Ultrahigh-resolution diffraction studies, molecular dynamics, quantum mechanics and molecular mechanic simulations are often adopted to gain insights in the PA mechanisms in enzymes. These methods require expertise and effort to setup and can be computationally intensive. We present a push button methodology--Proton abstraction Simulation (PRISM)--to enumerate the possible pathways of PA in a protein with known 3D structure based on the spatial and electrostatic properties of residues in the proximity of a given nucleophilic residue. Proton movements are evaluated in the vicinity of this nucleophilic residue based on distances, potential differences, spatial channels and characteristics of the individual residues (polarity, acidic, basic, etc). Modulating these parameters eliminates their empirical nature and also might reveal pathways that originate from conformational changes. We have validated our method using serine proteases and concurred with the dichotomy in PA in Class A β-lactamases, both of which are hydrolases. The PA mechanism in a transferase has also been corroborated. The source code is made available at www.sanchak.com/prism.

Testing high-level QM/MM methods for modeling enzyme reactions: acetyl-CoA deprotonation in citrate synthase

Combined quantum mechanics/molecular mechanics (QM/MM) calculations with high levels of correlated ab initio theory can now provide benchmarks for enzyme-catalyzed reactions. Here, we use such methods to test various QM/MM methods and the sensitivity of the results to details of the models for an important enzyme reaction, proton abstraction from acetyl-coenzyme A in citrate synthase. We calculate multiple QM/MM potential energy surfaces up to the local coupled cluster theory (LCCSD(T0)) level, with structures optimized at hybrid density functional theory and Hartree-Fock levels. The influence of QM methods, basis sets, and QM region size is shown to be significant. Correlated ab initio QM/MM calculations give barriers in agreement with experiment for formation of the acetyl-CoA enolate intermediate. In contrast, B3LYP fails to identify the enolate as an intermediate, whereas BH&HLYP does. The results indicate that QM/MM methods and setup should be tested, ideally using high-level calculations, to draw reliable mechanistic conclusions.