ADMET Evaluation in Drug Discovery. 12. Development of Binary

Department of Biochemistry, The University of Texas Southwestern Medical Center, 5323 Harry Hines Boulevard, Dallas, Texas 75390, United States. Mol. ...
1 downloads 11 Views 742KB Size
Article pubs.acs.org/molecularpharmaceutics

ADMET Evaluation in Drug Discovery. 12. Development of Binary Classification Models for Prediction of hERG Potassium Channel Blockage Sichao Wang,† Youyong Li,† Junmei Wang,‡ Lei Chen,† Liling Zhang,† Huidong Yu,† and Tingjun Hou*,† †

Institute of Functional Nano & Soft Materials (FUNSOM) and Jiangsu Key Laboratory for Carbon-Based Functional Materials & Devices, Soochow University, Suzhou, Jiangsu 215123, China ‡ Department of Biochemistry, The University of Texas Southwestern Medical Center, 5323 Harry Hines Boulevard, Dallas, Texas 75390, United States S Supporting Information *

ABSTRACT: Inhibition of the human ether-a-go-go related gene (hERG) potassium channel may result in QT interval prolongation, which causes severe cardiac side effects and is a major problem in clinical studies of drug candidates. The development of in silico tools to filter out potential hERG potassium channel blockers in early stages of the drug discovery process is of considerable interest. Here, a diverse set of 806 compounds with hERG inhibition data was assembled, and the binary hERG classification models using naive Bayesian classification and recursive partitioning (RP) techniques were established and evaluated. The naive Bayesian classifier based on molecular properties and the ECFP_8 fingerprints yielded 84.8% accuracy for the training set using the leave-one-out (LOO) cross-validation procedure and 85% accuracy for the test set of 120 molecules. For the two additional test sets, the model achieved 89.4% accuracy for the WOMBAT-PK test set, and 86.1% accuracy for the PubChem test set. The naive Bayesian classifiers gave better predictions than the RP classifiers. Moreover, the Bayesian classifier, employing molecular fingerprints, highlights the important structural fragments favorable or unfavorable for hERG potassium channel blockage, which offers extra valuable information for the design of compounds avoiding undesirable hERG activity. KEYWORDS: hERG, ADMET, QSAR, naive Bayesian classification, recursive partitioning



INTRODUCTION The human ether-a-go-go related gene (hERG) protein is a tetrameric potassium channel, which plays an important role in cardiac action potential.1−3 Blockage of hERG potassium ion channel is believed to be the major cause of drug-induced QT syndrome, which can lead to sudden death.4 Additionally, it may cause acquired long QT syndrome, leading to torsades de points (TdP), a severe cardiac side effect which represents a major problem in clinical studies of drug candidates.5 Many major drugs, such as terfenadine, cisapride, sertindole, thioridazine, and grepa-floxacin, were withdrawn from the market as a result of their undesirable hERG-related cardiotoxicity. Therefore, it is important to assess the hERG channel binding of lead compounds early in the preclinical phase of drug discovery. The hERG IC50 value can be measured experimentally through in vitro or in vivo assays. There are many techniques for early cardiotoxicity assessment, such as rubidium-flux assays, radioligand binding assays, in vitro electrophysiology measurements, and fluorescence-based assays.6 However, these in vitro experimental assays for evaluating hERG binding are time-consuming, © 2012 American Chemical Society

expensive, and labor-intensive. Therefore, it is valuable to develop in silico models that provide rapid and cheap screening platforms for identifying hERG inhibitors in the early stage of drug design and optimization. Many computational models have been developed by structure-based and ligand-based approaches to discriminate hERG blockers from nonblockers. Since the crystal structure of hERG potassium channel is not available, homology modeling techniques, in combination with molecular docking and free energy calculations, were used to predict the binding capabilities of drug candidates to hERG potassium channel. For example, Farid and co-workers docked five sertindole analogues into the active site of the homology model of hERG channel using KvAP as a template,7 and the predicted binding affinities show good relationship with the experimental values. Received: Revised: Accepted: Published: 996

January 12, 2012 February 17, 2012 March 1, 2012 March 1, 2012 dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

based on molecular properties and fingerprints were developed and validated by the external test sets. The comparison studies show that naive Bayesian classification outperformed RP in respect to prediction accuracy. Moreover, the naive Bayesian classifiers identified the key structural features for differentiating hERG active and inactive.

Rajamani et al. docked 27 known hERG channel binders into the active sites of the multiple homology models of hERG channel,8 and they found that the flexibility was essential to rank the binding affinities of a set of diverse ligands correctly. However, the structure-based studies were only based on very limited analogue data sets, and the prediction capability based on the homology model has still not been validated by large data sets with high structural diversity. It is still difficult to successfully utilize hERG potassium channel homology models for prospective molecular docking studies. Due to the lack of a crystal structure of hERG potassium channel, the ligand-based prediction models are still the mainstream for the predictions of hERG potassium channel blockers. In 2002, Ekins et al. first built a preliminary hERG pharmacophore model based on a training set of 11 molecules.9 This pharmacophore model with four hydrophobic and one positive ionizable points could successfully rank 15 molecules in a test set. In 2002, Cavalli et al. established a pharmacophore model from the CoMFA analysis of a set of 31 QT-prolonging drugs,10 and this pharmacophore model contained a positively charged tertiary amine flanked by three aromatic or hydrophobic centers. In 2003, Pearlstein and co-workers performed CoMSIA analysis on a data set of 10 structurally diverse hERG blockers and 22 sertindole analogues,11 and the results given by the CoMSIA analysis were consistent with the structural information given by a homology model of hERG. In 2005, Cianchetta et al. developed correlation models using pharmacophore-based GRIND descriptors for the data sets with and without a basic nitrogen, and the models showed satisfactory predictions for the training and test sets.12 In 2006, Sun developed a naive Bayesian classifier based on a large training set of 1979 compounds, and this model exhibited an ROC accuracy of 0.87 and could correctly classify 58 molecules from 66 tested molecules.13 In 2008, Thai and co-workers developed binary classification models based on P_VSA descriptors, and the prediction accuracy of the models was 82−88%.14 In 2006, Song and Clark calculated 2D fragmentbased descriptors for a training set of 71 compounds and a test set of 19 compounds,15 resulting in a hERG QSAR model with q2 values of 0.912 and 0.848 for the training and test sets, respectively, developed by linear support vector regression. In 2008, Li et al. developed a hERG classification model based on GRIND descriptors and support vector machine (SVM),16 and this model achieved good prediction accuracy for the training set (up to 92%); however, it could not give satisfactory prediction for the test set (global accuracy is ∼60%). Recently, Su and coworkers assembled a data set of 250 structurally diverse compounds with hERG activity, and they developed a regression model and a binary classification model based on 204 traditional 2D descriptors and 76 3D VolSurf-like descriptors.17 This binary model achieved 91% accuracy for the training set, but it could not give satisfactory prediction for the test sets. In summary, most reported data sets have a limited number of compounds (less or close to 300). Although Sun used a large data set of 1979 molecules in modeling,13 the data set is not an open data source. The broad multispecificity of hERG and the lack of an extensive data set of hERG activities become two barriers for establishing accurate prediction models. Here we present a large data set of 806 molecules that are categorized into the blocker and nonblocker classes, and we examined the impact of eight important molecular properties widely used in ADME prediction on hERG blocking. Then, utilizing the new data set, the naive Bayesian and recursive partitioning classifiers



METHODS AND MATERIALS 1. Preparation of Data Set. The data set has 806 molecules in total, and about 60% of the data (495) was collected by Li and co-workers.16 The second data source is an external test set of 66 compounds with hERG information extracted from the WOMBAT-PK database. In addition, we updated the data set with an extra set of 245 molecules collected from recent publications.6,14,15,18−37 The hERG blocking activity was determined primarily by the IC50 measurements using mammalian cell lines, such as HEK, CHO, and COS. When mammalian cell line data were not available, IC50 measurements from a nonmammalian cell line, XO (Xenopus laevis oocytes), were included. An ideal training set needs to be extensive, consistent, and diverse. But since the purpose of our work is to create a qualitative rather than a quantitative model for predicting hERG blockage, the variation of the activities among different cell lines may be tolerated.38 The distributions of molecular weight and IC50 shown in Figure S1 in the Supporting Information indicate that the data set is relatively diverse and will not lead to data-biased prediction models. The molecules in the data set were minimized in MOE by using molecular mechanics (MM) with the MMFF94 force field.39 The data set is available from the supporting Web site: http://cadd.suda.edu.cn/admet. We used different thresholds of 1, 5, 10, 20, 30, and 40 μM, respectively, to find the most suitable threshold to identify the hERG blockers and nonblockers in our data set. The threshold up to 50 μM was not considered because it was not typically relevant based upon previous investigation.16 First, the whole data set was split into three parts: a training set of 620 molecules randomly selected from the data set, a test set (test set I) of 120 molecules randomly selected from the data set, and another test set (test set II) with 66 molecules from the WOMBAT-PK database. The classification models were developed based on the training set and validated by the two test sets. Moreover, in order to give a more stringent validation, the data set of 1953 molecules with hERG binding activity from the PubChem bioassay database was used. 2. Molecular Descriptors. Fourteen molecular descriptors widely used in ADME predictions40,41 were used in our study. The descriptors include octanol−water partitioning coefficient (A log P) based on the Ghose and Crippen method,42 apparent partition coefficient at pH = 7.4 (log D) based on the Csizmadia method,43 molecular solubility (log S) based on the multiple linear regression model developed by Tetko et al.,44 molecular weight (MW), the number of hydrogen bond donors (nHBD), the number of hydrogen bond acceptors (nHBA), the number of rotatable bonds (nrot), the number of rings (nR), the number of aromatic rings (nAR), the sum of oxygen and nitrogen atoms (nO+N), polar surface area (PSA), molecular fractional polar surface area (MFPSA), and molecular surface area (MSA). All descriptors can be divided into two classes: physiochemical properties (A log P, log D, log S, MW, nHBD, nHBA, nR, nAR and nO+N), and geometry-related descriptors (PSA, MFPSA and MSA).45 All the descriptors were calculated by using Discovery Studio molecular simulation package.46 997

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

3. Calculation of Molecular Fingerprints. The SciTegic extended-connectivity fingerprints (ECFP, FCFP and LCFP) and Daylight-style path-based fingerprints (EPFP, FPFP and LPFP) were generated using a variant of the Morgan algorithm.47 The initial assignment of atom identifiers has different rules (F, E or L): F represents atom’s functional role code, E represents properties used in the Daylight atomic invariants rule, and L represents A log P atom type code. The number of connections to an atom, the element types, the charge, and the atomic mass form the atom type code. In order to denote the type of fingerprint, the second character, C or P, in the fingerprint name was used. Character C represents extended-connectivity fingerprints, and P represents path-based fingerprints. A fingerprint class is followed by an underscore and a number, which designates the effective diameter of the largest feature. Here for each fingerprint class, three diameters, 4, 6 and 8, were considered. The smaller diameter, 2, was not used because the structural fragments based on the diameter of 2 are so small and too general. The detailed descriptions of these fingerprints were reported in the literature.48 The fingerprints we used here are significantly different from the substructures in the prediction of ADME properties used previously.49−52 First, the fingerprints used here represent a much larger set of features than a set of predefined substructures. Furthermore, these fingerprints are generated directly from molecules, and they do not need to be preselected or predefined. Therefore, novel molecular classes are handled as easily as the common classes. Here, we generated the structural fingerprints by using Discovery Studio molecular simulation package.46 4. Naive Bayesian Classifiers. The naive Bayesian classification technique was used to develop classifiers to discriminate between hERG blockers and nonblockers. Bayes’s theorem relates the conditional and marginal probabilities of events A and B, provided that the probability of B does not equal zero. Its simple form is shown by eq 1: P(A|B) =

P(B|A)P(A) P(B)

probability of the given descriptors that will occur in the training set. Then, we get an assumption that each feature, Ai, is conditionally independent of every other feature Aj. The mathematical procedure to train a naive Bayesian classifier was described previously.53 An advantage of the naive Bayesian classification is that it requires a small amount of training data to estimate the parameters (means and variances of the variables) necessary for classification. Moreover, Bayesian classification can process large amounts of data, learn fast, and be tolerant of random noise. The naive Bayesian classifiers were developed in Discovery Studio molecular simulation package.46 5. Recursive Partitioning Classifiers. RP is a statistical method for multivariable analysis that can create a decision tree to classify the members in a data set based on predefined variables. Models are constructed by successively splitting a data set into increasingly homogeneous subsets until it is infeasible to continue, based on a set of “stopping rules”. At each splitting point or node, the RP algorithm searches a pool of independent variables (i.e., descriptors) and identifies a single variable and the corresponding splitting value that best purifies the group of compounds entering the node. The splitting process continues until either no further improvement can be achieved or the number of compounds in each purified group is too small to justify further splitting. 5-fold crossvalidation was used to determine the degree of pruning required for the best predictive performance. The minimum number of samples at each node was set to 7, and the maximum tree depth was changed from 2 to 10 systematically in order to evaluate the effect of the tree depth on the prediction. The decision trees were created in Discovery Studio molecular simulation package46 based on the training set and validated by three different test sets. 6. Validating the Prediction Accuracy of the RP and Bayesian Models. The quality of the Bayesian and RP classifiers was measured by the quantity of true positives (TP), true negatives (TN), false positives (FP), false negatives (FN), sensitivity (SE), specificity (SP), the prediction accuracy for active (or blockers) (Q+), the prediction accuracy for nonactive (or nonblockers) (Q−), the global accuracy (GA), and the Matthews correlation coefficient (C):

(1)

A naive Bayesian classifier is a simple probabilistic classifier based on applying Bayes’s theorem with strong (naive) independence assumptions. A more descriptive term for the underlying probability model is “independent feature model”. Abstractly, the probability model for a classifier is a conditional model P(C|A1,...,An) over a dependent class variable C with a small number of outcomes or classes, conditional on feature variables A1 through An. In our study, C represents the hERG activity class of a molecule: blocker class (+) or nonblocker class (−). A1 to An represent the calculated values for the feature variables (molecular properties and fingerprints). When Bayes’s theorem is used, we get P(A1, ..., A n|+ )P( +) P( +|A1, ..., A n) = P(A1, ..., A n)

SE =

TP TP + FN

(4)

SP =

TN TN + FP

(5)

PRE1 =

TP TP + FP

(6)

PRE2 =

TN TN + FN

(7)

GA =

(2)

In plain English, the above equation is written as likelihood × prior posterior = evidence

C= (3)

TP + TN TP + TN + FP + FN

(8)

TP × TN − FN × FP (TP + FN)(TP + FP)(TN + FN) (TN + FP) (9)

So in eq 2, P(A1,...,An|+) is the conditional probability of a particular compound being classified as hERG avtive; P(+) is the prior probability, a probability induced from a set of compounds in the training set; P(A1,...,An) is the marginal

The values of GA and C are two important indicators for the classification accuracy of models. The above quantities were calculated for both training and test sets. 998

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

Figure 1. Distributions of eight molecular properties, including A log P, log D, log S, MW, PSA, nrot, nHBD and nHDA, for the active and inactive classes.



RESULTS AND DISCUSSIONS

used in ADME predictions. Eight molecular properties include MW, log S, log D7.4, A log P, PSA, nHBA, nHBD and nrot. The distributions of these eight molecular properties for the active and inactive classes at the threshold of 30 μM are shown in Figure 1. We evaluated the significance of the difference between the means by Student’s t-test. As a complementary test, the linear correlations between each of these eight molecular properties and the IC50 values of 557 compounds are shown in Figure 2.

1. The Impact of Eight Molecular Properties on hERG Binding. For the predictions of ADME properties,41,50,52,54−67 many molecular descriptors have been proven to be helpful, and they were used to describe a lot of molecular properties, such as lipophilicity, hydrogen bonding ability, molecular flexibility, molecular weight, etc. Here, we systematically examined the relationships between eight molecular properties that are widely 999

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

Figure 2. Correlations between eight molecular properties, including A log P, log D, log S, MW, PSA, nrot, nHBD and nHDA, and pIC50 values.

In the eight molecular properties we studied here, A log P, log D and log S are related to hydrophobicity of a molecule. A log P is distributed between −5.339 and 11.512, with a mean of 1.856, log D between −5.122 and 11.512, with a mean of 2.186, and log S between −14.751 and 2.218, with a mean of −5.226. The mean values of A log P for the 369 inactive versus the 344 active are 1.334 and 2.416, respectively; it suggests that, with the rise of hydrophobicity, hERG binding capability increases. The p-value, which measures the difference in the mean of two distributations, is 3.26e−18 for A log P at the 95% confidence

level, indicating that the distributions of actives and inactives are significantly different. Compared with A log P, log D shows better classification capability because of a smaller p-value of 5.14e−23 as shown in Figure 1. It is obvious that the actives tend to be more lipophilic than the inactive agents. Besides hERG channel, it is a general phenomenon that increasing lipophilicity is usually favorable for the binding of the studied molecules to protein receptors.21 The mean values of log S for the 369 inactives versus the 344 actives are −4.01 and −6.53, respectively. It is interesting that log S shows even better 1000

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

0.696 0.664 0.634 0.777 0.728 0.695 0.689 0.665 0.691 0.684 0.697 0.667 0.700 0.703 0.684 0.685 0.668 0.653 0.883 0.850 0.850 0.892 0.858 0.858 0.842 0.833 0.850 0.842 0.842 0.833 0.850 0.850 0.842 0.842 0.833 0.825

GA PRE2

0.954 0.975 0.952 0.958 0.969 0.907 0.898 0.859 0.929 0.867 0.923 0.852 0.860 0.887 0.857 0.870 0.855 0.797 0.697 0.610 0.622 0.796 0.727 0.778 0.787 0.804 0.740 0.817 0.779 0.814 0.841 0.821 0.828 0.818 0.815 0.857 0.892 0.828 0.849 0.872 0.808 0.872 0.803 0.833 0.833 0.825 0.762 0.825 0.831 0.797 0.814 0.797 0.797 0.864 0.852 0.926 0.852 0.929 0.952 0.833 0.889 0.833 0.881 0.860 0.930 0.842 0.869 0.902 0.869 0.885 0.869 0.787 83 77 79 68 63 68 53 55 65 52 48 52 49 47 48 47 47 51 10 16 14 10 15 10 13 11 13 11 15 11 10 12 11 12 12 8 4 2 4 3 2 7 6 9 5 8 4 9 8 6 8 7 8 13 23 25 23 39 40 35 48 45 37 49 53 48 53 55 53 54 53 48 0.733 0.622 0.657 0.590 0.617 0.616 0.561 0.544 0.599 0.663 0.653 0.652 0.696 0.671 0.679 0.623 0.659 0.678 0.915 0.880 0.875 0.823 0.826 0.830 0.775 0.783 0.796 0.830 0.818 0.827 0.848 0.833 0.838 0.811 0.829 0.836 0.966 0.945 0.971 0.884 0.911 0.903 0.889 0.839 0.899 0.883 0.911 0.860 0.864 0.878 0.877 0.827 0.848 0.795 0.722 0.638 0.603 0.696 0.679 0.692 0.654 0.699 0.682 0.776 0.741 0.789 0.832 0.793 0.803 0.796 0.813 0.889 0.929 0.907 0.874 0.859 0.832 0.847 0.730 0.807 0.759 0.800 0.745 0.821 0.842 0.790 0.803 0.783 0.800 0.903 0.850 0.757 0.879 0.742 0.813 0.791 0.849 0.742 0.858 0.867 0.909 0.833 0.855 0.880 0.876 0.838 0.858 0.769 458 447 431 358 347 353 273 302 284 264 246 271 261 245 249 235 240 271 35 46 62 59 70 64 101 72 90 66 84 59 49 65 61 65 60 29 16 26 13 47 34 38 34 58 32 35 24 44 41 34 35 49 43 70

PRE1

classification capability than log P and log D because its p-value is 3.85e−49. In Figure 2, we also observed similar results that log S showed a better linear correlation (r = −0.480) with IC50 than log P (r = 0.28) and log D (r = 0.338). Based on the above analysis, we can draw a conclusion that log S is better to predict hERG blockage than log D and A log P. However, the two distributions of log S for the actives and inactives are still strongly overlapped. Our results are interesting because it appears that hERG blockage does not have direct connection with solubility. However, solubility is also an important indicator of hydrophobicity.52 Therefore, it is quite possible that log S gives better description for the hydrophobicity of the studied molecules. As seen from Figure 1, molecular weight also demonstrates a high impact on hERG binding. The mean values for the actives and inactives are 405.77 and 313.65, respectively. At the 95% confidence level, the p-value of the difference of MW for the two classes is 5.47e−27. It means that MW also has great impact on hERG activity, but its ability to discriminate the actives from the inactives is worse than that of log S. Moreover, in Figure 1, the descriptor, nrot, is a less desirable classifier itself on hERG binding because the p-value associated with the difference in the mean nrot values of the two classes is 1.74e−14, even less than that of log D (5.20e−21). nrot was used as a descriptor in many cases to characterize the bulkiness of a molecule because a larger molecule usually has more rotatable bonds.68 But in this study, nrot does not show an obvious correlation with molecular weight (r = −0.13). According to the former analysis, if molecular weight of a molecule is less than 250, it is less likely to be a hERG blocker,69 which is consistent with the analysis shown in Figure 1. Almost all hERG blockers have molecular weight larger than 250; however, most nonblockers also have molecular weight larger than 250. Therefore, molecular weight of 250 is not a good indicator for hERG binding. The other three descriptors (PSA, nHBD and nHDA) to characterize electrostatic or H-bond features show different capability to discriminate the actives from the inactives. The pvalues associated with the difference in the mean values of the two classes are 7.836e−6, 0.0041 and 2.69e−4 for PSA, nHBD and nHDA, respectively. In Figure 2, very low correlations between IC50 values and PSA (r = −0.157), nHBD (r = 0.0378) and nHDA (r = −0.164) were observed. Our findings are not consistent with the reported pharmacophore models of hERG blockers because all these pharmacophore models have H-bond acceptor features. These three descriptors, PSA, nHBD and nHDA, are very general and have low capacity to characterize the spatial constraint of pharmacophore features imposed by a pharmacophore model. 2. Naive Bayesian Classifiers. The analysis in the previous section shows that a single molecular property is not a good classification criterion for hERG binding. For this reason, a statistical technique, naive Bayesian classification, popularly used in bioinformatic analysis, was applied to develop classification models. We divided our data set into three parts, a training set of 620 molecules, a test set of 120 molecules (test set I) and an external test set of 66 molecules from the WOMBAT-PK database (test set II). For the threshold to class hERG blocker or nonblocker, different values of IC50 have been proposed. Aronov and Goldman adopted a threshold of 40 μM;38 Sun selected a threshold value of 30 μM;13 Tobite et al. used IC50 thresholds of 1 and 40 μM;19 Roche et al. opted for threshold values of 1 and 10 μM.21 Here, 6 levels of thresholds were used,

91 81 94 135 148 144 191 167 193 229 240 220 242 249 248 254 260 233 40

30

20

10

BC-1 BC-2 BC-3 BC-4 BC-5 BC-6 BC-7 BC-8 BC-9 BC-10 BC-11 BC-12 BC-13 BC-14 BC-15 BC-16 BC-17 BC-18 1

5

test

SP SE TN FP FN TP C GA PRE2 PRE1

training

SP SE TN FP FN TP fingerprints model thresholds (μM)

Table 1. The Performance of the Bayesian Classifiers for the Training and Test Sets Based on Different Thresholds

Article

LPFC_8 ECFP_4 EPFC_8 FCFP_8 LCFC_6 SEFP_4 FCFP_6 LCFP_8 FCFC_6 FCFP_8 LCFP_6 FCFP_4 ECFP_8 LCFP_6 FCFP_6 LCFP_4 LCFC_4 ECFC_8

C

Molecular Pharmaceutics

1001

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

Table 2. Mean and Standard Deviation (SD) of Global Accuracy (GA) and C Values for the Baeysian Classifiers Based on Different Thresholds training set

test set C

GA

C

GA

threshold (μM)

av

SD

av

SD

av

SD

av

SD

1 5 10 20 30 40

0.829 0.806 0.776 0.826 0.831 0.819

0.0470 0.0368 0.0340 0.0286 0.0254 0.0260

0.558 0.576 0.549 0.653 0.663 0.643

0.0751 0.0694 0.0650 0.0575 0.0507 0.0524

0.800 0.809 0.774 0.785 0.794 0.792

0.0431 0.0471 0.0405 0.0407 0.0424 0.0405

0.572 0.618 0.553 0.573 0.591 0.588

0.0655 0.1050 0.0858 0.0828 0.0836 0.0795

Figure 3. The comparison of global accuracy of the Bayesian classsifiers using different thresholds based on molecular properties and fingerprints for the (a) training and (b) test sets.

Figure 4. The distributions of the Bayesian scores given by the Bayesian classifier based on molecular properties and the ECFP_8 fingerprint set for the active and inactive classes for the (a) training and (b) test sets. The Bayesian scores for the training set were obtained by using the LOO cross-validation process.

and they are 1, 5, 10, 20, 30, and 40 μM, respectively. It must be noted that some compounds were removed from the training set as those compounds are ungroupable for certain thresholds. Take dexrazoxane (IC50 > 30) as an example, and it is not included in the training set when a threshold of 40 μM is used. To find out the best threshold and fingerprints for hERG classification, 36 classifiers based on the 14 molecular properties and one type of structural fingerprint set were developed. The 36 types of structural fingerprint sets include

ECFC, ECFP, EPFC, EPFP, FCFC, FCFP, FPFC, FPFP, LCFC, LCFP, LPFC and LPFP with three different diameters, 4, 6, and 8. The statistical results for these Bayesian classifiers are summarized in Table 1. For the training set, among all the fingerprint sets, the best Bayesian classifier based on the 14 molecular properties and the LPFC_8 fingerprint set has a sensitivity of 85.0%, a specificity of 92.9%, a prediction accuracy 1002

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

0.792 0.817 0.817 0.792 0.808 0.808 0.800 0.808 0.808

C GA PRE2

0.783 0.849 0.849 0.870 0.875 0.860 0.843 0.860 0.860 0.800 0.791 0.791 0.743 0.764 0.771 0.768 0.771 0.771

PRE1 SP

0.797 0.763 0.763 0.678 0.712 0.729 0.729 0.729 0.729 0.787 0.869 0.869 0.902 0.902 0.885 0.869 0.885 0.885

SE TN

47 45 45 40 42 43 43 43 43 12 14 14 19 17 16 16 16 16

FP FN

13 8 8 6 6 7 8 7 7 48 53 53 55 55 54 53 54 54

TP C

0.547 0.636 0.649 0.698 0.730 0.738 0.738 0.738 0.738 0.774 0.818 0.825 0.847 0.863 0.868 0.868 0.868 0.868

GA PRE2

0.772 0.839 0.837 0.895 0.902 0.892 0.892 0.892 0.892 0.777 0.797 0.812 0.804 0.828 0.845 0.845 0.845 0.845

PRE1 SP

0.806 0.806 0.826 0.800 0.829 0.852 0.852 0.852 0.852 0.739 0.830 0.823 0.898 0.901 0.887 0.887 0.887 0.887

SE TN

250 250 256 248 257 264 264 264 264 60 60 54 62 53 46 46 46 46

FP FN model

RP-1 RP-2 RP-3 RP-4 RP-5 RP-6 RP-7 RP-8 RP-9

tree depth

2 3 4 5 6 7 8 9 10

TP

74 48 50 29 28 32 32 32 32

test set training set

Table 3. The Performance of the RP Models Based on the ECFP_8 Fingerprint Set for the Training and Test Sets Based on Different Tree Depths 1003

209 235 233 254 255 251 251 251 251

of 72.2% for the active class, a prediction accuracy of 96.6% for the inactive class and a GA of 91.5% based on a leave-one-out (LOO) cross-validation when the threshold is 1 μM. But at different thresholds the fingerprint sets that achieve the highest global accuracies are not always the same. For example, the fingerprint sets in the best five Bayesian classifiers based on the other five thresholds are EPFC_8 (5 μM, GA = 87.1%), EPFC_8 (10 μM, GA = 84.3%), LPFC_8 (20 μM, GA = 90.4%), LPFC_6 (30 μM, GA = 89.0%), and LPFC_8 (40 μM, GA = 88.7), respectively. All the Bayesian classifiers were then validated by the predictions on the test set I of 120 molecules. As shown in Table 1, the classifier that gives the best prediction for this test set is also based on the threshold of 1 μM. In addition, we find that if we choose GAtraining ≥ GAtest ≥ 0.85 as a criterion for good prediction, only the classifier based on the threshold of 1 μM can satisfy this criterion. The best prediction of the Bayesian classifier based on the threshold of 1 μM is understandable. 1 μM is a relatively low threshold for hERG binding; therefore, the number of the nonblockers in both the training and test sets are much larger than that of the blockers. The classifiers based on the highly biased data set tend to give better predictions for the nonblocker class and yield better global accuracy. At each threshold, the average and standard deviation of the GA and C values of the 36 classifiers are shown in Table 2. For the training set, the classifiers based on the threshold of 30 μM perform best judged by the mean value and standard deviation of GA. For the test set, the classifiers based on the threshold of 5 μM give the highest mean GA, but the standard deviation of GA also reaches the highest. The 5 μM models have the best mean C value, but their standard deviation is much higher than those of the models based on the other five thresholds, so they may be not stable models for classification. In our study, taking consideration of both GA and C of the training and test sets, the threshold of 30 μM may be the best choice. For the training and test sets, we compared the performance of the classifiers based on different fingerprint sets with or without 14 molecular properties, and the result is shown in Figure 3. Obviously, compared with the classifiers only based on molecular properties, the addition of structural fingerprints can improve the classification accuracy significantly. The prediction accuracy of the Bayesian classifier based on the 14 molecular properties and the ECFP_8 fingerprint set to distinguish blockers from nonblockers was evaluated with two bimodal histograms of the training and test sets (Figure 4). The histograms show that the blockers tend to be more positive while the nonblockers tend to be more negative. Compared with the histograms of the molecular properties shown in Figure 1, the Bayesian classifier separates two classes substantially better. The overlapped area is much smaller than that shown in Figure 1. For the training set, the compounds in both classes overlap between −10 and 10. So the region between −10 and 0 can be defined as the “uncertain zone”. When the Bayesian score of a molecule is in the uncertain zone, the prediction for this molecule is not reliable. For the 120 molecules in the test set, 28 of them are located in the uncertain zone. If these 28 molecules without accurate predictions are eliminated from the test set, the sensitivity and specificity are improved from 86.9% and 83.0% to 88.9% and to 92.1%. 3. Recursive Partitioning Models. Compared with the Bayesian classifiers, the decision trees given by RP are more “visible”, because a RP model can be converted to simple

0.583 0.636 0.636 0.596 0.626 0.623 0.604 0.623 0.623

Article

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

the RP model becomes much better. The GA values of the training and test sets are 86.8% and 80.8%, respectively. The depth of decision tree plays a very important role in RP analysis, and it indicates the complexity of a decision tree. So the tree depth should be carefully calibrated based on the predictions on the external test set. Here, the tree depth was changed from 2 to 10 and the corresponding performance of the models for the training and test sets was evaluated (Table 3). The change of the GA values versus the tree depth is shown in Figure 5. For the training set the value of GA increases with the increase of the tree depth; however for the test set the value of C does not always increase. For the test set, the tree depth of 3 or 4 reaches the best prediction. So according to our analysis, the depth of 4 is the best choice. The best RP model gives a sensitivity of 0.869, a specificity of 0.763, a classification accuracy of 0.791 for the blocker class, a classification accuracy of 0.849 for the nonblocker class, and C of 0.636 for the test set. Based on the same descriptors and data set, the naive Bayesian classifier performs much better than the RP classifier. 4. Model Validation Using the Additional Test Sets. To evaluate and validate the actual prediction capability of the Bayesian classifier, two additional data sets were used. The test set II from the WOMBAT-PK database has 66 molecules. The six best Bayesian classifiers based on six different thresholds were used to evaluate the actual prediction accuracy for this test set (Table 4). Here we built the Bayesian classifiers based on the data set of 740 molecules (the training set plus the test set I), not only the separated training set of 620 molecules. For test set II, using the threshold of 30 μM, the Bayesian classifier reaches a GA of 89.4%, a specificity of 53.8% and a sensitivity of 98.1%, respectively; while the Bayesian classifier reaches a GA of 89.4%, a specificity of 72.7% and a sensitivity of 92.7%, respectively, using the threshold of 40 μM. Among all the classifiers shown in Table 4, five of them achieve high prediction performance for the tested molecules, indicated

Figure 5. The change of GA versus the tree depth for the training and test sets. The RP model was constructed based on molecular properties and the ECFP_8 fingerprint set.

hierarchical rules and it is easier to be interpreted. Here the threshold of 30 μM was used. First, a decision tree based on the 14 molecular properties without fingerprints was constructed and evaluated, and the tree depth was set to 7. The decision tree with a depth of 7 is shown in Figure 6. For the training and test sets, the GA and C values are 84.5% and 79.2%, and 69.0% and 58.5%, respectively, indicating that the prediction accuracy for the test set is much worse than that for the training set. Then, molecular fingerprints, together with molecular properties, were used simultaneously as the descriptors in RP analysis. Here the fingerprint ECFP_8 set was added into RP modeling. Obviously, after the addition of fingerprints, the performance of

Figure 6. Decision tree to classify compounds into the blocker and nonblocker classes created by RP (tree depth is set to 7). 1004

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

1 5 10 20 30 40

ECFC_8 FCFP_8 LCFC_6 LCFC_4 ECFC_4 ECFC_4

fingerprints

TP 106 180 240 278 279 289

28 44 39 43 65 75

FN 54 76 101 79 48 40

FP

TN 532 419 339 314 321 319

SE 0.791 0.804 0.860 0.866 0.811 0.794

SP 0.908 0.846 0.770 0.799 0.870 0.889

0.663 0.703 0.704 0.779 0.853 0.878

PRE1 0.950 0.905 0.897 0.880 0.832 0.810

PRE2

GA 0.886 0.833 0.805 0.829 0.842 0.841

C 0.654 0.629 0.615 0.662 0.683 0.685

13 31 42 46 52 51

TP 6 1 3 1 1 4

FN 7 19 10 9 6 3

FP 40 15 11 10 7 8

TN

SE 0.684 0.969 0.933 0.979 0.981 0.927

0.851 0.441 0.524 0.526 0.538 0.727

SP

test set 0.650 0.620 0.808 0.836 0.897 0.944

PRE1 0.870 0.938 0.786 0.909 0.875 0.667

PRE2

GA 0.803 0.697 0.803 0.848 0.894 0.894

1005

BC-25 BC-20 BC-26 BC-27 BC-28 BC-29

BC-25 BC-20 BC-26 BC-27 BC-28 BC-29

1 5 10 20 30 40

model

1 5 10 20 30 40

thresholds (μM)

LPFC_8 FCFP_8 FCFP_6 FCFP_8 ECFP_8 LCFP_4

LPFC_8 FCFP_8 FCFP_6 FCFP_8 ECFP_8 LCFP_4

fingerprints

117 180 247 281 281 292

117 180 247 281 281 292

TP

17 44 32 40 63 72

17 44 32 40 63 72

FN

36 76 125 73 63 58

36 76 125 73 63 58

FP

550 419 315 320 306 301

550 419 315 320 306 301

TN

0.873 0.804 0.885 0.875 0.817 0.802

0.873 0.804 0.885 0.875 0.817 0.802

SE

0.939 0.846 0.716 0.814 0.829 0.838

0.939 0.846 0.716 0.814 0.829 0.838

SP

training set

0.765 0.703 0.664 0.794 0.817 0.834

0.765 0.703 0.664 0.794 0.817 0.834

PRE1

PRE2

GA

Percentage of 20% 0.970 0.926 0.905 0.833 0.908 0.782 0.889 0.842 0.829 0.823 0.807 0.820 Percentage of 30% 0.970 0.926 0.905 0.833 0.908 0.782 0.889 0.842 0.829 0.823 0.807 0.820

0.772 0.629 0.586 0.686 0.646 0.641

0.772 0.629 0.586 0.686 0.646 0.641

C

42 149 254 331 345 178

32 87 135 124 115 92

TP

518 411 306 229 215 382

218 163 115 126 135 158

FN

43 169 337 288 289 221

53 231 456 393 389 307

FP

1350 1224 1056 1105 1104 1172

1650 1472 1247 1310 1314 1396

TN

0.075 0.266 0.454 0.591 0.616 0.318

0.128 0.348 0.540 0.496 0.460 0.368

SE

0.969 0.879 0.758 0.793 0.793 0.841

0.969 0.864 0.732 0.769 0.772 0.820

SP

test set

0.494 0.469 0.430 0.535 0.544 0.446

0.376 0.274 0.228 0.240 0.228 0.231

PRE1

0.723 0.749 0.775 0.828 0.837 0.754

0.883 0.900 0.916 0.912 0.907 0.898

PRE2

0.713 0.703 0.671 0.735 0.742 0.691

0.861 0.798 0.708 0.734 0.732 0.762

GA

Table 5. Performance of the Best Bayesian Classifiers for the PubChem Bioassay Test Set of hERG Blockage Using the Percentage of 20% or 30% as the Threshold

model

BC-19 BC-20 BC-21 BC-22 BC-23 BC-24

thresholds (μM)

training set

Table 4. Performance of the Bayesian Classifiers for the WOMBAT-PK Test Set C

0.098 0.177 0.208 0.374 0.395 0.179

0.159 0.192 0.198 0.201 0.177 0.156

C

0.527 0.478 0.521 0.614 0.633 0.632

Molecular Pharmaceutics Article

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

Figure 7. (a) The 15 good and (b) 15 bad fragments for hERG identified by the Bayesian classifier based on molecular properties and the ECFP_8 fingerprint set.

them. In the PubChem bioassay, for the positive control, the value obtained using 10 μM terfenadine is defined as 100% hERG blockade. To find an appropriate threshold to use for validating the PubChem test set, we tested the prediction performance of the best six classifiers shown in Table 1, and the results are summarized in Table 5. For the same PubChem data, Su and co-workers reported a rather low overall accuracy of 52% when using the threshold of 40 μM.17 Using the threshold of 40 μM, the predictions given by the Bayesian classifier developed here are much better than Su’s results, but the predictions are still not satisfactory, indicated by a relatively low GA (GA = 70.0%). Besides the threshold of 40 μM, Su et al. also tested the threshold of 10 μM, and they found that the classification model based on the threshold of 10 μM gave an overall accuracy of 65% for the PubChem data set. As shown in Table 5, the Bayesian classifier (GA = 71%) still outperforms

by the GA values larger than 80%. Compared with Li’s results, our predictions for the 66 tested molecules are much better. Based on the threshold of 40 μM, Li’s model only achieved the best overall accuracy of 72%, 85% (47/55) for the blockers and 36% (4/11) for the nonblockers. Therefore, we believe that the Bayesian classifiers based on molecular properties and structural fingerprints are highly reliable and functionally useful for the prediction of hERG toxicity across a wide range of chemical space. We also downloaded a larger test set of hERG bioassay data from PubChem. The blocking activities provided by PubChem are represented by the percentage of hERG blockade but not an actual IC50 value. The molecules in the PubChem data set were classified as blockers or nonblockers, using a threshold value of 20% hERG blockage. The IC50 values collected in our data set cannot directly match the percentage of hERG blockade in PubChem. So we investigated the internal relationship between 1006

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

Figure 8. (a) The 10 nonblockers and (b) 8 blockers misclassified by the Bayesian classifier based on molecular properties and the ECFP_8 fingerprint set.

Su’s model when using the threshold of 10 μM. Among the six thresholds, the classifier based on the threshold of 1 μM yields the best GA value (86.1%). The PubChem test set is highly biased, the inactive number is almost seven times the active one, and its biased degree is quite similar to that of the training set when the threshold of 1 μM is used. 5. Analysis of the Important Fragments Given by Naive Bayesian Classifier. The important fragments given by the Bayesian classifier may be useful for experimental scientists when designing molecules avoid hERG binding. The 15 good and 15 bad fragments favorable and unfavorable for hERG binding ranked by the Bayesian scores are shown in Figure 7. By analyzing the fingerprints with positive contributions to hERG binding shown in Figure 7a, we observe that most fragments (14 fragments) have nitrogen atoms, and four fragments among the top five have nitrogen atoms with positive charges. Acxcording to previous studies, it was believed that a

positively charged nitrogen in general increases the likelihood of hERG binding.69 The basic nitrogen is involved in the cation−π interaction with Tyr652. Moreover, we found that in some fragments, such as fragments 1, 2, 3, and 4, the tertiary amine group was linked by a hydrophobic tail and the hydrophobic part in these fragments may form strong van der Waals or hydrophobic interactions with some residues in hERG, such as Phe656. We also found that some fragments in Figure 7a had hydrophobic or aromatic rings, such as fragments 2, 5, 10, 11, 13, and 15. It is believed that these hydrophobes may form favorable hydrophobic or π−π interaction with Phe656 and Tyr652. The top 15 fingerprints unfavorable for hERG binding are shown in Figure 7b. Analysis of these fragments indicates that the top 5 fragments have at least one oxygen atom or a carboxylic acid. As we discussed above, the positive ionizable nitrogen is important for ligand−hERG interaction. The oxygen 1007

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

naive Bayesian classifiers are efficient for the predictions of hERG blockage. Moreover, the important molecular fragments favorable or unfavorable for hERG blockage are highlighted by the Bayesian analysis, and they are very helpful for the design of new drugs avoiding unfavorable hERG blockage. The spatial arrangement of the features important for hERG binding may not be well characterized simply by the Bayesian classifier if they are not located in the same fingerprint. In the future, based on different principles, we can develop different prediction models that complement each other. Then, it is possible that the combination of two or more models based on different principles can give higher confidence for predicting hERG blockage.

atom is negatively charged and cannot form cation−π interaction with Tyr652; moreover, the fragments with a carboxylic acid group are usually hydrophilic and they are also unfavorable for the favorable ligand−hERG hydrophobic interactions. Another interesting finding is that two fragments (9 and 12) in Figure 7b have nitrogen atoms; however, the nitrogen in fragment 9 is not positively charged and that in fragment 12 is a quaternary amine, and they are also not favorable for forming cation−π interaction with Tyr652. It is obvious that not all nitrogen atoms are favorable for ligand binding. 6. Why Some Molecules Cannot Be Predicted Correctly. In the test set, 10 of the 59 nonblockers were misclassified as false positives and 8 of the 61 blockers were misclassified as false negatives by the naive Bayesian classifier based on the threshold of 30 μM (Figure 8). Among the 18 outliers, 14 have the Bayesian scores located in the uncertain zone from −10 to 10. We believe that the following reasons can account for the misclassification. First, the misclassification may be primarily caused by the intrinsic limitation of the Bayesian classifiers based on fingerprints. All of the 10 misclassified nonblockers have aromatic rings that are thought to form π−π interaction with Phe656 and Tyr652. Moreover, most of them contain one or even more important groups, e.g. amine, favorable for hERG blockage shown in Figure 7. It should be noted that the hERG binding capability of a molecule is not solely determined by the presence of the important fragments favorable for hERG binding. The spatial location and arrangement of these important fragments are also essential, but they cannot be well characterized by the 2-D fingerprints used in our study. Moreover, if a positively charged nitrogen is already involved in the interactions with Tyr652, the other positively charged nitrogen atoms cannot give positive contribution any more. For example, desvenlafaxine has seven nitrogen atoms, and this molecule may be highly overestimated. Second, IC50 is an important but not a direct measure of hERG inhibition. The threshold of IC50 to define the blocker and nonblocker classes is also arbitrary. Six misclassified blockers do not have precise IC50 values, suggesting that they may be weak hERG blockers and be misclassified. Third, the experimental data is collected from a variety of sources, and the diversity of the experimental data will increase the data uncertainty and prediction complexity. Moreover, it is quite possible that the data collected from the literature may have mistakes in experimental evaluation.



ASSOCIATED CONTENT

S Supporting Information *

The distributions of molecular weight and IC50 (Figure S1). This material is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*Soochow University, Institute of Functional Nano & Soft Materials (FUNSOM), 199 Ren'ai Road, Suzhou, Jiangsu 215123, China. E-mail: [email protected] or tingjunhou@ hotmail.com. Phone: +86-512-65882039. Fax: +86-512-65882846. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This study was supported by the National Science Foundation of China (20973121 to T.H.), the National Basic Research Program of China (973 program, 2012CB932600 to T.H.), the NIH (R21GM097617 to J.W.) and the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD).



REFERENCES

(1) Sanguinetti, M. C.; Tristani-Firouzi, M. hERG potassium channels and cardiac arrhythmia. Nature 2006, 440, 463−469. (2) Recanatini, M.; Poluzzi, E.; Masetti, M.; Cavalli, A.; De Ponti, F. QT prolongation through hERG K+ channel blockade: Current knowledge and strategies for the early prediction during drug development. Med. Res. Rev. 2005, 25, 133−166. (3) Witchel, H. J. The hERG potassium channel as a therapeutic target. Expert Opin. Ther. Targets 2007, 11, 321−336. (4) Keating, M. T.; Sanguinetti, M. C. Molecular genetic insights into cardiovascular disease. Science 1996, 272, 681−685. (5) Hancox, J. C.; Mitcheson, J. S. Combined hERG channel inhibition and disruption of trafficking in drug-induced long QT syndrome by fluoxetine: a case-study in cardiac safety pharmacology. Br. J. Pharmacol. 2006, 149, 457−459. (6) Polak, S.; Wisniowska, B.; Brandys, J. Collation, assessment and analysis of literature in vitro data on hERG receptor blocking potency for subsequent modeling of drugs’ cardiotoxic properties. J. Appl. Toxicol. 2009, 29, 183−206. (7) Farid, R.; Day, T.; Friesner, R. A.; Pearlstein, R. A. New insights about HERG blockade obtained from protein modeling, potential energy mapping, and docking studies. Bioorg. Med. Chem. 2006, 14, 3160−3173. (8) Rajamani, R.; Tounge, B. A.; Li, J.; Reynolds, C. H. A two-state homology model of the hERG K+ channel: application to ligand binding. Bioorg. Med. Chem. Lett. 2005, 15, 1737−1741. (9) Ekins, S.; Crumb, W. J.; Sarazan, R. D.; Wikel, J. H.; Wrighton, S. A. Three-dimensional quantitative structure-activity relationship for



CONCLUSIONS In this study, based on an extensive hERG inhibition data set of 806 molecules, we first examined the relationships between eight important molecular properties and hERG binding capability. We found that solubility and molecular weight were more important for hERG binding than the other six molecular properties, but no single molecular property can be used to discriminate between hERG blockers and nonblockers efficiently. Then, the naive Bayesian classification and the RP technique were used to establish the classifiers to distinguish blockers from nonblockers. The Bayesian classifiers perform better than the RP models. Using the threshold of 30 μM, the best Bayesian classifier with the 14 molecular properties and the ECFP fingerprint set demonstrates good predictivity, indicated by the high prediction accuracies for the training set (GA = 84.8%), the test set I (GA = 85%), the test set II (GA = 89.4%), and the PubChem test set (GA = 75.3%), indicating that the 1008

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

inhibition of human ether-a-go-go-related gene potassium channel. J. Pharmacol. Exp. Ther. 2002, 301, 427−434. (10) Cavalli, A.; Poluzzi, E.; De Ponti, F.; Recanatini, M. Toward a pharmacophore for drugs inducing the long QT syndrome: Insights from a CoMFA study of HERG K+ channel blockers. J. Med. Chem. 2002, 45, 3844−3853. (11) Pearlstein, R. A.; Vaz, R. J.; Kang, J. S.; Chen, X. L.; Preobrazhenskaya, M.; Shchekotikhin, A. E.; Korolev, A. M.; Lysenkova, L. N.; Miroshnikova, O. V.; Hendrix, J.; Rampe, D. Characterization of HERG potassium channel inhibition using CoMSiA 3D QSAR and homology modeling approaches. Bioorg. Med. Chem. Lett. 2003, 13, 1829−1835. (12) Cianchetta, G.; Li, Y.; Kang, J. S.; Rampe, D.; Fravolini, A.; Cruciani, G.; Vaz, R. J. Predictive models for hERG potassium channel blockers. Bioorg. Med. Chem. Lett. 2005, 15, 3637−3642. (13) Sun, H. M. An accurate and interpretable Bayesian classification model for prediction of hERG liability. ChemMedChem 2006, 1, 315− 322. (14) Thai, K. M.; Ecker, G. F. Binary QSAR model for classification of hERG potassium channel blockers. Bioorg. Med. Chem. 2008, 16, 4107−4119. (15) Song, M. H.; Clark, M. Development and evaluation of an in silico model for hERG binding. J. Chem. Inf. Model. 2006, 46, 392− 400. (16) Li, Q. Y.; Jorgensen, F. S.; Oprea, T.; Brunak, S.; Taboureau, O. hERG classification model based on a combination of support vector machine method and GRIND descriptors. Mol. Pharmaceutics 2008, 5, 117−127. (17) Su, B. H.; Shen, M. Y.; Esposito, E. X.; Hopfinger, A. J.; Tseng, Y. J. In Silico Binary Classification QSAR Models Based on 4DFingerprints and MOE Descriptors for Prediction of hERG Blockage. J. Chem. Inf. Model. 2010, 50, 1304−1318. (18) Kongsamut, S.; Kang, J. S.; Chen, X. L.; Roehr, J.; Rampe, D. A comparison of the receptor binding and HERG channel affinities for a series of antipsychotic drugs. Eur. J. Pharmacol. 2002, 450, 37−41. (19) Tobita, M.; Nishikawa, T.; Nagashima, R. A discriminant model constructed by the support vector machine method for HERG potassium channel inhibitors. Bioorg. Med. Chem. Lett. 2005, 15, 2886−2890. (20) Du, L. P.; Li, M. Y.; You, Q. D.; Xia, L. A novel structure-based virtual screening model for the hERG channel blockers. Biochem. Biophys. Res. Commun. 2007, 355, 889−894. (21) Roche, O.; Trube, G.; Zuegge, J.; Pflimlin, P.; Alanine, A.; Schneider, G. A virtual screening method for prediction of the hERG potassium channel liability of compound libraries. ChemBioChem 2002, 3, 455−459. (22) Gintant, G. An evaluation of hERG current assay performance: Translating preclinical safety studies to clinical QT prolongation. Pharmacol. Ther. 2011, 129, 109−119. (23) De Bruin, M. L.; Pettersson, M.; Meyboom, R. H. B.; Hoes, A. W.; Leufkens, H. G. M. Anti-HERG activity and the risk of drug-induced arrhythmias and sudden death. Eur. Heart J. 2005, 26, 590−597. (24) Gavaghan, C. L.; Arnby, C. H.; Blomberg, N.; Strandlund, G.; Boyer, S. Development, interpretation and temporal evaluation of a global QSAR of hERG electrophysiology screening data. J. Comput.Aided Mol. Des. 2007, 21, 189−206. (25) Mittelstadt., S. W.; Hemenway., C. L.; Craig., M. P.; Hove., J. R. Evaluation of zebrafish embryos as a model for assessing inhibition of hERG. J. Pharmacol. Toxicol. Methods 2008, 57, 100−105. (26) Obrezanova, O.; Csanyi, G.; Gola, J. M. R.; Segall, M. D. Gaussian processes: A method for automatic QSAR Modeling of ADME properties. J. Chem. Inf. Model. 2007, 47, 1847−1857. (27) Obrezanova, O.; Segall, M. D. Gaussian Processes for Classification: QSAR Modeling of ADMET and Target Activity. J. Chem. Inf. Model. 2010, 50, 1053−1061. (28) Fraley, M. E.; Garbaccio, R. M.; Arrington, K. L.; Hoffman, W. F.; Tasber, E. S.; Coleman, P. J.; Buser, C. A.; Walsh, E. S.; Hamilton, K.; Fernandes, C.; Schaber, M. D.; Lobell, R. B.; Tao, W. K.; South, V. J.; Yan, Y. W.; Kuo, L. C.; Prueksaritanont, T.; Shu, C.; Torrent, M.;

Heimbrook, D. C.; Kohl, N. E.; Huber, H. E.; Hartman, G. D. Kinesin spindle protein (KSP) inhibitors. Part 2: The design, synthesis, and characterization of 2,4-diaryl-2,5-dihydropyrrole inhibitors of the mitotic kinesin KSP. Bioorg. Med. Chem. Lett. 2006, 16, 1775−1779. (29) Garbaccio, R. M.; Fraley, M. E.; Tasber, E. S.; Olson, C. M.; Hoffman, W. F.; Arrington, K. L.; Torrent, M.; Buser, C. A.; Walsh, E. S.; Hamilton, K.; Schaber, M. D.; Fernandes, C.; Lobell, R. B.; Tao, W. K.; South, V. J.; Yan, Y. W.; Kuo, L. C.; Prueksaritanont, T.; Slaughter, D. E.; Shu, C.; Heimbrook, D. C.; Kohl, N. E.; Huber, H. E.; Hartman, G. D. Kinesin spindle protein (KSP) inhibitors. Part 3: Synthesis and evaluation of phenolic 2,4-diaryl-2,5-dihydropyrroles with reduced hERG binding and employment of a phosphate prodrug strategy for aqueous solubility. Bioorg. Med. Chem. Lett. 2006, 16, 1780−1783. (30) McBriar, M. D.; Guzik, H.; Shapiro, S.; Xu, R.; Paruchova, J.; Clader, J. W.; O’Neill, K.; Hawes, B.; Sorota, S.; Margulis, M.; Tucker, K.; Weston, D. J.; Cox, K. Bicyclo[3.1.0]hexyl urea melanin concentrating hormone (MCH) receptor-1 antagonists: Impacting hERG liability via aryl modifications. Bioorg. Med. Chem. Lett. 2006, 16, 4262−4265. (31) Alberati, D.; Hainzl, D.; Jolidon, S.; Krafft, E. A.; Kurt, A.; Maier, A.; Pinard, E.; Thomas, A. W.; Zimmerli, D. Discovery of 4substituted-8-(2-hydroxy-2-phenyl-cyclohexyl)-2,8-diaza-spiro[4.5]decan-1-one as a novel class of highly selective GlyT1 inhibitors with improved metabolic stability. Bioorg. Med. Chem. Lett. 2006, 16, 4311− 4315. (32) Price, D. A.; Armour, D.; de Groot, M.; Leishman, D.; Napier, C.; Perros, M.; Stammen, B. L.; Wood, A. Overcoming HERG affinity in the discovery of the CCR5 antagonist maraviroc. Bioorg. Med. Chem. Lett. 2006, 16, 4633−4637. (33) Zhu, B. Y.; Jia, Z. J.; Zhang, P. L.; Su, T.; Huang, W. R.; Goldman, E.; Tumas, D.; Kadambi, V.; Eddy, P.; Sinha, U.; Scarborough, R. M.; Song, Y. H. Inhibitory effect of carboxylic acid group on hERG binding. Bioorg. Med. Chem. Lett. 2006, 16, 5507− 5512. (34) Fluxe, A.; Wu, S. D.; Sheffer, J. B.; Janusz, J. M.; Murawsky, M.; Fadayel, G. M.; Fang, B.; Hare, M.; Djandjighian, L. Discovery and synthesis of tetrahydroindolone-derived carbamates as Kv1.5 blockers. Bioorg. Med. Chem. Lett. 2006, 16, 5855−5858. (35) Wu, S. D.; Fluxe, A.; Janusz, J. M.; Sheffer, J. B.; Browning, G.; Blass, B.; Cobum, K.; Hedges, R.; Murawsky, M.; Fang, B.; Fadayel, G. M.; Hare, M.; Djandjighian, L. Discovery and synthesis of tetrahydroindolone derived semicarbazones as selective Kv1.5 blockers. Bioorg. Med. Chem. Lett. 2006, 16, 5859−5863. (36) Lynch, J. K.; Freeman, J. C.; Judd, A. S.; Iyengar, R.; Mulhern, M.; Zhao, G.; Napier, J. J.; Wodka, D.; Brodjian, S.; Dayton, B. D.; Falls, D.; Ogiela, C.; Reilly, R. M.; Campbell, T. J.; Polakowski, J. S.; Hernandez, L.; Marsh, K. C.; Shapiro, R.; Knourek-Segel, V.; Droz, B.; Bush, E.; Brune, M.; Preusser, L. C.; Fryer, R. M.; Reinhart, G. A.; Houseman, K.; Diaz, G.; Mikhail, A.; Limberis, J. T.; Sham, H. L.; Collins, C. A.; Kym, P. R. Optimization of chromone-2-carboxamide melanin concentrating hormone receptor 1 antagonists: Assessment of potency, efficacy, and cardiovascular safety. J. Med. Chem. 2006, 49, 6569−6584. (37) Nisius, B.; Goller, A. H.; Bajorath, J. Combining Cluster Analysis, Feature Selection and Multiple Support Vector Machine Models for the Identification of Human Ether-a-go-go Related Gene Channel Blocking Compounds. Chem. Biol. Drug Des. 2009, 73, 17−25. (38) Aronov, A. M.; Goldman, B. B. A model for identifying HERG K+ channel blockers. Bioorg. Med. Chem. 2004, 12, 2307−2315. (39) Halgren, T. A. Merck molecular force field 0.1. Basis, form, scope, parameterization, and performance of MMFF94. J. Comput. Chem. 1996, 17, 490−519. (40) Hou, T.; Wang, J. Structure - ADME relationship: still a long way to go? Expert Opin. Drug Metab. Toxicol. 2008, 4, 759−770. (41) Hou, T. J.; Li, Y. Y.; Zhang, W.; Wang, J. M. Recent Developments of In Silico Predictions of Intestinal Absorption and Oral Bioavailability. Comb. Chem. High Throughput Screening 2009, 12, 497−506. 1009

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010

Molecular Pharmaceutics

Article

(64) Norinder, U.; Haeberlein, M. Computational approaches to the prediction of the blood-brain distribution. Adv. Drug Delivery Rev. 2002, 54, 291−313. (65) Norinder, U.; Bergstrom, C. A. S. Prediction of ADMET properties. ChemMedChem 2006, 1, 920−937. (66) van de Waterbeemd, H.; Gifford, E. ADMET in silico modelling: Towards prediction paradise? Nat. Rev. Drug Discovery 2003, 2, 192−204. (67) Zhao, Y. H.; Le, J.; Abraham, M. H.; Hersey, A.; Eddershaw, P. J.; Luscombe, C. N.; Boutina, D.; Beck, G.; Sherborne, B.; Cooper, I.; Platts, J. A. Evaluation of human intestinal absorption data and subsequent derivation of a quantitative structure-activity relationship (QSAR) with the Abraham descriptors. J. Pharm. Sci. 2001, 90, 749− 784. (68) Hou, T. J.; Wang, J. M.; Zhang, W.; Xu, X. J. ADME evaluation in drug discovery. 6. Can oral bioavailability in humans be effectively predicted by simple molecular property-based rules? J. Chem. Inf. Model. 2007, 47, 460−463. (69) Aronov, M. M. Predictive in silico modeling for hERG channel blockers. Drug Discovery Today 2005, 10, 149−155.

(42) Ghose, A. K.; Viswanadhan, V. N.; Wendoloski, J. J. Prediction of hydrophobic (lipophilic) properties of small organic molecules using fragmental methods: An analysis of ALOGP and CLOGP methods. J. Phys. Chem. A 1998, 102, 3762−3772. (43) Csizmadia, F.; TsantiliKakoulidou, A.; Panderi, I.; Darvas, F. Prediction of distribution coefficient from structure 0.1. Estimation method. J. Pharm. Sci. 1997, 86, 865−871. (44) Tetko, I. V.; Tanchuk, V. Y.; Kasheva, T. N.; Villa, A. E. P. Estimation of aqueous solubility of chemical compounds using E-state indices. J. Chem. Inf. Comput. Sci. 2001, 41, 1488−1493. (45) Xue, Y.; Li, Z. R.; Yap, C. W.; Sun, L. Z.; Chen, X.; Chen, Y. Z. Effect of molecular descriptor feature selection in support vector machine classification of pharmacokinetic and toxicological properties of chemical agents. J. Chem. Inf. Comput. Sci. 2004, 44, 1630−1638. (46) Discovery Studio 2.5 Guide, Accelrys Inc., San Diego, http:// www.accelrys.com. 2009. (47) Rogers, D.; Hahn, M. Extended-Connectivity Fingerprints. J. Chem. Inf. Model. 2010, 50, 742−754. (48) Rogers, D.; Brown, R. D.; Hahn, M. Using extendedconnectivity fingerprints with Laplacian-modified Bayesian analysis in high-throughput screening follow-up. J. Biomol. Screening 2005, 10, 682−686. (49) Yoshida, F.; Topliss, J. G. QSAR model for drug human oral bioavailability. J. Med. Chem. 2000, 43, 2575−2585. (50) Hou, T. J.; Xia, K.; Zhang, W.; Xu, X. J. ADME evaluation in drug discovery. 4. Prediction of aqueous solubility based on atom contribution approach. J. Chem. Inf. Comput. Sci. 2004, 44, 266−275. (51) Hou, T. J.; Xu, X. J. ADME evaluation in drug discovery. 2. Prediction of partition coefficient by atom-additive approach based on atom-weighted solvent accessible surface areas (vol 43, pg 1058, 2003). J. Chem. Inf. Comput. Sci. 2004, 44, 1516−1516. (52) Wang, J. M.; Hou, T. J.; Xu, X. J. Aqueous Solubility Prediction Based on Weighted Atom Type Counts and Solvent Accessible Surface Areas. J. Chem. Inf. Model. 2009, 49, 571−581. (53) Chen, L.; Li, Y. Y.; Zhao, Q.; Peng, H.; Hou, T. J. ADME Evaluation in Drug Discovery. 10. Predictions of P-Glycoprotein Inhibitors Using Recursive Partitioning and Naive Bayesian Classification Techniques. Mol. Pharmaceutics 2011, 8, 889−900. (54) Hou, T. J.; Wang, J. M.; Li, Y. Y. ADME evaluation in drug discovery. 8. The prediction of human intestinal absorption by a support vector machine. J. Chem. Inf. Model. 2007, 47, 2408−2415. (55) Beresford, A. P.; Selick, H. E.; Tarbit, M. H. The emerging importance of predictive ADME simulation in drug discovery. Drug Discovery Today 2002, 7, 109−116. (56) Cruciani, C.; Crivori, P.; Carrupt, P. A.; Testa, B. Molecular fields in quantitative structure-permeation relationships: the VolSurf approach. J. Mol. Struct.: THEOCHEM 2000, 503, 17−30. (57) Gola, J.; Obrezanova, O.; Champness, E.; Segall, M. ADMET property prediction: The state of the art and current challenges. QSAR Combi. Sci. 2006, 25, 1172−1180. (58) Hou, T. J.; Wang, J. M.; Zhang, W.; Wang, W.; Xu, X. Recent advances in computational prediction of drug absorption and permeability in drug discovery. Curr. Med. Chem. 2006, 13, 2653− 2667. (59) Hou, T. J.; Wang, J. M.; Zhang, W.; Xu, X. J. ADME evaluation in drug discovery. 7. Prediction of oral absorption by correlation and classification. J. Chem. Inf. Model. 2007, 47, 208−218. (60) Huuskonen, J. Estimation of water solubility from atom-type electrotopological state indices. Environ. Toxicol. Chem. 2001, 20, 491− 497. (61) Dearden, J. C. In silico prediction of ADMET properties How far have we come? Expert Opin. Drug Metab. Toxicol. 2007, 3, 635− 639. (62) Jolivette, L. J.; Ekins, S. Methods for predicting human drug metabolism. Adv. Clin. Chem. 2007, 43, 131−176. (63) Liu, R. F.; Sun, H. M.; So, S. S. Development of quantitative structure-property relationship models for early ADME evaluation in drug discovery. 2. Blood-brain barrier penetration. J. Chem. Inf. Comput. Sci. 2001, 41, 1623−1632. 1010

dx.doi.org/10.1021/mp300023x | Mol. Pharmaceutics 2012, 9, 996−1010