Abstract
G protein–coupled receptors exist in a whole spectrum of conformations that are stabilized by the binding of ligands with different efficacy or intracellular effector proteins. Here, we investigate whether three-dimensional structures of receptor conformations in different states of activation can be used to enrich ligands with agonist behavior in prospective docking calculations. We focused on the β2-adrenergic receptor, as it is currently the receptor with the highest number of active-state crystal structures. Comparative docking calculations to distinct conformations of the receptor were used for the in silico prediction of ligands with agonist efficacy. The pharmacology of molecules selected based on these predictions was characterized experimentally, resulting in a hit rate of 37% ligands, all of which were agonists. The ligands furthermore contain a pyrazole moiety that has previously not been described for β2-adrenergic receptor ligands, and one of them shows an intrinsic efficacy comparable to salbutamol.
SIGNIFICANCE STATEMENT Structure-based ligand design for G protein–coupled receptors crucially depends on receptor conformation and, hence, their activation state. We explored the influence of using multiple active-conformation X-ray structures on the hit rate of docking calculations to find novel agonists, and how to predict the most fruitful strategy to apply. The results suggest that aggregating the ranks of molecules across docking calculations to more than one active-state structure exclusively yields agonists.
Introduction
G protein–coupled receptors (GPCRs) are highly flexible signal transduction proteins that are embedded in the outer membranes of eucaryotic cells. Structurally, they are relatively loose bundles of seven transmembrane (TM) domains that exist in a spectrum of different conformations (Latorraca et al., 2017). Binding of molecules that enhance signaling, i.e., agonists, to the orthosteric binding pocket (which is located within the transmembrane core for family A GPCRs) shifts the equilibrium between the various conformations such that binding of intracellular effector proteins (e.g., G proteins or arrestins) becomes more likely than in the basal state of the receptor. The receptor thus becomes “activated” and—depending on the effector—a signaling cascade is induced.
At a structural level, the various existing X-ray structures have shown that active conformations are characterized by an outward movement of the intracellular halves of TM V and VI relative to the basal state. This is accompanied by a comparatively small contraction of the orthosteric binding pocket (Rasmussen et al., 2011a). However, although binding of an agonist makes receptor-effector coupling more likely, there is also evidence that G protein binding increases the affinity of the receptor for agonists (De Lean et al., 1980). Recently, a “reverse pharmacology” study demonstrated that extracellular agonists preferentially bind to active receptor conformations (Pardon et al., 2018). In this study, the β2-adrenoceptor (β2AR) was “locked” in an active or inactive conformation through fusion with a G protein–mimicking (Nb80) or an irrelevant (Nbirr) nanobody, respectively. Agonists bound with higher affinity to the “active” β2AR-Nb80 fusion than to the “inactive” β2AR-Nbirr receptor fusion. For inverse agonists, the opposite binding preference was detected. Furthermore, those ligands that displayed higher affinity for the “active” conformation were confirmed to act as agonists in subsequent cellular assays. Thus, it was possible to predict ligand efficacy from the relative affinity to the two receptor-nanobody constructs, which correspond to different conformational ensembles.
Here, we reproduced the nanobody study in silico. It has been shown previously that such comparative docking calculations can increase the chance of finding molecules with agonist properties in large compound libraries (Weiss et al., 2013). At the same time, we have demonstrated that already small changes in receptor structure lead to completely different ligand sets (Kolb et al., 2012). The question thus remained whether the additional structures of the β2AR in active conformations that have become available in the past decade would contain novel information and lead to different agonists. Moreover, we investigated what the best way of aggregating the results of multiple docking calculations is. All molecules emerging from the screen were characterized pharmacologically. To investigate the selectivity of the ligands’ interactions, we performed all experiments on both the β2AR and the β1AR.
Here, we present the results of this study, which not only had a remarkable overall ligand hit rate of 37%, but all hits were agonists. Of note, several of the newly discovered ligands feature a pyrazole moiety that has never been described for ligands of the β2AR before. Somewhat in contrast to earlier work (Weiss et al., 2013) and the nanobody-based assay (Pardon et al., 2018), we found that the comparison of ligand ranks between active and inactive conformations did not lead to the highest enrichment of agonists. We therefore analyzed the retrospective enrichments of known β2AR ligands in the various conformations to determine whether this behavior can be predicted, thus allowing for appropriate choices in prospective screens.
Materials and Methods
Primary Docking Screen.
The crystal structures with Protein Data Bank (PDB) identifiers (IDs) 3SN6 (Rasmussen et al., 2011b), 4LDL (Ring et al., 2013), and 3NY9 (Wacker et al., 2010) were prepared for docking by protonation and subsequent minimization of all hydrogen atoms with CHARMM and the CHARMm22 force field (Momany and Rone, 1992). In addition, the binding pocket of the 3SN6 structure was relaxed by energy-minimizing residues W1093.28 [numbers in superscript are according to the Ballesteros-Weinstein enumeration scheme (Ballesteros and Weinstein, 1995)], T1103.29, D1133.32, V1173.36, F19345.52, N2936.55, and N3127.39 in the presence of the ligand visible in the crystal structure (BI-167107) to ameliorate clashes between ligand and protein, while ensuring that the polar hydrogen atoms of S2075.46 and N2936.55 were pointing toward the bound ligand, thus enabling the residues to act as hydrogen bond donors. The spheres that are used in DOCK to translate and rotate molecules into the binding pocket were moved to optimize the enrichment of a set of 124 enantiomers of known β2AR ligands over decoys generated with DUD-E (Mysinger et al., 2012). After this optimization step, the leads now subset of the ZINC12 library (3,687,621 molecules) (Irwin et al., 2012) was docked to both crystal structures of the β2AR in active conformations (PDB IDs 3SN6 and 4LDL; β2ARactive) using DOCK3.6 (Kuntz et al., 1982; Meng et al., 1992; Shoichet and Kuntz, 1993; Shoichet et al., 1999; Mysinger and Shoichet, 2010). In addition to the individual ranked lists (lists 1 and 3; Fig. 1) determined for each docking calculation, we also generated another ranked list from the two β2ARactive dockings to enrich molecules that ranked highly in the docking calculations to both structures (“dual reranking,” list 2 in Fig. 1). Furthermore, the rank-ordered list of the docking to the β2ARactive structure 3SN6 was compared with the list of docking to a structure of the β2AR in an inactive conformation (PDB ID 3NY9; β2ARinactive) to identify molecules ranked highly in an active conformation, but poorly in the inactive conformation used in this study (“selective reranking,” list 4 in Fig. 1). Evaluation of the top 500 molecule poses in each of the four lists (4LDL-based, 3SN6-based, dual reranking, selective reranking) was done visually to manually remove molecules with artificially inflated scores due to any one of the known deficiencies of scoring functions (Kolb et al., 2012). For molecules that showed favorable poses but were too big to fit the binding pocket, smaller derivatives were searched in ZINC12 and docked to 3SN6. A docking calculation to the β2ARinactive structure PDB ID 2RH1 (Cherezov et al., 2007) was performed to compare with our previous docking study (Kolb et al., 2009). This protein structure was prepared and the docking calculation conducted in the same manner as described earlier. The molecules with the highest numbers of favorable interactions in their poses from all considered ranking lists were pooled and the final list of molecules selected (for IDs and vendors, see Supplemental Table 1).
Reranking.
Two ranking lists containing the same molecules docked to two different receptor structures were compared with each other to find the highest-ranking molecules from two lists (dual reranking) or those that ranked highly in one but poorly in the other (selective reranking). Dual reranking was done as described previously (Schmidt et al., 2015). In brief, a relative rank (Rrel) for each molecule is calculated for both ranking lists using the rank r of each molecule and the total number of molecules (m) in the list (Eq 1):(1)where i is the indicator of each docking. With these relative ranks, a new score (D) is calculated for each molecule (Eq 2):(2)Molecules are then ranked according to the new score D. For further explanation on calculations see (Schmidt et al., 2015). Dual reranking was applied to compare the docking calculations to the two structures in an active conformation, i.e., PDB IDs 3SN6 and 4LDL.
For the selective reranking, the new score S was simply calculated as the ratio of the two ranks in the individual dockings (Eq 3):(3)where is the molecule rank in the list of the docking calculations against 3SN6, and is the molecule rank in the list of the docking calculations against 3NY9.
Similarity-Based Search for Derivatives and Structure-Based Evaluation.
Molecules similar to ligands causing the highest accumulation of cAMP in the primary cellular assay were searched in a library of 5,626,190 molecules from ZINC15 (Sterling and Irwin, 2015) using the extended-connectivity fingerprint 4 (ECFP4) and a Tanimoto coefficient cutoff of . The resulting set of molecules was then docked to the structure with PDB ID 3SN6 using DOCK3.6, and molecules were selected to be tested experimentally after visual inspection of the molecule poses as described for the primary screen (IDs and vendors in Supplemental Table 1).
Evaluation of Compound Novelty.
All molecules that have been tested in any way against the β2AR were downloaded from ChEMBL (Gaulton et al., 2012) and filtered for active molecules (2374 unique entries). All molecules that were tested within the present screen were then compared with this ChEMBL-derived data set by their Tanimoto similarity using ECFP4 fingerprints (Wawer and Bajorath, 2010). To evaluate the similarity of our ligands to any adrenergic receptor, a second data set was prepared by downloading the bioactivity data of all molecules targeting any one of the human adrenergic receptors from ChEMBL and filtering this data set for actives (7396 unique entries). The similarity evaluation was again done using ECFP4 fingerprints and the Tanimoto coefficient.
Cell Culture.
Stable cell lines of CHO cells expressing either the β2AR or β1AR and the cAMP response element (CRE)–secreted placental alkaline phosphatase (SPAP) reporter gene were grown in Dulbecco’s modified Eagle’s medium nutrient mix F12 containing 10% fetal calf serum and 2 mM l-glutamine at 37°C in a humidified 5% CO2:95% air atmosphere.
3H-CGP12177 Whole-Cell Binding.
Cells were plated to white 96-well plates and grown to confluence overnight. Medium was then removed from the cells and replaced by 100 μl of serum-free medium (sfm; Dulbecco’s modified Eagle’s medium nutrient mix F12 containing 2 mM l-glutamine) or compound (at twice the final concentration in sfm), followed immediately by the addition of 100 μl of 3H-CGP12177 in sfm (1:2 dilution in wells). The cells were then incubated at 37°C in a humidified 5% CO2:95% air atmosphere for 2 hours. After 2 hours, everything was removed from all wells, and the cells were washed twice with 200 μl of 4°C phosphate-buffered saline. Microscint 20 (100 μl) was added to each well, and the plates were left at room temperature in the dark for several hours before being counted on a TopCount (Packard BioScience/ PerkinElmer) (Baker, 2005b). Propranolol (10 μM) was used to define nonspecific binding in all plates, and the final concentration of 3H-CGP12177 was 0.5–0.9 nM. All compounds were examined in seven-point concentration-response curves, with each condition repeated in triplicate in each experiment. A sigmoidal curve was fitted to the data using GraphPad Prism 7 (GraphPad Software) and an IC50 determined from the following Eq 4:(4)where [A] is the concentration of the competing ligand, and IC50 is the concentration at which half of the specific binding of 3H-CGP12177 has been inhibited.
From the IC50 value and the known concentration of radioligand, [3H-CGP12177], a KD value (concentration at which half the receptors are bound) was calculated using the Cheng-Prusoff equation (Eq 5):(5)The KD of 3H-CGP12177 in these cells was 0.42 nM for the β1AR and 0.17 nM for the β2AR (Baker, 2005b).
CRE-SPAP Production.
Cells were plated to clear 96-well plates and grown to confluence overnight. Medium was then removed from all wells and replaced with 100 μl of sfm for 24 hours (i.e., cells were serum starved). The following day, medium was again removed and replaced by 100 μl of fresh sfm. Where used, antagonists diluted in sfm were used instead of pure sfm. Compound (10 μl diluted in sfm) was then added to the wells and the cells incubated for 5 hours at 37°C in a humidified 5% CO2:95% air atmosphere. After 5 hours, everything was removed from all wells and replaced by 40 μl of sfm. Cells were then incubated for 1 hour at 37°C in a humidified 5% CO2:95% air atmosphere before being placed in a 65°C oven for 30 minutes. The plates were cooled to room temperature, and 100 μl of 5 mM para-nitrophenylphosphate in Diethanolamine buffer (1 M diethanolamine, 0.28 M NaCl, 0.5 mM MgCl2, pH 9.85) was added to each well. Once the yellow color had developed, plates were read on an MRX plate reader (Dyntaec Scientific Labs) through a 405-nm filter (Baker et al., 2014). Isoprenaline (10 μM) was used as a positive control in all plates. All compounds were examined as a seven-point concentration-response curve with triplicate condition for each concentration of ligand. A sigmoidal concentration-response curve was fitted to the data using GraphPad Prism 7 (eq. 6):(6)where Emax is the maximum response, [A] is the agonist concentration, and EC50 is the concentration of agonist that produces 50% of the maximal response.)
To prove that the agonist responses were indeed occurring via the transfected βAR, the affinity (log KD value) of CGP20712A or ICI118551 was calculated from the rightward shift of the agonist concentration responses in the presence of a fixed concentration of antagonist using the Gaddum equation (Eq 7):(7)where DR (dose ratio) is the ratio of the agonist concentration required to stimulate an identical response in the presence and absence of the fixed concentration of antagonist [B].
Efficacy Ratios.
As the affinity (KD) and potency of the agonist response (EC50) for each ligand were determined in the same cell lines, and ligands were all examined in parallel experiments, an indicator of intrinsic efficacy (ability of a compound to stimulate a response) can be determined from the efficacy ratio at each receptor. This therefore takes into account the affinity of the ligand, and although the score cannot be compared across cell lines, it can be used to rank ligands in order of intrinsic efficacy at each receptor. Thus, cimaterol, with a β2AR affinity (KD) of 81 nM and a β2AR EC50 value of 0.21 nM, is very efficacious, as it hardly needs to occupy any receptors to stimulate a maximum response (efficacy ratio of 386). At the same receptor, salmeterol, although of higher affinity (KD 0.81 nM, EC50 0.012 nM), once bound, is less efficacious (efficacy ratio 67). CGP12177 is a partial agonist (as can be seen from the 37% stimulation in relation to isoprenaline). As a partial agonist, it needs to occupy all of the available receptors to stimulate its maximum 37% response, and the KD value (0.28 nM) and EC50 value (0.19 nM) are similar, giving a low efficacy ratio of 1.47.
Results
Primary Docking Screen.
The leads now subset of the ZINC12 library (Irwin et al., 2012), consisting of 3,687,621 molecules, was docked to two β2ARactive (PDB IDs 3SN6 and 4LDL). As described in the Materials and Methods, four different schemes were used to rank the molecules from the docking calculations (detailed in Fig. 1). Molecules were ranked individually in both conformations [list 1 (red) and list 3 (yellow) in Fig. 1] and reranked to identify consistently favorably ranked molecules [list 2 (orange) in Fig. 1] as well as those that had a large difference in ranks when compared with the docking to a β2ARinactive [PDB ID 3NY9; list 4 (green) in Fig. 1]. In each of these four ranking schemes, the top 500 molecules were visually inspected to remove molecules that are ranked highly because of insufficient penalties for unfavorable interactions due to the known deficiencies of current force fields. From the first three rankings [docking calculations to the β2ARactive; lists 1, 2, and 3 (red, orange, and yellow, respectively) in Fig. 1], 18 molecules were selected and purchased for testing. Five additional molecules showed favorable poses in the orthosteric pocket but had several parts protruding beyond it, and these parts consequently did not entertain favorable interactions. Smaller derivatives of these molecules were therefore searched in ZINC12, with the aim to keep the key interacting parts constant and only change the bulky noninteracting moieties. The nine resulting derivative molecules were docked to the 3SN6 structure, poses were evaluated visually, and three of the smaller derivatives were selected to be tested (Fig. 1, pink box). The parent molecules of these three derivatives emerged from the dual reranking (3SN6 vs. 4LDL; one derivative) and the selective reranking (3SN6 vs. 3NY9; two derivatives), respectively. Interestingly, ranking 4, contrasting the results from a docking calculation to a β2ARinactive and a β2ARactive, did not yield any molecules that were deemed worth testing after visual inspection.
One additional molecule was selected from a separate docking calculation to the β2ARinactive structure with PDB ID 2RH1 [list 6 (darker blue) in Fig. 1], which we had done to compare with our earlier study (Kolb et al., 2009). Molecule 8 was chosen because of its pose in the β2ARinactive, which featured interactions with D1133.32. Interestingly, this molecule did not feature the same docking pose in the various β2ARactive docking calculations because the bulky substituent of 8 forced a flipped pose in the pocket. Nonetheless, we decided to evaluate it as a possible agonist because of its significantly smaller size, which is generally considered to be one of the hallmarks of β2AR agonists.
Similarity-Based Search for Derivatives of the Initial Hits and Structure-Based Evaluation.
An initial evaluation of the compounds in a cAMP-accumulation-based assay in HEK293T cells (Supplemental Fig. 1; Supplemental Results) revealed two ligands with a higher-than-average percentage of activation, 1 and 2. These two molecules are close analogs of each other. A similarity search with these molecules as query retrieved 62 additional molecules. Since the parent molecules had the highest ranks in the docking calculation using β2ARactive structure 3SN6, the derivatives were docked to the same structure to evaluate their steric fit. A complementary docking to the β2ARinactive structure 3NY9 was used to compare the docking poses and validate them. Based on their favorable poses, five molecules were selected for detailed pharmacological analysis (salmon-colored box in Fig. 1).
Pharmacological Characterization.
A total of 27 molecules (22 from the primary docking screen and 5 from the similarity-based search) were examined in CHO cells stably expressing either the human β2AR (CHO-β2) or the human β1AR (CHO-β1).
Ligand affinity was determined from a whole-cell 3H-CGP12177 radioligand displacement assay. CGP20712A, a known β1AR antagonist, bound to the β2AR with much lower affinity than the β1AR (log KD −5.8 and −8.6, respectively), whereas ICI118551 bound to the β2AR with higher affinity (log KD −9.3 at the β2AR, −6.8 at the β1AR), thus demonstrating the presence of the two receptors in the respective cell lines (Table 1). Of the 27 ligands examined, 10 were found to have measurable affinity, eight from the primary screen and two from the similarity-based search. A further six molecules had very low affinity such that the calculation of KD values was not possible. The compound with the highest affinity was 1, with an affinity for both the β2AR and the β1AR of around 520 nM (Fig. 2; Table 1).
The ability of compounds to stimulate a functional response was then examined in the same cells using the CRE-SPAP reporter stably expressed in both cell lines. The benefit of this system is that it is an amplified downstream readout, thus maximizing the chances of detecting any agonist properties (Baker et al., 2004). Given that partial agonists may be amplified to appear more like full agonists in this system, several well known partial agonists were included in the study for comparative purposes. Cimaterol, a relatively nonselective βAR agonist, stimulated a full agonist response with respect to isoprenaline in both cell lines (log EC50 −9.7 at the β2AR, log EC50 −8.6 at the β1AR; Table 2). As expected, denopamine, a known β1AR-selective agonist, stimulated a full agonist response that was more potent in the β1AR cell line, whereas salbutamol, a known β2AR-selective agonist, was more potent in the β2AR cell line (Baker, 2010) (Table 2). As anticipated, the 17 compounds with no affinity in the binding assay were also not able to stimulate any significant agonist response in either cell line. These compounds, therefore, did not interact with either the β2AR or β1AR.
Of the 10 compounds that did have measureable affinity, all showed agonist efficacy (Fig. 2; Table 2), and 1 was the most potent ligand at both receptors. However, it was important to be sure that these agonist responses were indeed occurring via the transfected receptors. Dose responses were therefore examined in the parent CHO-CRE-SPAP cells that contained the stably transfected reporter gene but no transfected receptor. No agonist responses were elicited by any of the 10 novel agonist compounds, nor by cimaterol, denopamine, salbutamol, salmeterol, or CGP12177, in CHO-CRE-SPAP cells (n = 3 for each ligand).
Finally, the agonist responses were inhibited by the selective antagonists CGP20712A and ICI118551 in their respective cell lines. In the CHO-β2 cells, agonist responses to cimaterol, salbutamol, salmeterol, denopamine, and CGP12177 inhibited by ICI118551 yielded very similar log KD values [log KD approximately −9.6, similar to the log KD value of ICI118551 (−9.3) determined from the binding assay; Tables 1 and 2]. Similar high-affinity values for ICI118551 were obtained from its ability to inhibit the agonist responses to the 10 novel compounds (Table 2). This suggests that all of these ligands were acting via the transfected β2AR, its orthosteric pocket, and the same receptor conformation. Compound 8 stimulated a partial agonist response in the CHO-β2 cells, and although it was inhibited by ICI118551, a rightward shift and flattening of the curve was observed. As such, the dose-response curve did not reach the same maximum in the presence of ICI118551, and thus, a log KD could not be calculated using the Gaddum equation.
The β1AR exists in at least two conformations: a high-affinity catecholamine conformation (where responses are readily inhibited by antagonists) and a low-affinity secondary conformation (where inhibition of agonist responses requires much higher concentrations of antagonist), where CGP12177 is an agonist (Konkar et al., 2000; Granneman, 2001). In the CHO-β1 cells, agonist responses to cimaterol, salbutamol, salmeterol, and denopamine were inhibited by CGP20712A to yield high-affinity values for CGP20712A, again similar to those seen from the binding assay (Tables 1 and 2). The agonist response to CGP12177, however, required significantly higher concentrations of CGP20712A to achieve a rightward shift and thus yielded a log KD value for CGP20712A much lower than that of the other compounds (Table 2). This is because the agonist response to CGP12177 is occurring via the secondary conformation (Baker, 2005a). The agonist responses to the novel compounds were all inhibited by CGP20712A to yield values very similar to those for the literature compounds and similar to those obtained in the binding assay. All of these compounds are therefore exerting their agonist action through the orthosteric pocket in the primary catecholamine conformation of the β1AR—perhaps not surprisingly, as this is the conformation most similar to that of the β2AR.
Hit Rates.
In total, we tested 27 molecules, 22 emanating from the primary docking screen and 5 from the similarity-based search. The competition binding assay confirmed 10 new ligands of these 27 tested molecules, corresponding to a hit rate of 37% (see Supplemental Table 2; Table 1). All of the 10 compounds showed agonistic activity, i.e., again a 37% overall hit rate or 100% of the discovered ligands (Table 2). Two of these novel agonists resulted from the similarity-based search (i.e., two of five molecules, a 40% hit rate). Only 1 of the 22 molecules from the primary screen was retrieved from the docking calculation to the inactive receptor 2RH1 (5%). The majority of discovered agonists, seven, originated from the dual reranking (7/22 molecules; 32%). The overall hit rate of 37% is at the upper end compared with other docking studies using the β2AR (Sabio et al., 2008; Kolb et al., 2009; Weiss et al., 2013; Schmidt et al., 2017; Chevillard et al., 2019) or related aminergic GPCRs (Carlsson et al., 2011; Kruse et al., 2013).
Correlation between Retrospective Ligand Enrichments and Hit Rates.
The rather large variation in the number of ligands retrieved from the five docking and reranking schemes made us wonder whether this could have been predicted. To that end, we applied the ranking schemes to the set of 124 known active enantiomers and their decoys that had originally been used to optimize the docking procedure. For each of the resulting ranking lists, we calculated receiver operating characteristic plots. These plots allow assessment of the enrichment of agonists and antagonists in the top ranks versus the decoy molecules in each of the docking schemes. Agonists are enriched to a similar extent in each of the docking calculations to a β2ARactive as well as in both reranked lists, and better than antagonists (Fig. 3; Supplemental Fig. 2). Furthermore, selective reranking leads to below-random enrichment of antagonists, albeit with an early enrichment similar to the docking calculations to the two β2ARactive structures.
Interestingly, although overall enrichment of agonists and antagonists in the ranking lists of the dual reranking does not differ noticeably from the other three ranking schemes, antagonists are only enriched after roughly 6% of the decoys have been found. For a set of 3.6 million molecules (such as the one we used here), this would mean that no antagonists appear in the top 200,000 molecules of the dual reranking list, and therefore, all hits from the top 200,000 should be agonists. Although this is certainly overestimating the results, this reasoning is in line with the hit rate achieved during the large library docking screen.
Ligand Novelty.
Novelty of the discovered ligands was evaluated by comparison with all compounds binding to or acting at the β2AR or any adrenergic receptor according to ChEMBL (Gaulton et al., 2012). Although the most active novel agonist discovered during this screen, compound 1, is similar to known βAR agonists in that it contains a β-hydroxy-amine, it also contains a pyrazole moiety, a feature it shares with six of the novel agonists discovered in the primary docking screen. Molecules featuring this moiety have so far never been described as ligands of the β2AR. Because of this unprecedented ring, the ligands have relatively low Tanimoto similarity values compared with known β2AR ligands in ChEMBL (Gaulton et al., 2012). The highest similarity value was 0.47 (8 to CHEMBL599896), and only four molecules showed a similarity (see Supplemental Table 3). None of the most similar compounds contained the same structural motifs as the query molecules.
The more general comparison against ligands of any adrenergic receptor revealed a similar picture. The highest similarity was 0.49 (low-affinity 14 to CHEMBL15303) or 0.47 (8 to CHEMBL599896) within this search. Matches with a similarity were found for only seven of the molecules (see Supplemental Table 4). Casting an even wider net, all compounds containing the basic substructure of phenyl, β-hydroxy-amine, and pyrazole where retrieved from ChEMBL. None of them was associated with bioactivity data for the β2AR. We note that compound 7 has previously been described as a weak κ-opioid receptor antagonist (pK) (Zheng et al., 2017) but never as a β2AR ligand (Wawer and Bajorath, 2010).
Predicted Binding Poses.
In the predicted docking pose, the β-hydroxy-amine of 1 interacts with D1133.32 and N3127.39 (Fig. 4A). The fluorine substituent in ortho position on the aromatic moiety might contribute interactions with the protein (potentially N2936.55 or Y3087.35 in the β2AR) that explain its affinity and agonistic activity.
Compound 8 stands out from the other tested molecules since it is the only one that was discovered from a docking calculation to a β2ARinactive. This molecule has a bulkier structure than the other ligands, with an aliphatic chain and an additional aromatic moiety attached to the more standard benzene ring following the β-hydroxy-amine. Of note, other than the weak H-bond acceptor ether oxygen, there is no possibility to interact with any of the Ser residues in TM V. In the larger binding pocket of the β2ARinactive, compound 8 is predicted to adopt a pose with the alkyl-aromatic moiety pointing toward extracellular loop 2 (Fig. 4B). Unexpectedly, this molecule pose could not be reproduced in a docking calculation to the narrower β2ARactive binding pocket. It is interesting to mention that this molecule has only slightly lower affinity (1 μM at the β2AR) than the highest-affinity agonist 1. One potential reason for its counterintuitive binding preference might be that it is only a partial agonist at the β2AR (with a response of 35% that of isoprenaline and very similar to that of CGP12177 at 37%). As expected for a partial agonist, the efficacy ratio is low. Despite the large difference in affinity of 8 and CGP12177 (Table 1), once bound, both compounds have a similar ability to activate the receptor (i.e., have similar intrinsic efficacy, as demonstrated by the efficacy ratio). Given that the CRE-SPAP response is an amplified system (Baker et al., 2004), 8 can be regarded, just as CGP12177, as a weak partial agonist in the overall ranking of efficacy. Molecule poses of all molecules from the docking calculations to at least one of the used structures can be found in the Supplemental Materials section.
Discussion
This study investigated whether comparative docking to active (and inactive) conformations of the β2AR could preferentially detect novel compounds with βAR agonist activity, rather than just antagonist compounds. Moreover, we tried to determine the impact of using multiple β2ARactive structures and how to best aggregate the individual rankings.
Molecule Selection Based on Docking Ranks.
Since different approaches were used to retrieve potential ligands with agonistic action from the ranked lists of each docking calculation in the present work, we can evaluate each strategy in terms of the number of agonists retrieved. The initial hypothesis was that agonists would be ranked higher in docking calculations against the β2ARactive than the β2ARinactive. This assumption held true in retrospective studies. Somewhat unexpectedly, however, the selective reranking carried out in this work did not yield molecules that we considered worth testing in the prospective docking. Despite this low abundance of favorably interacting compounds, we still chose two smaller derivatives of compounds emerging from this ranking. All additional compounds tested here were selected either from the individual ranked lists of the docking calculations to the two β2ARactive (PDB IDs 3SN6 and 4LDL) or from the dual reranking taking into account the ranked lists of both of these docking calculations (cf. eq. 2).
In our hands, this dual reranking resulted in not only a very high hit rate but also a high number of favorably interacting molecule poses in general. Our results therefore suggest that the most fruitful approach to retrieve agonists is docking calculations to more than one structure in an active conformation, aggregating the lists to identify those molecules that are ranked highly in both docking calculations. Turning this around, and somewhat in contrast to our initial assumptions, a high agonist hit rate does not require the comparison of docking calculations against β2ARactive and β2ARinactive structures.
The comparison with the retrospective docking calculations of known β2AR ligands suggests that our prospective findings are echoed in the relative enrichments of agonists and antagonists over computer-matched decoys, and that such retrospective calculations could be used to indicate a priori which reranking scheme might work best.
Closer evaluation of the molecule poses in the β2ARactive and the β2ARinactive shows that orientations are similar, but that the molecule poses in the β2ARactive result, on average, in a higher number of favorable interactions between receptor and molecule. We suggest that this fact can be used as an additional criterion during pose evaluation: While poses of agonists obtained from docking calculations to the β2ARactive should be reproducible in terms of general orientation in docking calculations to the β2ARinactive (since the bigger pocket of the β2ARinactive should leave enough space), the interactions can allow differentiation (see Fig. 5 for an overlay of all four crystal structures used).
We note that one agonist with significantly lower efficacy (8) from the primary screen was chosen from the ranked list of a docking calculation to the inactive conformation structure 2RH1. The fact that 8 does not find a favorable pose in the β2ARactive while it is a partial agonist leads to the conclusion that, in certain cases, the discovery of bulky ligands from a docking calculation to the β2ARactive might be hampered by the rigid conformation used during the calculation. A complementary docking calculation to the β2ARinactive can be helpful here but also includes the risk of incorrect predictions, because statistically, a favorable rank in the β2ARinactive, but not in the β2ARactive, would be taken to indicate an antagonist or inverse agonist. That the discovery of this agonist from a β2ARinactive structure is an exception is supported by the fact that none of the ligands discovered in our earlier studies of docking to a β2ARinactive were agonists (Kolb et al., 2009; Schmidt et al., 2017).
Pharmacological Evaluation.
Of the 27 compounds chosen from the docking screens for pharmacological evaluation, 10 (37%) were found to be agonists of the β2AR, and most were also of the β1AR, which we used to evaluate the selectivity of the compounds. Importantly, all of these responses were demonstrated to occur at the β2AR and β1AR, respectively, because all responses were inhibited by selective antagonists, and there were no responses in the parent CHO-CRE-SPAP cells without the transfected receptors. As for many agonists, the affinity of these compounds for the β2AR and β1AR was relatively low (Baker, 2010). However, for all compounds that were found to have measurable binding, we found agonist activity, too.
Although it is easy to rank ligands in order of the direct measure of affinity, establishing a ligand’s intrinsic efficacy is substantially harder, as there is no single direct pharmacological measure of intrinsic efficacy. One simple way to be able to compare the intrinsic efficacy of compounds is to use an efficacy ratio—comparing the EC50 with the KD for the same compound, when all compounds are examined in parallel under identical conditions (Baker, 2010). Here, it can be seen that cimaterol is the compound with the highest intrinsic efficacy at both receptors (Table 2). As expected, denopamine had higher intrinsic efficacy at the β1AR than the β2AR, and salbutamol at the β2AR than the β1AR (Baker, 2010). Of the novel agonists, 1 was the compound with the highest intrinsic efficacy at both the β1AR and β2AR, followed by 2. Of note, none of the derivatives from the analog search was more active than the parent compounds. By the same measure, compound 1, the most active molecule, has a similar intrinsic efficacy as salmeterol, a long-acting β2AR agonist widely used in the treatment of asthma and chronic obstructive pulmonary disease. Thus, our approach of docking compounds to different active structures to find novel compounds with agonism was able to discover novel compounds with sufficient agonism to be potentially clinically useful. However, interestingly, this approach detected compounds of medium efficacy, and not novel compounds with very high efficacy [such as catecholamines, fenoterol, or even cimaterol or salbutamol (Baker, 2010)].
Structure-Activity Relationship.
Several of the agonist compounds (1–7, 9, 11, 12) have rather similar chemical structures and sometimes differ only by one substituent. This allows structure-activity relationship conclusions to be drawn. Most compounds feature fluorine atoms at different positions on the aromatic moiety closest to the hydroxyl group of the β-hydroxy-amine. An additional fluorine atom in para position (6 compared with 5) had essentially no effect on affinity and intrinsic efficacy. In meta position, an additional fluorine atom (1 compared with 2), however, leads to increased binding affinity. As expected, the potency (EC50) is also increased, leading to an almost unchanged intrinsic efficacy. The strongest effect on the efficacy can be seen for the ortho position. Moving the ortho-fluorine atom (1) to the para position (12) had no effect on affinity but did reduce the efficacy ratio (especially for the β2AR; from 1.76 to 1.05; Table 2). Exchanging the fluorine for a chlorine in ortho position (3 compared with 4) results in almost unchanged affinity and pEC50 values. Of note, however, this exchange to a chlorine atom results in a higher percentage of maximum response for the β2AR (79% for 4 to 95% for 3), whereas it results in a lower percentage of maximum response for the β1AR (60% for 4 to 33% for 3). This increase can likely be explained by the ability of chlorine to form stronger halogen bonds than fluorine. This could result in a stabilization of a more active conformation of the receptor. Considering the difference in responses at the β2AR and the β1AR, an interaction of the substituent in ortho position with the protein seems likely. The most plausible candidate is Y3087.35 in TM VII (distances of 2.5–4 Å between the ortho substituents of the various ligands and the oxygen of the hydroxyl group of Y3087.35 in the respective energy-minimized docking poses), which is only present in the β2AR (it is a phenylalanine in the β1AR).
Finally, the comparison of 2 and 14 indicates that removing the hydroxyl group and introducing an ether in the alkylic chain between the aromatic moiety and secondary amine results in abolishment of binding.
In summary, multiconformation docking screens appear as a productive strategy to identify novel molecules with agonism similar to that of clinically used drugs. Although the consideration of the ranks of each molecule in the docking calculation against the β2ARinactive was not essential, comparison of the poses oftentimes helped in deciding whether a particular molecule was likely an agonist or not. The specific conformations elucidated in the crystallographically determined receptor structures appear to be sufficient to enrich agonists, i.e., that “function follows form.” Our study also demonstrated that information is added with new structures, and that a dual reranking considering ranking lists of docking calculations to β2ARactive leads to the largest number of agonists.
Last, but not least, we identified novel agonists for the β2AR featuring a previously undescribed pyrazole moiety and intrinsic efficacies on par with clinically used drugs. Several derivatives were explored during this project, revealing insights into their structure-activity relationship. Thus, docking calculations once again yielded a novel chemotype (Sabio et al., 2008; Kolb et al., 2009; Weiss et al., 2013; Schmidt et al., 2017) particularly remarkable for a target as well explored as the β2AR.
Although our findings should be readily transferable to other class A GPCRs, with similarly successful results, this approach could, at the moment, be more difficult when targeting GPCRs from other classes due to the smaller number of available crystal structures. Yet, with structures appearing at the rate that they currently do, this is an issue that will rapidly disappear in the future.
Acknowledgments
We thank Florent Chevillard for contributing the prepared structure of PDB 4LDL and PDB 2RH1 and the script for similarity searches based on the Tanimoto coefficient. Niek van Hilten contributed the method and the script for selective reranking. We thank Richard Proudman for technical help, especially with cell culture. We also acknowledge Felix Terwesten, Corey Taylor, and Tobias Hüfner for maintaining the computer network and software environment used in this study.
Authorship Contributions
Participated in research design: Scharf, Kolb.
Conducted experiments: Scharf, Baker.
Performed data analysis: Scharf, Baker, Kolb.
Wrote or contributed to the writing of the manuscript: Scharf, Bünemann, Baker, Kolb.
Footnotes
- Received June 26, 2019.
- Accepted October 8, 2019.
M.M.S. was supported by a travel grant to Nottingham by the “MID: Marburg International Doctorate–IPID4all” program of Philipps-University Marburg. This work was supported by the Heisenberg professorship [KO4095/4-1] and a project grant [Grant KO4095/3-1] of the German Research Foundation DFG granted to P.K. J.G.B. and P.K. participated in COST Action CM1207 “GLISTEN.”
↵This article has supplemental material available at molpharm.aspetjournals.org.
Abbreviations
- β2AR
- β
- CRE
- cAMP response element
- ECFP
- extended-connectivity fingerprint
- GPCR
- G protein–coupled receptor
- ID
- identifier
- PDB
- Protein Data Bank
- sfm
- serum-free medium
- SPAP
- secreted placental alkaline phosphatase
- TM
- transmembrane
- Copyright © 2019 by The American Society for Pharmacology and Experimental Therapeutics