Application of Shape Similarity in Pose Selection and Virtual

Aug 6, 2015 - Fax: +81-45-503-9559. ... To evaluate the applicability of shape similarity in docking-based pose selection and virtual screening, we pa...
0 downloads 4 Views 1MB Size
Article pubs.acs.org/jcim

Application of Shape Similarity in Pose Selection and Virtual Screening in CSARdock2014 Exercise Ashutosh Kumar and Kam Y. J. Zhang* Structural Bioinformatics Team, Center for Life Science Technologies, RIKEN, 1-7-22 Suehiro, Tsurumi, Yokohama, Kanagawa 230-0045, Japan S Supporting Information *

ABSTRACT: To evaluate the applicability of shape similarity in dockingbased pose selection and virtual screening, we participated in the CSARdock2014 benchmark exercise for identifying the correct docking pose of inhibitors targeting factor XA, spleen tyrosine kinase, and tRNA methyltransferase. This exercise provides a valuable opportunity for researchers to test their docking programs, methods, and protocols in a blind testing environment. In the CSARdock2014 benchmark exercise, we have implemented an approach that uses ligand 3D shape similarity to facilitate docking-based pose selection and virtual screening. We showed here that ligand 3D shape similarity between bound poses could be used to identify the native-like pose from an ensemble of docking-generated poses. Our method correctly identified the native pose as the top-ranking pose for 73% of test cases in a blind testing environment. Moreover, the pose selection results also revealed an excellent correlation between ligand 3D shape similarity scores and RMSD to X-ray crystal structure ligand. In the virtual screening exercise, the average RMSD for our pose prediction was found to be 1.02 Å, and it was one of the top performances achieved in CSARdock2014 benchmark exercise. Furthermore, the inclusion of shape similarity improved virtual screening performance of docking-based scoring and ranking. The coefficient of determination (r2) between experimental activities and docking scores for 276 spleen tyrosine kinase inhibitors was found to be 0.365 but reached 0.614 when the ligand 3D shape similarity was included.



INTRODUCTION

Molecular docking is one of the most common virtual screening approaches and plays a central role in identification of initial hits or optimization of hit compounds in many inhibitor discovery studies.17−19 Molecular docking predicts the conformation and orientation of ligands within the target protein binding site. In the process of ligand binding prediction, hundreds to thousands of poses are generated that are subsequently ranked by scoring functions. A reliable scoring function would rank near native binding modes higher than others and could also prioritize active ligands over inactives. Over the years, several scoring functions were developed for molecular docking, which can be broadly classified into empirical, force field-based, and knowledge-based.17−19,26,27 However, the performance of these scoring functions varies, and some fared better than others in benchmark studies.28−31 Consensus scoring methodologies were also suggested,32−36 leading to the development of consensus scoring functions. Evaluation of ligand 3D shape similarity is one of the commonly used ligand-based virtual screening approaches.37−43 Ligand 3D shape similarity could be used as a virtual screening approach to identify initial hits,44−46 to hop from one chemical scaffold to another,47,48 and as a means to derive structure activity relationships for inhibitor optimization.49,50 Apart from

In the last two decades, virtual screening of chemical libraries has become an important method to identify initial hits in modern drug discovery.1−6 These initial hits are later optimized to improve potency, selectivity, bioavailability, and metabolic stability. Several virtual screening methodologies are in common use now and can be broadly categorized into ligand and structure-based virtual screening approaches.7−10 Ligandbased virtual screening approaches rely only on the knowledge derived from at least one active ligand. Ligand-based methods can work in the absence of structural information for the target protein. These approaches include two-dimensional or threedimensional (2D or 3D) similarity searches,11,12 ligand-based pharmacophore modeling,8,13 quantitative structure activity relationships (QSAR),8,14 machine learning,15,16 etc. Structure-based virtual screening methods require structural knowledge of target proteins. Molecular docking,17−19 structurebased pharmacophore modeling,20,21 and de novo design22,23 are most commonly used structure-based virtual screening methods. Although considerable success has been achieved in ligand- and structure-based virtual screening, shortcomings and problems still exist with both approaches.24,25 To overcome the shortcomings of individual ligand-based or structure-based methods, integrated approaches have been developed by combining these two methods in either parallel or hierarchical manner.1,3,7 © XXXX American Chemical Society

Special Issue: Community Structure Activity Resource (CSAR) Received: May 14, 2015

A

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

Pose Selection Methodology. The correct binding pose of 22 ligands (3 for FXA, 5 for SYK, and 14 for TRMD) from a set of docking-generated decoys was selected using ligand 3D shape similarity. All available inhibitor-bound crystal structures were retrieved for these three target proteins from the RCSB protein data bank.57 Small molecule inhibitor-bound crystal structures for FXA, SYK, and TRMD were downloaded and prepared by removing all the heteroatoms except the bound inhibitor (Supporting Information, Table S1). The inhibitorbound crystal structures were then superimposed on FXA, SYK, and TRMD proteins from the CSARdock2014 benchmark exercise data set. ROCS42 was then used to calculate in-place shape similarity between each docking-generated decoy and a known crystal structure ligand. The calculations were performed using the “scoreonly” option of ROCS, which turns off the alignment and optimization and scores the in-place poses only. TanimotoCombo was used to measure 3D shape similarity between the two ligands. The ranking of 200 dockinggenerated decoys was prepared by selecting the best and average 3D shape similarity scores between each crystal structure ligand and docking decoy. Pose Prediction and Virtual Screening Methodology. In this phase of the CSARdock2014 benchmark exercise (phase 2), rank-ordering of ligands in the phase 2 data set was performed using a consensus approach by combining the results of docking with shape similarity. Here, ligand 3D shape similarity was only used to rank-order the compounds and did not play any role in pose selection or prediction. For pose prediction, docking of ligands in the phase 2 data set was performed using multiple crystal structures. The 22 crystal structures released after the completion of phase 1 were utilized (3 for FXA, 5 for SYK, and 14 for TRMD for respective ligands) for docking. The protein preparation utility of Schrodinger’s Maestro58 was employed to prepare receptor structures for molecular docking. Receptor structures were prepared by adding hydrogens, assigning bond orders, and by determining the protonation state of charged residues. LigPrep59 was employed to prepare ligands for molecular docking. Ligands were prepared by converting them to 3D format, and subsequently, hydrogens were added and ionization states and tautomers were generated. All structures were then subjected to minimization using OPLS-2005 force field.60,61 The binding poses of ligands in the phase 2 data set were predicted using the Glide program62−65 in extra precision mode. Docking was performed using multiple crystal structures, and only a single top-scoring pose from Glide-XP docking was retained. The pose with the best Glide-XP score from multiple receptor docking runs was used to perform shape similarity calculations between each of the data set ligands and 22 crystal structure ligands (3 for FXA, 5 for SYK, and 14 for TRMD) as described in the previous section. To rank-order phase 2 data set ligands, docking and shape similarity scores were combined using the sum score data fusion ranking method.66,67 In sum score data fusion methodology, the relative docking and ligand 3D shape similarity scores for each ligand was calculated by dividing them by the highest score achieved by any compound. The new scores were then combined to give a hybrid score, which was used to prepare rankings for CSARdock phase 2 data set ligands.

these applications, shape similarity can also be used as a scoring scheme in molecular docking to rank/score sampled poses if at least one active ligand is available. One such approach is implemented in the HYBRID program,51 where ligand 3D shape and chemical similarity with receptor-bound ligands is utilized to dock input ligands. An improved methodology is available as POSIT program52 that uses both protein and known ligand structural information to predict ligand poses. Although the use of ligand 3D shape similarity as a scoring scheme seems to be an attractive option, we are unaware of its application in a real world scenario. Moreover, unbiased evaluation of any computational approach would be in either real life applications or a blind testing environment. The CSARdock exercise is an annual event conducted by the Community Structure Activity Resource (CSAR, http://www. csardock.org) in association with various academic groups and pharmaceutical companies. The CSARdock exercise provides a platform for researchers to test their docking and scoring methods in a blind environment. So far, four rounds of CSARdock exercises have been completed, 53−55 which generated immense interest in the scientific community to share and learn from their experiences for the testing and development of docking and scoring methods. The CSARdock2014 benchmark exercise is based on the unreleased crystal structures and binding affinity data donated by GlaxoSmithKline (GSK) and is staged in two phases: pose selection and virtual screening. In the first phase, participants have to identify the correct binding mode from a set of 200 docking generated decoys for each ligand in the data set. The second phase involves the prediction of docking poses and ranking of binding affinities of several series of inhibitors for three target proteins.



MATERIALS AND METHODS Data Set. In the CSARdock2014 benchmark exercise, two data sets were released to participants: (1) CSARdock2014 benchmark exercise phase 1 data set and (2) CSARdock2014 benchmark exercise phase 2 data set. To maintain the blind nature of the exercise, the second data set was released only after the completion of phase 1. (1). CSARdock2014 Benchmark Exercise Phase 1 Data Set. For the pose selection challenge, 200 docking generated decoys for each of the 22 ligand−protein pairs (3 for factor Xa (FXA), 5 for spleen tyrosine kinase (SYK), and 14 for tRNA methyltransferase (TRMD)) were released to the participants. The decoys were generated using the Dock56 program and contained at least one correct pose. Twenty-two co-crystal structures from three proteins were released only after the completion of phase 1 to facilitate the participants’ own analyses. (2). CSARdock2014 Benchmark Exercise Phase 2 Data Set. The phase 2 data set consists of three congeneric series of inhibitors for three target proteins (163 for FXA, 276 for SYK, and 31 for TRMD). The participants received these ligands as SMILES strings. Docking and ranking/scoring was required using any of the corresponding crystal structures. The biological activity and 20 co-crystal structures of the phase 2 data set ligands from three proteins (1 for FXA, 5 for SYK, and 14 for TRMD) were released to participants after the closing of this phase of competition to facilitate their analyses of pose prediction and virtual screening results. These 20 co-crystal structures were different from those released after phase 1.



RESULTS AND DISCUSSION In this CSARdock2014 exercise, we attempted to extend the utility of ligand 3D shape similarity evaluations to docking pose B

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

Figure 1. Scatterplots of ROCS shape similarity score (TanimotoCombo) vs RMSD to crystal structure ligand of 200 molecular docking-generated ligand poses for (A) FXA, (B) SYK, and (C) TRMD.

Figure 2. Scatterplots of ROCS shape similarity score (TanimotoCombo) vs RMSD to crystal structure ligand of ligand poses that are 5 Å away for (A) FXA, (B) SYK, and (C) TRMD.

Figure 3. Ranking of ligands in CSARdock2014 pose selection using the best and average ligand 3D shape similarity scores.

near-native binding poses higher than poses far away from the native. Through the participation of the CSARdock2014 exercise, we sought to answer the question as to whether ligand 3D shape similarity alone can be used to accurately rankorder binding poses. We hypothesized that near native binding poses will share high 3D shape similarity with crystal structure ligand than the far native ones. To address our question, a retrospective analysis was carried out utilizing the 22 protein− ligand crystal structures released after the completion of phase 1. The ligand 3D shape similarity of binding poses in the data set was calculated with target protein co-crystal structures (Supporting Information, Table S1) as described in the

selection and virtual screening. We found that ligand 3D shape similarity matching is a feasible pose selection approach when at least one ligand-bound crystal structure is available. We also observed improved virtual screening/affinity prediction performance when docking-based scoring was complemented with ligand 3D shape similarity in suitable cases. Ligand 3D Shape Similarity Can Be Used for Pose Selection. Ranking of docking-generated poses is a challenging problem as sampling algorithms can generate binding poses close to the correct solution, but identifying them as the top or within the top few using scoring functions has been difficult.17,18 A reliable scoring function would be able to rank C

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

one or two ligands. Specifically, TRMD-451 was a difficult case to predict, and the majority of approaches including ours failed to identify the native pose. Our approach failed in this case because a part of TRMD-451 extends to a previously unexplored subpocket by crystal structure ligands and shape overlap could not be obtained. Similarly, in the case of the SYK protein, our approach could only identify two out of five correct poses as the top solution (although all native poses were identified in top five) because the indazole ring in SYK ligands extended toward a subpocket formed by Lys402, Glu420, Leu446, and Met448 amino acid residues. Ligand 3D shape similarity score due to this pocket was lower, and ranking of the native pose was slightly curtailed. To elucidate the reason for the tremendous performance by molecular docking scoring functions and other methods, we analyzed the RMSD of the best pose with crystal structure ligands. As expected, RMSD was very low (0.2 to 0.77 Å, average = 0.38 Å), and conformation of the best pose is very near to that of the crystal structure ligand (Supporting Information, Figure S2). Furthermore, the next best pose was more than 2 Å away, and median RSMD to crystal structure ligand was more than 8 Å for a majority of the data set ligands. It is well known that docking programs and scoring function could recapitulate crystal structure conformation very well when the same starting conformation is used for docking. However, their performance drastically decreases when either apo structure or a different crystal structure is used for docking. Ligand 3D shape similarity-based ranking, however, only indirectly uses the protein−ligand interaction information. Hence, the ranking/ scoring of docking poses is free from errors associated with the charge of ligands and proteins, protonation state of residues, metal ions, etc. The only requirement for the ligand 3D shape similarity-based ranking approach is the availability of at least one ligand-bound crystal structure targeting the same binding pocket to facilitate the evaluation of its 3D shape similarity with the ligand poses to be scored. Despite the above-mentioned advantages, the use of ligand 3D shape similarity for pose selection is limited to cases where sufficient ligand conformational diversity exists. It is not suitable for small fragment-like ligands or other ligands with few rotatable bonds where 3D shape differences are limited. Volume-based clustering of ligands could be easily used to identify the applicability domain of our approach. Furthermore, as our approach relies on 3D shape of crystal ligands, predicting poses of ligands that extend toward previously unexplored subpockets in a binding site would not be reliable due to the absence of shape overlap for this region. Ligand 3D Shape Similarity Contributes to Virtual Screening Performance. Another property of a good scoring/ranking approach is that it should be able to rankorder biological affinities of ligands. Furthermore, it should be able to discriminate active from inactive ligands. However, prediction of biological affinities by docking scoring functions is a challenging task, and very little success has been achieved in this area.17−19 One important prerequisite for accurate rank ordering of biological affinities is the ability of molecular docking to generate ligand-binding modes as close as possible to the native pose. The failure to generate reliable binding modes would cause inaccurate rank ordering of biological activities. To analyze the accuracy of binding modes of FXA, SYK, and TRMD ligands in the CSARdock2014 virtual screening challenge, the reproducibility of crystal structure binding modes was evaluated. The CSARdock2014 exercise

previous section. The ligand 3D shape similarity score (here TanimotoCombo from ROCS program) was plotted against root-mean-square deviation (RSMD) to crystal structure ligand (Figure 1). Although we have used TanimotoCombo, a combined similarity score for shape and chemical feature overlap, the contribution from chemical feature similarity was negligible for most poses. As ligand 3D shape similarity was calculated “in place” without aligning two ligands, chemical feature similarity provided notable contribution only for poses close to native ones. The 3D ligand shape similarity score correlates very well with RMSD to crystal structure ligand with average Spearman correlation coefficients (ρ) of 0.882, 0.862, and 0.917 for FXA, SYK, and TRMD ligands, respectively (Figure 1 and Supporting Information, Table S2). A similar trend was observed with other correlation coefficients as well (Supporting Information, Table S2). Although the performance of ligand 3D shape similarity was excellent, this was no surprise since most of the poses were far way from the true solution (32.1%, 45.8%, and 60.5% beyond 10 Å for FXA, SYK, and TRMD, respectively, Figure 1). To rule out data distribution related bias, only poses with less than 5 Å RMSD were plotted, and as shown in Figure 2, good correlation was also observed with average ρ values of 0.682, 0.510, and 0.732 for FXA, SYK, and TRMD ligands, respectively (Supporting Information, Table S2). The other correlation coefficients also follow a similar trend (Supporting Information, Table S2). It seems reasonable to expect that ligand 3D shape similarity could be used to discriminate near native binding poses from incorrect ones, but it is uncertain whether this approach will work in a blind testing environment. Our result in the CSARdock2014 benchmark exercise pose selection challenge has addressed this question. We employed the best and average ligand 3D shape similarity scores to rank-order dockinggenerated binding poses, and the results for our submission are presented in Figure 3. The usage of the best ligand 3D shape similarity score outperformed the average ligand 3D shape similarity score with 16 out of 22 ligands (72.7%) correctly identified as the top solution and another 5 in the top 5 solutions (Figure 3). Overall, the correct pose was identified within the top 5 rank in 21 out of 22 ligands (95.4%), while missing the correct pose for only one ligand (TRMD-451, ninth in rank). However, when the ligand 3D shape similarity scores were averaged, the performance was worse, and the correct pose was identified in only 4 out of 22 cases (18.1%). There was no significant improvement even when the top five ranks were considered, and the correct pose for only nine more ligands could be identified (Figure 3). This result reveals that if several co-crystal structures are available, then the crystal ligand with the highest 3D shape similarity to the test ligands should be used since the averaging of 3D shape similarities led to reduced performance. We next checked how our approach fared as a binding pose scoring/ranking scheme when compared with other participants/methods. Seven out of 52 methods performed extremely well and correctly identified the native pose from docking generated decoys in all the cases.68 Another eight methods correctly identified the closest RMSD pose in the top three solutions. Our own comparison with other force field-based, empirical-, and knowledge-based scoring functions used in molecular docking (Autodock-Vina,69 DrugScore,70 Chemgauss3,71 PLP score,72 Chemscore73,74) revealed extremely good results (Supporting Information, Figure S1). Most scoring functions accurately identified the correct pose for all except D

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling Table 1. Overview of Our Pose Prediction Performance Using CSARdock2014 Benchmark Exercise Phase 2 Data Set CSAR quality crystal structures

all crystal structures

target protein

total number of ligands docked

number of crystal structures

median RMSD

mean RMSD

number of crystal structures

median RMSD

mean RMSD

FXA SYK TRMD total/averagea

163 276 31 470

1 5 14 20

0.496 0.293 0.769 0.636

0.496 0.877 1.114 1.024

5 9 31 45

2.684 0.557 0.947 1.062

1.995 0.770 1.417 1.351

a

Median and mean RMSD values for each protein family were provided by the CSARdock2014 organizers. The average RMSD value for all protein families was the sum of the values for each protein family weighted by the number of crystal structures divided by the total number of crystal structures.

Figure 4. Comparison of the pose prediction performance of CSARdock2014 blind exercise participants for CSARdock2014 phase 2 data set ligands. Twenty unreleased crystal structures (1, 5, and 14 for FXA, SYK, and TRMD, respectively) were used to assess the prediction accuracy. Group-B represents our prediction.

Table 2. Correlation between Experimental Activities and Our Scoring Using CSARdock2014 Benchmark Exercise Phase 2 Data Set molecular docking only target proteins

number of inhibitors

coefficient of determination (r2)

FXA SYK TRMD

163 276 31

0 0.366 0.138

hybrid scoring

Pearson r Spearman ρ 0.021 0.605 0.372

−0.097 0.617 0.472

Kendall τ

coefficient of determination (r2)

Pearson r

Spearman ρ

Kendall τ

−0.061 0.436 0.325

0.018 0.614 0.032

0.135 0.784 0.179

0.049 0.782 0.339

0.028 0.590 0.266

Figure 5. Correlation plots between the biological activities and (A) docking score only and (B) hybrid score for 276 inhibitors of SYK protein of the CSARdock2014 blind exercise data set. Docking scores represent scores obtained from Glide-XP docking, while hybrid scores are the concatenation of sum scores of docking and ROCS shape similarity TanimotoCombo.

organizers provided 20 previously unreleased CSAR quality crystal structures (1, 5, and 14 from FXA, SYK, and TRMD, respectively) after the closing of the virtual screening phase. The CSAR quality crystal structures were with exceptional

ligand density. These 20 crystal structures were used to evaluate the performance of our docking in generating reliable binding modes for virtual screening. Additional ligand-bound crystal structures with regular ligand density were also made available. E

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling The performance was measured by the RMSD of the ligandbinding pose with the corresponding crystal structure pose. The performance of our docking in reproducing the binding modes of 1, 5, and 14 CSAR quality crystal structures from FXA, SYK, and TRMD, respectively, is presented in Table 1. Overall, an exceptional performance was achieved with median RMSD values of 0.496, 0.293, and 0.769 Å for FXA, SYK, and TRMD, respectively. The average median and mean RMSD values were found to be 0.636 and 1.024, respectively. The performance decreased when those low quality crystal structures were also included in the analysis. Although slightly higher average median and mean RMSD values of 1.062 and 1.351 were observed, these were still in a reasonable range for dockinggenerated ligand poses. To facilitate comparison with other submissions, CSARdock2014 exercise organizers also shared the data from various participants of virtual screening. As some submissions included the binding mode prediction for only one protein, for fair comparison, we only included the submission where binding modes for ligands belonging to all three proteins were predicted. Our docking performed very well (Group B) and was among the top performers with RMSD less than 1 Å for 70% of the ligands (Figure 4) (100%, 80%, and 64.2% for FXA, SYK, and TRMD, respectively). The median RMSD value of our prediction was 0.636 Å, which is second only to “GroupS” among all the participants. Here, in phase 2 of the exercise, we sought to assess the effect of ligand 3D shape similarity in complementing docking scoring functions for virtual screening. Therefore, 163 FXA, 276 SYK, and 31 TRMD inhibitors were rank ordered using molecular docking-based scoring and hybrid scoring (a combination of docking and ligand 3D shape similarity scores using the sum score data fusion approach) as discussed previously and submitted as “Groups B-1 and B-2”. Correlation was calculated between experimental activities (−logIC50) and our molecular docking-based scoring and hybrid scoring (Table 2). The most impressive result for our scoring strategy involving ligand 3D shape similarity was obtained for 276 inhibitors of SYK protein (Table 2 and Figure 5). Our molecular docking has produced admirable binding modes for 276 SYK inhibitors, and as a result, a reasonable r2 value of 0.366 was obtained when only molecular docking scores were employed for rank ordering. However, the r2 value improved almost 2-fold from 0.366 to 0.614 when hybrid scoring was used by including the ligand 3D shape similarity term (Table 2 and Figure 5). When the performance of our hybrid approach was compared with other participants/methods, our scoring scheme was the top performer among all the participants in predicting the biological activity of SYK inhibitors (Figure 6). One reason for the excellent performance of hybrid scoring may be the use of five co-crystal structures of SYK protein (released after phase 1) for ligand 3D shape similarity calculations with phase 2 data set ligands. As both crystal structure ligands and phase 2 data set ligands belong to the same congeneric series, higher shape overlap could be obtained. Another reason for the excellent performance of hybrid scoring may be the identification of near native docking pose as the top pose by Glide, which is exemplified by 0.293 Å median RMSD for five SYK ligands for which co-crystal structures were available. Further, the availability of sufficient diversity in 3D shape of SYK ligands assured top performance. Moderate correlation was observed for 31 inhibitors of TRMD protein when using only molecular docking-based scoring (Table 2). However, no improvement in the correlation was observed upon the inclusion of ligand 3D

Figure 6. Comparison of the virtual screening performance (activity prediction) of CSARdock2014 blind exercise participants for 276 inhibitors of SYK protein. GroupB-1 and GroupB-2 represent our prediction, while employing docking scores and hybrid scores respectively for rank ordering of 276 ligands.

shape similarity with molecular docking-based scoring. Although CSAR provided 14 co-crystal structures for the TRMD ligand, the performance of 3D shape similarity-based ranking was not up to the task for the virtual screening challenge. Most of the TRMD ligands (including the crystal structure ligand) were small fragment-like ligands. The 3D shapes of these ligands were very similar to each other. As our approach relies on 3D shape diversity to rank ligands, it could not perform well. In comparison with TRMD ligands, 276 SYK ligands presented sufficient 3D shape diversity to differentiate between active and inactive ligands, and hence, the performance was significantly better. Though several groups were successful in accurately predicting the biological activities (Supporting Information, Figure S3), our method is not suitable for ligands with very low ligand 3D shape diversity. We were surprised to see absolutely no correlation between biological activities and hybrid scoring for 163 inhibitors of FXA protein. A similar trend was observed for ranking based on only docking scores. In fact, all of the participants performed equally poor in predicting the activities of FXA protein (Supporting Information, Figure S4). It seems that FXA is a challenging system for scoring ligands by virtue of a number of missing domains in FXA crystal structures. These missing domains may be involved in ligand binding and make the scoring of ligands problematic. The difficulty in scoring FXA ligands has been previously observed in the CSARdock2010 exercise.53



CONCLUSION Through the participation in the CSARdock2014 exercise, we have shown that ligand 3D shape similarities could be used in docking-based pose selection and virtual screening. Ligand 3D similarity matching is a valid scoring approach and can be used to score ligand poses generated by a docking program or any other approach. Although ligand 3D shape within a binding pocket is a property of the binding pocket shape and protein− ligand interaction information, the ligand 3D shape similaritybased scoring approach does not involve any direct calculation of protein−ligand interactions. Hence, the scoring of poses is not influenced by factors like charge, protonation state, metal ions, etc. However, at least one ligand-bound crystal structure targeting the same binding pocket is required to facilitate the evaluation of its 3D shape similarity with the ligand poses to be scored. The bound small molecule does not have to be chemically similar. As long as the 3D shape of bound molecule F

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling

(8) Sliwoski, G.; Kothiwale, S.; Meiler, J.; Lowe, E. W. Computational Methods in Drug Discovery. Pharmacol. Rev. 2013, 66, 334−395. (9) Sukumar, N.; Das, S. Current Trends in Virtual High Throughput Screening Using Ligand-Based and Structure-Based Methods. Comb. Chem. High Throughput Screening 2011, 14, 872−888. (10) Fukunishi, Y. Structure-Based Drug Screening and Ligand-Based Drug Screening with Machine Learning. Comb. Chem. High Throughput Screening 2009, 12, 397−408. (11) Maggiora, G.; Vogt, M.; Stumpfe, D.; Bajorath, J. Molecular Similarity in Medicinal Chemistry. J. Med. Chem. 2014, 57, 3186− 3204. (12) Stumpfe, D.; Bajorath, J. Similarity Searching. Wiley Interdisciplinary Reviews: Computational Molecular Science 2011, 1, 260−282. (13) Caporuscio, F.; Tafi, A. Pharmacophore Modelling: A Forty Year Old Approach and Its Modern Synergies. Curr. Med. Chem. 2011, 18, 2543−2553. (14) Guha, R. On Exploring Structure-Activity Relationships. Methods Mol. Biol. 2013, 993, 81−94. (15) Melville, J. L.; Burke, E. K.; Hirst, J. D. Machine Learning in Virtual Screening. Comb. Chem. High Throughput Screening 2009, 12, 332−343. (16) Gertrudes, J. C.; Maltarollo, V. G.; Silva, R. A.; Oliveira, P. R.; Honorio, K. M.; da Silva, A. B. Machine Learning Techniques and Drug Design. Curr. Med. Chem. 2012, 19, 4289−4297. (17) Yuriev, E.; Agostino, M.; Ramsland, P. A. Challenges and Advances in Computational Docking: 2009 in Review. J. Mol. Recognit. 2011, 24, 149−164. (18) Yuriev, E.; Ramsland, P. A. Latest Developments in Molecular Docking: 2010−2011 in Review. J. Mol. Recognit. 2013, 26, 215−239. (19) Yuriev, E.; Holien, J.; Ramsland, P. A. Improvements, Trends, and New Ideas in Molecular Docking: 2012−2013 in Review. J. Mol. Recognit. 2015, n/a. (20) Pirhadi, S.; Shiri, F.; Ghasemi, J. B. Methods and Applications of Structure Based Pharmacophores in Drug Discovery. Curr. Top. Med. Chem. 2013, 13, 1036−1047. (21) Sanders, M. P. A.; McGuire, R.; Roumen, L.; de Esch, I. J. P.; de Vlieg, J.; Klomp, J. P. G.; de Graaf, C. From the Protein’s Perspective: The Benefits and Challenges of Protein Structure-Based Pharmacophore Modeling. MedChemComm 2012, 3, 28−38. (22) Kutchukian, P. S.; Shakhnovich, E. I. De Novo Design: Balancing Novelty and Confined Chemical Space. Expert Opin. Drug Discovery 2010, 5, 789−812. (23) Loving, K.; Alberts, I.; Sherman, W. Computational Approaches for Fragment-Based and De Novo Design. Curr. Top. Med. Chem. 2010, 10, 14−32. (24) Scior, T.; Bender, A.; Tresadern, G.; Medina-Franco, J. L.; Martínez-Mayorga, K.; Langer, T.; Cuanalo-Contreras, K.; Agrafiotis, D. K. Recognizing Pitfalls in Virtual Screening: A Critical Review. J. Chem. Inf. Model. 2012, 52, 867−881. (25) Schneider, G. Virtual Screening: An Endless Staircase? Nat. Rev. Drug Discovery 2010, 9, 273−276. (26) Huang, S.-Y.; Grinter, S. Z.; Zou, X. Scoring Functions and Their Evaluation Methods for Protein-Ligand Docking: Recent Advances and Future Directions. Phys. Chem. Chem. Phys. 2010, 12, 12899−12908. (27) Wang, J. C.; Lin, J. H. Scoring Functions for Prediction of Protein-Ligand Interactions. Curr. Pharm. Des. 2013, 19, 2174−2182. (28) Li, Y.; Han, L.; Liu, Z.; Wang, R. Comparative Assessment of Scoring Functions on an Updated Benchmark: 2. Evaluation Methods and General Results. J. Chem. Inf. Model. 2014, 54, 1717−1736. (29) Li, Y.; Liu, Z.; Li, J.; Han, L.; Liu, J.; Zhao, Z.; Wang, R. Comparative Assessment of Scoring Functions on an Updated Benchmark: 1. Compilation of the Test Set. J. Chem. Inf. Model. 2014, 54, 1700−1716. (30) Cheng, T.; Li, X.; Li, Y.; Liu, Z.; Wang, R. Comparative Assessment of Scoring Functions on a Diverse Test Set. J. Chem. Inf. Model. 2009, 49, 1079−1093.

and compound to be scored are similar, the ligand 3D shape similarity-based scoring is going to produce good results. Ligand 3D shape similarities can also complement molecular docking and could be used as a virtual screening tool. Although ligand 3D shape similarity is a simple and fast pose prediction approach, it is limited to cases where sufficient ligand conformational diversity exists. It is not suitable for small fragment-like ligands or ligands whose activities are driven by specific functional groups such as cysteine/serine protease inhibitors.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jcim.5b00279. Four additional figures: (1) Comparison of ligand 3D shape similarity-based ranking with docking-based scoring functions for CSARdock2014 benchmark exercise phase 1 data set. (2) Distribution of RMSD of docking poses CSARdock2014 benchmark exercise phase 1 ligands with their respective crystal structure ligand. (3) Comparison of the virtual screening performance (activity prediction) of CSARdock2014 blind exercise participants for 163 inhibitors of FXA protein. (4) Comparison of the virtual screening performance (activity prediction) of CSARdock2014 blind exercise participants for 31 inhibitors of TRMD protein. Two additional tables: (1) Ligand-bound crystal structures used for ligand 3D shape similarity calculations. (2) Correlation between ligand 3D shape similarity and RMSD to crystal structure ligand. (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +81-45-503-9560. Fax: +81-45-503-9559. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We thank RIKEN Integrated Cluster of Clusters (RICC) at RIKEN for the supercomputing resources used for this study. We thank members of our lab for help and discussions.



REFERENCES

(1) Tanrikulu, Y.; Krüger, B.; Proschak, E. The Holistic Integration of Virtual Screening in Drug Discovery. Drug Discovery Today 2013, 18, 358−364. (2) Walters, W. P.; Stahl, M. T.; Murcko, M. A. Virtual Screening an Overview. Drug Discovery Today 1998, 3, 160−178. (3) Kumar, A.; Zhang, K. Y. J. Hierarchical Virtual Screening Approaches in Small Molecule Drug Discovery. Methods 2015, 71, 26−37. (4) Lavecchia, A.; Di Giovanni, C. Virtual Screening Strategies in Drug Discovery: A Critical Review. Curr. Med. Chem. 2013, 20, 2839− 2860. (5) Muegge, I. Synergies of Virtual Screening Approaches. Mini-Rev. Med. Chem. 2008, 8, 927−933. (6) Muegge, I.; Oloff, S. Advances in Virtual Screening. Drug Discovery Today: Technol. 2006, 3, 405−411. (7) Drwal, M. N.; Griffith, R. Combination of Ligand- and StructureBased Methods in Virtual Screening. Drug Discovery Today: Technol. 2013, 10, e395−e401. G

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling (31) Warren, G. L.; Andrews, C. W.; Capelli, A. M.; Clarke, B.; LaLonde, J.; Lambert, M. H.; Lindvall, M.; Nevins, N.; Semus, S. F.; Senger, S.; Tedesco, G.; Wall, I. D.; Woolven, J. M.; Peishoff, C. E.; Head, M. S. A Critical Assessment of Docking Programs and Scoring Functions. J. Med. Chem. 2006, 49, 5912−5931. (32) Liu, S.; Fu, R.; Zhou, L. H.; Chen, S. P. Application of Consensus Scoring and Principal Component Analysis for Virtual Screening against Beta-Secretase (Bace-1). PLoS One 2012, 7, e38086. (33) Zhong, S.; Zhang, Y.; Xiu, Z. Rescoring Ligand Docking Poses. Curr. Opin. Drug Discovery Dev. 2010, 13, 326−334. (34) Bar-Haim, S.; Aharon, A.; Ben-Moshe, T.; Marantz, Y.; Senderowitz, H. Selex-Cs: A New Consensus Scoring Algorithm for Hit Discovery and Lead Optimization. J. Chem. Inf. Model. 2009, 49, 623−633. (35) Teramoto, R.; Fukunishi, H. Supervised Consensus Scoring for Docking and Virtual Screening. J. Chem. Inf. Model. 2007, 47, 526− 534. (36) Oda, A.; Tsuchida, K.; Takakura, T.; Yamaotsu, N.; Hirono, S. Comparison of Consensus Scoring Strategies for Evaluating Computational Models of Protein-Ligand Complexes. J. Chem. Inf. Model. 2006, 46, 380−391. (37) Kalaszi, A.; Szisz, D.; Imre, G.; Polgar, T. Screen3d: A Novel Fully Flexible High-Throughput Shape-Similarity Search Method. J. Chem. Inf. Model. 2014, 54, 1036−1049. (38) Sato, T.; Yuki, H.; Takaya, D.; Sasaki, S.; Tanaka, A.; Honma, T. Application of Support Vector Machine to Three-Dimensional ShapeBased Virtual Screening Using Comprehensive Three-Dimensional Molecular Shape Overlay with Known Inhibitors. J. Chem. Inf. Model. 2012, 52, 1015−1026. (39) Liu, X.; Jiang, H.; Li, H. Shafts: A Hybrid Approach for 3d Molecular Similarity Calculation. 1. Method and Assessment of Virtual Screening. J. Chem. Inf. Model. 2011, 51, 2372−2385. (40) Hawkins, P. C. D.; Skillman, A. G.; Nicholls, A. Comparison of Shape-Matching and Docking as Virtual Screening Tools. J. Med. Chem. 2007, 50, 74−82. (41) Wei, N.-N.; Hamza, A. Sabre: Ligand/Structure-Based Virtual Screening Approach Using Consensus Molecular-Shape Pattern Recognition. J. Chem. Inf. Model. 2014, 54, 338−346. (42) Rocs 3.2.0.4. Openeye Scientific Software. http://www. eyesopen.com (accessed July 2015). (43) Vainio, M. J.; Puranen, J. S.; Johnson, M. S. Shaep: Molecular Overlay Based on Shape and Electrostatic Potential. J. Chem. Inf. Model. 2009, 49, 492−502. (44) Kumar, A.; Ito, A.; Takemoto, M.; Yoshida, M.; Zhang, K. Y. J. Identification of 1,2,5-Oxadiazoles as a New Class of Senp2 Inhibitors Using Structure Based Virtual Screening. J. Chem. Inf. Model. 2014, 54, 870−880. (45) Hevener, K. E.; Mehboob, S.; Su, P.-C.; Truong, K.; Boci, T.; Deng, J.; Ghassemi, M.; Cook, J. L.; Johnson, M. E. Discovery of a Novel and Potent Class of F. Tularensis Enoyl-Reductase (Fabi) Inhibitors by Molecular Shape and Electrostatic Matching. J. Med. Chem. 2012, 55, 268−279. (46) Naylor, E.; Arredouani, A.; Vasudevan, S. R.; Lewis, A. M.; Parkesh, R.; Mizote, A.; Rosen, D.; Thomas, J. M.; Izumi, M.; Ganesan, A.; Galione, A.; Churchill, G. C. Identification of a Chemical Probe for Naadp by Virtual Screening. Nat. Chem. Biol. 2009, 5, 220−226. (47) Vasudevan, S. R.; Singh, N.; Churchill, G. C. Scaffold Hopping with Virtual Screening from Ip3 to a Drug-Like Partial Agonist of the Inositol Trisphosphate Receptor. ChemBioChem 2014, 15, 2774− 2782. (48) Rush, T. S., 3rd; Grant, J. A.; Mosyak, L.; Nicholls, A. A ShapeBased 3-D Scaffold Hopping Method and Its Application to a Bacterial Protein-Protein Interaction. J. Med. Chem. 2005, 48, 1489−1495. (49) Pan, X.; Tan, N.; Zeng, G.; Huang, H.; Yan, H. 3d Qsar Studies on Ketoamides of Human Cathepsin K Inhibitors Based on Two Different Alignment Methods. Eur. J. Med. Chem. 2010, 45, 667−681. (50) Crisan, L.; Pacureanu, L.; Avram, S.; Bora, A.; Avram, S.; Kurunczi, L. Pls and Shape-Based Similarity Analysis of Maleimides– Gsk-3 Inhibitors. J. Enzyme Inhib. Med. Chem. 2014, 29, 599−610.

(51) Hybrid 3.0.1. Openeye Scientific Software. http://www. eyesopen.com (accessed July 2015). (52) Posit 3.1.0.5. Openeye Scientific Software. http://www. eyesopen.com (accessed July 2015).. (53) Smith, R. D.; Dunbar, J. B.; Ung, P. M.-U.; Esposito, E. X.; Yang, C.-Y.; Wang, S.; Carlson, H. A. Csar Benchmark Exercise of 2010: Combined Evaluation across All Submitted Scoring Functions. J. Chem. Inf. Model. 2011, 51, 2115−2131. (54) Damm-Ganamet, K. L.; Smith, R. D.; Dunbar, J. B., Jr.; Stuckey, J. A.; Carlson, H. A. Csar Benchmark Exercise 2011−2012: Evaluation of Results from Docking and Relative Ranking of Blinded Congeneric Series. J. Chem. Inf. Model. 2013, 53, 1853−1870. (55) Dunbar, J. B., Jr.; Smith, R. D.; Damm-Ganamet, K. L.; Ahmed, A.; Esposito, E. X.; Delproposto, J.; Chinnaswamy, K.; Kang, Y. N.; Kubish, G.; Gestwicki, J. E.; Stuckey, J. A.; Carlson, H. A. Csar Data Set Release 2012: Ligands, Affinities, Complexes, and Docking Decoys. J. Chem. Inf. Model. 2013, 53, 1842−1852. (56) Ewing, T. A.; Makino, S.; Skillman, A. G.; Kuntz, I. Dock 4.0: Search Strategies for Automated Molecular Docking of Flexible Molecule Databases. J. Comput.-Aided Mol. Des. 2001, 15, 411−428. (57) Bernstein, F. C.; Koetzle, T. F.; Williams, G. J.; Meyer, E. F., Jr.; Brice, M. D.; Rodgers, J. R.; Kennard, O.; Shimanouchi, T.; Tasumi, M. The Protein Data Bank: A Computer-Based Archival File for Macromolecular Structures. J. Mol. Biol. 1977, 112, 535−542. (58) Maestro, Version 9.2; Schrödinger, LLC., New York, 2011. (59) Ligprep, Version 2.5; Schrödinger, LLC., New York, 2011. (60) Jorgensen, W. L.; Tirado-Rives, J. The Opls [Optimized Potentials for Liquid Simulations] Potential Functions for Proteins, Energy Minimizations for Crystals of Cyclic Peptides and Crambin. J. Am. Chem. Soc. 1988, 110, 1657−1666. (61) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the Opls All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (62) Friesner, R. A.; Banks, J. L.; Murphy, R. B.; Halgren, T. A.; Klicic, J. J.; Mainz, D. T.; Repasky, M. P.; Knoll, E. H.; Shelley, M.; Perry, J. K.; Shaw, D. E.; Francis, P.; Shenkin, P. S. Glide: A New Approach for Rapid, Accurate Docking and Scoring. 1. Method and Assessment of Docking Accuracy. J. Med. Chem. 2004, 47, 1739−1749. (63) Friesner, R. A.; Murphy, R. B.; Repasky, M. P.; Frye, L. L.; Greenwood, J. R.; Halgren, T. A.; Sanschagrin, P. C.; Mainz, D. T. Extra Precision Glide: Docking and Scoring Incorporating a Model of Hydrophobic Enclosure for Protein-Ligand Complexes. J. Med. Chem. 2006, 49, 6177−6196. (64) Halgren, T. A.; Murphy, R. B.; Friesner, R. A.; Beard, H. S.; Frye, L. L.; Pollard, W. T.; Banks, J. L. Glide: A New Approach for Rapid, Accurate Docking and Scoring. 2. Enrichment Factors in Database Screening. J. Med. Chem. 2004, 47, 1750−1759. (65) Glide, Version 5.7; Schrödinger, LLC, New York, 2011. (66) Seal, A.; Yogeeswari, P.; Sriram, D.; Consortium, O.; Wild, D. J. Enhanced Ranking of Pknb Inhibitors Using Data Fusion Methods. J. Cheminf. 2013, 5, 2. (67) Sastry, G. M.; Inakollu, V. S.; Sherman, W. Boosting Virtual Screening Enrichments with Data Fusion: Coalescing Hits from TwoDimensional Fingerprints, Shape, and Docking. J. Chem. Inf. Model. 2013, 53, 1531−1542. (68) Community Structure-Activity Resource. http://CSARdock.org (accessed August 2015). (69) Trott, O.; Olson, A. J. Autodock Vina: Improving the Speed and Accuracy of Docking with a New Scoring Function, Efficient Optimization, and Multithreading. J. Comput. Chem. 2009, 31, 455− 461. (70) Neudert, G.; Klebe, G. Dsx: A Knowledge-Based Scoring Function for the Assessment of Protein−Ligand Complexes. J. Chem. Inf. Model. 2011, 51, 2731−2745. (71) McGann, M. R.; Almond, H. R.; Nicholls, A.; Grant, J. A.; Brown, F. K. Gaussian Docking Functions. Biopolymers 2003, 68, 76− 90. H

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Information and Modeling (72) Gehlhaar, D. K.; Verkhivker, G. M.; Rejto, P. A.; Sherman, C. J.; Fogel, D. R.; Fogel, L. J.; Freer, S. T. Molecular Recognition of the Inhibitor Ag-1343 by Hiv-1 Protease: Conformationally Flexible Docking by Evolutionary Programming. Chem. Biol. 1995, 2, 317−324. (73) Baxter, C. A.; Murray, C. W.; Clark, D. E.; Westhead, D. R.; Eldridge, M. D. Flexible Docking Using Tabu Search and an Empirical Estimate of Binding Affinity. Proteins: Struct., Funct., Genet. 1998, 33, 367−382. (74) Eldridge, M.; Murray, C.; Auton, T.; Paolini, G.; Mee, R. Empirical Scoring Functions: I. The Development of a Fast Empirical Scoring Function to Estimate the Binding Affinity of Ligands in Receptor Complexes. J. Comput.-Aided Mol. Des. 1997, 11, 425−445.

I

DOI: 10.1021/acs.jcim.5b00279 J. Chem. Inf. Model. XXXX, XXX, XXX−XXX