In silico Study on the Interaction between P-glycoprotein and Its Inhibitors at the Drug Binding Pocket

Abstract

P-glycoprotein (P-gp) is a member of the ATP-Binding Cassette transporter superfamily and mediates transmembrane efflux of many drugs. Since it is involved in multi-drug resistance activity in various cancer cells, the development of P-gp inhibitor is one of the major concerns in anticancer therapy. Human P-gp protein has at least two "functional" drug binding sites that are called "H" site and "R" site, hence it has multi-binding-specificities. Though the amino acid residues that constitute in drug binding pockets have been proposed by previous experimental evidences, the shapes and the binding poses are not revealed clearly yet. In this study, human P-gp structure was built by homology modeling with available crystal structure of mouse P-gp as a template and docking simulations were performed with inhibitors such as verapamil, hoechst33342, and rhodamine123 to construct the interaction between human P-gp and its inhibitors. The docking simulations were performed 500 times for each inhibitor, and then the interaction frequency of the amino acids at the binding poses was analyzed. With the analysis results, we proposed highly contributing residues that constitute binding pockets of the human P-gp for the inhibitors. Using the highly contributing residues, we proposed the locations and the shapes of verapamil binding site and "R" site, and suggested the possible position of "H" site.

Introduction

P-glycoprotein (P-gp), also known as ABCB1 and MDR1, is a member of the ATP-Binding Cassette (ABC) transporter superfamily.1 P-gp is the first identified eukaryotic ABC transporter.2 Since P-gp mediates the transmembrane transport of many drugs that are unrelated in structure and function, it is involved in the multi-drug resistance (MDR) in human cancer.3

The drug resistance of cancers to anticancer drugs originates in the decrease of absorption, increase of the rate of metabolism, or increase of drug excretion.4 When some cancer cells show resistance to structurally unrelated drugs, this phenomenon is called MDR.4 Generally, MDR is the result of expression of efflux pumps like P-gp that has large spectrum of substrates.4 P-gp reduces the concentration of anticancer agents in various cancer cells by pumping out treated drugs and this results in MDR.5 Since up to 40％ of human cancers express MDR phenomenon, a lot of anticancer agents become ineffective to cancer patients.2,3 So, Pgp inhibitor development has remained as one of the main concerns in anticancer therapy, even though several inhibitors for P-gp have been developed.2,3

P-gp is 170kDa membrane protein with 1280 amino acids6 and is widely distributed in several human cells such as intestinal enterocyte, kidney proximal tubule, hepatocyte and brain endothelia.7 Like most of ABC transporters, it has two homologous halves each contains one nucleotide binding domain(NBD) which has ATP-binding site and 6 transmembrane( TM) helices which form multiple drug-binding sites.2,8 While the NBDs of P-gp show high sequence homology with most ABC transporters, the sequence homology of TMDs are comparatively low.2 It is known that the drug transport by P-gp is energy dependent process, however the mechanism of the coupling between ATP hydrolysis and drug transport machinery is not clearly explained.2,6 To understand the P-gp transport mechanism, it is necessary to define drug binding pockets and find the protein-drug interaction modes.

To determine the drug binding sites of the P-gp, various methodologies have been introduced and based on the experimental results, it was proposed that P-gp has at least two drug binding sites that are named “H” site and “R” site, referring to hoechst33342 and rhodamine123, respectively.9 Hoechst33342 has been used to several experiments as a P-gp ligand because it emits fluorescence in the lipid and loses fluorescence in aqueous solution.10 Similarly, rhodamine123 is widely used P-gp ligand with fluorescence emittion.10 In 1997, Shapiro and Ling experimentally confirmed these two distinct sites, “H” site and “R” site, by using fluorescence monitoring measurement.10 In the study, they showed that these two sites have a positively cooperative interaction.10 Furthermore, Martin et al. suggested that there are at least four distinct drug interaction sites on P-gp.3 The possibility of these multiple drug binding sites is thought to be a cause of remarkably large spectrum of P-gp substrates.3 Many researchers have been trying to determine these binding sites with several methods. Through drug binding site mapping studies using a Fluorescence Resonance Energy Transfer (FRET) the locations of “H” and “R” site were illuminated.9,11 And Loo et al. suggested interaction residues of “R” site and verapamil binding site by mutation studies.12-18 However, the information of drug binding sites and the residues involved in the interaction have not been clearly revealed yet, even though the interactions between P-gp and substrates or inhibitors have been widely studied.

In this study, the human P-gp structure will be built by homology modeling and the docking simulations will be performed using the modeled structure with verapamil, hoechst33342, and rhodamine123 to construct the interaction map between human P-gp and its inhibitors. With the highly contributing residues, we will propose the shapes and binding poses of verapamil binding site and “R” site, and suggest the possible position of “H” site.

 

Method

Homology Modeling. Due to the absence of crystal structure of human P-gp(hP-gp), the 3D structure of hP-gp was built by homology modeling. The protein sequence of hP-gp was retrieved from NCBI protein database (www.ncbi. nlm.nih.gov/protein). The template protein for homology modeling was selected by Blast search (blast.ncbi.nlm.nih. gov/). Sequence alignment of the reference protein and hPgp was performed using ClustalX(2.1).19 Subsequently, homology modeling was performed with MODELLER(9.10).20 The implicit membrane model for hP-gp was generated on Discovery Studio(3.5).21 After that, the loop refinement, side chain refinement and energy minimization of generated model were performed on Discovery Studio. Profiles-3D protocol in Discovery Studio and MODELLER were used for the final validation of modeled structure.

Docking Study. Three hP-gp inhibitors (verapamil, rhodamine123, and hoechst33342) were selected as docking ligands for the protein-ligand docking. The chemical structures of the inhibitors were obtained from the Pubchem project (pubchem.ncbi.nlm.nih.gov/). The structures of inhibitors are shown in Figure 1. Rhodamine123 and hoechst33342 represent the ligand for “R” site and “H” site, respectively, and verapamil is one of the well-known P-gp ligand.

Figure 1.The structures and molecular weights of the P-gp inhibitors verapamil, rhodamine123, and hoechst33342.

Verapamil was selected as inhibitor to be investigated through docking study because its binding residues were proposed using mutation studies by Loo et al.13-16,22,23 So it is necessary to confirm the accuracy of our docking results by comparing with experimentally proposed verapamil binding residues, and to compare docking results of verapamil with those of rhodamine123 and hoechst33342.

We determined the inhibitor binding sites of hP-gp at the TM region which was surrounded by 12 helices. The docking simulations of hP-gp with inhibitors were performed using LibDock protocol in Discovery Studio. Up to 500 binding poses were generated for each inhibitor. Subsequently, the amino acids within 3.5 Å from docked inhibitor were taken from each pose and the frequencies of taken amino acids were calculated for each inhibitor. Based on the results of these frequency calculations, the ligand binding sites of hP-gp will be defined and compared with experimentally proposed binding sites.

 

Results and Discussion

Homology Modeling. In order to select the template structure to construct the homology model, protein sequence homology search was performed using hP-gp sequence with the ABC transporter proteins of which crystal structures were known such as mouse P-gp (mP-gp) (PDB ID: 3G5U),24 Sav1866 (PDB ID: 2HYD),25 MsbA (PDB ID: 3B60),26 and P-gp of C. elegans (PDB ID: 4F4C).27 The crystal structure of P-gp of Mus musculus (mP-gp)24 was selected as a template protein because mP-gp had the highest sequence homology score with hP-gp by blast search. Three X-ray structures of mP-gp were obtained by Aller et al., one of them was apostructure and the others were protein-ligand complex structures.24 The superimpositions of apo-structure (PDB ID: 3G5U) and complex structures (PDB ID: 3G60 and 3G61) showed that their structures had no significant difference (RMSD: 0.65 Å and 0.67 Å, respectively). Hence, we finally selected apo-structure as the template for homology modeling.

ClustalX2 alignment result of hP-gp and mP-gp showed 87％ of sequence identity and 94.2％ of similarity which were very satisfactory to construct the homology model. The alignment result was further used to make hP-gp homology model by the modeling tool of MODELLER. The calculation yielded 10 models for hP-gp. Among the generated 10 models, we selected the best model by considering DOPE assessment score (calculated by MODELLER program) and validation results (calculated by MODELLER and Profiles- 3D). The Ramachandran plot of homology modeled hP-gp (Fig. 2(a)) shows that most of the amino acids satisfy reasonable Ф-Ψ angles, and the plot of transmembrane region (Fig. 2(b)) indicates that most of the amino acids in the membrane make proper α-helix structure.

Figure 2.The Ramachandran plots of the human P-glycoprotein homology model. (a) The Ramachandran plot for overall structure of the hP-gp homology model. (b) The Ramachandran plot for transmembrane region of the hP-gp homology model.

It was argued that the mP-gp crystal structure was incompatible in the orientations of TM3, TM4, and TM5 with the results of arginine-scanning mutagenesis for identifying drug translocation pathway residues by Loo et al..27,28 Furthermore, it was suggested that the crystal structure of C. elegans had more reliable orientation in those helices.27 However, we suspected that the probable errors driven by some residue orientation incompatibilities in the mP-gp crystal structure were not large enough to disrupt the inhibitor docking simulations for following reasons; The identity and similarity between mP-gp and hP-gp are 87％ and 94.2％, respectively, while those between C. elegans P-gp and hP-gp are 43.8％ and 65.2％, respectively. Although the values between C. elegans P-gp and hP-gp are not improper for homology modeling, the identity and similarity values between mP-gp and hP-gp are dominantly high. Furthermore, we have also generated homology model using C. elegans P-gp crystal structure(PDB ID: 4F4C) as a template, but the orientations of residues those were proposed as binding residues in the mutagenesis by Loo et al. had no significant difference between the model by mP-gp and the model by C. elegans P-gp (the result not shown).

The crystal structure of mP-gp is heterodimer that has asymmetrical monomers and presents inward-facing conformation.24 In the structure, 12 TM helices form a large cavity which is accessible from the cytosol and inner-membrane leaflet, and the volume of the cavity is under 6000 Å3.24 The bound inhibitors are located between TM6 and TM12 helices and surrounded by the residues on TM1, 2, 5, 6, 7, 8, 9, 11, and 12.24 The number of residues which compose the binding pocket of the inhibitors in mP-gp crystal structure are 21 in total, among them only one residue is different from that of hP-gp.24 When comparing the structure of the mP-gp with the hP-gp focusing on those 21 residues, the orientations of the residues in the mP-gp crystal structure were equivalent to corresponding amino acid residues in modeled hP-gp structure (data not shown). These highly conserved residues between the hP-gp and the mP-gp at ligand binding site confirmed the reliability of generated hPgp homology model. Therefore, we performed docking simulation with the modeled hP-gp structure.

Molecular Docking. It was suggested that the P-gp ligands approach to the binding site of the P-gp from the inner membrane leaflet, not from the cytosol29 and both neutral and positively-charged hydrophobic drugs are probably approach to the TM regions of the P-gp from the membrane lipid.4 Several researchers proved that the drug binding sites are located on TM region of P-gp8, and some of the important binding or binding-related residues of verapamil and rhodamine were proposed.12 We will include all the possible residues in the pocket to define the site for docking simulations and will compare the results with the experimentally proposed verapamil and rhodamine binding residues. Qu and Sharom proposed the location of “H”site as cytoplasmic leaflet of cell membrane based on FRET study though the explicit interaction residues with hoechst33342 were not indicated.11 Considering the proposed “H”site by Qu and Sharom, the docking simulations were performed within large sphere as described in Figure 3(c) in order to include experimentally proposed “H”site and ligand binding residues in the pocket altogether. The ligand binding residues proposed by Loo et al. are summarized in Table 1 and shown in Figure 3(a) and 3(b).

Figure 3.(a), (b) Homology model of hP-gp structure. Proposed verapamil and rhodamine binding residues in mutation studies by Loo et al. are colored by yellow. (c) Docking region for docking simulation. (a) is bottom-up view. The region between the green plane and the blue plane refers to membrane, and cytoplasm is below the blue plane in (b) and (c).

Table 1.The a.a. that contribute to the binding pocket of verapamil and “R” site of the hP-gp in the studies by Loo et al. are summarized

Since the binding pocket of the hP-gp is large enough to have multi-binding sites and multi-drug specificities, we performed statistical analysis with 500 docking poses rather than the best docking pose to find the residues that constitute binding pocket. Verapamil, rhodamine123, and hoechst33342 were docked to the modeled hP-gp structure 500 times. In the docking simulations for rhodamine123 and hoechst33342, docking poses were generated less than 500 times. The residues which had the frequency above 60％ were regarded as the interaction residues in our results. The frequencies for each binding residue in each inhibitor are summarized by the graph in Figure 4. In order to confirm the hypothesis that verapamil binding site, “R” site, and “H” site are clearly distinct sites, the residues significantly contributed to the interaction with the inhibitors were summarized by Table 2 and were clustered by Figure 5.

Figure 4.The results of docking analysis. The residues which have conspicuously low frequency are not presented. (a) The results of verapamil docking. The residues which have the frequency over 300 times in overall 500 binding poses are checked by black circle. The residues checked by gray arrow are the verapamil binding residues proposed by mutation studies (Loo et al.). (b) The results of rhodamine123 docking. The residues which have the frequency over 170 times in overall 285 binding poses are checked by black circle. The residues checked by gray arrow are the rhodamine binding residues proposed by mutation studies (Loo et al.). (c) The results of hoechst33342 docking. The residues which have the frequency over 215 times in overall 360 binding poses are checked by black circle.

Table 2.The a.a. that contribute to the binding pocket of each inhibitor of the hP-gp in docking simulation are summarized. The cells with gray color refer to proposed residues by both Loo et al. and our study

Figure 5.The a.a. residues that highly contribute to the binding pocket of hP-gp are clustered into V, R, and H referring to verapamil, rhodamine123, and hoechst33342 binding pocket, respectively. Some a.a. residues are shared by all three binding pockets, some a.a. are shared by two pockets, and some a.a. are involved only one pocket.

Verapamil Binding Site. The location of the verapamil binding residues proposed by Loo et al. and those obtained from our docking simulation are shown in Figure 6(a) and 6(b), respectively. Although both Loo et al. and our proposed pocket residue members included 14 residues respectively (Tables 1 and 2), the components were not identical. Only four residues, Phe728, Leu975, Val982, and Ala985, were included in both Loo et al. and our proposed residues. Although the frequency of Ala985 was less than 300 in our results, it was selected as the interaction residue because of relatively high frequency. From the comparison, it seems that there are considerable discrepancies between the verapamil binding residues proposed by Loo et al. and our study.

Figure 6.Verapamil binding residues (bottom-up view). (a) Verapamil binding residues proposed by mutation studies (Loo et al.) (b) Verapamil binding residues included in our results. Gray residues are the residues which have the frequency over 200 times, but less than 300 times. The residues included in both the results of Loo et al. and our results are indicated with circles in (a) and (b), and colored with sky-blue in (b).

In the modeled hP-gp structure, some verapamil binding residues proposed by Loo et al. (Ser222, Ile306, Gly984, Ala342, and so on) may not have proper location for the inhibitors to contact with (Fig. 6(a)). Notably, the above mentioned improperly located residues are found in crystal structure of mP-gp. It might be due to the absence of membrane bilayer in crystal structure or to the induced-fit binding mechanism. Hence, it is reasonable that those kinds of residues have no interaction with inhibitors in docking simulations. In addition, the residues proposed by Loo et al. are dispersed widely in the modeled hP-gp structure (Fig. 6(a)). Excluding the probability that verapamil could bind to all the region of hP-gp in which the experimentally proposed residues are located, it could be interpreted that experimentally proposed residues might include the residues which do not bind to verapamil directly, but are related to binding by other mechanisms such as making hydrophobic region or inducing conformational changes when inhibitor binds to hP-gp. Furthermore, it could not be neglected that point mutation might influence on the results of experiments by changing the structures or properties of binding site in the mutation study. Considering discussed above, the overall shape and key residues of verapamil binding site are thought to be shared by both the results of Loo et al. and our docking simulations. Hence, we could illuminate the binding site for verapamil from the results of experimental studies and our docking simulations.

Among the 4 residues those are included in both Loo et al. and our proposed residue members, 3 residues are located in TM12 helix (Fig. 6). And the rest one residue, Phe728, is neighboring residue of the 3 residues of TM12 helix. These results indicate that TM12 helix has important role in verapamil binding and that Phe728, Leu975, Val982, and Ala985 altogether could make hydrophobic interactions with verapamil. Taken together, we assumed that verapamil binding site is the inner surface of TM12 and TM7 helices. The potential verapamil binding pose is shown in Figure 7. In the proposed binding pose, verapamil makes hydrogen bond with Tyr307 and π-π interaction with Tyr953 and Phe978. The remaining residues such as Phe728, Leu975, and Val982 could make hydrophobic environments.

Figure 7.Verapamil binding pose: Pink molecule is verapamil. The residues colored by sky-blue indicate those included in both the results of Loo et al. and our results. The residues colored by yellow are those have interactions with docked verapamil among the proposed residues by docking simulation. (a) Bottom-up view. The helix that includes Phe728 is TM7 and the helix that includes Leu975, Val982, and Ala985 is TM12. (b) The helix in the center is TM12.

“R” Site. Although the exact location of the “R” site has not been clarified, Loo et al. performed mutation studies using rhodamine B.17,18 We hypothesized that the proposed residues from above stated studies could be possible components of “R” site. The comparison of “R” site residues proposed by mutation studies and docking simulations is shown in Figure 8. Six binding residues were proposed by the studies of Loo et al. (Table 1) whereas 12 residues were selected by docking simulations (Table 2). Particularly, only two residues were presented in both the results of Loo et al. and our docking study. This inconsistency is analogous to the case of verapamil binding site. Therefore, it could be explained similarly with the reasons discussed in the results of verapamil; Improper residue orientation for inhibitor binding, the probability of existence of residues which do not interact with inhibitor directly, and the influence of mutation to the structures or properties. Furthermore, the possibility that there could be other binding residues should be considered since only 6 members of “R” site were proposed by experimental studies.

Figure 8.“R” site residues (bottom-up view). (a) “R” site residues proposed by mutation studies (Loo et al.) (b) “R” site residues included in our results. The residues included in both the results of Loo et al. and our results are indicated with circles in (a) and (b), and colored with sky-blue in (b).

Among the 6 residues proposed by Loo et al., 2 residues are in TM6 helix and 3 residues including common residues in the proposition of Loo et al. and our study, Leu975 and Val982, are in TM12 helix (Fig. 8). From these results, we could assume that TM12 helix is also important to rhodamine binding in the hP-gp and “R” site region estimated by the experimental studies is the region between TM12 helix and TM6 helix. This estimated region seems to have considerable consensus with the overall “R” site region proposed by docking simulations in the location (Fig. 8). Thus we could conclude that the shape and location of the “R” site are shared by both experimental studies and our study, even though the inconsistency of interaction residues. The proposed docking pose of rhodamine123 by our study is shown in Figure 9. In the proposed binding pose, rhodamine123 makes π-π interactions with Phe336 and Phe978, and makes hydrogen bond with Tyr953. The remaining hydrophobic residues could make hydrophobic environments.

Figure 9.Rhodamine123 binding pose: Purple molecule is rhodamine123. The residues colored by sky-blue indicate the residues proposed by mutation studies (Loo et al.). The residues colored by yellow are those have interactions with docked rhodamine123 among the proposed residues by docking simulation. (a) Bottom-up view. (b) The helix in the center is TM12. (c) The left helix is TM6 and right one is TM12.

From our proposition of verapamil binding site and “R” site, TM12 helix is thought to have important role in the drug binding. Indeed, TM12 helix has been suggested to play a key-role in the process of drug binding and transport. 30,31 From the docking results, we could also conclude that the verapamil binding site and “R” site share a lot of common region since all the binding residues for rhodamine123 are included in those for verapamil. Furthermore, in the result of mutation studies by Loo et al., the overall shape and location of verapamil binding site and “R” site are similar despite the residues of each site are not consistent. The recent docking study, however, proposed that the verapamil binding site was different from “R” site, not only from “H” site.32 The result of this docking study is controversial comparing to our result of docking study that verapamil binding site is similar to “R” site. Verapamil was known as the ligand which preferentially binds to “H” site in 1990’s.10 However, verapamil may not bind only one site and it is classified as a drug which bind to both “R” site and “H” site.8 Therefore, although it is not well clarified if verapamil binds to both sites, we might conclude from the results of Loo et al. and our study that verapamil binding site and “R” site are not obviously distinguishable and they are overlapped in some parts.

“H” Site. According to our docking results, only the half of the binding residues contributed for rhodamine123 were included in those for hoechst33342 whereas all the binding residues contacted with rhodamine123 presented in verapamil binding site. Furthermore, hoechst33342 had its own three binding residues which had interaction with only hoechst33342, even if hoechst33342 shared many of binding residues with verapamil (Fig. 9). As we assumed, these results could be interpreted that “H” site and “R” site occupy considerably distinct regions each other and “H” site includes unique region as binding site although it overlaps with the verapamil binding site and “R” site in part.

There is no experimental information about the binding residues of hoechst33342. But it has been known that “R” site and “H” site interact each other with positively cooperative manner, and the ligands of each site bind to P-gp noncompetitively.10 Therefore, “H” site is thought to be located at the site that could be discriminated against “R” site. In the FRET study by Qu and Sharom, “H” site was suggested to be in the cytoplasmic leaflet of membrane.11 In our results, the site which accesses those conditions was the neighborhood of the residues Leu304(TM5), Gln725(TM7), and Leu762(TM8) (Fig. 10(a)). This site did not show entire consensus with the region illuminated by the FRET study. Nevertheless, considering that the FRET study measures the distance instead of the location and that the results could be influenced by vicinal residues, the region mentioned above might be a proper candidate of “H” site. Interestingly, the three residues, Leu304, Gln725, and Leu762, were also included in the ‘lower’ binding site of QZ59-SSS in the mPgp crystal structure.24 Hence, we suggested that TM5, 7, and 8 helices form “H” site at the region around these three residues.

Figure 10.Suggested hoechst33342 binding pose. (a) The green and blue lines indicate membrane. (c) Presented four helices are TM5, TM8, TM7 and TM12 from the left.

To confirm if other “H” site ligands also bound to this region, we performed additional docking simulations using quercetin and colchicine which have been known to bind to “H” site. In these docking results, we could identify that both compounds bound the region we suggested to be “H” site (the results not shown). So we concluded that the “H” site is likely to be in the region including residues Leu304, Gln725, and Leu762 between TM5, TM7, TM8 and TM12 helices. The suggested binding pose of hoechst33342 makes hydrogen bond with Gln725, and two leucine residues seem to make hydrophobic environment with surrounding residues (Fig. 10). The suggested “H” site overlaps with verapamil binding site partially, but it is closer to inner leaflet than the verapamil binding site. When compared with “R” site, “H” site is located in the region distinguishable from the location of “R” site (Fig. 11).

Figure 11.Predicted binding pose comparison of (a) verapamil, (b) rhodamine123, and (c) hoechst33342. The green and blue lines indicate membrane.

 

Conclusion

In this article, we suggested the binding sites of three inhibitors of the hP-gp by homology modeling and docking study. The TM12 helix is thought to be important for inhibitor binding, and the verapamil binding site could be inner surface of TM12 and TM7 in the outer-leaflet of membrane. We proposed that “R” site is between TM12 and TM6 and shares common region with the verapamil binding site. On the other hand, “H” site might be in the region between TM5, 7, 8, 12 and has its own region although it seems to share the binding site with verapamil binding site partially. The shapes and locations of proposed binding sites by docking simulations are consistent with the results of experimental studies, despite the fact that not all residues show the consensus between the results of experimental studies and docking study.

Unfortunately, there are not abundant quantitative data such as the activities and binding affinities of ligands yet. This makes it difficult to validate and concretize the ligand binding sites and interactions. The absence of available crystal structure of hP-gp is another limitation of in silico approaches. Nevertheless, these approaches such as docking and simulation method could be helpful to elucidate the interactions of the hP-gp with inhibitors. Furthermore, using more data accessible by further experimental studies we could improve the accuracy of in silico approaches.