Accurate Binding Free Energy Predictions in Fragment Optimization

Peptide Drug Design Using Alchemical Free Energy Calculation: An Application and Validation on Agonists of Ghrelin Receptor

With recent large-scale applications and validations, the relative binding free energy (RBFE) calculated using alchemical free energy methods has been proven to be an accurate measure to probe the binding of small-molecule drug candidates. On the other hand, given the flexibility of peptides, it is of great interest to find out whether sufficient sampling could be achieved within the typical time scale of such calculation, and a similar level of accuracy could be reached for peptide drugs. However, the systematic evaluation of such calculations on protein-peptide systems has been less reported. Most reported studies of peptides were restricted to a limited number of data points or lacking experimental support. To demonstrate the applicability of the alchemical free energy method for protein-peptide systems in a typical real-world drug discovery project, we report an application of the thermodynamic integration (TI) method to the RBFE calculation of ghrelin receptor and its peptide agonists. Along with the calculation, the synthesis and in vitro EC50 activity of relamorelin and 17 new peptide derivatives were also reported. A cost-effective criterion to determine the data collection time was proposed for peptides in the TI simulation. The average of three TI repeats yielded a mean absolute error of 0.98 kcal/mol and Pearson's correlation coefficient (R) of 0.77 against the experimental free energy derived from the in vitro EC50 activity, showing good repeatability of the proposed method and a slightly better agreement than the results obtained from the arbitrary time frames up to 20 ns. Although it is limited by having one target and a deduced binding pose, we hope that this study can add some insights into alchemical free energy calculation of protein-peptide systems, providing theoretical assistance to the development of peptide drugs.

Explorations of Agonist Selectivity for the α9* nAChR with Novel Substituted Carbamoyl/Amido/Heteroaryl Dialkylpiperazinium Salts and Their Therapeutic Implications in Pain and Inflammation

There is an urgent need for nonopioid treatments for chronic and neuropathic pain to provide effective alternatives amid the escalating opioid crisis. This study introduces novel compounds targeting the α9 nicotinic acetylcholine receptor (nAChR) subunit, which is crucial for pain regulation, inflammation, and inner ear functions. Specifically, it identifies novel substituted carbamoyl/amido/heteroaryl dialkylpiperazinium iodides as potent agonists selective for human α9 and α9α10 over α7 nAChRs, particularly compounds 3f, 3h, and 3j. Compound 3h (GAT2711) demonstrated a 230 nM potency as a full agonist at α9 nAChRs, being 340-fold selective over α7. Compound 3c was 10-fold selective for α9α10 over α9 nAChR. Compounds 2, 3f, and 3h inhibited ATP-induced interleukin-1β release in THP-1 cells. The analgesic activity of 3h was fully retained in α7 knockout mice, suggesting that analgesic effects were potentially mediated through α9* nAChRs. Our findings provide a blueprint for developing α9*-specific therapeutics for pain.

Toward physics-based precision medicine: Exploiting protein dynamics to design new therapeutics and interpret variants

The goal of precision medicine is to utilize our knowledge of the molecular causes of disease to better diagnose and treat patients. However, there is a substantial mismatch between the small number of food and drug administration (FDA)-approved drugs and annotated coding variants compared to the needs of precision medicine. This review introduces the concept of physics-based precision medicine, a scalable framework that promises to improve our understanding of sequence-function relationships and accelerate drug discovery. We show that accounting for the ensemble of structures a protein adopts in solution with computer simulations overcomes many of the limitations imposed by assuming a single protein structure. We highlight studies of protein dynamics and recent methods for the analysis of structural ensembles. These studies demonstrate that differences in conformational distributions predict functional differences within protein families and between variants. Thanks to new computational tools that are providing unprecedented access to protein structural ensembles, this insight may enable accurate predictions of variant pathogenicity for entire libraries of variants. We further show that explicitly accounting for protein ensembles, with methods like alchemical free energy calculations or docking to Markov state models, can uncover novel lead compounds. To conclude, we demonstrate that cryptic pockets, or cavities absent in experimental structures, provide an avenue to target proteins that are currently considered undruggable. Taken together, our review provides a roadmap for the field of protein science to accelerate precision medicine.

Accurate calculation of affinity changes to the close state of influenza A M2 transmembrane domain in response to subtle structural changes of adamantyl amines using free energy perturbation methods in different lipid bilayers

Experimental binding free energies of 27 adamantyl amines against the influenza M2(22-46) WT tetramer, in its closed form at pH 8, were measured by ITC in DPC micelles. The measured Kd's range is ~44 while the antiviral potencies (IC50) range is ~750 with a good correlation between binding free energies computed with Kd and IC50 values (r = 0.76). We explored with MD simulations (ff19sb, CHARMM36m) the binding profile of complexes with strong, moderate and weak binders embedded in DMPC, DPPC, POPC or a viral mimetic membrane and using different experimental starting structures of M2. To predict accurately differences in binding free energy in response to subtle changes in the structure of the ligands, we performed 18 alchemical perturbative single topology FEP/MD NPT simulations (OPLS2005) using the BAR estimator (Desmond software) and 20 dual topology calculations TI/MD NVT simulations (ff19sb) using the MBAR estimator (Amber software) for adamantyl amines in complex with M2(22-46) WT in DMPC, DPPC, POPC. We observed that both methods with all lipids show a very good correlation between the experimental and calculated relative binding free energies (r = 0.77-0.87, mue = 0.36-0.92 kcal mol-1) with the highest performance achieved with TI/MBAR and lowest performance with FEP/BAR in DMPC bilayers. When antiviral potencies are used instead of the Kd values for computing the experimental binding free energies we obtained also good performance with both FEP/BAR (r = 0.83, mue = 0.75 kcal mol-1) and TI/MBAR (r = 0.69, mue = 0.77 kcal mol-1).

Decision tree-based identification of important molecular fragments for protein-ligand binding

Fragment-based drug design is an emerging technology in pharmaceutical research and development. One of the key aspects of this technology is the identification and quantitative characterization of molecular fragments. This study presents a strategy for identifying important molecular fragments based on molecular fingerprints and decision tree algorithms and verifies its feasibility in predicting protein-ligand binding affinity. Specifically, the three-dimensional (3D) structures of protein-ligand complexes are encoded using extended-connectivity fingerprints (ECFP), and three decision tree models, namely Random Forest, XGBoost, and LightGBM, are used to quantitatively characterize the feature importance, thereby extracting important molecular fragments with high reliability. Few-shot learning reveals that the extracted molecular fragments contribute significantly and consistently to the binding affinity even with a small sample size. Despite the absence of location and distance information for molecular fragments in ECFP, 3D visualization, in combination with the reverse ECFP process, shows that the majority of the extracted fragments are located at the binding interface of the protein and the ligand. This alignment with the distance constraints critical for binding affinity further supports the reliability of the strategy for identifying important molecular fragments.

Fast free energy estimates from λ-dynamics with bias-updated Gibbs sampling

Relative binding free energy calculations have become an integral computational tool for lead optimization in structure-based drug design. Classical alchemical methods, including free energy perturbation or thermodynamic integration, compute relative free energy differences by transforming one molecule into another. However, these methods have high operational costs due to the need to perform many pairwise perturbations independently. To reduce costs and accelerate molecular design workflows, we present a method called λ-dynamics with bias-updated Gibbs sampling. This method uses dynamic biases to continuously sample between multiple ligand analogues collectively within a single simulation. We show that many relative binding free energies can be determined quickly with this approach without compromising accuracy. For five benchmark systems, agreement to experiment is high, with root mean square errors near or below 1.0 kcal mol-1. Free energy results are consistent with other computational approaches and within statistical noise of both methods (0.4 kcal mol-1 or less). Notably, large efficiency gains over thermodynamic integration of 18-66-fold for small perturbations and 100-200-fold for whole aromatic ring substitutions are observed. The rapid determination of relative binding free energies will enable larger chemical spaces to be more readily explored and structure-based drug design to be accelerated.

Surface Engineering of Migratory Corrosion Inhibitors: Controlling the Wettability of Calcium Silicate Hydrate in the Nanoscale

Migratory corrosion inhibitors (MCIs) are regarded as effective additives to prevent harmful ion transmission and improve concrete durability, but their behavior in the porosity of concrete is still unclarified. This paper proposes a unique perspective to evaluate the effects of surfactant-like MCIs in calcium silicate hydrate (C-S-H) nanoporosity through molecular and electronic structural information. Advanced enhanced sampling methods and perturbation theory methods were applied to evaluate the role of different MCIs. The reduced density gradient of MCI molecules was obtained by using quantum chemical calculations. This calculation is instrumental in elucidating the intensity of interactions among distinct MCI molecule head groups and the C-S-H matrix. It is found that MCIs can effectively improve the interfacial tension (IFT) between C-S-H and water, which corresponds to the inhibitory ability of transmission. Free energy indicates that the MCI has the properties of strong adsorption and weak dissolution, facilitating the improvement of IFT. The relationship between the MCI functional group and the ability of adsorption and dissolution is revealed. This study suggests that MCIs work as surface controllers of C-S-H pores and that their properties can be assessed on the nanoscale.

The maximal and current accuracy of rigorous protein-ligand binding free energy calculations

Computational techniques can speed up the identification of hits and accelerate the development of candidate molecules for drug discovery. Among techniques for predicting relative binding affinities, the most consistently accurate is free energy perturbation (FEP), a class of rigorous physics-based methods. However, uncertainty remains about how accurate FEP is and can ever be. Here, we present what we believe to be the largest publicly available dataset of proteins and congeneric series of small molecules, and assess the accuracy of the leading FEP workflow. To ascertain the limit of achievable accuracy, we also survey the reproducibility of experimental relative affinity measurements. We find a wide variability in experimental accuracy and a correspondence between binding and functional assays. When careful preparation of protein and ligand structures is undertaken, FEP can achieve accuracy comparable to experimental reproducibility. Throughout, we highlight reliable protocols that can help maximize the accuracy of FEP in prospective studies.

FEP Protocol Builder: Optimization of Free Energy Perturbation Protocols Using Active Learning

Significant improvements have been made in the past decade to methods that rapidly and accurately predict binding affinity through free energy perturbation (FEP) calculations. This has been driven by recent advances in small-molecule force fields and sampling algorithms combined with the availability of low-cost parallel computing. Predictive accuracies of ∼1 kcal mol-1 have been regularly achieved, which are sufficient to drive potency optimization in modern drug discovery campaigns. Despite the robustness of these FEP approaches across multiple target classes, there are invariably target systems that do not display expected performance with default FEP settings. Traditionally, these systems required labor-intensive manual protocol development to arrive at parameter settings that produce a predictive FEP model. Due to the (a) relatively large parameter space to be explored, (b) significant compute requirements, and (c) limited understanding of how combinations of parameters can affect FEP performance, manual FEP protocol optimization can take weeks to months to complete, and often does not involve rigorous train-test set splits, resulting in potential overfitting. These manual FEP protocol development timelines do not coincide with tight drug discovery project timelines, essentially preventing the use of FEP calculations for these target systems. Here, we describe an automated workflow termed FEP Protocol Builder (FEP-PB) to rapidly generate accurate FEP protocols for systems that do not perform well with default settings. FEP-PB uses an active-learning workflow to iteratively search the protocol parameter space to develop accurate FEP protocols. To validate this approach, we applied it to pharmaceutically relevant systems where default FEP settings could not produce predictive models. We demonstrate that FEP-PB can rapidly generate accurate FEP protocols for the previously challenging MCL1 system with limited human intervention. We also apply FEP-PB in a real-world drug discovery setting to generate an accurate FEP protocol for the p97 system. FEP-PB is able to generate a more accurate protocol than the expert user, rapidly validating p97 as amenable to free energy calculations. Additionally, through the active-learning workflow, we are able to gain insight into which parameters are most important for a given system. These results suggest that FEP-PB is a robust tool that can aid in rapidly developing accurate FEP protocols and increasing the number of targets that are amenable to the technology.

New Insights into the Cooperativity and Dynamics of Dimeric Enzymes

A survey of protein databases indicates that the majority of enzymes exist in oligomeric forms, with about half of those found in the UniProt database being homodimeric. Understanding why many enzymes are in their dimeric form is imperative. Recent developments in experimental and computational techniques have allowed for a deeper comprehension of the cooperative interactions between the subunits of dimeric enzymes. This review aims to succinctly summarize these recent advancements by providing an overview of experimental and theoretical methods, as well as an understanding of cooperativity in substrate binding and the molecular mechanisms of cooperative catalysis within homodimeric enzymes. Focus is set upon the beneficial effects of dimerization and cooperative catalysis. These advancements not only provide essential case studies and theoretical support for comprehending dimeric enzyme catalysis but also serve as a foundation for designing highly efficient catalysts, such as dimeric organic catalysts. Moreover, these developments have significant implications for drug design, as exemplified by Paxlovid, which was designed for the homodimeric main protease of SARS-CoV-2.

A mutation-induced drug resistance database (MdrDB)

Mutation-induced drug resistance is a significant challenge to the clinical treatment of many diseases, as structural changes in proteins can diminish drug efficacy. Understanding how mutations affect protein-ligand binding affinities is crucial for developing new drugs and therapies. However, the lack of a large-scale and high-quality database has hindered the research progresses in this area. To address this issue, we have developed MdrDB, a database that integrates data from seven publicly available datasets, which is the largest database of its kind. By integrating information on drug sensitivity and cell line mutations from Genomics of Drug Sensitivity in Cancer and DepMap, MdrDB has substantially expanded the existing drug resistance data. MdrDB is comprised of 100,537 samples of 240 proteins (which encompass 5119 total PDB structures), 2503 mutations, and 440 drugs. Each sample brings together 3D structures of wild type and mutant protein-ligand complexes, binding affinity changes upon mutation (ΔΔG), and biochemical features. Experimental results with MdrDB demonstrate its effectiveness in significantly enhancing the performance of commonly used machine learning models when predicting ΔΔG in three standard benchmarking scenarios. In conclusion, MdrDB is a comprehensive database that can advance the understanding of mutation-induced drug resistance, and accelerate the discovery of novel chemicals.

Deciphering antidiarrheal effects of Meda pata (Litsea glutinosa (Lour.) C.B.Rob.) leaf extract in chemical-induced models of albino rats

ETHNOPHARMACOLOGICAL RELEVANCE

Diarrhea is one of the leading causes of preventable death in developing countries, mainly caused by bacterial infections and traditional therapies are very common in diarrheal incidences. Meda Pata (Litsea glutionsa) has a long history of use as traditional medicine for diarrhea, dysentery, and spasm in Bangladesh, India, and some other Asian countries.

AIM OF THE STUDY

This research reports the antidiarrheal effects of Meda Pata (Litsea glutinosa leaf extract, LGLEx) in animal models. The work has been supported by in silico molecular docking study to verify the effects.

MATERIALS AND METHODS

The antidiarrheal effect of LGLEx was investigated in castor oil-induced diarrhea, magnesium sulfate-induced diarrhea, and gastrointestinal motility test models. Antidiarrheal effects were supported by a molecular docking study through an interaction between LGLEx's GC-MS analyzed imidazole-containing compounds and muscarinic acetylcholine receptor (PDB: 4U14) and 5-HT3 receptor (PDB: 5AIN).

RESULTS

LGLEx potentially reduced the diarrheal incidences in in vivo assays reducing gastrointestinal motility. The maximum diarrheal inhibition was obtained in the castor oil-induced model (62.63%) and and BaSO₄-induced model (73.14%), which were statistically significant (P < 0.05) when compared to the reference drug loperamide. In the castor-oil and BaSO4-induced models, peristaltic movement was reduced by 25.96% and 32.17%, respectively. Biochemical markers particularly IgE, C-reactive proteins, and serum electrolytes were significantly (P < 0.0) restored in treated groups. A Molecular docking analysis revealed that two compounds (1-Ethyl-2-hydroxymethylimidazole and 1,6-Anhydro-beta-D-glucofuranose demonstrated the highest binding affinity with target receptors muscarinic acetylcholine receptor (PDB: 4U14) and 5-HT3 receptor (PDB: 5AIN) confirming their drug likeliness. The findings indicate a high potential antidiarrheal impact that warrants further investigation for its therapeutic application.

Combining QM/MM Calculations with Classical Mining Minima to Predict Protein-Ligand Binding Free Energy

We developed an effective binding free energy prediction protocol which incorporates quantum mechanical/molecular mechanical (QM/MM) calculations to substitute the specified atomic charges of force fields with quantum-mechanically recalculated ones at a proposed pose using a mining minima approach with the VeraChem mining minima engine. We tested this protocol using seven well-known targets with 147 different ligands and compared it with classical mining minima and the most popular binding free energy (BFE) methods using different metrics. Our new protocol, dubbed Qcharge-VM2, yielded an overall Pearson correlation of 0.86, which was better than all the methods examined. Qcharge-VM2 performed significantly better than implicit solvent-based methods, such as MM-GBSA and MM-PBSA, but not as good as explicit water-based free energy perturbation methods, such as FEP+, in terms of root-mean-square error, RMSE (1.75 kcal/mol) and mean unsigned error, MUE (1.39 kcal/mol) on a limited set of targets. However, our protocol is substantially less computationally demanding compared with FEP+. The combined accuracy and efficiency of our method can be valuable in drug discovery campaigns.

Past, Present, and Future Perspectives on Computer-Aided Drug Design Methodologies

The application of computational approaches in drug discovery has been consolidated in the last decades. These families of techniques are usually grouped under the common name of "computer-aided drug design" (CADD), and they now constitute one of the pillars in the pharmaceutical discovery pipelines in many academic and industrial environments. Their implementation has been demonstrated to tremendously improve the speed of the early discovery steps, allowing for the proficient and rational choice of proper compounds for a desired therapeutic need among the extreme vastness of the drug-like chemical space. Moreover, the application of CADD approaches allows the rationalization of biochemical and interactive processes of pharmaceutical interest at the molecular level. Because of this, computational tools are now extensively used also in the field of rational 3D design and optimization of chemical entities starting from the structural information of the targets, which can be experimentally resolved or can also be obtained with other computer-based techniques. In this work, we revised the state-of-the-art computer-aided drug design methods, focusing on their application in different scenarios of pharmaceutical and biological interest, not only highlighting their great potential and their benefits, but also discussing their actual limitations and eventual weaknesses. This work can be considered a brief overview of computational methods for drug discovery.

SophosQM: Accurate Binding Affinity Prediction in Compound Optimization

The optimization of compounds' binding affinity for a biological target is a crucial aspect of the drug development process. Being able to accurately predict binding energies in advance of synthesizing compounds would have a massive impact on the speed of the drug discovery process. The ideal binding affinity prediction method should combine accuracy, reliability, and speed. In this paper, we present SophosQM, a quantum mechanics (QM)-based approach, which can accurately predict the binding affinities of compounds to proteins. The binding affinity predictive models generated by SophosQM are based on the fragment molecular orbital (FMO) method to estimate the enthalpic component of the binding free energy, and a macroscopic descriptor, clog P, is used as an approximation of the entropic component. The affinity prediction is performed using multilinear regression, fitting the experimental values against the FMO-computed enthalpic term and clog P. The quality of the prediction can be assessed in terms of the correlation coefficient between experimental and predicted values. In this work, the method's reliability and accuracy are exemplified by applying SophosQM to 70 compounds binding to six different targets of pharmaceutical relevance. Overall, the results show a very satisfactory performance with a global correlation coefficient in the order of 0.9. Our predictions also show a satisfactory performance compared to data based on free energy perturbation. Finally, SophosQM can also be applied in high-throughput mode by using semiempirical QM methods to evaluate large portions of chemical space, while retaining a good level of accuracy, but decreasing the computing time to just a few seconds per compound.

Ligand binding free energy evaluation by Monte Carlo Recursion

The correct evaluation of ligand binding free energies by computational methods is still a very challenging active area of research. The most employed methods for these calculations can be roughly classified into four groups: (i) the fastest and less accurate methods, such as molecular docking, designed to sample a large number of molecules and rapidly rank them according to the potential binding energy; (ii) the second class of methods use a thermodynamic ensemble, typically generated by molecular dynamics, to analyze the endpoints of the thermodynamic cycle for binding and extract differences, in the so-called 'end-point' methods; (iii) the third class of methods is based on the Zwanzig relationship and computes the free energy difference after a chemical change of the system (alchemical methods); and (iv) methods based on biased simulations, such as metadynamics, for example. These methods require increased computational power and as expected, result in increased accuracy for the determination of the strength of binding. Here, we describe an intermediate approach, based on the Monte Carlo Recursion (MCR) method first developed by Harold Scheraga. In this method, the system is sampled at increasing effective temperatures, and the free energy of the system is assessed from a series of terms W(b,T), computed from Monte Carlo (MC) averages at each iteration. We show the application of the MCR for ligand binding with datasets of guest-hosts systems (N = 75) and we observed that a good correlation is obtained between experimental data and the binding energies computed with MCR. We also compared the experimental data with an end-point calculation from equilibrium Monte Carlo calculations that allowed us to conclude that the lower-energy (lower-temperature) terms in the calculation are the most relevant to the estimation of the binding energies, resulting in similar correlations between MCR and MC data and the experimental values. On the other hand, the MCR method provides a reasonable view of the binding energy funnel, with possible connections with the ligand binding kinetics, as well. The codes developed for this analysis are publicly available on GitHub as a part of the LiBELa/MCLiBELa project (https://github.com/alessandronascimento/LiBELa).

QSAR via Multisite λ-Dynamics in the Orphaned TSSK1B Kinase

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. This article is protected by copyright. All rights reserved.

Recent Advances in Alchemical Binding Free Energy Calculations for Drug Discovery

Rigorous physics-based methods to calculate binding free energies of protein-ligand complexes have become a valued component of structure-based drug design. Relative and absolute binding free energy calculations have been deployed prospectively in support of solving diverse drug discovery challenges. Here we review recent applications of binding free energy calculations to fragment growing and linking, scaffold hopping, binding pose validation, virtual screening, covalent enzyme inhibition, and positional analogue scanning. Furthermore, we discuss the merits of using protein models and highlight recent efforts to replace costly binding free energy calculations with predictions from machine learning models trained on a limited number of free energy perturbation or thermodynamic integration calculations thereby allowing for extended chemical space exploration.

Thirty years of molecular dynamics simulations on posttranslational modifications of proteins

Posttranslational modifications (PTMs) are an integral component to how cells respond to perturbation. While experimental advances have enabled improved PTM identification capabilities, the same throughput for characterizing how structural changes caused by PTMs equate to altered physiological function has not been maintained. In this Perspective, we cover the history of computational modeling and molecular dynamics simulations which have characterized the structural implications of PTMs. We distinguish results from different molecular dynamics studies based upon the timescales simulated and analysis approaches used for PTM characterization. Lastly, we offer insights into how opportunities for modern research efforts on in silico PTM characterization may proceed given current state-of-the-art computing capabilities and methodological advancements.

Application of Computational Biology and Artificial Intelligence in Drug Design

Traditional drug design requires a great amount of research time and developmental expense. Booming computational approaches, including computational biology, computer-aided drug design, and artificial intelligence, have the potential to expedite the efficiency of drug discovery by minimizing the time and financial cost. In recent years, computational approaches are being widely used to improve the efficacy and effectiveness of drug discovery and pipeline, leading to the approval of plenty of new drugs for marketing. The present review emphasizes on the applications of these indispensable computational approaches in aiding target identification, lead discovery, and lead optimization. Some challenges of using these approaches for drug design are also discussed. Moreover, we propose a methodology for integrating various computational techniques into new drug discovery and design.

Predicting antibiotic resistance in complex protein targets using alchemical free energy methods

Drug resistant Mycobacterium tuberculosis, which mostly results from single nucleotide polymorphisms in antibiotic target genes, poses a major threat to tuberculosis treatment outcomes. Relative binding free energy (RBFE) calculations can rapidly predict the effects of mutations, but this approach has not been tested on large, complex proteins. We use RBFE calculations to predict the effects of M. tuberculosis RNA polymerase and DNA gyrase mutations on rifampicin and moxifloxacin susceptibility respectively. These mutations encompass a range of amino acid substitutions with known effects and include large steric perturbations and charged moieties. We find that moderate numbers (n = 3-15) of short RBFE calculations can predict resistance in cases where the mutation results in a large change in the binding free energy. We show that the method lacks discrimination in cases with either a small change in energy or that involve charged amino acids, and we investigate how these calculation errors may be decreased.

Are Deep Learning Structural Models Sufficiently Accurate for Free-Energy Calculations? Application of FEP+ to AlphaFold2-Predicted Structures

The availability of AlphaFold2 has led to great excitement in the scientific community─particularly among drug hunters─due to the ability of the algorithm to predict protein structures with high accuracy. However, beyond globally accurate protein structure prediction, it remains to be determined whether ligand binding sites are predicted with sufficient accuracy in these structures to be useful in supporting computationally driven drug discovery programs. We explored this question by performing free-energy perturbation (FEP) calculations on a set of well-studied protein-ligand complexes, where AlphaFold2 predictions were performed by removing all templates with >30% identity to the target protein from the training set. We observed that in most cases, the ΔΔG values for ligand transformations calculated with FEP, using these prospective AlphaFold2 structures, were comparable in accuracy to the corresponding calculations previously carried out using crystal structures. We conclude that under the right circumstances, AlphaFold2-modeled structures are accurate enough to be used by physics-based methods such as FEP in typical lead optimization stages of a drug discovery program.

Evaluating the use of absolute binding free energy in the fragment optimisation process

Key to the fragment optimisation process within drug design is the need to accurately capture the changes in affinity that are associated with a given set of chemical modifications. Due to the weakly binding nature of fragments, this has proven to be a challenging task, despite recent advancements in leveraging experimental and computational methods. In this work, we evaluate the use of Absolute Binding Free Energy (ABFE) calculations in guiding fragment optimisation decisions, retrospectively calculating binding free energies for 59 ligands across 4 fragment elaboration campaigns. We first demonstrate that ABFEs can be used to accurately rank fragment-sized binders with an overall Spearman's r of 0.89 and a Kendall τ of 0.67, although often deviating from experiment in absolute free energy values with an RMSE of 2.75 kcal/mol. We then also show that in several cases, retrospective fragment optimisation decisions can be supported by the ABFE calculations. Comparing against cheaper endpoint methods, namely N_wat-MM/GBSA, we find that ABFEs offer better ranking power and correlation metrics. Our results indicate that ABFE calculations can usefully guide fragment elaborations to maximise affinity.

Drug Discovery by Automated Adaptation of Chemical Structure and Identity

Computer-aided drug design offers the potential to dramatically reduce the cost and effort required for drug discovery. While screening-based methods are valuable in the early stages of hit identification, they are frequently succeeded by iterative, hypothesis-driven computations that require recurrent investment of human time and intuition. To increase automation, we introduce a computational method for lead refinement that combines concerted dynamics of the ligand/protein complex via molecular dynamics simulations with integrated Monte Carlo-based changes in the chemical formula of the ligand. This approach, which we refer to as ligand-exchange Monte Carlo molecular dynamics, accounts for solvent- and entropy-based contributions to competitive binding free energies by coupling the energetics of bound and unbound states during the ligand-exchange attempt. Quantitative comparison of relative binding free energies to reference values from free energy perturbation, conducted in vacuum, indicates that ligand-exchange Monte Carlo molecular dynamics simulations sample relevant conformational ensembles and are capable of identifying strongly binding compounds. Additional simulations demonstrate the use of an implicit solvent model. We speculate that the use of chemical graphs in which exchanges are only permitted between ligands with sufficient similarity may enable an automated search to capture some of the benefits provided by human intuition during hypothesis-guided lead refinement.

Exploration of novel ligands to target C-C Motif Chemokine Receptor 2 (CCR2) as a promising pharmacological treatment against traumatic brain injury

It is widely reported that the overexpression of the C-C Motif Chemokine Receptor 2 (CCR2) has negative implications in neuroinflammatory diseases such as traumatic brain injury (TBI), although promising drugs to tackle this have been less forthcoming. As of 2016, only 2 drugs specifically targeting this receptor have made their way to market, with unsuccessful outcome unfortunately, suggesting that the search for more specific and precise ligands is utterly necessary. In this paper we hypothesized that by targeting Glu291, Met295, Trp98, Leu45 and Val189 amino acids, essential in the binding of CCR2 with C-C Motif Chemokine Ligand 2 (CCL2), the endogenous substrate, mitigates its activity in TBI. We used a pharmacophore model to screen for suitable ligands that may bind to CCR2, which returned 871 ligands. Docking and molecular dynamics results uncovered that two ligands (A102) and (A435) contained several of those important residues and showed a stability and compactness when in complex with CCR2, with these results confirmed by MMGBSA calculations with A102 recording a better interaction compared to A435. Finally, a PPI network was built to explore downstream signaling being regulated by both ligands in TBI, showing amyloid precursor protein (APP) as a key target and neuroactive-ligand receptor interaction (1.80E-27) the top functional annotated category. In conclusion, for the first time we report novel ligands A102 and A435 targeting CCR2 as a potential new pharmacological approach to target inflammation post-TBI.

Alchemical Free Energy Estimators and Molecular Dynamics Engines: Accuracy, Precision, and Reproducibility

The binding free energy between a ligand and its target protein is an essential quantity to know at all stages of the drug discovery pipeline. Assessing this value computationally can offer insight into where efforts should be focused in the pursuit of effective therapeutics to treat a myriad of diseases. In this work, we examine the computation of alchemical relative binding free energies with an eye for assessing reproducibility across popular molecular dynamics packages and free energy estimators. The focus of this work is on 54 ligand transformations from a diverse set of protein targets: MCL1, PTP1B, TYK2, CDK2, and thrombin. These targets are studied with three popular molecular dynamics packages: OpenMM, NAMD2, and NAMD3 alpha. Trajectories collected with these packages are used to compare relative binding free energies calculated with thermodynamic integration and free energy perturbation methods. The resulting binding free energies show good agreement between molecular dynamics packages with an average mean unsigned error between them of 0.50 kcal/mol. The correlation between packages is very good, with the lowest Spearman's, Pearson's and Kendall's tau correlation coefficients being 0.92, 0.91, and 0.76, respectively. Agreement between thermodynamic integration and free energy perturbation is shown to be very good when using ensemble averaging.

In Silico Structure-Based Approach for Group Efficiency Estimation in Fragment-Based Drug Design Using Evaluation of Fragment Contributions

The notion of a contribution of a specific group in an organic molecule’s property and/or activity is both common in our thinking and is still not strictly correct due to the inherent non-additivity of free energy with respect to molecular fragments composing a molecule. The fragment- based drug discovery (FBDD) approach has proven to be fruitful in addressing the above notions. The main difficulty of the FBDD, however, is in its reliance on the low throughput and expensive experimental means of determining the fragment-sized molecules binding. In this article we propose a way to enhance the throughput and availability of the FBDD methods by judiciously using an in silico means of assessing the contribution to ligand-receptor binding energy of fragments of a molecule under question using a previously developed in silico Reverse Fragment Based Drug Discovery (R-FBDD) approach. It has been shown that the proposed structure-based drug discovery (SBDD) type of approach fills in the vacant niche among the existing in silico approaches, which mainly stem from the ligand-based drug discovery (LBDD) counterparts. In order to illustrate the applicability of the approach, our work retrospectively repeats the findings of the use case of an FBDD hit-to-lead project devoted to the experimentally based determination of additive group efficiency (GE)—an analog of ligand efficiency (LE) for a group in the molecule—using the Free-Wilson (FW) decomposition. It is shown that in using our in silico approach to evaluate fragment contributions of a ligand and to estimate GE one can arrive at similar decisions as those made using the experimentally determined activity-based FW decomposition. It is also shown that the approach is rather robust to the choice of the scoring function, provided the latter demonstrates a decent scoring power. We argue that the proposed approach of in silico assessment of GE has a wider applicability domain and expect that it will be widely applicable to enhance the net throughput of drug discovery based on the FBDD paradigm.

SPLDExtraTrees: robust machine learning approach for predicting kinase inhibitor resistance

Drug resistance is a major threat to the global health and a significant concern throughout the clinical treatment of diseases and drug development. The mutation in proteins that is related to drug binding is a common cause for adaptive drug resistance. Therefore, quantitative estimations of how mutations would affect the interaction between a drug and the target protein would be of vital significance for the drug development and the clinical practice. Computational methods that rely on molecular dynamics simulations, Rosetta protocols, as well as machine learning methods have been proven to be capable of predicting ligand affinity changes upon protein mutation. However, the severely limited sample size and heavy noise induced overfitting and generalization issues have impeded wide adoption of machine learning for studying drug resistance. In this paper, we propose a robust machine learning method, termed SPLDExtraTrees, which can accurately predict ligand binding affinity changes upon protein mutation and identify resistance-causing mutations. Especially, the proposed method ranks training data following a specific scheme that starts with easy-to-learn samples and gradually incorporates harder and diverse samples into the training, and then iterates between sample weight recalculations and model updates. In addition, we calculate additional physics-based structural features to provide the machine learning model with the valuable domain knowledge on proteins for these data-limited predictive tasks. The experiments substantiate the capability of the proposed method for predicting kinase inhibitor resistance under three scenarios and achieve predictive accuracy comparable with that of molecular dynamics and Rosetta methods with much less computational costs.

Unbinding ligands from SARS-CoV-2 Mpro via umbrella sampling simulations

The umbrella sampling (US) simulation is demonstrated to be an efficient approach for determining the unbinding pathway and binding affinity to the SARS-CoV-2 Mpro of small molecule inhibitors. The accuracy of US is in the same range as the linear interaction energy (LIE) and fast pulling of ligand (FPL) methods. In detail, the correlation coefficient between US and experiments does not differ from FPL and is slightly smaller than LIE. The root mean square error of US simulations is smaller than that of LIE. Moreover, US is better than FPL and poorer than LIE in classifying SARS-CoV-2 Mpro inhibitors owing to the reciever operating characteristic-area under the curve analysis. Furthermore, the US simulations also provide detailed insights on unbinding pathways of ligands from the binding cleft of SARS-CoV-2 Mpro. The residues Cys44, Thr45, Ser46, Leu141, Asn142, Gly143, Glu166, Leu167, Pro168, Ala191, Gln192 and Ala193 probably play an important role in the ligand dissociation. Therefore, substitutions at these points may change the mechanism of binding of inhibitors to SARS-CoV-2 Mpro.

Efficient Hit-to-Lead Searching of Kinase Inhibitor Chemical Space via Computational Fragment Merging

In early-stage drug discovery, the hit-to-lead optimization (or "hit expansion") stage entails starting from a newly identified active compound and improving its potency or other properties. Traditionally, this process relies on synthesizing and evaluating a series of analogues to build up structure-activity relationships. Here, we describe a computational strategy focused on kinase inhibitors, intended to expedite the process of identifying analogues with improved potency. Our protocol begins from an inhibitor of the target kinase and generalizes the synthetic route used to access it. By searching for commercially available replacements for the individual building blocks used to make the parent inhibitor, we compile an enumerated library of compounds that can be accessed using the same chemical transformations; these huge libraries can exceed many millions─or billions─of compounds. Because the resulting libraries are much too large for explicit virtual screening, we instead consider alternate approaches to identify the top-scoring compounds. We find that contributions from individual substituents are well described by a pairwise additivity approximation, provided that the corresponding fragments position their shared core in precisely the same way relative to the binding site. This key insight allows us to determine which fragments are suitable for merging into single new compounds and which are not. Further, the use of pairwise approximation allows interaction energies to be assigned to each compound in the library without the need for any further structure-based modeling: interaction energies instead can be reliably estimated from the energies of the component fragments, and the reduced computational requirements allow for flexible energy minimizations that allow the kinase to respond to each substitution. We demonstrate this protocol using libraries built from six representative kinase inhibitors drawn from the literature, which target five different kinases: CDK9, CHK1, CDK2, EGFR^T790M, and ACK1. In each example, the enumerated library includes additional analogues reported by the original study to have activity, and these analogues are successfully prioritized within the library. We envision that the insights from this work can facilitate the rapid assembly and screening of increasingly large libraries for focused hit-to-lead optimization. To enable adoption of these methods and to encourage further analyses, we disseminate the computational tools needed to deploy this protocol.

Progress in Free Energy Perturbation: Options for Evolving Fragments

One of the remaining bottlenecks in fragment-based drug design (FBDD) is the initial exploration and optimization of the identified hit fragments. There is a growing interest in computational approaches that can guide these efforts by predicting the binding affinity of newly designed analogues. Among others, alchemical free energy (AFE) calculations promise high accuracy at a computational cost that allows their application during lead optimization campaigns. In this review, we discuss how AFE could have a strong impact in fragment evolution, and we raise awareness on the challenges that could be encountered applying this methodology in FBDD studies.

PremPLI: a machine learning model for predicting the effects of missense mutations on protein-ligand interactions

Resistance to small-molecule drugs is the main cause of the failure of therapeutic drugs in clinical practice. Missense mutations altering the binding of ligands to proteins are one of the critical mechanisms that result in genetic disease and drug resistance. Computational methods have made a lot of progress for predicting binding affinity changes and identifying resistance mutations, but their prediction accuracy and speed are still not satisfied and need to be further improved. To address these issues, we introduce a structure-based machine learning method for quantitatively estimating the effects of single mutations on ligand binding affinity changes (named as PremPLI). A comprehensive comparison of the predictive performance of PremPLI with other available methods on two benchmark datasets confirms that our approach performs robustly and presents similar or even higher predictive accuracy than the approaches relying on first-principle statistical mechanics and mixed physics- and knowledge-based potentials while requires much less computational resources. PremPLI can be used for guiding the design of ligand-binding proteins, identifying and understanding disease driver mutations, and finding potential resistance mutations for different drugs. PremPLI is freely available at https://lilab.jysw.suda.edu.cn/research/PremPLI/ and allows to do large-scale mutational scanning.

Sun et al. present PremPLI, a machine learning approach and web tool to predict how missense mutations in a drug’s target protein will affect the drug’s binding affinity. PremPLI can be applied to identify drug resistance mechanisms in cancer cells and microorganisms and develop drugs to combat drug resistance.

Insights on peptide topology in the computational design of protein ligands: the example of lysozyme binding peptides

Herein, we compared the ability of linear and cyclic peptides generated in silico to target different protein sites: internal pockets and solvent-exposed sites. We selected human lysozyme (HuL) as a model target protein combined with the computational evolution of linear and cyclic peptides. The sequence evolution of these peptides was based on the PARCE algorithm. The generated peptides were screened based on their aqueous solubility and HuL binding affinity. The latter was evaluated by means of scoring functions and atomistic molecular dynamics (MD) trajectories in water, which allowed prediction of the structural features of the protein-peptide complexes. The computational results demonstrated that cyclic peptides constitute the optimal choice for solvent exposed sites, while both linear and cyclic peptides are capable of targeting the HuL pocket effectively. The most promising binders found in silico were investigated experimentally by surface plasmon resonance (SPR), nuclear magnetic resonance (NMR), and electrospray ionization mass spectrometry (ESI-MS) techniques. All tested peptides displayed dissociation constants in the micromolar range, as assessed by SPR; however, both NMR and ESI-MS suggested multiple binding modes, at least for the pocket binding peptides. A detailed NMR analysis confirmed that both linear and cyclic pocket peptides correctly target the binding site they were designed for.

Systematic Development of Peptide Inhibitors Targeting the CXCL12/HMGB1 Interaction

During inflammatory reactions, the production and release of chemotactic factors guide the recruitment of selective leukocyte subpopulations. The alarmin HMGB1 and the chemokine CXCL12, both released in the microenvironment, can form a heterocomplex, which exclusively acts on the chemokine receptor CXCR4, enhancing cell migration, and in some pathological conditions such as rheumatoid arthritis exacerbates the immune response. An excessive cell influx at the inflammatory site can be diminished by disrupting the heterocomplex. Here, we report the computationally driven identification of the first peptide (HBP08) binding HMGB1 and selectively inhibiting the activity of the CXCL12/HMGB1 heterocomplex. Furthermore, HBP08 binds HMGB1 with the highest affinity reported so far (K_d of 0.8 ± 0.4 μM). The identification of this peptide represents an important step toward the development of innovative pharmacological tools for the treatment of severe chronic inflammatory conditions characterized by an uncontrolled immune response.

Fast Prediction of Binding Affinities of the SARS-CoV-2 Spike Protein Mutant N501Y (UK Variant) with ACE2 and Miniprotein Drug Candidates

A recently identified variant of SARS-CoV-2 virus, known as the United Kingdom (UK) variant (lineage B.1.1.7), has an N501Y mutation on its spike protein. SARS-CoV-2 spike protein binds with angiotensin-converting enzyme 2 (ACE2), a key protein for the viral entry into the host cells. Here, we report an efficient computational approach, including the simple energy minimizations and binding free energy calculations, starting from an experimental structure of the binding complex along with experimental calibration of the calculated binding free energies, to rapidly and reliably predict the binding affinities of the N501Y mutant with human ACE2 (hACE2) and recently reported miniprotein and hACE2 decoy (CTC-445.2) drug candidates. It has been demonstrated that the N501Y mutation markedly increases the ACE2-spike protein binding affinity (Kd) from 22 to 0.44 nM, which could partially explain why the UK variant is more infectious. The miniproteins are predicted to have ∼10,000- to 100,000-fold diminished binding affinities with the N501Y mutant, creating a need for design of novel therapeutic candidates to overcome the N501Y mutation-induced drug resistance. The N501Y mutation is also predicted to decrease the binding affinity of a hACE2 decoy (CTC-445.2) binding with the spike protein by ∼200-fold. This convenient computational approach along with experimental calibration may be similarly used in the future to predict the binding affinities of potential new variants of the spike protein.

Potent Antimalarials with Development Potential Identified by Structure-Guided Computational Optimization of a Pyrrole-Based Dihydroorotate Dehydrogenase Inhibitor Series

Dihydroorotate dehydrogenase (DHODH) has been clinically validated as a target for the development of new antimalarials. Experience with clinical candidate triazolopyrimidine DSM265 (1) suggested that DHODH inhibitors have great potential for use in prophylaxis, which represents an unmet need in the malaria drug discovery portfolio for endemic countries, particularly in areas of high transmission in Africa. We describe a structure-based computationally driven lead optimization program of a pyrrole-based series of DHODH inhibitors, leading to the discovery of two candidates for potential advancement to preclinical development. These compounds have improved physicochemical properties over prior series frontrunners and they show no time-dependent CYP inhibition, characteristic of earlier compounds. Frontrunners have potent antimalarial activity in vitro against blood and liver schizont stages and show good efficacy in Plasmodium falciparum SCID mouse models. They are equally active against P. falciparum and Plasmodium vivax field isolates and are selective for Plasmodium DHODHs versus mammalian enzymes.

Variational Method for Networkwide Analysis of Relative Ligand Binding Free Energies with Loop Closure and Experimental Constraints

We describe an efficient method for the simultaneous solution of all free energies within a relative binding free-energy (RBFE) network with cycle closure and experimental/reference constraint conditions using Bennett Acceptance Ratio (BAR) and Multistate BAR (MBAR) analysis. Rather than solving the BAR or MBAR equations for each transformation independently, the simultaneous solution of all transformations are obtained by performing a constrained minimization of a global objective function. The nonlinear optimization of the objective function is subjected to affine linear constraints that couple the free energies between the network edges. The constraints are used to enforce the closure of thermodynamic cycles within the RBFE network, and to enforce an additional set of linear constraint conditions demonstrated here to be subsets of (1 or 2) experimental values. We describe details of the practical implementation of the network BAR/MBAR procedure, including use of generalized coordinates in the minimization of the free-energy objective function, propagation of bootstrap errors from those coordinates, and performance and memory optimization. In some cases it is found that use of restraints in the optimization is more practical than use of generalized coordinates for enforcing constraint conditions. The fast BARnet and MBARnet methods are used to analyze the RBFEs of six prototypical protein-ligand systems, and it is shown that enforcement of cycle closure conditions reduces the error in the predictions only modestly, and further reduction in errors can be achieved when one or two experimental RBFEs are included in the optimization procedure. These methods have been implemented into FE-ToolKit, a new free-energy analysis toolkit. The BARnet/MBARnet framework presented here opens the door to new, more efficient and robust free-energy analysis with enhanced predictive capability for drug discovery applications.

Automation of absolute protein-ligand binding free energy calculations for docking refinement and compound evaluation

Absolute binding free energy calculations with explicit solvent molecular simulations can provide estimates of protein-ligand affinities, and thus reduce the time and costs needed to find new drug candidates. However, these calculations can be complex to implement and perform. Here, we introduce the software BAT.py, a Python tool that invokes the AMBER simulation package to automate the calculation of binding free energies for a protein with a series of ligands. The software supports the attach-pull-release (APR) and double decoupling (DD) binding free energy methods, as well as the simultaneous decoupling-recoupling (SDR) method, a variant of double decoupling that avoids numerical artifacts associated with charged ligands. We report encouraging initial test applications of this software both to re-rank docked poses and to estimate overall binding free energies. We also show that it is practical to carry out these calculations cheaply by using graphical processing units in common machines that can be built for this purpose. The combination of automation and low cost positions this procedure to be applied in a relatively high-throughput mode and thus stands to enable new applications in early-stage drug discovery.

Binding Thermodynamics and Interaction Patterns of Inhibitor-Major Urinary Protein-I Binding from Extensive Free-Energy Calculations: Benchmarking AMBER Force Fields

Mouse major urinary protein (MUP) plays a key role in the pheromone communication system. The one-end-closed β-barrel of MUP-I forms a small, deep, and hydrophobic central cavity, which could accommodate structurally diverse ligands. Previous computational studies employed old protein force fields and short simulation times to determine the binding thermodynamics or investigated only a small number of structurally similar ligands, which resulted in sampled regions far from the experimental structure, nonconverged sampling outcomes, and limited understanding of the possible interaction patterns that the cavity could produce. In this work, extensive end-point and alchemical free-energy calculations with advanced protein force fields were performed to determine the binding thermodynamics of a series of MUP-inhibitor systems and investigate the inter- and intramolecular interaction patterns. Three series of inhibitors with a total of 14 ligands were simulated. We independently simulated the MUP-inhibitor complexes under two advanced AMBER force fields. Our benchmark test showed that the advanced AMBER force fields including AMBER19SB and AMBER14SB provided better descriptions of the system, and the backbone root-mean-square deviation (RMSD) was significantly lowered compared with previous computational studies with old protein force fields. Surprisingly, although the latest AMBER force field AMBER19SB provided better descriptions of various observables, it neither improved the binding thermodynamics nor lowered the backbone RMSD compared with the previously proposed and widely used AMBER14SB. The older but widely used AMBER14SB actually achieved better performance in the prediction of binding affinities from the alchemical and end-point free-energy calculations. We further analyzed the protein-ligand interaction networks to identify important residues stabilizing the bound structure. Six residues including PHE38, LEU40, PHE90, ALA103, LEU105, and TYR120 were found to contribute the most significant part of protein-ligand interactions, and 10 residues were found to provide favorable interactions stabilizing the bound state. The two AMBER force fields gave extremely similar interaction networks, and the secondary structures also showed similar behavior. Thus, the intra- and intermolecular interaction networks described with the two AMBER force fields are similar. Therefore, AMBER14SB could still be the default option in free-energy calculations to achieve highly accurate binding thermodynamics and interaction patterns.

Rapid, accurate, precise and reproducible ligand-protein binding free energy prediction

A central quantity of interest in molecular biology and medicine is the free energy of binding of a molecule to a target biomacromolecule. Until recently, the accurate prediction of binding affinity had been widely regarded as out of reach of theoretical methods owing to the lack of reproducibility of the available methods, not to mention their complexity, computational cost and time-consuming procedures. The lack of reproducibility stems primarily from the chaotic nature of classical molecular dynamics (MD) and the associated extreme sensitivity of trajectories to their initial conditions. Here, we review computational approaches for both relative and absolute binding free energy calculations, and illustrate their application to a diverse set of ligands bound to a range of proteins with immediate relevance in a number of medical domains. We focus on ensemble-based methods which are essential in order to compute statistically robust results, including two we have recently developed, namely thermodynamic integration with enhanced sampling and enhanced sampling of MD with an approximation of continuum solvent. Together, these form a set of rapid, accurate, precise and reproducible free energy methods. They can be used in real-world problems such as hit-to-lead and lead optimization stages in drug discovery, and in personalized medicine. These applications show that individual binding affinities equipped with uncertainty quantification may be computed in a few hours on a massive scale given access to suitable high-end computing resources and workflow automation. A high level of accuracy can be achieved using these approaches.

Alchemical Binding Free Energy Calculations in AMBER20: Advances and Best Practices for Drug Discovery

Predicting protein-ligand binding affinities and the associated thermodynamics of biomolecular recognition is a primary objective of structure-based drug design. Alchemical free energy simulations offer a highly accurate and computationally efficient route to achieving this goal. While the AMBER molecular dynamics package has successfully been used for alchemical free energy simulations in academic research groups for decades, widespread impact in industrial drug discovery settings has been minimal because of the previous limitations within the AMBER alchemical code, coupled with challenges in system setup and postprocessing workflows. Through a close academia-industry collaboration we have addressed many of the previous limitations with an aim to improve accuracy, efficiency, and robustness of alchemical binding free energy simulations in industrial drug discovery applications. Here, we highlight some of the recent advances in AMBER20 with a focus on alchemical binding free energy (BFE) calculations, which are less computationally intensive than alternative binding free energy methods where full binding/unbinding paths are explored. In addition to scientific and technical advances in AMBER20, we also describe the essential practical aspects associated with running relative alchemical BFE calculations, along with recommendations for best practices, highlighting the importance not only of the alchemical simulation code but also the auxiliary functionalities and expertise required to obtain accurate and reliable results. This work is intended to provide a contemporary overview of the scientific, technical, and practical issues associated with running relative BFE simulations in AMBER20, with a focus on real-world drug discovery applications.

Estimating the ligand-binding affinity via λ-dependent umbrella sampling simulations

The umbrella sampling (US) approach has been demonstrated to be a very efficient method for estimating the ligand-binding affinity. However, most of the calculated values overestimate experimental ones that are probably caused by the inaccurate representation of the interaction between the ligand and the surrounding molecules. The issue can be resolved via the implementation aspects of λ-alteration simulation into the US approach, which we call the λ-dependent umbrella sampling (λUS) scheme. In particular, the electrostatic and van der Waals interactions were simultaneously changed by using the coupling parameter λ during λUS simulations. The mean value of obtained results, ∆ G US λ = 0.20 = - 11.59 ± 1.51 kcal mol^-1 , is in good fitting to the mean value of respective experiments, ∆G_EXP =  - 11.26 ± 0.89 kcal mol^-1 . Moreover, the correlation between the proposed approach and experiment is quite good with a value of R US λ = 0.20 = 0.82 ± 0.10 . The λUS scheme significantly enhances the calculated accuracy since the RMSE of the proposed scheme is smaller than traditional US simulations, RMSE US λ = 0.20 = 2.99 ± 0.82 kcal mol^-1 versus RMSE US λ = 0.00 = 5.48 ± 0.81 kcal mol^-1 . Furthermore, the precision is increased since the computed error via λUS approach, δ US λ = 0.20 = 1.51 kcal mol^-1 , was smaller than those of the US simulation, δ US λ = 0.00 = 1.78 kcal mol^-1 . Overall, the proposed approach perhaps provides an efficient way to accurately and precisely estimate the ligand-binding free energy.

Identify potent SARS-CoV-2 main protease inhibitors via accelerated free energy perturbation-based virtual screening of existing drugs

Drug repurposing effort for treatment of a new disease, such as COVID-19, usually starts from a virtual screening of existing drugs, followed by experimental validation, but the actual hit rate is generally rather low with traditional computational methods. It has been demonstrated that a new virtual screening approach with accelerated free energy perturbation-based absolute binding free energy (FEP-ABFE) predictions can reach an unprecedentedly high hit rate, leading to successful identification of 15 potent inhibitors of SARS-CoV-2 main protease (M^pro) from 25 computationally selected drugs under a threshold of K_i = 4 μM. The outcomes of this study are valuable for not only drug repurposing to treat COVID-19 but also demonstrating the promising potential of the FEP-ABFE prediction-based virtual screening approach.

The COVID-19 pandemic caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) has become a global crisis. There is no therapeutic treatment specific for COVID-19. It is highly desirable to identify potential antiviral agents against SARS-CoV-2 from existing drugs available for other diseases and thus repurpose them for treatment of COVID-19. In general, a drug repurposing effort for treatment of a new disease, such as COVID-19, usually starts from a virtual screening of existing drugs, followed by experimental validation, but the actual hit rate is generally rather low with traditional computational methods. Here we report a virtual screening approach with accelerated free energy perturbation-based absolute binding free energy (FEP-ABFE) predictions and its use in identifying drugs targeting SARS-CoV-2 main protease (M^pro). The accurate FEP-ABFE predictions were based on the use of a restraint energy distribution (RED) function, making the practical FEP-ABFE−based virtual screening of the existing drug library possible. As a result, out of 25 drugs predicted, 15 were confirmed as potent inhibitors of SARS-CoV-2 M^pro. The most potent one is dipyridamole (inhibitory constant K_i = 0.04 µM) which has shown promising therapeutic effects in subsequently conducted clinical studies for treatment of patients with COVID-19. Additionally, hydroxychloroquine (K_i = 0.36 µM) and chloroquine (K_i = 0.56 µM) were also found to potently inhibit SARS-CoV-2 M^pro. We anticipate that the FEP-ABFE prediction-based virtual screening approach will be useful in many other drug repurposing or discovery efforts.

X-Ray Crystallography and Free Energy Calculations Reveal the Binding Mechanism of A2A Adenosine Receptor Antagonists

We present a robust protocol based on iterations of free energy perturbation (FEP) calculations, chemical synthesis, biophysical mapping and X-ray crystallography to reveal the binding mode of an antagonist series to the A2A adenosine receptor (AR). Eight A2A AR binding site mutations from biophysical mapping experiments were initially analyzed with sidechain FEP simulations, performed on alternate binding modes. The results distinctively supported one binding mode, which was subsequently used to design new chromone derivatives. Their affinities for the A2A AR were experimentally determined and investigated through a cycle of ligand-FEP calculations, validating the binding orientation of the different chemical substituents proposed. Subsequent X-ray crystallography of the A2A AR with a low and a high affinity chromone derivative confirmed the predicted binding orientation. The new molecules and structures here reported were driven by free energy calculations, and provide new insights on antagonist binding to the A2A AR, an emerging target in immuno-oncology.

The rise of molecular simulations in fragment-based drug design (FBDD): an overview

Fragment-based drug discovery (FBDD) is an innovative approach, progressively more applied in the academic and industrial context, to enhance hit identification for previously considered undruggable biological targets. In particular, FBDD discovers low-molecular-weight (LMW) ligands (<300Da) able to bind to therapeutically relevant macromolecules in an affinity range from the micromolar (μM) to millimolar (mM). X-ray crystallography (XRC) and nuclear magnetic resonance (NMR) spectroscopy are commonly the methods of choice to obtain 3D information about the bound ligand-protein complex, but this can occasionally be problematic, mainly for early, low-affinity fragments. The recent development of computational fragment-based approaches provides a further strategy for improving the identification of fragment hits. In this review, we summarize the state of the art of molecular dynamics simulations approaches used in FBDD, and discuss limitations and future perspectives for these approaches.

FragPELE: Dynamic Ligand Growing within a Binding Site. A Novel Tool for Hit-To-Lead Drug Design

The early stages of drug discovery rely on hit-to-lead programs, where initial hits undergo partial optimization to improve binding affinities for their biological target. This is an expensive and time-consuming process, requiring multiple iterations of trial and error designs, an ideal scenario for applying computer simulation. However, most state-of-the-art modeling techniques fail to provide a fast and reliable answer to the Induced-Fit protein-ligand problem. To aid in this matter, we present FragPELE, a new tool for in silico hit-to-lead drug design, capable of growing a fragment from a bound core while exploring the protein-ligand conformational space. We tested the ability of FragPELE to predict crystallographic data, even in cases where cryptic sub-pockets open because of the presence of particular R-groups. Additionally, we evaluated the potential of the software on growing and scoring five congeneric series from the 2015 FEP+ dataset, comparing them to FEP+, SP and Induced-Fit Glide, and MMGBSA simulations. Results show that FragPELE could be useful not only for finding new cavities and novel binding modes in cases where standard docking tools cannot but also to rank ligand activities in a reasonable amount of time and with acceptable precision.

Best Practices for Alchemical Free Energy Calculations [Article v1.0]

Alchemical free energy calculations are a useful tool for predicting free energy differences associated with the transfer of molecules from one environment to another. The hallmark of these methods is the use of "bridging" potential energy functions representing alchemical intermediate states that cannot exist as real chemical species. The data collected from these bridging alchemical thermodynamic states allows the efficient computation of transfer free energies (or differences in transfer free energies) with orders of magnitude less simulation time than simulating the transfer process directly. While these methods are highly flexible, care must be taken in avoiding common pitfalls to ensure that computed free energy differences can be robust and reproducible for the chosen force field, and that appropriate corrections are included to permit direct comparison with experimental data. In this paper, we review current best practices for several popular application domains of alchemical free energy calculations performed with equilibrium simulations, in particular relative and absolute small molecule binding free energy calculations to biomolecular targets.

What Next for Quantum Mechanics in Structure-Based Drug Discovery?

There is significant potential for electronic structure methods to improve the quality of the predictions furnished by the tools of computer-aided drug design, which typically rely on empirically derived functions. In this perspective, we consider some recent examples of how quantum mechanics has been applied in predicting protein-ligand geometries, protein-ligand binding affinities and ligand strain on binding. We then outline several significant developments in quantum mechanics methodology likely to influence these approaches: in particular, we note the advent of more computationally expedient ab initio quantum mechanical methods that can provide chemical accuracy for larger molecular systems than hitherto possible. We highlight the emergence of increasingly accurate semiempirical quantum mechanical methods and the associated role of machine learning and molecular databases in their development. Indeed, the convergence of improved algorithms for solving and analyzing electronic structure, modern machine learning methods, and increasingly comprehensive benchmark data sets of molecular geometries and energies provides a context in which the potential of quantum mechanics will be increasingly realized in driving future developments and applications in structure-based drug discovery.