Search

The Lipari-Szabo generalized order parameter probes the picosecond to nanosecond time scale motions of a protein and is useful for rationalizing a multitude of biological processes such as protein recognition and ligand binding. Although these fast motions are an important and intrinsic property of proteins, it remains unclear what simulation conditions are most suitable to reproduce methyl symmetry axis side chain order parameter data ( S axis 2 ) from molecular dynamics simulations. In this study, we show that, while S axis 2 tends to converge within tens of nanoseconds, it is essential to run 10 to 20 replicas starting from configurations close to the experimental structure to obtain the best agreement with experimental S axis 2 values. Additionally, in a comparison of force fields, AMBER ff14SB outperforms CHARMM36m in accurately capturing these fast time scale motions, and we suggest that the origin of this performance gap is likely attributed to differences in side chain torsional parametrization and not due to differences in the global protein conformations sampled by the force fields. This study provides insight into obtaining accurate and reproducible S axis 2 values from molecular simulations and underscores the necessity of using replica simulations to compute equilibrium properties.

Biocatalysis can be powerful in organic synthesis but is often limited by enzymes' substrate scope and selectivity. Developing a biocatalytic step involves identifying an initial enzyme for the target reaction followed by optimization through rational design, directed evolution, or both. These steps are time consuming, resource-intensive, and require expertise beyond typical organic chemistry. Thus, an effective strategy for streamlining the process from enzyme identification to implementation is essential to expanding biocatalysis. Here, we present a strategy combining bioinformatics-guided enzyme mining and ancestral sequence reconstruction (ASR) to resurrect enzymes for biocatalytic synthesis. Specifically, we achieve an enantioselective synthesis of azaphilone natural products using two ancestral enzymes: a flavin-dependent monooxygenase (FDMO) for stereodivergent oxidative dearomatization and a substrate-selective acyltransferase (AT) for the acylation of the enzymatically installed hydroxyl group. This cascade, stereocomplementary to established chemoenzymatic routes, expands access to enantiomeric linear tricyclic azaphilones. By leveraging the co-occurrence and coevolution of FDMO and AT in azaphilone biosynthetic pathways, we identified an AT candidate, CazE, and addressed its low solubility and stability through ASR, obtaining a more soluble, stable, promiscuous, and reactive ancestral AT (AncAT). Sequence analysis revealed AncAT as a chimeric composition of its descendants with enhanced reactivity likely due to ancestral promiscuity. Flexible receptor docking and molecular dynamics simulations showed that the most reactive AncAT promotes a reactive geometry between substrates. We anticipate that our bioinformatics-guided, ASR-based approach can be broadly applied in target-oriented synthesis, reducing the time required to develop biocatalytic steps and efficiently access superior biocatalysts.

Since its inception nearly a half century ago, CHARMM has been playing a central role in computational biochemistry and biophysics. Commensurate with the developments in experimental research and advances in computer hardware, the range of methods and applicability of CHARMM have also grown. This review summarizes major developments that occurred after 2009 when the last review of CHARMM was published. They include the following: new faster simulation engines, accessible user interfaces for convenient workflows, and a vast array of simulation and analysis methods that encompass quantum mechanical, atomistic, and coarse-grained levels, as well as extensive coverage of force fields. In addition to providing the current snapshot of the CHARMM development, this review may serve as a starting point for exploring relevant theories and computational methods for tackling contemporary and emerging problems in biomolecular systems. CHARMM is freely available for academic and nonprofit research at https://academiccharmm.org/program.

While large library docking has discovered potent ligands for multiple targets, as the libraries have grown the hit lists can become dominated by rare artifacts that cheat our scoring functions. Here, we investigate rescoring top-ranked docked molecules with orthogonal methods to identify these artifacts, exploring implicit solvent models and absolute binding free energy perturbation as cross-filters. In retrospective studies, this approach deprioritized high-ranking nonbinders for nine targets while leaving true ligands relatively unaffected. We tested the method prospectively against hits from docking against AmpC β-lactamase. We prioritized 128 high-ranking molecules for synthesis and testing, a mixture of 39 molecules flagged as likely cheaters and 89 that were plausible inhibitors. None of the predicted cheating compounds inhibited AmpC detectably, while 57% of the 89 plausible compounds did so. As our libraries continue to grow, deprioritizing docking artifacts by rescoring with orthogonal methods may find wide use.

Alchemical free energy methods are useful in computer-aided drug design and computational protein design because they provide rigorous statistical mechanics-based estimates of free energy differences from molecular dynamics simulations. λ dynamics is a free energy method with the ability to characterize combinatorial chemical spaces spanning thousands of related systems within a single simulation, which gives it a distinct advantage over other alchemical free energy methods that are mostly limited to pairwise comparisons. Recently developed methods have improved the scalability of λ dynamics to perturbations at many sites; however, the size of chemical space that can be explored at each individual site has previously been limited to fewer than ten substituents. As the number of substituents increases, the volume of alchemical space corresponding to nonphysical alchemical intermediates grows exponentially relative to the size corresponding to the physical states of interest. Beyond nine substituents, λ dynamics simulations become lost in an alchemical morass of intermediate states. In this work, we introduce new biasing potentials that circumvent excessive sampling of intermediate states by favoring sampling of physical end points relative to alchemical intermediates. Additionally, we present a more scalable adaptive landscape flattening algorithm for these larger alchemical spaces. Finally, we show that this potential enables more efficient sampling in both protein and drug design test systems with up to 24 substituents per site, enabling, for the first time, simultaneous simulation of all 20 amino acids.

While large library docking has discovered potent ligands for multiple targets, as the libraries have grown, the very top of the hit-lists can become populated with artifacts that cheat our scoring functions. Though these cheating molecules are rare, they become ever-more dominant with library growth. Here, we investigate rescoring top-ranked molecules from docking screens with orthogonal methods to identify these artifacts, exploring implicit solvent models and absolute binding free energy perturbation (AB-FEP) as cross-filters. In retrospective studies, this approach deprioritized high-ranking non-binders for nine targets while leaving true ligands relatively unaffected. We tested the method prospectively against results from large library docking AmpC β -lactamase. From the very top of the docking hit lists, we prioritized 128 molecules for synthesis and experimental testing, a mixture of 39 molecules that rescoring flagged as likely cheaters and another 89 that were plausible true actives. None of the 39 predicted cheating compounds inhibited AmpC up to 200 μ M in enzyme assays, while 57% of the 89 plausible true actives did do so, with 19 of them inhibiting the enzyme with apparent K i values better than 50 μ M . As our libraries continue to grow, a strategy of catching docking artifacts by rescoring with orthogonal methods may find wide use in the field.

The AMP transferase, FICD, is an emerging drug target finetuning stress signaling in the endoplasmic reticulum (ER). FICD is a bi-functional enzyme, catalyzing both AMP addition (AMPylation) and removal (deAMPylation) from the ER resident chaperone BiP/GRP78. Despite increasing evidence linking excessive BiP/GRP78 AMPylation to human diseases, small molecules to inhibit pathogenic FICD variants are lacking. Using an in-vitro high-throughput screen, we identify two small-molecule FICD inhibitors, C22 and C73. Both molecules significantly inhibit FICD-mediated BiP/GRP78 AMPylation in intact cells while only weakly inhibiting BiP/GRP78 deAMPylation. C22 and C73 also efficiently inhibit pathogenic FICD variants and improve proinsulin processing in β cells. Our study identifies and validates FICD inhibitors, highlighting a novel therapeutic avenue against pathologic protein AMPylation., Competing Interests: Declaration of Interests The authors declare no competing interests.

Accurate force field parameters, potential energy functions, and receptor-ligand models are essential for modeling the solvation and binding of drug-like molecules to a receptor. A large and ever-growing chemical space of medicinally relevant scaffolds has also required these factors, especially force field parameters, to be highly transferable. Generalized force fields such as the CHARMM General Force Field (CGenFF) and the generalized AMBER force field (GAFF) have accomplished this feat along with other contemporaneous ones like OPLS. Here, we analyze the limits in the parametrization of drug-like small molecules by CGenFF and GAFF in terms of the various functional groups represented within them. Specifically, we link the presence of specific functional groups to the error in the absolute hydration free energy of over 600 small molecules, predicted by alchemical free energy methods implemented in the CHARMM program. Our investigation reveals that molecules with (i) a nitro group in CGenFF and GAFF are, respectively, over- or undersolubilized in aqueous medium, (ii) amine groups are undersolubilized more so in CGenFF than in GAFF, and (iii) carboxyl groups are more oversolubilized in GAFF than in CGenFF. We present our analyses of the potential factors underlying these trends. We also showcase the use of a machine-learning-based approach combined with the SHapley Additive exPlanations framework to attribute these trends to specific functional groups, which can be easily adopted to explore the limits of other general force fields.

Using a combination of multisite λ-dynamics (MSλD) together with in vitro IC 50 assays, we evaluated the polypharmacological potential of a scaffold currently in clinical trials for inhibition of human neutrophil elastase (HNE), targeting cardiopulmonary disease, for efficacious inhibition of Proteinase 3 (PR3), a related neutrophil serine proteinase. The affinities we observe suggest that the dihydropyrimidinone scaffold can serve as a suitable starting point for the establishment of polypharmacologically targeting both enzymes and enhancing the potential for treatments addressing diseases like chronic obstructive pulmonary disease.

The relevance of motions in biological macromolecules has been clear since the early structural analyses of proteins by X-ray crystallography. Computer simulations have been applied to provide a deeper understanding of the dynamics of biological macromolecules since 1976, and are now a standard tool in many labs working on the structure and function of biomolecules. In this mini-review we highlight some areas of current interest and active development for simulations, in particular all-atom molecular dynamics simulations., Competing Interests: Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Charles L. Brooks III reports financial support was provided by NIH grant GM130587. Alexander D MacKerell Jr reports financial support was provided by NIH grant GM131710. Carol B. Post reports financial support was provided by NIH grant RO1GM039478., (Copyright © 2023 The Authors. Published by Elsevier B.V. All rights reserved.)

3-Hydroxyindolenines can be used to access several structural motifs that are featured in natural products and pharmaceutical compounds, yet the chemical synthesis of 3-hydroxyindolenines is complicated by overoxidation, rearrangements, and complex product mixtures. The selectivity possible in enzymatic reactions can overcome these challenges and deliver enantioenriched products. Herein, we present the development of an asymmetric biocatalytic oxidation of 2-arylindole substrates aided by a curated library of flavin-dependent monooxygenases (FDMOs) sampled from an ancestral sequence space, a sequence similarity network, and a deep-learning-based latent space model. From this library of FDMOs, a previously uncharacterized enzyme, Champase, from the Valley fever fungus, Coccidioides immitis strain RS, was found to stereoselectively catalyze the oxidation of a variety of substituted indole substrates. The promiscuity of this enzyme is showcased by the oxidation of a wide variety of substituted 2-arylindoles to afford the respective 3-hydroxyindolenine products in moderate to excellent yields and up to 95:5 er.

Understanding natural protein evolution and designing novel proteins are motivating interest in development of high-throughput methods to explore large sequence spaces. In this work, we demonstrate the application of multisite λ dynamics (MSλD), a rigorous free energy simulation method, and chemical denaturation experiments to quantify evolutionary selection pressure from sequence-stability relationships and to address questions of design. This study examines a mesophilic phylogenetic clade of ribonuclease H (RNase H), furthering its extensive characterization in earlier studies, focusing on E. coli RNase H (ecRNH) and a more stable consensus sequence (AncCcons) differing at 15 positions. The stabilities of 32,768 chimeras between these two sequences were computed using the MSλD framework. The most stable and least stable chimeras were predicted and tested along with several other sequences, revealing a designed chimera with approximately the same stability increase as AncCcons, but requiring only half the mutations. Comparing the computed stabilities with experiment for 12 sequences reveals a Pearson correlation of 0.86 and root mean squared error of 1.18 kcal/mol, an unprecedented level of accuracy well beyond less rigorous computational design methods. We then quantified selection pressure using a simple evolutionary model in which sequences are selected according to the Boltzmann factor of their stability. Selection temperatures from 110 to 168 K are estimated in three ways by comparing experimental and computational results to evolutionary models. These estimates indicate selection pressure is high, which has implications for evolutionary dynamics and for the accuracy required for design, and suggests accurate high-throughput computational methods like MSλD may enable more effective protein design., Competing Interests: Competing interests statement:The authors declare no competing interest.

The COVID-19 pandemic continues to pose a substantial threat to human lives and is likely to do so for years to come. Despite the availability of vaccines, searching for efficient small-molecule drugs that are widely available, including in low- and middle-income countries, is an ongoing challenge. In this work, we report the results of an open science community effort, the "Billion molecules against COVID-19 challenge", to identify small-molecule inhibitors against SARS-CoV-2 or relevant human receptors. Participating teams used a wide variety of computational methods to screen a minimum of 1 billion virtual molecules against 6 protein targets. Overall, 31 teams participated, and they suggested a total of 639,024 molecules, which were subsequently ranked to find 'consensus compounds'. The organizing team coordinated with various contract research organizations (CROs) and collaborating institutions to synthesize and test 878 compounds for biological activity against proteases (Nsp5, Nsp3, TMPRSS2), nucleocapsid N, RdRP (only the Nsp12 domain), and (alpha) spike protein S. Overall, 27 compounds with weak inhibition/binding were experimentally identified by binding-, cleavage-, and/or viral suppression assays and are presented here. Open science approaches such as the one presented here contribute to the knowledge base of future drug discovery efforts in finding better SARS-CoV-2 treatments., (© 2023 The Authors. Molecular Informatics published by Wiley-VCH GmbH.)

Allostery is involved in innumerable biological processes and plays a fundamental role in human disease. Thus, the exploration of allosteric modulation is crucial for research on biological mechanisms and in the development of novel therapeutics. The development of small-molecule allosteric effectors can be used as tools to probe biological mechanisms of interest. One of the main limitations in targeting allosteric sites is the difficulty in uncovering them for specific receptors. Furthermore, upon discovery of novel allosteric modulation, early lead generation is made more difficult as compared to that at orthosteric sites because there is likely no information about the types of molecules that can bind at the site. In the work described here, we present a novel drug discovery pipeline, FASTDock, which allows one to uncover ligandable sites as well as small molecules that target the given site without requiring pre-existing knowledge of ligands that can bind in the targeted site. By using a hierarchical screening strategy, this method has the potential to enable high-throughput screens of an exceptionally large database of targeted ligand space.

CHARMM is rich in methodology and functionality as one of the first programs addressing problems of molecular dynamics and modeling of biological macromolecules and their partners, e.g., small molecule ligands. When combined with the highly developed CHARMM parameters for proteins, nucleic acids, small molecules, lipids, sugars, and other biologically relevant building blocks, and the versatile CHARMM scripting language, CHARMM has been a trendsetting platform for modeling studies of biological macromolecules. To further enhance the utility of accessing and using CHARMM functionality in increasingly complex workflows associated with modeling biological systems, we introduce pyCHARMM, Python bindings, functions, and modules to complement and extend the extensive set of modeling tools and methods already available in CHARMM. These include access to CHARMM function-generated variables associated with the system (psf), coordinates, velocities and forces, atom selection variables, and force field related parameters. The ability to augment CHARMM forces and energies with energy terms or methods derived from machine learning or other sources, written in Python, CUDA, or OpenCL and expressed as Python callable routines is introduced together with analogous functions callable during dynamics calculations. Integration of Python-based graphical engines for visualization of simulation models and results is also accessible. Loosely coupled parallelism is available for workflows such as free energy calculations, using MBAR/TI approaches or high-throughput multisite λ-dynamics (MSλD) free energy methods, string path optimization calculations, replica exchange, and molecular docking with a new Python-based CDOCKER module. CHARMM accelerated platform kernels through the CHARMM/OpenMM API, CHARMM/DOMDEC, and CHARMM/BLaDE API are also readily integrated into this Python framework. We anticipate that pyCHARMM will be a robust platform for the development of comprehensive and complex workflows utilizing Python and its extensive functionality as well as an optimal platform for users to learn molecular modeling methods and practices within a Python-friendly environment such as Jupyter Notebooks.

The COVID-19 pandemic has highlighted the need for new antiviral approaches because many of the currently approved drugs have proven ineffective against mitigating SARS-CoV-2 infections. The host transmembrane serine protease TMPRSS2 is a promising antiviral target because it plays a role in priming the spike protein before viral entry occurs for the most virulent variants. Further, TMPRSS2 has no established physiological role, thereby increasing its attractiveness as a target for antiviral agents. Here, we utilize virtual screening to curate large libraries into a focused collection of potential inhibitors. Optimization of a recombinant expression and purification protocol for the TMPRSS2 peptidase domain facilitates subsequent biochemical screening and characterization of selected compounds from the curated collection in a kinetic assay. In doing so, we identify new noncovalent TMPRSS2 inhibitors that block SARS-CoV-2 infectivity in a cellular model. One such inhibitor, debrisoquine, has high ligand efficiency, and an initial structure-activity relationship study demonstrates that debrisoquine is a tractable hit compound for TMPRSS2., Competing Interests: The authors declare no competing financial interest., (© 2023 American Chemical Society.)

Controlling the selectivity of a reaction is critical for target-oriented synthesis. Accessing complementary selectivity profiles enables divergent synthetic strategies, but is challenging to achieve in biocatalytic reactions given enzymes' innate preferences of a single selectivity. Thus, it is critical to understand the structural features that control selectivity in biocatalytic reactions to achieve tunable selectivity. Here, we investigate the structural features that control the stereoselectivity in an oxidative dearomatization reaction that is key to making azaphilone natural products. Crystal structures of enantiocomplementary biocatalysts guided the development of multiple hypotheses centered on the structural features that control the stereochemical outcome of the reaction; however, in many cases, direct substitutions of active site residues in natural proteins led to inactive enzymes. Ancestral sequence reconstruction (ASR) and resurrection were employed as an alternative strategy to probe the impact of each residue on the stereochemical outcome of the dearomatization reaction. These studies suggest that two mechanisms are active in controlling the stereochemical outcome of the oxidative dearomatization reaction: one involving multiple active site residues in AzaH and the other dominated by a single Phe to Tyr switch in TropB and AfoD. Moreover, this study suggests that the flavin-dependent monooxygenases (FDMOs) adopt simple and flexible strategies to control stereoselectivity, which has led to stereocomplementary azaphilone natural products produced by fungi. This paradigm of combining ASR and resurrection with mutational and computational studies showcases sets of tools for understanding enzyme mechanisms and provides a solid foundation for future protein engineering efforts.

Multisite λ-dynamics (MSλD) is a novel method for the calculation of relative free energies of binding for ligands to their targeted receptors. It can be readily used to examine a large number of molecules with multiple functional groups at multiple sites around a common core. This makes MSλD a powerful tool in structure-based drug design. In the present study, MSλD is applied to calculate the relative binding free energies of 1296 inhibitors to the testis specific serine kinase 1B (TSSK1B), a validated target for male contraception. For this system, MSλD requires significantly fewer computational resources compared to traditional free energy methods like free energy perturbation or thermodynamic integration. From MSλD simulations, we examined whether modifications of a ligand at two different sites are coupled or not. Based on our calculations, we established a quantitative structure-activity relationship (QSAR) for this set of molecules and identified a site in the ligand where further modification, such as adding more polar groups, may lead to increased binding affinity., (© 2023 The Authors. Protein Science published by Wiley Periodicals LLC on behalf of The Protein Society.)

A recently proposed lipid-chaperone hypothesis suggests that free lipid molecules, not bound to membranes, affect the aggregation of amyloidogenic peptides such as amyloid-β (Aβ) peptides, whose aggregates are the hallmarks of Alzheimer's disease. Here, we combine experiments with all-atom molecular dynamics simulations in explicit solvent to explore the effects of neuronal ganglioside GM1, abundant in mammalian brains, on the aggregation of two principal isoforms of Aβ, Aβ40 and Aβ42. Our simulations show that free GM1 forms stable, highly water-soluble complexes with both isoforms, and nuclear magnetic resonance experiments support the formation of well-ordered, structurally compact GM1+Aβ complexes. By simulation, we also show that Aβ40 monomers display a preference for binding to GM1-containing hetero-oligomers over GM1-lacking homo-oligomers, while Aβ42 monomers have the opposite preference. These observations explain why GM1 dose-dependently inhibits Aβ40 aggregation but has no effect on Aβ42 aggregation, as assessed by thioflavin T fluorescence.

Allostery is the phenomenon of coupling between distal binding sites in a protein. Such coupling is at the crux of protein function and regulation in a myriad of scenarios, yet determining the molecular mechanisms of coupling networks in proteins remains a major challenge. Here, we report mechanisms governing pH-dependent myristoyl switching in monomeric hisactophilin, whereby the myristoyl moves between a sequestered state, i.e., buried within the core of the protein, to an accessible state, in which the myristoyl has increased accessibility for membrane binding. Measurements of the pH and temperature dependence of amide chemical shifts reveal protein local structural stability and conformational heterogeneity that accompany switching. An analysis of these measurements using a thermodynamic cycle framework shows that myristoyl-proton coupling at the single-residue level exists in a fine balance and extends throughout the protein. Strikingly, small changes in the stereochemistry or size of core and surface hydrophobic residues by point mutations readily break, restore, or tune myristoyl switch energetics. Synthesizing the experimental results with those of molecular dynamics simulations illuminates atomistic details of coupling throughout the protein, featuring a large network of hydrophobic interactions that work in concert with key electrostatic interactions. The simulations were critical for discerning which of the many ionizable residues in hisactophilin are important for switching and identifying the contributions of nonnative interactions in switching. The strategy of using temperature-dependent NMR presented here offers a powerful, widely applicable way to elucidate the molecular mechanisms of allostery in proteins at high resolution.

The coactivator KIX of CBP uses two binding surfaces to recognize multiple activators and exhibits allostery in ternary complex formation. Activator•coactivator interactions are central to transcriptional regulation, yet the microscopic origins of allostery in dynamic proteins like KIX are largely unknown. Here, we investigate the molecular recognition and allosteric manifestations involved in two KIX ternary systems c-Myb•KIX•MLL and pKID•KIX•MLL. Exploring the hypothesis that binary complex formation prepays an entropic cost for positive cooperativity, we utilize molecular dynamics simulations, side chain methyl order parameters, and differential scanning fluorimetry (DSF) to explore conformational entropy changes in KIX. The protein's configurational micro-states from structural clustering highlight the utility of protein plasticity in molecular recognition and allostery. We find that apo KIX occupies a wide distribution of lowly-populated configurational states. Each binding partner has its own suite of KIX states that it selects, building a model of molecular recognition fingerprints. Allostery is maximized with MLL pre-binding, which corresponds to the observation of a significant reduction in KIX micro-states observed when MLL binds. With all binding partners, the changes in KIX conformational entropy arise predominantly from changes in the most flexible loop. Likewise, we find that a small molecule and mutations allosterically inhibit/enhance activator binding by tuning loop dynamics, suggesting that loop-targeting chemical probes could be developed to alter KIX•activator interactions. Experimentally capturing KIX stabilization is challenging, particularly because of the disordered nature of particular activators. However, DSF melting curves allow for inference of relative entropic changes that occur across complexes, which we compare to our computed entropy changes using simulation methyl order parameters., Competing Interests: The authors have declared that no competing interests exist.

Alchemical free energy methods are playing a growing role in molecular design, both for computer-aided drug design of small molecules and for computational protein design. Multisite λ dynamics (MSλD) is a uniquely scalable alchemical free energy method that enables more efficient exploration of combinatorial alchemical spaces encountered in molecular design, but simulations have typically been limited to a few hundred ligands or sequences. Here, we focus on coupling between sites to enable scaling to larger alchemical spaces. We first discuss updates to the biasing potentials that facilitate MSλD sampling to include coupling terms and show that this can provide more thorough sampling of alchemical states. We then harness coupling between sites by developing a new free energy estimator based on the Potts models underlying direct coupling analysis, a method for predicting contacts from sequence coevolution, and find it yields more accurate free energies than previous estimators. The sampling requirements of the Potts model estimator scale with the square of the number of sites, a substantial improvement over the exponential scaling of the standard estimator. This opens up exploration of much larger alchemical spaces with MSλD for molecular design.

With the ability to sample combinations of alchemical perturbations at multiple sites off a small molecule core, multisite λ-dynamics (MSλD) has become an attractive alternative to conventional alchemical free energy methods for exploring large combinatorial chemical spaces. However, current software implementations dictate that combinatorial sampling with MSλD must be performed with a multiple topology model (MTM), which is nontrivial to create by hand, especially for a series of ligand analogues which may have diverse functional groups attached. This work introduces an automated workflow, referred to as msld_py_prep , to assist in the creation of a MTM for use with MSλD. One approach for partitioning partial atomic charges between ligands to create a MTM, called charge renormalization, is also presented and rigorously evaluated. We find that msld_py_prep greatly accelerates the preparation of MSλD ready-to-use files and that charge renormalization can provide a successful approach for MTM generation, as long as bookending calculations are applied to correct small differences introduced by charge renormalization. Charge renormalization also facilitates the use of many different force field parameters with MSλD, broadening the applicability of MSλD for computer-aided drug design.

Accurate and rapid predictions of the binding affinity of a compound to a target are one of the ultimate goals of computer aided drug design. Alchemical approaches to free energy estimations follow the path from an initial state of the system to the final state through alchemical changes of the energy function during a molecular dynamics simulation. Herein, we explore the accuracy and efficiency of two such techniques: relative free energy perturbation (FEP) and multisite lambda dynamics (MSλD). These are applied to a series of inhibitors for the bromodomain-containing protein 4 (BRD4). We demonstrate a procedure for obtaining accurate relative binding free energies using MSλD when dealing with a change in the net charge of the ligand. This resulted in an impressive comparison with experiment, with an average difference of 0.4 ± 0.4 kcal mol -1 . In a benchmarking study for the relative FEP calculations, we found that using 20 lambda windows with 0.5 ns of equilibration and 1 ns of data collection for each window gave the optimal compromise between accuracy and speed. Overall, relative FEP and MSλD predicted binding free energies with comparable accuracy, an average of 0.6 kcal mol -1 for each method. However, MSλD makes predictions for a larger molecular space over a much shorter time scale than relative FEP, with MSλD requiring a factor of 18 times less simulation time for the entire molecule space.

The binding of small-molecule ligands to protein or nucleic acid targets is important to numerous biological processes. Accurate prediction of the binding modes between a ligand and a macromolecule is of fundamental importance in structure-based structure-function exploration. When multiple ligands with different sizes are docked to a target receptor, it is reasonable to assume that the residues in the binding pocket may adopt alternative conformations upon interacting with the different ligands. In addition, it has been suggested that the entropic contribution to binding can be important. However, only a few attempts to include the side chain conformational entropy upon binding within the application of flexible receptor docking methodology exist. Here, we propose a new physics-based scoring function that includes both enthalpic and entropic contributions upon binding by considering the conformational variability of the flexible side chains within the ensemble of docked poses. We also describe a novel hybrid searching algorithm that combines both molecular dynamics (MD)-based simulated annealing and genetic algorithm crossovers to address the enhanced sampling of the increased search space. We demonstrate improved accuracy in flexible cross-docking experiments compared with rigid cross-docking. We test our developments by considering five protein targets, thrombin, dihydrofolate reductase(DHFR), T4 L99A, T4 L99A/M102Q, and PDE10A, which belong to different enzyme classes with different binding pocket environments, as a representative set of diverse ligands and receptors. Each target contains dozens of different ligands bound to the same binding pocket. We also demonstrate that this flexible docking algorithm may be applicable to RNA docking with a representative riboswitch example. Our findings show significant improvements in top ranking accuracy across this set, with the largest improvement relative to rigid, 23.64%, occurring for ligands binding to DHFR. We then evaluate the ability to identify lead compounds among a large chemical space for the proposed flexible receptor docking algorithm using a subset of the DUD-E containing receptor targets MCR, GCR, and ANDR. We demonstrate that our new algorithms show improved performance in modeling flexible binding site residues compared to DOCK. Finally, we select the T4 L99A and T4 L99A/M102Q decoy sets, containing dozens of binders and experimentally validated nonbinders, to test our approach in distinguishing binders from nonbinders. We illustrate that our new algorithms for searching and scoring have superior performance to rigid receptor CDOCKER as well as AutoDock Vina. Finally, we suggest that flexible CDOCKER is sufficiently fast to be utilized in high-throughput docking screens in the context of hierarchical approaches.

There is an accelerating interest in practical applications of alchemical free energy methods to problems in protein design, constant pH simulations, and especially computer-aided drug design. In the present paper, we describe a basic lambda dynamics engine (BLaDE) that enables alchemical free energy simulations, including multisite λ dynamics (MSλD) simulations, on graphical processor units (GPUs). We find that BLaDE is 5 to 8 times faster than the current GPU implementation of MSλD-based free energy calculations in CHARMM. We also demonstrate that BLaDE running standard molecular dynamics attains a performance competitive with and sometimes exceeding that of the highly optimized OpenMM GPU code. BLaDE is available as a standalone program and through an API in CHARMM.

In this work, the discrete λ variant of the Gibbs sampler-based λ-dynamics ( d -GSλD) method is developed to enable multiple functional group perturbations to be investigated at one or more sites of substitution off a common ligand core. The theoretical framework and special considerations for constructing discrete λ states for multisite d -GSλD are presented. The precision and accuracy of the d -GSλD method is evaluated with three test cases of increasing complexity. Specifically, methyl → methyl symmetric perturbations in water, 1,4-benzene hydration free energies and protein-ligand binding affinities for an example HIV-1 reverse transcriptase inhibitor series are computed with d -GSλD. Complementary MSλD calculations were also performed to compare with d -GSλD's performance. Excellent agreement between d -GSλD and MSλD is observed, with mean unsigned errors of 0.12 and 0.22 kcal/mol for computed hydration and binding free energy test cases, respectively. Good agreement with experiment is also observed, with errors of 0.5-0.7 kcal/mol. These findings support the applicability of the d -GSλD free energy method for a variety of molecular design problems, including structure-based drug design. Finally, a discussion of d -GSλD versus MSλD approaches is presented to compare and contrast features of both methods.

Inhibitors of transcriptional protein-protein interactions (PPIs) have high value both as tools and for therapeutic applications. The PPI network mediated by the transcriptional coactivator Med25, for example, regulates stress-response and motility pathways, and dysregulation of the PPI networks contributes to oncogenesis and metastasis. The canonical transcription factor binding sites within Med25 are large (∼900 Å 2 ) and have little topology, and thus, they do not present an array of attractive small-molecule binding sites for inhibitor discovery. Here we demonstrate that the depsidone natural product norstictic acid functions through an alternative binding site to block Med25-transcriptional activator PPIs in vitro and in cell culture. Norstictic acid targets a binding site comprising a highly dynamic loop flanking one canonical binding surface, and in doing so, it both orthosterically and allosterically alters Med25-driven transcription in a patient-derived model of triple-negative breast cancer. These results highlight the potential of Med25 as a therapeutic target as well as the inhibitor discovery opportunities presented by structurally dynamic loops within otherwise challenging proteins.

Computation of the thermodynamic consequences of protein mutations holds great promise in protein biophysics and design. Alchemical free energy methods can give improved estimates of mutational free energies, and are already widely used in calculations of relative and absolute binding free energies in small molecule design problems. In principle, alchemical methods can address any amino acid mutation with an appropriate alchemical pathway, but identifying a strategy that produces such a path for proline and glycine mutations is an ongoing challenge. Most current strategies perturb only side chain atoms, while proline and glycine mutations also alter the backbone parameters and backbone ring topology. Some strategies also perturb backbone parameters and enable glycine mutations. This work presents a strategy that enables both proline and glycine mutations and comprises two key elements: a dual backbone with restraints and scaling of bonded terms, facilitating backbone parameter changes, and a soft bond in the proline ring, enabling ring topology changes in proline mutations. These elements also have utility for core hopping and macrocycle studies in computer-aided drug design. This new strategy shows slight improvements over an alternative side chain perturbation strategy for a set T4 lysozyme mutations lacking proline and glycine, and yields good agreement with experiment for a set of T4 lysozyme proline and glycine mutations not previously studied. To our knowledge this is the first report comparing alchemical predictions of proline mutations with experiment. With this strategy in hand, alchemical methods now have access to the full palette of amino acid mutations., (© 2021 Wiley Periodicals LLC.)

Acid-base catalysis, which involves one or more proton transfer reactions, is a chemical mechanism commonly employed by many enzymes. The molecular basis for catalysis is often derived from structures determined at the optimal pH for enzyme activity. However, direct observation of protons from experimental structures is quite difficult; thus, a complete mechanistic description for most enzymes remains lacking. Dihydrofolate reductase (DHFR) exemplifies general acid-base catalysis, requiring hydride transfer and protonation of its substrate, DHF, to form the product, tetrahydrofolate (THF). Previous X-ray and neutron crystal structures coupled with theoretical calculations have proposed that solvent mediates the protonation step. However, visualization of a proton transfer has been elusive. Based on a 2.1 Å resolution neutron structure of a pseudo-Michaelis complex of E. coli DHFR determined at acidic pH, we report the direct observation of the catalytic proton and its parent solvent molecule. Comparison of X-ray and neutron structures elucidated at acidic and neutral pH reveals dampened dynamics at acidic pH, even for the regulatory Met20 loop. Guided by the structures and calculations, we propose a mechanism where dynamics are crucial for solvent entry and protonation of substrate. This mechanism invokes the release of a sole proton from a hydronium (H 3 O + ) ion, its pathway through a narrow channel that sterically hinders the passage of water, and the ultimate protonation of DHF at the N5 atom., Competing Interests: The authors declare no competing financial interest., (© 2021 The Authors. Published by American Chemical Society.)

Accurate predictions of changes to protein-ligand binding affinity in response to chemical modifications are of utility in small-molecule lead optimization. Relative free-energy perturbation (FEP) approaches are one of the most widely utilized for this goal but involve significant computational cost, thus limiting their application to small sets of compounds. Lambda dynamics, also rigorously based on the principles of statistical mechanics, provides a more efficient alternative. In this paper, we describe the development of a workflow to set up, execute, and analyze multisite lambda dynamics (MSLD) calculations run on GPUs with CHARMM implemented in BIOVIA Discovery Studio and Pipeline Pilot. The workflow establishes a framework for setting up simulation systems for exploratory screening of modifications to a lead compound, enabling the calculation of relative binding affinities of combinatorial libraries. To validate the workflow, a diverse data set of congeneric ligands for seven proteins with experimental binding affinity data is examined. A protocol to automatically tailor fit biasing potentials iteratively to flatten the free-energy landscape of any MSLD system is developed, which enhances sampling and allows for efficient estimation of free-energy differences. The protocol is first validated on a large number of ligand subsets that model diverse substituents, which shows accurate and reliable performance. The scalability of the workflow is also tested to screen more than 100 ligands modeled in a single system, which also resulted in accurate predictions. With a cumulative sampling time of 150 ns or less, the method results in average unsigned errors of under 1 kcal/mol in most cases for both small and large combinatorial libraries. For the multisite systems examined, the method is estimated to be more than an order of magnitude more efficient than contemporary FEP applications. The results thus demonstrate the utility of the presented MSLD workflow to efficiently screen combinatorial libraries and explore the chemical space around a lead compound and thus are of utility in lead optimization.

When the electrostatic environment surrounding binding partners changes between unbound and bound states, the net uptake or release of a proton is possible by either binding partner. This process is pH-dependent in that the free energy required to uptake or release the proton varies with pH. This pH-dependence is typically not considered in conventional free energy methods where the use of fixed protonation states is the norm. In the present paper, we apply a simple two-step approach to calculate the pH-dependent binding free energy of a model cucubit[7]uril host/guest system. By use of λ-dynamics with an enhanced sampling protocol, adaptive landscape flattening, p K a shifts and reference binding free energies upon complexation were determined. This information enables the construction of pH-dependent binding profiles that accurately capture the p K a shifts and reproduce binding free energies at the different pH conditions that were observed experimentally. Our calculations illustrate a general framework for computing pH-dependent binding free energies but also point to some issues in modeling the molecular charge distributions within this series of molecules with CGenFF. However, by introducing some minor charge modifications to the CGenFF force field, we saw significant improvement in accuracy of the calculated p K a shifts.

Fast Fourier transform (FFT)-based protein ligand docking together with parallel simulated annealing for both rigid and flexible receptor docking are implemented on graphical processing unit (GPU) accelerated platforms to significantly enhance the throughput of the CDOCKER and flexible CDOCKER - the docking algorithms in the CHARMM program for biomolecule modeling. The FFT-based approach for docking, first applied in protein-protein docking to efficiently search for the binding position and orientation of proteins, is adapted here to search ligand translational and rotational spaces given a ligand conformation in protein-ligand docking. Running on GPUs, our implementation of FFT docking in CDOCKER achieves a 15 000 fold speedup in the ligand translational and rotational space search in protein-ligand docking problems. With this significant speedup it becomes practical to exhaustively search ligand translational and rotational space when docking a rigid ligand into a protein receptor. We demonstrate in this paper that this provides an efficient way to calculate an upper bound for docking accuracy in the assessment of scoring functions for protein-ligand docking, which can be useful for improving scoring functions. The parallel molecular dynamics (MD) simulated annealing, also running on GPUs, aims to accelerate the search algorithm in CDOCKER by running MD simulated annealing in parallel on GPUs. When utilized as part of the general CDOCKER docking protocol, acceleration in excess of 20 times is achieved. With this acceleration, we demonstrate that the performance of CDOCKER for redocking is significantly improved compared with three other popular protein-ligand docking programs on two widely used protein ligand complex data sets: the Astex diverse set and the SB2012 test set. The flexible CDOCKER is similarly improved by the parallel MD simulated annealing on GPUs. Based on the results presented here, we suggest that the accelerated CDOCKER platform provides a highly competitive docking engine for both rigid-receptor and flexible-receptor docking studies and will further facilitate continued improvement in the physics-based scoring function employed in CDOCKER docking studies.

The generalized Born with molecular volume and solvent accessible surface area (GBMV2/SA) implicit solvent model provides an accurate description of molecular volume and has the potential to accurately describe the conformational equilibria of structured and disordered proteins. However, its broader application has been limited by the computational cost and poor scaling in parallel computing. Here, we report an efficient implementation of both the electrostatic and nonpolar components of GBMV2/SA on graphics processing unit (GPU) within the CHARMM/OpenMM module. The GPU-GBMV2/SA is numerically equivalent to the original CPU-GBMV2/SA. The GPU acceleration offers ~60- to 70-fold speedup on a single NVIDIA TITAN X (Pascal) graphics card for molecular dynamic simulations of both folded and unstructured proteins of various sizes. The current implementation can be further optimized to achieve even greater acceleration with minimal reduction on the numerical accuracy. The successful development of GPU-GBMV2/SA greatly facilitates its application to biomolecular simulations and paves the way for further development of the implicit solvent methodology. © 2019 Wiley Periodicals, Inc., (© 2019 Wiley Periodicals, Inc.)

The matrix-2 (M2) protein from influenza A virus is a tetrameric, integral transmembrane (TM) protein that plays a vital role in viral replication by proton flux into the virus. The His37 tetrad is a pH sensor in the center of the M2 TM helix that activates the channel in response to the low endosomal pH. M2 consists of different regions that are believed to be involved in membrane targeting, packaging, nucleocapsid binding, and proton transport. Although M2 has been the target of many experimental and theoretical studies that have led to significant insights into its structure and function under differing conditions, the main mechanism of proton transport, its conformational dynamics, and the role of the amphipathic helices (AHs) on proton conductance remain elusive. To this end, we have applied explicit solvent constant pH molecular dynamics using the multisite λ-dynamics approach (CpHMD MSλD ) to investigate the buried ionizable residues comprehensively and to elucidate their effect on the conformational transition. Our model recapitulates the pH-dependent conformational transition of M2 from closed to open state when the AH domain is included in the M2 construct, revealing the role of the amphipathic helices on this transition and shedding light on the proton-transport mechanism. This work demonstrates the importance of including the amphipathic helices in future experimental and theoretical studies of ion channels. Finally, our work shows that explicit solvent CpHMD MSλD provides a realistic pH-dependent model for membrane proteins., Competing Interests: The authors declare no competing interest.

Allosteric regulation by intrinsically disordered proteins (IDPs) is an important class of cellular processes, including transcription. Molecular dynamics (MD) simulation is a promising approach to unravel the complex molecular interactions involved in the allosteric regulation by IDPs. While allosteric regulation is often characterized by the effect of a ligand on the binding affinity of a distal ligand, the binding affinity is often challenging to calculate by MD simulations because of insufficient sampling of the rare events in this binding-unbinding process. In the current work, we present a new sampling approach based on Hamiltonian replica exchange that allows accurate and efficient calculation of binding affinities using a native-centric coarse-grained model. We also demonstrate the utility of the new method by studying the positive allostery in the kinase-inducible domain interacting domain of the CREB binding protein, in which a prebound ligand enhances the binding of the second ligand.

Class A flavin-dependent hydroxylases (FdHs) catalyze the hydroxylation of organic compounds in a site- and stereoselective manner. In stark contrast, conventional synthetic routes require environmentally hazardous reagents and give modest yields. Thus, understanding the detailed mechanism of this class of enzymes is essential to their rational manipulation for applications in green chemistry and pharmaceutical production. Both electrophilic substitution and radical intermediate mechanisms have been proposed as interpretations of FdH hydroxylation rates and optical spectra. While radical mechanistic steps are often difficult to examine directly, modern quantum chemistry calculations combined with statistical mechanical approaches can yield detailed mechanistic models providing insights that can be used to differentiate reaction pathways. In the current work, we report quantum mechanical/molecular mechanical (QM/MM) calculations on the fungal TropB enzyme that shows an alternative reaction pathway in which hydroxylation through a hydroxyl radical-coupled electron-transfer mechanism is significantly favored over electrophilic substitution. Furthermore, QM/MM calculations on several modified flavins provide a more consistent interpretation of the experimental trends in the reaction rates seen experimentally for a related enzyme, para -hydroxybenzoate hydroxylase. These calculations should guide future enzyme and substrate design strategies and broaden the scope of biological spin chemistry.

Yersinia Protein Tyrosine Phosphatase (YopH) is the most efficient enzyme among all known PTPases and relies on its catalytic loop movements for substrate binding and catalysis. Fluorescence, NMR, and UV resonance Raman (UVRR) techniques have been used to study the thermodynamic and dynamic properties of the loop motions. In this study, a computational approach based on the pathway refinement methods nudged elastic band (NEB) and harmonic Fourier beads (HFB) has been developed to provide structural interpretations for the experimentally observed kinetic processes. In this approach, the minimum potential energy pathways for the loop open/closure conformational changes were determined by NEB using a one-dimensional global coordinate. Two dimensional data analyses of the NEB results were performed as an efficient method to qualitatively evaluate the energetics of transitions along several specific physical coordinates. The free energy barriers for these transitions were then determined more precisely using the HFB method. Kinetic parameters were estimated from the energy barriers using transition state theory and compared against experimentally determined kinetic parameters. When the calculated energy barriers are calibrated by a simple "scaling factor", as have been done in our previous vibrational frequency calculations to explain the ligand frequency shift upon its binding to protein, it is possible to make structural interpretations of several observed enzyme dynamic rates. For example, the nanosecond kinetics observed by fluorescence anisotropy may be assigned to the translational motion of the catalytic loop and microsecond kinetics observed in fluorescence T-jump can be assigned to the loop backbone dihedral angle flipping. Furthermore, we can predict that a Trp354 conformational conversion associated with the loop movements would occur on the tens of nanoseconds time scale, to be verified by future UVRR T-jump studies.

Alchemical free energy calculations have made a dramatic impact upon the field of structure-based drug design by allowing functional group modifications to be explored computationally prior to experimental synthesis and assay evaluation, thereby informing and directing synthetic strategies. In furthering the advancement of this area, a series of 21 β-secretase 1 (BACE1) inhibitors developed by Janssen Pharmaceuticals were examined to evaluate the ability to explore large substituent perturbations, some of which contain scaffold modifications, with multisite λ-dynamics (MSλD), an innovative alchemical free energy framework. Our findings indicate that MSλD is able to efficiently explore all structurally diverse ligand end-states simultaneously within a single MD simulation with a high degree of precision and with reduced computational costs compared to the widely used approach TI/MBAR. Furthermore, computational predictions were shown to be accurate to within 0.5-0.8 kcal/mol when CM1A partial atomic charges were combined with CHARMM or OPLS-AA-based force fields, demonstrating that MSλD is force field independent and a viable alternative to FEP or TI approaches for drug design.

Triosephosphate isomerase (TIM) barrel proteins have not only a conserved architecture that supports a myriad of enzymatic functions, but also a conserved folding mechanism that involves on- and off-pathway intermediates. Although experiments have proven to be invaluable in defining the folding free-energy surface, they provide only a limited understanding of the structures of the partially folded states that appear during folding. Coarse-grained simulations employing native centric models are capable of sampling the entire energy landscape of TIM barrels and offer the possibility of a molecular-level understanding of the readout from sequence to structure. We have combined sequence-sensitive native centric simulations with small-angle X-ray scattering and time-resolved Förster resonance energy transfer to monitor the formation of structure in an intermediate in the Sulfolobus solfataricus indole-3-glycerol phosphate synthase TIM barrel that appears within 50 μs and must at least partially unfold to achieve productive folding. Simulations reveal the presence of a major and 2 minor folding channels not detected in experiments. Frustration in folding, i.e., backtracking in native contacts, is observed in the major channel at the initial stage of folding, as well as late in folding in a minor channel before the appearance of the native conformation. Similarities in global and pairwise dimensions of the early intermediate, the formation of structure in the central region that spreads progressively toward each terminus, and a similar rate-limiting step in the closing of the β-barrel underscore the value of combining simulation and experiment to unravel complex folding mechanisms at the molecular level., Competing Interests: The authors declare no conflict of interest.

Modeling the pH dependence of protein and peptide chemical shifts outside the range of physiological values (6.5-7) is key to understanding structure-function relationships of these systems. These capabilities are largely not available in current chemical shift prediction software. In this study, we utilize a combination of molecular dynamics and quantum mechanics to investigate the through-space and through-bond contributions of protonation-dependent chemical shift perturbations (CSPs) in model tripeptides. By altering the protonation state of the titratable group in the tripeptides, we observe a notable difference in the conformational ensembles and attendantly compute significant CSPs for all nuclei near the site of protonation. We thus demonstrate the ability to recapitulate experimental pH-dependent CSPs with good agreement (R = 0.85, 0.99, and 0.98 for 13 C, 15 N, and 1 H, respectively). Broadly, we provide the groundwork for incorporating pH effects into empirical and semiempirical chemical shift predictors., (Copyright © 2019 Biophysical Society. Published by Elsevier Inc. All rights reserved.)

Biocatalysts have the potential to perform reactions with exceptional selectivity and high catalytic efficiency while utilizing safe and sustainable reagents. Despite these positive attributes, the utility of a biocatalyst can be limited by the breadth of substrates that can be accommodated in the active site in a reactive pose. Proven strategies exist for optimizing the performance of a biocatalyst toward unnatural substrates, including protein engineering; however, these methods can be time intensive and require specialized equipment that renders these approaches inaccessible to synthetic chemists. Strategies accessible to chemists for the expansion of a natural enzyme's substrate scope, while maintaining high levels of site- and stereoselectivity, remain elusive. Here, we employ a computationally guided substrate engineering strategy to expand the synthetic utility of a flavin-dependent monooxygenase. Specifically, experimental observations and computational modeling led to the identification of a critical interaction between the substrate and protein which is responsible for orienting the substrate in a pose productive for catalysis. The fundamental hypothesis for this positioning group strategy is supported by binding and kinetic assays as well as computational studies with a panel of compounds. Further, incorporation of this positioning group into substrates through a cleavable ester linkage transformed compounds not oxidized by the biocatalyst SorbC into substrates efficiently oxidatively dearomatized by the wild-type enzyme with the highest levels of site- and stereoselectivity known for this transformation., Competing Interests: The authors declare no competing financial interest.

Biocatalytic reactions embody many features of ideal chemical transformations, including the potential for impeccable selectivity, high catalytic efficiency, mild reaction conditions and the use of environmentally benign reagents. These advantages have created a demand for biocatalysts that expand the portfolio of complexity-generating reactions available to synthetic chemists. However, the tradeoff that often exists between the substrate scope of a biocatalyst and its selectivity limits the application of enzymes in synthesis. We recently demonstrated that a flavin-dependent monooxygenase, TropB, maintains high levels of site- and stereoselectivity across a range of structurally diverse substrates. Herein, we disclose the structural basis for substrate binding in TropB, which performs a synthetically challenging asymmetric oxidative dearomatization reaction with exquisite site- and stereoselectivity across a range of phenol substrates, providing a foundation for future protein engineering and reaction development efforts. Our hypothesis for substrate binding is informed by a crystal structure of TropB and molecular dynamics simulations with the corresponding computational TropB model and is supported by experimental data. In contrast to canonical class A FAD-dependent monooxygenases in which substrates bind in a protonated form, our data indicate that the phenolate form of the substrate binds in the active site. Furthermore, the substrate position is controlled through twopoint binding of the phenolate oxygen to Arg206 and Tyr239, which are shown to have distinct and essential roles in catalysis. Arg206 is involved in the reduction of the flavin cofactor, suggesting a role in flavin dynamics. Further, QM/MM simulations reveal the interactions that govern the facial selectivity that leads to a highly enantioselective transformation. Thus, the structural origins of the high levels of site-and stereoselectivity observed in reactions of TropB across a range of substrates are elucidated, providing a foundation for future protein engineering and reaction development efforts.

The multistate Bennett acceptance ratio method (MBAR) and unbinned weighted histogram analysis method (UWHAM) are widely employed approaches to calculate relative free energies of multiple thermodynamic states that gain statistical precision by employing free energy contributions from configurations sampled at each of the simulated λ states. With the increasing availability of high throughput computing resources, a large number of configurations can be sampled from hundreds or even thousands of states. Combining sampled configurations from all states to calculate relative free energies requires the iterative solution of large scale MBAR/UWHAM equations. In the current work, we describe the development of a fast solver to iteratively solve these large scale MBAR/UWHAM equations utilizing our previous findings that the MBAR/UWHAM equations can be derived as a Rao-Blackwell estimator. The solver is implemented and distributed as a Python module called FastMBAR. Our benchmark results show that FastMBAR is more than 2 times faster than the widely used solver pymbar, when it runs on a central processing unit (CPU) and more than 100 times faster than pymbar when it runs on a graphical processing unit (GPU). The significant speedup achieved by FastMBAR running on a GPU is useful not only for solving large scale MBAR/UWHAM equations but also for estimating uncertainty of calculated free energies using bootstrapping where the MBAR/UWHAM equations need to be solved multiple times.

Immobilization on solid supports provides an effective way to improve enzyme stability and simplify downstream processing for biotechnological applications, which has been widely used in research and in applications. However, surface immobilization may disrupt enzyme structure due to interactions between the enzyme and the supporting substrate, leading to a loss of the enzyme catalytic efficiency and stability. Here, we use a model enzyme, nitroreductase (NfsB), to demonstrate that engineered variants with two strategically positioned surface-tethering sites exhibit improved enzyme stability when covalently immobilized onto a surface. Tethering sites were designed based on molecular dynamics (MD) simulations, and enzyme variants containing cysteinyl residues at these positions were expressed, purified, and immobilized on maleimide-terminated self-assembled monolayer (SAM) surfaces. Sum frequency generation (SFG) vibrational spectroscopy and attenuated total reflection Fourier transform infrared (ATR-FTIR) spectroscopy were used to deduce the NfsB enzyme orientations, which were found to be consistent with those predicted from the MD simulations. Thermal stability analyses demonstrated that NfsB variants immobilized through two tethering sites exhibited generally improved thermal stability compared with enzymes tethered at only one position. For example, NfsB enzyme chemically immobilized via positions 423 and 111 exhibits at least 60% stability increase compared to chemically immobilized NfsB mutant via a single site. This research develops a generally applicable and systematic approach using a combination of simulation and experimental methods to rationally select protein immobilization sites for the optimization of surface-immobilized enzyme activity and stability.

The estimation of changes in free energy upon mutation is central to the problem of protein design. Modern protein design methods have had remarkable success over a wide range of design targets, but are reaching their limits in ligand binding and enzyme design due to insufficient accuracy in mutational free energies. Alchemical free energy calculations have the potential to supplement modern design methods through more accurate molecular dynamics based prediction of free energy changes, but suffer from high computational cost. Multisite λ dynamics (MSλD) is a particularly efficient and scalable free energy method with potential to explore combinatorially large sequence spaces inaccessible with other free energy methods. This work aims to quantify the accuracy of MSλD and demonstrate its scalability. We apply MSλD to the classic problem of calculating folding free energies in T4 lysozyme, a system with a wealth of experimental measurements. Single site mutants considering 32 mutations show remarkable agreement with experiment with a Pearson correlation of 0.914 and mean unsigned error of 1.19 kcal/mol. Multisite mutants in systems with up to five concurrent mutations spanning 240 different sequences show comparable agreement with experiment. These results demonstrate the promise of MSλD in exploring large sequence spaces for protein design., (© 2018 The Protein Society.)

Numerous classes of riboswitches have been found to regulate bacterial gene expression in response to physiological cues, offering new paths to antibacterial drugs. As common studies of isolated riboswitches lack the functional context of the transcription machinery, we here combine single-molecule, biochemical, and simulation approaches to investigate the coupling between co-transcriptional folding of the pseudoknot-structured preQ 1 riboswitch and RNA polymerase (RNAP) pausing. We show that pausing at a site immediately downstream of the riboswitch requires a ligand-free pseudoknot in the nascent RNA, a precisely spaced sequence resembling the pause consensus, and electrostatic and steric interactions with the RNAP exit channel. While interactions with RNAP stabilize the native fold of the riboswitch, binding of the ligand signals RNAP release from the pause. Our results demonstrate that the nascent riboswitch and its ligand actively modulate the function of RNAP and vice versa, a paradigm likely to apply to other cellular RNA transcripts., (Copyright © 2018 Elsevier Inc. All rights reserved.)