Lessons learned from participating in D3R 2016 Grand Challenge 2: compounds targeting the farnesoid X receptor

D3R 2016 Grand Challenge 2 focused on predictions of binding modes and affinities for 102 compounds against the farnesoid X receptor (FXR). In this challenge, two distinct methods, a docking-based method and a template-based method, were employed by our team for the binding mode prediction. For the new template-based method, 3D ligand similarities were calculated for each query compound against the ligands in the co-crystal structures of FXR available in Protein Data Bank. The binding mode was predicted based on the co-crystal protein structure containing the ligand with the best ligand similarity score against the query compound. For the FXR dataset, the template-based method achieved a better performance than the docking-based method on the binding mode prediction. For the binding affinity prediction, an in-house knowledge-based scoring function ITScore2 and MM/PBSA approach were employed. Good performance was achieved for MM/PBSA, whereas the performance of ITScore2 was sensitive to ligand composition, e.g. the percentage of carbon atoms in the compounds. The sensitivity to ligand composition could be a clue for the further improvement of our knowledge-based scoring function.     

Binding affinities of the farnesoid X receptor in the D3R Grand Challenge 2 estimated by free-energy perturbation and docking

We have studied the binding of 102 ligands to the farnesoid X receptor within the D3R Grand Challenge 2016 blind-prediction competition. First, we employed docking with five different docking software and scoring functions. The selected docked poses gave an average root-mean-squared deviation of 4.2 Å. Consensus scoring gave decent results with a Kendall’s τ of 0.26?±?0.06 and a Spearman’s ρ of 0.41?±?0.08. For a subset of 33 ligands, we calculated relative binding free energies with free-energy perturbation. Five transformations between the ligands involved a change of the net charge and we implemented and benchmarked a semi-analytic correction (Rocklin et al., J Chem Phys 139:184103, 2013) for artifacts caused by the periodic boundary conditions and Ewald summation. The results gave a mean absolute deviation of 7.5 kJ/mol compared to the experimental estimates and a correlation coefficient of R ²?=?0.1. These results were among the four best in this competition out of 22 submissions. The charge corrections were significant (7–8 kJ/mol) and always improved the results. By employing 23 intermediate states in the free-energy perturbation, there was a proper overlap between all states and the precision was 0.1–0.7 kJ/mol. However, thermodynamic cycles indicate that the sampling was insufficient in some of the perturbations.     

Molecular docking performance evaluated on the D3R Grand Challenge 2015 drug-like ligand datasets

The D3R Grand Challenge 2015 was focused on two protein targets: Heat Shock Protein 90 (HSP90) and Mitogen-Activated Protein Kinase Kinase Kinase Kinase 4 (MAP4K4). We used a protocol involving a preliminary analysis of the available data in PDB and PubChem BioAssay, and then a docking/scoring step using more computationally demanding parameters that were required to provide more reliable predictions. We could evidence that different docking software and scoring functions can behave differently on individual ligand datasets, and that the flexibility of specific binding site residues is a crucial element to provide good predictions.     

Performance of multiple docking and refinement methods in the pose prediction D3R prospective Grand Challenge 2016

We describe the performance of multiple pose prediction methods for the D3R 2016 Grand Challenge. The pose prediction challenge includes 36 ligands, which represent 4 chemotypes and some miscellaneous structures against the FXR ligand binding domain. In this study we use a mix of fully automated methods as well as human-guided methods with considerations of both the challenge data and publicly available data. The methods include ensemble docking, colony entropy pose prediction, target selection by molecular similarity, molecular dynamics guided pose refinement, and pose selection by visual inspection. We evaluated the success of our predictions by method, chemotype, and relevance of publicly available data. For the overall data set, ensemble docking, visual inspection, and molecular dynamics guided pose prediction performed the best with overall mean RMSDs of 2.4, 2.2, and 2.2 Å respectively. For several individual challenge molecules, the best performing method is evaluated in light of that particular ligand. We also describe the protein, ligand, and public information data preparations that are typical of our binding mode prediction workflow.     

Workflows and performances in the ranking prediction of 2016 D3R Grand Challenge 2: lessons learned from a collaborative effort

The 2016 D3R Grand Challenge 2 includes both pose and affinity or ranking predictions. This article is focused exclusively on affinity predictions submitted to the D3R challenge from a collaborative effort of the modeling and informatics group. Our submissions include ranking of 102 ligands covering 4 different chemotypes against the FXR ligand binding domain structure, and the relative binding affinity predictions of the two designated free energy subsets of 15 and 18 compounds. Using all the complex structures prepared in the same way allowed us to cover many types of workflows and compare their performances effectively. We evaluated typical workflows used in our daily structure-based design modeling support, which include docking scores, force field-based scores, QM/MM, MMGBSA, MD-MMGBSA, and MacroModel interaction energy estimations. The best performing methods for the two free energy subsets are discussed. Our results suggest that affinity ranking still remains very challenging; that the knowledge of more structural information does not necessarily yield more accurate predictions; and that visual inspection and human intervention are considerably important for ranking. Knowledge of the mode of action and protein flexibility along with visualization tools that depict polar and hydrophobic maps are very useful for visual inspection. QM/MM-based workflows were found to be powerful in affinity ranking and are encouraged to be applied more often. The standardized input and output enable systematic analysis and support methodology development and improvement for high level blinded predictions.     

Sensitivity of molecular docking to induced fit effects in influenza virus neuraminidase

Many proteins undergo small side chain or even backbone movements on binding of different ligands into the same protein structure. This is known as induced fit and is potentially problematic for virtual screening of databases against protein targets. In this report we investigate the limits of the rigid protein approximation used by the docking program, GOLD, through cross-docking using protein structures of influenza neuraminidase. Neuraminidase is known to exhibit small but significant induced fit effects on ligand binding. Some neuraminidase crystal structures caused concern due to the bound ligand conformation and GOLD performed poorly on these complexes. A `clean' set, which contained unique, unambiguous complexes, was defined. For this set, the lowest energy structure was correctly docked (i.e. RMSD < 1.5 Å away from the crystal reference structure) in 84% of proteins, and the most promiscuous protein (1mwe) was able to dock all 15 ligands accurately including those that normally required an induced fit movement. This is considerably better than the 70% success rate seen with GOLD against general validation sets. Inclusion of specific water molecules involved in water-mediated hydrogen bonds did not significantly improve the docking performance for ligands that formed water-mediated contacts but it did prevent docking of ligands that displaced these waters. Our data supports the use of a single protein structure for virtual screening with GOLD in some applications involving induced fit effects, although care must be taken to identify the protein structure that performs best against a wide variety of ligands. The performance of GOLD was significantly better than the GOLD implementation of ChemScore and the reasons for this are discussed. Overall, GOLD has shown itself to be an extremely good, robust docking program for this system.     

Sampling and refinement protocols for template-based macrocycle docking: 2018 D3R Grand Challenge 4

We describe a new template-based method for docking flexible ligands such as macrocycles to proteins. It combines Monte-Carlo energy minimization on the manifold, a fast manifold search method, with BRIKARD for complex flexible ligand searching, and with the MELD accelerator of Replica-Exchange Molecular Dynamics simulations for atomistic degrees of freedom. Here we test the method in the Drug Design Data Resource blind Grand Challenge competition. This method was among the best performers in the competition, giving sub-angstrom prediction quality for the majority of the targets.     

Exploring hierarchical refinement techniques for induced fit docking with protein and ligand flexibility

We present a series of molecular‐mechanics‐based protein refinement methods, including two novel ones, applied as part of an induced fit docking procedure. The methods used include minimization; protein and ligand sidechain prediction; a hierarchical ligand placement procedure similar to a‐priori protein loop predictions; and a minimized Monte Carlo approach using normal mode analysis as a move step. The results clearly indicate the importance of a proper opening of the active site backbone, which might not be accomplished when the ligand degrees of freedom are prioritized. The most accurate method consisted of the minimized Monte Carlo procedure designed to open the active site followed by a hierarchical optimization of the sidechain packing around a mobile flexible ligand. The methods have been used on a series of 88 protein‐ligand complexes including both cross‐docking and apo‐docking members resulting in complex conformations determined to within 2.0 Å heavy‐atom RMSD in 75% of cases where the protein backbone rearrangement upon binding is less than 1.0 Å α‐carbon RMSD. We also demonstrate that physics‐based all‐atom potentials can be more accurate than docking‐style potentials when complexes are sufficiently refined. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Docking of small molecules to farnesoid X receptors using AutoDock Vina with the Convex-PL potential: lessons learned from D3R Grand Challenge 2

The 2016 D3R Grand Challenge 2 provided an opportunity to test multiple protein–ligand docking protocols on a set of ligands bound to farnesoid X receptor that has many available experimental structures. We participated in the Stage 1 of the Challenge devoted to the docking pose predictions, with the mean RMSD value of our submission poses of 2.9 Å. Here we present a thorough analysis of our docking predictions made with AutoDock Vina and the Convex-PL rescoring potential by reproducing our submission protocol and running a series of additional molecular docking experiments. We conclude that a correct receptor structure, or more precisely, the structure of the binding pocket, plays the crucial role in the success of our docking studies. We have also noticed the important role of a local ligand geometry, which seems to be not well discussed in literature. We succeed to improve our results up to the mean RMSD value of 2.15–2.33 Å  dependent on the models of the ligands, if docking these to all available homologous receptors. Overall, for docking of ligands of diverse chemical series we suggest to perform docking of each of the ligands to a set of multiple receptors that are homologous to the target.     

CDOCKER and $$\lambda$$-dynamics for prospective prediction in D3R Grand Challenge 2

The opportunity to prospectively predict ligand bound poses and free energies of binding to the Farnesoid X Receptor in the D3R Grand Challenge 2 provided a useful exercise to evaluate CHARMM based docking (CDOCKER) and \(\lambda\)-dynamics methodologies for use in “real-world” applications in computer aided drug design. In addition to measuring their current performance, several recent methodological developments have been analyzed retrospectively to highlight best procedural practices in future applications. For pose prediction with CDOCKER, when the protein structure used for rigid receptor docking was close to the crystallographic holo structure, reliable poses were obtained. Benzimidazoles, with a known holo receptor structure, were successfully docked with an average RMSD of 0.97 \(\AA\). Other non-benzimidazole ligands displayed less accuracy largely because the receptor structures we chose for docking were too different from the experimental holo structures. However, retrospective analysis has shown that when these ligands were re-docked into their holo structures, the average RMSD dropped to 1.18 \(\AA\) for all ligands. When sulfonamides and spiros were docked with the apo structure, which agrees more with their holo structure than the structures we chose, five out of six ligands were correctly docked. These docking results emphasize the need for flexible receptor docking approaches. For \(\lambda\)-dynamics techniques, including multisite \(\lambda\)-dynamics (MS\(\lambda\)D), reasonable agreement with experiment was observed for the 33 ligands investigated; root mean square errors of 2.08 and 1.67 kcal/mol were obtained for free energy sets 1 and 2, respectively. Retrospectively, soft-core potentials, adaptive landscape flattening, and biasing potential replica exchange (BP-REX) algorithms were critical to model large substituent perturbations with sufficient precision and within restrictive timeframes, such as was required with participation in Grand Challenge 2. These developments, their associated benefits, and proposed procedures for their use in future applications are discussed.     

Binding pose and affinity prediction in the 2016 D3R Grand Challenge 2 using the Wilma-SIE method

The Farnesoid X receptor (FXR) exhibits significant backbone movement in response to the binding of various ligands and can be a challenge for pose prediction algorithms. As part of the D3R Grand Challenge 2, we tested Wilma-SIE, a rigid-protein docking method, on a set of 36 FXR ligands for which the crystal structures had originally been blinded. These ligands covered several classes of compounds. To overcome the rigid protein limitations of the method, we used an ensemble of publicly available structures for FXR from the PDB. The use of the ensemble allowed Wilma-SIE to predict poses with average and median RMSDs of 2.3 and 1.4 Å, respectively. It was quite clear, however, that had we used a single structure for the receptor the success rate would have been much lower. The most successful predictions were obtained on chemical classes for which one or more crystal structures of the receptor bound to a molecule of the same class was available. In the absence of a crystal structure for the class, observing a consensus binding mode for the ligands of the class using one or more receptor structures of other classes seemed to be indicative of a reasonable pose prediction. Affinity prediction proved to be more challenging with generally poor correlation with experimental IC50s (Kendall tau?~?0.3). Even when the 36 crystal structures were used the accuracy of the predicted affinities was not appreciably improved. A possible cause of difficulty is the internal energy strain arising from conformational differences in the receptor across complexes, which may need to be properly estimated and incorporated into the SIE scoring function.     

The sensitivity of the results of molecular docking to induced fit effects: Application to thrombin, thermolysin and neuraminidase

This paper describes the application of PRO_–LEADS to the flexible docking of ligands into crystallographically derived enzyme structures that are assumed to be rigid. PRO_–LEADS uses a Tabu search methodology to perform the flexible search and an empirically derived estimate of the binding affinity to drive the docking process. The paper tests the extent to which the assumption of a rigid enzyme compromises the accuracy of the results. All-pairs docking experiments are performed for three enzymes (thrombin, thermolysin and influenza virus neuraminidase) based on six or more ligand-enzyme crystal structures for each enzyme. In 76% of the cases, PRO_–LEADS can successfully identify the correct ligand conformation as the lowest energy configuration when the enzyme structure is derived from that ligand's crystal structure, but the methodology only docks 49% of the cases successfully when the ligand is docked against enzyme crystal structures derived from other ligands. Small movements in the enzyme structure lead to an under-prediction in the energy of the correct binding mode by up to 14 kJ/mol and in some cases this under-prediction can lead to the native mode not being recognised as the lowest energy solution. The type of movements responsible for mis-docking are: the movement of sidechains as a result of changes in C position; the movement of sidechains without changes in C position; the movement of flexible portions of main chains to facilitate the formation of hydrogen bonds; and the movement of metal atoms bound to the enzyme active site. The work illustrates that the assumption of a rigid active site can lead to errors in identification of the correct binding mode and the assessment of binding affinity, even for enzymes which show relatively small shift in atomic positions from one ligand to the next. A good docking code, such as PRO_–LEADS, can usually dock successfully if there is induced fit in relatively rigid enzymes but there remains the need to develop improved strategies for dealing with enzyme flexibility. The work implies that treatments of enzyme flexibility which focus only on sidechain rotations will not deal with the critical shifts responsible for mis-docking of ligands in thrombin, thermolysin and neuraminidase. The paper demonstrates the utility of all pairs docking experiments as a method of assessing the effectiveness of docking methodologies in dealing with enzyme flexibility.     

Targeting the RdRp of Emerging RNA Viruses: The Structure-Based Drug Design Challenge

The RNA-dependent RNA polymerase (RdRp) is an essential enzyme for the viral replication process, catalyzing the viral RNA synthesis using a metal ion-dependent mechanism. In recent years, RdRp has emerged as an optimal target for the development of antiviral drugs, as demonstrated by recent approvals of sofosbuvir and remdesivir against Hepatitis C virus (HCV) and severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), respectively. In this work, we overview the main sequence and structural features of the RdRp of emerging RNA viruses such as Coronaviruses, Flaviviruses, and HCV, as well as inhibition strategies implemented so far. While analyzing the structural information available on the RdRp of emerging RNA viruses, we provide examples of success stories such as for HCV and SARS-CoV-2. In contrast, Flaviviruses’ story has raised attention about how the lack of structural details on catalytically-competent or ligand-bound RdRp strongly hampers the application of structure-based drug design, either in repurposing and conventional approaches.     

Binding-affinity predictions of HSP90 in the D3R Grand Challenge 2015 with docking,MM/GBSA,QM/MM,and free-energy simulations

We have estimated the binding affinity of three sets of ligands of the heat-shock protein 90 in the D3R grand challenge blind test competition. We have employed four different methods, based on five different crystal structures: first, we docked the ligands to the proteins with induced-fit docking with the Glide software and calculated binding affinities with three energy functions. Second, the docked structures were minimised in a continuum solvent and binding affinities were calculated with the MM/GBSA method (molecular mechanics combined with generalised Born and solvent-accessible surface area solvation). Third, the docked structures were re-optimised by combined quantum mechanics and molecular mechanics (QM/MM) calculations. Then, interaction energies were calculated with quantum mechanical calculations employing 970–1160 atoms in a continuum solvent, combined with energy corrections for dispersion, zero-point energy and entropy, ligand distortion, ligand solvation, and an increase of the basis set to quadruple-zeta quality. Fourth, relative binding affinities were estimated by free-energy simulations, using the multi-state Bennett acceptance-ratio approach. Unfortunately, the results were varying and rather poor, with only one calculation giving a correlation to the experimental affinities larger than 0.7, and with no consistent difference in the quality of the predictions from the various methods. For one set of ligands, the results could be strongly improved (after experimental data were revealed) if it was recognised that one of the ligands displaced one or two water molecules. For the other two sets, the problem is probably that the ligands bind in different modes than in the crystal structures employed or that the conformation of the ligand-binding site or the whole protein changes.     

Combining self- and cross-docking as benchmark tools: the performance of DockBench in the D3R Grand Challenge 2

Molecular docking is a powerful tool in the field of computer-aided molecular design. In particular, it is the technique of choice for the prediction of a ligand pose within its target binding site. A multitude of docking methods is available nowadays, whose performance may vary depending on the data set. Therefore, some non-trivial choices should be made before starting a docking simulation. In the same framework, the selection of the target structure to use could be challenging, since the number of available experimental structures is increasing. Both issues have been explored within this work. The pose prediction of a pool of 36 compounds provided by D3R Grand Challenge 2 organizers was preceded by a pipeline to choose the best protein/docking-method couple for each blind ligand. An integrated benchmark approach including ligand shape comparison and cross-docking evaluations was implemented inside our DockBench software. The results are encouraging and show that bringing attention to the choice of the docking simulation fundamental components improves the results of the binding mode predictions.     

Exhaustive docking and solvated interaction energy scoring: lessons learned from the SAMPL4 challenge

We continued prospective assessments of the Wilma–solvated interaction energy (SIE) platform for pose prediction, binding affinity prediction, and virtual screening on the challenging SAMPL4 data sets including the HIV-integrase inhibitor and two host–guest systems. New features of the docking algorithm and scoring function are tested here prospectively for the first time. Wilma–SIE provides good correlations with actual binding affinities over a wide range of binding affinities that includes strong binders as in the case of SAMPL4 host–guest systems. Absolute binding affinities are also reproduced with appropriate training of the scoring function on available data sets or from comparative estimation of the change in target’s vibrational entropy. Even when binding modes are known, SIE predictions lack correlation with experimental affinities within dynamic ranges below 2 kcal/mol as in the case of HIV-integrase ligands, but they correctly signaled the narrowness of the dynamic range. Using a common protein structure for all ligands can reduce the noise, while incorporating a more sophisticated solvation treatment improves absolute predictions. The HIV-integrase virtual screening data set consists of promiscuous weak binders with relatively high flexibility and thus it falls outside of the applicability domain of the Wilma–SIE docking platform. Despite these difficulties, unbiased docking around three known binding sites of the enzyme resulted in over a third of ligands being docked within 2 Å from their actual poses and over half of the ligands docked in the correct site, leading to better-than-random virtual screening results.     

Lessons learned from the design of chemical space networks and opportunities for new applications

The concept of chemical space is of fundamental relevance in chemical informatics and computer-aided drug discovery. In a series of articles published in the Journal of Computer-Aided Molecular Design, principles of chemical space design were evaluated, molecular networks proposed as an alternative to conventional coordinate-based chemical reference spaces, and different types of chemical space networks (CSNs) constructed and analyzed. Central to the generation of CSNs was the way in which molecular similarity relationships were assessed and a primary focal point was the network-based representation of biologically relevant chemical space. The design and comparison of CSNs based upon alternative similarity measures can be viewed as an evolutionary path with interesting lessons learned along the way. CSN design has matured to the point that such chemical space representations can be used in practice. In this contribution, highlights from the sequence of CSN design efforts are discussed in context, providing a perspective for future practical applications.