Prediction and Identification of GPCRs Targeting for Drug Repurposing in Osteosarcoma

Abstract

Background:

Osteosarcoma (OS) is a malignant bone tumor common in children and adolescents. The 5-year survival rate is only 67-69% and there is an urgent need to explore novel drugs effective for the OS. G protein-coupled receptors (GPCRs) are the common drug targets and have been found to be associated with the OS, but have been seldom used in OS.

Methods:

The GPCRs were obtained from GPCRdb, and the GPCRs expression profile of the OS was downloaded from the UCSC Xena platform including clinical data. 10-GPCRs model signatures related to OS risk were identified by risk model analysis with R software. The predictive ability and pathological association of the signatures in OS were explored by bio-informatics analysis. The therapeutic effect of the target was investigated, followed by the investigation of the targeting drug by the colony formation experiment were.

Results:

We screened out 10 representative GPCRs from 50 GPCRs related to OS risk and established a 10-GPCRs prognostic model (with CCR4, HCRTR2, DRD2, HTR1A, GPR158, and GPR3 as protective factors, and HTR1E, OPN3, GRM4, and GPR144 as risk factors). We found that the low-risk group of the model was significantly associated with the higher survival probability, with the area under the curve (AUC) of the ROC greater than 0.9, conforming with the model. Moreover, both risk-score and metastasis were the independent risk factor of the OS, and the risk score was positively associated with the metastatic. Importantly, the CD8 T-cells were more aggregated in the low-risk group, in line with the predict survival rate of the model. Finally, we found that DRD2 was a novel target with approved drugs (cabergoline and bromocriptine), and preliminarily proved the therapeutic effects of the drugs on OS. These novel findings might facilitate the development of OS drugs.

Conclusion:

This study offers a satisfactory 10-GPCRs model signature to predict the OS prognostic, and based on the model signature, candidate targets with approved drugs were provided.

Introduction

Osteosarcoma is a malignant bone tumor that is common in children and adolescents. Surprisingly, there has been little improvement in the survival rate of patients since the 1970s, with a 5-year survival rate of only 67-69%, which is one of the highest mortality rates among childhood and adolescent cancers (1). Patients with metastases and relapses can hardly be cured by conventional methods like surgery and radiotherapy (2). As an alternative, chemotherapy is also limited by cytotoxicity, drug resistance, etc., with unsatisfactory survival rate and sometimes serious side effects (3, 4), for example, methotrexate-induced acute encephalopathy (5). This may be due to the lack of tissue specificity of chemotherapeutic. Therefore, there is an urgent need to explore novel drugs specifically targeting OS.

Recently, the critical role of GPCRs in bone development, remodeling, and disease has been verified, with 92 GPCRs associated with bone disease and dysfunction (6). The GPCRs are the largest family of transmembrane proteins involved in multiple biological processes, including bone development and remodeling (7, 8), inflammation/immune response (9), tumor growth, and metastasis (10), etc. A growing body of evidence suggests that GPCRs serve as pro-tumor (GPR56 (11), GPR110 (12), PTHR1 (13), PAR1 (14), S1PR3 (15)) or anti-tumor (GPER1 (16), A3AR (17)) in OS. However, most studies only focused on one type of GPCRs, which may not reflect the real condition since the receptors synergistically exert the pro-tumor or anti-tumor effects by highly integrated interactions of numerous receptor signals. Thus, it is imperative to integrate all the GPCRs to investigate their correlation with OS-risk for the screening of the representative GPCRs to predict the OS-risk. Specifically, GPCRs have easy-to-target ligand-binding domains to bind with a variety of chemical regulators, which constitute the most important class of drug targets (18, 19). This may benefit the selection of the appropriate drug for OS therapy by repurposing the known drugs for new indications (20). In this study, in an attempt to find the suitable known-drug repurposing for OS, we investigated the association of the GPCRs with OS progression and determined the GPCRs that serve as drug targets by risk evaluation. Firstly, we analyzed the 85 samples of OS from the TARGET database and assessed the risk correlation of all GPCRs in the GPCRdb with OS. GPCRdb is a database containing all human non-olfactory GPCRs, including over 2000 approved drugs and in-trial agents and nearly 200,000 ligands with activity and availability data (21). Second, to avoid excessive intervention for low-risk patients, we established a risk-score model by using risk-GPCRs and explored the predictive ability of the model and its components on OS risk. We also investigated the correlation between the model and the infiltration abundance of immune cells. Finally, we preliminarily verified the effects of the GPCRs-targeting drugs on OS by using two types of OS cell lines. This research will provide a novel insight into OS therapy.

Materials and Methods

Data Source and Preprocessing

Human non-olfactory GPCRs were obtained from GPCRdb (https://gpcrdb.org/), and a total of 395 GPCRs were obtained. Eighty-five OS samples containing 395 GPCRs expression data and clinical information were downloaded from the GDC TARGET on the UCSC XENA platform (22) (https://xena.ucsc.edu/). The single-cell RNA-seq dataset (GSE152048) was downloaded from the GEO database (www.ncbi.nlm.nih.gov/geo).

We preprocessed the TARGET-OS dataset to filter low-expressed GPCRs, the deletion measure was to remove those GPCRs with expression level < 1 and accounting for more than 50% of all samples. Due to the low expression, only 333 GPCRs were finally included in the analysis. The details of these GPCRs and the clinical characteristics of the samples were shown in Figure 2A.

Identification of Key GPCRs and Construction of Prognostic Models

Univariate Cox analysis of the cohort. Perform univariate Cox proportional hazard regression for each GPCRs in the cohort, and performed the “coxph” function of the “survival” package in R with P < 0.05 as the filter threshold.

LASSO-Cox analysis. In order to decrease the number of GPCRs in the risk model, we performed Least Absolute Shrinkage and Selection Operator (LASSO) regression on the GPCRs obtained by univariate Cox analysis to reduce the over-fitting phenomenon. This was realized by performed “glmnet” and “cv.glmnet” functions of the “glmnet” package in R.

Multivariate Cox analysis of risk GPCRs. In order to screen out GPCRs that can collaboratively predict OS risk, and to further simplify the number of GPCRs in the risk model, multivariate Cox analysis was performed on the GPCRs obtained from the LASSO-Cox analysis, and the “Akaike Information Criteria (AIC)” was used for stepwise reduction of the variables. The GPCRs combination with the lowest AIC score was used to further construct a risk model. The receiver operator characteristic (ROC) curve was used to evaluate the accuracy of the model.

Predict Survival Rate

Kaplan Meier Plotter (23) (http://kmplot.com/analysis/index.php?p=service&cancer=custom_plot) was used to predict the survival rate of the risk model and its genes. Auto-select best cutoff.

Calculation of TME Immune Cell Infiltration Abundance

The infiltration of the tumor microenvironment (TME) immune cells may important to the patient’s outcome. The CIBERSORTx algorithm (https://cibersortx.stanford.edu/) was used to quantify the infiltration abundance of 22 types of immune cells in OS. The CIBERSORTx parameters were as follows: the input matrix was the RNA-Seq data containing 85 samples and 33464 genes, the 22 immune cell types from Newman et al. (24) were input as the reference of genes signature, 500 times for permutation test, and RNA-seq data without quantile normalization.

Calculating the EMT Score

We obtained 77 marker genes of epithelial-to-mesenchymal transition (EMT) from Mak (25). The EMT score for each sample was calculated as , like the described in Chen’s study (26), where M represents the expression of the mesenchymal genes, E represents the expression of the epithelial genes, N and n represents the number of mesenchymal genes and epithelial genes respectively.

Mapping Risk GPCRs to the OS Single-Cell Atlas

In order to detect and evaluate the expression of risk GPCRs in the OS micro-environment, we used the OS single-cell atlas. The “Seurat” package (27) in R was used for the single-cell atlas analysis of the GSE152048 dataset (28). We preprocessed the GSE152048 dataset to filter the low-quality cells and genes, the inclusion criteria for genes were expressed in at least 5 cells, the inclusion criteria for cells were that at least 300 genes were expressed in the cell, and the proportion of mitochondrial genes and hemoglobin genes in the cell were less than 5%. The biomarkers provided in the original article (28) were used for the annotation of the cell population.

Mining the Targeted Drugs of Risk-Relation-GPCRs

We searched the targeted drugs of GPCRs in the drug&ligands menu of GPCRdb. The screening criteria of the targeted drugs were those have reported in tumor treatment.

Cell Culture

The osteosarcoma cell line (143B, ATCC: CRL-8303, and HOS, ATCC: HTB-96TM, were purchased from the American Type Culture Collection) were cultured in DMEM/F-12 medium (DMEM, Gibco, Shanghai, China) containing 10% fetal bovine serum (Gibco, American), and 100 µg/mL streptomycin and 100 U/mL penicillin (Solarbio, Beijing, China). All the cells were placed in an incubator containing humidified air with 5% carbon dioxide at 37°C, and the medium was replaced with fresh medium once every 3 days.

CCK8-Kit to Detect Cell Viability

143B and HOS cells suspension were plated into 96-well plates, the number of cells in each well was about 8×10 ³, and they were treated with cabergoline and bromocriptine (0, 1, 10, 50, 100, 150, 200 μmol/L), after 24 hours of incubation, remove the drug, add 100 μl of culture medium and 10μl of CCK-8 detection reagent, and then incubated in an incubator for two hours. Detect the absorbance of each well of the plate in the microplate reader, and set the wavelength to 450 nm. The experiments were repeated three times under the same conditions.

Colony Formation Experiment

Cells (2×10³) of 143B and HOS cells were plated into 6-well plate. After 24 hours, 143B cells 100 μM cabergoline and bromocriptine, and HOS cells were treated with 50-100 μM cabergoline and bromocriptine, 7days later, it was washed twice with phosphate buffer solution and fixed with 4% paraformaldehyde for 10 minutes. After washing with deionized water, it was stained with crystal violet dye for 10 minutes. After washing again, the cell colonies were counted with Image J.

Statistical Analysis

In this study, R software (version 4.0.3) was applied to the statistical analysis process, unless specifically stated, otherwise, only when p-value < 0.05 was considered statistically significant (main packages include: survival, glmnet, Survminer, survivorROC, Seurat, limma, ggplot2, clusterProfiler, org.Hs.eg. db).

Result

A protocol was designed to construct a multivariate model for predicting the OS prognosis and to develop targets and targeted drugs. The analysis process was carried out by following protocol (Figure 1).

Figure 1

Identification and Functional Annotation of Risk-Related GPCRs

The information of the GPCRs and the OS samples in this study were shown in (Figure 2A). 50 GPCRs related to OS risk were obtained through univariate Cox proportional hazard regression analysis (P < 0.05). Among them, the high expression of 45 GPCRs was associated with the OS-low risk (protective factors), while 5 GPCRs’ high expression was associated with the OS-high risk (risk factors) (Figure 2B). specified, the GPR158 was the most significant protective factor, and MTNR1B was the most significant risk factor. GO (gene ontology) and KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis was used to annotate the biological processes, molecular functions, and signaling pathways of the 50 GPCRs. The results showed (Figures 2C, D) that, the most enriched GO term were those who participated in the regulation of cytosolic calcium ion concentration, G protein-coupled peptide receptor activity, and immune receptor activity. In terms of KEGG, neuroactive ligand−receptor interaction, calcium signaling pathway, and chemokine signaling pathway were the most enrichment signaling pathways for risk-related GPCRs.

Figure 2

Construction of Risk Prognosis Model Based on GPCRs

The 50 GPCRs were considered to be numerous prognostic factors and the workload was still heavy in screening the targets. Therefore, we performed LASSO-Cox regression analysis to decrease the numbers of GPCRs, meanwhile reducing the over-fitting phenomenon. The results show that with the increase of Lambda, the coefficient of the independent variables was decreasing, and the number of independent variables whose coefficient’ absolute value was greater than 0 was decreasing (Figure 3A). When using 10-fold cross-validation to construct the model, we could see when the cross-validation was the smallest, only 19 GPCRs remained as candidate GPCRs for subsequent analysis (Figure 3B). To further decrease the number of GPCRs, 19 GPCRs were performed multivariate Cox regression analysis and AIC stepwise regression analysis. There is a consensus that when the AIC is the smallest, the model is the best, it gave the model sufficient fit with fewest variables. Finally, the fewest of GPCRs were obtained to construct the model (Figure 3C), they were: CCR4, DRD2, HCRTR2, HTR1A, HTR1E, OPN3, GPR3, GRM4, GPR158, GPR144. Together, these 10 GPCRs would be used to build the final prognostic model.

Figure 3

The 10 GPCRs were integrated as a risk prognosis model to evaluate the prognostic value of the OS patients. First, integrate these GPCRs into the multivariate Cox model, and the “predict” function of the “survival” package in R was used to calculate the risk score of each OS patient. Then, according to the median risk score, the OS patients were divided into high-risk and low-risk groups. The risk curves and scatter plots was illustrated the relation of the model’s risk score and the vital status of OS patients (Figures 4A, B), we could see that the higher of the risk score, the more of the dead patients. The heat map showed the relationship between the expression level of 10 GPCRs and the risk score (Figure 4C). Among them, 6 GPCRs with HR < 1 (CCR4, HCRTR2, GPR3, DRD2, HTR1A, GPR158), and their high expression was related to the low-risk. The 4 GPCRs with HR>1 (GRM4, OPN3, GPR144, HTR1E), and their high expression was related to the high-risk. The KM survival curve of the patient suggested that the survival probability of the OS patients in the low risk-score group was significant prolonged (p=8.5e-9) than those in the high risk-score group, which indicated that the risk score of the model has a satisfactory prognostic value (Figure 4D). In addition, the R package “survivalROC” was used to calculate the area under the ROC curve to quantify the prediction accuracy of the model (Figure 4E). The results showed a satisfactory prognostic value of the model in the existing samples, it could accurately predict the survival rate of patients in 2, 3, and 5 years (with AUC > 0.9).

Figure 4

The KM survival curves of 10 GPCRs were shown in (Figure 5A), and all of 10 GPCRs have a significant prognostic value (P <0.05). The grouping of samples was performed by Kaplan-Meier plotter according to GPCRs expression auto select best cutoff (Figure 5B).

Figure 5

Univariate and Multivariate Cox Analysis of the Model and Clinical Features

The independence of 10-GPCRs signature models in clinical applications was tested by using univariate and multivariate Cox regression analyses (Table 1). Using the binary classification method, the samples were divided into two groups for analysis according to metastatic (yes VS. no), gender (male VS. female), age (>14 VS. <14), and risk score (high VS. low). The univariate analysis results suggested that the risk-score of the model was significant associated with the hazard ratio (p = 9.81e-06). At the same time, the results of the multivariate analysis indicated that the risk score could be used as an independent factor because it was not interfered with by other factors (hazard ratio=15.627, 95% confidence interval = 4.364 - 55.956, p = 2.40e-05). Similarly, here we found that metastatic was also an independent prognostic factor.

The Correlation Between the 10 GPCRs and Metastasis

Metastasis is generally considered to be an intractable feature of the OS. Here the KM survival curve showed that the survival rate of non-metastatic patients were significantly prolonged than those of metastatic patients (P = 1.5e-05), confirming the metastasis has a useful prognostic value for OS (Figure 6A). Based on this result, we also explored the Pearson correlation between the expression of 10-GPCRs and OS-metastasis and found that DRD2 was the GPCRs that most negatively correlated with OS-metastasis, indicating that the higher expression of DRD2, the lower rate of metastasis (Figure 6B). To further reveal the expression of 10 GPCRs in metastatic-OS and non-metastatic-OS (Figure 6C), we analyzed the expression of 10 GPCRs in both type OS, and the results showed that the expression of DDR2 in non-metastatic-OS was significantly higher than the metastatic-OS (P = 0.0031), while the expression of HTR1E in metastatic-OS was significantly higher than the non-metastatic-OS (P = 0.018). In addition, it’s may note that the expression of GPR158 in non-metastatic-OS was seemed higher than the metastatic-OS (P = 0.057, close to significance). There was no significant difference in the expression of other GPCRs in metastatic-OS and non-metastatic-OS.

Figure 6

The Association Between 10-GPCRs Model and TME-Infiltrating Immune Cells

As mentioned above, the GPCRs related to OS-risk were highly enriched in immune receptor activity and Chemokine signaling pathway. Therefore, we attempted to investigate whether the risk grouping of the 10-GPCRs model was associated with the infiltration of immune cells. We used the CIBERSORTx to calculate the abundance of immune cells in each sample, by removing one outlier sample, and then, “limma” was used to analyze the difference infiltration of immune cells between the high-risk group and the low-risk group. We observed that CD8 T cells (P=0.00191) and monocytes (P=0.01735) were more aggregated in the low-risk group (Figure 7A). The heatmap illustrates the relationship between the risk grouping and the patient’s status, metastasis, EMT score, and the infiltration abundance of CD8 T-cells and monocytes (Figure 7B). Statistical analysis indicated the EMT score was not a significant difference between the low-risk group and high-risk group (Figure 7C), and the infiltration of CD8 T-cell and monocytes was significantly more aggregated in the low-risk group (Figure 7D). By survival analysis, we found that EMTscore could not predict the survival rate of OS patients (Figure 7E), while the infiltration abundance of CD8 T-cells could effectively predict the survival rate of patients (Figure 7F), which conforms to the logic that CD8 T-cells could kill tumor cells for treatment, but unfortunately monocytes could not effectively predict the survival rate (Figure 7G). All in all, monocytes and CD8 T-cells related to immune activation were more aggregated in the low-risk group, this might be why the risk score of the model could predict patients’ survival.

Figure 7

The Expression of 10 GPCRs in OS Microenvironment

By analyzing the OS single-cell RNA-seq dataset, we obtained the OS atlas (Figures 8A–C), the atlas contains a total of 82,567 cells and 26,163 genes. The OS atlas was annotated into 11 cell subgroups and divided into immune cell subgroups and non-immune cell subgroups, among which the immune cell subgroups were the cluster 2, 5, 11, the non-immune cell subgroups were the cluster 0, 1, 3, 4, 6, 7, 8, 9, 10. Then mapping the expression of 10 GPCRs in the OS atlas (Figure 8D), the results showed that CCR4 was mainly expressed in tumor-infiltrating lymphocytes (TILs), and its positive-cell proportion (PCP) in the tumor micro-environment was about 0.046%. HCRTR2 was mainly expressed in non-immune cells, and its PCP was about 0.13%. GPR3 was mainly expressed in non-immune cells, and its PCP was about 1.36%. DRD2 was mainly expressed in non-immune cells, and its PCP was about 0.11%. GPR158 was mainly expressed in non-immune cells, and its PCP was about 1.05%. GRM4 was mainly expressed in non-immune cells, and its PCP was about 1.05%. OPN3 was generally expressed in various types of cells in the OS micro-environment, and its PCP was about 9.87%. HTR1E was mainly expressed in non-immune cells, and its PCP was about 0.02%. Among them, the expression of HTR1A and GPR144 could not be detected in this OS atlas, it might be caused by the low RNA capture rate at single-cell level. This result simply revealed the cellular expression context of the GPCRs.

Figure 8

OS Targeted Drugs Based on the Model Signature

The prognostic value of 10-GPCRs in the OS has been demonstrated above. Therefore, aimed at the 10 GPCRs as therapeutic targets may improve the treatment effect of the OS. Consistent with our research expectations, we found 3 targets with known drugs in 10GPCRs, they were CCR4, DRD2 and GRM4. Among them, CCR4 has a drug in the trial, DRD2 has 66 approved drugs and 32 trial drugs, GRM4 has 2 drugs in trials (21). Thus, DRD2 was the easiest target for drug repurposing, and due to the DRD2 was a protective factor, we hypothesize its agonist will give a better treatment effect on OS. Therefore, we used the agonists of DRD2 (cabergoline and bromocriptine) to test our hypothesis, and we found that cabergoline and bromocriptine did indeed have an inhibiting effect on OS by CCK8 cell viability experiment, where the half-inhibitory concentration (IC50) of cabergoline and bromocriptine in metastatic-OS (143B) were about 100 μM, the IC50 value of cabergoline in non-metastatic-OS (HOS) was about 100 μM, the IC50 value of bromocriptine in the HOS cells was about 150 μM, (Figure 9A). The colony formation experiment confirmed that cabergoline and bromocriptine could significantly inhibit the colony formation of 143B and HOS, which was characterized by smaller colony spots and fewer colonies (Figure 9B), the quantitative comparison between the experimental group and the control group showed that the drugs inhibitory effects were significant (Figure 9C). Here we provided an effective novel target with known drugs, this might increase new drug options for the treatment of OS.

Figure 9

Discussion

According to cancer statistics in 2020, the 5-year survival rate of the OS is still only 67-69% (1). The conventional treatment for OS patients is surgical resection and neoadjuvant chemotherapy, but the clinical outcomes have been little improvement since the 1970s. Therefore, there is an urgent need for potential bio-targets to assess the risk of OS patients.

In this study, we obtained the expression profiles of GPCRs from the TARGET-OS database on the UCSC Xena platform. By using univariate Cox regression analysis, Lasso-Cox regression analysis, and multivariate Cox regression analysis, 10 GPCRs that were significantly related to OS risk were obtained. A prognostic model of OS based on the above 10 GPCRs was established. This model had high reliability (Kaplan–Meier curve’s p-value < 0.001, and AUC > 0.9), and high discriminatory ability in predicting survival rate. In the model, 6 GPCRs (CCR4, HCRTR2, GPR3, DRD2, HTR1A, GPR158) were highly associated with a low risk of OS (protective factor), and 4 GPCRs (GRM4, OPN3, GPR144, HTR1E) were highly associated with a high risk of OS (risk factors). Furthermore, we determined that the model and the metastasis could be used as independent risk prognostic factors through the joint analysis of univariate and multivariate. In addition, based on the Pearson correlation analysis between the metastasis and the expression level of 10-GPCRs, it was believed that the expression of DRD2 (P = 0.0031) was negatively correlated with the OS-metastasis, and the expression of HTR1E (P = 0.018) was positively correlated with OS-metastasis. We also found that the model could distinguish the samples with different immune cells, for example, the CD8 T-cells and monocytes were more aggregated in the low-risk group, and the infiltration ratio of CD8 T-cells was associated with a high survival rate of OS. The single-cell atlas was used to detect the cellular identity of 10-GPCRs in the OS micro-environment. Finally, in 10-GPCRs, we found that DRD2 was a target with approved drugs. By targeting DRD2, the agonists (bromocriptine and cabergoline) have the effect of inhibiting OS. Thus, we found targets with known drugs based on the risk model signature by using total GPCRs modeling for cancer prediction.

For protective factors, CCR4 was expressed in multiple subtypes of T cells, including effector CD8 T cells, Chemotaxis assay indicated that CCR4⁺ CD8 T cells could be recruited by CCL22 treatment (29), which was consistent with the finding that CCR4 was higher expressed in the low-risk group. Meanwhile, CD8 T cells were relatively higher infiltrated in the low-risk group. Both were correlated with a prolonger OS survival. In various cancer cell lines, HCRTR2 expression was down-regulated, which may initiate promoter hypermethylation (30, 31). After orexin receptor HCRTR2 was activated by orexins, apoptosis would be induced (32), indicating that the HCRTR2 had similar epigenetic characteristics to tumor suppressors. For GPR3, its agonists may render the breast cancer cells susceptible to cytotoxicity induced by cationic amphiphilic drugs, leading to lysosomal-dependent cell death (33), which showed the tumor-suppressive effect (34). By using the specific agonist, HTR1A can inhibit DNA synthesis (35), with lower expression in invasive cancer cells compared with non-invasive cancer cells (36). GPR158 has higher methylation level in esophageal squamous cell carcinoma compared with mucosa (37), which can induce apoptosis and is related to the high survival rate of glioma patients (38). For DRD2, it was considered as a tumor suppressor (39). It was reported that the DRD2 on cell membrane could exert anti-tumor effects by down-regulating eEF1A2 (40), while the eEF1A2 could promote the proliferation, migration, and invasion of OS cells by activating the Akt/mTOR signaling pathway (41), which implies that DRD2 has the anti-OS effect. It was also reported that the DRD2 agonists could exert anti-tumor effects through ROCK-mediated inactivation of Cofilin-1 or inhibiting EGFR/AKT/MMP-13 pathway (42, 43). It’s worth noting that the positive expression of Cofilin-1 in OS was related to clinical stage, distant metastasis, and tumor grade (44), and the high expression of EGFR in OS was related to high proliferative activity, metastatic potential, and poor prognosis (45), suggesting that DRD2 could be a drug target of OS. To summarize, the reports about the function of the protective factors were consistent with the view that high expression of the protective factors was associated with a higher survival rate of OS in this study.

For risk factors, our previous study showed that the high expression of GRM4 was associated with the lower survival rate of OS patients (46). For OPN3, its high expression was associated with a lower survival rate in lung adenocarcinoma and acral melanoma (47, 48). GPR144 (ADGRD2) was an important mediator in the hypoxic response of glioblastoma and had a significant tumor-promoting effect (49), associated with rheumatoid factor and osteosarcoma proliferation and invasion (50). For HTR1E, the activation of HTR1E could regulate the release of cytokines IL-6 and CXCL8 (51), which mediate osteosarcoma-lung interactions crucial to metastasis (52). The above report about the risk factors was consistent with our statistical results.

Bromocriptine and cabergoline are known drugs targeting DRD2 and are clinically used to treat pituitary adenomas and Parkinson’s disease. Treated for 24 h, the IC50 value of bromocriptine and cabergoline was about 100μM for some pituitary tumor cell lines (53), which is a similar lethal dose. Notably, bromocriptine can induct the expression level of DRD2 (54), which may upregulate the tumor suppressor (39). Besides, cabergoline and bromocriptine may also target other receptors, such as cabergoline can target 5HT(2A-C), D₃, and 5HT(1A) (also call HTR1A, as a protective factor in this study) (55). Bromocriptine was reported to regulate amino acid biosynthesis and metabolism (56). This study may provide a novel insight into the drug therapy of OS.

Conclusion

In summary, based on GDC TARGET and GPCRdb, we constructed an effective 10-GPCRS risk-score model (including CCR4, HCRTR2, GPR3, DRD2, HTR1A, GPR158, GRM4, OPN3, GPR144, HTR1E) to predict OS prognosis. The stability and accuracy of the model were evaluated by ROC curve. In addition, we found 3 GPCRs targets with known drugs in 10-GPCRs model, among which one target (DRD2) and its drugs (bromocriptine and cabergoline) were confirmed to inhibit the OS cells in this study. This study may provide reference for OS therapy in clinic.

Funding

This study was financially supported by the National Key R&D Program of China (2018YFC1105900), National Natural Science Fund of China (Grant No. 81972120), the Guangxi Science and Technology Base and Talent Special Project (Grant No. GuikeAD19254003), the Distinguished Young Scholars Program of Guangxi Medical University, the Youth Science Foundation of Guangxi Medical University(GXMUYSF202107), and the Scientific Research Project for Young Teachers of Guangxi-Collaborative Innovation Center for Biomedical (02406221008C).

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

• 1

  SiegelRLMillerKDJemalA. Cancer Statistics, 2020. CA: Cancer J Clin (2020) 70:7–30. doi: 10.3322/caac.21590

  • CrossRef
  • Google Scholar

• 2

  WhelanJSDavisLE. Osteosarcoma, Chondrosarcoma, and Chordoma. J Clin Oncol (2018) 36:188–93. doi: 10.1200/JCO.2017.75.1743

  • CrossRef
  • Google Scholar

• 3

  LagmayJPKrailoMDDangHKimAHawkinsDSBeatyOIIIet al. Outcome of Patients With Recurrent Osteosarcoma Enrolled in Seven Phase II Trials Through Children’s Cancer Group, Pediatric Oncology Group, and Children’s Oncology Group: Learning From the Past to Move Forward. J Clin Oncol (2016) 34:3031. doi: 10.1200/JCO.2015.65.5381

  • CrossRef
  • Google Scholar

• 4

  PalmeriniEJonesRLMarchesiEPaioliACesariMLonghiAet al. Gemcitabine and Docetaxel in Relapsed and Unresectable High-Grade Osteosarcoma and Spindle Cell Sarcoma of Bone. BMC Cancer (2016) 16:1–8. doi: 10.1186/s12885-016-2312-3

  • CrossRef
  • Google Scholar

• 5

  InabaHKhanRLaninghamFCrewsKPuiC-HDawNC. Clinical and Radiological Characteristics of Methotrexate-Induced Acute Encephalopathy in Pediatric Patients With Cancer. Ann Oncol (2008) 19:178–84. doi: 10.1093/annonc/mdm466

  • CrossRef
  • Google Scholar

• 6

  LuoJSunPSiwkoSLiuMXiaoJ. The Role of GPCRs in Bone Diseases and Dysfunctions. Bone Res (2019) 7:1–19. doi: 10.1038/s41413-019-0059-6

  • CrossRef
  • Google Scholar

• 7

  LuoJZhouWZhouXLiDWengJYiZet al. Regulation of Bone Formation and Remodeling by G-Protein-Coupled Receptor 48. Development (2009) 136:2747–56. doi: 10.1242/dev.033571

  • CrossRef
  • Google Scholar

• 8

  LuoJYangZMaYYueZLinHQuGet al. LGR4 is a Receptor for RANKL and Negatively Regulates Osteoclast Differentiation and Bone Resorption. Nat Med (2016) 22:539–46. doi: 10.1038/nm.4076

  • CrossRef
  • Google Scholar

• 9

  SharmaNAkhadeASQadriA. Sphingosine-1-Phosphate Suppresses TLR-Induced CXCL8 Secretion From Human T Cells. J Leukocyte Biol (2013) 93:521–8. doi: 10.1189/jlb.0712328

  • CrossRef
  • Google Scholar

• 10

  DorsamRTGutkindJS. G-Protein-Coupled Receptors and Cancer. Nat Rev Cancer (2007) 7:79–94. doi: 10.1038/nrc2069

  • CrossRef
  • Google Scholar

• 11

  ChenZGaoPLiZ. Expression of G Protein-Coupled Receptor 56 is an Unfavorable Prognostic Factor in Osteosarcoma Patients. Tohoku J Exp Med (2016) 239:203–11. doi: 10.1620/tjem.239.203

  • CrossRef
  • Google Scholar

• 12

  LiuZZhangGZhaoCLiJ. Clinical Significance of G Protein-Coupled Receptor 110 (GPR110) as a Novel Prognostic Biomarker in Osteosarcoma. Med Sci monitor: Int Med J Exp Clin Res (2018) 24:5216. doi: 10.12659/MSM.909555

  • CrossRef
  • Google Scholar

• 13

  SunCLiS. PTHR1 in Osteosarcoma: Specific Molecular Mechanisms and Comprehensive Functional Perspective. J Cell Mol Med (2021) 25:3175–81. doi: 10.1111/jcmm.16420

  • CrossRef
  • Google Scholar

• 14

  McAuleyJRBaileyKMEkambaramPKleiLRKangHHuDet al. MALT1 is a Critical Mediator of PAR1-Driven NF-κb Activation and Metastasis in Multiple Tumor Types. Oncogene (2019) 38:7384–98. doi: 10.1038/s41388-019-0958-4

  • CrossRef
  • Google Scholar

• 15

  ShenYZhaoSWangSPanXZhangYXuJet al. S1P/S1PR3 Axis Promotes Aerobic Glycolysis by YAP/c-MYC/PGAM1 Axis in Osteosarcoma. EBioMedicine (2019) 40:210–23. doi: 10.1016/j.ebiom.2018.12.038

  • CrossRef
  • Google Scholar

• 16

  WangZChenXZhaoYJinYZhengJ. G-Protein-Coupled Estrogen Receptor Suppresses the Migration of Osteosarcoma Cells via Post-Translational Regulation of Snail. J Cancer Res Clin Oncol (2019) 145:87–96. doi: 10.1007/s00432-018-2768-4

  • CrossRef
  • Google Scholar

• 17

  IyerSVRanjanAEliasHKParralesASasakiHRoyBCet al. Genome-Wide RNAi Screening Identifies TMIGD3 Isoform1 as a Suppressor of NF-κb and Osteosarcoma Progression. Nat Commun (2016) 7:1–13. doi: 10.1038/ncomms13561

  • CrossRef
  • Google Scholar

• 18

  SantosRUrsuOGaultonABentoAPDonadiRSBologaCGet al. A Comprehensive Map of Molecular Drug Targets. Nat Rev Drug Discovery (2017) 16:19–34. doi: 10.1038/nrd.2016.230

  • CrossRef
  • Google Scholar

• 19

  HauserASAttwoodMMRask-AndersenMSchiöthHBGloriamDE. Trends in GPCR Drug Discovery: New Agents, Targets and Indications. Nat Rev Drug Discovery (2017) 16:829–42. doi: 10.1038/nrd.2017.178

  • CrossRef
  • Google Scholar

• 20

  PushpakomSIorioFEyersPAEscottKJHopperSWellsAet al. Drug Repurposing: Progress, Challenges and Recommendations. Nat Rev Drug Discovery (2019) 18:41–58. doi: 10.1038/nrd.2018.168

  • CrossRef
  • Google Scholar

• 21

  KooistraAJMordalskiSPándy-SzekeresGEsguerraMMamyrbekovAMunkCet al. GPCRdb in 2021: Integrating GPCR Sequence, Structure and Function. Nucleic Acids Res (2021) 49:D335–43. doi: 10.1093/nar/gkaa1080

  • CrossRef
  • Google Scholar

• 22

  GoldmanMJCraftBHastieMRepečkaKMcDadeFKamathAet al. Visualizing and Interpreting Cancer Genomics Data via the Xena Platform. Nat Biotechnol (2020) 38:675–8. doi: 10.1038/s41587-020-0546-8

  • CrossRef
  • Google Scholar

• 23

  NagyÁMunkácsyGGyőrffyB. Pancancer Survival Analysis of Cancer Hallmark Genes. Sci Rep (2021) 11:1–10. doi: 10.1038/s41598-021-84787-5

  • CrossRef
  • Google Scholar

• 24

  NewmanAMLiuCLGreenMRGentlesAJFengWXuYet al. Robust Enumeration of Cell Subsets From Tissue Expression Profiles. Nat Methods (2015) 12:453–7. doi: 10.1038/nmeth.3337

  • CrossRef
  • Google Scholar

• 25

  MakMPTongPDiaoLCardnellRJGibbonsDLWilliamWNet al. Patient-Derived, Pan-Cancer EMT Signature Identifies Global Molecular Alterations and Immune Target Enrichment Following Epithelial-to-Mesenchymal Transition. Clin Cancer Res (2016) 22:609–20. doi: 10.1158/1078-0432.CCR-15-0876

  • CrossRef
  • Google Scholar

• 26

  ChenHYaoJBaoRDongYZhangTDuYet al. Cross-Talk of Four Types of RNA Modification Writers Defines Tumor Microenvironment and Pharmacogenomic Landscape in Colorectal Cancer. Mol Cancer (2021) 20:1–21. doi: 10.1186/s12943-021-01322-w

  • CrossRef
  • Google Scholar

• 27

  StuartTButlerAHoffmanPHafemeisterCPapalexiEMauckWMIIIet al. Comprehensive Integration of Single-Cell Data. Cell (2019) 177:1888–902.e1821. doi: 10.1016/j.cell.2019.05.031

  • CrossRef
  • Google Scholar

• 28

  ZhouYYangDYangQLvXHuangWZhouZet al. Single-Cell RNA Landscape of Intratumoral Heterogeneity and Immunosuppressive Microenvironment in Advanced Osteosarcoma. Nat Commun (2020) 11:1–17. doi: 10.1038/s41467-020-20059-6

  • CrossRef
  • Google Scholar

• 29

  LiHChenXZengWZhouWZhouQWangZet al. Radiation-Enhanced Expression of CCL22 in Nasopharyngeal Carcinoma is Associated With CCR4+ CD8 T Cell Recruitment. Int J Radiat Oncol Biol Phys (2020) 108:126–39. doi: 10.1016/j.ijrobp.2020.05.001

  • CrossRef
  • Google Scholar

• 30

  HuSNiuJZhangRLiXLuoMSangTet al. Orexin A Associates With Inflammation by Interacting With OX1R/OX2R Receptor and Activating Prepro-Orexin in Cancer Tissues of Gastric Cancer Patients. Gastroenterología y Hepatología (English Edition) (2020) 43:240–7. doi: 10.1016/j.gastre.2019.10.009

  • CrossRef
  • Google Scholar

• 31

  MisawaKImaiAMochizukiDMisawaYEndoSHosokawaSet al. Genes Encoding Neuropeptide Receptors are Epigenetic Markers in Patients With Head and Neck Cancer: A Site-Specific Analysis. Oncotarget (2017) 8:76318. doi: 10.18632/oncotarget.19356

  • CrossRef
  • Google Scholar

• 32

  LaburtheMVoisinTEl FirarA. Orexins/hypocretins and Orexin Receptors in Apoptosis: A Mini-Review. Acta Physiol (2010) 198:393–402. doi: 10.1111/j.1748-1716.2009.02035.x

  • CrossRef
  • Google Scholar

• 33

  AnandALiuBGiacobiniJDMaedaKRohdeMJäätteläM. Cell Death Induced by Cationic Amphiphilic Drugs Depends on Lysosomal Ca2+ Release and Cyclic AMP. Mol Cancer Ther (2019) 18:1602–14. doi: 10.1158/1535-7163.MCT-18-1406

  • CrossRef
  • Google Scholar

• 34

  MehlmannLMSaekiYTanakaSBrennanTJEvsikovAVPendolaFLet al. The Gs-Linked Receptor GPR3 Maintains Meiotic Arrest in Mammalian Oocytes. Science (2004) 306:1947–50. doi: 10.1126/science.1103974

  • CrossRef
  • Google Scholar

• 35

  SulaimanPJosephBKaimalSBPauloseCS. Decreased Hepatic 5-HT 1A Receptors During Liver Regeneration and Neoplasia in Rats. Neurochem Res (2008) 33:444–9. doi: 10.1007/s11064-007-9452-4

  • CrossRef
  • Google Scholar

• 36

  HenriksenRDizeyiNAbrahamssonP-A. Expression of Serotonin Receptors 5-HT1A, 5-HT1B, 5-HT2B and 5-HT4 in Ovary and in Ovarian Tumours. Anticancer Res (2012) 32:1361–6.

  • Google Scholar

• 37

  OkaDYamashitaSTomiokaTNakanishiYKatoHKaminishiMet al. The Presence of Aberrant DNA Methylation in Noncancerous Esophageal Mucosae in Association With Smoking History: A Target for Risk Diagnosis and Prevention of Esophageal Cancers. Cancer: Interdiscip Int J Am Cancer Soc (2009) 115:3412–26. doi: 10.1002/cncr.24394

  • CrossRef
  • Google Scholar

• 38

  LiNZhangYSidlauskasKEllisMEvansIFrankelPet al. Inhibition of GPR158 by microRNA-449a Suppresses Neural Lineage of Glioma Stem/Progenitor Cells and Correlates With Higher Glioma Grades. Oncogene (2018) 37:4313–33. doi: 10.1038/s41388-018-0277-1

  • CrossRef
  • Google Scholar

• 39

  TanYSunRLiuLYangDXiangQLiLet al. Tumor Suppressor DRD2 Facilitates M1 Macrophages and Restricts NF-κb Signaling to Trigger Pyroptosis in Breast Cancer. Theranostics (2021) 11:5214. doi: 10.7150/thno.58322

  • CrossRef
  • Google Scholar

• 40

  LeeJM. The Role of Protein Elongation Factor Eef1a2 in Ovarian Cancer. Reprod Biol Endocrinol (2003) 1:1–5. doi: 10.21236/ADA431692

  • CrossRef
  • Google Scholar

• 41

  YangJTangJLiJCenYChenJDaiG. Effect of Activation of the Akt/mTOR Signaling Pathway by EEF1A2 on the Biological Behavior of Osteosarcoma. Ann Trans Med (2021) 9:158. doi: 10.21037/atm-20-7974

  • CrossRef
  • Google Scholar

• 42

  PeverelliEGiardinoETreppiediDLocatelliMVairaVFerreroSet al. Dopamine Receptor Type 2 (DRD2) Inhibits Migration and Invasion of Human Tumorous Pituitary Cells Through ROCK-Mediated Cofilin Inactivation. Cancer Lett (2016) 381:279–86. doi: 10.1016/j.canlet.2016.08.005

  • CrossRef
  • Google Scholar

• 43

  HuangHWuKMaJDuYCaoCNieY. Dopamine D2 Receptor Suppresses Gastric Cancer Cell Invasion and Migration via Inhibition of EGFR/AKT/MMP-13 Pathway. Int Immunopharmacol (2016) 39:113–20. doi: 10.1016/j.intimp.2016.07.002

  • CrossRef
  • Google Scholar

• 44

  LiZFYaoYDZhaoYYLiuYLiuZHHuPet al. Effects of PAK4/LIMK1/Cofilin-1 Signaling Pathway on Proliferation, Invasion, and Migration of Human Osteosarcoma Cells. J Clin Lab Anal (2020) 34:e23362. doi: 10.1002/jcla.23362

  • CrossRef
  • Google Scholar

• 45

  LinderMGlitznerESrivatsaSBakiriLMatsuokaKShahrouziPet al. EGFR is Required for FOS-Dependent Bone Tumor Development via RSK2/CREB Signaling. EMBO Mol Med (2018) 10:e9408. doi: 10.15252/emmm.201809408

  • CrossRef
  • Google Scholar

• 46

  YangWMaolinHJinminZZheW. High Expression of Metabotropic Glutamate Receptor 4: Correlation With Clinicopathologic Characteristics and Prognosis of Osteosarcoma. J Cancer Res Clin Oncol (2014) 140:419–26. doi: 10.1007/s00432-013-1581-3

  • CrossRef
  • Google Scholar

• 47

  XuCWangRYangYXuTLiYXuJet al. Expression of OPN3 in Lung Adenocarcinoma Promotes Epithelial-Mesenchymal Transition and Tumor Metastasis. Thorac Cancer (2020) 11:286–94. doi: 10.1111/1759-7714.13254

  • CrossRef
  • Google Scholar

• 48

  ZengWZhangWFengJHeXLuH. Expression of OPN3 in Acral Lentiginous Melanoma and its Associated With Clinicohistopathologic Features and Prognosis. Immunity Inflammation Dis (2021) 9:840–50. doi: 10.1002/iid3.438

  • CrossRef
  • Google Scholar

• 49

  BayinNFrensterJKaneJRubensteinJModrekABaitalmalRet al. GPR133 (ADGRD1), an Adhesion G-Protein-Coupled Receptor, is Necessary for Glioblastoma Growth. Oncogenesis (2016) 5:e263–3. doi: 10.1038/oncsis.2016.63

  • CrossRef
  • Google Scholar

• 50

  TsengW-YWuY-JJYangT-YChiangN-YTsaiW-PGordonSet al. High Levels of Soluble GPR56/ADGRG1 are Associated With Positive Rheumatoid Factor and Elevated Tumor Necrosis Factor in Patients With Rheumatoid Arthritis. J Microbiol Immunol Infect (2018) 51:485–91. doi: 10.1016/j.jmii.2016.11.010

  • CrossRef
  • Google Scholar

• 51

  BayerHMullerTMyrtekDSorichterSZiegenhagenMNorgauerJet al. Serotoninergic Receptors on Human Airway Epithelial Cells. Am J Respir Cell Mol Biol (2007) 36:85–93. doi: 10.1165/rcmb.2006-0151OC

  • CrossRef
  • Google Scholar

• 52

  GrossACCamHPhelpsDASarafAJBidHKCamMet al. IL-6 and CXCL8 Mediate Osteosarcoma-Lung Interactions Critical to Metastasis. JCI Insight (2018) 3:e99791. doi: 10.1172/jci.insight.99791

  • CrossRef
  • Google Scholar

• 53

  TangCSunRWenGZhongCYangJZhuJet al. Bromocriptine and Cabergoline Induce Cell Death in Prolactinoma Cells via the ERK/EGR1 and AKT/mTOR Pathway Respectively. Cell Death Dis (2019) 10:1–14. doi: 10.1038/s41419-019-1526-0

  • CrossRef
  • Google Scholar

• 54

  YangYMamouniKLiXChenYKavuriSDuYet al. Repositioning Dopamine D2 Receptor Agonist Bromocriptine to Enhance Docetaxel Chemotherapy and Treat Bone Metastatic Prostate Cancer. Mol Cancer Ther (2018) 17:1859–70. doi: 10.1158/1535-7163.MCT-17-1176

  • CrossRef
  • Google Scholar

• 55

  PlataniaCBMLeggioGMDragoFSalomoneSBucoloC. Regulation of Intraocular Pressure in Mice: Structural Analysis of Dopaminergic and Serotonergic Systems in Response to Cabergoline. Biochem Pharmacol (2013) 86:1347–56. doi: 10.1016/j.bcp.2013.08.010

  • CrossRef
  • Google Scholar

• 56

  EricsonEGebbiaMHeislerLEWildenhainJTyersMGiaeverGet al. Off-Target Effects of Psychoactive Drugs Revealed by Genome-Wide Assays in Yeast. PloS Genet (2008) 4:e1000151. doi: 10.1371/journal.pgen.1000151

  • CrossRef
  • Google Scholar