Combinatorial Pharmacophore Modeling of ... - ACS Publications

Organic cation transporter 2 (OCT2) is responsible for the entry step of many drugs in renal elimination, of which the changing activity may cause unw...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/molecularpharmaceutics

Combinatorial Pharmacophore Modeling of Organic Cation Transporter 2 (OCT2) Inhibitors: Insights into Multiple Inhibitory Mechanisms Yuan Xu,† Xian Liu,† Shanshan Li,† Nannan Zhou,‡ Likun Gong,† Cheng Luo,† Xiaomin Luo,*,† Mingyue Zheng,*,† Hualiang Jiang,†,‡ and Kaixian Chen† †

State Key Laboratory of Drug Research, Shanghai Institute of Materia Medica, Chinese Academy of Sciences, 555 Zuchongzhi Road, Shanghai 201203, China ‡ State Key Laboratory of Bioreactor Engineering and Shanghai Key Laboratory of Chemical Biology, School of Pharmacy, East China University of Science and Technology, Shanghai 200237, China S Supporting Information *

ABSTRACT: Organic cation transporter 2 (OCT2) is responsible for the entry step of many drugs in renal elimination, of which the changing activity may cause unwanted drug−drug interactions (DDIs). To develop drugs with favorable safety profile and provide instruction for rational clinical drug administration, it is of great interest to investigate the multiple mechanisms of OCT2 inhibition. In this study, we designed a combinatorial scheme to screen the optimum combination of pharmacophores from a pool of hypotheses established based on 162 OCT2 inhibitors. Among them, one single pharmacophore hypothesis represents a potential binding mode that may account for one unique inhibitory mechanism, and the obtained pharmacophore combination describes the multimechanisms of OCT2 inhibition. The final model consists of four individual pharmacophores, i.e., DHPR18, APR2, PRR5 and HHR4. Given a query ligand, it is considered as an inhibitor if it matches at least one of the hypotheses, or a noninhibitor if it fails to match any of four hypotheses. Our combinatorial pharmacophore model performs reasonably well to discriminate inhibitors and noninhibitors, yielding an overall accuracy around 0.70 for a test set containing 81 OCT2 inhibitors and 218 noninhibitors. Intriguingly, we found that the number of matched hypotheses was positively correlated with inhibition rate, which coincides with the pharmacophore modeling result of P-gp substrate binding. Further analysis suggested that the hypothesis PRR5 was responsible for competitive inhibition of OCT2, and other hypotheses were important for interaction between the inhibitor and OCT2. In light of the results, a hypothetical model for inhibiting transporting mediated by OCT2 was proposed. KEYWORDS: organic cation transporter 2 (OCT2), drug clearance, drug−drug interaction (DDI), pharmacophore, multiple inhibition mechanism

1. INTRODUCTION Since P-glycoprotein (P-gp, MDR1) was first detected and cloned more than 25 years ago,1 considerable interest has focused on the impact of transporters on drug absorption, distribution, metabolism, excretion and toxicity (ADME/T).2,3 Transporters affect pharmacokinetic properties of drugs and consequently may induce drug−drug interactions (DDIs). During the last decades, increasing members of transporters have been reported, which can be divided into two major families: ATP-binding cassette (ABC) transporters and solute carrier (SLC) transporters. Except for P-gp, very limited © 2013 American Chemical Society

information has been so far obtained for other members of transporters. Renal elimination is a major route of clearance for many drugs and drug metabolites. Transporters located on the membrane of the renal proximal tubule cells facilitate the flux of drugs from blood to the tubular lumen for excretion. Organic Received: Revised: Accepted: Published: 4611

July 22, 2013 September 30, 2013 October 23, 2013 October 23, 2013 dx.doi.org/10.1021/mp400423g | Mol. Pharmaceutics 2013, 10, 4611−4619

Molecular Pharmaceutics

Article

validations. The further clustering results showed that compounds from distinct clusters differ in molecular properties and their activities, suggesting that the clusters may represent distinct inhibitory mechanisms. However, since each cluster contains both inhibitors and noninhibitors, the molecular properties of different inhibition mechanisms were not carefully analyzed. Overall, previous computational studies focused on discerning OCT2 inhibitors from noninhibitors, but very little is known about structural determinants underlying the multimechanisms of OCT2 inhibition. Pharmacophore modeling is commonly used for understanding mechanism of action for a collection of active ligands. A single pharmacophore hypothesis explains how structurally diverse ligands can bind to a common site, which represents a specific mode of action. Since OCT2 inhibitors may interact with OCT2 through different patterns, it is difficult to establish a single model explaining the multiple inhibitory mechanisms. Likewise, since a single pharmacophore represents only one mechanism and matches a portion of inhibitors, it will result in a poor recall rate (i.e., sensitivity or true positive rate, defined in eq 3) in searching OCT2 inhibitors. To address this issue, here we developed a combinatorial phamacophore (CP) model of OCT2 inhibitors that consists of multiple pharmacophore hypotheses. Given an input compound, it will be predicted as an inhibitor if it matches at least one hypothesis of a CP model, or a noninhibitor if it is unable to match any of the hypotheses. From the view of classification, a true positive prediction means an experimentally observed OCT2 inhibitor (i.e., a positive sample) is predicted as an inhibitor. This strategy requires that a single hypothesis should match as few noninhibitors as possible, which controls the false positive rate that will otherwise significantly increase after combining different pharmacophores. The approach can be outlined as follows. First, a training set containing 162 experimentally verified OCT2 inhibitors was collected to establish three-, four- and five-point pharmacophore hypotheses. Altogether 30 hypotheses that meet predefined criteria were selected as candidates, including 5 three-point, 23 four-point and 2 five-point hypotheses. Then, any 3 to 5 members of the 30 candidate hypotheses were combined without repetition to comprise a CP model. The obtained CP models were assessed by internal prediction accuracy measured by F score. The CP model with the highest F score and the fewer hypotheses was further analyzed and validated by a test set including 299 compounds. By comparing the pharmacophore hypotheses constituting the CP model, we aim to provide insights into the multimechanisms of OCT2 inhibition.

cation transporter 2 (OCT2) is the key transporter for cation influx in the renal epithelial cells. It involves in the first step of renal elimination that uptakes drug molecules from blood, across the basolateral membrane, into the proximal tubule cell.4 OCT2 belongs to SLC family, consisting of 12 transmembrane domains, a large extracellular loop between TMD1 and TMD2, as well as an intracellular loop between TMD6 and TMD7. The main function of OCT2 is detoxification, which is responsible for the transporting of many clinically important drugs such as metformin,5 oxaliplatin,6 some antiretroviral drugs,7 and toxic substances [e.g., 1-methyl-4-phenylpyridinium (MPP+)]. It was reported that the renal secretion of tetraethylammonium (TEA) in OCT1/2−/− mice was completely abolished.8 Inhibition of OCT2 will lead to higher circulating plasma concentrations of its substrates, and cause dose and concentration-dependent adverse effects, even toxication.9 Therefore, to avoid drug−drug interactions (DDIs) in clinical trials, it is urgent to determine what common structures of OCT2 inhibitors are and how these inhibitors affect the transport activity of OCT2. Comprehensive understanding of how ligands interact with OCT2 is still lacking. It has been implied that there are at least three main patterns of OCT2 inhibition:10−12 (1) Competitive inhibitors interact with the substrate-binding site of OCT2. (2) Noncompetitive inhibitors cause a short-range allosteric effect on the substrate binding site.10 (3) Noncompetitive inhibitors occlude the substrate binding site and lock the transition of the transporter from outward-facing to inward-facing.11 Currently, there are over hundreds of inhibitors for hOCT2 that have been reported,13−19 but the common structural features required to inhibit the transporting are poorly understood. Computational modeling is frequently used to elicit such structural−activity relationships (SARs). Suhre et al. reported a two-dimensional quantitative structure−activity relationship (2D-QSAR) model based on 34 OCT2 inhibitors that inhibit the TEA transporting.20 They found that the determinants in influencing the activity of transporting include hydrophobicity represented by the descriptors such as AlogP98 and CLogP, ionization described by relative negative charge, together with a steric factor Shadow YZ and relative polar surface area (PSA). Besides, a phamacophore model generated with six selective hOCT2 inhibitors suggested a positive charge feature 5.72 Å away from a hydrogen bond acceptor feature and a 92.76° angle between the two features are required for their binding to OCT2. Zolk et al. reported a computational analysis of 26 commonly used drugs for inhibition of MPP+ uptake, and observed a significant correlation between their topological polar surface area (TPSA) and activity on MPP+ inhibition.21 A two-point pharmacophore was generated with ten most potent OCT2 inhibitors, showing a pattern of an ion-pair interaction site and a hydrophobic aromatic site separated by 5.0 Å. Both of the pharmacophore models revealed some structural determinants for inhibitor interactions with OCT2. However, the data sets used to construct the models are small. It is not clear to what extent the models can be applied to explain the inhibitory mechanisms. Recently, Kido et al. screened a drug library consisting of 910 compounds, among which 244 inhibited OCT2 transporting of 4-(4-(dimethylamino)styryl)-N-methylpyridinium (ASP+).12 Computational analyses identified the molecular charge as one key property for differentiating inhibitors from noninhibitors. Discriminant structure−activity models were derived using partial least-squares (PLS), which yielded reasonably good results on 100 times of 5-fold cross

2. MATERIALS AND METHODS Data Set. Altogether 907 compounds tested for OCT2 ASP+ transporting were collected from the literature.12 Among them, two compounds with invalid (PubChem ID: Quinolone) or duplicated (PubChem ID: 5355) records were eliminated. Additionally, a compound (PubChem ID: 5538) possessing fewer than 3 pharmacophore features was excluded as it failed in generating common pharmacophores. Altogether, after removing these three compounds, 904 ligands including 243 inhibitors were used in this study. These compound structures with their inhibition rates were then randomly divided into 2/3 for training and 1/3 for test data. The training set contains 162 inhibitors among 605 compounds; the test set contains 81 inhibitors among 299 compounds. 4612

dx.doi.org/10.1021/mp400423g | Mol. Pharmaceutics 2013, 10, 4611−4619

Molecular Pharmaceutics

Article

Figure 1. Work flow for developing combinatorial pharmacophore model. Hypo: Hypothesis. Variant: Any one of a series of pharmacopores that possess same pharmacophore features but have different spatial arrangement of these features.

Pharmacophore Development. The Phase (Schrödinger, 9.0.211, LLC, New York, NY, 2009)22 was used to prepare ligand structures and generate pharmacophore hypotheses with default parameters, unless otherwise noted. Four steps were followed to develop pharmacophore hypotheses with inhibitors in the training set. (1) Prepare ligands. First, the chemical structures downloaded from the PubChem Web site were prepared with the LigPrep module to produce low-energy three-dimensional (3D) structures with specified chiralities. During the process, the tautomerization and ionization states were adjusted at a pH of 7.4, and a maximum 1000 of conformations per structure were generated. (2) Create sites. Phase module was used to build hypotheses with six pharmacophore features including the following: hydrogen bond acceptor (A), donor (D), hydrophobic (H), negative charge (N), positive charge (P), and aromatic ring (R). (3) Find common pharmacophores. Pharmacophores that contain identical sets of features with very similar spatial arrangements were grouped together. A tree-based search was performed to find common pharmacophores with 3 to 5 sites, and each identified common pharmacophore hypothesis should be matched by at least 10 active ligands. (4) Score hypothesis. Each hypothesis generated in the previous step was scored and ranked by default with the “survival score”. The hypothesis with RMSD below 1.2 Å and vector score above 0.5 can pass the score process, and the top 10% ranked hypotheses were kept for further study. Combinatorial Pharmacophore Modeling. Two rounds of filtering were performed for the generated single pharmacophore hypotheses prior to the CP modeling. The first filtering is selectivity score, which is an empirical estimate of the rarity of a hypothesis. Selectivity score measures what fraction of molecules is likely to match the hypothesis, regardless of their activity toward the receptor. Higher selectivity means that the hypothesis is more likely to be unique to the active-set ligands, which is useful to control the false positive rate. Here, the hypotheses with “selectivity score” greater than or equal to 1.0 were chosen to search for matches in the training set. The second filtering is the number of

noninhibitors matching the hypothesis. For 4-point and 5-point pharmacophore hypotheses, the threshold for the number of matching noninhibitors was 50; for 3-point hypotheses, the threshold was set to 100. An ideal hypothesis will match more inhibitors but fewer noninhibitors. In the end, one hypothesis was selected for each series of hypotheses possessing the same set of features but different spatial arrangements. A combinatorial strategy was used to yield a complete enumeration of all possible combinations of N (N = 3 to 5) different pharmacophore hypotheses; each is hereafter called a CP model. All the generated CP models were then screened to select the one showing the best internal prediction performance on the training set. A CP model will predict a ligand as an inhibitor if it can match any of the hypotheses in the model. Our data set contains many more negative samples than positive samples. When learning from the imbalanced data, the performance measure used for model selection plays a vital role. Here, the ability for prediction of these CP models was evaluated with F score, which is the harmonic mean of precision (PR) and recall (RE). These parameters are defined below: F=2×

PR·RE PR + RE

(1)

PR =

TP TP + FP

(2)

RE =

TP TP + FN

(3)

where TP is true positive, FP is false positive, and FN is false negative prediction. In this study, PR means the fraction of positive predictions that are “true” (experimentally verified OCT2 inhibitors), and RE means the fraction of the “true” inhibitors that can be recognized (predicted as positive). Both PR and RE are therefore based on an understanding and measure of a model’s ability to identify true OCT2 inhibitors. We use F score here to achieve a balance point of PR and RE, and a higher F score means a better performance on discriminating inhibitors and noninhibitors based on an overall 4613

dx.doi.org/10.1021/mp400423g | Mol. Pharmaceutics 2013, 10, 4611−4619

Molecular Pharmaceutics

Article

consideration. In the end, the CP model which yields the highest F score and with fewer numbers of component hypotheses was kept for further analysis.

3. RESULTS Figure 1 shows the flowchart to develop combinatorial pharmacophore model. As depicted, altogether 30 individual pharmacophore hypotheses were retained after two rounds of filtering, including 2 hypotheses with 5-point, 23 hypotheses with 4-point and 5 hypotheses with 3-point features. Details about these selected hypotheses including their scores and ligand matching results can be found in Supporting Information Table S1. Altogether 173,971 CP models were generated by exhaustively enumerating all combinations of 3 to 5 hypotheses obtained above. The top ranked candidate CP models are listed in Supporting Information Table S2. One model consisting of four hypotheses, i.e., DHPR18, APR2, PRR5 and HHR4, was selected (Table 1), which showed balanced performance on

Figure 2. Four pharmacophore hypotheses for the combinatorial pharmacophore model: (A) DHPR18 with norfloxacin; (B) APR2 with orphenadrine; (C) PRR5 with pyrilamine; (D) HHR4 with mycophenolate mofetil.

Table 1. Prediction Performance of Single Pharmacophore Hypotheses and Combinatorial Pharmacophore Model on Both Training and Test Sets training set

Table 2. The Site Distances of Single Pharmacophore Hypotheses

test set

hypothesis

SEa,d

SPb,d

ACCc,d

SE

SP

ACC

DHPR18 APR2 PRR5 HHR4 CP model

0.07 0.37 0.34 0.25 0.67

0.94 0.84 0.91 0.78 0.66

0.51 0.60 0.63 0.51 0.66

0.09 0.28 0.31 0.30 0.60

0.97 0.83 0.94 0.85 0.72

0.53 0.56 0.62 0.57 0.69

a SE = RE = TP/(TP + FN). bSP = TN/(TN + FP). cACC = (TP + TN)/(TP + TN + FP + FN). dTP is true positive, FN is false negative, TN is true negative and FP is false positive.

hypothesis

site1

site2

distance (Å)

DHPR18

D2 D2 D2 H4 H4 P7 A1 A1 P3 P4 P4 R5 H12 H12 H11

H4 P7 R9 P7 R9 R9 P3 R4 R4 R5 R6 R6 H11 R15 R15

6.21 1.01 6.21 5.30 2.77 5.45 3.09 3.69 5.30 5.77 5.63 5.43 5.78 2.95 3.89

APR2

discriminating inhibitors and noninhibitors. For the training set, SE is 0.67 and SP is 0.66. For the test set, SE and SP are 0.60 and 0.72, respectively. As listed in Table 1, compared with any of the constituent hypotheses, the CP model showed a remarkably improved recall rate, suggesting that the proposed scheme is appropriate for investigating the multimechanism system. Figure 2 displays pharmacophore hypotheses with their reference ligands. The site spatial measurements of every hypothesis are listed in Table 2 and Table 3. Of note is that an aromatic ring vector is present in all four hypotheses, which highlighted the role of hydrophobic aromatic interaction or π−π interaction in ligand-OCT2 recognition. Positive charge is another distinguished feature for OCT2 inhibitors, which was found in three hypotheses. Models developed earlier revealed some important features for inhibitor binding. As previously reported for OCT1 pharmacophore models23,24 and other OCT2 computational models based on ligands inhibiting different substrates,20,21 hydrophobic and positive charge are two main properties that contribute to their inhibition. These findings implied that these two features are required for getting close to or into the pocket, while other features in different hypotheses represented the diversity of inhibitors that may act via multiple mechanisms. The pharmacophore model developed for MPP+ transportation inhibiting consisted of a positive interaction site and a hydrophobic aromatic interaction site with 5.0 Å.21 This mode resembles the positive charge (P) and aromatic ring (R)

PRR5

HHR4

features in our DHPR18, APR2, PRR5 pharmacophore hypotheses (5.45 Å, 5.30 Å, 5.63 Å respectively, as shown in Table 3). Another pharmacophore model was generated with inhibitors tested using TEA as substrate probe, which consisted of a positive charge feature and a hydrogen bond donor at a distance of 5.89 Å, and an angle of 129.97° between the positive charge and the hydrogen bond donor vector.20 These structural patterns were not observed in any hypotheses of our CP model, implying that the pharmacophore model of OCT2 inhibitors depends on the substrate employed in measuring the activities of the inhibitors. It also suggests that the mode of action for inhibiting ASP+ transporting shares more similarities with the mode for inhibiting MPP+ transporting, while it is very different from the mode for inhibiting TEA transporting. Another observation of interest is the relationship between the activity of an inhibitor and the number of pharmacophore hypotheses that the inhibitor can match. As shown in Table 4, if all the matched ligands were further divided into three groups based on their inhibition rates, namely, strong inhibitors (inhibition rate >75%), weak inhibitors (inhibition rate 4614

dx.doi.org/10.1021/mp400423g | Mol. Pharmaceutics 2013, 10, 4611−4619

Molecular Pharmaceutics

Article

inhibitors for multiple mechanisms were collected. (1) Competitive inhibitors: tripelennamine, phenyltoloxamine and diphenidol. (2) Occluding inhibitors: carvedilol, cisapride and perospirone. (3) Allosteric inhibitors: adrenosterone, medroxyprogesterone and spironolactone.12 The matching results of these inhibitors to different pharmacophore hypotheses are shown in Table 5 (the matching results of the whole data set are provided in Supporting Information Table S3). It can be noticed that competitive inhibitors only match PRR5 and APR2, and three inhibitors match PRR5. Three occlusion inhibitors, in contrast, match a variety of hypotheses, and all of them match APR2. Particularly, DHPR18 and HHR4 were only found in the category of occlusion inhibitors, suggesting that these two hypotheses may account for the noncompetitive inhibition by occlusion. Compared with the former two mechanisms, the allosteric inhibition is significantly different because it involves the binding at an “allosteric” site and the rearrangement of the active site. The difference in inhibition mechanism is also reflected in hypothesis matching. As shown, the allosteric inhibitors match none of the hypotheses, indicating that the current CP model mainly addresses the effect on active site binding of OCT2 inhibitors. A comparison of molecular weights for inhibitors and noninhibitors (Supporting Information, Figure S1) revealed that OCT2 inhibitors possess a relatively narrower range of molecular weight. In order to understand the pharmacophore hypotheses for different mechanisms, we further analyzed molecular weight of inhibitors for different hypotheses. Here we focused on the inhibitors that only matched one pharmacophore hypothesis. Since there were only three inhibitors specifically matched to DHPR18, this hypothesis was not included in analysis. PRR5 and APR2 have similar distribution that has the largest distribution in the range of 250−300 Da (Figure 3). The main difference between them is that PRR5 is more centralized in 250−480, while APR2 is widely distributed, showing the presence around 200 and 650. In contrast, the overall distribution of HHR4 is significantly different, which is flat and more normal-like. It has the largest distribution of molecular weight in the range of 350−400, which is much larger than those of PRR5 or APR2. Previous study demonstrated that competitive inhibitors of SLC6 bind to a relatively small pocket of the big transporting chamber, which confined the competitive inhibitors as small molecules. The noncompetitive inhibitors occluding substrates, instead, interacted with a different site of the chamber that allows big molecules obstructing the substrate to enter the transportation site.11 It is consistent with a computational

Table 3. The site angles of single pharmacophore hypotheses hypothesis

site1

site2

site3

angle (deg)

DHPR18

H4 H4 P7 D2 D2 P7 D2 D2 H4 D2 D2 H4 P3 A1 A1 R5 P4 P4 H11 H12 H12

D2 D2 D2 H4 H4 H4 P7 P7 P7 R9 R9 R9 A1 P3 R4 P4 R5 R6 H12 H11 R15

P7 R9 R9 P7 R9 R9 H4 R9 R9 H4 P7 P7 R4 R4 P3 R6 R6 R5 R15 R15 H11

24.1 25.7 37.7 4.5 77.2 78.0 151.4 135.7 29.8 77.1 6.5 72.2 102.4 42.9 34.7 56.8 60.3 62.8 37.7 27.6 114.7

APR2

PRR5

HHR4

Table 4. The Relationship between Inhibition Rate and the Number of Matching Hypotheses no. of matching hypotheses inhibn rate (%)

1

2

3

4

ratio of matching multiple hypotheses

>75 25−75