The organic anion transporting polypeptides (OATPs) are key membrane bound transporters expressed in many organs including intestine, liver, lung, choroid plexus, blood-brain barrier, and other organs (Tamai et al., 2000). This transporter family mediates sodium-independent transport of a diverse array of molecules that are mostly anions as well as organic cations, steroid conjugates, organic anions, and xenobiotics (Bossuyt et al., 1996; Hagenbuch and Meier, 2004) by coupling uptake with the efflux of bicarbonate (Satlin et al., 1997) or glutathione (Li et al., 1998). The OATPs share some substrate overlapping specificity with other promiscuous efflux transporters such as P-gp and multidrug resistance-associated protein 2, indicative of a degree of coordination. As a result, OATPs have been implicated in drug-drug interactions (Kim, 2003) as exemplified by the interactions between cerivastatin and cyclosporine A as well as cerivastatin, gemfibrozil, and its glucuronide metabolite (Shitara et al., 2003, 2004). Thirty-six mammalian OATPs have been identified, but only a few of these have been characterized in any detail. Currently, 11 human OATPs have been identified, and a new species-independent nomenclature system has been proposed for all of the OATPs (SLCO) (Hagenbuch and Meier, 2004). OATP1B1 (previous names OATP-C, LST-1, OATP2, and SLC21A6), consisting of 670 amino acids with 12 putative membrane spanning domains, represents the most studied human OATP to date (Meier and Stieger, 2000) and is expressed on the basolateral plasma membrane of hepatocytes. This transporter has been well characterized when transfected in Xenopus laevis oocyte and other expression systems such as CHO, HeLa, and Hek-293 cells. Many diverse molecules with a range of K_m values are substrates of this protein, including pravastatin, dehydroepiandrosterone sulfate (DHEAS), and bromosulfophthalein (BSP) (Table 1).

An equally well studied rat SLCO, namely Oatp1a1 (Oatp1, Oatp, Slc21a1) that is expressed in the liver, kidney, and choroid plexus appears to share many of the same substrates as OATP1B1 (67% identity) as well as organic cations (Table 2). The considerable substrate overlap of these rat and human transporters, although by no means identical (Meier et al., 1997), could suggest some degree of structural homology in key substrate recognition areas of the proteins that may require assessment using site-directed mutagenesis or other methods. As no high-res-olution structures are presently available for these trans-porters, it is necessary to gain structural insight from alternative methods. If sufficient binding data are generated in vitro, these can be used to build three-dimensional quantitative structure-activity relationship (3D-QSAR) models (Ekins and Swaan, 2004). A recent report has briefly discussed 3D-QSAR studies of Oatp/OATP substrates which apparently produced a pharmacophore containing two hydrogen bond acceptors, one hydrogen bond donor, and two hydrophobic regions, although no model was shown or further details provided (Hagenbuch and Meier, 2004). It is likely that computational modeling techniques can provide structural insights into the molecular features of substrates that are common to these OATP transporters; alternatively, they may indicate those structural features that differentiates one OATP from the other. In the present study, we generated pharmacophore models for rat Oatp1a1 and human OATP1B1 with published K_m values from different experimental systems. Meta-analysis of the individual models allowed us to provide an understanding of the key molecular features for substrate-OATP transporter interactions.

Materials and Methods

Literature Search. A complete literature search was performed up to September 2004 to retrieve substrates for each OATP evaluated in this study. Various databases were consulted for this task including PubMed (http://www.ncbi.nlm.nih.gov/entrez/query.fcgi), the TP-search transporter database (http://www.ilab.rise.waseda.ac.jp/tp-search/), and ISI Web of Knowledge (http://isiwebofknowledge.com/). Details from the publications were noted, such as the cell system used for transfection, the K_m value, and the species. Considerable care was taken to assign molecules to the correct OATP/Oatp. Comparisons between substrates with K_m data shared between Oatp1a1 and OATP1B1 were used to derive a correlation (Table 3). In cases where multiple values for the same substrate could be found for each cell line, the mean value was then computed.

OATP Pharmacophore Development with Catalyst. The computational molecular modeling studies were carried out using a Silicon Graphics (Palo Alto, CA) Octane workstation. Briefly, models were constructed using Catalyst version 4.9 (Accelrys, San Diego, CA) to generate hypotheses. Molecules for each OATP were initially sketched in ChemDraw version 7.0.1 (CambridgeSoft, Cambridge, MA), saved in mol file format, and then imported into Catalyst. The 3-D molecular structures were produced using up to 255 conformers with the fast conformer generation method, allowing a maximum energy difference of 20 kcal/mol. Ten hypotheses were generated using these conformers for each of the molecules in the training sets and the K_m values after selection of the following features: hydrophobic, hydrogen bond acceptor, hydrogen bond donor, and the positive and negative ionizable features for the substrates in the algorithm within Catalyst. After assessing all 10 hypotheses generated, the lowest energy cost hypothesis was considered the best as this possessed features representative of all the hypotheses and had the lowest total cost.

The quality of the structure-activity correlation between the estimated and observed activity values was estimated by means of an r value. Statistical significance of the retrieved hypothesis was verified by permuting (randomizing) the response variable ten times, i.e., the activities and structures of the training set compounds were mixed ten times (so that each value was no longer assigned to the original molecule) and the Catalyst hypothesis generation procedure was repeated. The total energy cost of the generated pharmacophores can be calculated from the deviation between the estimated activity and the observed activity, combined with the complexity of the hypothesis (i.e., the number of pharmacophore features). A null hypothesis can also be calculated which presumes that there is no relationship in the data, and the experimental activities are normally distributed about their mean. Hence, the greater the difference between the energy cost of the generated hypothesis and the energy cost of the null hypothesis, the less likely it is that the hypothesis reflects a chance correlation.

Validation of the Catalyst Models. The test sets contained molecules not included in the initial training sets as described previously. These test set molecules were fit by the fast-fit algorithm to the respective Catalyst models to predict a value as previously described for P450s (Ekins et al., 1999). Fast fit refers to the method of finding the optimum fit of the substrate to the hypothesis among all the conformers of the molecule without performing an energy minimization on the conformers of the molecule.

Merging Catalyst Models. Pharmacophore models generated for different cell types for each transporter were merged in Catalyst using a pairwise comparison. Final pharmacophore models with the interfeature distances labeled for clarity were illustrated using MOLEKEL 4.3 (Swiss Center for Scientific Computing, Manno, Switzerland) (Portmann and Luthi, 2000).

Statistical Evaluation of Test Set Predictions. Observed and predicted data were graphed and fit using KaleidaGraph (Synergy Software, Reading, PA) to generate an r (correlation) value. The nonparametric Spearman rho values were calculated as a measure of rank ordering in the test data using JMP version 5.1 (SAS Institute Inc., Cary, NC), where a value of 1 is optimal.

Results

Correlation between Substrate Data for Rat Oatp1a1 and Human OATP1B1. Overlapping substrate affinity for OATPs has been previously suggested due to the relatively close phylogenetic proximity of members of the OATP family. We have therefore made several comparisons using literature data published by many groups. Comparisons between substrates with K_m data shared between rat Oatp1a1 and human OATP1B1 (Table 3) were used to derive a correlation (log K_{m, oatp1a1} = 0.55 log K_{m, OATP1B1} + 0.68). The analysis of all published data available to date across both of these transporters suggests eight molecules in common which provides a reasonable correlation for the K_m values (r² of 0.64).

Pharmacophores for Human OATP1B1. Using the substrates described in the literature (Table 1), sufficient data were available from Hek-293 cells and X. laevis oocytes expressing OATP1B1 to generate two separate pharmacophores. The OATP1B1-oocyte data set consisted of 12 molecules (K_m range 0.0076–17 μM) in the training set. The pharmacophore generated in Catalyst contained two hydrogen bond acceptors and three hydrophobes (Fig. 1A) with an observed versus predicted correlation (r = 0.97) but minimal difference between the null and total cost values (Supplemental Table 1). Similarly, after scrambling the molecules and K_m data, insignificant difference was found in the statistics. A second pharmacophore, OATP1B1-Hek, consisted of 12 molecules (K_m range 0.1–268 μM) in the training set. The pharmacophore generated in Catalyst contained two hydrogen bond acceptors and one hydrophobe (Fig. 1B) with an observed versus predicted correlation (r = 0.91) but minimal difference between the null and total cost values. However, after scrambling the molecules and K_m data, the correlation decreased. Both pharmacophores compared favorably when merged, with good overlap of the hydrogen bond acceptor features and central hydrophobic features (Fig. 1F).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 1.

Pharmacophores generated from substrate data for human OATP1B1 expressed in oocytes (showing bilirubin mapped to features) (A), Hek cells (showing bilirubin monoglucuronide mapped to features) (B), rat Oatp1a1 expressed in oocytes (showing aldosterone mapped to features) (C), CHO cells (showing BSP mapped to features) (D), HeLa cells (showing taurohyodeoxycholate mapped to features) (E), merged OATP1B1 model using pharmacophores described in A and B (F), meta-analysis model using all cell type molecule data for human OATP1B1 (showing bilirubin mapped to features) (G), and merged Oatp1a1 model using pharmacophores described in C, D, and E (H), with showing aldosterone mapped to features (I). Pharmacophore features include hydrophobes (cyan), negative ionizable (blue), and hydrogen bond acceptors (green).

Pharmacophores for Rat Oatp1a1. Among the data available on substrates described in the literature on oocytes, CHO and HeLa cells expressing Oatp1a1 (Table 2), it was apparent that three separate pharmacophores could be generated. The Oatp1a1-oocyte data set consisted of 14 molecules (K_m range 0.015–3300 μM) in the training set (S-dini-trophenyl-glutathione was excluded as it was an outlier). The pharmacophore generated in Catalyst contained two hydrogen bond acceptors and three hydrophobes (Fig. 1C) with an observed versus predicted correlation (r = 0.92; Fig. 2) and a 30-point difference between the null and total cost values (Supplemental Table 1). Similarly, a decrease in correlation and increase in the total cost were found after scrambling the molecules and K_m data. A second pharmacophore, Oatp1a1-CHO, consisted of 12 molecules (K_m range 3–3000 μM) in the training set. The pharmacophore generated in Catalyst contained one hydrogen bond acceptor and two hydrophobes (Fig. 1D) with an observed versus predicted correlation (r = 0.90) but minimal difference between the null and total cost values. However, the correlation was found to decrease with no change in the total cost upon scrambling the molecules and K_m data. A third pharmacophore, Oatp1a1-HeLa, consisted of nine molecules (K_m range 3.1–214 μM) in the training set. The pharmacophore generated in Catalyst contained one hydrogen bond acceptor, one hydrophobe, and one negative ionizable feature (Fig. 1E) with an observed versus predicted correlation (r = 0.95) and larger total than null cost value. However, no change was found after scrambling the molecules and K_m data. The three pharmacophores compared favorably when merged with good overlap of the hydrogen bond acceptor features and central hydrophobic feature, whereas the negative ionizable feature was only present in one of three models which was generated with a narrower K_m range (Fig. 1H).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 2.

Test set predictions for the rat Oatp1a1-CHO pharmacophore. Fifteen molecules were used with K_m data generated in other cell types.

Meta-Pharmacophores for Rat Oatp1a1 and Human OATP1B1. Due to the limited size of the training sets generated for each model with individual cell types, we attempted to increase the training set scope by combining data from multiple cell types. Combining data for the rat Oatp1a1 resulted in a training set of 26 molecules (K_m range 0.015–3300 μM). The meta-pharmacophore generated in Catalyst contained two hydrogen bond acceptors and three hydrophobes (Fig. 1I) with a good observed versus predicted correlation (r = 0.90) and over a 30-point difference between the null and total cost values (Supplemental Table 1). After scrambling the molecules and K_m data, the correlation decreased and the total cost increased to close to the null value. Combining data for the human OATP1B1 resulted in a training set of 18 molecules (K_m range 0.0076–268 μM). The meta-pharmacophore generated in Catalyst contained two hydrogen bond acceptors and three hydrophobes (Fig. 1G) with an observed versus predicted correlation (r = 0.92) and nearly a 20-point difference between the null and total cost values (Supplemental Table 1). After scrambling the molecules and K_m data, the correlation decreased and the total cost increased.

Conformational Features of the Meta-Pharmacophores. The molecules BSP, cholate, DHEAS, estrone-3-sulfate, and taurocholate exist in the training sets for meta-OATP1B1 and meta-Oatp1a1. To illustrate that these molecules would map to the meta-pharmacophores, their relative positioning is shown in Fig. 3. DHEAS shows medium affinity to Oatp1a1 and low affinity to OATP1B1. This is validated in our OATP1B1 meta-pharmacophore models because the molecule fails to occupy a hydrophobic feature and hydrogen bond acceptor (Fig. 3, A and C). DHEAS fits reasonably well to all Oatp1a1 pharmacophore features as expected from its relatively high affinity to Oatp1a1 (Fig. 3B). The relative molecular dimensions of the pharmacophoric feature points is illustrated in Fig. 3, D–F. It should be noted that the OATP1B1 pharmacophore that includes bilirubin occupies twice the dimensions of the oatp1a1 and OATP1B1 (without bilirubin) models.

External Test Sets for Human OATP1B1 and Rat Oatp1a1 Pharmacophores. Additional molecules not included in the individual training sets for each cell type were used as test sets for the other respective model(s) for each transporter. The test set correlations using each single cell type OATP pharmacophore were as follows: OATP1B1-Hek (r = 0.214, n = 7; Spearman rho –0.14, p = 0.76), OATP1B1-oocyte (r = 0.658, n = 7; Spearman rho 0.57, p = 0.18), Oatp1a1-oocyte (r = 0.412, n = 13; Spearman rho 0.39, p = 0.18), Oatp1a1-HeLa (r = 0.224, n = 18; Spearman rho –0.26, p = 0.28), and Oatp1a1-CHO (r = 0.912, n = 15; Spearman rho 0.92, p = <0.0001, Fig. 2).

External Test Set Validation. We have previously shown (Pang et al., 1998) that BSP-glutathione is transported by rat Oatp1a1, whereas enalaprilat, hippuric acid, benzoic acid, and harmol sulfate were not. These molecules form an excellent external test for our Oatp1a1 pharmacophore model. Hippuric acid and benzoic acid failed to fit to the model (indicative of them not being substrates), and enalaprilat and harmol sulfate resulted in modest to high K_m values (69 and 550 μM, respectively). BSP-Glutathione, the glutathione adduct of BSP, was predicted to have a K_m value of 23 μM, which indicates it is a good substrate, further validating the Oatp1a1 meta-pharmacophore model. Recently, troglitazone, its metabolites, and similar compounds were shown to produce statistically significant inhibition of estrone-3-sulfate uptake in OATP1B1 (Nozawa et al., 2004). Indeed, when these molecules are fitted to the OATP1B1 meta model (without bilirubin), they are all predicted to have affinity for this transporter: troglitazone (predicted K_m 7.9 μM), troglitazone glucuronide-M2 (1.6 μM), troglitazone sulfate-M1 (5.6 μM), troglitazone quinone-M3 (6.7 μM), pioglitazone (14 μM), and rosiglitazone (9.9 μM). The OATP1B1 inhibitor indocyanine green (Cui et al., 2001) was also fitted to this same pharmacophore showing an overlap to the two hydrophobes and one hydrogen bond acceptor (predicted K_m 6.3 μM).

Discussion

In the absence of crystal structures for many of the membrane bound proteins involved in interactions with xenobiotics and endobiotics such as enzymes and transporters, computational approaches have been extremely useful in gaining insight into the ligand-protein interaction(s). However, the quality and consistency of datasets have been a determining factor in the overall predictive value of the QSAR models to date. It has been particularly challenging to assimilate and model data acquired across species and experimental cell systems. As a result, most QSAR studies have focused on datasets gathered from one species, cell type, and frequently, one laboratory setting. In the present study, we aimed to overcome these issues by a combination of pharmacophore building and meta-analysis.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 3.

Mapping of DHEAS to three pharmacophore models. Human OATP1B1 meta (A), rat Oatp1a1 meta (B), human OATP1B1 meta (no bilirubin) (C) model. Pharmacophore features include hydrophobes (cyan) and hydrogen bond acceptors (green). Illustration of the interatomic distances within the pharmacophore models for human OATP1B1 meta (D), rat Oatp1a1 meta (E), and human OATP1B1 meta (no bilirubin) (F). Red spheres represent hydrogen bond acceptor features, whereas green spheres represent hydrophobic features. The labeled interfeature distances are shown in Ångstroms. G–N, mapping of taurocholate (G, K), BSP (H, L), cholate (I, M), and estradiol-3-sulfate (J, N) to the oatp1a1 (G–J) and OATP1B1 (K–N) pharmacophores.

To our knowledge, the application of QSAR models for OATP or Oatp has been limited largely due to the absence of consistent datasets. A notable exception is a recent study by Yarim and coworkers on rat oatp1a5 (Yarim et al., 2005) who used comparative molecular field analysis on 18 substrates. An improved understanding of the structural requirements of the OATPs may explain the mechanisms underlying the reported drug-drug interactions as due to transporter inhibition (Kim, 2003). Because some studies have described multiple inhibitors of uptake with EC₅₀ values (Tirona et al., 2003), it may be possible to generate similar pharmacophores for inhibitors of the respective transporter. However, the difficulty in interpretation of whether these molecules are interacting with the same site or even sites responsible for transport is a disadvantage compared with modeling substrate K_m data. It is therefore important to determine the critical features of OATP substrates, in particular, the extensively studied rat Oatp1a1 and human OATP1B1.

Using the spectrum of substrates collated from the literature (Table 1), we built two pharmacophores with Hek-293 cells and X. laevis oocytes expressing OATP1B1 (Fig. 1, A and B). Both pharmacophores had high correlations and contained multiple hydrogen bond acceptors and hydrophobic features, which overlapped when merged, suggesting overlap (Fig. 1F). Small differences between the null and total cost values (Supplemental Table 1, online) before and after scrambling the molecules and K_m data indicated these models may be limited. Ideally, upon scrambling one would expect the training set correlation to diminish as the affinity data are randomly assigned with molecular structure and the total cost for the hypothesis will be similar to the null hypothesis cost. This simple test can give some idea of whether the model derived originally is meaningful or likely to be similar to one generated with random data. Literature data for the Oatp1a1 substrates (Table 2) yielded three pharmacophores (Fig. 1, C–E). Interestingly, all models were subtly different but shared hydrogen bond acceptors and hydrophobes, resulting in good overlaps in the merged pharmacophore. One Oatp1a1 pharmacophore contained a negative ionizable feature that likely relates to the most active molecules in the training set. Most of these models showed decreases in the training correlation and increases in the total cost after scrambling, indicating these models are of some utility. It is important to note that relatively small literature training sets (i.e., <20 compounds) can dramatically impact the model cost statistics and their behavior upon scrambling.

Test sets were generated from the molecules not initially included in the individual cell type pharmacophores, but results were variable due to the limited nature of the test sets. By combining data from different cell types, we increased the training set scope for each transporter model both in terms of the number of molecules and the activity range. These so-called meta-pharmacophores for rat Oatp1a1 and human OATP1B1 both contained two hydrogen bond acceptors and three hydrophobes (Fig. 1, G–I), although with different 3-D coordinates (Supplemental Table 2). Both models had larger differences between null and total cost values with similar observed versus predicted correlations compared with the individual models generated with separate cell lines. Also, after scrambling the molecules and K_m data for these meta-pharmacophores, the correlations decreased, and the total cost increased for both models, indicating that they are statistically significant. It seems, therefore, that the degree of correlation between substrates for these two transporters observed earlier could, in part, be due to the possession of similar molecular features necessary for interaction with each transporter, whereas the exact 3-D coordinates differ to some extent. The key pharmacophore features for both of these transporters appear to be two hydrogen bond acceptors at either end of a large hydrophobic area. These features will correspond with hydrogen bond donors and a large planar hydrophobic recognition site on the respective transporters.

Five molecules (BSP, cholate, DHEAS, estrone-3-sulfate, and taurocholate) exist in both training sets for meta-OATP1B1 and meta-Oatp1a1. All but DHEAS have similar relative affinity to both transporters and map to both pharmacophores (Fig. 3, G–N). The fact that DHEAS shows medium affinity to Oatp1a1 and low affinity to OATP1B1 makes it a good substrate to compare both pharmacophore models. Hence, DHEAS was fitted to both pharmacophore models (Fig. 3). A hydrophobic feature and hydrogen bond acceptor was missed and explained the low affinity of DHEAS to OATP1B1 (Fig. 3A). DHEAS fits reasonably well to all Oatp1a1 pharmacophore features and correlates with its relatively high affinity to Oatp1a1 (Fig. 3B). Our models could, therefore, successfully distinguish the affinity of DHEAS.

There is some controversy in the literature as to whether unbound bilirubin is a substrate for OATP1B1, and its inclusion in the meta-pharmacophore may necessitate further evaluation. Although it was shown to be a substrate in both Hek-293 cells and oocytes (Table 1), no transport was observed in HeLa and Hek-293 cells by Wang and colleagues (Wang et al., 2003). We have therefore generated a meta-pharmacophore without bilirubin and found that this pharmacophore no longer has three hydrophobic features and has less optimal training set statistics due to the narrowing of the K_m range. The mapping of DHEAS to this modified meta-model is consistent with the initial model as it maps to both hydrophobic features and one hydrogen bond acceptor feature (Fig. 3C), suggesting that the hydrophobic feature removed after exclusion of bilirubin is not critical to the overall model. Our data are unlikely to end the controversy over whether or not bilirubin is a substrate for OATP1B1.

The approaches taken in this study indicate that relatively small amounts of literature data for drug transporters can be used to provide qualitative models of the key features involved in the ligand-protein interaction. By combining data from different experimental cell systems expressing the transporter, meta-models can be generated which comply with most criteria for robust pharmacophores. Namely, these models have a sizeable difference between the null and total cost of the hypothesis and show a decrease in the average model correlation following scrambling of the molecule structures and the K_m values. External testing validated the predictive ability of the models. These meta-pharmacophores are grossly similar to merged models of the multiple pharmacophores derived with individual datasets. This approach of merging individual models for the same transporter had previously been taken with P-gp to show the concordance of five inhibitor pharmacophores and the tight clustering of hydrophobic features at the extremities and central hydrogen bonding features (Ekins et al., 2002). The pharmacophores for OATP1B1 and Oatp1a1 (Fig. 3, D–F) are distinctly different to those previously generated for P-gp, in that the hydrophobic features are now centrally located pharmacophore features with the hydrogen bond acceptor features at the extremities.

In conclusion, we have shown for the first time that using limited data sets from different laboratories, cell types, and species can be used successfully to derive robust computational pharmacophore models describing the key features for substrate interaction with rat Oatp1a1 and human OATP1B1 transporters. Our ultimate meta-pharmacophore approach combines data from different cell types and produces models which suggest a degree of similarity consistent with the biological data generated for each transporter to date (Fig. 3, D–F). We believe this meta-pharmacophore approach can be used with other transporters for which in vitro data are scattered across multiple experimental systems and species. Ultimately, this approach may aid in predicting potential drug-transporter or drug-drug interactions occurring at the key human transport proteins and overcome some of the variability due to the unique qualities of individual cell types and experimental systems used.