Implementation of the SCC-DFTB method for hybrid QM/MM simulations within the amber molecular dynamics package

Theoretical Insights into the Catalytic Oxidation of Phenols and Arylamines by Laccases via the Proton-Coupled Electron Transfer Mechanism

Laccases play a vital role in the degradation of toxic phenolic and aromatic amine compounds, generating considerable attention in ecological pollution remediation. However, the distinct mechanism of the laccase-catalyzed oxidation of phenols and arylamines remains unclear. Here, we examined the catalytic oxidation mechanisms of phenols and arylamines by Trametes versicolor (TvL) and Melanocarpus albomyces (MaL) laccases using molecular docking, quantum mechanics (QM), and QM/molecular mechanics (QM/MM) calculations. We docked four phenolic substrates, including 1,2-benzenediol, 2-propenylphenol, 2-methoxyhydroquinone, and 2-aminophenol, to TvL and identified their favorable reaction conformations, in which Asp206 of TvL plays an important role in binding substrates to promote the catalytic reactions. Based on the docking conformations, the QM and QM/MM calculations revealed that the oxidation reactions take place via a proton-coupled electron transfer mechanism, with proton transfer (PT) from the hydroxyl groups of substrates to the side chain of Asp206 and synchronous electron hopping from the aromatic ring of substrates to the type one copper (T1Cu) of TvL. For the MaL and 2,6-dimethoxyphenol interacting system, the oxidation reactions occur through a concerted double-proton-coupled electron transfer mechanism with a water-mediated indirect PT from the hydroxyl group of substrates to the conserved Glu235 and electron hopping from the substrate to T1Cu at the same time. The corresponding energy barriers change from 0.7 to 18.4 kcal/mol, indicating the different degradation rates of the phenols and arylamines by laccases. These findings provide insights into the oxidation mechanism of phenols and arylamines by laccases and may extend the applications of laccases.

PyDFT-QMMM: A modular, extensible software framework for DFT-based QM/MM molecular dynamics

PyDFT-QMMM is a Python-based package for performing hybrid quantum mechanics/molecular mechanics (QM/MM) simulations at the density functional level of theory. The program is designed to treat short-range and long-range interactions through user-specified combinations of electrostatic and mechanical embedding procedures within periodic simulation domains, providing necessary interfaces to external quantum chemistry and molecular dynamics software. To enable direct embedding of long-range electrostatics in periodic systems, we have derived and implemented force terms for our previously described QM/MM/PME approach [Pederson and McDaniel, J. Chem. Phys. 156, 174105 (2022)]. Communication with external software packages Psi4 and OpenMM is facilitated through Python application programming interfaces (APIs). The core library contains basic utilities for running QM/MM molecular dynamics simulations, and plug-in entry-points are provided for users to implement custom energy/force calculation and integration routines, within an extensible architecture. The user interacts with PyDFT-QMMM primarily through its Python API, allowing for complex workflow development with Python scripting, for example, interfacing with PLUMED for free energy simulations. We provide benchmarks of forces and energy conservation for the QM/MM/PME and alternative QM/MM electrostatic embedding approaches. We further demonstrate a simple example use case for water solute in a water solvent system, for which radial distribution functions are computed from 100 ps QM/MM simulations; in this example, we highlight how the solvation structure is sensitive to different basis-set choices due to under- or over-polarization of the QM water molecule's electron density.

Mechanistic insights into the conversion of flavin adenine dinucleotide (FAD) to 8-formyl FAD in formate oxidase: a combined experimental and in-silico study

Formate oxidase (FOx), which contains 8-formyl flavin adenine dinucleotide (FAD), exhibits a distinct advantage in utilizing ambient oxygen molecules for the oxidation of formic acid compared to other glucose-methanol-choline (GMC) oxidoreductase enzymes that contain only the standard FAD cofactor. The FOx-mediated conversion of FAD to 8-formyl FAD results in an approximate 10-fold increase in formate oxidase activity. However, the mechanistic details underlying the autocatalytic formation of 8-formyl FAD are still not well understood, which impedes further utilization of FOx. In this study, we employ molecular dynamics simulation, QM/MM umbrella sampling simulation, enzyme activity assay, site-directed mutagenesis, and spectroscopic analysis to elucidate the oxidation mechanism of FAD to 8-formyl FAD. Our results reveal that a catalytic water molecule, rather than any catalytic amino acids, serves as a general base to deprotonate the C8 methyl group on FAD, thus facilitating the formation of a quinone-methide tautomer intermediate. An oxygen molecule subsequently oxidizes this intermediate, resulting in a C8 methyl hydroperoxide anion that is protonated and dissociated to form OHC-RP and OH-. During the oxidation of FAD to 8-formyl FAD, the energy barrier for the rate-limiting step is calculated to be 22.8 kcal/mol, which corresponds to the required 14-hour transformation time observed experimentally. Further, the elucidated oxidation mechanism reveals that the autocatalytic formation of 8-formyl FAD depends on the proximal arginine and serine residues, R87 and S94, respectively. Enzymatic activity assay validates that the mutation of R87 to lysine reduces the kcat value to 75% of the wild-type, while the mutation to histidine results in a complete loss of activity. Similarly, the mutant S94I also leads to the deactivation of enzyme. This dependency arises because the nucleophilic OH- group and the quinone-methide tautomer intermediate are stabilized through the noncovalent interaction provided by R87 and S94. These findings not only explain the mechanistic details of each reaction step but also clarify the functional role of R87 and S94 during the oxidative maturation of 8-formyl FAD, thereby providing crucial theoretical support for the development of novel flavoenzymes with enhanced redox properties.

Condensed Matter Systems Exposed to Radiation: Multiscale Theory, Simulations, and Experiment

This roadmap reviews the new, highly interdisciplinary research field studying the behavior of condensed matter systems exposed to radiation. The Review highlights several recent advances in the field and provides a roadmap for the development of the field over the next decade. Condensed matter systems exposed to radiation can be inorganic, organic, or biological, finite or infinite, composed of different molecular species or materials, exist in different phases, and operate under different thermodynamic conditions. Many of the key phenomena related to the behavior of irradiated systems are very similar and can be understood based on the same fundamental theoretical principles and computational approaches. The multiscale nature of such phenomena requires the quantitative description of the radiation-induced effects occurring at different spatial and temporal scales, ranging from the atomic to the macroscopic, and the interlinks between such descriptions. The multiscale nature of the effects and the similarity of their manifestation in systems of different origins necessarily bring together different disciplines, such as physics, chemistry, biology, materials science, nanoscience, and biomedical research, demonstrating the numerous interlinks and commonalities between them. This research field is highly relevant to many novel and emerging technologies and medical applications.

CHARMM-GUI QM/MM Interfacer for a Quantum Mechanical and Molecular Mechanical (QM/MM) Simulation Setup: 1. Semiempirical Methods

Quantum mechanical (QM) treatments, when combined with molecular mechanical (MM) force fields, can effectively handle enzyme-catalyzed reactions without significantly increasing the computational cost. In this context, we present CHARMM-GUI QM/MM Interfacer, a web-based cyberinfrastructure designed to streamline the preparation of various QM/MM simulation inputs with ligand modification. The development of QM/MM Interfacer has been achieved through integration with existing CHARMM-GUI modules, such as PDB Reader and Manipulator, Solution Builder, and Membrane Builder. In addition, new functionalities have been developed to facilitate the one-stop preparation of QM/MM systems and enable interactive and intuitive ligand modifications and QM atom selections. QM/MM Interfacer offers support for a range of semiempirical QM methods, including AM1(+/d), PM3(+/PDDG), MNDO(+/d, +/PDDG), PM6, RM1, and SCC-DFTB, tailored for both AMBER and CHARMM. A nontrivial setup related to ligand modification, link-atom insertion, and charge distribution is automatized through intuitive user interfaces. To illustrate the robustness of QM/MM Interfacer, we conducted QM/MM simulations of three enzyme-substrate systems: dihydrofolate reductase, insulin receptor kinase, and oligosaccharyltransferase. In addition, we have created three tutorial videos about building these systems, which can be found at https://www.charmm-gui.org/demo/qmi. QM/MM Interfacer is expected to be a valuable and accessible web-based tool that simplifies and accelerates the setup process for hybrid QM/MM simulations.

The development of the QM/MM interface and its application for the on-the-fly QM/MM nonadiabatic dynamics in JADE package: Theory, implementation, and applications

Understanding the nonadiabatic dynamics of complex systems is a challenging task in computational photochemistry. Herein, we present an efficient and user-friendly quantum mechanics/molecular mechanics (QM/MM) interface to run on-the-fly nonadiabatic dynamics. Currently, this interface consists of an independent set of codes designed for general-purpose use. Herein, we demonstrate the ability and feasibility of the QM/MM interface by integrating it with our long-term developed JADE package. Tailored to handle nonadiabatic processes in various complex systems, especially condensed phases and protein environments, we delve into the theories, implementations, and applications of on-the-fly QM/MM nonadiabatic dynamics. The QM/MM approach is established within the framework of the additive QM/MM scheme, employing electrostatic embedding, link-atom inclusion, and charge-redistribution schemes to treat the QM/MM boundary. Trajectory surface-hopping dynamics are facilitated using the fewest switches algorithm, encompassing classical and quantum treatments for nuclear and electronic motions, respectively. Finally, we report simulations of nonadiabatic dynamics for two typical systems: azomethane in water and the retinal chromophore PSB3 in a protein environment. Our results not only illustrate the power of the QM/MM program but also reveal the important roles of environmental factors in nonadiabatic processes.

Generating High-Precision Force Fields for Molecular Dynamics Simulations to Study Chemical Reaction Mechanisms Using Molecular Configuration Transformer

Theoretical studies on chemical reaction mechanisms have been crucial in organic chemistry. Traditionally, calculating the manually constructed molecular conformations of transition states for chemical reactions using quantum chemical calculations is the most commonly used method. However, this way is heavily dependent on individual experience and chemical intuition. In our previous study, we proposed a research paradigm that used enhanced sampling in molecular dynamics simulations to study chemical reactions. This approach can directly simulate the entire process of a chemical reaction. However, the computational speed limited the use of high-precision potential energy functions for simulations. To address this issue, we presented a scheme for training high-precision force fields for molecular modeling using a previously developed graph-neural-network-based molecular model, molecular configuration transformer. This potential energy function allowed for highly accurate simulations at a low computational cost, leading to more precise calculations of the mechanism of chemical reactions. We applied this approach to study a Claisen rearrangement reaction and a carbonyl insertion reaction catalyzed by manganese.

Analytic Gradients for the Electrostatic Embedding QM/MM Model in Periodic Boundary Conditions Using Particle-Mesh Ewald Sums and Electrostatic Potential Fitted Charge Operators

Long-range electrostatic effects are fundamental for describing chemical reactivity in the condensed phase. Here, we present the methodology of an efficient quantum mechanical/molecular mechanical (QM/MM) model in periodic boundary conditions (PBC) compatible with QM/MM boundaries at chemical bonds. The method combines electrostatic potential fitted charge operators and electrostatic potentials derived from the smooth particle-mesh Ewald (PME) sum approach. The total energy and its analytic first derivatives with respect to QM, MM, and lattice vectors allow QM/MM molecular dynamics (MD) in the most common thermodynamic ensembles. We demonstrate the robustness of the method by performing a QM/MM MD equilibration of methanol in water. We simulate the cis/trans isomerization free-energy profiles in water of proline amino acid and a proline-containing oligopeptide, showing a correct description of the reaction barrier. Our PBC-compatible QM/MM model can efficiently be used to study the chemical reactivity in the condensed phase and enzymatic catalysis.

Three-metal ion mechanism of cross-linked and uncross-linked DNA polymerase β: A theoretical study

In our recent publication, we have proposed a revised base excision repair pathway in which DNA polymerase β (Polβ) catalyzes Schiff base formation prior to the gap-filling DNA synthesis followed by β-elimination. In addition, the polymerase activity of Polβ employs the "three-metal ion mechanism" instead of the long-standing "two-metal ion mechanism" to catalyze phosphodiester bond formation based on the fact derived from time-resolved x-ray crystallography that a third Mg2+ was captured in the polymerase active site after the chemical reaction was initiated. In this study, we develop the models of the uncross-linked and cross-linked Polβ complexes and investigate the "three-metal ion mechanism" vs the "two-metal ion mechanism" by using the quantum mechanics/molecular mechanics molecular dynamics simulations. Our results suggest that the presence of the third Mg2+ ion stabilizes the reaction-state structures, strengthens correct nucleotide binding, and accelerates phosphodiester bond formation. The improved understanding of Polβ's catalytic mechanism provides valuable insights into DNA replication and damage repair.

Unraveling the anti-corrosion mechanisms of a novel hydrazone derivative on steel in contaminated concrete pore solutions: An integrated study

INTRODUCTION

Corrosion-induced deterioration of infrastructure is a growing global concern. The development and application of corrosion inhibitors are one of the most effective approaches to protect steel rebar from corrosion. Hence, this study focuses on a novel hydrazone derivative, (E)-N'-(4-(dimethylamino)benzylidene)-2-(5-methoxy-2-methyl-1H-indol-3-yl)aceto-hydrazide (HIND), and its potential application to mitigate corrosion in steel rebar exposed to chloride-contaminated concrete pore solutions (ClSCPS).

OBJECTIVES

The research aims to evaluate the anti-corrosion capabilities of HIND on steel rebar within a simulated corrosive environment, focusing on the mechanisms of its inhibitory effect.

METHODS

The corrosion of steel rebar exposed to the ClSCPS was studied through weight loss and electrochemical methods. The surface morphology of steel rebar surface was characterized by FE-SEM-EDS, AFM; oxidation states of the steel rebar and crystal structures were examined using X-ray photoelectron spectroscopy (XPS) and X-ray diffraction (XRD) methods. Further, experimental findings were complemented by theoretical studies using self-consistent-charge density-functional tight-binding (SCC-DFTB) simulations. The performance of HIND was monitored at an optimal concentration over a period of 30 days.

RESULTS

The results indicated a significant reduction in steel rebar corrosion upon introducing HIND. The inhibitor molecules adhered to the steel surface, preventing further deterioration and achieving an inhibition efficiency of 88.4% at 0.5 mmol/L concentration. The surface morphology analysis confirmed the positive effect of HIND on the rebar surface, showing a decrease in the surface roughness of the steel rebar from 183.5 in uninhibited to 50 nm in inhibited solutions. Furthermore, SCC-DFTB simulations revealed the presence of coordination between iron atoms and HIND active sites.

CONCLUSION

The findings demonstrate the potential of HIND as an effective anti-corrosion agent in chloride-contaminated environments. Its primary adsorption mechanism involves charge transfer from the inhibitor molecules to iron atoms. Therefore, applying HIND could be an effective strategy to address corrosion-related challenges in reinforced infrastructure.

The Effect of 8,5'-Cyclo 2'-deoxyadenosine on the Activity of 10-23 DNAzyme: Experimental and Theoretical Study

The in vivo effectiveness of DNAzymes 10-23 (Dz10-23) is limited due to the low concentration of divalent cations. Modifications of the catalytic loop are being sought to increase the activity of Dz10-23 in physiological conditions. We investigated the effect of 5'S or 5'R 5',8-cyclo-2'deoxyadenosine (cdA) on the activity of Dz10-23. The activity of Dz10-23 was measured in a cleavage assay using radiolabeled RNA. The Density Functional Tight Binding methodology with the self-consistent redistribution of Mulliken charge modification was used to explain different activities of DNAzymes. The substitution of 2'-deoxyadenosine with cdA in the catalytic loop decreased the activity of DNAzymes. Inhibition was dependent on the position of cdA and its absolute configuration. The order of activity of DNAzymes was as follows: wt-Dz > ScdA5-Dz ≈ RcdA15-Dz ≈ ScdA15-Dz > RcdA5-Dz. Theoretical studies revealed that the distance between phosphate groups at position 5 in RcdA5-Dz was significantly increased compared to wt-Dz, while the distance between O4 of dT4 and nonbonding oxygen of PO2 attached to 3'O of dG2 was much shorter. The strong inhibitory effect of RcdA5 may result from hampering the flexibility of the catalytic loop (increased rigidity), which is required for the proper positioning of Me2+ and optimal activity.

Addressing electronic and dynamical evolution of molecules and molecular clusters: DFTB simulations of energy relaxation in polycyclic aromatic hydrocarbons

We present a review of the capabilities of the density functional based Tight Binding (DFTB) scheme to address the electronic relaxation and dynamical evolution of molecules and molecular clusters following energy deposition via either collision or photoabsorption. The basics and extensions of DFTB for addressing these systems and in particular their electronic states and their dynamical evolution are reviewed. Applications to PAH molecules and clusters, carbonaceous systems of major interest in astrochemical/astrophysical context, are reported. A variety of processes are examined and discussed such as collisional hydrogenation, fast collisional processes and induced electronic and charge dynamics, collision-induced fragmentation, photo-induced fragmentation, relaxation in high electronic states, electronic-to-vibrational energy conversion and statistical versus non-statistical fragmentation. This review illustrates how simulations may help to unravel different relaxation mechanisms depending on various factors such as the system size, specific electronic structure or excitation conditions, in close connection with experiments.

Computational Exploration of Phenolic Compounds in Corrosion Inhibition: A Case Study of Hydroxytyrosol and Tyrosol

The corrosion of materials remains a critical challenge with significant economic and infrastructural impacts. A comprehensive understanding of adsorption characteristics of phytochemicals can facilitate the effective design of high-performing environmentally friendly inhibitors. This study conducted a computational exploration of hydroxytyrosol (HTR) and tyrosol (TRS) (potent phenolic compounds found in olive leaf extracts), focusing on their adsorption and reactivity on iron surfaces. Utilizing self-consistent-charge density-functional tight-binding (SCC-DFTB) simulations, molecular dynamics (MD) simulations, and quantum chemical calculations (QCCs), we investigated the molecules' structural and electronic attributes and interactions with iron surfaces. The SCC-DFTB results highlighted that HTR and TRS coordinated with iron atoms when adsorbed individually, but only HTR maintained bonding when adsorbed alongside TRS. At their individual adsorption, HTR and TRS had interaction energies of -1.874 and -1.598 eV, which became more negative when put together (-1.976 eV). The MD simulations revealed parallel adsorption under aqueous and vacuum conditions, with HTR demonstrating higher adsorption energy. The analysis of quantum chemical parameters, including global and local reactivity descriptors, offered crucial insights into molecular reactivity, stability, and interaction-prone atomic sites. QCCs revealed that the fraction of transferred electron ∆N aligned with SCC-DFTB results, while other parameters of purely isolated molecules failed to predict the same. These findings pave the way for potential advancements in anticorrosion strategies leveraging phenolic compounds.

QM/MM Modeling Aided Enzyme Engineering in Natural Products Biosynthesis

Natural products and their derivatives are widely used across various industries, particularly pharmaceuticals. Modern engineered biosynthesis provides an alternative way of producing and meeting the growing need for diverse natural products. Natural enzymes, on the other hand, often exhibit unsatisfactory catalytic characteristics and necessitate further enzyme engineering modifications. QM/MM, as a powerful and extensively used computational tool in the field of enzyme catalysis, has been increasingly applied in rational enzyme engineering over the past decade. In this review, we summarize recent advances in QM/MM computational investigation on enzyme catalysis and enzyme engineering for natural product biosynthesis. The challenges and perspectives for future QM/MM applications aided enzyme engineering in natural product biosynthesis will also be discussed.

Unraveling the catalytic mechanism of SARS-CoV-2 papain-like protease with allosteric modulation of C270 mutation using multiscale computational approaches

Papain-like protease (PL^pro) is a promising therapeutic target against SARS-CoV-2, but its restricted S1/S2 subsites pose an obstacle in developing active site-directed inhibitors. We have recently identified C270 as a novel covalent allosteric site for SARS-CoV-2 PL^pro inhibitors. Here we present a theoretical investigation of the proteolysis reaction catalyzed by the wild-type SARS-CoV-2 PL^pro as well as the C270R mutant. Enhanced sampling MD simulations were first performed to explore the influence of C270R mutation on the protease dynamics, and sampled thermodynamically favorable conformations were then submitted to MM/PBSA and QM/MM MD simulations for thorough characterization of the protease-substrate binding and covalent reactions. The disclosed proteolysis mechanism of PL^pro, as characterized by the occurrence of proton transfer from the catalytic C111 to H272 prior to the substrate binding and with deacylation being the rate-determining step of the whole proteolysis process, is not completely identical to that of the 3C-like protease, another key cysteine protease of coronaviruses. The C270R mutation alters the structural dynamics of the BL2 loop that indirectly impairs the catalytic function of H272 and reduces the binding of the substrate with the protease, ultimately showing an inhibitory effect on PL^pro. Together, these results provide a comprehensive understanding at the atomic level of the key aspects of SARS-CoV-2 PL^pro proteolysis, including the catalytic activity allosterically regulated by C270 modification, which is crucial to the follow-up inhibitor design and development.

Unraveling Bonding Mechanisms and Electronic Structure of Pyridine Oximes on Fe(110) Surface: Deeper Insights from DFT, Molecular Dynamics and SCC-DFT Tight Binding Simulations

The development of corrosion inhibitors with outstanding performance is a never-ending and complex process engaged in by researchers, engineers and practitioners. The computational assessment of organic corrosion inhibitors' performance is a crucial step towards the design of new task-specific materials. Herein, the electronic features, adsorption characteristics and bonding mechanisms of two pyridine oximes, namely 2-pyridylaldoxime (2POH) and 3-pyridylaldoxime (3POH), with the iron surface were investigated using molecular dynamics (MD), and self-consistent-charge density-functional tight-binding (SCC-DFTB) simulations. SCC-DFTB simulations revealed that the 3POH molecule can form covalent bonds with iron atoms in its neutral and protonated states, while the 2POH molecule can only bond with iron through its protonated form, resulting in interaction energies of -2.534, -2.007, -1.897, and -0.007 eV for 3POH, 3POH⁺, 2POH⁺, and 2POH, respectively. Projected density of states (PDOSs) analysis of pyridines-Fe(110) interactions indicated that pyridine molecules were chemically adsorbed on the iron surface. Quantum chemical calculations (QCCs) revealed that the energy gap and Hard and Soft Acids and Bases (HSAB) principles were efficient in predicting the bonding trend of the molecules investigated with an iron surface. 3POH had the lowest energy gap of 1.706 eV, followed by 3POH⁺ (2.806 eV), 2POH⁺ (3.121 eV), and 2POH (3.431 eV). In the presence of a simulated solution, MD simulation showed that the neutral and protonated forms of molecules exhibited a parallel adsorption mode on an iron surface. The excellent adsorption properties and corrosion inhibition performance of 3POH may be attributed to its low stability compared to 2POH molecules.

Mechanism of Chiral-Selective Aminoacylation of an RNA Minihelix Explored by QM/MM Free-Energy Simulations

Aminoacylation of a primordial RNA minihelix composed of D-ribose shows L-amino acid preference over D-amino acid without any ribozymes or enzymes. This preference in the amino acylation reaction likely plays an important role in the establishment of homochirality in L-amino acid in modern proteins. However, molecular mechanisms of the chiral selective reaction remain unsolved mainly because of difficulty in direct observation of the reaction at the molecular scale by experiments. For seeking a possible mechanism of the chiral selectivity, quantum mechanics/molecular mechanics (QM/MM) umbrella sampling molecular dynamics (MD) simulations of the aminoacylation reactions in a modeled RNA were performed to investigate differences in their free-energy profiles along the reactions for L- and D-alanine and its physicochemical origin. The reaction is initiated by approaching a 3′-oxygen of the RNA minihelix to the carbonyl carbon of an aminoacyl phosphate oligonucleotide. The QM/MM umbrella sampling MD calculations showed that the height of the free-energy barrier for L-alanine aminoacylation reaction was 17 kcal/mol, which was 9 kcal/mol lower than that for the D-alanine system. At the transition state, the distance between the negatively charged 3′-oxygen and the positively charged amino group of L-alanine was shorter than that of D-alanine, which was caused by the chirality difference of the amino acid. These results indicate that the transition state for L-alanine is more electrostatically stabilized than that for D-alanine, which would be a plausible mechanism previously unexplained for chiral selectivity in the RNA minihelix aminoacylation.

On the recognition of Yersinia protein tyrosine phosphatase by carboxylic acid derivatives

Some members of Yersinia (Y), a genus of bacteria in the family Yersiniaceae, are pathogenic in humans, causing a range of health problems, from gastrointestinal syndromes to the plague. The Y protein tyrosine phosphatase (PTP) YopH is a crucial virulence determinant, considering the vital roles of PTPs in the intracellular signal transduction pathways and cell cycle control. The structural understanding of YopH as a cellular target in pathogenic conditions caused by Y infection is a prerequisite for designing potent and selective YopH inhibitors. Thus, by using molecular docking simulations, the open and closed conformations of the so-called 'WPD loop' (352-Gly-Asn-Trp-Pro-Asp-Gln-Thr-Ala-Val-Ser-361), located nearby the active site (403-Cys-Arg-Ala-Gly-Val-Gly-Arg-Thr-410) in YopH structure, are shown to be relevant for recognition by carboxylic acid derivatives, and the closed conformation is a more preferable receptor in terms of the quantitative correlation with experimental data. In both cases, aurintricarboxylic acid (ATA) has the greatest affinity to YopH. Consequently, a quantum mechanics/molecular mechanics (QM/MM) molecular model is derived to see into the extent of the ATA-induced open-closed conformational change. Active site residues and the WPD loop, as well as ATA are treated using SCC-DFTB-D (QM level), while the rest of the complex is treated using AMBER force field (MM level). The active/inactive functional behavior of YopH is explored by observing the interaction mode of ATA with the wild-type (wt)/Cys403Ser receptor and evaluating the competitive inhibition parameters. Implications of the present study for experimental research are discussed. Communicated by Ramaswamy H. Sarma.

QDπ: A Quantum Deep Potential Interaction Model for Drug Discovery

We report QDπ-v1.0 for modeling the internal energy of drug molecules containing H, C, N, and O atoms. The QDπ model is in the form of a quantum mechanical/machine learning potential correction (QM/Δ-MLP) that uses a fast third-order self-consistent density-functional tight-binding (DFTB3/3OB) model that is corrected to a quantitatively high-level of accuracy through a deep-learning potential (DeepPot-SE). The model has the advantage that it is able to properly treat electrostatic interactions and handle changes in charge/protonation states. The model is trained against reference data computed at the ωB97X/6-31G* level (as in the ANI-1x data set) and compared to several other approximate semiempirical and machine learning potentials (ANI-1x, ANI-2x, DFTB3, MNDO/d, AM1, PM6, GFN1-xTB, and GFN2-xTB). The QDπ model is demonstrated to be accurate for a wide range of intra- and intermolecular interactions (despite its intended use as an internal energy model) and has shown to perform exceptionally well for relative protonation/deprotonation energies and tautomers. An example application to model reactions involved in RNA strand cleavage catalyzed by protein and nucleic acid enzymes illustrates QDπ has average errors less than 0.5 kcal/mol, whereas the other models compared have errors over an order of magnitude greater. Taken together, this makes QDπ highly attractive as a potential force field model for drug discovery.

Structural Insights into the Antiparallel G-Quadruplex in the Presence of K+ and Mg2+ Ions

G-Quadruplex (GQ) is a secondary structural unit of DNA, formed at the telomere region of the chromosome with a high guanine content. It is reported that the GQs can hinder many biological processes. Thus, research thrives to explore the structural stability of GQs. Studies based on circular dichroism (CD) and nuclear magnetic resonance (NMR) experiments established the vital role of cations such as K⁺ and Mg²⁺ in the stability of antiparallel G-quadruplexes (AGQs). However, there is a need to understand how stability in AGQ is attained in the presence of cations. Here, we employed molecular dynamics (MD) simulations, steered MD (SMD) simulations, and QM/MM calculations to understand the biophysical and electronic bases of the stability imparted to AGQs via cation binding. Our results showed that Mg²⁺ prefers to bind in the plane with the guanine tetrad, whereas K⁺ binds in between the AGQ tetrads. Thus, three Mg²⁺ cations or two K⁺ ions are needed to stabilize an AGQ molecule, where each and every tetrad binds to Mg²⁺ more robustly with a higher binding affinity. SMD revealed that the traversal of K⁺ through the AGQ central channel required less force than that of Mg²⁺, illustrating the presence of more strong interactions between Mg²⁺ and AGQ tetrads compared to K⁺. The stabilization in the AGQ tetrads due to cation binding is reassessed by employing ab initio simulations. Mixed QM/MM calculations confirmed that Mg²⁺ binds strongly to AGQ compared to K⁺, and it induces higher interactions between the guanine tetrads. However, K⁺ binding to AGQ induces a higher stabilization energy than Mg²⁺ binding to AGQ tetrads. Despite the higher binding energy, Mg²⁺ binding imparts lower stabilization to AGQ due to its unfavorable fermionic quantum energy.

Challenges in Protein QM/MM Simulations with Intra-Backbone Link Atoms

Hybrid quantum mechanical/molecular mechanical (QM/MM) simulations fuel discoveries in many fields of science including computational biochemistry and enzymology. Development of more convenient tools leads to an increase in the number of works in which mechanical insights into enzymes' mode of operation are obtained. Most commonly, these tools feature hydrogen-capping (link atom) approach to provide coupling between QM and MM subsystems across a covalent bond. Extensive studies were conducted to provide a solid foundation for the correctness of such an approach when a bond to a nonpolar MM atom is considered. However, not every task may be accomplished this way. Certain scenarios of using QM/MM in computational enzymology encourage or even necessitate the incorporation of backbone atoms into the QM region. Two out of three backbone atoms are polar, and in QM/MM with electrostatic embedding, a neighboring link atom will be hyperpolarized. Several schemes to mitigate this effect were previously proposed alongside a rigorous assessment of quantitative effects on model systems. However, it was not clear whether they may translate into qualitatively different results and how link atom hyperpolarization may manifest itself in a real-life enzymological scenario. Here, we show that the consequences of such an artifact may be severe and may completely overturn the conclusions drawn from the simulations. Our case advocates for the use of charge redistribution schemes whenever intra-backbone QM/MM boundaries are considered. Moreover, we addressed how different boundary types and charge redistribution schemes influence backbone dynamics. We showed that the results are heavily dependent on which boundary MM terms are retained, with charge alteration being of secondary importance. In the worst case, only three intra-backbone boundaries may be used with relative confidence in the adequacy of resulting simulations, irrespective of the hyperpolarization mitigation scheme. Thus, advances in the field are certainly needed to fuel new discoveries. As of now, we believe that issues raised in this work might encourage authors in the field to report what boundaries, boundary MM terms, and charge redistribution schemes they are using, so their results may be correctly interpreted.

An efficient electrostatic embedding QM/MM method using periodic boundary conditions based on particle-mesh Ewald sums and electrostatic potential fitted charge operators

Hybrid quantum mechanics/molecular mechanics (QM/MM) models are successful at describing the properties and reactivity of biological macromolecules. Combining ab initio QM/MM methods and periodic boundary conditions (PBC) is currently the optimal approach for modeling chemical processes in an infinite environment, but frequently, these models are too time-consuming for general applicability to biological systems in a solution. Here, we define a simple and efficient electrostatic embedding QM/MM model in PBC, combining the benefits of electrostatic potential fitted atomic charges and particle-mesh Ewald sums, which can efficiently treat systems of an arbitrary size at a reasonable computational cost. To illustrate this, we apply our scheme to extract the lowest singlet excitation energies from a model for Arabidopsis thaliana cryptochrome 1 containing circa 93 000 atoms, accurately reproducing the experimental absorption maximum.

Conformational transitions induced by NADH binding promote reduction half-reaction in 2-hydroxybiphenyl-3-monooxygenase catalytic cycle

2-Hydroxybiphenyl-3-monoxygenase from Pseudomonas azelaica is an effective catalyst of the regiospecific conversions of various aromatic compounds. A comprehensive understanding of the complete catalytic cycle, including the as yet unclear details of NADH binding, NADH/FAD interaction as well as related conformational changes could facilitate the rational design of improved enzyme variants for practical applications. Induced fit formation of a specific pocket for the nicotinamide ring at NADH binding has been revealed using advanced molecular simulation methods including metadynamics and QM/MM modeling. The resulting triple stacking interaction of the nicotinamide as well as isoalloxazine rings and evolutionarily correlated amino acid residues of the active site greatly contributes to the stabilization of the charge-transfer complex and determines the Pro-S stereospecificity of the hydride transfer and the low energy barrier 11 kcal/mol. Then the resulting FADH^- anion undergoes the consequent conformational transition of the FAD isoalloxazine ring from the open out to the closed in position which is followed by the binding of an oxygen molecule what is crucial for the next step of substrate oxidation and the completion of the catalytic cycle.

On the inhibition of cytochrome P450 3A4 by structurally diversified flavonoids

Cytochrome P450 3A4 (CYP3A4) is the most versatile enzyme involved in drug metabolism. The time-dependent inhibition of CYP3A4 by acacetin, apigenin, chrysin, and pinocembrin was experimentally detected, but not entirely elaborated so far. Thus, a two-level QM/MM (Quantum Mechanics/Molecular Mechanics) model is developed to yield insights into the receptor-flavonoid recognition at the molecular scale. Active site residues and the flavonoid are modelled using SCC-DFTB-D (QM level), while the rest of the complex is treated using AMBER force field (MM level). QM/MM binding free energies are well correlated with experimental data, indicating the largest inhibitory effect of chrysin on enzyme activity at a submicromolar concentration. Consequently, quercetin (QUE) and flavopiridol (FLP) are observed as representative examples of structurally different flavonoids. The inhibition parameters for QUE and FLP are evaluated using the well-calibrated QM/MM strategy, thereby aiding to quantitatively conceive the functional behavior of the whole family of flavonoids. A kinetic threshold for further assessment of the drug-drug interactions underlying the time-dependent inhibition of CYP3A4 by flavonoids is explored.Communicated by Ramaswamy H. Sarma.

An understanding of coronavirus and exploring the molecular dynamics simulations to find promising candidates against the Mpro of nCoV to combat the COVID-19: A systematic review

The first infection case of new coronavirus was reported at the end of 2019 and after then, the cases are reported in all nations across the world in a very short period. Further, the regular news of mutations in the virus has made life restricted with appropriate behavior. To date, a new strain (Omicron and its new subvariant Omicron XE) has brought fear amongst us due to a higher trajectory of increase in the number of cases. The researchers thus started giving attention to this viral infection and discovering drug-like candidates to cure the infections. Finding a drug for any viral infection is not an easy task and takes plenty of time. Therefore, computational chemistry/bioinformatics is followed to get promising molecules against viral infection. Molecular dynamics (MD) simulations are being explored to get drug candidates in a short period. The molecules are screened via molecular docking, which provides preliminary information which can be further verified by molecular dynamics (MD) simulations. To understand the change in structure, MD simulations generated several trajectories such as root mean square deviation (RMSD), root mean square fluctuation (RMSF), hydrogen bonding, and radius of gyration for the main protease (Mpro) of the new coronavirus (nCoV) in the presence of small molecules. Additionally, change in free energy for the formation of complex of Mpro of nCoV with the small molecule can be determined by applying molecular mechanics with generalized born and surface area solvation (MM-GBSA). Thus, the promising molecules can be further explored for clinical trials to combat coronavirus disease-19 (COVID-19).

Glycosyltransferase engineering and multi-glycosylation routes development facilitating synthesis of high-intensity sweetener mogrosides

Mogrosides are widely served as natural zero-calorie sweeteners. To date, the biosynthesis of high-intensity sweetness mogrosides V from mogrol has not been achieved because of inefficient and uncontrollable multi-glycosylation process. To address this challenge, we reported three UDP-glycosyltransferases (UGTs) catalyzing the primary and branched glycosylation of mogrosides and increased the catalytic efficiency by 74-400-folds toward branched glycosylation using an activity-based sequence conservative analysis engineering strategy. The computational studies provided insights into the origin of improved catalytic activity. By virtue of UGT mutants, we provided regio- and bond-controllable multi-glycosylation routes, successfully facilitating sequential glycosylation of mogrol to three kinds of mogroside V in excellent yield of 91-99%. Meanwhile, the feasibility of the routes was confirmed in engineered yeasts. It suggested that the multi-glycosylation routes would be combined with mogrol synthetic pathway to de novo produce mogrosides from glucose by aid of metabolic engineering and synthetic biology strategies in the future.

Enhanced Sampling Simulation Reveals How Solvent Influences Chirogenesis of the Intra-Molecular Diels-Alder Reaction

The timescale involved in chemical reactions is quite often beyond that of normal molecular dynamics simulations. Here, we combine metadynamics with selective integrated tempering sampling to simulate an intra-molecular Diels-Alder reaction in explicit solvents. Based on a one-dimensional collective variable obtained from harmonic linear discriminant analysis, four chiral isomers of products were observed in the simulation. Analyses of reactive trajectories showed that this reaction follows a concerted mechanism in all four solvents. In addition, the hydrogen bond between the reactant and water solvent plays an important role in the water-accelerated reaction mechanism. The dynamics of chirality formation varies significantly with solvents. The chirality of products forms significantly before the transition state, especially in ionic liquid.

QM/MM investigation of the catalytic mechanism of processive endoglucanase Cel9G from Clostridium cellulovorans

Carbohydrate degradation catalyzed by glucoside hydrolases (GHs) is a major mechanism in biomass conversion. GH family 9 endoglucanase (Cel9G) from Clostridium cellulovorans, a typical multimodular enzyme, contains a catalytic domain closely linked to a family 3c carbohydrate-binding module (CBM3c). Unlike the conventional behavior proposed for other carbohydrate-binding modules, CBM3c has a direct impact on catalytic activity. In this work, extensive molecular dynamics (MD) simulations were employed to clarify the functional role of CBM3c. Furthermore, the detailed catalytic mechanism of Cel9G was investigated at the atomistic level using the combined quantum mechanical and molecular mechanical (QM/MM) method. Based on these simulations, owing to the rigidity of the peptide linker, CBM3c may affect the enzymatic activity via direct interactions with alpha helix 4 of GH9, especially with the K123 and H125 residues. In addition, using cellohexaose as a substrate, the QM/MM MD simulations confirmed that this enzyme can cleave the β-1,4-glycosidic linkage via an inverting mechanism. An oxocarbenium ion-like transition state was located with a barrier height of 19.6 kcal mol^-1. Furthermore, the G(-1) pyranose unit preferentially adopted a distorted ¹S₅/⁴H₅ conformer in the enzyme-substrate complex. For the cleavage of the glycosidic bond, we were able to identify a plausible route (¹S₅/⁴H₅ → [⁴H₅/⁴E]^# → ⁴C₁) from the reactant to the product at the G(-1) site.

Deciphering the Spectral Tuning Mechanism in Proteorhodopsin: The Dominant Role of Electrostatics Instead of Chromophore Geometry

Proteorhodopsin (PR) is a photoactive proton pump found in marine bacteria. There are two phenotypes of PR exhibiting an environmental adaptation to the ocean's depth which tunes their maximum absorption: blue‐absorbing proteorhodopsin (BPR) and green‐absorbing proteorhodopsin (GPR). This blue/green color‐shift is controlled by a glutamine to leucine substitution at position 105 which accounts for a 20 nm shift. Typically, spectral tuning in rhodopsins is rationalized by the external point charge model but the Q105L mutation is charge neutral. To study this tuning mechanism, we employed the hybrid QM/MM method with sampling from molecular dynamics. Our results reveal that the positive partial charge of glutamine near the C₁₄−C₁₅ bond of retinal shortens the effective conjugation length of the chromophore compared to the leucine residue. The derived mechanism can be applied to explain the color regulation in other retinal proteins and can serve as a guideline for rational design of spectral shifts.

Path to Actinorhodin: Regio- and Stereoselective Ketone Reduction by a Type II Polyketide Ketoreductase Revealed in Atomistic Detail

In type II polyketide synthases (PKSs), which typically biosynthesize several antibiotic and antitumor compounds, the substrate is a growing polyketide chain, shuttled between individual PKS enzymes, while covalently tethered to an acyl carrier protein (ACP): this requires the ACP interacting with a series of different enzymes in succession. During biosynthesis of the antibiotic actinorhodin, produced by Streptomyces coelicolor, one such key binding event is between an ACP carrying a 16-carbon octaketide chain (actACP) and a ketoreductase (actKR). Once the octaketide is bound inside actKR, it is likely cyclized between C7 and C12 and regioselective reduction of the ketone at C9 occurs: how these elegant chemical and conformational changes are controlled is not yet known. Here, we perform protein-protein docking, protein NMR, and extensive molecular dynamics simulations to reveal a probable mode of association between actACP and actKR; we obtain and analyze a detailed model of the C7-C12-cyclized octaketide within the actKR active site; and we confirm this model through multiscale (QM/MM) reaction simulations of the key ketoreduction step. Molecular dynamics simulations show that the most thermodynamically stable cyclized octaketide isomer (7R,12R) also gives rise to the most reaction competent conformations for ketoreduction. Subsequent reaction simulations show that ketoreduction is stereoselective as well as regioselective, resulting in an S-alcohol. Our simulations further indicate several conserved residues that may be involved in selectivity of C7-12 cyclization and C9 ketoreduction. Detailed insights obtained on ACP-based substrate presentation in type II PKSs can help design ACP-ketoreductase systems with altered regio- or stereoselectivity.

Multiscale Simulations Identify Origins of Differential Carbapenem Hydrolysis by the OXA-48 β-Lactamase

OXA-48 β-lactamases are frequently encountered in bacterial infections caused by carbapenem-resistant Gram-negative bacteria. Due to the importance of carbapenems in the treatment of healthcare-associated infections and the increasingly wide dissemination of OXA-48-like enzymes on plasmids, these β-lactamases are of high clinical significance. Notably, OXA-48 hydrolyzes imipenem more efficiently than other commonly used carbapenems, such as meropenem. Here, we use extensive multiscale simulations of imipenem and meropenem hydrolysis by OXA-48 to dissect the dynamics and to explore differences in the reactivity of the possible conformational substates of the respective acylenzymes. Quantum mechanics/molecular mechanics (QM/MM) simulations of the deacylation reaction for both substrates demonstrate that deacylation is favored when the 6α-hydroxyethyl group is able to hydrogen bond to the water molecule responsible for deacylation but disfavored by the increasing hydration of either oxygen of the carboxylated Lys73 general base. Differences in free energy barriers calculated from the QM/MM simulations correlate well with the experimentally observed differences in hydrolytic efficiency between meropenem and imipenem. We conclude that the impaired breakdown of meropenem, compared to imipenem, which arises from a subtle change in the hydrogen bonding pattern between the deacylating water molecule and the antibiotic, is most likely induced by the meropenem 1β-methyl group. In addition to increased insights into carbapenem breakdown by OXA β-lactamases, which may aid in future efforts to design antibiotics or inhibitors, our approach exemplifies the combined use of atomistic simulations in determining the possible different enzyme-substrate substates and their influence on enzyme reaction kinetics.

DFT-based QM/MM with Particle-Mesh Ewald for Direct, Long-Range Electrostatic Embedding

We present a DFT-based, QM/MM implementation with long-range electrostatic embedding achieved by direct real-space integration of the particle mesh Ewald (PME) computed electrostatic potential. The key transformation is the interpolation of the electrostatic potential from the PME grid to the DFT quadrature grid, from which integrals are easily evaluated utilizing standard DFT machinery. We provide benchmarks of the numerical accuracy with choice of grid size and real-space corrections, and demonstrate that good convergence is achieved while introducing nominal computational overhead. Furthermore, the approach requires only small modification to existing software packages, as is demonstrated with our implementation in the OpenMM and Psi4 software. After presenting convergence benchmarks, we evaluate the importance of long-range electrostatic embedding in three solute/solvent systems modeled with QM/MM. Water and BMIM/BF4 ionic liquid were considered as ``simple' and ``complex' solvents respectively, with water and p-phenylenediamine (PPD) solute molecules treated at QM level of theory. While electrostatic embedding with standard real-space truncation may introduce negligible error for simple systems such as water solute in water solvent, errors become more significant when QM/MM is applied to complex solvents such as ionic liquids. An extreme example is the electrostatic embedding energy for oxidized PPD in BMIM/BF4 for which real-space truncation produces severe error even at 2-3 nm cutoff distances. This latter example illustrates that utilization of QM/MM to compute redox potentials within concentrated electrolytes/ionic media requires carefully chosen long-range electrostatic embedding algorithms, with our presented algorithm providing a general and robust approach.

Specificity and Reactivity of Mycobacterium tuberculosis Serine/Threonine Kinases PknG and PknB

Mycobacterium tuberculosis (Mtb), the causative agent of Tuberculosis, has 11 eukaryotic-like serine/threonine protein kinases, which play essential roles in cell growth, signal transduction, and pathogenesis. Protein kinase G (PknG) regulates the carbon and nitrogen metabolism by phosphorylation of the glycogen accumulation regulator (GarA) protein at Thr21. Protein kinase B (PknB) is involved in cell wall synthesis and cell shape, as well as phosphorylates GarA but at Thr22. While PknG seems to be constitutively activated and recognition of GarA requires phosphorylation in its unstructured tail, PknB activation is triggered by phosphorylation of its activation loop, which allows binding of the forkhead-associated domain of GarA. In the present work, we used molecular dynamics and quantum-mechanics/molecular mechanics simulations of the catalytically competent complex and kinase activity assays to understand PknG/PknB specificity and reactivity toward GarA. Two hydrophobic residues in GarA, Val24 and Phe25, seem essential for PknG binding and allow specificity for Thr21 phosphorylation. On the other hand, phosphorylated residues in PknB bind Arg26 in GarA and regulate its specificity for Thr22. We also provide a detailed analysis of the free energy profile for the phospho-transfer reaction and show why PknG has a constitutively active conformation not requiring priming phosphorylation in contrast to PknB. Our results provide new insights into these two key enzymes relevant for Mtb and the mechanisms of serine/threonine phosphorylation in bacteria.

Solvation Effects on the Thermal Helix Inversion of Molecular Motors from QM/MM Calculations

Molecular motors convert light and thermal energies into mechanical work, offering good opportunities to design novel molecular devices. Among them, molecular motors alternate a photoisomerization and a thermal helix inversion to achieve unidirectional rotation. The rotational speed is limited by the helix inversion step, which in turn is governed by a barrier in the electronic ground state. In this work, we systematically study the solvation effect on the thermal process of selected molecular motors, comparing reaction barriers obtained from both density functional theory (DFT) in the isolated system and umbrella sampling within a hybrid quantum mechanics/molecular mechanics (QM/MM) model in solution. We find more prominent solvation effects on those molecular motors with larger dipole moments. The results could provide insight into how to functionalize molecular motors to speed up their rotation.

Modeling the syn-cycle in the light activated opening of the channelrhodopsin-2 ion channel

The ion channel of channelrhodopsin-2 (ChR2) is activated by absorbing light. The light stimulates retinal to isomerize to start the photocycle. There are two pathways for photocycles, which are caused by isomerization of the retinal from all-trans, 15-anti to 13-cis, 15-anti in the dark-adapted state (anti-cycle) and from 13-cis, 15-syn to all-trans, 15-syn in the light-adapted state (syn-cycle). In this work, the structure of the syn-cycle intermediate and mechanism of channel opening were studied by molecular dynamics (MD) and steered molecular dynamics (SMD) simulations. Due to the lack of crystal structure of intermediates in the syn-cycle of ChR2, the intermediate models were constructed from the homologous intermediates in the anti-cycle. The isomerization of retinal was shown to cause the central gate (CG) hydrogen bond network to rearrange, cutting the link between TM2 and TM7. TM2 is moved by the intrahelical hydrogen bond of E90 and K93, and induced the intracellular gate (ICG) to expand. The ion penetration pathway between TM1, TM2, TM3 and TM7 in the P500* state was observed by MD simulations. However, this channel is not fully opened compared with the homologous P500 state in the anti-cycle. In addition, the protons on Schiff bases were found to be unable to form hydrogen bonds with the counter residues (E123 and D253) in the P500* state, preventing an evolution of the P500* state to a P390-like state in the syn-cycle.

Modelling the syn-cycle is a series of operations on the ChR2 crystal structure (PDB ID: 6EID). By replacement and isomerization, we obtained P500* and P480 intermediates. A feasible explanation that no P390* was observed in experiment was inferred.

Exploring the Catalytic Mechanism of the RNA Cap Modification by nsp16-nsp10 Complex of SARS-CoV-2 through a QM/MM Approach

The inhibition of key enzymes that may contain the viral replication of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) have assumed central importance in drug discovery projects. Nonstructural proteins (nsps) are essential for RNA capping and coronavirus replication since it protects the virus from host innate immune restriction. In particular, nonstructural protein 16 (nsp16) in complex with nsp10 is a Cap-0 binding enzyme. The heterodimer formed by nsp16-nsp10 methylates the 5′-end of virally encoded mRNAs to mimic cellular mRNAs and thus it is one of the enzymes that is a potential target for antiviral therapy. In this study, we have evaluated the mechanism of the 2′-O methylation of the viral mRNA cap using hybrid quantum mechanics/molecular mechanics (QM/MM) approach. It was found that the calculated free energy barriers obtained at M062X/6-31+G(d,p) is in agreement with experimental observations. Overall, we provide a detailed molecular analysis of the catalytic mechanism involving the 2′-O methylation of the viral mRNA cap and, as expected, the results demonstrate that the TS stabilization is critical for the catalysis.

Crystallography and QM/MM Simulations Identify Preferential Binding of Hydrolyzed Carbapenem and Penem Antibiotics to the L1 Metallo-β-Lactamase in the Imine Form

Widespread bacterial resistance to carbapenem antibiotics is an increasing global health concern. Resistance has emerged due to carbapenem-hydrolyzing enzymes, including metallo-β-lactamases (MβLs), but despite their prevalence and clinical importance, MβL mechanisms are still not fully understood. Carbapenem hydrolysis by MβLs can yield alternative product tautomers with the potential to access different binding modes. Here, we show that a combined approach employing crystallography and quantum mechanics/molecular mechanics (QM/MM) simulations allow tautomer assignment in MβL:hydrolyzed antibiotic complexes. Molecular simulations also examine (meta)stable species of alternative protonation and tautomeric states, providing mechanistic insights into β-lactam hydrolysis. We report the crystal structure of the hydrolyzed carbapenem ertapenem bound to the L1 MβL from Stenotrophomonas maltophilia and model alternative tautomeric and protonation states of both hydrolyzed ertapenem and faropenem (a related penem antibiotic), which display different binding modes with L1. We show how the structures of both complexed β-lactams are best described as the (2S)-imine tautomer with the carboxylate formed after β-lactam ring cleavage deprotonated. Simulations show that enamine tautomer complexes are significantly less stable (e.g., showing partial loss of interactions with the L1 binuclear zinc center) and not consistent with experimental data. Strong interactions of Tyr32 and one zinc ion (Zn1) with ertapenem prevent a C6 group rotation, explaining the different binding modes of the two β-lactams. Our findings establish the relative stability of different hydrolyzed (carba)penem forms in the L1 active site and identify interactions important to stable complex formation, information that should assist inhibitor design for this important antibiotic resistance determinant.

The Impact of Retinal Configuration on the Protein-Chromophore Interactions in Bistable Jumping Spider Rhodopsin-1

Bistable rhodopsins have two stable forms that can be interconverted by light. Due to their ability to act as photoswitches, these proteins are considered as ideal candidates for applications such as optogenetics. In this work, we analyze a recently crystalized bistable rhodopsin, namely the jumping spider rhodopsin-1 (JSR1). This rhodopsin exhibits identical absorption maxima for the parent and the photoproduct form, which impedes its broad application. We performed hybrid QM/MM simulations to study three isomers of the retinal chromophore: the 9-cis, 11-cis and all-trans configurations. The main aim was to gain insight into the specific interactions of each isomer and their impact on the absorption maximum in JSR1. The absorption spectra were computed using sampled snapshots from QM/MM molecular dynamics trajectories and compared to their experimental counterparts. The chromophore-protein interactions were analyzed by visualizing the electrostatic potential of the protein and projecting it onto the chromophore. It was found that the distance between a nearby tyrosine (Y126) residue plays a larger role in the predicted absorption maximum than the primary counterion (E194). Geometric differences between the isomers were also noted, including a structural change in the polyene chain of the chromophore, as well as changes in the nearby hydrogen bonding network.

Exploring Dihydroflavonol-4-Reductase Reactivity and Selectivity by QM/MM-MD Simulations

Flavonoids are secondary metabolites ubiquitously found in plants. Their antioxidant properties make them highly interesting natural compounds for use in pharmacology. Therefore, unravelling the mechanisms of flavonoid biosynthesis is an important challenge. Among all the enzymes involved in this biosynthetic pathway, dihydroflavonol-4-reductase (DFR) plays a key role in the production of anthocyanins and proanthocyanidins. Here, we provide new information on the mechanism of action of this enzyme by using QM/MM-MD simulations applied to both dihydroquercetin (DHQ) and dihydrokaempferol (DHK) substrates. The consideration of these very similar compounds shed light on the major role played by the enzyme on the stabilization of the transition state but also on the activation of the substrate before the reaction through near-attack conformer effects.

2-Mercaptomethyl Thiazolidines (MMTZs) Inhibit All Metallo-β-Lactamase Classes by Maintaining a Conserved Binding Mode

Metallo-β-lactamase (MBL) production in Gram-negative bacteria is an important contributor to β-lactam antibiotic resistance. Combining β-lactams with β-lactamase inhibitors (BLIs) is a validated route to overcoming resistance, but MBL inhibitors are not available in the clinic. On the basis of zinc utilization and sequence, MBLs are divided into three subclasses, B1, B2, and B3, whose differing active-site architectures hinder development of BLIs capable of "cross-class" MBL inhibition. We previously described 2-mercaptomethyl thiazolidines (MMTZs) as B1 MBL inhibitors (e.g., NDM-1) and here show that inhibition extends to the clinically relevant B2 (Sfh-I) and B3 (L1) enzymes. MMTZs inhibit purified MBLs in vitro (e.g., Sfh-I, K_i 0.16 μM) and potentiate β-lactam activity against producer strains. X-ray crystallography reveals that inhibition involves direct interaction of the MMTZ thiol with the mono- or dizinc centers of Sfh-I/L1, respectively. This is further enhanced by sulfur-π interactions with a conserved active site tryptophan. Computational studies reveal that the stereochemistry at chiral centers is critical, showing less potent MMTZ stereoisomers (up to 800-fold) as unable to replicate sulfur-π interactions in Sfh-I, largely through steric constraints in a compact active site. Furthermore, in silico replacement of the thiazolidine sulfur with oxygen (forming an oxazolidine) resulted in less favorable aromatic interactions with B2 MBLs, though the effect is less than that previously observed for the subclass B1 enzyme NDM-1. In the B3 enzyme L1, these effects are offset by additional MMTZ interactions with the protein main chain. MMTZs can therefore inhibit all MBL classes by maintaining conserved binding modes through different routes.

DNA repair glycosylase hNEIL1 triages damaged bases via competing interaction modes

DNA glycosylases must distinguish the sparse damaged sites from the vast expanse of normal DNA bases. However, our understanding of the nature of nucleobase interrogation is still limited. Here, we show that hNEIL1 (human endonuclease VIII-like 1) captures base lesions via two competing states of interaction: an activated state that commits catalysis and base excision repair, and a quarantine state that temporarily separates and protects the flipped base via auto-inhibition. The relative dominance of the two states depends on key residues of hNEIL1 and chemical properties (e.g. aromaticity and hydrophilicity) of flipped bases. Such a DNA repair mechanism allows hNEIL1 to recognize a broad spectrum of DNA damage while keeps potential gratuitous repair in check. We further reveal the molecular basis of hNEIL1 activity regulation mediated by post-transcriptional modifications and provide an example of how exquisite structural dynamics serves for orchestrated enzyme functions.

hNEIL1 (human endonuclease VIII-like 1) is a broadly specific DNA glycosylase for base excision repair. Here, the authors show that hNEIL1 can assume activated or triage conformations: the structural basis for the mechanism that enables broad specificity and reduces futile repair of normal bases.

Finite-temperature effect in the O-acylation of (R,S)-propranolol catalyzed by Candida antarctica lipase B

CalB is a triacylglycerol hydrolase (E.C.3.1.1.3) used in the O-acylation of the beta-adrenergic blocking agent (R,S)-propranolol. The catalytic mechanism involves two steps: enzyme acylation and enzyme deacylation. The enantioselectivity of the O-acylation of (R,S)-propranolol originates from the second step, where the acyl-enzyme transfers the acyl group to the racemic substrate. This step occurs via an initial Michaelis complex (MCC) and a tetrahedral intermediate (TI-2). To gain more insight into the molecular basis of this reaction, we performed an exhaustive conformational sampling along the reaction coordinate of the enantioselective step of the reaction (MCC→TI-2→EPC) applying a QM/MM MD protocol (SCC-DFTB/CHARMM) in combination with umbrella sampling and the weighted histogram analysis method. To identify finite temperature effects we compare the PMF and the potential energy pathway. It is found that the effect of the finite temperature in this reaction is a destabilization of the tetrahedral intermediate and an increase of the barrier height of its formation. This increase is higher for the S-enantiomer.

Side-Chain Polarity Modulates the Intrinsic Conformational Landscape of Model Dipeptides

The intrinsic conformational preferences of small peptides may provide additional insight into the thermodynamics and kinetics of protein folding. In this study, we explore the underlying energy landscapes of two model peptides, namely, Ac-Ala-NH₂ and Ac-Ser-NH₂, using geometry-optimization-based tools developed within the context of energy landscape theory. We analyze not only how side-chain polarity influences the structural preferences of the dipeptides, but also other emergent properties of the landscape, including heat capacity profiles, and kinetics of conformational rearrangements. The contrasting topographies of the free energy landscape agree with recent results from Fourier transform microwave spectroscopy experiments, where Ac-Ala-NH₂ was found to exist as a mixture of two conformers, while Ac-Ser-NH₂ remained structurally locked, despite exhibiting an apparently rich conformational landscape.

QM/MM Simulations of Enzymatic Hydrolysis of Cellulose: Probing the Viability of an Endocyclic Mechanism for an Inverting Cellulase

Glycoside hydrolases (GH) cleave carbohydrate glycosidic bonds and play pivotal roles in living organisms and in many industrial processes. Unlike acid-catalyzed hydrolysis of carbohydrates in solution, which can occur either via cyclic or acyclic oxocarbenium-like transition states, it is widely accepted that GH-catalyzed hydrolysis proceeds via a general acid mechanism involving a cyclic oxocarbenium-like transition state with protonation of the glycosidic oxygen. The GH45 subfamily C inverting endoglucanase from Phanerochaete chrysosporium (PcCel45A) defies the classical inverting mechanism as its crystal structure conspicuously lacks a general Asp or Glu base residue. Instead, PcCel45A has an Asn residue, a notoriously weak base in solution, as one of its catalytic residues at position 92. Moreover, unlike other inverting GHs, the relative position of the catalytic residues in PcCel45A impairs the proton abstraction from the nucleophilic water that attacks the anomeric carbon, a key step in the classical mechanism. Here, we investigate the viability of an endocyclic mechanism for PcCel45A using hybrid quantum mechanics/molecular mechanics (QM/MM) simulations, with the QM region treated with the self-consistent-charge density-functional tight-binding level of theory. In this mechanism, an acyclic oxocarbenium-like transition state is stabilized leading to the opening of the glucopyranose ring and formation of an unstable acyclic hemiacetal that can be readily decomposed into hydrolysis product. In silico characterization of the Michaelis complex shows that PcCel45A significantly restrains the sugar ring to the ⁴C₁ chair conformation at the −1 subsite of the substrate binding cleft, in contrast to the classical exocyclic mechanism in which ring puckering is critical. We also show that PcCel45A provides an environment where the catalytic Asn92 residue in its standard amide form participates in a cooperative hydrogen bond network resulting in its increased nucleophilicity due to an increased negative charge on the oxygen atom. Our results for PcCel45A suggest that carbohydrate hydrolysis catalyzed by GHs may take an alternative route from the classical mechanism.

Active site dynamics and catalytic mechanism in arabinan hydrolysis catalyzed by GH43 endo-arabinanase from QM/MM molecular dynamics simulation and potential energy surface

The endo-1,5-α-L-arabinanases, belonging to glycoside hydrolase family 43 (GH43), catalyse the hydrolysis of α-1,5-arabinofuranosidic bonds in arabinose-containing polysaccharides. These enzymes are proposed targets for industrial and medical applications. Here, molecular dynamics (MD), potential energy surface and free energy (potential of mean force) simulations are undertaken using hybrid quantum mechanical/molecular mechanical (QM/MM) potentials to understand the active site dynamics, catalytic mechanism and the electrostatic influence of active site residues of the GH43 endo-arabinanase from G. stearothermophilus. The calculated results give support to the single-displacement mechanism proposed for the inverting GH43 enzymes: first a proton is transferred from the general acid E201 to the substrate, followed by a nucleophilic attack by water, activated by the general base D27, on the anomer carbon. A conformational change (²E ↔E₃ ↔ ⁴E) in the -1 sugar ring is observed involving a transition state featuring an oxocarbenium ion character. Residues D87, K106, H271 are highlighted as potential targets for future mutation experiments in order to increase the efficiency of the reaction. To our knowledge, this is the first QM/MM study providing molecular insights into the glycosidic bond hydrolysis of a furanoside substrate by an inverting GH in solution.Communicated by Ramaswamy H. Sarma.

Transition Path Sampling Study of the Feruloyl Esterase Mechanism

Serine hydrolases cleave peptide and ester bonds and are ubiquitous in nature, with applications in biotechnology, in materials, and as drug targets. The serine hydrolase two-step mechanism employs a serine-histidine-aspartate/glutamate catalytic triad, where the histidine residue acts as a base to activate poor nucleophiles (a serine residue or a water molecule) and as an acid to allow the dissociation of poor leaving groups. This mechanism has been the subject of debate regarding how histidine shuttles the proton from the nucleophile to the leaving group. To elucidate the reaction mechanism of serine hydrolases, we employ quantum mechanics/molecular mechanics-based transition path sampling to obtain the reaction coordinate using the Aspergillus niger feruloyl esterase A (AnFaeA) as a model enzyme. The optimal reaction coordinates include terms involving nucleophilic attack on the carbonyl carbon and proton transfer to, and dissociation of, the leaving group. During the reaction, the histidine residue undergoes a reorientation on the time scale of hundreds of femtoseconds that supports the "moving histidine" mechanism, thus calling into question the "ring flip" mechanism. We find a concerted mechanism, where the transition state coincides with the tetrahedral intermediate with the histidine residue pointed between the nucleophile and the leaving group. Moreover, motions of the catalytic aspartate toward the histidine occur concertedly with proton abstraction by the catalytic histidine and help stabilize the transition state, thus partially explaining how serine hydrolases enable poor nucleophiles to attack the substrate carbonyl carbon. Rate calculations indicate that the second step (deacylation) is rate-determining, with a calculated rate constant of 66 s^-1. Overall, these results reveal the pivotal role of active-site dynamics in the catalytic mechanism of AnFaeA, which is likely similar in other serine hydrolases.

The effect on ion channel of different protonation states of E90 in channelrhodopsin-2: a molecular dynamics simulation

Channelrhodopsin-2 (ChR2) is a cationic channel protein that has been extensively studied in optogenetics. The ion channel is opened via a series of proton transfers and H-bond changes during the photocycle but the detailed mechanism is still unknown. Molecular dynamics (MD) simulations with enhanced sampling were performed on the dark-adapted state (i.e., D470) and two photocycle intermediates (P₁⁵⁰⁰ and P₂³⁹⁰) to study the proton transfer path of the Schiff base and the subsequent conformational changes. The results suggest there are two possible proton transfer pathways from the Schiff base to proton acceptors (i.e., E123 or D253), depending on the protonation of E90. If E90 is protonated in the P₁⁵⁰⁰ state, the proton on the Schiff base will transfer to E123. The polyene chain of 13-cis retinal tilts and opens the channel that detours the blocking central gate (CG) and forms a narrow channel through the transmembrane helices (TM) 2, 3, 6 and 7. In contrast, if E90 deprotonates after retinal isomerization, the primary proton acceptor is D253, and an almost-open channel through TM1, 2, 3 and 7 is generated. The channel diameter is very close to the experimental value. The potential mean force (PMF) suggests that the free energy is extremely low for ions passing through this channel.

With E90 protonated, the proton acceptor of RSBH⁺ is E123 with a narrow channel along TM3; while with E90 deprotonated, proton transfer from RSBH⁺ to D253 generates an approximately open channel along TM2.