Skip to main content
Download PDF

Immunological role and clinical prognostic significance of P2RY6 in lung adenocarcinoma: a multi-omics studies and single-cell sequencing analysis

World Journal of Surgical Oncology volume 21, Article number: 341 (2023) Cite this article

Abstract

Background

There is increasing evidence that recombinant human P2Y purinoceptor 6 (P2RY6) may be involved in inflammatory responses. However, the role of P2RY6 in lung adenocarcinoma (LUAD) remains unknown.

Methods

We used transcriptomic, genomic, single-cell transcriptomic, and methylation sequencing data from The Cancer Genome Atlas database to analyze the aberrant status and prognostic value of P2RY6 in a variety of tumors. The LUAD single-cell sequencing dataset was used to explore the effect of P2RY6 on the tumor microenvironment. Cell-type identification by estimating relative subsets of RNA transcripts (CIBERSORT) was used to quantify immune cells in the tumor microenvironment. We also analyzed the correlation of P2RY6 with immune checkpoints and immune regulation-related genes. The correlation of between tumor mutation burden (TMB), microsatellite instability (MSI), and P2RY6 expression was also analyzed simultaneously. Tissue microarray and immunohistochemistry were employed to assess the expression of P2RY6 in internal tumor samples.

Results

Our findings indicate that P2RY6 exhibits significantly higher expression levels in various cancer tissues, particularly in LUAD. High expression of P2RY6 was closely associated with a poor prognosis for patients, and it plays a role in regulating immune-related pathways, such as cytokine-cytokine receptor interaction. Notably, P2RY6 expression is closely linked to the abundance of CD8 + T cells. Furthermore, we have developed a P2RY6-related inflammation prediction model that demonstrates promising results in predicting the prognosis of LUAD patients, with an AUC (area under the curve) value of 0.83. This performance is significantly better than the traditional TNM staging system. Through single-cell transcriptome sequencing analysis, we observed that high P2RY6 expression is associated with increased intercellular communication. Additionally, pathway enrichment analysis revealed that P2RY6 influences antigen presentation and processing pathways within the LUAD microenvironment.

Conclusions

This study suggests that P2RY6 would be a new target for immunotherapy in LUAD.

Introduction

Lung cancer is currently the most prevalent form of malignancy worldwide, with the highest incidence and mortality rates [1]. It can be classified into two subtypes: small cell lung cancer (SCLC) and non-small cell lung cancer (NSCLC). NSCLC accounts for 85% of all cases [2]. Among NSCLC, lung adenocarcinoma (LUAD) is the most prevalent histological type of primary lung tumor [3]. Current treatment options for lung cancer include targeted therapy and immunotherapy, in combination with traditional surgery, radiotherapy, and chemotherapy. These advancements have significantly improved the survival rate of patients [4, 5]. However, due to diagnostic limitations, many LUAD patients are diagnosed at advanced stages of the disease, often accompanied by metastasis, resulting in a poor prognosis [6]. Despite extensive research on the pathogenesis of LUAD and the development of new treatment approaches, it unfortunately remains one of the most aggressive and deadly types of tumors, with a 5-year survival rate of approximately 23% [7]. In recent years, immune checkpoint inhibitors have demonstrated strong efficacy in many types of cancer [8,9,10]. The tumor microenvironment of the lung is known to be infiltrated by a diverse range of immune cells [11, 12], including natural killer cells, B cells, macrophages, dendritic cells, and T lymphocytes. As such, it has been suggested that LUAD may be a “hot” tumor type, and could be particularly sensitive to immunotherapy [13, 14]. Exploring new immunotherapy targets of LUAD and exploring the immune escape mechanism of LUAD is of great significance for improving the prognosis of LUAD patients.

P2RY6, also known as P2Y6, belongs to the P2 receptor family and is activated by extracellular nucleotides. It undergoes variable splicing, resulting in multiple transcript variants that encode different protein isoforms [15]. Recent evidence suggests that P2RY6 signaling plays a role in regulating inflammation [16, 17]. Specifically, it has been found to have significant pro-inflammatory and pro-atherogenic roles in macrophages. Under inflammatory conditions, P2RY6 is involved in regulating neutrophil migration by controlling TLR2-induced IL-8 release from human monocytes [18, 19]. Furthermore, studies have demonstrated that P2RY6 receptors can contribute to host defense against bacterial infections, suggesting that they may act as novel mediators in enhancing the innate immune response to invading pathogens [20]. While previous research has investigated the role of P2RY6 in pancreatic cancer [21] and colorectal cancer [22], there is currently no corresponding article exploring its role in LUAD.

It is noteworthy to mention that this study is the first extensive analysis of P2RY6 expression and its correlation with various malignancies, including LUAD. By further investigating the association between P2RY6 expression and the prognosis of LUAD, the study has confirmed its findings through immunohistochemistry. The results of this study provide novel insights into the function of P2RY6 in LUAD and elucidate the intrinsic mechanisms by which P2RY6 influences the tumor microenvironment. Furthermore, these findings underscore the potential utility of P2RY6 in tumor immunotherapy applications.

Methods

Clinical specimens

In this study, twenty patients with LUAD were included, and samples of cancer tissue and adjacent normal paracancerous tissue were collected. It is crucial to note that these patients were diagnosed with LUAD based on pathological examination and had not undergone any prior radiotherapy or chemotherapy treatments. Prior informed consent was obtained from all patients, indicating their voluntary participation in the study. To ensure the ethical conduct of the research, the study protocol underwent review and approval by the Ethics and Anthropology Committee. All experiments and methods were performed in accordance with applicable standards and regulations, thereby guaranteeing the integrity and reliability of the study results.

Data collection

We downloaded The Cancer Genome Atlas (TCGA) [23], genotypic tissue expression (GTEx) [24], RNA expression information, and clinical data from the UCSC Xena database [25] (https://xenabrowser.net/datapages/), which included more than 1000 tumor samples from 33 human cancers. We downloaded DNA copy number and methylation data from the cBioPortal database [26] (https://www.cbioportal.org/), and we ensured that the acquisition and application methods were consistent with relevant guidelines and regulations.

Single-cell transcriptome sequencing data acquisition and analysis

The single-cell sequencing datasets were derived from GEO database (GSE131907). We obtained processed and cell annotated data of 11 LUAD primary tumor. All samples were sequenced using the Hiseq X10 (Illumina, San Diego, CA) with standard parameters. R (version 3.5.2) and Seurat R package (version 3.1.1) were used for quality control (QC) and analysis.

Differential expression analysis of P2RY6

We used SangerBox (http://vip.sangerbox.com/home.html) software to compare the difference between P2RY6 in cancer tissues and normal tissues base on TCGA data, and the normal samples from the GTEx database were added to compare the differences of P2RY6 expression in cancer samples and normal paracancerous samples.

Analysis of genetic alterations

Gene alteration characteristics of P2RY6 were queried by using the cBioPortal network (https://www.cbioportal.org/). After selecting the “Quick Select” as well as “TCGA Pan-Cancer Atlas Studies” options, information about the frequency of alterations, copy number alteration, and their mutation type were described in the “Cancer Types Summary” section of the TCGA database.

Immunomodulators

To elucidate tumor-immune system interactions, we searched the online comprehensive database TISIDB (http://cis.hku.hk/TISIDB/) for immunosuppressants and immunostimulants significantly associated with P2RY6-related gene expression (Spearman correlation test, p < 0.05).

Survival analysis

Kaplan–Meier analysis was performed to show the survival of patients in the high and low P2RY6 expression groups and to assess the potential of P2RY6 as a prognostic marker. Survivor, glmnet [27], survivor ROC, and survminer were the R packages used in the operations. To explore the relationship between P2RY6 expression and disease-specific survival (DSS), disease-free survival (DFS), overall survival (OS), and progression-free survival (PFS) in patients with each cancer type, we performed Cox regression analysis using TCGA in an R environment.

Immune infiltration analysis

We explored the correlation between P2RY6 expression and immune infiltration by using the TIMER2 database and the CIBERSORT method. P2RY6 expression was also analyzed at the following cellular abundances, including B cells, CD4 + T cells, CD8 + T cells, neutrophils, dendritic cells (DCs), and macrophages. We compared their expression with P2RY6 after obtaining stromal score, ESTIMATE score, and immune score for each tumor type sample in TCGA by using ESTIMATE and summarized their relationship.

Genome enrichment analysis

We identified the pathways associated with P2RY6 by gene set enrichment analysis (GSEA) by the cluster Profiler R package. Fold changes in mean gene expression between patients with high and low P2RY6 expression were sorted, and the sorted gene list was then represented in the input file, and assessed pathways by using HALLMARK pathways and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways.

Construction of nomogram

After combining the clinical characteristics and risk scores of the patients, we made a nomogram of the cancer prognosis, which was created by the rms package of the R software. In order to predict the deviation of the probability from the actual situation, we display it by applying the bootstrap method (1000 repetitions) and using a calibration curve.

Statistical analysis

P2RY6 high and low expression groups were defined according to TPM expression values. We used Spearman correlation test to mine the relationship between P2RY6 expression and targets such as tumor mutation burden and microsatellite instability. We calculated log-rank p values and HR values in survival analysis by COX regression analysis. p values were considered statistically significant when they were less than 0.05.

Immunohistochemical staining

Tissue microarrays were prepared using tumor and paracancerous samples from surgical specimens. An experienced pathologist selected typical wax blocks of cancer and paraneoplastic tissues and marked representative lesions according to the sections. Holes were then punched into the wax blocks using a 1 mm diameter and 4 mm long tissue microarrayer. The tissue from the perforated wax blocks was accurately placed into small holes of empty white wax blocks. The resulting tissue microarray was subjected to baking and bleaching using a tissue section machine until all tissue specimens were planted in the empty white wax blocks. Finally, the tissue microarrays were baked at 70 °C, dewaxed in xylene, and hydrated in gradient ethanol. The routinely dewaxed tissue microarrays were subjected to high-temperature antigen repair at 100 °C for 2.5 min, with an additional 15 min of boiling while maintaining the fire. Immunohistochemical staining was performed according to the kit instructions (CK0062-20 T, Signalway Antibody, USA). Subsequently, the sections were incubated with 30% hydrogen peroxide for 15 min and immersed in 5% bovine serum albumin for 30 min at 37 °C to inhibit and block endogenous peroxidase activity. The tissue microarrays were then incubated overnight at 4 °C with anti-P2RY6 (dilution 1:300, Signalway Antibody, USA). After re-warming for 30 min, the sections were rinsed three times in phosphate-buffered solution (PBS), titrated with murine rabbit general secondary antibody, and stained with DAB. The sections were then restained with hematoxylin and sealed.

Cell communication analysis

The identifying and illustrating alterations in the intercellular signaling network (iTALK) R package is a novel tool for intercellular communication analysis based on scRNA-seq [28], which captures high levels of down- or upregulated ligand-receptor gene pairs. We applied iTALK to analyze the differences in communication between cells with high versus low P2RY6 expression (cut-off TPM > 1). Differential analysis was performed using MAST [29]. Ligands were classified into four groups, namely cytokines/chemokines, immune checkpoint genes, and growth factors.

Weighted gene co-expression network analysis

WGCNA is a method used to analyze the expression patterns of genes in multiple samples, which focuses on grouping genes with highly similar expression patterns. It focuses on grouping genes with highly similar expression patterns into modules and further analyzes the intrinsic associations between the module or key genes in the module and the features [30]. Genes with a mean expression value (TPM) greater than 1 are selected to construct the WGCNA network. The weighted gene co-expression network requires the selection of a soft thresholding power to make the relationships between genes in the co-expression network conform to the scale-free network distribution, making the correlation values more consistent with scale-free network characteristics and more biologically meaningful; the “pickSoftThreshold” function was applied for the soft threshold beta screening process. After obtaining the expression matrix and the best beta value, the “blockwiseModules” function of the WGCNA package was used for co-expression matrix construction. The association between gene modules and clinical features was judged by module significance (MS), which is the mean value of the significance of all genes in a module, and the higher the value, the stronger the association between the module and clinical features [31]. The higher the MS value, the stronger the association between the module and the clinical features, which can be used as a significant module for in-depth analysis.

Results

Increased expression of P2RY6 in pan-cancer

Based on the data obtained from the pan-cancer database, we conducted an evaluation of the mRNA expression level of P2RY6. By comparing the expression level of P2RY6 in cancer and paracancerous tissues, we observed a significant difference between the two. The expression of P2RY6 mRNA was found to be significantly higher in most tumor tissues compared to normal tissues in various types of cancers including BLCA (bladder urothelial carcinoma), ESCA (esophageal carcinoma), CHOL (cholangiocarcinoma), BRCA (breast invasive carcinoma), LIHC (liver hepatocellular carcinoma), HNSC (head and neck squamous cell carcinoma), KIRC (kidney renal clear cell carcinoma), LUAD (lung adenocarcinoma), LUSC (lung squamous cell carcinoma), READ (rectum adenocarcinoma), THCA (thyroid carcinoma), STAD (stomach adenocarcinoma), and UCEC (uterine corpus endometrial carcinoma). The only exception was KICH, where P2RY6 expression showed a significant decrease. Considering the limited number of normal samples in the TCGA database, we decided to further analyze the expression level of P2RY6 in different types of cancers using TCGA’s tumor tissue data and the normal tissue data from the GTEx database. Our analysis included 27 types of tumors. The results revealed that apart from the aforementioned tumors, P2RY6 was also highly expressed in CESC (cervical squamous cell carcinoma and endocervical adenocarcinoma), GBM (glioblastoma multiforme), COAD (colon adenocarcinoma), LAML (acute myeloid leukemia), PAAD (pancreatic adenocarcinoma), OV (ovarian serous cystadenocarcinoma), LGG (brain lower grade glioma), PRAD (prostate adenocarcinoma), SKCM (skin cutaneous melanoma), and TGCT (testicular germ cell tumors) (Fig. 1B, all p < 0.01). These findings suggest that the differential expression of P2RY6 may be associated with the occurrence and development of various types of tumors, and it could potentially promote tumor growth, playing a significant role in tumorigenesis.

Fig. 1
figure 1

P2RY6 expression level in different human cancers. A The pan-cancer expression level of P2RY6 from TCGA database. Tumor tissues are represented by yellow spindles, and normal tissues are represented by blue spindles (*, p < 0.05; **, p < 0.01; and ***, p < 0.001). B The pan-cancer expression level of P2RY6 from TCGA and GTEx databases. Tumor tissues are represented by yellow spindles, and normal tissues are represented by blue spindles (*, p < 0.05; **, p < 0.01; and ***, p < 0.001)

Full size image

The expression of P2RY6 is proportional to the stage of tumor

Using the “pathological staging” module of the World Health Organization, we conducted an analysis to examine the correlation between the expression of P2RY6 and the pathological stage of cancer. Our findings revealed a significant increase in the expression level of P2RY6 in stages III–IV of BLCA (bladder urothelial carcinoma) and THCA (thyroid carcinoma) (Supplementary Fig. 1). Interestingly, although there was an increasing trend in the expression of P2RY6 with higher tumor stages, such as KICH (kidney chromophobe), KIRC (kidney renal clear cell carcinoma), LUSC (lung squamous cell carcinoma), KIRP (kidney renal papillary cell carcinoma), MESO (mesothelioma), LUAD (lung adenocarcinoma), SKCM (skin cutaneous melanoma), TGCT (testicular germ cell tumors), and UVM (uveal melanoma), STAD (stomach adenocarcinoma), the differences were not statistically significant (Supplementary Fig. 1). Although there appears to be a correlation between P2RY6 expression and tumor stage, further evidence is required to clarify the relationship between P2RY6 and tumor progression.

P2RY6 has clinical prognostic significance

We utilized the Kaplan–Meier plotter database and conducted univariate Cox regression analysis to investigate the impact of P2RY6 expression on cancer patient survival. The results from Kaplan–Meier analysis revealed a significant correlation between high P2RY6 expression and a decline in overall survival (OS) for various cancer types, including CESC, BRCA, GBM, COAD, ESCA, KIRC, HNSC, KIRP, LIHC, LGG, LUAD, PAAD, PCPG (pheochromocytoma and paraganglioma), TGCT, THYM (thymoma), STAD, UCEC, and UVM (Supplementary Fig. 2, all p < 0.05). Additionally, we observed a significant association between decreased progression-free interval (PFI) and high P2RY6 expression in BLCA, DLBC, ESCA, GBM, HNSC, KIRC, LUSC, LGG, LUAD, PAAD, PRAD, TYM, STAD, UCEC, and UVM (Supplementary Fig. 3, all p < 0.05). Furthermore, high P2RY6 expression was found to be significantly associated with disease-specific survival (DSS) in BRCA, BLCA, COAD, HNSC, GBM, ESCA, KIRC, LUSC, KIRP, LGG, LIHC, LUAD, MESO, PAAD, TCGT, STAD, THYM, UCEC, and UVM (Supplementary Fig. 4, all p < 0.05). The results of Cox regression analysis indicated that elevated P2RY6 expression was closely related to poor OS, PFI, and DSS, suggesting that it may serve as a risk factor for most tumor patients (Fig. 2, all p < 0.05). Based on the aforementioned findings, we speculate that P2RY6 expression is an important factor influencing cancer patient survival, and the high expression of P2RY6 may be closely associated with unfavorable prognosis in various cancers.

Fig. 2
figure 2

Forest plot for univariate Cox regression analysis of P2RY6. The univariate cox regression analysis results of P2RY6 on OS of TGCA pan-cancer patients are shown in A. The univariate cox regression analysis results of P2RY6 on the PFI of TGCA pan-cancer patients are shown in B. The results of univariate cox regression analysis of P2RY6 on DSS of TGCA pan-cancer patients are shown in C. Red represents a significant result

Full size image

Gene set enrichment analysis of P2RY6

To investigate the potential molecular functions of P2RY6 in various cancers, we performed gene cluster enrichment analysis. The results of the KEGG signaling pathway enrichment analysis indicated that the high expression group of P2RY6 was significantly enriched in pathways such as the IL-17 signaling pathway, cytokine-cytokine receptor interaction, inflammatory bowel disease, toll-like receptor signaling pathway, and NF-kappa B signaling pathway (Fig. 3A–B, all p < 0.01). Conversely, the low expression group was significantly enriched in pathways related to bile secretion, drug metabolism-cytochrome P450, cholesterol metabolism, steroid hormone biosynthesis, and retinol metabolism. These findings suggest that high P2RY6 expression is closely associated with various immune-related signaling pathways. Additionally, we verified these results using HALLMARK pathway analysis. It was observed that the high expression group of P2RY6 was significantly enriched in pathways such as allograft rejection, TNF-α signaling via NF-kB, interferon-alpha response, inflammatory response, and IL6/JAK/STAT3 signaling (Fig. 3C–D, all p < 0.01). On the other hand, the low expression group was significantly enriched in pathways related to adipogenesis, fatty acid metabolism, androgen response, xenobiotic metabolism, bile acid metabolism, and others. These findings further support the association between P2RY6 and immune-related signaling pathways, suggesting that P2RY6 may play a role in the immune response during the development and progression of different cancers.

Fig. 3
figure 3

The results of gene set enrichment analysis of P2RY6 high expression group and low expression group through different ways. A The enrichment analysis results of the P2RY6 high expression group in the KEGG pathway. B The enrichment analysis results of the P2RY6 low expression group in the KEGG pathway. C The enrichment analysis results of the P2RY6 high expression group in the HALLMARK pathway. D The enrichment analysis results of the P2RY6 low expression group in the HALLMARK pathway. p < 0.05 was considered significant

Full size image

Analysis of immune cell infiltration and related immune regulatory genes and chemokines in P2RY6

Considering the growing recognition of the correlation between the tumor microenvironment, especially immunological characteristics, and the prognosis of various cancers, we conducted an analysis using the “CIBERSORT” method to assess the correlation of 22 immune cells. Our findings revealed a significant relationship between the expression of P2RY6 and the infiltration levels of CD4 memory resting T cells, Treg cells, and M2 macrophages in these immune cells across most tumors (Supplementary Fig. 5, all p < 0.01). These findings indicate that P2RY6 might be involved in the immune response of different cancers by promoting the infiltration of specific immune cells, which could potentially explain its cancer-promoting effect in most tumor types. Additionally, we explored the immune effects of P2RY6 in cancer from other perspectives. We analyzed the relationship between P2RY6 and immunosuppressive genes, immune-activating genes, chemokines, and chemokine receptors. Our analysis revealed a positive association between the expression of P2RY6 and immunosuppressive genes such as CD86, CD80, IL2RA, TGFB1, and LGALS9 in pan-cancer (Supplementary Fig. 6). Moreover, P2RY6 expression showed a positive correlation with chemokines and chemokine receptors such as CCL5, CCL8, and their corresponding receptors CCR5 and CCR8 (Supplementary Fig. 6 C–D). Furthermore, we observed a strong correlation between the expression of P2RY6 and StromalScore, ImmuneScore, and ESTIMATEScore in most tumors, demonstrating a significant positive association (Supplementary Fig. 7, 8 and 9, all p < 0.05). Based on these results, we speculate that P2RY6 might play a role in the chemotactic recruitment of immune cells and the regulation of immune signaling pathways in the tumor immune response.

Correlation between P2RY6 and immune checkpoint genes, TMB, and MSI in human cancer

Immune checkpoint genes play a crucial role in tumor immunotherapy. To further investigate the potential of P2RY6 in clinical immunotherapy, we analyzed its correlation with immune checkpoint genes. Our analysis encompassed 47 common immune checkpoint genes, and we observed a strong correlation between the expression of P2RY6 and these genes in various types of cancers (Supplementary Fig. 10). This indicates that P2RY6 might coordinate the activity of immune checkpoint genes in different signal transduction pathways, making it a promising target for immunotherapy. Furthermore, we conducted a Spearman correlation test to explore the relationship between P2RY6, tumor mutational burden (TMB), and microsatellite instability (MSI). The results demonstrated a positive correlation between P2RY6 and TMB in BLCA, BRCA, GBM, ESCA, KIRP, LGG, LUAD, SARC, and THYM. Additionally, P2RY6 exhibited a positive correlation with MSI in BLCA, PAAD, COAD, BRCA, DLBC, KIRC, READ, HNSC, SKCM, STAD, and TGCT (Supplementary Fig. 11). These findings shed light on the role of P2RY6 in the immune response and immune mechanisms within the tumor microenvironment. The positive associations with TMB and MSI suggest that P2RY6 may influence the immune response in these specific cancer types. These results deepen our understanding of the potential implications of P2RY6 in the context of tumor immunotherapy.

Expression level of P2RY6 in LUAD and its relationship with immune cells

The data analysis revealed a significant difference in the expression level of P2RY6 in LUAD compared to other cancer tissues, with LUAD showing notably higher expression. This higher expression of P2RY6 in LUAD suggests its potential involvement in promoting the development of LUAD and its association with prognosis. Kaplan–Meier analysis and univariate Cox regression analysis further supported this observation, demonstrating a significant correlation between high P2RY6 expression and worsened overall survival (OS), progression-free interval (PFI), and disease-specific survival (DSS) in LUAD patients. To delve deeper into the role of P2RY6 in LUAD, we examined its protein expression levels using immunohistochemical staining. The staining of LUAD tissues showed that P2RY6 was predominantly located in the cell membrane and cytoplasm, with some cells displaying nuclear staining. Additionally, P2RY6 staining appeared darker in tumors compared to adjacent para-cancerous tissues (Fig. 4A). Moreover, based on corresponding immunohistochemical scores, the expression level of P2RY6 in LUAD was significantly higher than that in para-cancerous tissues (Fig. 4B, p = 0.0081). These experimental results align with the findings from the database analysis, reinforcing the notion that P2RY6 plays a pivotal role in LUAD. Further investigation into the signal pathways regulated by P2RY6 in the LUAD immune response revealed a range of immunomodulators, including IDO1, LGALS9, IL10, TIGIT, IL10RB, TNFRSF18, CD40, and IL2RA, which exhibited significant positive correlations with P2RY6 in LUAD (Supplementary Fig. 12, p < 0.01). Additionally, by utilizing the TCGA database, we assessed the level of immune cell infiltration in LUAD according to P2RY6 expression. Notably, the high P2RY6 expression group displayed lower levels of resting CD4 T cell infiltration (Supplementary Fig. 13A, p < 0.01), which is associated with poor OS in LUAD. Furthermore, the abundance of Treg cells was significantly increased in the high P2RY6 expression group (Supplementary Fig. 13B, p < 0.05). These findings collectively indicate that P2RY6 plays a crucial role in immune cell infiltration in LUAD and exerts varying degrees of influence on tumor occurrence and progression.

Fig. 4
figure 4

Analysis of the protein expression level of P2RY6. A Representative images of immunohistochemical detection of P2RY6 expression in LUAD. B Box plot of immunohistochemical scores related to LUAD and adjacent tissues

Full size image

Prognostic value and risk model of P2RY6-related immunomodulators in LUAD

The study employed univariate Cox regression analysis to investigate the prognostic significance of P2RY6-related immunomodulators in LUAD. Several genes, including PVRL2, TGFBR1, CD276, CD70, MICB, NT5E, and RAET1E, exhibited significant associations with P2RY6 (Supplementary Fig. 14A, all p < 0.05). Subsequently, multivariate Cox regression analysis was conducted on the genes that demonstrated statistical significance in the univariate analysis. This analysis identified several genes, namely TGFBR1, CD70, CD276, NT5E, MICB, and RAET1E, which were closely correlated with poor prognosis in LUAD (Supplementary Fig. 14B, all p < 0.05). Based on these findings, a related prognostic risk model was constructed by calculating the gene expression values and coefficients, which were then used to calculate the risk scores for patients with LUAD. The correlation between these genes and overall survival (OS) was further evaluated using the Kaplan–Meier plotter database. The results revealed that patients with a high-risk score exhibited significantly lower survival time and worse OS (Supplementary Fig. 14C, p < 0.01). Moreover, the corresponding receiver operating characteristic (ROC) curve analysis demonstrated that the area under the curve (AUC) for the risk score was 0.807, while the AUC for the TNM staging was 0.683. Importantly, when combining the two, the AUC value increased to 0.83, indicating that the P2RY6-based prognostic model possessed superior discrimination capabilities compared to TNM staging (Supplementary Fig. 14D).

Prognostic value of P2RY6-related risk models in LUAD cohort

To explore the accuracy and prognostic value of the risk model, the study utilized correlation heat maps and conducted Cox regression analysis. Figure 5A depicted the distribution of risk scores, survival time, and gene expression profiles in LUAD. The chart clearly demonstrated that higher risk scores were associated with shorter survival time and higher mortality rates. Notably, it was evident from the figure that the expression of the P2RY6 gene increased as the risk score increased. In both univariate and multivariate Cox regression analyses shown in Fig. 5B (all p < 0.05), the risk score was found to be significantly correlated with the survival rate of LUAD patients. It emerged as an independent prognostic factor, exerting a stronger impact on the survival rate compared to age and stage. Therefore, the risk score based on the P2RY6 gene exhibited favorable prognostic performance and could serve as an independent prognostic index for LUAD.

Fig. 5
figure 5

Analysis of the prognostic value of the risk scoring model in LUAD patients. A The risk score distribution, patient’s survival status, and related heat map of gene expression profile for LUAD. B The results of univariate and multivariate Cox regression analysis on the overall survival risk score of LUAD

Full size image

Construction of P2RY6-related nomogram

Finally, we constructed a prognostic nomogram for LUAD (Fig. 6A), which could be used to predict the 3-year and 5-year survival probability of related individuals. From the calibration curve, it could be seen that the predicted survival probability of the nomogram for 3 and 5 years was in good agreement with the ideal reference line, which was very close to the ideal 45 degree (Fig. 6B–C).

Fig. 6
figure 6

The establishment of a nomogram related to the risk score is used for the prognosis of LUAD patients. A A nomogram that can be used to predict the 3-year and 5-year survival rates of LUAD patients. B A calibration curve used to calibrate the 3-year survival rate of the nomogram related to LUAD patients. C A calibration curve used to calibrate the 5-year survival rate of the nomogram related to LUAD patients

Full size image

Single-cell transcriptome analysis of the P2RY6 in the LUAD tumor microenvironment

LUAD scRNA-seq data from 11 primary foci was obtained from the GEO database. Cell clustering analysis based on the t-distributed stochastic neighbor embedding (tSNE) algorithm showed that the above cells could be classified into eight clusters, namely T lymphocytes, endothelial cells, B lymphocytes, epithelial cells, NK cells, MAST cells, myeloid cells, and fibroblasts (Fig. 7A). Additionally, we found that all eight of these cells were present in samples from different sources (Fig. 7B). The expression of P2RY6 in the LUAD tumor microenvironment is cell-specific, being highest in myeloid cells and endothelial cells and almost absent in NK cells and MAST cells (Fig. 7C and D). It should be noted that this phenomenon can also be caused by bias in the inclusion of cells in scRNA-seq or by the lack of sensitivity of the sequencing method. Overall, the P2RY6 differed significantly among the different LUAD cell clusters, suggesting that P2RY6 may be a characteristic of LUAD cells.

Fig. 7
figure 7

Single-cell transcriptomic atlas of LUAD. A A tSEN plot of LUAD samples with eight distinct cell types. B A tSEN plot of LUAD from 11 different samples. C A tSEN plot of P2RY6 in different cell types. D Comparison of P2RY6 in different tumor microenvironment cells. The blue horizontal line on the violin plot indicates the median of P2RY6 expression. The letters at the top indicate whether there is a statistical difference between the two cells. Different letters indicate that the difference is statistically significant

Full size image

To further explore the mechanisms of P2RY6 in the LUAD tumor microenvironment, we performed a cellular communication analysis. We first depict the top 100 receptor-ligand pairs with the highest abundance. We found close cellular communication between T lymphocytes, NK cells, and fibroblasts (Fig. 8A). CD3D was the main receptor responsible for cellular communication in T lymphocytes. KLRD1 was the main receptor responsible for cellular communication in NK cells. The B2M and HLA gene families were important ligands responsible for cellular communication in the cells mentioned above (Fig. 8B). Our analysis of cellular communication between high and low P2RY6 expressing cells identified a total of 102 statistically significant receptor-ligand pairs (Fig. 8B). Interestingly, we found that most of the ligands were significantly upregulated in the P2RY6 high expression group. For example, the ligand COL6A1, which is predominantly upregulated in epithelial cells (malignant cells), has been reported to promote bone metastasis in lung cancer [32]. The above results partly explain the potential mechanism for the poor prognosis of LUAD patients with high P2RY6 expression.

Fig. 8
figure 8

Cell communication analysis and WGCNA analysis of LUAD. A Network plot showing the number of ligand-receptor interactions detected in LUAD. B Circos plots showing top 100 highly expressed ligand-receptor interactions (upper panel). Circos plot showing significantly different ligand-receptor interactions between cells high and low P2RY6 expression (lower panel). C Clustering dendrograms of genes. D Gene module–trait associations. Each row corresponds to an ME. The corresponding correlation and p-value are shown in the cells. The cells are color coded by the correlation according to the color legend. E Significantly enriched KEGG pathway in model blue

Full size image

To clarify the gene modules regulated by P2RY6, we used WGCNA to identify P2RY6-related gene modules. The gene network was constructed using the “pickSoftThreshold” function to select β = 4 at “SFT.R.sq” of 0.90. The gene clustering results were further sheared according to the criteria of hybrid dynamic shearing tree by setting the lower limit of the number of gene modules to 30 to obtain different gene modules. After determining the gene modules, the eigenvector value ME was then calculated for each module. The ME represents the average level of gene expression values within each module. All modules were then re-clustered using the “cutreeDynamic” function to merge similar modules into new ones. The modules with high similarity were merged to obtain a total of 8 gene modules (Fig. 8C). The black module contained 42 genes, the blue module contained 1635 genes, the brown module contained 462 genes, the green module contained 294 genes, the grey module contained 1212 genes, the red module contained 72 genes, the turquoise module contained 1736 genes, and the yellow module contained 392 genes. It is important to note that the grey module is a mixed module; i.e., the genes contained in it cannot be combined into one module. We next calculated the gene modules associated with P2RY6 expression using Pearson correlation analysis of the correlation coefficients and corresponding p-values between the eigenvector gene ME and clinical features in the modules to identify significant modules. The results showed that the blue module was significantly positively correlated with P2RY6 expression (p = 8 × 10−32, Fig. 8D). We performed a pathway enrichment analysis of the genes in the blue module and the results showed that the blue module was significantly enriched in the pathways “Antigen processing and presentation” and Th1 and Th2 cell differentiation. The blue module was significantly enriched in “Antigen processing and presentation,” “Th1 and Th2 cell differentiation,” and other pathways (Fig. 8E). The above pathways are closely related to tumor immunity, suggesting that P2RY6 may serve as a potential target.

Discussion

In our study, we first analyzed the expression level of P2RY6 from the perspective of pan-cancer, and found that compared with normal tissues, the expression level of P2RY6 in a variety of cancer tissues was significantly higher, and the results of Kaplan–Meier analysis and univariate Cox regression analysis showed that the poor prognosis of patients was closely correlated with the high expression of P2RY6, so we concluded that the high expression of P2RY6 in most tumors may not be conducive to the survival of patients with related tumors. Previous studies have pointed out that the progression of the disease is closely correlated with the expression of P2RY6, and there is a strong correlation between different levels of P2RY6 expression and different prognosis. For example, Wan et al. pointed out that P2RY6 is closely correlated with the prognosis and treatment of gastric cancer [33]. Qing et al. pointed out that the high expression of P2RY6 can be used as a prognostic indicator of pancreatic cancer [34]. Therefore, after analyzing several databases including GEPIA, Timer, Kaplan–Meier plotter, and TCGA, we speculate that P2RY6 may be a cancer-promoting gene, and its expression may be a valuable prognostic biomarker in cancer.

To further investigate the potential mechanisms of P2RY6 in various cancers, we analyzed the enrichment of P2RY6-related signaling pathways using the TCGA database and found significant correlations with multiple immune-related pathways. Our work also demonstrated a close relationship between P2RY6 gene expression and the degree of immune cell infiltration. Previous studies have shown that tumor-infiltrating lymphocytes can serve as a novel prognostic index [35] and be utilized for the treatment of solid tumors [36]. Interestingly, we found that P2RY6 is closely related to the immunity of many cancers, especially LUAD. Combined with the previously observed low infiltration of CD4 memory resting T cells and high infiltration of Treg cells in LUAD, we speculate that immune cell infiltration may play an extraordinary role in the formation and treatment of LUAD. It has been proved that these cells play an important role in innate immunity and participate in the function of immune cells [37,38,39]. Richard Kennedy et al. pointed out that CD4 + T memory cells can enhance anti-tumor response by enhancing clonal expansion of tumor site, preventing activation-induced cell death, and giving priority to the generation of immune memory cells by CTL as antigen presenting cells [40], which is consistent with the low infiltration level of CD4 memory resting T cells caused by high expression of P2RY6 in our study. In the process of tumorigenesis and development, Treg cells can inhibit the anti-tumor response, so a large number of Treg cells infiltrate into tumor tissues, which is often closely related to poor prognosis [38]. From these, we can also explain why the increased level of Treg cell infiltration in tumor microenvironment often leads to worse prognosis of patients with high expression of P2RY6. In addition, it is worth noting that when we explored the potential relationship between P2RY6 and immunomodulators, we found that the expression of this gene was positively related to a variety of immune activators and immunosuppressants, which means that P2RY6 may have a synergistic effect with other checkpoint members, thus playing a specific role in the regulation of tumor immune microenvironment.

Tumor microenvironment is an important part of tumor tissue, and its composition varies with tumor types, but its specific characteristics include blood vessels, stromal cells, immune cells, and extracellular matrix [41,42,43,44]. With the rapid development of tumor immunotherapy, there is a lot of evidence that tumor microenvironment has high application value in clinical treatment [44,45,46,47]. and now immunotherapy has been clinically selected for the treatment of a variety of tumors (such as renal cell carcinoma [48] and melanoma [49]). After analyzing the 47 immune checkpoint genes collected, we thought that P2RY6 may be a new target for immunotherapy, and after analyzing the relationship between P2RY6, TMB, and MSI, we understand that P2RY6 was positively related to TMB and MSI in most tumors. These results may provide a theoretical basis for combined molecular targeted immunotherapy in the future. Our understanding of the immunological characteristics of P2RY6 microenvironment makes us confident to infer that LUAD may affect tumor immunity and thus affect the prognosis of LUAD patients. Some genes of the same family with P2RY6 have developed corresponding immunotherapy [50]. Inspired by this, we think that P2RY6 may have a good prospect in the immune regulation of tumor microenvironment.

Furthermore, based on the immunomodulators genes related to P2RY6, we successfully constructed the relevant risk model by using Cox regression analysis. Finally, the predicted value of the risk model is verified by the nomogram. Our results show that the risk score is closely correlated with the survival rate of patients with LUAD, and may become an independent prognostic factor with high prognostic value. As a result, patients with early LUAD can be divided into high-risk and low-risk groups with significantly different overall survival rates. The risk scores obtained from P2RY6-related immunomodulators can distinguish the risk groups defined by different gene differential expression, which may be helpful for the further study of LUAD-related prognostic signals [51, 52]. Combining these results, we can find that the risk score is closely related to the immune microenvironment of LUAD, so we are more confident that P2RY6 will be a new target of LUAD immunotherapy and bring results different from traditional treatments in clinical treatment [53, 54].

Our research still has some limitations. Firstly, despite adhering to strict standards, individual variations among different cancer patients may introduce bias to our findings. Secondly, our analysis of P2RY6 using Cox regression may be influenced by other potential prognostic factors such as patients’ socioeconomic status and treatment plans, which are inherent disadvantages of database research. Importantly, this study is based on bioinformatics analysis, and the limitations of public databases may result in a lack of clinical data analysis for LUAD patients. Therefore, before applying our model to clinical practice, further in vivo and in vitro studies are necessary to confirm the interaction between P2RY6 and the LUAD tumor microenvironment.

Conclusions

In summary, our findings suggest that P2RY6 holds promise as a novel prognostic marker for LUAD. We have also explored the potential role of P2RY6 in tumor immunology and its prognostic value in LUAD. Additionally, we developed a nomogram integrating independent clinical factors, which can accurately predict prognosis. Based on our results, we propose that P2RY6 may serve as a target for personalized immunotherapy in LUAD patients. We hope that these findings will provide valuable references for future laboratory research and clinical applications of P2RY6 in cancer prognosis and immunotherapy.

Availability of data and materials

The dataset supporting the conclusions of this article is included within the article. Xena database [23] (https://xenabrowser.net/datapages/), cBioPortal database [24] (https://www.cbioportal.org/), cBioPortal network (https://www.cbioportal.org/), and TISIDB (http://cis.hku.hk/TISIDB/).

Abbreviations

LUAD:

Lung adenocarcinoma

NSCLC:

Non-small cell lung cancer

TCGA:

The Cancer Genome Atlas

OS:

Overall survival

DFS:

Disease-free survival

DSS:

Disease-specific survival

References

  1. Siegel RL, Miller KD, Fuchs HE, Jemal A. Cancer Statistics, 2021. CA Cancer J Clin. 2021;71:7–33.

    Article  PubMed  Google Scholar 

  2. Sher T, Dy GK, Adjei AA. Small cell lung cancer. Mayo Clin Proc. 2008;83:355–67.

    Article  CAS  PubMed  Google Scholar 

  3. Shi J, Hua X, Zhu B, Ravichandran S, Wang M, Nguyen C, Brodie SA, Palleschi A, Alloisio M, Pariscenti G, Jones K, Zhou W, Bouk AJ, Boland J, Hicks B, Risch A, Bennett H, Luke BT, Song L, Duan J, Liu P, Kohno T, Chen Q, Meerzaman D, Marconett C, Laird-Offringa I, Mills I, Caporaso NE, Gail MH, Pesatori AC, Consonni D, Bertazzi PA, Chanock SJ, Landi MT. Somatic genomics and clinical features of lung adenocarcinoma: a retrospective study. PLoS Med. 2016;13:e1002162.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Tayoun T, Faugeroux V, Oulhen M, Déas O, Michels J, Brulle-Soumare L, Cairo S, Scoazec JY, Marty V, Aberlenc A, Planchard D, Remon J, Ponce S, Besse B, Kannouche PL, Judde JG, Pawlikowska P, Farace F. Targeting genome integrity dysfunctions impedes metastatic potency in non-small cell lung cancer circulating tumor cell-derived explants. JCI insight. 2022;7:e155804.

    Article  PubMed  PubMed Central  Google Scholar 

  5. Shi J, Chen Y, Peng C, Kuang L, Zhang Z, Li Y, et al. Advances in targeted therapy against driver mutations and epigenetic alterations in non-small cell lung cancer. Oncologie. 2022;24:613–48.

    Article  Google Scholar 

  6. Zappa C, Mousa SA. Non-small cell lung cancer: current treatment and future advances. Transl Lung Cancer Res. 2016;5:288–300.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Miller KD, Nogueira L, Mariotto AB, Rowland JH, Yabroff KR, Alfano CM, Jemal A, Kramer JL, Siegel RL. Cancer treatment and survivorship statistics, 2019. CA Cancer J Clin. 2019;69:363–85.

    Article  PubMed  Google Scholar 

  8. Calles A, Liao X, Sholl LM, Rodig SJ, Freeman GJ, Butaney M, Lydon C, Dahlberg SE, Hodi FS, Oxnard GR, Jackman DM, Jänne PA. Expression of PD-1 and its ligands, PD-L1 and PD-L2, in smokers and never smokers with KRAS-mutant lung cancer. J Thorac Oncol. 2015;10:1726–35.

    Article  CAS  PubMed  Google Scholar 

  9. Zhou C, Chen G, Huang Y, Zhou J, Lin L, Feng J, Wang Z, Shu Y, Shi J, Hu Y, Wang Q, Cheng Y, Wu F, Chen J, Lin X, Wang Y, Huang J, Cui J, Cao L, Liu Y, Zhang Y, Pan Y, Zhao J, Wang L, Chang J, Chen Q, Ren X, Zhang W, Fan Y, He Z, Fang J, Gu K, Dong X, Zhang T, Shi W, Zou J. Camrelizumab plus carboplatin and pemetrexed versus chemotherapy alone in chemotherapy-naive patients with advanced non-squamous non-small-cell lung cancer (CameL): a randomised, open-label, multicentre, phase 3 trial. Lancet Respir Med. 2021;9:305–14.

    Article  CAS  PubMed  Google Scholar 

  10. Wang H, Fang J, Wang Y, Li S, Wang Z, Li S, et al. Gene editing in non-small cell lung cancer: current application and future perspective. Oncologie. 2022;24:65–83.

    Article  Google Scholar 

  11. Zheng X, Hu Y, Yao C. The paradoxical role of tumor-infiltrating immune cells in lung cancer. Intractable Rare Dis Res. 2017;6:234–41.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Bremnes RM, Busund LT, Kilvær TL, Andersen S, Richardsen E, Paulsen EE, Hald S, Khanehkenari MR, Cooper WA, Kao SC, Dønnem T. The role of tumor-infiltrating lymphocytes in development, progression, and prognosis of non-small cell lung cancer. J Thorac Oncol. 2016;11:789–800.

    Article  PubMed  Google Scholar 

  13. Zuo S, Wei M, Wang S, Dong J, Wei J. Pan-cancer analysis of immune cell infiltration identifies a prognostic immune-cell characteristic score (ICCS) in lung adenocarcinoma. Front Immunol. 2020;11:1218.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Chen Y, Wu H, Jiao A, Tong J, Zhu J, Zhang M, et al. Chinese herbal prescription QYSL prevents progression of lung cancer by targeting tumor microenvironment. Oncologie. 2022;24:295–307.

    Article  Google Scholar 

  15. Molle CM, Arguin G, Jemfer C, Placet M, Dagenais Bellefeuille S, Gendron FP. The expression of the P2Y(6) receptor is regulated at the transcriptional level by p53. Biochem Biophys Res Commun. 2020;524:798–802.

    Article  CAS  PubMed  Google Scholar 

  16. Garcia RA, Yan M, Search D, Zhang R, Carson NL, Ryan CS, Smith-Monroy C, Zheng J, Chen J, Kong Y, Tang H, Hellings SE, Wardwell-Swanson J, Dinchuk JE, Psaltis GC, Gordon DA, Glunz PW, Gargalovic PS. P2Y6 receptor potentiates pro-inflammatory responses in macrophages and exhibits differential roles in atherosclerotic lesion development. PLoS ONE. 2014;9:e111385.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Bar I, Guns PJ, Metallo J, Cammarata D, Wilkin F, Boeynams JM, Bult H, Robaye B. Knockout mice reveal a role for P2Y6 receptor in macrophages, endothelial cells, and vascular smooth muscle cells. Mol Pharmacol. 2008;74:777–84.

    Article  CAS  PubMed  Google Scholar 

  18. Salem M, El Azreq MA, Pelletier J, Robaye B, Aoudjit F, Sévigny J. Exacerbated intestinal inflammation in P2Y(6) deficient mice is associated with Th17 activation. Biochim Biophys Acta Mol Basis Dis. 2019;1865:2595–605.

    Article  CAS  PubMed  Google Scholar 

  19. Ben Yebdri F, Kukulski F, Tremblay A, Sévigny J. Concomitant activation of P2Y(2) and P2Y(6) receptors on monocytes is required for TLR1/2-induced neutrophil migration by regulating IL-8 secretion. Eur J Immunol. 2009;39:2885–94.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Zhang Z, Wang Z, Ren H, Yue M, Huang K, Gu H, Liu M, Du B, Qian M. P2Y(6) agonist uridine 5′-diphosphate promotes host defense against bacterial infection via monocyte chemoattractant protein-1-mediated monocytes/macrophages recruitment. J Immunol (Baltimore, Md: 1950). 2011;186:5376–87.

    Article  CAS  Google Scholar 

  21. He QL, Jiang HX, Zhang XL, Qin SY. Relationship between a 7-mRNA signature of the pancreatic adenocarcinoma microenvironment and patient prognosis (a STROBE-compliant article). Medicine. 2020;99:e21287.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Placet M, Arguin G, Molle CM, Babeu JP, Jones C, Carrier JC, Robaye B, Geha S, Boudreau F, Gendron FP. The G protein-coupled P2Y6 receptor promotes colorectal cancer tumorigenesis by inhibiting apoptosis. Biochim Biophys Acta Mol Basis Dis. 2018;1864:1539–51.

  23. Wang Z, Jensen MA, Zenklusen JC. A practical guide to The Cancer Genome Atlas (TCGA). Methods Mol Biol (Clifton, NJ). 2016;1418:111–41.

    Article  Google Scholar 

  24. GTEx Consortium. The Genotype-Tissue Expression (GTEx) project. Nat Genet. 2013;45:580–5.

    Article  Google Scholar 

  25. Goldman MJ, Craft B, Hastie M, Repečka K, McDade F, Kamath A, Banerjee A, Luo Y, Rogers D, Brooks AN, Zhu J, Haussler D. Visualizing and interpreting cancer genomics data via the Xena platform. Nat Biotechnol. 2020;38:675–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Gao J, Aksoy BA, Dogrusoz U, Dresdner G, Gross B, Sumer SO, Sun Y, Jacobsen A, Sinha R, Larsson E, Cerami E, Sander C, Schultz N. Integrative analysis of complex cancer genomics and clinical profiles using the cBioPortal. Sci Sig. 2013;6:pl1.

    Google Scholar 

  27. Engebretsen S, Bohlin J. Statistical predictions with glmnet. Clin Epigenetics. 2019;11:123.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Wang Y, Wang R, Zhang S, Song S, Jiang C, Han G, Wang M, Ajani J, Futreal A, Wang LJB. iTALK: an R package to characterize and illustrate intercellular communication. 2019. p. 507871.

    Google Scholar 

  29. Finak G, McDavid A, Yajima M, Deng J, Gersuk V, Shalek AK, Slichter CK, Miller HW, McElrath MJ, Prlic M, Linsley PS, Gottardo R. MAST: a flexible statistical framework for assessing transcriptional changes and characterizing heterogeneity in single-cell RNA sequencing data. Genome Biol. 2015;16:278.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Tang X, Huang X, Wang D, Yan R, Lu F, Cheng C, Li Y, Xu J. Identifying gene modules of thyroid cancer associated with pathological stage by weighted gene co-expression network analysis. Gene. 2019;704:142–8.

    Article  CAS  PubMed  Google Scholar 

  32. Li N, Liu M, Cao X, Li W, Li Y, Zhao Z. Identification of differentially expressed genes using microarray analysis and COL6A1 induction of bone metastasis in non-small cell lung cancer. Oncol Lett. 2021;22:693.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Wan H, Xie R, Xu J, He J, Tang B, Liu Q, Wang S, Guo Y, Yang X, Dong TX, Carethers JM, Yang S, Dong H. Anti-proliferative effects of nucleotides on gastric cancer via a novel P2Y6/SOCE/Ca(2+)/β-catenin pathway. Sci Rep. 2017;7:2459.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Künzli BM, Berberat PO, Giese T, Csizmadia E, Kaczmarek E, Baker C, Halaceli I, Büchler MW, Friess H, Robson SC. Upregulation of CD39/NTPDases and P2 receptors in human pancreatic disease. Am J Physiol Gastrointest Liver Physiol. 2007;292:G223–30.

    Article  PubMed  Google Scholar 

  35. Reynders K, De Ruysscher D. Tumor infiltrating lymphocytes in lung cancer: a new prognostic parameter. J Thorac Dis. 2016;8:E833–5.

    Article  PubMed  PubMed Central  Google Scholar 

  36. Badalamenti G, Fanale D, Incorvaia L, Barraco N, Listì A, Maragliano R, Vincenzi B, Calò V, Iovanna JL, Bazan V, Russo A. Role of tumor-infiltrating lymphocytes in patients with solid tumors: can a drop dig a stone? Cell Immunol. 2019;343:103753.

    Article  CAS  PubMed  Google Scholar 

  37. Jaynes JM, Sable R, Ronzetti M, Bautista W, Knotts Z, Abisoye-Ogunniyan A, Li D, Calvo R, Dashnyam M, Singh A, Guerin T, White J, Ravichandran S, Kumar P, Talsania K, Chen V, Ghebremedhin A, Karanam B, Bin Salam A, Amin R, Odzorig T, Aiken T, Nguyen V, Bian Y, Zarif JC, de Groot AE, Mehta M, Fan L, Hu X, Simeonov A, Pate N, Abu-Asab M, Ferrer M, Southall N, Ock CY, Zhao Y, Lopez H, Kozlov S, de Val N, Yates CC, Baljinnyam B, Marugan J, Rudloff U. Mannose receptor (CD206) activation in tumor-associated macrophages enhances adaptive and innate antitumor immune responses. Sci Transl Med. 2020;12:eaax6337.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Tanaka A, Sakaguchi S. Regulatory T cells in cancer immunotherapy. Cell Res. 2017;27:109–18.

    Article  CAS  PubMed  Google Scholar 

  39. Kim HJ, Cantor H. CD4 T-cell subsets and tumor immunity: the helpful and the not-so-helpful. Cancer Immunol Res. 2014;2:91–8.

    Article  CAS  PubMed  Google Scholar 

  40. Kennedy R, Celis E. Multiple roles for CD4+ T cells in anti-tumor immune responses. Immunol Rev. 2008;222:129–44.

    Article  CAS  PubMed  Google Scholar 

  41. Meurette O, Mehlen P. Notch signaling in the tumor microenvironment. Cancer Cell. 2018;34:536–48.

    Article  CAS  PubMed  Google Scholar 

  42. Kalluri R. The biology and function of fibroblasts in cancer. Nat Rev Cancer. 2016;16:582–98.

    Article  CAS  PubMed  Google Scholar 

  43. Hui L, Chen Y. Tumor microenvironment: sanctuary of the devil. Cancer Lett. 2015;368:7–13.

    Article  CAS  PubMed  Google Scholar 

  44. Hinshaw DC, Shevde LA. The tumor microenvironment innately modulates cancer progression. Can Res. 2019;79:4557–66.

    Article  CAS  Google Scholar 

  45. Locy H, de Mey S, de Mey W, De Ridder M, Thielemans K, Maenhout SK. Immunomodulation of the tumor microenvironment: turn foe into friend. Front Immunol. 2018;9:2909.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Jarosz-Biej M, Smolarczyk R, Cichoń T, Kułach N. Tumor microenvironment as a “game changer” in cancer radiotherapy. Int J Mol Sci. 2019;20:3212.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Bader JE, Voss K, Rathmell JC. Targeting metabolism to improve the tumor microenvironment for cancer immunotherapy. Mol Cell. 2020;78:1019–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Lai Y, Zeng T, Liang X, Wu W, Zhong F, Wu W. Cell death-related molecules and biomarkers for renal cell carcinoma targeted therapy. Cancer Cell Int. 2019;19:221.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Namikawa K, Yamazaki N. Targeted therapy and immunotherapy for melanoma in Japan. Curr Treat Options Oncol. 2019;20:7.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Matsuo K, Nishiuma S, Hasegawa Y, Kawabata F, Kitahata K, Nakayama T. Vaccination with antigen combined with αβ-ATP as a vaccine adjuvant enhances antigen-specific antibody production via dendritic cell activation. Biol Pharm Bull. 2016;39:1073–6.

    Article  CAS  PubMed  Google Scholar 

  51. Song C, Guo Z, Yu D, Wang Y, Wang Q, Dong Z, Hu W. A prognostic nomogram combining immune-related gene signature and clinical factors predicts survival in patients with lung adenocarcinoma. Front Oncol. 2020;10:1300.

    Article  PubMed  PubMed Central  Google Scholar 

  52. Zuo S, Wei M, Zhang H, Chen A, Wu J, Wei J, Dong J. A robust six-gene prognostic signature for prediction of both disease-free and overall survival in non-small cell lung cancer. J Transl Med. 2019;17:152.

    Article  PubMed  PubMed Central  Google Scholar 

  53. Shedden K, Taylor JM, Enkemann SA, Tsao MS, Yeatman TJ, Gerald WL, Eschrich S, Jurisica I, Giordano TJ, Misek DE, Chang AC, Zhu CQ, Strumpf D, Hanash S, Shepherd FA, Ding K, Seymour L, Naoki K, Pennell N, Weir B, Verhaak R, Ladd-Acosta C, Golub T, Gruidl M, Sharma A, Szoke J, Zakowski M, Rusch V, Kris M, Viale A, Motoi N, Travis W, Conley B, Seshan VE, Meyerson M, Kuick R, Dobbin KK, Lively T, Jacobson JW, Beer DG. Gene expression-based survival prediction in lung adenocarcinoma: a multi-site, blinded validation study. Nat Med. 2008;14:822–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Sun Z, Wigle DA, Yang P. Non-overlapping and non-cell-type-specific gene expression signatures predict lung cancer survival. J Clin Oncol. 2008;26:877–83.

    Article  PubMed  Google Scholar 

Download references

Acknowledgements

None.

Funding

None.

Author information

Authors and Affiliations

  1. Department of Pneumology, Yiwu Central Hospital, Zhejiang, China

    Hong Wu & Xuhui Dong

Authors
  1. Hong Wu
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Xuhui Dong
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Conceived and designed the experiments: HW, XD; Performed the data collection: HW, XD; Analyzed the data: HW, XD; Contributed reagents/materials/analysis tools: HW, XD; Contributed to the writing of the manuscript: HW, XD; All authors reviewed the manuscript.

Corresponding author

Correspondence to Hong Wu.

Ethics declarations

Ethics approval and consent to participate

This study was approved by the Ethics and Human Subject Committee of Yiwu Central Hospital.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: 

Supplementary Figure 1. The expression level of P2RY6 in different stages of tumor. (A–B). The logging scale uses Log2 (TPM+1), and analyzes the expression of P2RY6 gene through the TCGA database according to the main pathological stages (I-II, III-IV).

Additional file 2: Supplementary Figure 2.

Kaplan-Meier survival curve of OS in difference of various tumor.

Additional file 3: Supplementary Figure 3.

Kaplan-Meier survival curve of PFI in difference of various tumor types.

Additional file 4: Supplementary Figure 4

. Kaplan-Meier survival curve of DSS in difference of various tumor types.

Additional file 5: Supplementary Figure 5

. The results of correlation analysis between the expression of P2RY6 in pan-cancer and the impact of immune microenvironment are related to the level of immune cell infiltration. Positive correlation is represented by yellow, and negative correlation is represented by blue, and the darker the color, the stronger the correlation. *, p <0.05, **, p <0.01, and ***, p <0.001.

Additional file 6: Supplementary Figure 6.

Correlation analysis between P2RY6 expression and immune regulation-related genes. (A) Heat map of the correlation between P2RY6 expression in pan-cancer and immune activation genes. (B) The heat map of the correlation between the expression of P2RY6 in pan-cancer and the immunosuppressive state-related genes. (C) Heat map of the correlation between the expression of P2RY6 in pan-cancer and chemokine genes. (D) Heat map of the correlation between the expression of P2RY6 in pan-carcinoma and chemokine receptor genes. Positive correlation is represented by yellow, and negative correlation is represented by blue, and the darker the color, the stronger the correlation. *, p <0.05, **, p <0.01, and ***, p<0.001.

Additional file 7: Supplementary Figure 7.

Results of correlation analysis between P2RY6 expression and Stromal Score in pan-cancer.

Additional file 8: Supplementary Figure 8.

Results of correlation analysis between P2RY6 expression and Immune Score in pan-cancer.

Additional file 9: Supplementary Figure 9.

Results of correlation analysis between P2RY6 expression and ESTIMATE Score in pan-cancer.

Additional file 10: Supplementary Figure 10.

Correlation analysis between P2RY6 expression and pan-cancer immune marker set. *p<0.05; **p <0.01; *** p <0.001.

Additional file 11: Supplementary Figure 11.

The relationship between P2RY6 expression and TMB and MSI in pan-cancer. (A) TMB. (B) MSI. TMB, tumor mutation burden; MSI, microsatellite instability. The radar chart shows Spearman's correlation coefficient and p-value.

Additional file 12: Supplementary Figure 12.

Correlation analysis between P2RY6 gene and immunomodulator.

Additional file 13: Supplementary Figure 13.

Relationship between P2RY6 expression level and tumor infiltrating immune cells in LUAD. (A) The level of resting T cell infiltration in CD4. (B) The infiltration level of Treg cells. *p <0.05; **p <0.01; *** p <0.001.

Additional file 14: 

Supplementary Figure 14. The construction of P2RY6 related risk model. (A) Forest plot results of univariate Cox regression analysis for patients with LUAD. (B) Forest plot results of multivariate Cox regression analysis for patients with LUAD. (C) Kaplan-Meier curve for predicting OS in patients with LUAD based on the TCGA data set. (D) Receiver operating characteristic curve (ROC) analysis results of patients with LUAD.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Wu, H., Dong, X. Immunological role and clinical prognostic significance of P2RY6 in lung adenocarcinoma: a multi-omics studies and single-cell sequencing analysis. World J Surg Onc 21, 341 (2023). https://doi.org/10.1186/s12957-023-03216-1

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12957-023-03216-1

Keywords

World Journal of Surgical Oncology

ISSN: 1477-7819

Contact us