GPCR Modulation of Thieno[2,3-b]pyridine Anti-Proliferative Agents

A panel of docking scaffolds was developed for the known molecular targets of the anticancer agents, thieno[2,3-b]pyridines, in order to glean insight into their mechanism of action. The reported targets are the copper-trafficking antioxidant 1 protein, tyrosyl DNA phosphodiesterase 1, the colchicine binding site in tubulin, adenosine A2A receptor, and, finally, phospholipase C-δ1. According to the panel, the A2A receptor showed the strongest binding, inferring it to be the most plausible target, closely followed by tubulin. To investigate whether the thieno[2,3-b]pyridines modulate G protein-coupled receptors (GPCRs) other than A2A, a screen against 168 GPCRs was conducted. According to the results, ligand 1 modulates five receptors in the low µM region, four as an antagonist; CRL-RAMP3 (IC50—11.9 µM), NPSR1B (IC50—1.0 µM), PRLHR (IC50—9.3 µM), and CXCR4 (IC50—6.9 µM). Finally, one agonist, GPRR35, was found (EC50 of 7.5 µM). Molecular modelling showed good binding to all of the receptors investigated; however, none of these surpass the A2A receptor. Furthermore, the newly-identified receptors are relatively modestly expressed in the cancer cell lines most affected by the thieno[2,3-b]pyridines, making them less likely to be the main targets of the mechanism of action for this compound class. Nevertheless, new modulators against GPCRs are of an interest as potential hits for further drug development.


Introduction
Reynisson et al. [1] have established that the class of thieno [2,3-b]pyridines (TPs) have potent anticancer activity against a variety of tumour cell lines. Their activity was discovered by the virtual high throughput screen (vHTS) against the phospholipase C-γ2 (PLC-γ2) isoform, a potential molecular target for anticancer therapy [2]. Testing against the National Cancer Institute's human tumour cell line panel (NCI-60) [3] revealed that the most potent analogues had low nano-molar growth inhibition, particularly against breast and melanoma tumour cell lines, but also against other cancer types [4][5][6][7].
TPs are reported to modulate a number of molecular targets that are implicated in cancer progression: (i) phospholipase C-δ1/3 (PLC), deduced by the same cellular behaviour of the MDA-MB-231 breast cancer cell line upon administration of TP and the knockout (silencing) of the PLC-δ1/3 genes; [8,9] (ii) copper-trafficking antioxidant 1 (ATOX1) protein, the inhibition of which reduces the proliferation of cancer cells [10]; (iii) tyrosyl DNA phosphodiesterase 1 (TDP 1), a phospholipase D enzyme, involved in repairing DNA damage [11]; (iv) the colchicine binding site in tubulin [12,13], an established target for anticancer drugs; and (v) adenosine A 2A receptor an established target for anticancer drugs; and (v) adenosine A2A receptor (A2AAR) [14], a G proteincoupled receptor (GPCR) implicated in the immune response to tumours [15]. The reported targets are shown in Figure 1, along with the most potent measured inhibition value and the corresponding TP molecular structure. To complicate the interpretation of this data, different assay formats are used for each of the targets. Also, various derivatives of the TPs give the best response, i.e., 1 and 2 give excellent overall potency in cell-based assays pointing to PLC-δ inhibition, 3 gives the most favourable IC50 for TDP1, 4 has the best tubulin inhibition profile, whereas 5 has the lowest reported Kd against ATOX1, and, finally, 6 is the most active analogue for A2AAR. Since these derivatives were discovered by different research groups, not using the same compound collections and synthetic strategies, it is not surprising that various derivatives of TP were discovered and developed. Furthermore, each molecular target should have a unique structure-activity relationship, resulting in different TP derivatives being tested. When this is considered in terms of chemical space, Figure 2 can be produced to map locations of the active derivatives, in which each ellipse represents the active volume of chemical space for each target. Taking into account the available data, an extensive overlap between the aforementioned targets can be envisioned. It is currently not known which receptor/enzyme's modulation produces the pronounced anticancer potential of the TPs or indeed whether a synergistic effect is responsible for the efficacy. For the further development of the TPs as potential clinical candidates, it is vital to elucidate the molecular mechanism responsible for its potency. The aim of this study is to glean insight into the molecular mechanism responsible for the anticancer potential of the TPs by conducting molecular modelling against the reported targets. The To complicate the interpretation of this data, different assay formats are used for each of the targets. Also, various derivatives of the TPs give the best response, i.e., 1 and 2 give excellent overall potency in cell-based assays pointing to PLC-δ inhibition, 3 gives the most favourable IC 50 for TDP1, 4 has the best tubulin inhibition profile, whereas 5 has the lowest reported K d against ATOX1, and, finally, 6 is the most active analogue for A 2A AR. Since these derivatives were discovered by different research groups, not using the same compound collections and synthetic strategies, it is not surprising that various derivatives of TP were discovered and developed. Furthermore, each molecular target should have a unique structure-activity relationship, resulting in different TP derivatives being tested. When this is considered in terms of chemical space, Figure 2 can be produced to map locations of the active derivatives, in which each ellipse represents the active volume of chemical space for each target. Taking into account the available data, an extensive overlap between the aforementioned targets can be envisioned. It is currently not known which receptor/enzyme's modulation produces the pronounced anticancer potential of the TPs or indeed whether a synergistic effect is responsible for the efficacy. For the further development of the TPs as potential clinical candidates, it is vital to elucidate the molecular mechanism responsible for its potency. an established target for anticancer drugs; and (v) adenosine A2A receptor (A2AAR) [14], a G proteincoupled receptor (GPCR) implicated in the immune response to tumours [15]. The reported targets are shown in Figure 1, along with the most potent measured inhibition value and the corresponding TP molecular structure. To complicate the interpretation of this data, different assay formats are used for each of the targets. Also, various derivatives of the TPs give the best response, i.e., 1 and 2 give excellent overall potency in cell-based assays pointing to PLC-δ inhibition, 3 gives the most favourable IC50 for TDP1, 4 has the best tubulin inhibition profile, whereas 5 has the lowest reported Kd against ATOX1, and, finally, 6 is the most active analogue for A2AAR. Since these derivatives were discovered by different research groups, not using the same compound collections and synthetic strategies, it is not surprising that various derivatives of TP were discovered and developed. Furthermore, each molecular target should have a unique structure-activity relationship, resulting in different TP derivatives being tested. When this is considered in terms of chemical space, Figure 2 can be produced to map locations of the active derivatives, in which each ellipse represents the active volume of chemical space for each target. Taking into account the available data, an extensive overlap between the aforementioned targets can be envisioned. It is currently not known which receptor/enzyme's modulation produces the pronounced anticancer potential of the TPs or indeed whether a synergistic effect is responsible for the efficacy. For the further development of the TPs as potential clinical candidates, it is vital to elucidate the molecular mechanism responsible for its potency. The aim of this study is to glean insight into the molecular mechanism responsible for the anticancer potential of the TPs by conducting molecular modelling against the reported targets. The The aim of this study is to glean insight into the molecular mechanism responsible for the anticancer potential of the TPs by conducting molecular modelling against the reported targets.
The predicted binding affinities should give an indication as to which experiments are most likely to shed light on the mode of action of this class of anticancer compounds.

Results
A docking panel of the putative bio-molecular targets of the TPs was developed consisting of PLC-δ1, ATOX1, TDP1, tubulin's three binding sites, and the A 2A AR GPCR. The aim was to compare the scores of the docked TPs derivatives in order to establish which target is most likely to be modulated, thus leading to the observed anticancer effects.

PLC-δ1
There is evidence for PLC-δ1 supporting the growth and migration of neoplastic mammary epithelial cells, linking this target to cancer growth [16]. Upon docking to PLC-δ1, derivatives 1-6 were predicted to have the same/similar binding modes previously reported using both GOLD and Glide [1,4,5,11]. Similar poses were generated with consistent interactions and strong affinity, as well as good overlap with the co-crystallised inositol 1,4,5-triphosphate. Hydrogen bonds with Asn312, Lys438, and His311 were observed as well as the side chains of Glu341, Tyr551, and Gly554, as shown in Figure 3. predicted binding affinities should give an indication as to which experiments are most likely to shed light on the mode of action of this class of anticancer compounds.

Results
A docking panel of the putative bio-molecular targets of the TPs was developed consisting of PLC-δ1, ATOX1, TDP1, tubulin's three binding sites, and the A2AAR GPCR. The aim was to compare the scores of the docked TPs derivatives in order to establish which target is most likely to be modulated, thus leading to the observed anticancer effects.

PLC-δ1
There is evidence for PLC-δ1 supporting the growth and migration of neoplastic mammary epithelial cells, linking this target to cancer growth [16]. Upon docking to PLC-δ1, derivatives 1-6 were predicted to have the same/similar binding modes previously reported using both GOLD and Glide [1,4,5,11]. Similar poses were generated with consistent interactions and strong affinity, as well as good overlap with the co-crystallised inositol 1,4,5-triphosphate. Hydrogen bonds with Asn312, Lys438, and His311 were observed as well as the side chains of Glu341, Tyr551, and Gly554, as shown in Figure 3.

ATOX1
The ATOX1 copper trafficking protein is linked to reduced proliferation of cancer cells demonstrated by knockout experiments of its gene [10]. Compound 6 is reported to form hydrogen bonds with the side chains of Glu17, Arg21, and Lys 60 amino acid residues, blocking the proteinprotein interface of ATOX1 for copper ion delivery [10]. In the case for 1 and 2 using GOLD, similar poses were predicted as those reported. Compound 1 blocked the copper-trafficking interface by forming a hydrogen bond with Arg21 and lipophilic contacts with Ala18 and Lys60, as shown in Figure 4, whilst 2 displayed hydrogen bonding interactions with Cys15 and Thr58, as well as a lipophilic contact between the phenyl and Lys60. Glide predicted a binding mode for 6 in which its thieno [2,3-b]pyridine moiety is positioned between Glu17, Arg21, and Lys60 side chains, with a hydrogen bond with the Glu17 side chain via its amino group and with the Cys12 side chain with the amide NH moiety. Cys12 and Cys15 are two of the four cysteine residues, which coordinate with a copper ion at the interface of two protein chains of ATOX1, thus blocking access to these residues inhibiting ion transfer. In the case of 1, however, a different binding pose was observed with the phenyl ring close to Lys60 and formed a hydrogen bond with Arg21.

ATOX1
The ATOX1 copper trafficking protein is linked to reduced proliferation of cancer cells demonstrated by knockout experiments of its gene [10]. Compound 6 is reported to form hydrogen bonds with the side chains of Glu17, Arg21, and Lys 60 amino acid residues, blocking the protein-protein interface of ATOX1 for copper ion delivery [10]. In the case for 1 and 2 using GOLD, similar poses were predicted as those reported. Compound 1 blocked the copper-trafficking interface by forming a hydrogen bond with Arg21 and lipophilic contacts with Ala18 and Lys60, as shown in Figure 4, whilst 2 displayed hydrogen bonding interactions with Cys15 and Thr58, as well as a lipophilic contact between the phenyl and Lys60. Glide predicted a binding mode for 6 in which its thieno [2,3-b]pyridine moiety is positioned between Glu17, Arg21, and Lys60 side chains, with a hydrogen bond with the Glu17 side chain via its amino group and with the Cys12 side chain with the amide NH moiety. Cys12 and Cys15 are two of the four cysteine residues, which coordinate with a copper ion at the interface of two protein chains of ATOX1, thus blocking access to these residues inhibiting ion transfer. In the case of 1, however, a different binding pose was observed with the phenyl ring close to Lys60 and formed a hydrogen bond with Arg21.

TDP1
Camptothecin and its analogues (e.g., topotecan and irinotecan) constitute an effective and widely used class of cancer therapeutics. They work by inhibiting the enzyme DNA topoisomerase I (topo I), preventing DNA replication. In human cells, TDP1 is important to the mechanism for removing DNA lesions [17,18]. In cancer cells, however, TDP1 counters the action of camptothecin by repairing drug-created stalled topo I-DNA complexes, thus making TDP1 a potential anticancer target [19,20]. For TDP1, both His263 and His493 are proposed to play a key role in its biological function and, as such, the binding pocket is defined at these amino acids [21][22][23]. When docked to TDP1 using GOLD, all the ligands fitted into the pocket in a similar way as those depicted for 1 in Figure 5. These results are in line with previously published modelling data for the TPs [11]. Glide predicted a good fit for all the ligands, despite some variations in hydrogen bonding patterns, as for the GOLD predictions. The main difference between the Glide-and GOLD-predicted poses was the reverse orientations of the phenyl and thieno[2,3-b]pyridine rings.

TDP1
Camptothecin and its analogues (e.g., topotecan and irinotecan) constitute an effective and widely used class of cancer therapeutics. They work by inhibiting the enzyme DNA topoisomerase I (topo I), preventing DNA replication. In human cells, TDP1 is important to the mechanism for removing DNA lesions [17,18]. In cancer cells, however, TDP1 counters the action of camptothecin by repairing drug-created stalled topo I-DNA complexes, thus making TDP1 a potential anticancer target [19,20]. For TDP1, both His263 and His493 are proposed to play a key role in its biological function and, as such, the binding pocket is defined at these amino acids [21][22][23]. When docked to TDP1 using GOLD, all the ligands fitted into the pocket in a similar way as those depicted for 1 in Figure 5. These results are in line with previously published modelling data for the TPs [11]. Glide predicted a good fit for all the ligands, despite some variations in hydrogen bonding patterns, as for the GOLD predictions. The main difference between the Glide-and GOLD-predicted poses was the reverse orientations of the phenyl and thieno[2,3-b]pyridine rings.

TDP1
Camptothecin and its analogues (e.g., topotecan and irinotecan) constitute an effective and widely used class of cancer therapeutics. They work by inhibiting the enzyme DNA topoisomerase I (topo I), preventing DNA replication. In human cells, TDP1 is important to the mechanism for removing DNA lesions [17,18]. In cancer cells, however, TDP1 counters the action of camptothecin by repairing drug-created stalled topo I-DNA complexes, thus making TDP1 a potential anticancer target [19,20]. For TDP1, both His263 and His493 are proposed to play a key role in its biological function and, as such, the binding pocket is defined at these amino acids [21][22][23]. When docked to TDP1 using GOLD, all the ligands fitted into the pocket in a similar way as those depicted for 1 in Figure 5. These results are in line with previously published modelling data for the TPs [11]. Glide predicted a good fit for all the ligands, despite some variations in hydrogen bonding patterns, as for the GOLD predictions. The main difference between the Glide-and GOLD-predicted poses was the reverse orientations of the phenyl and thieno[2,3-b]pyridine rings.

Tubulin
Tubulin-targeting agents are widely used for the treatment of cancer. They disrupt/inhibit tubulin dynamics by acting as stabilisers or destabilisers. Three different binding sites have been identified for targeting microtubules, namely colchicine-, taxol-, and Vinca alkaloid-binding pockets [24].
In order to establish the most plausible binding pocket for the TPs, 1-6 were docked to the three identified sites and the results are given in Table 1. A clear trend is seen for the colchicine site predicted for the GS and Glide values, with the ChemPLP results also leaning towards this site. CS and ASP gave similar scores for all three sites. Furthermore, the results are in line with earlier findings where TP derivatives have been predicted to bind to the colchicine site [12,13]. According to the GOLD docking results, compounds 1-6 have a similar predicted configuration in the colchicine-binding site and showed stacking/polar interactions with tubulin ( Figure 6), sitting deep in the hydrophobic pocket (Ser178, lipophilic contacts with Ala180 and Val181). Docking results with Glide were highly consistent for TPs 1-6, yet somewhat different from those produced by GOLD. The binding affinities of the Glide predictions were high, and a common hydrogen bond with Val181 backbone NH, also a feature of co-crystallised colchicine, was observed.

Tubulin
Tubulin-targeting agents are widely used for the treatment of cancer. They disrupt/inhibit tubulin dynamics by acting as stabilisers or destabilisers. Three different binding sites have been identified for targeting microtubules, namely colchicine-, taxol-, and Vinca alkaloid-binding pockets [24].
In order to establish the most plausible binding pocket for the TPs, 1-6 were docked to the three identified sites and the results are given in Table 1. A clear trend is seen for the colchicine site predicted for the GS and Glide values, with the ChemPLP results also leaning towards this site. CS and ASP gave similar scores for all three sites. Furthermore, the results are in line with earlier findings where TP derivatives have been predicted to bind to the colchicine site [12,13]. According to the GOLD docking results, compounds 1-6 have a similar predicted configuration in the colchicine-binding site and showed stacking/polar interactions with tubulin ( Figure 6), sitting deep in the hydrophobic pocket (Ser178, lipophilic contacts with Ala180 and Val181). Docking results with Glide were highly consistent for TPs 1-6, yet somewhat different from those produced by GOLD. The binding affinities of the Glide predictions were high, and a common hydrogen bond with Val181 backbone NH, also a feature of co-crystallised colchicine, was observed.

Testing for GPCRs Activity
One of the reported targets for the TPs is the A2AAR GPCR receptor with nano-molar activity [14]. It is estimated that 30-50% of all drugs on the market target GPCR receptors [25]. It was therefore decided to investigate whether the TPs modulate other GPCR receptors.
Derivative 1 was screened against the gpcrMAX Panel containing 168 receptors in both agonist and antagonist formats at 10 µM by the DiscoverX company [26]. The analysis of the data revealed five potential antagonist receptors and one agnostic receptor. In brief, ligands where considered as potential antagonists with inhibition >35% and agonists with activity >30% (see Tables S8-S10 in Supplementary Materials). These six receptors were tested for dose-response to quantify their effect and five targets (four antagonists/one agonist) were active; the results are shown in Table 2. The

Testing for GPCRs Activity
One of the reported targets for the TPs is the A 2A AR GPCR receptor with nano-molar activity [14]. It is estimated that 30-50% of all drugs on the market target GPCR receptors [25]. It was therefore decided to investigate whether the TPs modulate other GPCR receptors.
Derivative 1 was screened against the gpcrMAX Panel containing 168 receptors in both agonist and antagonist formats at 10 µM by the DiscoverX company [26]. The analysis of the data revealed five potential antagonist receptors and one agnostic receptor. In brief, ligands where considered as potential antagonists with inhibition >35% and agonists with activity >30% (see Tables S8-S10 in Supplementary  Materials). These six receptors were tested for dose-response to quantify their effect and five targets (four antagonists/one agonist) were active; the results are shown in Table 2. The receptor, HTR1B, showed activity in the primary screen but had IC 50 > 50 µM for the dose-response and was therefore not considered to be a hit. The dose response curves are given in the Supplementary Information  pages 15 and 16. Furthermore, the complete technical reports from the company DiscoverX containing controls are available in the Supplementary Information.
The five receptors are modulated by TP 1 in the low µ-molar region and are therefore less potent than A 2A AR (ADORA2A), the known GPCR, which was unfortunately not in the screening panel. Interestingly, a related Adenosine A3 (ADORA3) receptor was in the panel and some activity was observed (25%) in the antagonist format, but not enough to be taken forward for quantitative testing. Recently, a panel of TPs was tested for in vitro anti-platelet activity and showed both platelet activation and aggregation [27]. The authors concluded that the TPs are potential inhibitors for the P2Y 12 receptor, the target for clopidogrel, the standard clinical therapeutic. The P2Y 12 receptor (P2RY12) had 13% agonistic activity according to the assay, but no antagonist modulation was observed.
The activity of 1 against the identified targets was relatively modest, in the low µ-molar range (see Table 2), which implies that the low nano-molar activity seen against some tumour cell lines cannot be explained by the modulation of these receptors. Nevertheless, identifying activity against any GPCRs is helpful in elucidating their mechanisms of action and in the development of potent small molecules against specific receptors. Molecular modelling was performed in order to analyse and compare the binding mode of the TPs against A 2A AR and the identified receptors. A 2A AR and CXCR4 both have published crystal structures [28,29], and CRL-RAMP and GPR35 have published homologies models [30,31]. Protein structures for NPSR1B and PRLHP were obtained from the GPCR database [32].

A 2A AR
The A 2A AR receptor is linked to the immune response of tumours by T regulatory cell function, suppression and modulation of natural killer cell cytotoxicity, and tumour-specific CD4+/CD8+ activity [15]. In the case of A 2A AR antagonists, a common binding motif involves π-π stacking between aromatic moieties of the ligands and the Phe168 side chain, as well as hydrogen bonding interactions with Asn253 and Glu169 side chains [14].
GOLD and Glide predicted very similar binding poses for compounds 1-6. Key interactions such as hydrogen bonding with Asn253 and Glu169 side chains and π-π stacking with the Phe168 aromatic side chain were common among the derivatives, as were high binding affinities (Figure 7).

CXCR4
CXCR4 is a chemokine receptor involved in cancer metastasis and inflammatory diseases, as well as HIV-1 infection [29]. Interestingly, modulation of this receptor with small organic ligands has been reported for both metastasis and inflammation [33].
GOLD predicted hydrogen bonding between Try94 and Asp97 side chains via the amino group of the TPs as seen for 1 in Figure 8. Also, π-π stacking with Trp94 side chain was observed for some compounds with the thieno[2,3-b]pyridine ring. Glide docking was very consistent in terms of binding poses, which were in agreement with the GOLD predictions and the co-crystallised conformation of ITD ((6,6-dimethyl-5,6-dihydroimidazo [2,1-b] [1,3]thiazol-3-yl)methyl N,N′dicyclohexylimidothiocarbamate).

CRL-RAMP3
Calcitonin receptor-like activity modifying protein 3 (CRL-RAMP3), is a member of the RAMP family of receptors, single transmembrane-domain proteins. CRL-RAMP3 regulates calcitonin and is a linear polypeptide hormone involved in the calcium ion regulation in blood [34]. The CRL-RAMP3

CXCR4
CXCR4 is a chemokine receptor involved in cancer metastasis and inflammatory diseases, as well as HIV-1 infection [29]. Interestingly, modulation of this receptor with small organic ligands has been reported for both metastasis and inflammation [33].
GOLD predicted hydrogen bonding between Try94 and Asp97 side chains via the amino group of the TPs as seen for 1 in Figure 8. Also, π-π stacking with Trp94 side chain was observed for some compounds with the thieno[2,3-b]pyridine ring. Glide docking was very consistent in terms of binding poses, which were in agreement with the GOLD predictions and the co-

CXCR4
CXCR4 is a chemokine receptor involved in cancer metastasis and inflammatory diseases, as well as HIV-1 infection [29]. Interestingly, modulation of this receptor with small organic ligands has been reported for both metastasis and inflammation [33].
GOLD predicted hydrogen bonding between Try94 and Asp97 side chains via the amino group of the TPs as seen for 1 in Figure 8. Also, π-π stacking with Trp94 side chain was observed for some compounds with the thieno[2,3-b]pyridine ring. Glide docking was very consistent in terms of binding poses, which were in agreement with the GOLD predictions and the co-crystallised conformation of ITD ((6,6-dimethyl-5,6-dihydroimidazo [2,1-b] [1,3]thiazol-3-yl)methyl N,N′dicyclohexylimidothiocarbamate).

CRL-RAMP3
Calcitonin receptor-like activity modifying protein 3 (CRL-RAMP3), is a member of the RAMP family of receptors, single transmembrane-domain proteins. CRL-RAMP3 regulates calcitonin and is a linear polypeptide hormone involved in the calcium ion regulation in blood [34]. The CRL-RAMP3

CRL-RAMP3
Calcitonin receptor-like activity modifying protein 3 (CRL-RAMP3), is a member of the RAMP family of receptors, single transmembrane-domain proteins. CRL-RAMP3 regulates calcitonin and is a linear polypeptide hormone involved in the calcium ion regulation in blood [34]. The CRL-RAMP3 complex has been indicated to be important to angiogenesis and therefore to tumour growth [35]. A published homology model of RAMP3 was used for the modelling [30]. The model of RAMP3 suggests that Tyr84 can contact the Try52 of adrenomedullin (AM) and calcitonin gene-related peptides (CGRP) Phe37 phenyl rings and is therefore expected to decrease the potency of both peptides at the AM2 receptor [36,37]. Furthermore, Glu74 and Tyr84 may interact with ligands that bind to the CRL-RAMP3 heterodimer.
The GOLD and Glide docking predicted similar poses for all ligands, which occupied the cleft between CRL and RAMP3. They were observed to π-π stack with Trp72 and Tyr124 aromatic side chains of CRL. The thieno[2,3-b]pyridine moiety was surrounded by RAMP3 residues Glu74, Trp84, and Arg38 (Figure 9). complex has been indicated to be important to angiogenesis and therefore to tumour growth [35]. A published homology model of RAMP3 was used for the modelling [30]. The model of RAMP3 suggests that Tyr84 can contact the Try52 of adrenomedullin (AM) and calcitonin gene-related peptides (CGRP) Phe37 phenyl rings and is therefore expected to decrease the potency of both peptides at the AM2 receptor [36,37]. Furthermore, Glu74 and Tyr84 may interact with ligands that bind to the CRL-RAMP3 heterodimer. The GOLD and Glide docking predicted similar poses for all ligands, which occupied the cleft between CRL and RAMP3. They were observed to π-π stack with Trp72 and Tyr124 aromatic side chains of CRL. The thieno[2,3-b]pyridine moiety was surrounded by RAMP3 residues Glu74, Trp84, and Arg38 (Figure 9). Docking of the TPs using GOLD yielded plausible binding poses. Hydrogen bonds were predicted with Arg239 side chain amide NH2 via the carbonyl on cycloalkane and with Arg148 via NH2 on the thiophene ring. The TPs are predicted to form π-π stacking with Phe147.
Glide docking predicted two different binding trends for 1-6. Compounds 1 (Figure 10), 3, and 6 aligned in a similar way as was observed with GOLD, forming hydrogen bonds with Asn230 and Arg239 side chains. Alternatively, 2, 4, and 5 formed two hydrogen bonds; one with the Cys146 backbone and one with the Ser246 side chain. Docking of the TPs using GOLD yielded plausible binding poses. Hydrogen bonds were predicted with Arg239 side chain amide NH 2 via the carbonyl on cycloalkane and with Arg148 via NH 2 on the thiophene ring. The TPs are predicted to form π-π stacking with Phe147.
Glide docking predicted two different binding trends for 1-6. Compounds 1 (Figure 10), 3, and 6 aligned in a similar way as was observed with GOLD, forming hydrogen bonds with Asn230 and Arg239 side chains. Alternatively, 2, 4, and 5 formed two hydrogen bonds; one with the Cys146 backbone and one with the Ser246 side chain.

NPSRb1
Neuropeptide S receptor b1 (NPSRb1) is linked to bronchial smooth muscles and asthma as well as inflammatory bowel disease [38]. Two residues in the second transmembrane region (TM2) of human NPSRb1, namely Asn107 and Asp105, are reported to be important for antagonist affinity [39,40]. Furthermore, three known antagonists of NPSRb1, with nano-molar activities, are known ( Figure 11) [41,42]. This information was used to locate the binding site of the receptor and to assess which of the amino acid residues are important for a successful binding. Docking the known antagonists in the vicinity of the target residues yielded plausible binding poses. The TPs almost perfectly overlapped with the predicted poses of known antagonists, SHA66 and SHA68. The phenyl ring was inserted into the cleft between two TM domains and parallel to helices, just like the p-fluorophenyl moiety of SHA68 and the TP moiety was in contact with the extracellular domain residues (Figure 12). In this conformation, the carbonyl group formed hydrogen bonds with the amino group of Ile106.

NPSRb1
Neuropeptide S receptor b1 (NPSRb1) is linked to bronchial smooth muscles and asthma as well as inflammatory bowel disease [38]. Two residues in the second transmembrane region (TM2) of human NPSRb1, namely Asn107 and Asp105, are reported to be important for antagonist affinity [39,40]. Furthermore, three known antagonists of NPSRb1, with nano-molar activities, are known ( Figure 11) [41,42]. This information was used to locate the binding site of the receptor and to assess which of the amino acid residues are important for a successful binding.

NPSRb1
Neuropeptide S receptor b1 (NPSRb1) is linked to bronchial smooth muscles and asthma as well as inflammatory bowel disease [38]. Two residues in the second transmembrane region (TM2) of human NPSRb1, namely Asn107 and Asp105, are reported to be important for antagonist affinity [39,40]. Furthermore, three known antagonists of NPSRb1, with nano-molar activities, are known ( Figure 11) [41,42]. This information was used to locate the binding site of the receptor and to assess which of the amino acid residues are important for a successful binding. Docking the known antagonists in the vicinity of the target residues yielded plausible binding poses. The TPs almost perfectly overlapped with the predicted poses of known antagonists, SHA66 and SHA68. The phenyl ring was inserted into the cleft between two TM domains and parallel to helices, just like the p-fluorophenyl moiety of SHA68 and the TP moiety was in contact with the extracellular domain residues (Figure 12). In this conformation, the carbonyl group formed hydrogen bonds with the amino group of Ile106. Docking the known antagonists in the vicinity of the target residues yielded plausible binding poses. The TPs almost perfectly overlapped with the predicted poses of known antagonists, SHA66 and SHA68. The phenyl ring was inserted into the cleft between two TM domains and parallel to helices, just like the p-fluorophenyl moiety of SHA68 and the TP moiety was in contact with the extracellular domain residues (Figure 12). In this conformation, the carbonyl group formed hydrogen bonds with the amino group of Ile106.

PRLHR
Prolactin releasing hormone receptor (PRLHR) has seven transmembrane domains. Prolactin is a protein hormone involved in lactation of female mammals. There is some evidence that it increases the risk of breast cancer [43]. To our knowledge, no biological or theoretical information regarding binding site for PRLHR ligands exists. Therefore, possible binding site cavities were determined using Site Map followed by docking of 1-6 to the identified cavities on the homology model. Docking results obtained from the cavity, which was detected at the centre of TM2-6 and close to extracellular loops, were reasonably consistent for both GOLD and Glide.
In the binding mode of 1 obtained from GOLD, a hydrogen bond with Asn298 was observed and π-π stacking with His321 ( Figure 13). Compound 2 formed hydrogen bonds with Asn298 and Thr117, and π-π stacking with Phe138 and Trp291. Glide predicted similar binding for 1, 2, 4, and 5, with consistent interactions with the putative binding site residues. They had hydrogen bonding with Asn298 side chain donor via their amide carbonyl and π-π stacking with Trp291 and His295 side chains with the phenyl ring of TS.

PRLHR
Prolactin releasing hormone receptor (PRLHR) has seven transmembrane domains. Prolactin is a protein hormone involved in lactation of female mammals. There is some evidence that it increases the risk of breast cancer [43]. To our knowledge, no biological or theoretical information regarding binding site for PRLHR ligands exists. Therefore, possible binding site cavities were determined using Site Map followed by docking of 1-6 to the identified cavities on the homology model. Docking results obtained from the cavity, which was detected at the centre of TM2-6 and close to extracellular loops, were reasonably consistent for both GOLD and Glide.
In the binding mode of 1 obtained from GOLD, a hydrogen bond with Asn298 was observed and π-π stacking with His321 ( Figure 13). Compound 2 formed hydrogen bonds with Asn298 and Thr117, and π-π stacking with Phe138 and Trp291. Glide predicted similar binding for 1, 2, 4, and 5, with consistent interactions with the putative binding site residues. They had hydrogen bonding with Asn298 side chain donor via their amide carbonyl and π-π stacking with Trp291 and His295 side chains with the phenyl ring of TS.

PRLHR
Prolactin releasing hormone receptor (PRLHR) has seven transmembrane domains. Prolactin is a protein hormone involved in lactation of female mammals. There is some evidence that it increases the risk of breast cancer [43]. To our knowledge, no biological or theoretical information regarding binding site for PRLHR ligands exists. Therefore, possible binding site cavities were determined using Site Map followed by docking of 1-6 to the identified cavities on the homology model. Docking results obtained from the cavity, which was detected at the centre of TM2-6 and close to extracellular loops, were reasonably consistent for both GOLD and Glide.
In the binding mode of 1 obtained from GOLD, a hydrogen bond with Asn298 was observed and π-π stacking with His321 ( Figure 13). Compound 2 formed hydrogen bonds with Asn298 and Thr117, and π-π stacking with Phe138 and Trp291. Glide predicted similar binding for 1, 2, 4, and 5, with consistent interactions with the putative binding site residues. They had hydrogen bonding with Asn298 side chain donor via their amide carbonyl and π-π stacking with Trp291 and His295 side chains with the phenyl ring of TS.

Discussion
As 1-6 were docked using five different scoring functions, an estimate was generated to identify which molecular targets are most likely to be modulated by the TPs. The scoring function results are given for each derivative 1-6 in the Supplementary Information. The results, for the six ligands, were averaged and the standard deviation was derived (Table 3). Table 3. Predicted average scores and standard deviations for derivatives 1-6 against the identified targets. The best scores are in bold and the second best in italics. Looking at the results, for GS A 2A AR has the highest score, closely followed by tubulin, CXCR4, and GPR35, while the other targets have clearly lower predicted scores. CS again predicts A 2A AR to have the highest score, with the second highest shared by the rest of the GPCRs, except NSPRb1. The ChemPLP scoring function put A 2A AR as the top target, followed by CXCR4 and CRL-RAMP3. ChemPLP has been shown to be the best, or one of the best, performing scoring functions available [44,45]. ASP predicted A 2A AR to be the best, followed by CXCR4 and GPR35. Finally, Glide predicted A 2A AR to have the highest affinity with the TPs, followed by tubulin. Based on these results, it is apparent that A 2A AR is the most likely target, with the highest predicted score for all of the scoring functions used.
The results for all the scoring functions were normalised to the highest score (A 2A AR), giving them an equal weight. They were then summed and divided to give the best score of 10 (Total column in Table 3) and standard deviations were derived. When these numbers were analysed, it was seen that the GPCRs were predicted to be modulated, except NSPRb1, with tubulin also obtaining a good score, followed by PLC-δ1. According to this analysis, the TPs mostly exert their anticancer effect by GPCR modulation.
As can be seen in Figure 1, the TPs derivatives do not have radically different chemical structures; however, when the scoring data is analysed, notable shifts in predicted affinities are observed, e.g., removing the substituents on the phenyl ring from 1 resulting in 2 increases the affinity of the latter towards PLC-δ1 for both the GS and PLP scoring functions (see Tables S2-S7 in the Supplementary Materials). This means that each modification of the TPs changes the affinity to many of the targets, rendering it very difficult to generate a reliable structure activity relationship (SAR) for cell-based data. Thus, it is a multidimensional problem to optimise the TPs' structure for anticancer activity. The changes in water solubility and, consequently, cell penetration, between the derivatives also needs to be considered, adding further complexity to the problem.

Transcript Levels in the Cancer Cells Affected by TPs
The NCI60 transcriptome database and Cellminer tools [46] (discover.nci.nih.gov/cellminer) were used to assess the expression of mRNA of the identified targets in the National Institute of Cancer human tumour cell lines panel (NCI-60) [3]. The z-score representations, a measure of mRNA transcript intensity [46], were collected for all the available targets and the results are given in the Supplementary Materials. Unfortunately, CRL-RAMP3 and PRLHR were not available in the database. The melanoma cell line MDA-MB-435 is the most susceptible tumour cell for the TPs in the NCI-60 panel and the z-scores for the targets are given in Table 4, as well as the minimum and maximum values given for the cell line panel. According to the values in Table 4, only PLC-δ1 and ATOX1 have positive values, i.e., significant mRNA is expressed in this cell line making these targets plausible. TDP1 is expressed to some level but TDP1's expression levels are expected to rise under stress, leading to DNA damage (such a chemotherapy) to aid the cell to survive. Also, tubulin is heavily expressed during mitosis, making it abundant during cell duplication in all of the cancer cell lines. In general, the transcript intensity of the GPCR receptors is low, in particular for CXCR4, but GPR35, A 2A AR, and NSPRb1 are expressed to some extent.
For further analysis, the expression levels of each target were considered and correlated with the GI 50 values of ligand 1, as shown in Tables S11-S20 in the Supplementary Materials.
In general, the data suggest that gene expression level is necessary, but not sufficient, to predict the efficacy of inhibition. mRNA abundance and consequent protein expression are known to correlate poorly due to complex biological and technical issues [47].
For Glide, the crystal structures were prepared using Protein Preparation Wizard of Maestro [64], where hydrogens were added, ionisation and tautomeric states were generated by Epik, and proton orientations were set by PROPKA [64]. The same coordinates as those used for GOLD were used for search the space centre in grid generation, and the van der Waals radii of non-polar receptor atoms (partial charge cut-off 0.25) were scaled down to 0.8 Å. Ligands were docked at extra precision, post-docking minimisation was enabled, Epikstate penalties were added to docking scores, and a maximum of 10 poses were recorded per ligand.
In order to verify the robustness of the docking protocol, the available co-crystallised ligands were re-docked and their RMSD (Root-Mean-Square Deviation) values generated. The results are shown in Table S1 in the Supplementary Materials, and an adequate to excellent overlap was obtained.

Homology Modelling
The homology models were prepared with Protein Preparation Wizard of Maestro [64]. For the models of NPRS1B and PRLHR, which were obtained from the GPCR data base [32], restrained force field minimisation was performed with convergence of heavy atoms to an RMSD of 0.5 Å to prevent clashes. The coordinates of N1 of Trp-72 (chain A) side chain for CALCRL-RAMP3, the centroid of docked bufrolin for GPR35, the midpoint of Asp105 and Asn107 for NPSR1B, and the centroid of the predicted cavity by Site Map (Schrödinger, LLC, New York, NY, USA, 2015) [64] for PRLHR were picked as the centres of search spaces for each model.

Biological Testing
Compound 1 was screened against the gpcrMAX Panel using 168 different GPCR receptors, in both agonist and antagonist formats. The assays were performed using the PathHunter beta-arrestin enzyme fragment complementation (EFC) technology at DiscoverX [26]. The full data set and the analysis protocol are given in the Supplementary Materials.

Conclusions
This work communicates the discovery of activity of TPs against GPCRs, albeit in the low µM range. Finding potent ligands against unexplored GPCRs is very helpful in the elucidation of their nature and is the first step in developing a specific ligand for modulation and possibly a clinical candidate. The modelling data indicates A 2A AR and tubulin to be the most likely targets of the TPs; however, when protein expression data is considered, other targets such as ATOX1 and PLC-δ1 are evidenced to play a role in the anticancer activity of this class of compounds. It is clear that the TPs are prolific inhibitors/modulators of various targets and can therefore be categorised as a privileged structure in drug discovery [65]. The challenge is to develop the TPs to be specific against the target of interest in order to develop them into viable drug candidates.