Evaluating AM1/d-CB1 for Chemical Glycobiology QM/MM Simulations

Multidimensional Free Energy and Accelerated Quantum Library Methods Provide a Gateway to Glycoenzyme Conformational, Electronic, and Reaction Mechanisms

Enzyme reactions are complex to simulate accurately, and none more so than glycoenzymes (glycosyltransferase and glycosidases). A rigorous sampling of the protein frame and the conformationally plural carbohydrate substrate coupled with an unbiased treatment of the electron dynamics is needed to discover the true reaction landscapes. Here, we demonstrate the effectiveness of two computational methods ported in libraries that we have developed. The first is a flat histogram free energy method called FEARCF capable of multidimensional sampling and rapidly converging to a complete coverage of phase space. The second, the Quantum Supercharger Library (QSL), is a method that accelerates the computation of the ab initio electronic wave function as well as the integral derivatives on graphical processing units (GPUs). These QSL accelerated computations form the core components needed for direct quantum dynamics and QM/MM dynamics when coupled with legacy codes such as GAMESS and NWCHEM, making state of the art hyper-parallel electronic computations in chemistry and chemical biology possible. The combination of QSL (acceleration of ab initio QM computation) and FEARCF (multidimensional hyper-parallel reaction dynamics) makes the simulation of ab initio QM/MM reaction dynamics of enzyme catalysis feasible. Enzymes that process carbohydrates pose an added challenge as their pyranose ring substrates span multidimensional conformational space whose sampling is an intimate function of the catalytic mechanism. Here, we use the pairing of FEARCF and QSL to simulate the catalytic effect of the glycoenzyme β-N-acetylglucosamine transferase (OGT). The reaction mechanism is discovered from a variable three bond reaction surface using SCCDFTB. The role of the OGT in distorting the pyranose ring of β-N-acetylglucosamine (GlcNAc) away from the equilibrium ⁴C₁ chair conformation toward the E₃ envelope needed for the transition state is discovered from its pucker free energy hypersurfaces (or free energy volume, FEV). A complete GlcNAc ring pucker HF 6-31g FEV is constructed from ab initio QM dynamics in vacuum and ab initio QM/MM dynamics in the OGT catalytic domain. The OGT is shown to clearly lower the pathway toward the transition state E₃ ring conformer as well as stabilize it by 1.63 kcal/mol. Illustrated here is the use of QSL accelerated ab initio QM/MM dynamics that thoroughly explores carbohydrate catalyzed reactions through a FEARCF multidimensional sampling of the interdependence between reaction and conformational space. This demonstrates how experimentally inaccessible molecular and electronic mechanisms that underpin enzyme catalysis can be discovered by directly modeling the dynamics of these complex reactions.

Three-Dimensional Structures of Carbohydrates and Where to Find Them

Analysis and systematization of accumulated data on carbohydrate structural diversity is a subject of great interest for structural glycobiology. Despite being a challenging task, development of computational methods for efficient treatment and management of spatial (3D) structural features of carbohydrates breaks new ground in modern glycoscience. This review is dedicated to approaches of chemo- and glyco-informatics towards 3D structural data generation, deposition and processing in regard to carbohydrates and their derivatives. Databases, molecular modeling and experimental data validation services, and structure visualization facilities developed for last five years are reviewed.

Computational analysis of carbohydrate recognition based on hybrid QM/MM modeling: a case study of norovirus capsid protein in complex with Lewis antigen

Norovirus is a major pathogen of nonbacterial acute gastroenteritis in humans and animals. Carbohydrate recognition between norovirus capsid proteins and Lewis antigens is considered to play a critical role in initiating infection of eukaryotic cells. In this article, we first report a detailed atomistic simulation study of the norovirus capsid protein in complex with the Lewis antigen based on ab initio QM/MM combined with MD-FEP simulations. To understand the mechanistic details of ligand binding, we analyzed and compared the carbohydrate recognition mechanism of the wild-type P domain protein with a mutant protein. Small structural differences between two capsid proteins are observed on the weak interaction site of residue 389, which is located on the solvent exposed surface of the P domain. To further clarify affinity differences in ligand binding, we directly evaluated free energy changes of the ligand binding process. Although the mutant protein loses its interaction energy with the Lewis antigen, this small amount of energy penalty is compensated for by an increase in the solvation stability, which is induced by structural reorganization at the ligand binding site on the protein surface. As a sum of these opposite energy components, the mutant P domain obtains a slightly enhanced binding affinity for the Lewis antigen. The present computational study clearly demonstrated that a detailed free energy balance of the interaction energy between the capsid protein and the surrounding aqueous solvent is the mechanistic basis of carbohydrate recognition in the norovirus capsid protein.

Producing DFT/MM enzyme reaction trajectories from SCC-DFTB/MM driving forces to probe the underlying electronics of a glycosyltransferase reaction

The SCC-DFTB/MIO/CHARMM free energy surface for a glycosyltransferase, TcTS, is benchmarked against a DFT/MM reaction trajectory using the same CHARMM MM force field ported to the NWChem package. The popular B3LYP functional, against which the MIO parameter set was parameterized is used to optimize TS structures and run DFT reaction dynamics. A novel approach was used to generate reaction forces from a SCC-DFTB/MIO/CHARMM reaction surface to drive B3LYP/6-31G/MM and B3LYP/6-31G(d)/MM reaction trajectories. Although TS structures compare favorably, differences stemming primarily from a minimal basis set approximation prevented a successful 6-31G(d) FEARCF reaction dynamics trajectory. None the less, the dynamic evolution of the B3LYP/6-31G/MM-computed electron density provided an opportunity to perform NBO analysis along the reaction trajectory. Here, we illustrate that a successful ab initio reaction trajectory is computationally accessible when the underlying potential energy function of the semi-empirical method used to produce driving forces is sufficiently close to the ab initio potential. © 2017 Wiley Periodicals, Inc.

Improved Sugar Puckering Profiles for Nicotinamide Ribonucleoside for Hybrid QM/MM Simulations

The coenzyme nicotinamide adenine dinucleotide (NAD⁺) and its reduced form (NADH) play ubiquitous roles as oxidizing and reducing agents in nature. The binding, and possibly the chemical redox step, of NAD⁺/NADH may be influenced by the cofactor conformational distribution and, in particular, by the ribose puckering of its nicotinamide-ribonucleoside (NR) moiety. In many hybrid quantum mechanics-molecular mechanics (QM/MM) studies of NAD⁺/NADH dependent enzymes, the QM region is treated by semiempirical (SE) methods. Recent work suggests that SE methods do not adequately describe the ring puckering in sugar molecules. In the present work we adopt an efficient and practical strategy to correct for this deficiency for NAD⁺/NADH. We have implemented a cost-effective correction to a SE Hamiltonian by adding a correction potential, which is defined as the difference between an accurate benchmark density functional theory (DFT) potential energy surface (PES) and the SE PES. In practice, this is implemented via a B-spline interpolation scheme for the grid-based potential energy difference surface. We find that the puckering population distributions obtained from free energy QM(SE)/MM simulations are in good agreement with DFT and in fair accord with experimental results. The corrected PES should facilitate a more accurate description of the ribose puckering in the NAD⁺/NADH cofactor in simulations of biological systems.

Revision of the GROMOS 56A6(CARBO) force field: Improving the description of ring-conformational equilibria in hexopyranose-based carbohydrates chains

This article describes a revised version 56A6(CARBO_R) of the GROMOS 56A6(CARBO) force field for hexopyranose-based carbohydrates. The simulated properties of unfunctionalized hexopyranoses are unaltered with respect to 56A6CARBO . In the context of both O1 -alkylated hexopyranoses and oligosaccharides, the revision stabilizes the regular (4) C1 chair for α-anomers, with the opposite effect for β-anomers. As a result, spurious ring inversions observed in α(1→4)-linked chains when using the original 56A6(CARBO) force field are alleviated. The (4) C1 chair is now the most stable conformation for all d-hexopyranose residues, irrespective of the linkage type and anomery, and of the position of the residue along the chain. The methylation of a d-hexopyranose leads to a systematic shift in the ring-inversion free energy ((4) C1 to (1) C4 ) by 7-8 kJ mol(-1), positive for the α-anomers and negative for the β-anomers, which is qualitatively compatible with the expected enhancement of the anomeric effect upon methylation at O1. The ring-inversion free energies for residues within chains are typically smaller in magnitude compared to those of the monomers, and correlate rather poorly with the latter. This suggests that the crowding of ring substituents upon chain formation alters the ring flexibility in a nonsystematic fashion. In general, the description of carbohydrate chains afforded by 56A6(CARBO_R) suggests a significant extent of ring flexibility, i.e., small but often non-negligible equilibrium populations of inverted chairs, and challenges the "textbook" picture of conformationally locked carbohydrate rings.

Enhanced semiempirical QM methods for biomolecular interactions

Recent successes and failures of the application of 'enhanced' semiempirical QM (SQM) methods are reviewed in the light of the benefits and backdraws of adding dispersion (D) and hydrogen-bond (H) correction terms. We find that the accuracy of SQM-DH methods for non-covalent interactions is very often reported to be comparable to dispersion-corrected density functional theory (DFT-D), while computation times are about three orders of magnitude lower. SQM-DH methods thus open up a possibility to simulate realistically large model systems for problems both in life and materials science with comparably high accuracy.

AM1/d-CB1: A Semiempirical Model for QM/MM Simulations of Chemical Glycobiology Systems

A semiempirical method based on the AM1/d Hamiltonian is introduced to model chemical glycobiological systems. We included in the parameter training set glycans and the chemical environment often found about them in glycoenzymes. Starting with RM1 and AM1/d-PhoT models we optimized H, C, N, O, and P atomic parameters targeting the best performing molecular properties that contribute to enzyme catalyzed glycan reaction mechanisms. The training set comprising glycans, amino acids, phosphates and small organic model systems was used to derive parameters that reproduce experimental data or high-level density functional results for carbohydrate, phosphate and amino acid heats of formation, amino acid proton affinities, amino acid and monosaccharide dipole moments, amino acid ionization potentials, water-phosphate interaction energies, and carbohydrate ring pucker relaxation times. The result is the AM1/d-Chemical Biology 1 or AM1/d-CB1 model that is considerably more accurate than existing NDDO methods modeling carbohydrates and the amino acids often present in the catalytic domains of glycoenzymes as well as the binding sites of lectins. Moreover, AM1/d-CB1 computed proton affinities, dipole moments, ionization potentials and heats of formation for transition state puckered carbohydrate ring conformations, observed along glycoenzyme catalyzed reaction paths, are close to values computed using DFT M06-2X. AM1/d-CB1 provides a platform from which to accurately model reactions important in chemical glycobiology.