ProSelection: A novel algorithm to select proper ... - ACS Publications

Subscriber access provided by Nanyang Technological Univ

Article

ProSelection: A novel algorithm to select proper protein structure subsets for in silico target identification and drug discovery research Nanyi Wang, Lirong Wang, and Xiang-Qun (Sean) Xie J. Chem. Inf. Model., Just Accepted Manuscript • DOI: 10.1021/acs.jcim.7b00277 • Publication Date (Web): 10 Oct 2017 Downloaded from http://pubs.acs.org on October 11, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Information and Modeling is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

ProSelection: A novel algorithm to select proper protein structure subsets for in silico target identification and drug discovery research Nanyi Wang †, Lirong Wang†, *, and Xiang-Qun Xie†, * †

Department of Pharmaceutical Sciences and Computational Chemical Genomics

Screening Center, School of Pharmacy; NIH National Center of Excellence for Computational Drug Abuse Research; Drug Discovery Institute; University of Pittsburgh, Pittsburgh, Pennsylvania 15260, United States *Corresponding Author: Xiang-Qun (Sean) Xie, MBA, Ph.D. Professor of Pharmaceutical Sciences/Drug Discovery Institute Director of CCGS and NIDA CDAR Centers 335 Sutherland Drive, 206 Pavilion, School of Pharmacy University of Pittsburgh Pittsburgh, PA15261, USA 412-383-5276 (Phone) 412-383-7436 (Fax) Email: [email protected] *Co-corresponding Author: Lirong Wang, Ph.D. Assistant Professor of Pharmaceutical Sciences, School of Pharmacy University of Pittsburgh 517 Salk Hall, 3501 Terrace Street Pittsburgh, PA15261, USA 412-624-8118 (Phone) Email: [email protected]

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ABSTRACT: Molecular docking is widely applied to computer-aided drug design and has become relatively mature in the recent decades. Application of docking in modeling varies from single lead compound optimization to large-scale virtual screening. The performance of molecular docking is highly dependent on the protein structures selected. It is especially challenging for large-scale target prediction research when multiple structures are available for a single target. Therefore, we have established ProSelection, a docking preferred-protein selection algorithm, in order to generate the proper structure subset(s). By ProSelection algorithm, protein structures of “weak selector” are filtered out whereas structures of “strong selector” are kept. Specifically, the structure which has a good statistical performance of distinguishing active ligands from inactive ligands is defined as a “strong selector”. In this study, 249 protein structures of 14 autophagy-related targets are investigated. Surflex-dock was used as the docking engine to distinguish active and inactive compounds against these protein structures. Both t-test and Mann-Whitney U test were used to distinguish the “strong selector” from “weak selector” based on the normality of the docking score distribution. The suggested docking score threshold for active ligands (SDA) was generated for each “strong selector” structure according to the receiver operating characteristic (ROC) curve. The performance of ProSelection was further validated by predicting the potential off-targets of 43 FDA-approved small molecule antineoplastic drugs. Overall, ProSelection will accelerate the computational work in protein structure selection and could be a useful tool for molecular docking, target prediction, and protein-chemical database establishment research.


Page 2 of 32

Page 3 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


INTRODUCTION Molecular docking, which is commonly applied to computer-aided drug design, has become relatively mature in the recent decades. Since the 1980s, pioneering works have been reported to achieve rigid docking of small molecules to protein structures1. Since then the number of algorithms available to evaluate ligand-protein interaction is large and ever increasing. Docking programs like Autodock2, Surflex-Dock3, GOLD4, Glide5, FlexX6, DOCK7, ICM8, and FRED9 are now widely used in receptor to small molecule docking. The number of users for relatively new algorithms such as Autodock Vina2, EADock DSS10, H-DOCK11, and ParaDockS12 is still growing. Application of docking in modeling varies from single lead compound optimization to large-scale virtual screening. Typical performance of prediction for widely used programs in docking has a success rate close to 20%13. This rate may further increase to over 50% when other conditions, like protein conformations from multiple complexes, are optimized. Computational methods for molecular docking rely heavily on the representation of proteins. However, multiple crystal structures are usually available for a single protein target. For instance, 19 protein targets are randomly selected from the PDB database and listed in the Supplementary table 1, while 18 of them have more than one crystal structure available. Multiple crystal structures can cause problems when docking an identical ligand that yields different results from these structures. One possible solution to this is incorporating multiple conformations of a single receptor into the docking algorithms14. Ensemble docking (ED) is one of the primary research focuses within the area of receptor flexibility, in which compounds are docked against multiple conformations of the target receptor15. This approach mimics the dynamic behavior of the protein by providing a structural degree of freedom, where each ligand may find a matching receptor16. Combined with ligands’ flexibility, ensemble docking based prediction will achieve a higher success rate, compared to the rigid binding17. Receptor ensembles for ED can be divided into two categories: (1) ensembles base on experimental structures56,57,58, and (2) ensembles base on in silico-generated structures. In the experimentally derived ensembles, the area under the ROC curve (AUC) is always employed as the enrichment metric. It reveals that the experimentally derived ensembles only yield better enrichment in some cases when



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

compared to individual members18. Abagyan’s group even reported an anticooperative behavior of the experimentally derived ensembles, in which the performances of ensembles actually decreased, due to the addition of some ensemble members19. Encouragingly, the performance of multiple conformations (MRCs) in virtual ligand screening can be further improved by using a proper subset of receptor structures, which is selected from among the originally available conformers19. Researchers then attempt to define a method to select proper subsets of protein structures for generating the receptor ensembles. Some recommended filtering out receptor subsets that provide only positive docking scores20. However, this filter protocol is still based on the empirical observation of the correspondence between the runs that yielded positive docking scores for all molecules in the dataset and those performing worse than a random selection. In this work, we introduce the Protein Structure Selection (ProSelection), a new computational and statistical algorithm for picking the proper protein structure subsets according to the user-defined docking protocol. ProSelection is tested using the currently widely used docking algorithm Surflex-dock and is designed for the virtual ligand screening against protein targets. Both t-test and Mann-Whitney U test are used as the enrichment metrics in ProSelection, which is different from the previous selection approaches of conformers using AUC20. This research cannot be considered as traditional virtual compound screening. Instead, our method is more like “structure screening”. That is, our approach is designed to find top structures with better capability for distinguishing compounds. It could be more useful when applying docking approach for target identification. That is, a small molecule will be used to search against thousands of PDB structures. In this situation, if strong selectors are used, we can reduce the prediction errors. ProSelection exhibits a significant power in protein structures selection, making it a promising option for the proper subsets selection for both the experimental structures based receptor ensembles and rigid docking. The research design of this work is presented in the flowchart (Figure 1).

METHODS AND DATA SETS Protein and Ligand Test-set Construction. ProSelection was tested in one of our in-house Autophagy Knowledgebase (http://www.cbligand.org/autophagy) which was


Page 4 of 32

Page 5 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


established using a procedure similar to AlzPlatform21, a chemogenomics knowledgebase for Alzheimer's disease. In the Autophagy knowledgebase, 102 autophagy-related protein targets were collected from literature and public databases, including

PubMed

(www.ncbi.nlm.nih.gov/pubmed),

(pubchem.ncbi.nlm.nih.gov/), (https://scifinder.cas.org/),

DrugBank

UniProt

PubChem

(http://www.drugbank.ca),

(http://www.uniprot.org/),

and

SciFinder ChEMBL

(https://www.ebi.ac.uk/chembl/). Protein structures corresponding to the targets in the Autophagy Knowledgebase were downloaded from RCSB Protein Data Bank22. The number of known active and inactive chemical ligands for each autophagy-target was calculated. Only the targets which have more than 40 known ligands (20 active ligands and 20 inactive ligands) were kept in this study as well as their associated protein structures. These targets (14 in total) and their protein structures (249 in total) were then used to generate the protein test-set and structure test-set (Figure 1), respectively. Active and inactive chemical ligands, which are associated with each of those autophagic targets, are acquired from the ChEMBL database to build the ligand testset. Specifically, compounds with IC50 values lower than 5 nM against a target are considered as the active ligands for that target. While compounds with IC50 values higher than 50 μM are considered as the inactive ligands. Compounds between 5 nM and 50 μM were excluded in this research, because our aim is to find structures with better capability to discriminate well-separated actives and inactives. We use a restrict IC50 threshold to separate apart structures (5 nM for actives). If a PDB structure can’t even distinguish more distinctive actives/inactives, we will not consider it as a better selector. However, in the validation section, we consider that usually compounds with inhibition IC50 of 10 uM will be categorized as actives in practice, so we use it as more relaxed threshold. All the co-crystalized inhibitors in the original protein structures are excluded in the ligand test-set. Reasonable molecular weight ranges from 200 Da to 600 Da23. Well defined active/inactives with 1:1 ratio was used in this study. None of the original structure co-crystalized compounds are included in the ligand test set. The searching and generation of compounds test-set are conducted automatically using Python script. As of November 2016, 14 protein targets in our Autophagy Knowledgebase each have at least 20 active and 20 inactive ligands recorded in ChEMBL database (Figure 2). These targets and their 249 associated



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

crystal structures are then used to build the protein targets test-set and the structure test-set in this research. Validation Compound Set Establishment. Compound set for testing the performance of “strong selector” structures is generated from the list of FDAapproved small molecular antineoplastic drugs since 2005 (www.fda.gov). Although 200-600 was used to select the ligand test set, the range of molecular weight for chemicals within the validation compound set was set as from 200 Da to 800 Da (4 FDA drugs are between 600-800 among all the 43 FDA drugs we used). We do this because the original structure selection is supposed to be more “strict” while the validation can be “loose” (and can include more compounds) since both 600 or 800 are reasonable Mw in the molecular docking. With this criterion, 43 qualified compounds were included in the validation compound set. Docking Protocol. Docking is conducted for all ligands in the ligand test-set following the same procedure using the Tripos molecular modeling packages SYBYL-X 8.024. All molecules for docking are built using standard bond lengths and angles from SYBYL-X 8.0 and are then optimized using the Tripos force field. The flexible docking method, named Surflex-dock, is used with all the parameters setting to default. Before docking, co-crystallized ligands and water molecules are removed from the original PDB structures, and hydrogen atoms are added. Chain A of each structure will be selected for docking by default if no other evidence is available. The active binding pocket of each protein is defined by the residues near the cocrystallized active ligands (most are inhibitors) or is automatically generated by the SYBYL-X utility scripts. Most of the active binding pockets generated for protein kinases in the Autophagy Knowledgebase are pockets which contain the ATP binding domain. This is because most of the co-crystalized ligands we used are ATP competitive inhibitors for these structures. Therefore, most of the ligand predictions in our test of ProSelection are predicting the potential ligands which bind to the ATP domain. The scoring function in Surflex-dock, which contains hydrophobic, polar, repulsive, entropic, and solvation terms, was trained to estimate the dissociation constant (Kd) expressed in log(Kd)25. For large-scale compound screening, it is usually not applicable to check each docking pose of individual compound, even those on the top of the list. Therefore, only the best docking score of a compound against each protein structure is kept for assessing and ranking the potential protein


Page 6 of 32

Page 7 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


partners or targets. Selection Metrics. For a protein structure, molecular docking generates higher docking scores against active ligands compared to inactive ligands, statistically. Therefore, the performance of virtual screening based on molecular docking can be quantified by the p-value. This p-value either from the t-test of two independent samples

26

or the Mann-Whitney U test27, which generally measures the statistical

mean difference between the docking scores in the active ligands set and scores in the inactive ligands set. Choosing of the statistical method was determined by the normality of docking score distribution. Specifically, t-test for two independent samples is used in this work only if the docking scores are normally distributed in both active and inactive ligands set. Otherwise, the Mann-Whitney U test will be used to calculate the p-value instead, due to the required normality assumption of the t-test. All statistical tests in this research are one-tailed. A widely used open source python package Scipy (https://www.scipy.org/) is used to conduct the statistical analysis in this research. The detailed selection method is described below: First, the normality of docking score distributions is determined by the D’Agostino and Pearson’s test28, 29, which combines skew and kurtosis to produce an omnibus test of normality (function name in Python: scipy.stats.mstats.normaltest). Protein structures which have normal-distributed docking scores for both its active and inactive ligand sets will be selected to run the t-test for two independent samples with unequal variance at 95% significance level30, (1) Compute the test statistic 𝑡𝑡 =

𝑥𝑥̅1 − 𝑥𝑥̅2

�

𝑠𝑠12 𝑠𝑠22 𝑛𝑛1 + 𝑛𝑛2

(2) Compute the approximate degrees of freedom d’, where

𝑑𝑑 ′ =

𝑠𝑠 2 𝑠𝑠 2 (𝑛𝑛1 + 𝑛𝑛2 )2

1 2 𝑠𝑠12

2

2

𝑠𝑠 2 �𝑛𝑛 � �𝑛𝑛2 � 1 2 𝑛𝑛1 − 1 + 𝑛𝑛2 − 1



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(3) Round 𝑑𝑑′ down to the nearest integer 𝑑𝑑′ ′

If t ≤ 0, then 𝑝𝑝 = 2 × (𝑎𝑎𝑎𝑎𝑎𝑎𝑎𝑎 𝑡𝑡𝑡𝑡 𝑡𝑡ℎ𝑒𝑒 𝑙𝑙𝑙𝑙𝑙𝑙𝑙𝑙 𝑜𝑜𝑜𝑜 𝑡𝑡 𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢 𝑎𝑎 𝑡𝑡𝑑𝑑′′ 𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑟𝑟𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖𝑖)

If t ≥ 0, then 𝑝𝑝 = 2 × (𝑎𝑎𝑎𝑎𝑎𝑎𝑎𝑎 𝑡𝑡𝑡𝑡 𝑡𝑡ℎ𝑒𝑒 𝑟𝑟𝑟𝑟𝑟𝑟ℎ𝑡𝑡 𝑜𝑜𝑜𝑜 𝑡𝑡 𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢𝑢 𝑎𝑎 𝑡𝑡𝑑𝑑′′ 𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑𝑑)

Where 𝑑𝑑′′ is the nearest integer to the approximate degrees of freedom, 𝑥𝑥̅1 and 𝑥𝑥̅2 are

the averages of docking scores for active ligands and inactive ligands, respectively.

Ligand numbers in each ligand set are represented by 𝑛𝑛1 and 𝑛𝑛2 . Variance of the two input samples is referred as 𝑠𝑠12 and 𝑠𝑠22 .

Structures without normally distributed docking scores against either of the ligand

sets will be chosen to run Mann-Whitney U test at 95% significance level31 for which the general calculations after the ranking procedure is shown below: If R1 ≠ n1(n1 + n2 + 1)/2 and there are no ties, then 𝑇𝑇 = ��𝑅𝑅1 −

𝑛𝑛1 𝑛𝑛2 𝑛𝑛1 (𝑛𝑛1 + 𝑛𝑛2 + 1) 1 � − �� (𝑛𝑛1 + 𝑛𝑛2 + 1) 12 2 2

If R1 ≠ n1(n1 + n2 + 1)/2 and there are ties, then 𝑇𝑇 = ��𝑅𝑅1 −

𝑔𝑔

∑𝑖𝑖=1 𝑡𝑡𝑖𝑖 (𝑡𝑡𝑖𝑖2 − 1) 𝑛𝑛1 𝑛𝑛2 𝑛𝑛1 (𝑛𝑛1 + 𝑛𝑛2 + 1) 1 � − �� 𝑛𝑛1 + 𝑛𝑛2 + 1 − (𝑛𝑛1 + 𝑛𝑛2 )(𝑛𝑛1 + 𝑛𝑛2 − 1) 12 2 2

If R1 = n1(n1 + n2 + 1)/2, then T = 0.

Where 𝑅𝑅1 refers to the rank sum in the first sample, 𝑡𝑡𝑖𝑖 refers to the number of observations with the same value in the ith tied group, g is the number of tied groups, and 𝑛𝑛1 , 𝑛𝑛2 represent the number of ligands in the active ligand set and inactive ligand set, respectively. Then the exact p-value can be computed by: 𝑝𝑝 = 2 × [1 −

(𝑇𝑇)]

(𝑇𝑇) represents the cumulative-distribution function (cdf) for a standard normal

distribution and is denoted by normal distribution ( N (0,1) ).

(𝑥𝑥) = Pr(𝑋𝑋 ≤ 𝑥𝑥) where X follows an standard

However, multiple structures are usually available for a single protein target in the PDB database. This number can reach over 50 in some cases. Therefore, Bonferroni


Page 8 of 32

Page 9 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


correction is applied after multiple comparisons in this research because of the relatively large number of comparisons32. In Bonferroni correction, the significance level α is adjusted to

𝛼𝛼

𝑚𝑚

, where m is the number of hypothesis tests33. The equation

used for converting the Mann-Whitney U test or two independent samples t-test pvalue to Bonferroni corrected p-value is shown below: (1 − 𝑝𝑝)𝑘𝑘 = 1 − 𝑝𝑝0

Where 𝑝𝑝 refers to the original p-value, 𝑝𝑝0 is the Bonferroni corrected p-value, k is the number of PDB structures for a target.

Selection of “strong selector” structures is based on the p-value after the Bonferroni correction. Specifically, PDB structures with a p-value larger than 0.05 will be considered as “weak selector” due to their incapability of distinguishing its active ligands from inactive ligands based on the selected docking protocol. In the opposite, structures which can significantly distinguish active ligands from inactive ligands (p < 0.05) will be labeled as “strong selector”. Particularly, most of the structures for kinases in this research are labeled based on their ability to distinguish ATP competitive inhibitors. This is because most of the docking pockets we previous defined for kinases in the docking protocol are ATP binding pockets (See docking protocol section). ProSelection is docking protocol sensitive, which further indicates that the “weak selector” labeled by ProSelection may not be a bad crystal structure if taking other docking pockets (e.g. Pocket other than ATP binding pocket) or other research

purposes

(e.g.

Predicting

a

non-ATP-competitive

inhibitor)

into

consideration. Individual Docking Score Criterion. When conducting the molecular docking on a bunch of protein targets and their associated chemical ligands, how to interpret the docking results may become a major issue. We raised a hypothesis that the distributions of docking scores are different among different protein structures. This assumption indicates that the docking score criterions for active ligands may be specific for each target and protein structure. A suggested docking score threshold for active ligands (SDA) is used in this research as an individual docking score criterion. SDA is generated based on the receiver operating characteristic (ROC) curves34. This calculation can be divided into four steps: 1) ROC curve establishment, 2)



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Screen through the ROC curve and measure the distance from points-on-the-curve to the top-left point on the ROC plot, and, 3) Choose the point with minimum distance, 4) Use the point corresponding docking score as SDA. The ROC curve is established by stepping sequentially through the ranked list of test-set compounds arranged in order of increasing docking scores for each target. Sensitivity and specificity are calculated at each ranking score. Plotting the sensitivity versus (1 - specificity) for all positions in the ranked list will generate the ROC plot. One interesting feature of the ROC curve is the capability of generating the optimized cutoff point to optimize the sensitivity and specificity. In this work, the optimal cutoff point is calculated for each ROC curve using Pathagoras' theorem, which minimizes the distance between the topleft point in the ROC plot and the selected points on the ROC59 (Supplementary Figure 1). However, in order to keep a relatively low false-positive-rate (FPR) during prediction, we require the 𝐹𝐹𝐹𝐹𝐹𝐹 3 ≤ 0.0001 for each point to be considered as a “valid

point” before using the Pathagoras' theorem. This is similar to require at least one

correct prediction in the top 3 predictions when calculating the optimal cutoff points (also means that at least one of the top 3 predictions is correct). The docking scores associated with these final optimized cutoff points were reported as SDAs. Any ligand that has a docking score higher than its corresponding SDA will be considered as an “active” ligand against its associated crystal structure.

RESULTS Distributions of Docking Scores. Figure 3 shows the normality of the docking scores distributions for 249 structures in the structure test-set, which belongs to the 14 different autophagic protein targets. For each single protein structure, the distributions of docking scores can be different between active ligands set and inactive ligand set (Figure 4). Generally, 68% of the structures have normally distributed docking scores against their corresponding inactive ligands. Only 47% of the structures have normally distributed docking scores against their corresponding active ligands. These numbers are not surprising because the inactive compounds usually have more “random” chemical structures when compared to active compounds. According to the distributions of docking scores, the Mann-Whitney U test was applied to 167 (67%) crystal structures while the t-test was applied to 82 structures (33%) structures in the


Page 10 of 32

Page 11 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


structure test-set. Using different statistical tests based on the statistical normality in this work is expected to yield a better statistical accuracy during the structure selection. Structure selection. Either the t-test for two independent samples or the MannWhitney U test is used to generate the statistical p-value, which depends on the normality of docking score distribution for each PDB structure. Structures with a pvalue lower than 0.05 are categorized as “strong selector” structures while others are classified as “weak selector” for the docking protocol we selected. Percentages of two categories of structures for each target in the test-set are shown in Figure 5. This result suggests that for the molecular docking protocol we choose (Sybyl 8.0 Surflex-dock with default parameters used), a significant portion of PDB structures are “weak selector” (Figure 5). These “weak selector” structures may result in unseparated docking scores against active ligands and inactive ligands during docking process. If no modifications are made either for these structures or the docking protocol, high score decoys can be expected, which can be misleading for further research approaches. However, these results of the structure selection, specifically “strong selector” and “weak selector” structures, are docking protocol sensitive. Either the docking parameters setting, pockets definition, active and inactive ligands choosing, or the docking mode selection can contribute to the results, independent of the protein structure itself. If other binding pockets, such as the allosteric binding pockets, are used, different results of selection can be expected. Therefore, the exact results of structural selection in this work can only be referred when the similar docking protocol is conducted. For a closer look at the “weak selector” structures, several examples are used below to show why these structures are “weak selector” in our docking protocol. PDB structures 3NAY and 3NAX. 3NAY and 3NAX are X-ray structures of 3phosphoinositide dependent protein kinase-1 (PDPK1) binding to compound 2 and compound 7 (Supplementary figure 3), which are PDPK1 inhibitors used by Jannik’s group35. In our structure selection, most of the known active ligands for this target are ATP competitive inhibitors. The typical active ATP binding site is used as the docking pocket. As a result, 3NAY is a “strong selector” structure (p = 3.9×10−3)

while 3NAX is labeled as a “weak selector” (p = 0.12). Compound 2 in 3NAY binds



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 32

to the active DFG-in conformation of PDPK1 ATP binding site, like most typical ATP competitive inhibitors (type I). Compound 7 in 3NAX binds in the ATP pocket and induces substantial conformational changes, causing the DFG loop to exist in an “out” conformation (DFG-out) that typifies inactive form kinase inhibitors (type II)35. Compound 7 also disrupts the αC helix which is a unique feature to this compound among inactive form kinase inhibitors35. The conformational change in the αC helix results in catalytic residue Glu-130 being displaced from the active site. It can be concluded that compound 7 binds to the inactive kinase conformation of PDPK1 and induces conformational changes. Therefore, compared to the structure of active ATP binding site, the PDPK1 protein structure co-crystalized with compound 7 is an inactive kinase conformation with some additional conformational changes. For verification,

we

use

another

known

PDPK1

ATP

competitive

inhibitor,

CHEMBL3640476, to dock both the structure 3NAY and 3NAX. In 3NAY, CHEMBL3640476 forms H-bonds to Ala162, Glu166, and Thr222 within the ATPsite of PDPK1 (Figure 6A). Glu130 is also in the active binding site. However, only one H-bond to Glu166 is formed and Glu130 is out of the active site in the 3NAX (Figure 6B). This indicates that the inactive kinase conformation, like 3NAX, is unsurprisingly being a “weak selector” structure for the typical ATP competitive inhibitors (type I) binding and prediction. Similar results can also be found for other pairs of kinase structures like 4L2Y and 4L23 (Supplementary Figure 4), 1MRY/1MRV and 2X39 (Supplementary Figure 5). PDB structures 2YXJ and 4HNJ. 2YXJ and 4HNJ are structures of the Bcl-2-like protein 1 (also named Bcl-xL). Bcl-xL is a potent non-kinase protein inhibitor of cell death and caspases activation. Structure 2YXJ is the Bcl-xL in complex with ABT737 (Supplementary figure 3), a known Bcl-xL inhibitor which belongs to the class of BH3 mimics60. While 4HNJ is the crystallographic structure of BCL-xL domainswapped dimer in complex with PUMA (p53-upregulated modulator of apoptosis), which induced partial unfolding of two α-helix within Bcl-xL BH3 binding pocket39. In our research, 2YXJ is a “strong selector” (p = 0.029) and 4HNJ is a “weak selector” (p = 0.49). This categorization is validated by the docking research where the known BH3 mimics inhibitor of Bcl-xL (CHEMBL2312484) forms two H-bonds to Phe150 in the BH3 binding pocket (Figure 7A). However, no H-bond is formed in the BH3 binding pocket in the 4HNJ (Figure 7B). Therefore, 4HNJ is unsurprisingly a “weak


Page 13 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


selector” for BH3 mimics because of its structural change in the BH3 binding domain induced by PUMA. Individual docking score criterion. Before validating the target-prediction using the structures selected by ProSelection, we raise the hypothesis that each protein structure of a protein target may have its unique distribution of docking scores for even the same set of ligands. This individualism is due to the existence of diverse protein/compound structures and docking protocols in real life. For instance, a docking score of 6.0 may be high enough for those compounds against one structure to be considered as “active” ligands. However, the same score 6.0 may not be sufficient for another structure to determine its active ligands. This phenomenon further indicates that if an individual docking score criterion can be established for each structure, a better understanding of docking score in the molecular docking can be expected to distinguish active ligands from inactive ligands. In this work, the suggested docking score threshold for active ligands (SDA) is used as a score criterion for individual protein structures base on the ROC curve with both a restriction of FPR and an optimized sensitivity versus specificity (see more details in method section). Supplementary figure 2 presents the ROC curves for 4 of the “strong selector” structures in our structure test-set. SDA is then calculated by finding the optimal cutoff points in the ROC curves for each crystal structure. Particularly, SDA is a new concept raised in this work which is designed for a better understanding of the docking scores in the molecular docking research. The SDAs for all the 142 “strong selector” structures in the structure test-set distribute from 4.0 to 10.0 (Supplementary Table 2 and Figure 8). It presents a nonnormal distribution with a mean value of 7.37 and standard deviation of 1.265. Although most of the SDAs fall in the range of 6 to 8, few reach an even higher score of 10.0 which was previously considered as “very high” in molecular docking. On the other hand, 4PPI, one of the PDB structures of the Bcl-2 family protein Bcl-xL, has a low SDA of 3.76. These results are particularly inspiring for active ligands prediction in the future because previously compounds with a docking score lower than 4.0 will tend to be recognized as “inactive ligands” no matter which target the score is against in the Sybyl-X. Validation of the target prediction for selected targets. Target identification of



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

small molecules is essential for unraveling the underlying mechanisms of their bioactivities or side effects. In cancer therapeutics, autophagy is considered as one of the most critical drug resistance mechanisms, due to its attempt at removing damaged cellular components in the drug administrated cancer cells40. As a validation procedure, we used our established ProSelection algorithm to predict the potential autophagy related off-targets for 43 FDA-approved small molecule antineoplastic agents. Only the “strong selector” structures are used. Previously, 13 targets in our protein test-set have at least one “strong selector” structure among all the 14 targets we explored according to the ProSelection (human ABL1 target has no “strong selector”). ChEMBL is then used as the source for searching the experimental bioactivity data (IC50) for the drug binding of these targets. Two autophagy related targets in the protein test-set, insulin receptor (INSR) and histone deacetylase 6 (HDAC6), are selected for validation due to their relatively sufficient data in the ChEMBL database. INSR is commonly considered as a drug target for anti-hyperglycemia and the treatment of diabetes. However, it also regulates the phagophore formation, which is one of the early stages in the macroautophagy. Enhanced insulin signaling will inhibit the autophagy in the initiation process by inducing the downstream PI3K/AKT/mTOR pathway41. A molecular docking study for 43 FDA-approved antineoplastic drugs against insulin receptor was conducted (Figure 9A). In this docking process, only the 18 “strong selector” structures of insulin receptor were used. Among these 43 drugs, six agents have IC50 data for binding in the ChEMBL database (Table 1). Their associated docking scores were summarized in Figure 9B. One of the commonly used IC50 threshold for active compounds for inhibition is 10μM42. Specifically, in this research, one compound can be considered as “active” once it has an IC50 lower than 10 μM (10000 nM). Osimertinib is an “active” compound for INSR in the experiment (Table 1) and it shows a docking score higher than SDA for 8 of total 18 INSR structures in the molecular docking after ProSelection. A similar result was also observed for the “active” drugs Crizotinib, Ceritinib, and Sunitinib. However, for “inactive” compounds like Lapatinib, scores against some structures are unexpectedly higher than SDA (Figure 9B). This incorrect prediction is acceptable since molecular docking has some extent of “randomness”43 and most of the inactive compounds in our test have a lower quantity of high scores


Page 14 of 32

Page 15 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


1

1

( < 3 all available structures) compared to active compounds (≥ 3 all available

structures). On the functional level, drugs that inhibit the INSR will further inhibit the downstream signaling of PI3K/AKT/mTOR pathway. An enhancement of autophagy function can be observed as the final result of this signaling, which is consistent to the autophagy-related drug resistant nature of many antineoplastic agents40. HDAC6 is another target which is closely related to autophagy. Autophagy acts as a compensatory degradation system when the ubiquitin-proteasome system (UPS) is impaired. Histone deacetylase 6 (HDAC6), a microtubule-associated deacetylase that interacts with polyubiquitinated proteins, is an essential mechanistic link in this compensatory interaction49. Although HDAC6 is not required for autophagy activation, it rather controls the fusion of autophagosomes to lysosomes during the autophagy process50. Molecular docking results of 43 FDA-approved antineoplastic drugs against HDAC6 were presented in Figure 10A. Again, only the “strong selector” structures are used. Belinostat, Panobinostat, and Bendamustine are known to bind HDAC6 with IC50 of 15, 11, and 6 nM respectively (Table 2). Belinostat and Panobinostat are successfully predicted by “strong selector” structures while Bendamustine is accidently predicted as “inactive” (Figure 10B). For known HDAC6 “inactive” drugs Lapatinib and Romidepsin, none of the “strong selector” structures show a docking score higher than SDA (Figure 10B). These results indicate that the protein structures selected by ProSelection work accurately during active ligands prediction.

CONCLUSION ProSelection, a computational statistical algorithm designed for generating the proper protein structure subset, is developed for docking-based target identification research. The virtual screening performances of different PDB structures were measured by using targets associated active and inactive ligands with Sybyl-X default Surflex-dock docking protocol. This study examined 249 crystal protein structures from 14 autophagy related targets. Either the t-test or the Mann-Whitney U test was used to distinguish “strong selector” structures from “weak selector” structures based on the normality of the docking score distributions for each protein structure. The



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 32

receiver operating characteristic (ROC) curve was made and suggested docking scores threshold for active ligands (SDA) were generated for each “strong selector” structure according to the ROC curve. The performance of the target prediction using only “strong selector” structures was validated by recently FDA approved small molecule antineoplastic drugs for predicting their potential off-targets among autophagy-related proteins. Protein selection is getting more and more important in the docking-based virtual screening. The structures come from the crystallization environments, which may yield different outcomes even when using the same compound for docking. One approach to incorporate multiple protein structures into the molecular docking is the using of an ensemble consisting of multiple conformations of the same protein17. This method can be categorized as the “ensemble docking” (ED) strategy, which achieves the protein “flexibility” compared to the rigid docking. In rigid docking, protein structure is not able to change its conformation during the docking process. By using ensemble docking, the docking performance can be increased54,

55

. Although ED

strategy can be relatively fast compared to the sequentially docking during large-scale database screening and even more accurate, it could still be “over-fitting” if some “weak selector” structures are taken into consideration during ED. It is common that different crystal structures are generated based on various research purpose. Thus, the “strong selector” and the “weak selector” protein structures may be distinct when facing diverse docking protocols and research aims. Therefore, we developed the ProSelection, an algorithm designed for protein structure selection based on the known ligands and docking protocols. After using the personalized ligands and protocols to run ProSelection, users can put the structures selected by ProSelection directly into their rigid docking or act as a new source for their further ensemble docking approaches. It is worth considering the weaknesses of this current study to explore the directions for further research. ProSelection depends heavily on the abundance of known active and inactive ligands. As shown in Figure 2, the number of active ligands is usually much larger than inactive ligands in the public databases. This ratio can even approach 50 for some targets. One solution to rebalance the number of ligands is replacing the inactive ligands with random ligands. This is expected to slightly reduce the performance of the structure selection because of the existence of active ligands in


Page 17 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


a random compound-set. However, this modification is still acceptable because a “strong selector” structure is also expected to distinguish active ligands from random ligands. More importantly, the number of random compounds is not limited. More protein targets will be available for ProSelection if the active-ligands set and the random-ligands set are used. Specifically, every target, with the size of active ligand test-set larger than 20, is sufficient for ProSelection if a random-ligands set is used. However, for targets which do not have enough active ligands, statistical approaches like ProSelection can’t be used. Therefore, machine learning algorithms can be an alternative solution to predict the proper protein structure subset based on the similar logic of designing ProSelection. Using of the DUD dataset as source is another consideration, however, it will become a new study design and only a limited number of the targets included in the DUD dataset. This means the proteins can only be from those DUD-targets. This is not the original purpose of this study. A further consideration in this work is the collaboration of ProSelection with ensemble docking. Because of the backstage usage of popular docking tools like Sybyl-x, statistical methods like ProSelection are still categorized as “rigid” docking but with flexible ligands. If ensemble docking (which has taken receptor flexibility into consideration) can be integrated into ProSelection, a much better performance is expected because of the consideration of the flexibility in both the protein target and the ligand.

AUTHOR INFORMATION Corresponding Author Author to whom correspondence should be addressed: Dr. Xiang-Qun (Sean) Xie, email: [email protected]; Tel.: +1-412-383-5276; Fax: +1-412-383-7436. Or Dr. Lirong Wang, email: [email protected]; Tel: +1-412-624-8118. Notes The authors declare no competing financial interest. Acknowledgement



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The project is supported by funding to Xie laboratory from the NIH NIDA (P30 DA035778A1) and DOD (W81XWH-16-1-0490).

ASSOCIATED CONTENT The Supporting Information is available which contains the Suggested docking score thresholds for active ligands (SDA) of 142 “strong selector” protein structures, SDA generation, ROC examples of “strong selector” structures, structures of the chemical inhibitors, and the docking comparisons between “strong selector” and “weak selector” for known inhibitors.

REFERENCE 1.

Kuntz, I. D.; Blaney, J. M.; Oatley, S. J.; Langridge, R.; Ferrin, T. E., A geometric approach to macromolecule-ligand interactions. J Mol Biol 1982, 161, 269-88.

2.

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 2010, 31, 455-61.

3.

Jain, A. N., Surflex: fully automatic flexible molecular docking using a molecular similarity-based search engine. Journal of medicinal chemistry 2003, 46, 499-511.

4.

Jones, G.; Willett, P.; Glen, R. C.; Leach, A. R.; Taylor, R., Development and validation of a genetic algorithm for flexible docking. Journal of molecular biology 1997, 267, 727-748.

5.

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., Glide: a new approach for rapid, accurate docking and scoring. 1. Method and assessment of docking accuracy. Journal of medicinal chemistry 2004, 47, 1739-1749.

6.

Rarey, M.; Kramer, B.; Lengauer, T.; Klebe, G., A fast flexible docking method using an incremental construction algorithm. Journal of molecular biology 1996, 261, 470-489.

7.

Ewing, T. J.; Makino, S.; Skillman, A. G.; Kuntz, I. D., DOCK 4.0: search strategies for automated molecular docking of flexible molecule databases. Journal of computer-aided molecular design 2001, 15, 411-428.

8.

Brylinski, M., In; ACS Publications: 2013.

9.

McGann, M., FRED and HYBRID docking performance on standardized datasets. Journal of computer-aided molecular design 2012, 26, 897-906.

10. Grosdidier, A.; Zoete, V.; Michielin, O., Fast docking using the CHARMM force field with EADock DSS. Journal of computational chemistry 2011, 32, 2149-2159. 11. Luo, W.; Pei, J.; Zhu, Y., A fast protein-ligand docking algorithm based on hydrogen bond matching and surface shape complementarity. Journal of molecular modeling 2010, 16, 903-913.


Page 18 of 32

Page 19 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


12. Meier, R.; Pippel, M.; Brandt, F.; Sippl, W.; Baldauf, C., ParaDockS: a framework for molecular docking with population-based metaheuristics. Journal of chemical information and modeling 2010, 50, 879-889. 13. Jain, A. N., Effects of protein conformation in docking: improved pose prediction through protein pocket adaptation. Journal of computer-aided molecular design 2009, 23, 355-374. 14. Lill, M. A., Efficient incorporation of protein flexibility and dynamics into molecular docking simulations. Biochemistry 2011, 50, 6157-6169. 15. Totrov, M.; Abagyan, R., Flexible ligand docking to multiple receptor conformations: a practical alternative. Current opinion in structural biology 2008, 18, 178-184. 16. Yuriev, E.; Ramsland, P. A., Latest developments in molecular docking: 2010–2011 in review. Journal of Molecular Recognition 2013, 26, 215-239. 17. Huang, S. Y.; Zou, X., Ensemble docking of multiple protein structures: considering protein structural variations in molecular docking. Proteins: Structure, Function, and Bioinformatics 2007, 66, 399-421. 18. Craig, I. R.; Essex, J. W.; Spiegel, K., Ensemble docking into multiple crystallographically derived protein structures: an evaluation based on the statistical analysis of enrichments. Journal of chemical information and modeling 2010, 50, 511-524. 19. Rueda, M.; Bottegoni, G.; Abagyan, R., Recipes for the selection of experimental protein conformations for virtual screening. Journal of chemical information and modeling 2009, 50, 186193. 20. Bottegoni, G.; Rocchia, W.; Rueda, M.; Abagyan, R.; Cavalli, A., Systematic exploitation of multiple receptor conformations for virtual ligand screening. PLoS One 2011, 6, e18845. 21. Liu, H.; Wang, L.; Lv, M.; Pei, R.; Li, P.; Pei, Z.; Wang, Y.; Su, W.; Xie, X.-Q., AlzPlatform: an Alzheimer’s disease domain-specific chemogenomics knowledgebase for polypharmacology and target identification research. Journal of chemical information and modeling 2014, 54, 1050-1060. 22. Berman, H. M.; Battistuz, T.; Bhat, T.; Bluhm, W. F.; Bourne, P. E.; Burkhardt, K.; Feng, Z.; Gilliland, G. L.; Iype, L.; Jain, S., The protein data bank. Acta Crystallographica Section D: Biological Crystallography 2002, 58, 899-907. 23. Balius, T. E.; Mukherjee, S.; Rizzo, R. C., Implementation and evaluation of a docking‐rescoring method using molecular footprint comparisons. Journal of computational chemistry 2011, 32, 2273-2289. 24. Geldenhuys, W. J.; Nakamura, H., 3D-QSAR and docking studies on transforming growth factor (TGF)-β receptor 1 antagonists. Bioorganic & medicinal chemistry letters 2010, 20, 1918-1923. 25. Gu, S.-X.; Li, Z.-M.; Ma, X.-D.; Yang, S.-Q.; He, Q.-Q.; Chen, F.-E.; De Clercq, E.; Balzarini, J.; Pannecouque, C., Chiral resolution, absolute configuration assignment and biological activity of racemic diarylpyrimidine CH (OH)-DAPY as potent nonnucleoside HIV-1 reverse transcriptase inhibitors. European journal of medicinal chemistry 2012, 53, 229-234. 26. Paternoster, R.; Brame, R.; Mazerolle, P.; Piquero, A., Using the correct statistical test for the equality of regression coefficients. Criminology 1998, 36, 859-866. 27. Ruxton, G. D., The unequal variance t-test is an underused alternative to Student's t-test and the Mann–Whitney U test. Behavioral Ecology 2006, 17, 688-690. 28. d'Agostino, R. B., An omnibus test of normality for moderate and large size samples. Biometrika 1971, 341-348. 29. D'AGOSTINO, R.; Pearson, E. S., Tests for departure from normality. Empirical results for the distributions of b2 and√ b1. Biometrika 1973, 60, 613-622.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

30. Rosner, B. Two-Sample t Test for Independent Samples with Unequal Variance. In Fundamentals of Biostatistics (7th Edition); Cengage Learning: Boston, MA, 2010, pp 287-293. 31. Rosner, B. The Wilcoxon Rank-Sum Test. In Fundamentals of Biostatistics (7th Edition); Cengage Learning: Boston, MA, 2010; Chapter 339-344. 32. Curtin, F.; Schulz, P., Multiple correlations and Bonferroni’s correction. Biological psychiatry 1998, 44, 775-777. 33. Shaffer, J. P., Multiple hypothesis testing. Annual review of psychology 1995, 46, 561-584. 34. Fawcett, T., An introduction to ROC analysis. Pattern recognition letters 2006, 27, 861-874. 35. Nagashima, K.; Shumway, S. D.; Sathyanarayanan, S.; Chen, A. H.; Dolinski, B.; Xu, Y.; Keilhack, H.; Nguyen, T.; Wiznerowicz, M.; Li, L., Genetic and pharmacological inhibition of PDK1 in cancer cells characterization of a selective allosteric kinase inhibitor. Journal of Biological Chemistry 2011, 286, 6433-6448. 36. Zhao, Y.; Zhang, X.; Chen, Y.; Lu, S.; Peng, Y.; Wang, X.; Guo, C.; Zhou, A.; Zhang, J.; Luo, Y., Crystal structures of PI3Kα complexed with PI103 and its derivatives: new directions for inhibitors design. ACS medicinal chemistry letters 2013, 5, 138-142. 37. McHardy, T.; Caldwell, J. J.; Cheung, K.-M.; Hunter, L. J.; Taylor, K.; Rowlands, M.; Ruddle, R.; Henley, A.; de Haven Brandon, A.; Valenti, M., Discovery of 4-amino-1-(7 H-pyrrolo [2, 3-d] pyrimidin-4-yl) piperidine-4-carboxamides as selective, orally active inhibitors of protein kinase B (Akt). Journal of medicinal chemistry 2010, 53, 2239-2249. 38. Huang, X.; Begley, M.; Morgenstern, K. A.; Gu, Y.; Rose, P.; Zhao, H.; Zhu, X., Crystal structure of an inactive Akt2 kinase domain. Structure 2003, 11, 21-30. 39. Follis, A. V.; Chipuk, J. E.; Fisher, J. C.; Yun, M.-K.; Grace, C. R.; Nourse, A.; Baran, K.; Ou, L.; Min, L.; White, S. W., PUMA binding induces partial unfolding within BCL-xL to disrupt p53 binding and promote apoptosis. Nature chemical biology 2013, 9, 163-168. 40. Mohammad, R. M.; Muqbil, I.; Lowe, L.; Yedjou, C.; Hsu, H.-Y.; Lin, L.-T.; Siegelin, M. D.; Fimognari, C.; Kumar, N. B.; Dou, Q. P. Broad targeting of resistance to apoptosis in cancer. In Seminars in cancer biology, 2015; Elsevier: 2015; Vol. 35; pp S78-S103. 41. Mizushima, N.; Levine, B.; Cuervo, A. M.; Klionsky, D. J., Autophagy fights disease through cellular self-digestion. Nature 2008, 451, 1069-1075. 42. Hildmann, C.; Wegener, D.; Riester, D.; Hempel, R.; Schober, A.; Merana, J.; Giurato, L.; Guccione, S.; Nielsen, T.K.; Ficner, R. and Schwienhorst, A., Substrate and inhibitor specificity of class 1 and class 2 histone deacetylases. Journal of biotechnology 2006, 124, 258-270. 43. Erickson, J. A.; Jalaie, M.; Robertson, D. H.; Lewis, R. A.; Vieth, M., Lessons in molecular recognition: the effects of ligand and protein flexibility on molecular docking accuracy. Journal of medicinal chemistry 2004, 47, 45-55. 44. Marsilje, T. H.; Pei, W.; Chen, B.; Lu, W.; Uno, T.; Jin, Y.; Jiang, T.; Kim, S.; Li, N.; Warmuth, M., Synthesis, structure–activity relationships, and in vivo efficacy of the novel potent and selective anaplastic lymphoma kinase (ALK) inhibitor 5-Chloro-N 2-(2-isopropoxy-5-methyl-4(piperidin-4-yl) phenyl)-N 4-(2-(isopropylsulfonyl) phenyl) pyrimidine-2, 4-diamine (LDK378) currently in phase 1 and phase 2 clinical trials. Journal of medicinal chemistry 2013, 56, 56755690. 45. Cui, J. J.; Tran-Dubé, M.; Shen, H.; Nambu, M.; Kung, P.-P.; Pairish, M.; Jia, L.; Meng, J.; Funk, L.; Botrous, I., Structure based drug design of crizotinib (PF-02341066), a potent and selective dual inhibitor of mesenchymal–epithelial transition factor (c-MET) kinase and anaplastic lymphoma kinase (ALK). Journal of medicinal chemistry 2011, 54, 6342-6363. 46. Mahboobi, S.; Sellmer, A.; Winkler, M.; Eichhorn, E.; Pongratz, H.; Ciossek, T.; Baer, T.; Maier,


Page 20 of 32

Page 21 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


T.; Beckers, T., Novel chimeric histone deacetylase inhibitors: a series of lapatinib hybrides as potent inhibitors of epidermal growth factor receptor (EGFR), human epidermal growth factor receptor 2 (HER2), and histone deacetylase activity. Journal of medicinal chemistry 2010, 53, 8546-8555. 47. Finlay, M. R. V.; Anderton, M.; Ashton, S.; Ballard, P.; Bethel, P. A.; Box, M. R.; Bradbury, R. H.; Brown, S. J.; Butterworth, S.; Campbell, A., In; ACS Publications: 2014. 48. Apsel, B.; Blair, J. A.; Gonzalez, B.; Nazif, T. M.; Feldman, M. E.; Aizenstein, B.; Hoffman, R.; Williams, R. L.; Shokat, K. M.; Knight, Z. A., Targeted polypharmacology: discovery of dual inhibitors of tyrosine and phosphoinositide kinases. Nature chemical biology 2008, 4, 691-699. 49. Pandey, U. B.; Nie, Z.; Batlevi, Y.; McCray, B. A.; Ritson, G. P.; Nedelsky, N. B.; Schwartz, S. L.; DiProspero, N. A.; Knight, M. A.; Schuldiner, O., HDAC6 rescues neurodegeneration and provides an essential link between autophagy and the UPS. Nature 2007, 447, 860-864. 50. Lee, J. Y.; Koga, H.; Kawaguchi, Y.; Tang, W.; Wong, E.; Gao, Y. S.; Pandey, U. B.; Kaushik, S.; Tresse, E.; Lu, J., HDAC6 controls autophagosome maturation essential for ubiquitin‐selective quality‐control autophagy. The EMBO journal 2010, 29, 969-980. 51. Jones, P.; Altamura, S.; De Francesco, R.; Gallinari, P.; Lahm, A.; Neddermann, P.; Rowley, M.; Serafini, S.; Steinkühler, C., Probing the elusive catalytic activity of vertebrate class IIa histone deacetylases. Bioorganic & medicinal chemistry letters 2008, 18, 1814-1819. 52. Fredenhagen, A.; Kittelmann, M.; Oberer, L.; Kuhn, A.; Kühnöl, J.; Délémonté, T.; Aichholz, R.; Wang, P.; Atadja, P.; Shultz, M. D., Biocatalytic synthesis and structure elucidation of cyclized metabolites of the deacetylase inhibitor panobinostat (LBH589). Drug Metabolism and Disposition 2012, 40, 1041-1050. 53. Miller, T. A.; Witter, D. J.; Belvedere, S., Histone deacetylase inhibitors. Journal of medicinal chemistry 2003, 46, 5097-5116. 54. Evangelista, W.; Weir, R. L.; Ellingson, S. R.; Harris, J. B.; Kapoor, K.; Smith, J. C.; Baudry, J., Ensemble-based docking: From hit discovery to metabolism and toxicity predictions. Bioorganic & Medicinal Chemistry 2016, 24, 4928-4935. 55. Ellingson, S. R.; Miao, Y.; Baudry, J.; Smith, J. C., Multi-conformer ensemble docking to difficult protein targets. The Journal of Physical Chemistry B 2014, 119, 1026-1034. 56. Shen, M.; Tian, S.; Pan, P.; Sun, H.; Li, D.; Li, Y.; Zhou, H.; Li, C.; Lee, S. M.-Y.; Hou, T., Discovery of novel ROCK1 inhibitors via integrated virtual screening strategy and bioassays. Scientific reports 2015, 5 57. Tian, S.; Sun, H.; Li, Y.; Pan, P.; Li, D.; Hou, T., Development and evaluation of an integrated virtual screening strategy by combining molecular docking and pharmacophore searching based on multiple protein structures. Journal of chemical information and modeling 2013, 53, 2743-2756. 58. Tian, S.; Sun, H.; Pan, P.; Li, D.; Zhen, X.; Li, Y.; Hou, T., Assessing an ensemble docking-based virtual screening strategy for kinase targets by considering protein flexibility. Journal of chemical information and modeling 2014, 54, 2664-2679. 59. Cleophas, T. J.; Droogendijk, J.; van Ouwerkerk, B. M., Validating diagnostic tests, correct and incorrect methods, new developments. Current clinical pharmacology 2008, 3, 70-76. 60. Van Delft, MF.; Wei, AH.; Mason, KD.; Vandenberg, CJ.; Chen, L.; Czabotar, PE.; Willis, SN.; Scott, CL.; Day, CL.; Cory, S.; Adams, JM. The BH3 mimetic ABT-737 targets selective Bcl-2 proteins and efficiently induces apoptosis via Bak/Bax if Mcl-1 is neutralized. Cancer cell 2006, 10, 389-99.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

FIGURE LEGEND Figure 1. Method design of ProSelection. The number in parenthesis and bracket represents the number of protein targets or structures in the corresponding step, respectively.

Figure 2. Ligand information for the 14 targets in Autophagy Knowledgebase. Most of the targets in this protein test-set have more known active ligands than inactive ligands reported in the ChEMBL database. Abbreviations: PDPK1, 3phosphoinositide-dependent protein kinase 1; B2CL1, Bcl-2-like protein 1; CASP1, Caspase-1; CATB, Cathepsin B; CATD, Cathepsin D; HDAC6, Histone deacetylase 6; INSR, Insulin receptor; SIR1, NAD-dependent protein deacetylase sirtuin-1; SIR2, NAD-dependent protein deacetylase sirtuin-2; PK3CA, Phosphatidylinositol 4,5bisphosphate 3-kinase catalytic subunit alpha isoform; AKT1, RAC-alpha serine/threonine-protein kinase; AKT2, RAC-beta serine/threonine-protein kinase; MTOR, Serine/threonine-protein kinase mTOR; ABL1, Tyrosine-protein kinase ABL1.

Figure 3. Normality test for 249 structures in the structure test-set. P-value was calculated based on the D’Agostino and Pearson’s test (see section on evaluation metrics). Distributions of docking scores which have a p-value higher than 0.05 can be considered as normal distributions. Each dot represents the p value of the distribution of docking scores for one protein structure. The normtest_A represents the targets corresponding active ligands while normtest_N represents the corresponding inactive ligands.

Figure 4. Distributions of docking scores for PDB structure 4JT6. Distribution of docking scores and normal Q-Q plot for (A) active ligands (p < 0.05) and (B) inactive ligands (p = 0.83) against 4JT6, one of the PDB structures for mammalian target of rapamycin (mTOR). The p-value was calculated based on the D’Agostino and Pearson’s test for normality. P-values lower than 0.05 indicates a non-normal


Page 22 of 32

Page 23 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


distribution of the docking scores. For structure 4JT6, docking scores are normally distributed among negative ligands while non-normally distributed among active ligands.

Figure 5. Percentage of the “strong selector” structures among 14 targets in the protein test-set. Crystal structures which have a good performance in distinguishing active ligands from inactive ligands are labeled as “strong selector”. The number on the bar represents the quantity of “strong selector” or “weak selector” structures for each target. Totally 249 PDB structures in structure test-set are examined while 142 structures are “strong selector” and 107 structures are “weak selector”. Among all the 14 targets in this test, thirteen targets have at least one “strong selector” structure according to ProSelection. Abbreviations: PDPK1, 3-phosphoinositide-dependent protein kinase 1; B2CL1, Bcl-2-like protein 1; CASP1, Caspase-1; CATB, Cathepsin B; CATD, Cathepsin D; HDAC6, Histone deacetylase 6; INSR, Insulin receptor; SIR1, NAD-dependent protein deacetylase sirtuin-1; SIR2, NAD-dependent protein deacetylase sirtuin-2; PK3CA, Phosphatidylinositol 4,5-bisphosphate 3-kinase catalytic subunit alpha isoform; AKT1, RAC-alpha serine/threonine-protein kinase; AKT2, RAC-beta serine/threonine-protein kinase; MTOR, Serine/threonine-protein kinase mTOR; ABL1, Tyrosine-protein kinase ABL1.

Figure 6. Compound CHEMBL3640476 binds to the PDPK1 ATP binding site Structures of the compound CHEMBL3640476 binds to the PDPK1 ATP binding site in the kinase domain. (A) CHEMBL3640476 binds to the PDPK1 crystal structure 3NAY (Docking score 8.49). (B) CHEMBL3640476 binds to the PDPK1 crystal structure 3NAX (Docking score 5.35). (C) The original structure of 3NAX (PDPK1 co-crystalized with compound 7). Blue sticks represent the ATP competitive inhibitor of PDPK1, named CHEMBL3640476. Yellow sticks represent the important residues of PDPK1 ATP binding domain which are also labeled. H-bonds to the important residues are displayed as dashed red lines. In 3NAY, CHEMBL3640476 has H-bonds to Ala162, Glu166, and Thr222 within the ATP-site of PDPK1. Glu130 is also in the active binding site. However, only one H-bond to Glu166 is formed while Glu130 is



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

out of the active site in the 3NAX.

Figure 7. Structures of the compound CHEMBL2312484 bind to the BH3 binding domain within the Bcl-xL. (A) CHEMBL2312484 binds to the Bcl-xL crystal structure 2YXJ (Docking score 5.42). (B) CHEMBL2312484 binds to the BclxL crystal structure 4HNJ (Docking score 3.35). Blue sticks represent the inhibitor compound CHEMBL2312484. Yellow sticks represent the important residues in the BH3 binding domain which are also labeled. H-bonds to the important residues are displayed as dashed red lines. In 2YXJ, CHEMBL2312484 has two H-bonds to Phe105. In 4HNJ, however, no H-bond is formed between the compound and the important residues. One residue, Tyr195, is even “out” of the BH3 binding domain due to the protein structural change induced by PUMA in 4HNJ.

Figure 8. Suggested docking score threshold for active ligands (SDA) of 142 “strong selector” structures.

Figure 9. Molecular docking results for FDA-approved drugs against insulin receptor. (A) Docking scores for 43 FDA-approved antineoplastic drugs. (B) Docking scores for 6 FDA-approved drugs with experimental IC50 data of inhibition in ChEMBL. All docking scores are shown in folds compared to the corresponding suggested docking score threshold for active ligands (SDA) for each structure. Sybyl software generates negative scores for some compounds during the docking process. Scores lower than 0 are manipulated to 0 for figure generating. The folds are negative because the original docking scores are below zero. Sybyl 8.0 software sometimes generates negative scores for some compounds during the docking process. But the number of such compounds is small enough and can be ignored.

Figure 10. Molecular docking results for FDA-approved drugs against histone deacetylase 6. (A) Docking scores for 43 FDA-approved antineoplastic drugs. (B) Docking scores for 5 FDA-approved drugs with experimental inhibition data in


Page 24 of 32

Page 25 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ChEMBL. All docking scores are shown in folds compared to the corresponding SDA for each structure. Scores lower than 0 are manipulated to 0 for figure generating. Sybyl 8.0 software sometimes generates negative scores for some compounds during the docking process. But the number of such compounds is small enough and can be ignored.

TABLE Table 1. Experimental bioactivities of 6 drugs against human insulin receptor Drug

Ceritinib

Crizotinib

Lapatinib

Osimertinib

Sunitinib

Topotecan

IC50 (nM)

7

102

17000

912

3200

Not Active

Description

Inhibition

Inhibition

Inhibition

Inhibition

Inhibition

-

Ref

44

45

46

47

48

-

Table 2. Experimental bioactivities for 5 antineoplastic drugs against human histone deacetylase 6. Drug

Belinostat

Bendamustine

Lapatinib

Panobinostat

Romidepsin

IC50 (nM)

15

6

>10000

11

14000

Description

Inhibition

-

Inhibition

Inhibition

Inhibition

Ref

51

BindingDB: 7095

46

52

53



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

FIGURE Figure.1


Page 26 of 32

Page 27 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure.2

Figure.3



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure.4

Figure.5


Page 28 of 32

Page 29 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure.6



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure.7

Figure.8


Page 30 of 32

Page 31 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure.9



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure.10


Page 32 of 32

ProSelection: A novel algorithm to select proper ... - ACS Publications

Recommend Documents