Test MM-PB/SA on true conformational ensembles of protein-ligand complexes

Machine Learning Methods as a Cost-Effective Alternative to Physics-Based Binding Free Energy Calculations

The rank ordering of ligands remains one of the most attractive challenges in drug discovery. While physics-based in silico binding affinity methods dominate the field, they still have problems, which largely revolve around forcefield accuracy and sampling. Recent advances in machine learning have gained traction for protein-ligand binding affinity predictions in early drug discovery programs. In this article, we perform retrospective binding free energy evaluations for 172 compounds from our internal collection spread over four different protein targets and five congeneric ligand series. We compared multiple state-of-the-art free energy methods ranging from physics-based methods with different levels of complexity and conformational sampling to state-of-the-art machine-learning-based methods that were available to us. Overall, we found that physics-based methods behaved particularly well when the ligand perturbations were made in the solvation region, and they did not perform as well when accounting for large conformational changes in protein active sites. On the other end, machine-learning-based methods offer a good cost-effective alternative for binding free energy calculations, but the accuracy of their predictions is highly dependent on the experimental data available for training the model.

In Silico Investigations into the Selectivity of Psychoactive and New Psychoactive Substances in Monoamine Transporters

New psychoactive substances (NPS) are a group of compounds that mimic the effects of illicit substances. A range of NPS have been shown to interact with the three main classes of monoamine transporters (DAT, NET, and SERT) to differing extents, but it is unclear why these differences arise. To aid in understanding the differences in affinity between the classes of monoamine transporters, several in silico experiments were conducted. Docking experiments showed there was no direct correlation between a range of scoring functions and experimental activity, but Spearman ranking analysis showed a significant correlation (α = 0.1) for DAT, with the affinity ΔG (0.42), αHB (0.40), GoldScore (0.40), and PLP (0.41) scoring functions, and for DAT (0.38) and SERT (0.40) using a consensus scoring approach. Qualitative structure-activity relationship (QSAR) experiments resulted in the generation of robust and predictive three-descriptor models for SERT (r ² = 0.87, q ² = 0.8, and test set r ² = 0.74) and DAT (r ² = 0.68, q ² = 0.51, test set r ² = 0.63). Both QSAR models described similar characteristics for binding, i.e., rigid hydrophobic molecules with a biogenic amine moiety, and were not sufficient to facilitate a deeper understanding of differences in affinity between the monoamine transporters. This contextualizes the observed promiscuity for NPS between the isoforms and highlights the difficulty in the design and development of compounds that are isoform-selective.

Using AMBER18 for Relative Free Energy Calculations

With renewed interest in free energy methods in contemporary structure-based drug design, there is a pressing need to validate against multiple targets and force fields to assess the overall ability of these methods to accurately predict relative binding free energies. We computed relative binding free energies using graphics processing unit accelerated thermodynamic integration (GPU-TI) on a data set originally assembled by Schrödinger, Inc. Using their GPU free energy code (FEP+) and the OPLS2.1 force field combined with the REST2 enhanced sampling approach, these authors obtained an overall MUE of 0.9 kcal/mol and an overall RMSD of 1.14 kcal/mol. In our study using GPU-TI from AMBER with the AMBER14SB/GAFF1.8 force field but without enhanced sampling, we obtained an overall MUE of 1.17 kcal/mol and an overall RMSD of 1.50 kcal/mol for the 330 perturbations contained in this data set. A more detailed analysis of our results suggested that the observed differences between the two studies arise from differences in sampling protocols along with differences in the force fields employed. Future work should address the problem of establishing benchmark quality results with robust statistical error bars obtained through multiple independent runs and enhanced sampling, which is possible with the GPU-accelerated features in AMBER.

Detailed potential of mean force studies on host-guest systems from the SAMPL6 challenge

Accurately predicting receptor-ligand binding free energies is one of the holy grails of computational chemistry with many applications in chemistry and biology. Many successes have been reported, but issues relating to sampling and force field accuracy remain significant issues affecting our ability to reliably calculate binding free energies. In order to explore these issues in more detail we have examined a series of small host-guest complexes from the SAMPL6 blind challenge, namely octa-acids (OAs)-guest complexes and Curcurbit[8]uril (CB8)-guest complexes. Specifically, potential of mean force studies using umbrella sampling combined with the weighted histogram method were carried out on both systems with both known and unknown binding affinities. We find that using standard force fields and straightforward simulation protocols we are able to obtain satisfactory results, but that simply scaling our results allows us to significantly improve our predictive ability for the unknown test sets: the overall RMSD of the binding free energy versus experiment is reduced from 5.59 to 2.36 kcal/mol; for the CB8 test system, the RMSD goes from 8.04 to 3.51 kcal/mol, while for the OAs test system, the RSMD goes from 2.89 to 0.95 kcal/mol. The scaling approach was inspired by studies on structurally related known benchmark sets: by simply scaling, the RMSD was reduced from 6.23 to 1.19 kcal/mol and from 2.96 to 0.62 kcal/mol for the CB8 benchmark system and the OA benchmark system, respectively. We find this scaling procedure to correct absolute binding affinities to be highly effective especially when working across a "congeneric" series with similar charge states. It is less successful when applied to mixed ligands with varied charges and chemical characteristics, but improvement is still realized in the present case. This approach suggests that there are large systematic errors in absolute binding free energy calculations that can be straightforwardly accounted for using a scaling procedure. Random errors are still an issue, but near chemical accuracy can be obtained using the present strategy in select cases.

Discovery of Novel Adenosine Receptor Antagonists through a Combined Structure- and Ligand-Based Approach Followed by Molecular Dynamics Investigation of Ligand Binding Mode

An intense effort is made by pharmaceutical and academic research laboratories to identify and develop selective antagonists for each adenosine receptor (AR) subtype as potential clinical candidates for "soft" treatment of various diseases. Crystal structures of subtypes A_2A and A₁ARs offer exciting opportunities for structure-based drug design. In the first part of the present work, Maybridge HitFinder library of 14400 compounds was utilized to apply a combination of structure-based against the crystal structure of A_2AAR and ligand-based methodologies. The docking poses were rescored by CHARMM energy minimization and calculation of the desolvation energy using Poisson-Boltzmann equation electrostatics. Out of the eight selected and tested compounds, five were found positive hits (63% success). Although the project was initially focused on targeting A_2AAR, the identified antagonists exhibited low micromolar or micromolar affinity against A_2A/A₃, ARs, or A₃AR, respectively. Based on these results, 19 compounds characterized by novel chemotypes were purchased and tested. Sixteen of them were identified as AR antagonists with affinity toward combinations of the AR family isoforms (A_2A/A₃, A₁/A₃, A₁/A_2A/A₃, and A₃). The second part of this work involves the performance of hundreds of molecular dynamics (MD) simulations of complexes between the ARs and a total of 27 ligands to resolve the binding interactions of the active compounds, which were not achieved by docking calculations alone. This computational work allowed the prediction of stable and unstable complexes which agree with the experimental results of potent and inactive compounds, respectively. Of particular interest is that the 2-amino-thiophene-3-carboxamides, 3-acylamino-5-aryl-thiophene-2-carboxamides, and carbonyloxycarboximidamide derivatives were found to be selective and possess a micromolar to low micromolar affinity for the A₃ receptor.

Forging the Basis for Developing Protein-Ligand Interaction Scoring Functions

In structure-based drug design, scoring functions are widely used for fast evaluation of protein-ligand interactions. They are often applied in combination with molecular docking and de novo design methods. Since the early 1990s, a whole spectrum of protein-ligand interaction scoring functions have been developed. Regardless of their technical difference, scoring functions all need data sets combining protein-ligand complex structures and binding affinity data for parametrization and validation. However, data sets of this kind used to be rather limited in terms of size and quality. On the other hand, standard metrics for evaluating scoring function used to be ambiguous. Scoring functions are often tested in molecular docking or even virtual screening trials, which do not directly reflect the genuine quality of scoring functions. Collectively, these underlying obstacles have impeded the invention of more advanced scoring functions. In this Account, we describe our long-lasting efforts to overcome these obstacles, which involve two related projects. On the first project, we have created the PDBbind database. It is the first database that systematically annotates the protein-ligand complexes in the Protein Data Bank (PDB) with experimental binding data. This database has been updated annually since its first public release in 2004. The latest release (version 2016) provides binding data for 16 179 biomolecular complexes in PDB. Data sets provided by PDBbind have been applied to many computational and statistical studies on protein-ligand interaction and various subjects. In particular, it has become a major data resource for scoring function development. On the second project, we have established the Comparative Assessment of Scoring Functions (CASF) benchmark for scoring function evaluation. Our key idea is to decouple the "scoring" process from the "sampling" process, so scoring functions can be tested in a relatively pure context to reflect their quality. In our latest work on this track, i.e. CASF-2013, the performance of a scoring function was quantified in four aspects, including "scoring power", "ranking power", "docking power", and "screening power". All four performance tests were conducted on a test set containing 195 high-quality protein-ligand complexes selected from PDBbind. A panel of 20 standard scoring functions were tested as demonstration. Importantly, CASF is designed to be an open-access benchmark, with which scoring functions developed by different researchers can be compared on the same grounds. Indeed, it has become a popular choice for scoring function validation in recent years. Despite the considerable progress that has been made so far, the performance of today's scoring functions still does not meet people's expectations in many aspects. There is a constant demand for more advanced scoring functions. Our efforts have helped to overcome some obstacles underlying scoring function development so that the researchers in this field can move forward faster. We will continue to improve the PDBbind database and the CASF benchmark in the future to keep them as useful community resources.

Recent Advancements in Docking Methodologies

Nowadays molecular docking has become an important methodology in CADD (Computer-Aided Drug Design)-assisted drug discovery process. It is an important computational tool widely used to predict binding mode, binding affinity and binding free energy of a protein-ligand complex. The important factors responsible for accurate results in docking studies are correct binding site prediction, use of suitable small-molecule databases, consistent docking pose, high dock score with good MD (Molecular Dynamics), clarity whether the compound is an inhibitor or agonist, etc. However, still there are several limitations which make it difficult to obtain accurate results from docking studies. In this chapter, the main focus is on recent advancements in various aspects of molecular docking such as ligand sampling, protein flexibility, scoring functions, fragment docking, post-processing, docking into homology models and protein-protein docking.

The MM/PBSA and MM/GBSA methods to estimate ligand-binding affinities

INTRODUCTION

The molecular mechanics energies combined with the Poisson-Boltzmann or generalized Born and surface area continuum solvation (MM/PBSA and MM/GBSA) methods are popular approaches to estimate the free energy of the binding of small ligands to biological macromolecules. They are typically based on molecular dynamics simulations of the receptor-ligand complex and are therefore intermediate in both accuracy and computational effort between empirical scoring and strict alchemical perturbation methods. They have been applied to a large number of systems with varying success.

AREAS COVERED

The authors review the use of MM/PBSA and MM/GBSA methods to calculate ligand-binding affinities, with an emphasis on calibration, testing and validation, as well as attempts to improve the methods, rather than on specific applications.

EXPERT OPINION

MM/PBSA and MM/GBSA are attractive approaches owing to their modular nature and that they do not require calculations on a training set. They have been used successfully to reproduce and rationalize experimental findings and to improve the results of virtual screening and docking. However, they contain several crude and questionable approximations, for example, the lack of conformational entropy and information about the number and free energy of water molecules in the binding site. Moreover, there are many variants of the method and their performance varies strongly with the tested system. Likewise, most attempts to ameliorate the methods with more accurate approaches, for example, quantum-mechanical calculations, polarizable force fields or improved solvation have deteriorated the results.

Practical Considerations in Virtual Screening and Molecular Docking

Molecular docking has become an important common component of the drug discovery toolbox, and its relative low-cost implications and perceived simplicity of use has stimulated an everincreasing popularity within academic communities. The inherent “garbage-in-garbage-out” defect of molecular docking, however, leads a lot of researchers to dedicate countless hours to the identification of hit compounds that later prove to be inactive. Several considerations that can greatly improve the success and enrichment of true bioactive hit compounds are commonly overlooked at the initial stages of a molecular docking study. This chapter will cover several of these considerations, including protonation states, active site waters, separating actives from decoys, consensus docking and molecular mechanics generalized-Born/surface area (MM-GBSA) rescoring, and incorporation of pharmacophoric constraints, in an attempt to clarify what is, in fact, very complicated and inherent difficulties of a structure-based drug design study.

Bringing Clarity to the Prediction of Protein-Ligand Binding Free Energies via "Blurring"

We present a method to evaluate the free energies of ligand binding utilizing a Monte Carlo estimation of the configuration integrals concomitant with uncertainty quantification. Ensembles for integration are built through systematically perturbing an initial ligand conformation in a rigid binding pocket, which is optimized separately prior to incorporation of the ligand. We call the procedure producing the ensembles "blurring", and it is carried out using an in-house developed code. The Boltzmann factor contribution of each pose to the configuration integral is computed and from there the free energy is obtained. Potential function uncertainties are estimated using a fragment-based error propagation method. This method has been applied to a set of small aromatic ligands complexed with T4 Lysozyme L99A mutant. Microstate energies have been determined with the force fields ff99SB and ff94, and the semiempirical method PM6DH2 in conjunction with continuum solvation models including Generalized Born (GB), the Conductor-like Screening Model (COSMO), and SMD. Of the methods studied, PM6DH2-based scoring gave binding free energy estimates, which yielded a good correlation to the experimental binding affinities (R2 = 0.7). All methods overestimated the calculated binding affinities. We trace this to insufficient sampling, the single static protein structure, and inaccuracies in the solvent models we have used in this study.

Receptor-ligand molecular docking

Docking methodology aims to predict the experimental binding modes and affinities of small molecules within the binding site of particular receptor targets and is currently used as a standard computational tool in drug design for lead compound optimisation and in virtual screening studies to find novel biologically active molecules. The basic tools of a docking methodology include a search algorithm and an energy scoring function for generating and evaluating ligand poses. In this review, we present the search algorithms and scoring functions most commonly used in current molecular docking methods that focus on protein-ligand applications. We summarise the main topics and recent computational and methodological advances in protein-ligand docking. Protein flexibility, multiple ligand binding modes and the free-energy landscape profile for binding affinity prediction are important and interconnected challenges to be overcome by further methodological developments in the docking field.

Engineering protein therapeutics: predictive performances of a structure-based virtual affinity maturation protocol

The implementation of a structure based virtual affinity maturation protocol and evaluation of its predictivity are presented. The in silico protocol is based on conformational sampling of the interface residues (using the Dead End Elimination/A* algorithm), followed by the estimation of the change of free energy of binding due to a point mutation, applying MM/PBSA calculations. Several implementations of the protocol have been evaluated for 173 mutations in 7 different protein complexes for which experimental data were available: the use of the Boltzamnn averaged predictor based on the free energy of binding (ΔΔG(*)) combined with the one based on its polar component only (ΔΔE(pol*)) led to the proposal of a subset of mutations out of which 45% would have successfully enhanced the binding. When focusing on those mutations that are less likely to be introduced by natural in vivo maturation methods (99 mutations with at least two base changes in the codon), the success rate is increased to 63%. In another evaluation, focusing on 56 alanine scanning mutations, the in silico protocol was able to detect 89% of the hot-spots.

Use of virtual screening for discovering antiretroviral compounds interacting with the HIV-1 nucleocapsid protein

The HIV-1 nucleocapsid protein (NC) is considered as an emerging drug target for the therapy of AIDS. Several studies have highlighted the crucial role of NC within the viral replication cycle. However, although NC inhibition has provided in vitro and in vivo antiretroviral activity, drug-candidates which interfere with NC functions are still missing in the therapeutic arsenal against HIV. Based on previous studies, where the dynamic behavior of NC and its ligand binding properties have been investigated by means of computational methods, here we used a virtual screening protocol for discovering novel antiretroviral compounds which interact with NC. The antiretroviral activity of virtual hits was tested in vitro, whereas biophysical studies elucidated the direct interaction of most active compounds with NC(11-55), a peptide corresponding to the zinc finger domain of NC. Two novel antiretroviral small molecules capable of interacting with NC are presented here.

Free energy calculations on snake venom metalloproteinase BaP1

BaP1 is a snake venom metalloproteinase from the venom of Bothrops asper, showing high structural homology with the catalytic domain of human adamalysins and matrix metalloproteinases. It induces the release of cytokines, like interleukin-1 and tumor necrosis factor alpha. Recently, the high-resolution crystal structure of BaP1 with a bound inhibitor became available, representing an interesting model concerning inhibitor design for medicinally important metalloproteinases such as tumor necrosis factor alpha-converting enzyme and MMP13. We here use computational modeling to gain a better understanding about the binding properties of various ligands to BaP1, with a focus on computing ligand binding free energies. The obtained results should be of general significance for future research on medicinally important metalloproteinases. We have investigated the binding of the original inhibitor in detail and calculated its binding strength using MMP/GBSA free energy calculations. Additionally, the binding strengths of alternative ligands have been computed, and two of them are predicted and experimentally verified to strongly inhibit the enzyme. A suggestion for chemical modifications of BaP1 inhibitors could be made to guide future synthesis efforts. Furthermore, a contribution to the proteolytic reaction mechanism of metzincins is given. The pK value of the catalytically active glutamic acid residue 143 has been found to be significantly raised when compared with a free glutamate side chain. Calculations on other matrix metalloproteinases confirmed that this is not confined to BaP1, but seems to be a common feature of metzincins.

Solvent interaction energy calculations on molecular dynamics trajectories: increasing the efficiency using systematic frame selection

End-point methods such as linear interaction energy (LIE) analysis, molecular mechanics generalized Born solvent-accessible surface (MM/GBSA), and solvent interaction energy (SIE) analysis have become popular techniques to calculate the free energy associated with protein-ligand binding. Such methods typically use molecular dynamics (MD) simulations to generate an ensemble of protein structures that encompasses the bound and unbound states. The energy evaluation method (LIE, MM/GBSA, or SIE) is subsequently used to calculate the energy of each member of the ensemble, thus providing an estimate of the average free energy difference between the bound and unbound states. The workflow requiring both MD simulation and energy calculation for each frame and each trajectory proves to be computationally expensive. In an attempt to reduce the high computational cost associated with end-point methods, we study several methods by which frames may be intelligently selected from the MD simulation including clustering and address the question of how the number of selected frames influences the accuracy of the SIE calculations.

Prediction of ligand binding using an approach designed to accommodate diversity in protein-ligand interactions

Computational determination of protein-ligand interaction potential is important for many biological applications including virtual screening for therapeutic drugs. The novel internal consensus scoring strategy is an empirical approach with an extended set of 9 binding terms combined with a neural network capable of analysis of diverse complexes. Like conventional consensus methods, internal consensus is capable of maintaining multiple distinct representations of protein-ligand interactions. In a typical use the method was trained using ligand classification data (binding/no binding) for a single receptor. The internal consensus analyses successfully distinguished protein-ligand complexes from decoys (r², 0.895 for a series of typical proteins). Results are superior to other tested empirical methods. In virtual screening experiments, internal consensus analyses provide consistent enrichment as determined by ROC-AUC and pROC metrics.