- Email: [email protected]
Optimization of virtual screening protocols: FlexX based virtual screening for COX-2 inhibitors reveals the importance of tailoring screen parameters
Journal of Molecular Structure (Theochem) 676 (2004) 1–5 www.elsevier.com/locate/theochem
Optimization of virtual screening protocols: FlexX based virtual screening for COX-2 inhibitors reveals the importance of tailoring screen parameters Andra´s Bauera, Zolta´n Kova´ria,b, Gyo¨rgy M. Keseru˝a,b,* b
a Department of Computer Assisted Drug Discovery, Gedeon Richter Ltd, P.O. Box 27, H-1475 Budapest 10, Hungary Department of Chemical Information Technologies, Budapest University of Technology and Economics, H-1111, Budapest, Szt. Gelle´rt te´r 4, Hungary
Received 18 September 2003; revised 18 September 2003; accepted 16 January 2004
Abstract As molecular docking poses a challenging problem in rational drug design many workgroups have lavished their attention on finding and testing various docking techniques in order to maximize their benefits. One such experiment we report was conducted using the FlexX docking algorithm on a variety of data sets containing COX-2 inhibitors as well as randomly selected inactive molecules which were docked into the COX-2 crystal structure to test the predictive power of scoring functions available in either FlexX or CScore. Our goal was to determine both their abilities to separate active and inactive molecules, as well as to find which function shows the best correlation with the molecules’ IC50 values. We have shown, that F-score had the best ability to rank the COX-2 inhibitors according to their IC50 values, while ChemScore has proven itself to be the most adequate for finding possible hits in large databases. q 2004 Elsevier B.V. All rights reserved. Keywords: Virtual screening; Scoring functions; COX-2; Docking
1. Introduction The drawback of traditional NSAIDs have for long been apparent, creating an ever increasing problem in the treatment of certain pain types. With the discovery of the inducible cyclooxygenase enzyme [1,2], the possibility of creating selective, non-ulcerogenic anti-inflammatory drugs appeared. These drugs, i.e. selective COX-2 inhibitors are known to produce much milder side effects, and are considered safer and more effective than their non-selective counterparts. Throughout the years COX-2 inhibitors have been subjects to many virtual screening experiments [3,4]. Thus, it is not surprising that it acquired the status of a sort Abbreviations: COX-2: cyclooxygenase-2; RMS: root mean square; EF: enrichment factor, calculated by dividing the ratio of active molecules in the entire dataset by the ratio of active molecules in the top five or 10% of the aforementioned dataset; IC50: the concentration of inhibitors at which enzyme activity drops to half of its original activity; pIC50: the negative logarithm of IC50. * Corresponding author. Department of Computer Assisted Drug Discovery, Gedeon Richter Ltd, P.O. Box 27, H-1475 Budapest 10, Hungary. Tel.: þ36-1-4314605; fax: þ36-1-4326002. E-mail address: [email protected] (G.M. Keseru˝). 0166-1280/$ - see front matter q 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.theochem.2004.01.016
of golden standard, a benchmark for testing and validating new high throughput docking and virtual screening methods. Hence, our choice for conducting a virtual screening experiment seemed obvious. Not only have the development of new COX-2 selective drugs had a great medical importance, its abundance in related literature provided an excellent opportunity to test the effectiveness of FlexX [5] in virtual screening. Just as in in vitro HTS, virtual screening protocols—the virtual analogs of assays—also need refining and optimization in order to provide them with the versatility and robustness needed to perform acceptable. In the past, docking studies have generally concentrated on two methods to test the scoring functions for two distinct abilities. One method usually involved large databases which contained both active, and inactive molecules, and the software’s main job was, so to say, to separate the chaff from the wheat, to find the active molecules by ranking them as close to the top as possible [6]. The other type of testing generally had to do with the function’s ability to differentiate between the potency of different active molecules. In this latter case, the parameter used most frequently had been the correlation between the molecules rankings and their pIC50 values [7].
2
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5
What remained unclear is how to merge exactly these two types of data, or at least how to find common ground on which to interpret it. Another question raised was the efficacy of FlexX in both of these roles. However, as FlexX provides a set of solutions, these ‘poses’ could further be evaluated according to a plethora of scoring functions such as GoldScore [8], PMF-Score [9], ChemScore [10] or DScore [11] as well as the software’s own scoring function (F-Score).
2. Methods The protein structure used in this study was a crystal structure retrieved from the Protein Data Bank Bernstein et al. [12] (PDB ID:6COX) [13]. For the docking purposes protein atoms of chain A, the relating prosthetic group haem were used and the bound inhibitor SC-558 was used to define the binding pocket. The docking experiments were partly carried out with this structure, here referred to as ‘crystal structure’, and partly using a minimized version of this protein structure processed as follows, referred to as ‘relaxed structure’. The charges for both structures were computed with the MMFF94 Halgren [14] method. The structure was then minimized using the Powell algorithm for 20 steps with no initial optimization whatsoever, using the ˚ MMFF94 force field. The non-bonding cutoff was set to 9 A and a distance dependent dielectric constant was applied ðe ¼ 4rÞ: All atoms of the protein were treated as aggregates, with the exception of those being within the ˚ radius of the bound ligand. Fitting Ca atoms of the 8A crystal structure and relaxed structure resulted in an RMS ˚. distance of 0.0455 A A total of 54 known COX-2 inhibitors, from various sources [15 –17] were split into two sets. One of the sets was used to find the best docking method (test set), whilst the other was to serve as a means to validating the results (validate set). All molecules were subjected to SYBYL’s CONCORD clean-up option. Each molecule was preprocessed before the docking calculations during which it was minimized with 600 iterations of Powell algorithm using Tripos force field, and was given charges according to the Gasteiger-Hu¨ckel method. Both the cleaned-up and the minimized molecules were used in our experiments. To conduct the docking experiment we used the FlexX software, version 1.11.1. while the various scores were generated using the CScore software. The test set of 28 molecules was docked into both the crystal structure and the relaxed one. The FlexX Receptor Description Files were created from the structures in a similar manner. For both ˚ structures the active sites included all residues within an 8 A radius of the bound inhibitor. CScore calculations were enabled and set to serial mode. During CScore calculations both the relax structure and the score relaxed functions were enabled. The first option allowed the ligand more movement, while the second option performed CScore
calculations on the relaxed ligands. All other functions were set to default values. All results were retrieved by selecting the best docking solution for each ligand according to G-Score PMF-Score, D-Score, ChemScore and FlexX’s own score (F-Score). These poses were then postprocessed, as follows. Each ligand’s best pose was merged into the original structure, then the whole structure was minimized for 10, 30 or 50 iterations using MMFF94 force field and charges, by Powell algorithm. During this minimization, the whole protein was defined as an aggregate, so as only the atoms of the docked ligand were allowed to move. The minimized ligand was extracted, and rescored using the CScore functions. The results were then processed in the STATISTICA 6.0 (StatSoft, Inc. (2001)). Here, the ligands’ rankings according to each scoring function were compared to their COX-2 inhibiting potential, represented by their pIC50 values. The scoring functions’ predictive potentials were evaluated by calculating their Spearman rank-order correlation values using the 6.0 version of STATISTICA. The validating set had been both pre- and postprocessed. For the enrichment calculations two test databases were used, one of 25 (database 1), and another of 558 molecules (database 2). The sets contained five and 11 COX-2 selective inhibitor molecules taken from the test set, respectively. The rest of the molecules were taken from the WDI database (World Drug Index, Derwent Inc., (2001)), with no known COX-2 inhibiting activity. The latter database was subjected to docking both with and without preprocessing, whilst the database 1 was docked only after ligand minimization. For these experiments, we used the crystal structure of the COX-2 mentioned above. Docking settings were those used in the predictive ability testing experiments, with the relax ligand option disabled. The evaluation of these runs was conducted by calculating the sets’ Enrichment Factors (EF) for the top ranking 5 and 10% of databases. To validate the results achieved with the first two runs, a third database of 557 molecules was created (database 3), one that contained 10 COX-2 inhibiting molecules taken from the previous experiment’s validating set. This database was processed exactly the same manner as database 2.
3. Results and discussion As can be seen on Table 1, the rank-correlation values obtained from the FlexX runs on the test set clearly show that without any pre- or postprocessing, PMF score produced the best results in terms of predicting the ligands inhibitory potential (Table 1). It is also clear that even this result falls significantly short of correctly predicting the molecules ranks, only one, of the displayed correlation values, PMF score’s is significant. It is also noticeable that only 26 of the 28 ligands have been docked successfully by
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5 Table 1 R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 Scoring function
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.330415 0.390769 0.098621 0.206959 0.148769
0.099227 0.048401 0.631719 0.310379 0.468256
No. indicates the number of successfully docked ligands.
3
Table 4 Rank correlation values of different scoring functions with pre-docking minimization. R-values are calculated with Spearman rank-order correlation, with correlations being significant a p , 0:05 Scoring functions
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
28 28 28 28 28
0.427477 0.318555 0.144499 0.184455 0.195950
0.023268 0.098505 0.463180 0.347401 0.317637
No. indicates the number of successfully docked ligands. Results were obtained from calculations performed on the test set.
FlexX. Results with the minimized protein structure (Table 2) or with the use of FlexX’s relax ligand function (Table 3) are even more dismal. However, it has granted us great insight, for it has revealed that the relaxation of the protein structure did not improve the predictive value of any scoring function. It has also shown, that the FlexX option to relax the ligands before docking has not improved the results obtained (Table 3). Thus, the calculations showed that the docking of unrelaxed ligands into the crystal structure yielded the best results. Still better are the results shown on Table 4. Even though the relaxation of ligands have not helped achieve better results, the preprocessing discussed above not only improved the rank-order correlation but has also enabled FlexX to dock in those ligands it has previously failed to dock. It has also led us to the conclusion that the use of CONCORD clean-up procedure is not sufficient in dealing Table 2 The minimization of the protein structure led to worse correlation values than the use of the original structure. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 Scoring function
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.227115 20.219917 20.192562 20.043551 20.155516
1.14249 21.10441 20.96135 20.21356 20.77125
No. indicates the number of successfully docked ligands.
with the ligands properties. It is also clear, that by far, F-Score became the best performer with its correlation value soaring to 0.427. The postprocessing was, in light of evidence presented above, performed on those minimized ligands docked in the crystal structure. The optimal level of minimization was worked out by submitting the docking solutions to a series of minimizations, after which Spearman rank correlation was performed. The docking solutions receiving the best scores by FlexX’s score were minimized as described above for 10, 30 and 50 iterations, as it can be seen on Table 5. A minimization of 10 steps yielded the best results. Further minimization decreased the scoring function’s predictive capability. Spearman rank correlation after a 10 iteration minimization for the best solutions by F-Score demonstrate a staggering value of 0.564 (Table 6). Table 5 Impact of minimization on F-Score’s predictive value. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 No. of iterations
R
p-level
0 10 30 50
0.423156 0.564313 0.490421 0.227578
0.024856 0.001760 0.008061 0.244146
Data obtained from experiments with the test set.
Table 3 The enabling of FlexX’s relax ligand option greatly reduced correlation values. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Table 6 Rank correlation values after a 10 iteration postdocking minimization, performed on the test set. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Scoring function
No.
R
p-level
Scoring functions
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.268376 20.078974 20.007863 0.029065 0.184964
0.184964 0.701359 0.969589 0.887915 0.943177
F-Score PMF-Score D-Score G-Score ChemScore
28 28 28 28 28
0.564313 20.008758 0.163656 20.000547 0.220033
0.001760 0.964722 0.405331 0.997794 0.260556
No. indicates the number of successfully docked ligands.
No. indicates the number of successfully docked ligands.
4
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5
Table 7 Rank correlation values after a 10 iteration postdocking minimization. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Test set Validate set
No.
R
p-level
28 26
0.564313 0.501369
0.001760 0.009072
No. indicates the number of successfully docked ligands.
For validation of the result achieved with the test set, we used another set of 26 molecules (validate set) which was pre- and postprocessed in the same manner as the test set. Results of these calculations can be seen on Table 7 and show a Spearman rank correlation value of 0.501. The enrichment studies were conducted after the correlation experiments, to observe their efficiency in the enrichment studies. The database (database 2) containing 558 molecules was docked with and without preprocessing. Results of these experiments can be seen and compared on Table 8. The results of the docking experiments showed that contrary to the correlation studies, the scoring function yielding the best results was undoubtedly ChemScore. It is also clear that initial ligand minimization is also essential for obtaining better results. The first dataset consisting of 25 molecules (database 1) was submitted to both pre- and postprocessing. It was shown that postprocessing by minimization of the ligands in the protein structure did not significantly improve scores. To validate these results, another database of 557 molecules was created (database 3), and after preprocessing
F-Score D-Score PMF-Score G-Score ChemScore
4. Conclusions Eventually, FlexX has proven to be an effective tool for database screening, whilst its ability to correctly rank inhibitors according to their potency required some preand post processing in order to produce significant results. This processing involved minimization of ligands, and a post processing, in which the best scoring pose was minimized together with the protein’s active site, and rescored again. For the high throughput docking, ChemScore achieved results that show its excellent properties in database screening. It was also shown, however, that minimization of the ligand structures can and will improve results. The study also revealed a suite, with which, albeit using slightly different methods both types of studies can be conducted, thus improving efficiency, and providing a way of not only finding hits, but also assesing their potential.
Acknowledgements The financial support provided by the Grant OTKA (T 042933) of the Hungarian National Research Fund and the support of Tripos Inc. are gratefully acknowledged.
Table 8 Results of FlexX runs on database 2 Scoring functions
it was also docked in the crystal structure. Enrichment factors of this and the 558 molecule database (database 2) calculated based on results generated by ChemScore are shown for comparison in Table 9. The results of the validating run were much the same as those of the previous runs, with ChemScore producing the best enrichment in both the top 5 and the top 10%. In the validating set, 70% of the COX-2 active molecules made it in the top 10%.
Without preprocessing
With preprocessing
EF5%
EF10%
EF5%
EF10%
1.8117 1.8117 0 1.8117 3.6234
0.9058 1.8117 0.9058 0.9058 1.8117
0 9.0584 0 5.4351 9.0584
0.9058 4.5292 0.9058 2.7175 5.4351
Scoring functions indicate the method of extracting the best scoring solutions. EF stands for enrichment factor, calculated for the top 5 or 10 percent of the of the docked molecules.
Table 9 Best docking solutions were extracted according to ChemScore No. of dataset
EF5%
EF10%
Database 2 Database 3 (validating set)
9.0584 7.9714
5.4351 6.975
References [1] T. Hla, K. Neilson, Human cyclooxygenase-2 cDNA, Proc. Natl. Acad. Sci. USA 89 (16) (1992) 7384–7388. [2] D.A. Kujubu, B.S. Fletcher, B.C. Varnum, R.W. Lim, H.R. Herschman, TIS10, a phorbol ester tumor promoter-inducible mRNA from Swiss 3T3 cells, encodes a novel prostaglandin synthase/cyclooxygenase homologue, J. Biol. Chem. 266 (20) (1991) 12866–12872. [3] H.B. Broughton, A method for including protein flexibility in protein–ligand docking: improving tools for database mining and virtual screening, J. Mol. Graph. Model. 18 (3) (2000) 247– 257. see also pages 302– 304. [4] G.R. Desiraju, B. Gopalakrishnan, R.K. Jetti, A. Nagaraju, D. Raveendra, J.A. Sarma, M.E. Sobhia, R. Thilagavathi, Computeraided design of selective COX-2 inhibitors: comparative molecular field analysis, comparative molecular similarity indices analysis, and docking studies of some 1,2-diarylimidazole derivatives, J. Med. Chem. 45 (22) (2002) 4847–4857. [5] M. Rarey, B. Kramer, T. Lengauer, G. Klebe, A fast flexible docking method using an incremental construction algorithm, J. Mol. Biol. 261 (3) (1996) 470–489.
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5 [6] R. Wang, S. Wang, How does consensus scoring work for virtual library screening? An idealized computer experiment, J. Chem. Inf. Comput. Sci. 41 (5) (2001) 1422–1426. [7] S.S. Wesolowski, W.L. Jorgensen, Estimation of binding affinities for celecoxib analogues with COX-2 via Monte Carlo-extended linear response, Bioorg. Med. Chem. Lett. 12 (3) (2002) 267–270. [8] G. Jones, P. Willett, R.C. Glen, A.R. Leach, R. Taylor, Development and validation of a genetic algorithm for flexible docking, J. Mol. Biol. 267 (3) (1997) 727 –748. [9] I. Muegge, Y.C. Martin, A general and fast scoring function for protein– ligand interactions: a simplified potential approach, J. Med. Chem. 42 (5) (1999) 791–804. [10] M.D. Eldridge, C.W. Murray, T.R. Auton, G.V. Paolini, R.P. Mee, 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. 11 (5) (1997) 425–445. [11] I.D. Kuntz, J.M. Blaney, S.J. Oatley, R. Langridge, T.E. Ferrin, A geometric approach to macromolecule– ligand interactions, J. Mol. Biol. 161 (2) (1982) 269 –288. [12] F.C. Bernstein, T.F. Koetzle, G.J. Williams, E.F. Meyer Jr., M.D. Brice, J.R. Rodgers, O. Kennard, T. Shimanouchi, M. Tasumi, The Protein Data Bank: a computer-based archival file for macromolecular structures, J. Mol. Biol. 112 (3) (1977) 535–542. [13] R.G. Kurumbail, A.M. Stevens, J.K. Gierse, J.J. McDonald, R.A. Stegeman, J.Y. Pak, D. Gildehaus, J.M. Miyashiro, T.D. Penning, K. Seibert, P.C. Isakson, W.C. Stallings, Structural basis for selective inhibition of cyclooxygenase-2 by anti-inflammatory agents, Nature
[14]
[15]
[16]
[17]
5
384 (6610) (1996) 644 –648. Erratum in: Nature 1997 Feb 6; 385(6616):555. A.T. Halgren, Maximally diagonal force constants in dependent angle-bending coordinates. II. Implications for the design of empirical force fields, J. Am. Chem. Soc. 112 (12) (1990) 4710–4723. J.M. Janusz, P.A. Young, J.M. Ridgeway, M.W. Scherz, K. Enzweiler, L.I. Wu, L. Gan, R. Darolia, R.S. Matthews, D. Hennes, D.E. Kellstein, S.A. Green, J.L. Tulich, T. Rosario-Jansen, I.J. Magrisso, K.R. Wehmeyer, D.L. Kuhlenbeck, T.H. Eichhold, R.L. Dobson, S.P. Sirko, R.W. Farmer, New cyclooxygenase-2/5-lipoxygenase inhibitors. 1. 7-tert-buty1-2,3-dihydro-3,3-dimethylbenzofuran derivatives as gastrointestinal safe antiinflammatory and analgesic agents: discovery and variation of the 5-keto substituent, J. Med. Chem. 41 (7) (1998) 1112–1123. H.C. Huang, J.J. Li, D.J. Garland, T.S. Chamberlain, E.J. Reinhard, R.E. Manning, K. Seibert, C.M. Koboldt, S.A. Gregory, G.D. Anderson, A.W. Veenhuizen, Y. Zhang, W.E. Perkins, E.G. Burton, J.N. Cogburn, P.C. Isakson, D.B. Reitz, Diarylspiro[2.4]heptenes as orally active, highly selective cyclooxygenase-2 inhibitors: synthesis and structure-activity relationships, J. Med. Chem. 39 (1) (1996) 253–266. H. Liu, X. Huang, J. Shen, X. Luo, M. Li, B. Xiong, G. Chen, J. Shen, Y. Yang, H. Jiang, K. Chen, Inhibitory mode of 1,5-diarylpyrazole derivatives against cyclooxygenase-2 and cyclooxygenase-1: molecular docking and 3D QSAR analyses, J. Med. Chem. 45 (22) (2002) 4816– 4827.
Optimization of virtual screening protocols: FlexX based virtual screening for COX-2 inhibitors reveals the importance of tailoring screen parameters Andra´s Bauera, Zolta´n Kova´ria,b, Gyo¨rgy M. Keseru˝a,b,* b
a Department of Computer Assisted Drug Discovery, Gedeon Richter Ltd, P.O. Box 27, H-1475 Budapest 10, Hungary Department of Chemical Information Technologies, Budapest University of Technology and Economics, H-1111, Budapest, Szt. Gelle´rt te´r 4, Hungary
Received 18 September 2003; revised 18 September 2003; accepted 16 January 2004
Abstract As molecular docking poses a challenging problem in rational drug design many workgroups have lavished their attention on finding and testing various docking techniques in order to maximize their benefits. One such experiment we report was conducted using the FlexX docking algorithm on a variety of data sets containing COX-2 inhibitors as well as randomly selected inactive molecules which were docked into the COX-2 crystal structure to test the predictive power of scoring functions available in either FlexX or CScore. Our goal was to determine both their abilities to separate active and inactive molecules, as well as to find which function shows the best correlation with the molecules’ IC50 values. We have shown, that F-score had the best ability to rank the COX-2 inhibitors according to their IC50 values, while ChemScore has proven itself to be the most adequate for finding possible hits in large databases. q 2004 Elsevier B.V. All rights reserved. Keywords: Virtual screening; Scoring functions; COX-2; Docking
1. Introduction The drawback of traditional NSAIDs have for long been apparent, creating an ever increasing problem in the treatment of certain pain types. With the discovery of the inducible cyclooxygenase enzyme [1,2], the possibility of creating selective, non-ulcerogenic anti-inflammatory drugs appeared. These drugs, i.e. selective COX-2 inhibitors are known to produce much milder side effects, and are considered safer and more effective than their non-selective counterparts. Throughout the years COX-2 inhibitors have been subjects to many virtual screening experiments [3,4]. Thus, it is not surprising that it acquired the status of a sort Abbreviations: COX-2: cyclooxygenase-2; RMS: root mean square; EF: enrichment factor, calculated by dividing the ratio of active molecules in the entire dataset by the ratio of active molecules in the top five or 10% of the aforementioned dataset; IC50: the concentration of inhibitors at which enzyme activity drops to half of its original activity; pIC50: the negative logarithm of IC50. * Corresponding author. Department of Computer Assisted Drug Discovery, Gedeon Richter Ltd, P.O. Box 27, H-1475 Budapest 10, Hungary. Tel.: þ36-1-4314605; fax: þ36-1-4326002. E-mail address: [email protected] (G.M. Keseru˝). 0166-1280/$ - see front matter q 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.theochem.2004.01.016
of golden standard, a benchmark for testing and validating new high throughput docking and virtual screening methods. Hence, our choice for conducting a virtual screening experiment seemed obvious. Not only have the development of new COX-2 selective drugs had a great medical importance, its abundance in related literature provided an excellent opportunity to test the effectiveness of FlexX [5] in virtual screening. Just as in in vitro HTS, virtual screening protocols—the virtual analogs of assays—also need refining and optimization in order to provide them with the versatility and robustness needed to perform acceptable. In the past, docking studies have generally concentrated on two methods to test the scoring functions for two distinct abilities. One method usually involved large databases which contained both active, and inactive molecules, and the software’s main job was, so to say, to separate the chaff from the wheat, to find the active molecules by ranking them as close to the top as possible [6]. The other type of testing generally had to do with the function’s ability to differentiate between the potency of different active molecules. In this latter case, the parameter used most frequently had been the correlation between the molecules rankings and their pIC50 values [7].
2
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5
What remained unclear is how to merge exactly these two types of data, or at least how to find common ground on which to interpret it. Another question raised was the efficacy of FlexX in both of these roles. However, as FlexX provides a set of solutions, these ‘poses’ could further be evaluated according to a plethora of scoring functions such as GoldScore [8], PMF-Score [9], ChemScore [10] or DScore [11] as well as the software’s own scoring function (F-Score).
2. Methods The protein structure used in this study was a crystal structure retrieved from the Protein Data Bank Bernstein et al. [12] (PDB ID:6COX) [13]. For the docking purposes protein atoms of chain A, the relating prosthetic group haem were used and the bound inhibitor SC-558 was used to define the binding pocket. The docking experiments were partly carried out with this structure, here referred to as ‘crystal structure’, and partly using a minimized version of this protein structure processed as follows, referred to as ‘relaxed structure’. The charges for both structures were computed with the MMFF94 Halgren [14] method. The structure was then minimized using the Powell algorithm for 20 steps with no initial optimization whatsoever, using the ˚ MMFF94 force field. The non-bonding cutoff was set to 9 A and a distance dependent dielectric constant was applied ðe ¼ 4rÞ: All atoms of the protein were treated as aggregates, with the exception of those being within the ˚ radius of the bound ligand. Fitting Ca atoms of the 8A crystal structure and relaxed structure resulted in an RMS ˚. distance of 0.0455 A A total of 54 known COX-2 inhibitors, from various sources [15 –17] were split into two sets. One of the sets was used to find the best docking method (test set), whilst the other was to serve as a means to validating the results (validate set). All molecules were subjected to SYBYL’s CONCORD clean-up option. Each molecule was preprocessed before the docking calculations during which it was minimized with 600 iterations of Powell algorithm using Tripos force field, and was given charges according to the Gasteiger-Hu¨ckel method. Both the cleaned-up and the minimized molecules were used in our experiments. To conduct the docking experiment we used the FlexX software, version 1.11.1. while the various scores were generated using the CScore software. The test set of 28 molecules was docked into both the crystal structure and the relaxed one. The FlexX Receptor Description Files were created from the structures in a similar manner. For both ˚ structures the active sites included all residues within an 8 A radius of the bound inhibitor. CScore calculations were enabled and set to serial mode. During CScore calculations both the relax structure and the score relaxed functions were enabled. The first option allowed the ligand more movement, while the second option performed CScore
calculations on the relaxed ligands. All other functions were set to default values. All results were retrieved by selecting the best docking solution for each ligand according to G-Score PMF-Score, D-Score, ChemScore and FlexX’s own score (F-Score). These poses were then postprocessed, as follows. Each ligand’s best pose was merged into the original structure, then the whole structure was minimized for 10, 30 or 50 iterations using MMFF94 force field and charges, by Powell algorithm. During this minimization, the whole protein was defined as an aggregate, so as only the atoms of the docked ligand were allowed to move. The minimized ligand was extracted, and rescored using the CScore functions. The results were then processed in the STATISTICA 6.0 (StatSoft, Inc. (2001)). Here, the ligands’ rankings according to each scoring function were compared to their COX-2 inhibiting potential, represented by their pIC50 values. The scoring functions’ predictive potentials were evaluated by calculating their Spearman rank-order correlation values using the 6.0 version of STATISTICA. The validating set had been both pre- and postprocessed. For the enrichment calculations two test databases were used, one of 25 (database 1), and another of 558 molecules (database 2). The sets contained five and 11 COX-2 selective inhibitor molecules taken from the test set, respectively. The rest of the molecules were taken from the WDI database (World Drug Index, Derwent Inc., (2001)), with no known COX-2 inhibiting activity. The latter database was subjected to docking both with and without preprocessing, whilst the database 1 was docked only after ligand minimization. For these experiments, we used the crystal structure of the COX-2 mentioned above. Docking settings were those used in the predictive ability testing experiments, with the relax ligand option disabled. The evaluation of these runs was conducted by calculating the sets’ Enrichment Factors (EF) for the top ranking 5 and 10% of databases. To validate the results achieved with the first two runs, a third database of 557 molecules was created (database 3), one that contained 10 COX-2 inhibiting molecules taken from the previous experiment’s validating set. This database was processed exactly the same manner as database 2.
3. Results and discussion As can be seen on Table 1, the rank-correlation values obtained from the FlexX runs on the test set clearly show that without any pre- or postprocessing, PMF score produced the best results in terms of predicting the ligands inhibitory potential (Table 1). It is also clear that even this result falls significantly short of correctly predicting the molecules ranks, only one, of the displayed correlation values, PMF score’s is significant. It is also noticeable that only 26 of the 28 ligands have been docked successfully by
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5 Table 1 R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 Scoring function
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.330415 0.390769 0.098621 0.206959 0.148769
0.099227 0.048401 0.631719 0.310379 0.468256
No. indicates the number of successfully docked ligands.
3
Table 4 Rank correlation values of different scoring functions with pre-docking minimization. R-values are calculated with Spearman rank-order correlation, with correlations being significant a p , 0:05 Scoring functions
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
28 28 28 28 28
0.427477 0.318555 0.144499 0.184455 0.195950
0.023268 0.098505 0.463180 0.347401 0.317637
No. indicates the number of successfully docked ligands. Results were obtained from calculations performed on the test set.
FlexX. Results with the minimized protein structure (Table 2) or with the use of FlexX’s relax ligand function (Table 3) are even more dismal. However, it has granted us great insight, for it has revealed that the relaxation of the protein structure did not improve the predictive value of any scoring function. It has also shown, that the FlexX option to relax the ligands before docking has not improved the results obtained (Table 3). Thus, the calculations showed that the docking of unrelaxed ligands into the crystal structure yielded the best results. Still better are the results shown on Table 4. Even though the relaxation of ligands have not helped achieve better results, the preprocessing discussed above not only improved the rank-order correlation but has also enabled FlexX to dock in those ligands it has previously failed to dock. It has also led us to the conclusion that the use of CONCORD clean-up procedure is not sufficient in dealing Table 2 The minimization of the protein structure led to worse correlation values than the use of the original structure. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 Scoring function
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.227115 20.219917 20.192562 20.043551 20.155516
1.14249 21.10441 20.96135 20.21356 20.77125
No. indicates the number of successfully docked ligands.
with the ligands properties. It is also clear, that by far, F-Score became the best performer with its correlation value soaring to 0.427. The postprocessing was, in light of evidence presented above, performed on those minimized ligands docked in the crystal structure. The optimal level of minimization was worked out by submitting the docking solutions to a series of minimizations, after which Spearman rank correlation was performed. The docking solutions receiving the best scores by FlexX’s score were minimized as described above for 10, 30 and 50 iterations, as it can be seen on Table 5. A minimization of 10 steps yielded the best results. Further minimization decreased the scoring function’s predictive capability. Spearman rank correlation after a 10 iteration minimization for the best solutions by F-Score demonstrate a staggering value of 0.564 (Table 6). Table 5 Impact of minimization on F-Score’s predictive value. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05 No. of iterations
R
p-level
0 10 30 50
0.423156 0.564313 0.490421 0.227578
0.024856 0.001760 0.008061 0.244146
Data obtained from experiments with the test set.
Table 3 The enabling of FlexX’s relax ligand option greatly reduced correlation values. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Table 6 Rank correlation values after a 10 iteration postdocking minimization, performed on the test set. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Scoring function
No.
R
p-level
Scoring functions
No.
R
p-level
F-Score PMF-Score D-Score G-Score ChemScore
26 26 26 26 26
0.268376 20.078974 20.007863 0.029065 0.184964
0.184964 0.701359 0.969589 0.887915 0.943177
F-Score PMF-Score D-Score G-Score ChemScore
28 28 28 28 28
0.564313 20.008758 0.163656 20.000547 0.220033
0.001760 0.964722 0.405331 0.997794 0.260556
No. indicates the number of successfully docked ligands.
No. indicates the number of successfully docked ligands.
4
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5
Table 7 Rank correlation values after a 10 iteration postdocking minimization. R-values are calculated with Spearman rank-order correlation, with correlations being significant at p , 0:05
Test set Validate set
No.
R
p-level
28 26
0.564313 0.501369
0.001760 0.009072
No. indicates the number of successfully docked ligands.
For validation of the result achieved with the test set, we used another set of 26 molecules (validate set) which was pre- and postprocessed in the same manner as the test set. Results of these calculations can be seen on Table 7 and show a Spearman rank correlation value of 0.501. The enrichment studies were conducted after the correlation experiments, to observe their efficiency in the enrichment studies. The database (database 2) containing 558 molecules was docked with and without preprocessing. Results of these experiments can be seen and compared on Table 8. The results of the docking experiments showed that contrary to the correlation studies, the scoring function yielding the best results was undoubtedly ChemScore. It is also clear that initial ligand minimization is also essential for obtaining better results. The first dataset consisting of 25 molecules (database 1) was submitted to both pre- and postprocessing. It was shown that postprocessing by minimization of the ligands in the protein structure did not significantly improve scores. To validate these results, another database of 557 molecules was created (database 3), and after preprocessing
F-Score D-Score PMF-Score G-Score ChemScore
4. Conclusions Eventually, FlexX has proven to be an effective tool for database screening, whilst its ability to correctly rank inhibitors according to their potency required some preand post processing in order to produce significant results. This processing involved minimization of ligands, and a post processing, in which the best scoring pose was minimized together with the protein’s active site, and rescored again. For the high throughput docking, ChemScore achieved results that show its excellent properties in database screening. It was also shown, however, that minimization of the ligand structures can and will improve results. The study also revealed a suite, with which, albeit using slightly different methods both types of studies can be conducted, thus improving efficiency, and providing a way of not only finding hits, but also assesing their potential.
Acknowledgements The financial support provided by the Grant OTKA (T 042933) of the Hungarian National Research Fund and the support of Tripos Inc. are gratefully acknowledged.
Table 8 Results of FlexX runs on database 2 Scoring functions
it was also docked in the crystal structure. Enrichment factors of this and the 558 molecule database (database 2) calculated based on results generated by ChemScore are shown for comparison in Table 9. The results of the validating run were much the same as those of the previous runs, with ChemScore producing the best enrichment in both the top 5 and the top 10%. In the validating set, 70% of the COX-2 active molecules made it in the top 10%.
Without preprocessing
With preprocessing
EF5%
EF10%
EF5%
EF10%
1.8117 1.8117 0 1.8117 3.6234
0.9058 1.8117 0.9058 0.9058 1.8117
0 9.0584 0 5.4351 9.0584
0.9058 4.5292 0.9058 2.7175 5.4351
Scoring functions indicate the method of extracting the best scoring solutions. EF stands for enrichment factor, calculated for the top 5 or 10 percent of the of the docked molecules.
Table 9 Best docking solutions were extracted according to ChemScore No. of dataset
EF5%
EF10%
Database 2 Database 3 (validating set)
9.0584 7.9714
5.4351 6.975
References [1] T. Hla, K. Neilson, Human cyclooxygenase-2 cDNA, Proc. Natl. Acad. Sci. USA 89 (16) (1992) 7384–7388. [2] D.A. Kujubu, B.S. Fletcher, B.C. Varnum, R.W. Lim, H.R. Herschman, TIS10, a phorbol ester tumor promoter-inducible mRNA from Swiss 3T3 cells, encodes a novel prostaglandin synthase/cyclooxygenase homologue, J. Biol. Chem. 266 (20) (1991) 12866–12872. [3] H.B. Broughton, A method for including protein flexibility in protein–ligand docking: improving tools for database mining and virtual screening, J. Mol. Graph. Model. 18 (3) (2000) 247– 257. see also pages 302– 304. [4] G.R. Desiraju, B. Gopalakrishnan, R.K. Jetti, A. Nagaraju, D. Raveendra, J.A. Sarma, M.E. Sobhia, R. Thilagavathi, Computeraided design of selective COX-2 inhibitors: comparative molecular field analysis, comparative molecular similarity indices analysis, and docking studies of some 1,2-diarylimidazole derivatives, J. Med. Chem. 45 (22) (2002) 4847–4857. [5] M. Rarey, B. Kramer, T. Lengauer, G. Klebe, A fast flexible docking method using an incremental construction algorithm, J. Mol. Biol. 261 (3) (1996) 470–489.
A. Bauer et al. / Journal of Molecular Structure (Theochem) 676 (2004) 1–5 [6] R. Wang, S. Wang, How does consensus scoring work for virtual library screening? An idealized computer experiment, J. Chem. Inf. Comput. Sci. 41 (5) (2001) 1422–1426. [7] S.S. Wesolowski, W.L. Jorgensen, Estimation of binding affinities for celecoxib analogues with COX-2 via Monte Carlo-extended linear response, Bioorg. Med. Chem. Lett. 12 (3) (2002) 267–270. [8] G. Jones, P. Willett, R.C. Glen, A.R. Leach, R. Taylor, Development and validation of a genetic algorithm for flexible docking, J. Mol. Biol. 267 (3) (1997) 727 –748. [9] I. Muegge, Y.C. Martin, A general and fast scoring function for protein– ligand interactions: a simplified potential approach, J. Med. Chem. 42 (5) (1999) 791–804. [10] M.D. Eldridge, C.W. Murray, T.R. Auton, G.V. Paolini, R.P. Mee, 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. 11 (5) (1997) 425–445. [11] I.D. Kuntz, J.M. Blaney, S.J. Oatley, R. Langridge, T.E. Ferrin, A geometric approach to macromolecule– ligand interactions, J. Mol. Biol. 161 (2) (1982) 269 –288. [12] F.C. Bernstein, T.F. Koetzle, G.J. Williams, E.F. Meyer Jr., M.D. Brice, J.R. Rodgers, O. Kennard, T. Shimanouchi, M. Tasumi, The Protein Data Bank: a computer-based archival file for macromolecular structures, J. Mol. Biol. 112 (3) (1977) 535–542. [13] R.G. Kurumbail, A.M. Stevens, J.K. Gierse, J.J. McDonald, R.A. Stegeman, J.Y. Pak, D. Gildehaus, J.M. Miyashiro, T.D. Penning, K. Seibert, P.C. Isakson, W.C. Stallings, Structural basis for selective inhibition of cyclooxygenase-2 by anti-inflammatory agents, Nature
[14]
[15]
[16]
[17]
5
384 (6610) (1996) 644 –648. Erratum in: Nature 1997 Feb 6; 385(6616):555. A.T. Halgren, Maximally diagonal force constants in dependent angle-bending coordinates. II. Implications for the design of empirical force fields, J. Am. Chem. Soc. 112 (12) (1990) 4710–4723. J.M. Janusz, P.A. Young, J.M. Ridgeway, M.W. Scherz, K. Enzweiler, L.I. Wu, L. Gan, R. Darolia, R.S. Matthews, D. Hennes, D.E. Kellstein, S.A. Green, J.L. Tulich, T. Rosario-Jansen, I.J. Magrisso, K.R. Wehmeyer, D.L. Kuhlenbeck, T.H. Eichhold, R.L. Dobson, S.P. Sirko, R.W. Farmer, New cyclooxygenase-2/5-lipoxygenase inhibitors. 1. 7-tert-buty1-2,3-dihydro-3,3-dimethylbenzofuran derivatives as gastrointestinal safe antiinflammatory and analgesic agents: discovery and variation of the 5-keto substituent, J. Med. Chem. 41 (7) (1998) 1112–1123. H.C. Huang, J.J. Li, D.J. Garland, T.S. Chamberlain, E.J. Reinhard, R.E. Manning, K. Seibert, C.M. Koboldt, S.A. Gregory, G.D. Anderson, A.W. Veenhuizen, Y. Zhang, W.E. Perkins, E.G. Burton, J.N. Cogburn, P.C. Isakson, D.B. Reitz, Diarylspiro[2.4]heptenes as orally active, highly selective cyclooxygenase-2 inhibitors: synthesis and structure-activity relationships, J. Med. Chem. 39 (1) (1996) 253–266. H. Liu, X. Huang, J. Shen, X. Luo, M. Li, B. Xiong, G. Chen, J. Shen, Y. Yang, H. Jiang, K. Chen, Inhibitory mode of 1,5-diarylpyrazole derivatives against cyclooxygenase-2 and cyclooxygenase-1: molecular docking and 3D QSAR analyses, J. Med. Chem. 45 (22) (2002) 4816– 4827.