|Home | About | Journals | Submit | Contact Us | Français|
Homology modeling of the human A2A adenosine receptor (AR) based on bovine rhodopsin predicted a protein structure that was very similar to the recently determined crystallographic structure. The inaccuracy of previous antagonist docking is related to the loop structure of rhodopsin being carried over to the model of the A2A AR and was rectified when the β2-adrenergic receptor was used as a template is used for homology modeling. Docking of the triazolotriazine antagonist ligand ZM241385 1 was greatly improved by including water molecules of the X-ray structure or by using a constraint from mutagenesis. Automatic agonists docking to both a new homology modeled receptor and the A2A AR crystallographic structure produced similar results. Heterocyclic nitrogen atoms closely corresponded when the docked adenine moiety of agonists and 1 were overlayed. The cumulative mutagenesis data, which support the proposed mode of agonist docking, can be reexamined in light of the crystallographic structure. Thus, homology modeling of GPCRs remains a useful technique in probing the structure of the protein and predicting modes of ligand docking.
Crystallographic structural data are available today for four different GPCRs: bovine rhodopsin,1 human β2-adrenergic receptor,2 turkey β1-adrenergic receptor,3 and human A2A adenosine receptor (AR).4 All of these receptors are transmembrane proteins consisting of seven α-helices connected by three extracellular (ELs) and three intracellular loops (ILs). The general configuration of the transmembrane domains (TMs) is very similar for all GPCRs. In particular, the heavy atoms of TMs of the A2A AR and β2AR can be superimposed with a RMSD of 1.90Å, and the superimposition of α-helices of the A2A AR and rhodopsin provided a RMSD value of 2.16Å. The differences in the configuration of TMs of rhodpsin and the β2AR are represented by a RMSD of 1.90Å.
Since for a long time the only GPCR for which an experimental structure was available was bovine rhodopsin, rhodopsin was widely used as a template for homology modeling of other GPCRs. One of the first molecular models constructed for GPCRs was a model of the human A2A AR based on the electron density map of rhodopsin.5 During recent years, numerous homology models have been generated for various GPCRs, including the A2A and other subtypes of ARs.6-8 Many of the models proposed were successfully used for investigation of ligand-receptor interactions and the development of novel biologically active compounds, in particular, for the ARs.9 Now with an experimental structure of the A2A AR available it is possible to evaluate the quality of the proposed models and to refine hypotheses concerning the ligand binding modes. With this aim we compared our previously published rhodpsin-based model of the A2A AR6 (pdb code: 1UPE) with the X-ray structure of this receptor.
The docking orientation of the antagonist ligand 4-2-[7-amino-2-(2-furyl)-1,2,4-triazolo[1,5-a][1,3,5]triazin-5-yl-amino]ethylphenol 1 (ZM241385)10 in the human A2A AR was different from the antagonist docking modes typically predicted previously by modeling. A predicted antagonist binding site, e.g., for N-[9-chloro-2-(2-furanyl)[1,2,4]triazolo[1,5-c]quinazolin-5-amine 2 (CGS15943),6 corresponded more closely to the position of the retinal binding site in rhodopsin and the binding site of the inverse agonist carazolol in the β2-adrenergic receptor. In this study we have evaluated the use of molecular modeling of GPCRs and ligand docking in light of the newly reported crystallographic structures of the A2A AR and other GPCRs. Costanzi studied the β2-adrenergic receptor structure and its docked ligand to conclude that GPCR modeling is applicable to the design of site-directed mutagenesis experiments and to drug discovery.11 We have extended the analysis to the adenosine system.
All atoms of the α-helical TMs of a previously published rhodopsin-based homology model of the A2A AR6 and the X-ray structure of this receptor were aligned with an RMSD value for all TM atoms of 2.37Å. Not surprisingly, the configuration and orientation of the TMs of the theoretical model and experimental structure of the A2A AR were found to be very similar (Fig. 1). In contrast, the configurations of the ELs are significantly different in these two structures.
Previously, various residues located mostly in TMs 3, 5, 6, and 7 were predicted with modeling to be involved in ligand recognition (Supporting information, Table S1).6,12,13 In particular, it was suggested that Ile80 (3.28), Val84 (3.32), Leu85 (3.33), Thr88 (3.36), Gln89 (3.37), Ile135 (4.56), Leu167 (EL2), Phe168 (EL2), Asn181 (5.42), Phe182 (5.43), Val186 (5.47), Trp246 (6.48), Leu249 (6.51), His250 (6.52), Asn253 (6.55), Ile274 (7.39), Ser277 (7.42), His278 (7.43) are located in proximity to antagonist 2 and involved in interactions with the ligand6 (the numbers in parentheses correspond to the Ballesteros-Weinstein indexing system).14 The importance of residues located in the abovementioned TM regions and in EL2 was confirmed by site-directed mutagenesis data obtained not only for the A2A AR (Table S1), but also for other ARs and diverse GPCRs.15 In agreement with the results of molecular modeling and pharmacological studies, many of these critical residues are located in proximity to the antagonist ligand 1 in the recent crystallographic structure of the A2A AR. Examples of residues in proximity to the bound 1 in the X-ray structure that were predicted by modeling to be involved in antagonist recognition are: Leu85 (3.33), Phe168 (EL2), Asn181 (5.42), Trp246 (6.48), Leu249 (6.51), Asn253 (6.55), Ile274 (7.39) (Supporting information, Table S1).6 The differences in orientation of side chains of these amino acid residues between the homology model and X-ray structure were examined. It should be noted that the antagonist bound X-ray structure of the A2A AR was compared with the original rhodopsin-based model of the A2A AR.6 Since the inactive, retinal-bound structure of rhodopsin was used as a template for the modeling, the A2A AR model is more related to the ground state of the receptor and therefore can be compared with an antagonist occupied X-ray structure.
As shown in Fig. 2, homology modeling demonstrated an excellent ability to predict the correct orientation of residue side chains. The residues in their predicted orientation fit very well to the experimental structure. In addition, it was found that superimposition of Cα-atoms of selected residues in the ligand binding pocket, i.e., Thr88 (3.36), Asn181 (5.42), Trp246 (6.48), His250 (6.52), Asn253 (6.55), Ser277 (7.42), His278 (7.43), Ile274 (7.39), provided an RMSD value of 0.90Å. The proximity of the side chains of conserved His278 in TM7 and a Glu13 in TM1, previously identified through the combination of modeling and mutagenesis for both A2A and A3 ARs,16,17 was confirmed in the X-ray structure. Also, the H-bonding between the backbone carbonyl group of Ser281 (7.46) and the side chain amino group of Asn24 (1.50), and a H-bond between Ser281 (7.46) and Asp52 (2.50) side chains were previously predicted by the modeling. These interactions were observed in the recent A2A AR crystal structure. In addition, the proximity of Asn24 (1.50) to Asp52 (2.50) was suggested based on the modeling studies. In the A2A AR crystal structure the side chains of these two residues appear to be connected through a water molecule. Thus, the comparison performed clearly demonstrated the ability of a homology modeling approach to provide an accurate model of a GPCR protein structure, especially within its TM domains.
Interestingly, several residues shown by mutagenesis to be critical for antagonist binding are not involved in direct contact with 1 in the recent A2A AR crystal structure. Thus, the involvement might conceivably be through intervening amino acid residues or a conserved bound water molecule. In particular, it was demonstrated that replacement of Phe182 (5.43), Ser281(7.46) or Phe257 (5.59) with Ala resulted in loss of agonist and antagonist radioligand binding.5 However, the F182Y mutation resulted in only 2-fold decreased antagonist affinity. In the crystal structure Phe182 is involved in π-π interactions with His250 (6.52) and may play a role in helix-helix packing, as was suggested by Kim and co-workers.5 More importantly, His250 is directly involved in the interaction with ligands, and the replacement of Phe182 may result in the repositioning of the His residue. Ser281 (7.46) was shown to be critical for both agonist and antagonist binding.5 As mentioned above, in the X-ray structure this residue is involved in inter-helical interactions. Namely, Ser281 forms a H-bond with Asp52 (2.50), which is highly conserved among all GPCRs. In contrast, Phe257 does not form significant interactions with other residues in the crystal structure, and this residue is oriented outwardly from the binding cavity. Thus, the involvement of certain residues in ligand recognition in the A2A AR is now shown to be indirect.
To investigate an ability of molecular docking to reproduce an experimentally observed ligand binding mode, 1 was re-docked to the crystal structure of the A2A AR. Initially the Glide XP software was utilized to automatically dock 1 to the rigid binding site of the A2A AR without any constraints and without the presence of water molecules. The results obtained indicated that the orientation of the ligand inside the receptor was inverted along the long axis as compared with the crystal structure (Fig. 3A). The differences in the position and orientation of 1 docked to the crystal structure and its experimental position was characterized by an RMSD1 value of 10.27 Å.
The binding mode of 1 obtained after molecular docking performed with the MOE software was in better agreement with the experimental binding mode (Fig. 3B). However, due to the differences in the orientation of the phenylethylamine chain the RMSD1 value of 6.11 Å between predicted and experimental positions of 1 was relatively high.
Compound 1 was also docked to the crystal structure of the A2A AR containing the three water molecules that coordinated with the ligand in the X-ray structure. The ligand binding mode obtained was in an excellent agreement with the experimental binding mode (Fig. 3C), RMSD1 = 0.90Å. In practice, prediction of the accurate position of water molecules inside a ligand binding site is a challenge, which would hinder the general application of this finding to other GPCR modeling.
However, we found that data from site-directed mutagenesis to establish extra constraints in the model can greatly improve the accuracy of prediction. In particular, it was shown that Asn6.55 is critical for ligand binding at all four subtypes of ARs. Also, this residue formed a H-bond with the exocyclic amino group of 1. Molecular docking of 1 to the rigid A2A AR performed with the H-bond interaction between Asn253 (6.55) and the exocyclic amino group of exocyclic amino group of 1 used as a constraint provided a reasonable overlay between predicted and experimental binding modes of 1 (RMSD1 = 2.09 Å, Fig. 3D). In order to study the role of flexibility the receptor side chains in the docking studies, the InducedFit approach was also utilized to re-dock 1 to the crystal structure. The analysis of the binding mode of 1 obtained after InducedFit docking revealed that the core heterocycle of the ligand had the same position and orientation inside the receptor as appeared in the crystal structure.
In contrast, the predicted conformation of the flexible phenylethylamine chain of 1 was found different as compared with the X-ray structure. In the crystal structure, the OH-group of the phenyl ring is H-bonded with a water molecule. Since InducedFit docking was performed in the absence of water molecules, the interaction between the phenyl ring and water was obviously impossible. However, a H-bonded functional group is more favorable energetically than noncoordinated functional group. For this reason, a H-bond between the OH-group of 1 and the backbone oxygen atom of Ile67 was formed during the InducedFit docking resulting in different conformation of the phenylethylamine chain. The RMSD1 of 3.58 Å was found between the experimental position of 1 inside the A2A AR and its position obtained after InducedFit docking (Fig. 3E).
In the A2A AR X-ray structure, the antagonist binding site was shifted closer to TMs 6 and 7 and closer to the exofacial side than the bound ligands carazolol and cyanopindolol in the previous adrenergic receptor structures.2-4 Therefore, we examined the modeling of 1 binding in greater detail. Several molecular models of the ARs with various antagonists docked have been proposed. Recently, a putative binding mode of 1 was examined with the molecular modeling approach by Yuzlenco et al.23 However, most of the previous models concerned antagonists other than 1. In most of reported models, the amino acid residues predicted to be involved in antagonist recognition were in good agreement with the X-ray structure of the A2A AR. For example, Kim and co-workers suggested an involvement of Leu85 (3.33), Phe168 (EL2), Asn253 (6.55), Ile274 (7.39) in the binding of the antagonist 2.6 In particular, the model proposed that Asn253 (6.55) was H-bonded to the exocyclic amino group of 2, and Phe168 was found in proximity to its triazoloquinazoline core. Here, we modeled a complex of the A2A AR, derived from a new homology model, with its antagonist 1 to compare the results of molecular modeling with the crystal structure.
First of all, a homology model of the A2A AR to be utilized for studies of the docking of 1 was built based on the recently reported crystal structure of the β2AR. Since the configuration of EL2 is different in various GPCR crystal structures, an accurate prediction of the configuration of the entire loop was not attempted. However, the disulfide bridge formed between two conserved Cys residues located in TM3 and EL2 was observed in all available GPCR crystal structures. In view of this, in our model of the A2A AR only part of EL2 (Cys166 - Asn175) was constructed.
The obtained homology model was utilized for automatic docking of 1, and the ligand binding mode obtained was found to be in good agreement with the experimental structure (Fig. 7). In particular, in the resulting model Asn253 (6.55) formed a H-bond with the exocyclic amino group of 1, and Phe168 involved in π-π interactions with the ligand triazolotriazine ring. The furanyl ring was oriented inside the receptor and was found in proximity to Leu85 (3.33), Thr88 (3.36), Trp246 (6.48), Leu249 (6.51). The 4-hydroxyphenylethylamine chain of 1 was oriented toward the extracellular region of the A2A AR and was located between Ala265 (EL3) and Leu167 (EL2). Thus, the obtained binding mode of 1 corresponded very well to its binding mode observed in the crystal structure. This finding demonstrates the general applicability of the homology modeling approach to construct realistic models of GPCRs. The homology models of the A2A AR based on rhodopsin and on the β2AR are very similar in the region of the ligand binding pocket. Hence, the previous docking experiments, in which the ligands were not affected by the different loop structure, are expected to be valid.
The published crystal structure of the A2A AR with its antagonist 1 represents a conformation of the receptor resembling the inactive state. However, as in previous homology modeling studies based on the rhodopsin ground state, we assumed that this structure would be able to fit small agonists as well. With this assumption, the molecular docking of the principal native AR agonist, adenosine 3, and the synthetic agonist 5′-N-ethylcarboxamidoadenosine (4, NECA) to the A2A AR X-ray structure was performed automatically as described in the Methods section. It should be noted that these computational studies were conducted using standard, automatic techniques implemented in currently favored and widely-used software such as MOE and the Schrödinger Package.18,19 Standard procedures for protein refinement, such as the addition of missing hydrogens and the assignment of bond orders, were automatically performed during the protein preparation step. To take into account the flexibility of the receptor side chains, the InducedFit docking approach was utilized.20
The results of molecular docking revealed that both studied agonists fit well the X-ray structure of the A2A AR and independently adopted very similar binding modes (Fig. 4). In particular, the complexes obtained indicated that the N6-amino group of the ligands could form a H-bond with Asn253 (6.55) and Glu169 (EL2). The adenine ring was involved in π-π stacking with the aromatic ring of Phe168 (EL2), as was the heterocycle of 1 in the X-ray structure. The 2′- and 3′-OH groups were observed in proximity to Ser277 (7.42) and His728 (7.43) and can be involved in H-bonding with these residues. Also, in the case of adenosine, its 5′-OH group formed a H-bond with Ser277 (7.42). The NH-group of the ethylcarboxamide moiety of 4, which is important for the activation of ARs,21 was found at a distance of 4.5Å from the OH-group of Thr88 (3.36). Both, Ser277 (7.42) and Thr88 (3.36) were previously proposed with the site-directed mutagenesis to be critical for agonist but not for antagonist activity (Supporting information).
The orientation of side chains in proximity to the ligand obtained after InducedFit docking of NECA was compared with their orientation in the A2A AR crystal structure. The superimposition of all non-H atoms of amino acid residues located within 4 Å around 4 in the model obtained and the corresponding atoms of the receptor crystal structure provided an RMSD value of 0.65 Å. The side chains of most of the residues were found in their original conformations. However, the terminal amido group of Asn253 (6.55) was rotated by 68°, and the carboxylic group of Glu169 (EL2) was rotated by 67°. In addition, the side chain of Met270 (7.35) was shifted by 2.9 Å, and the indole ring of Trp246 (6.48) shifted approximately by 1 Å outward from 4.
Fig. 5 shows the general binding mode of the ARs agonists proposed previously based on molecular modeling and site-directed mutagenesis data. Thus, the binding mode obtained after the docking of adenosine and 4 to the experimental structure of the A2A AR is almost identical to their binding mode proposed in previous modeling studies and remains in good agreement with the mutagenesis data.
The previous studies of the binding modes of AR ligands suggested that the agonist adenine ring and an antagonist aromatic system should be accommodated in the same region within the putative binding site of the receptor.10 Our present results of molecular docking of adenosine and 4 have fully confirmed that hypothesis. The overlay of adenosine and 1 inside the A2A AR is shown in Fig. 6. The adenine moiety of adenosine and the triazolotriazine fragment of 1 have very similar positions inside the binding site and are involved in similar interactions with the receptor, namely with Phe168 (EL2), Glu169 (EL2), and Asn253 (6.55). In agreement with previous molecular modeling studies of ARs, the A2A AR agonists were found to be located between TMs 3, 5, 6 and 7. Moreover, in a comparison of corresponding atoms of the two docked structures both the ring N atoms and the exocyclic amines are in close proximity. The adenine C2-position of the agonists was oriented toward the extracellular part of the receptor. Therefore, C2 side chains of functionalized A2A AR agonists22 could occupy the same region as the phenylethylamino chain of 1. Previously, it was predicted with molecular modeling that 2-hydroxypropenyl moiety of 2-hydroxypropenyl-substituted derivatives of adenosine or 4 can interact with residues located in the second extracellular loop (EL2), namely with Cys166.7 Similarly, the phenylethylamine chain of 1 in the X-ray structure is also oriented toward the ELs. Thus, the orientation and contact points of the modeled agonist adenine ring inside the receptor were almost identical to those of the bicyclic triazolotriazine core of 1.
The crystallographic structure of the A2A AR was reasonably approximated using rhodopsin-based molecular modeling. Several key features of this receptor structure previously identified through the combination of modeling and mutagenesis were validated in the X-ray study.6,16 Overall, residues of the upper half of TMs 3, 5, 6, and 7 and of EL2 were predicted to be involved in ligand recognition, and this was the case in the X-ray structure. The crystallographic structure also revealed the direct involvement of EL2 in coordination of 1. The role of the ELs in ligand recognition has been well supported by modeling and mutagenesis.24,25 Thus, in cases where an X-ray structure is lacking, the homology modeled structure of a given family A GPCR could be considered a reasonable template for further modeling, with the possible exception of the extracellular and intracellular loops.
Molecular docking of 1 to its native binding site of the A2A AR crystal structure performed with different methods and different constraints resulted in different binding modes predicted. We found that the most accurate position of the ligand can be obtained by molecular docking performed in the presence of water molecules as they appeared in the ligand binding site of the crystal structure. However, if the receptor structure is unknown, the prediction of the position of the water molecules inside the binding site is problematic. Using the ligand-receptor interactions proposed based on the site-directed mutagenesis data also can significantly improve the accuracy of molecular docking. Our results demonstrated that the experimental binding mode of 1 can be reproduced using only one constraint, namely the H-bond between exocyclic amino group of the ligand and the side chain oxygen atom of Asn253 (6.55). Thus, our results demonstrated that molecular docking approach can be successfully applied to accurately reproduce an experimentally observed binding mode of a GPCR ligand. On the other hand, these data revealed the difficulties in the prediction of a ligand binding mode if the experimental data are limited. In addition and not surprisingly, even if the crystal structure of a receptor is available, the results of molecular docking are extremely dependent on the software and protocols used for the molecular docking. These are also the reasons for differences between the binding mode of 1 as appearing in the crystal structure and previously predicted binding modes of analogous AR antagonists.
The recently published X-ray structure of the A2A AR did not predict a docking mode for agonists. Therefore, a molecular modeling approach is still relevant, as it has been utilized in our recent drug design and receptor reengineering efforts that focused mainly on agonists. The fact that the agonist affinity of the T4-lysozyme construct of the A2A AR used in the X-ray study is even enhanced over the WT receptor4 would seem to justify this approach.
Several molecular models of ARs with their various agonists have been published; for a detailed review see Costanzi et al.12 Based on the predictions made with modeling combined with experimental data of site-directed mutagenesis, we proposed a typical binding mode for ARs agonists, shown in Fig. 4 for adenosine and the potent, nonselective agonist 4. It was suggested that the N6-amino group of the adenine ring of adenosine derivatives most likely forms a H-bond with Asn253 (6.55).6,26 Another hydrophilic residue, Asn181 (5.42), was predicted to be involved in H-bonding with the adenine ring N1-nitrogen atom. The ligand ribose ring was predicted to lie within a predominantly hydrophilic region. The 2′- and 3′-groups of the ribose ring were predicted to be located in proximity to His278 (7.43), which is conserved among the ARs. Experimental data indicate that the conserved His residue (7.43) is involved in H-bonding with the OH-groups of the ligand ribose ring. This proximity was indicated in a study of neoceptors, in which the ligand was modified in a fashion complementary to the single amino acid mutation of an AR.17,27,28 In the proposed molecular models, Ser277 (7.42) and Thr88 (3.36) were involved in H-bonding with the 5′-OH group of adenosine derivatives and with the ethylcarboxamide moiety of analogues of 4. Mutation of Thr88 in A2A AR neoceptors suggests a proximity to the terminal position of the 5′-uronamide.26 It was concluded that the ribose-binding region is the most important for receptor activation.
It is to be noted that the ligand-receptor interactions mentioned above were observed for two different conformations with respect to the glycosidic bond of adenosine and its derivatives. Molecular docking of various adenosine analogues to the rhodopsin-based models of all four ARs demonstrated the preferences for an anti-conformation of the glycosidic bond.7,37-39 However, computational studies of activation of the A2A AR performed by Kim and co-workers indicated the preferences of an intermediate conformation in the model obtained after rotation of the highly conserved Trp246 (6.48).6 In the present sstudy, both agonists, adenosine and 4, docked in the crystallographic structure of the A2A AR adopted an anti-conformation.
Docked agonist ligands display a close correspondence between the orientation of the adenine moiety and the heterocyclic moiety of the antagonist 1 in the crystallographic structure. Several sources of evidence support this mode of overlay of the agonist nucleobase and triazolotriazine moieties, including parallel SAR (structure activity relationship) between agonist and antagonist series.9,29,30 Also, the exocyclic amino group common to both chemical classes has been suggested to form a H-bond with Asn253. The importance of this critical residue for both agonist and antagonist binding has been supported with mutagenesis.
The overlay of AR agonists and antagonists has long been a concern in light of the structural similarities of nucleobase and core heterocycles of the antagonists. Three older hypotheses treated the mode of overlay of purine agonists and purine (e.g., xanthine) antagonists.31 In particular, it is possible to superimpose N1, N3, N7, and N9 nitrogen atoms of the purine ring of adenosine with the corresponding nitrogen atoms of theophylline. However, in this case the positive electrostatic potential of the exocyclic amino group of adenosine corresponds to the negative electrostatic potential of the carbonyl oxygen atom of theophylline. In view of this, van Galen et al. suggested the overlap of N1, N3, N7, and N9 atoms of the adenine ring with the C2, C6, N9, and N7 atoms, respectively, of the xanthine moiety.32 Since the binding of N6-substituted adenosine derivatives correlates well with the binding of C8-substituted xanthines it was proposed that substituents at these positions can occupy the same place inside the receptor.29 The superimposition of the N1, N3, N9 atoms of the adenine ring with the N9, N3, and N1 atoms of xanthine ring was suggested.33 Since we have not studied xanthine binding to the A2A AR in this study, we cannot distinguish these possibilities. However, our present study predicts that correspondence model in which N is mapped on N applies to antagonists such as 1, in which there is an exocyclic NH. This correspondence of amino groups was previously used to develop A3 selective antagonists through acylation of the amino group in the heterocyclic series.36 Other parallels in SAR with nucleoside agonists, i.e., at the 2 position, have been found in the series of adenine antagonists, which also have an exocyclic NH.34
The inaccuracy and uncertainty of the modeling approach applies more to the docking to the A2A AR, rather than its protein structural prediction, at least in the TM region. This problem is somewhat alleviated if there is a likely anchor point, such as an electrostatic bridge for the biogenic amine and nucleotide (P2Y) receptors.2,25 Moreover, the prediction of the correct binding mode of an antagonist seems to be a more difficult goal than modeling of the agonist-bound receptor. Small molecule agonists of a given Family A GPCR are perhaps more likely than antagonists to occupy a conserved region of the receptor, because of a commonality in the receptor activation mechanism.35 In theory, in order to block agonist action, antagonists need to overlap with the docked position of agonists only partially. Thus, the hypothesis that antagonists display a more diverse range of docking modes than agonists is reasonable. Although the orientation of the 1 antagonist in the X-ray structure was shifted toward the extracellular region in comparison to the previously obtained orientation of a distantly related antagonist 2, many the same amino acid residues predicted in the modeling are involved. Some receptor antagonists, such as the secondary propanolamines of the β2AR and the truncated nucleosides of the A3 AR21 are more likely to occupy the same binding mode as agonists due to structural similarity. Thus, in comparison to 1 the position of carazolol in the β2 adrenergic receptor is closer to the position of retinal in rhodopsin and, as predicted here, the position of adenosine in the A2A AR. Finally, in this modeling study we have predicted that the key heterocyclic pharmacophore of the antagonist 1 closely overlays the purine moiety of adenosine agonists.
The problem of flexibility of the loops has been discussed in the literature and, in general, differences in configuration between the computational models and the experimental structure were expected. For this reason, in some models of GPCRs the extracellular and intracellular loops were not constructed at all. The EL2 of the A2A AR crystallographic structure is more similar to the β-adrenergic receptors than to rhodopsin. Not surprisingly, it was possible to obtain a similar docked mode of 1 using as a template the β2 adrenergic receptor crystallographic structure, which has a more extended EL2 structure. Basing the previous models of ARs on rhodopsin possibly introduced some distortion of the docked ligand structure due to steric crowding in the region of EL2. This distortion of the homology model is particularly relevant to ligands such as 1, which have an extended chain. The presence of chains in functionalized congeners of both AR agonist and antagonist ligands anticipated that part of these extended ligands would protrude into the extracellular region.29 Also, Duong et al.30 demonstrated that mutagenesis of a charged residue in EL2 of another AR subtype, the A3 AR, selectively increased affinity in C2 derivatives of adenosine, in which a long chain was modified in a complementary fashion. The C2 substituent of adenosine agonists would correspond to the 2-phenylethylamino chain of 1, according to our agonist docking hypothesis.
To summarize, a comparison of the recently published crystal structure of the A2A AR and its structure predicted with the molecular modeling revealed similarities not only in the general configuration of TMs, but also in the position and orientation of side chains of amino acid residues involved in ligand recognition. Analysis of the ligand-receptor interactions observed in the X-ray structure fully confirmed much of the previous data obtained with rhodopsin-based homology modeling and site-directed mutagenesis. Molecular docking of AR agonists to the A2A AR crystal structure provided a similar agonist binding mode to that previously predicted based on the A2A AR homology models. The utilization of the β2AR-based model of the A2A AR receptor for the docking of 1 provided an accurate approximation of the experimentally observed binding mode of this A2A AR antagonist. Thus, our findings have demonstrated that homology modeling of GPCRs remains a useful technique in probing the structure of the protein and predicting modes of ligand docking. Further insights in ligand binding will be gained from crystallization of the ARs with chemically diverse ligands, including adenine derivatives.
The MOE software was utilized to construct a homology model of the human A2A AR. The recent crystal structure of the β2AR was used as a template. The sequence alignment of the A2A AR and β2AR was performed with MOE program taking into account positions of residues highly conserved for all GPCRs. Only a part (Cys166 - Asn175) of EL2 was modeled. A total of 25 independent models were built, and the highest scoring intermediate model was chosen automatically as the final model. The geometry optimization was performed with the AMBER99 force field. The Reaction Field salvation model was used. The “Fine” option was selected for both intermediate and final models refinement. The energy minimization of intermediate models was terminated when the Root Mean Square (RMS) of the energy gradient falls below 1.0, for the final model the value of RMS gradient was set to 0.5 kcal·mol-1·Å-1.
The InducedFit docking procedure of the Schrödinger Macromodel package was utilized to dock adenosine and 4 to the X-ray structure of the A2A AR. The receptor grid generation was performed for a box with the center in the centroid of 1. The box size was determined automatically. The default values of 0.50 were used for the receptor and ligand van der Waals scaling during initial Glide docking. The number of poses to remain after initial docking was set to 30. The residues located within 5.0Å of the docked ligand were subjected to refinement with the Prime. The optimization of side chains was enabled. The structures within 30.0 kcal·mol-1 of the best structure were redocked.
The A2A AR crystal structure was used for molecular docking of 1. Three different approaches were used: Schrodinger Glide XP and InducedFit docking, and MOE Tabu Seacrh docking. In all cases the size of the docking box was determined automatically. The center of the docking box had a center in the centorind of 1 inside the crystal structure of the A2A AR in its experimental position and conformation. For all, Glide XP, InducedFit docking, and MOE-dock, their default parameters were used.
The β2AR-based model of the A2A AR was used for molecular docking of 1. The automatic docking was performed with MOE-dock tool implemented in MOE software. A 3D docking box with the side of 24Å was defined with center placed in centroid of Thr88, Phe168, Asn181, Trp246, His250, Asn253, Ser277, and His278. The standard parameters of the Tabu Search (1000 steps/run, 10 attempts/step, 10 Tabu list length) protocol were applied. The default value of 25 docking runs was used. The calculations were performed in the MMFF94x force field.
To refine the orientation of the receptor side chains, the complex of the A2A AR with 1 docked was subjected to the MCMM calculations. The MCMM calculations were performed with the Schrödinger Macromodel software. 1 and all amino acid residues located within 4Å from the ligand were subjected to the MCMM search. A shell of constrained residues located within additional 2Å from 1 was used. The maximum number of steps was set to 100, and the energy window for saving structures was set to 100.0 kJ·mol-1. Water was used as an implicit solvent. Maximum of 500 interactions of Polak - Ribier conjugate gradient (PRCG) minimization method with a convergence threshold of 0.05 kJ·mol-1·Å-1 was used.
This research was supported by the Intramural Research Program of the NIH, National Institute of Diabetes and Digestive and Kidney Diseases. We thank Dr. Stefano Costanzi (NIDDK) and Dr. Soo-Kyung Kim (California Inst. of Technology) for helpful discussion.