The development of reliable ways of predicting the binding free energies of covalent inhibitors is a challenge for computer-aided drug design. Such development is important, for example, in the fight against the SARS-CoV-2 virus, in which covalent inhibitors can provide a promising tool for blocking Mpro, the main protease of the virus. This work develops a reliable and practical protocol for evaluating the binding free energy of covalent inhibitors. Our protocol presents a major advance over other approaches that do not consider the chemical contribution of the binding free energy. Our strategy combines the empirical valence bond method for evaluating the reaction energy profile and the PDLD/S-LRA/β method for evaluating the noncovalent part of the binding process. This protocol has been used in the calculations of the binding free energy of an α-ketoamide inhibitor of Mpro. Encouragingly, our approach reproduces the observed binding free energy. Our study of covalent inhibitors of cysteine proteases indicates that in the choice of an effective warhead it is crucial to focus on the exothermicity of the point on the free energy surface of a peptide cleavage that connects the acylation and deacylation steps. Overall, we believe that our approach should provide a powerful and effective method for in silico design of covalent drugs.

Understanding the reaction mechanism of Cas9 is crucial for the application of programmable gene editing. Despite the availability of the structures of Cas9 in apoand substrate-bound forms, the catalytically active structure is still unclear. Our first attempt to explore the catalytic mechanism of Cas9 HNH domain has been based on the reasonable assumption that we are dealing with the same mechanism as endonuclease VII, including the assumption that the catalytic water is in the first shell of the Mg 2+ . Trying this mechanism with the cryo-EM structure forced us to induce significant structural change driven by the movement of K848 (or other positively charged residue) close to the active site to facilitate the proton transfer step. In the present study, we explore a second reaction mechanism where the catalytic water is in the second shell of the Mg 2+ and assume that the cryo-EM structure by itself is a suitable representation of a catalytic-ready structure. The alternative mechanism indicates that if the active water is from the second shell then the calculate reaction barrier is This article is protected by copyright. All rights reserved. This is the author manuscript accepted for publication and has undergone full peer review but has not been through the copyediting, typesetting, pagination and proofreading process, which may lead to differences between this version and the Version of Record. Please cite this article as

The proteasome is a key protease in the eukaryotic cells which is responsible for various important cellular processes such as the control of the cell cycle, immune responses, protein homeostasis, inflammation, apoptosis, and the response to proteotoxic stress. Acting as a major molecular machine for protein degradation, proteasome first identifies damaged or obsolete regulatory proteins by attaching ubiquitin chains and subsequently utilizes conserved pore loops of the heterohexameric ring of AAA+ (ATPases associated with diverse cellular activities) to pull and mechanically unfold and translocate the misfolded protein to the active site for proteolysis. A detailed knowledge of the reaction mechanism for this proteasomal proteolysis is of central importance, both for fundamental understanding and for drug discovery. The present study investigates the mechanism of the proteolysis by the proteasome with full consideration of the protein’s flexibility and its impact on the reaction free energy. Major attention is paid to the role of the protein electrostatics in determining the activation barriers. The reaction mechanism is studied by considering a small artificial fluorogenic peptide substrate (Suc-LLVY-AMC) and evaluating the activation barriers and reaction free energies for the acylation and deacylation steps, by using the empirical valence bond method. Our results shed light on the proteolysis mechanism and thus should be important for further studies of the proteasome action.

Many cellular processes are controlled by GTPases, and gaining quantitative understanding of the activation of such processes has been a major challenge. In particular, it is crucial to obtain reliable free-energy surfaces for the relevant reaction paths both in solution and in GTPases active sites. Here, we revisit the energetics of the activation of EF-G and EF-Tu by the ribosome and explore the nature of the catalysis of the GTPase reaction. The comparison of EF-Tu to EF-G allows us to explore the impact of possible problems with the available structure of EF-Tu. Additionally, mutational effects are used for a careful validation of the emerging conclusions. It is found that the reaction may proceed by both a two-water mechanism and a one-water (GTP as a base) mechanism. However, in both cases, the activation involves a structural allosteric effect, which is likely to be a general-activation mechanism for all GTPases.

Computer-aided enzyme design is a field of great potential importance for biotechnological applications, medical advances, and a fundamental understanding of enzyme action. However, reaching a predictive ability in this direction is extremely challenging. It requires both the ability to predict quantitatively the activation barriers in cases where the structure and sequence are known and the ability to predict the effect of different mutations. In this work, we propose a protocol for predicting reasonable starting structures of mutants of proteins with known structures and for calculating the activation barriers of the generated mutants. Our approach also allows us to use the predicted structures of the generated mutant to predict structures and activation barriers for subsequent set of mutations. This protocol is used to examine the reliability of the in silico directed evolution of Kemp eliminase and haloalkane dehalogenase. We also used the results of single and double mutations as a base for predicting the effect of transition-state stabilization by multiple concurrent mutations. This strategy seems to be useful in creating an activity funnel that provides a qualitative ranking of the catalytic power of different mutants.

Increasing the accuracy of the evaluation of ligand-binding energies is one of the most important tasks of current computational biology. Here we explore the accuracy of free energy perturbation (FEP) approaches by comparing the performance of a "regular" FEP method to the one using replica exchange to enhance the sampling on a well-defined benchmark. The examination was limited to the so-called alchemical perturbations which are restricted to a fragment of the drug, and therefore, the calculation is a relative one rather than the absolute binding energy of the drug. Overall, our calculations reach the 1 kcal/mol accuracy limit. It is also shown that the accurate prediction of the position of water molecules around the binding pocket is important for FEP calculations. Interestingly, the replica exchange method does not significantly improve the accuracy of binding energies, suggesting that we reach the limit where the force field quality is a critical factor for accurate calculations.

Understanding the activation mechanism of the μ-opioid receptor (μ-OR) and its selective coupling to the inhibitory G protein (Gi) is vital for pharmaceutical research aimed at finding treatments for the opioid overdose crisis. Many attempts have been made to understand the mechanism of the μ-OR activation, following the elucidation of new crystal structures such as the antagonist- and agonist-bound μ-OR. However, the focus has not been placed on the underlying energetics and specificity of the activation process. An energy-based picture would not only help to explain this coupling but also help to explore why other possible options are not common. For example, one would like to understand why μ-OR is more selective to Gi than a stimulatory G protein (Gs). Our study used homology modeling and a coarse-grained model to generate all of the possible “end states” of the thermodynamic cycle of the activation of μ-OR. The end points were further used to generate reasonable intermediate structures of the receptor and the Gi to calculate two-dimensional free energy landscapes. The results of the landscape calculations helped to propose a plausible sequence of conformational changes in the μ-OR and Gi system and for exploring the path that leads to its activation. Furthermore, in silico alanine scanning calculations of the last 21 residues of the C terminals of Gi and Gs were performed to shed light on the selective binding of Gi to μ-OR. Overall, the present work appears to demonstrate the potential of multiscale modeling in exploring the action of G protein-coupled receptors.

Polyaromatic dye molecules employed in photovoltaic and electronic applications are often processed in organic solvents. The aggregation of these dyes is key to their applications, but a fundamental molecular understanding of how the solvent environment controls the stacking of polyaromatics is unclear. This study reports initial results from Monte Carlo simulations of how various acene molecule dimers stack when they are dissolved in different solvents. Free energies computed using full dispersion interactions versus those with sterics only suggest that solvent entropy alone accounts for the majority of the stacking free energy in solvents with compact molecular geometries such as carbon tetrachloride. However, in contrast with carbon tetrachloride, we also observe significant variations in the stacking free energies of naphthalene, anthracene, and tetracene across other solvents such as toluene and cyclohexane. The weak attractive dispersion interactions between the acene solutes and planar and near-planar solvent molecules enable them to intercalate between the acene monomers, inducing extra stability beyond what solvent entropic driving force alone could predict. In all three solvents studied (carbon tetrachloride, cyclohexane, toluene) the solvent environment helps facilitate stacking of all three acenes studied (naphthalene, anthracene, tetracene), inducing a significant stabilization free energy between −4 and −8 kcal/mol. Extensive free energy umbrella sampling along the other orthogonal directions allows us to accurately calculate the dimerization equilibrium constants of all three acenes, which vary over several orders of magnitude in a way that depends intricately on the solvent they are in. Given the prevalence of solution-based processing techniques for organic electronic and photonic devices, these results provide useful insights into the critical role that solvent structure and characteristics play in the solution-based aggregation of organic dyes.