CCR2 and PTPRC are regulators of tumor microenvironment and potential prognostic biomarkers of lung adenocarcinoma.

BACKGROUND: Tumor microenvironment (TME) plays an essential role in lung adenocarcinoma (LUAD) development and metastasis. With the development of TME research, it has been proved that differences in tumor-infiltrating immune cells (TICs) and gene expression profile are related to the prognosis of cancer. The aim of our study was to identify key genes affecting immune state in TME of LUAD.
METHODS: The RNA-seq data and clinical characteristics of 594 LUAD patients were downloaded from the TCGA database. ImmuneScore, StromalScore and ESTIMATEScore of each LUAD sample were calculated using ESTIMATE algorithm. Based on the median of different scores, LUAD samples were divided into high and low score groups. Differentially expressed genes (DEGs) between groups were obtained, and univariate Cox regression analysis and protein-protein interaction (PPI) network were used to screen the shared DEGs generating in the intersection analysis. Finally, the CIBORSORT algorithm was performed to calculate the relative contents of TICs for each LUAD sample, and the correlation analysis between TICs and key genes was used to determine the influence of key genes to the TME.
RESULTS: In the presented study, we found that three different scores were positively correlated with the prognosis of LUAD patients, and correlation analysis showed the different scores were closely related to tumor progression and metastasis. After performing the intersection analysis, a total of 585 up-regulated and 107 down-regulated DEGs between the high and low score groups were obtained, all of which were enriched in immune-related functions. Having used univariate COX regression analysis and PPI network, the key genes, CCR2 and PTPRC, affecting the immune status of TME and the prognosis of LUAD were acquired. Analysis based on the CIBERSORT algorithm suggested that CCR2 and PTPRC were correlated with a variety of TICs, and closely related to the clinical characteristics of the LUAD patients.
CONCLUSIONS: Our research showed that CCR2 and PTPRC may be potential prognostic markers in LUAD, which may affect the function of γδT cells and other immune cells by participating in the regulation of TME immune state. 2021 Annals of Translational Medicine. All rights reserved.

Introduction

Lung cancer is one of the malignant tumors that threaten human health (1), among which adenocarcinoma is the most common histological type of non-small cell lung cancer (NSCLC), accounting for approximately 60% of NSCLC (2). Given the early symptoms of lung cancer are nonspecific, a vast majority of patients are diagnosed at advanced stages of the disease (3). Although great progress has been made in the early diagnosis and treatment of lung cancer in recent years, the exploration of prognostic biomarker and novel therapeutic targets in lung cancer is still unsatisfactory (4,5). Currently, immunotherapy have brought a new direction to the clinical treatment of lung cancer, and a variety of immune checkpoint inhibitors (ICIs) have shown encouraging efficacy in clinical trials (6,7). However, some patients have not benefited from immunotherapy (8) or developed drug resistance during medication (9), and its side effects and adverse reactions are still a worrying issue. Accordingly, it is still necessary to explore and ascertain the pathogenic mechanism of lung adenocarcinoma (LUAD) to break through the current bottleneck.

Tumor microenvironment (TME) is a complex ecosystem composed of various types of cells, which can be divided into immune components dominated by immune cells and matrix components dominated by fibroblasts (10). Most of the cancer characteristics are activated and sustained by varying degrees through contributions from the distinctive cells of TME. They mediate the recruitment, activation, programming and persistence of tumor cells in a variety of ways (11). A growing body of studies have elucidated that tumor-infiltrating immune cells (TICs) can be used as a predictor of cancer prognosis and treatment target (12,13). Nevertheless, Further studies have found that TME components have great heterogeneity among different tumor stages and individuals, thereby leading to the differences in prognosis (14,15). For instance, in liver cancer, tumors with higher transcriptome diversity have lower T cell cytolytic activity, which was related to the poor prognosis (16). Likewise, in colorectal cancer increased M1 macrophage levels indicated shorter overall survival (OS) of patients (17). This heterogeneity not only results in differences in patient outcome, but also increases the difficulty for clinicians to use targeted drugs (18). Therefore, the extensive search for novel and effective TME targeting sites may provide a new window of opportunity for combined immunotherapy (19).

In this study, we used ESTIMATE algorithm to calculate the TME components of LUAD samples in the Cancer Genome Atlas (TCGA) database, and CIBERSORT algorithms were used to calculate the relative contents of various TICs, combining with univariate Cox regression analysis and protein-protein interaction (PPI) network to screen differentially expressed genes (DEGs) and finally obtained predictive biomarkers CCR2 and PTPRC. CCR2 is the receptor for monocyte chemotactic protein-1 (MCP-1), which were expressed by a variety of cell types, including monocytes, dendritic cells (DC), endothelial cells and cancer cells. The upregulation of MCP-1 was related to formation, metastasis and recurrence of multiple types of cancers (20,21). Nevertheless, some studies have indicated that the high expression of CCR2 inhibited the development of small cell lung cancer (SCLC) (22). The protein tyrosine phosphatase CD45, encoded by the PTPRC gene, which has thought to be a regulator of B cell and T cell receptor signaling, it has been reported that PTPRC expression was significantly down-regulated in acute lymphoblastic leukemia and Parkinson's disease (23,24). In patients with SCLC, its expression was positively correlated with the patients’ survival (25). Here, our study revealed that the CCR2 and PTPRC might be a potential indicator for the alteration of TME status in LUAD. We present the following article in accordance with the REMARK reporting checklist (available at https://dx.doi.org/10.21037/atm-21-3301).

Methods

Extraction of data

Through the TCGA database (https://portal.gdc.cancer.gov/), we obtained RNA-seq data and clinical characteristics of 594 LUAD samples (535 tumor samples and 59 healthy samples). The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).

Analyzing of ImmuneScore, StromalScore, and ESTIMATEScore

We analyzed the immune and matrix components in the TME using R language version 4.0.3 loaded with estimate package, presented by ImmuneScore, StromalScore, and ESTIMATEScore. The higher the score reflects the higher corresponding component in TME. ESTIMATEScore represents the sum of immune and matrix components.

Survival analysis

The survival and survminer R package were used for survival analysis. Kaplan-Meier analysis was used to perform the survival curve as well, P<0.05 was considered significant.

Analysis of correlation between clinical characteristics and scores

We analyzed the relationship between clinical characteristics and different scores of LUAD patients downloaded from the TCGA via R language, and used the Wilcoxon rank sum test or the Kruskal-Wallis rank sum test as a test of significance.

Screening of DEGs

LUAD samples were marked as high or low score group based on the median of ImmuneScore and StromalScore. We used the package “limma” to analyze the DEGs between the high and low score groups. DEGs with |log2FC|>1 and false discovery rate (FDR) <0.05 were statistically significant. R language with package pheatmap was used to construct the heatmaps of DEGs.

Enrichment function analysis

To explore the functions of the DEGs, we used clusterProfiler, enrichplot, and ggplot2 R package for the analysis. Terms with P value and q value <0.05 were considered significant.

PPI network construction

The STRING public database (https://string-db.org/) was used to generate the PPI network, and results were reconstructed with Cytoscape version 3.8.1. The nodes with the confidence of interactive association more than 0.7 were selected to generate the PPI network.

Univariate COX regression analysis

R language combined with “survival” package was used in univariate Cox regression analysis. DEGs with a value of P<0.05 were selected for following analysis.

TICs analysis

We used the CIBERSORT algorithm on 535 tumor samples to calculate the relative contents of 22 types of TICs in each LUAD sample, samples with P<0.05 were used for subsequent analysis.

Statistical analysis

Statistical analysis was performed using with R (version 4.0.3, https://www.r-project.org/). P value <0.05 was considered to be statistically significant.

Results

Survival analysis of LUAD patients in three different scores

We combined ESTIMATE algorithms and Kaplan-Meier survival analysis to profile the relationship between different scores and patients’ prognosis. ImmuneScore, StromalScore and ESTIMATEScore showed a positive correlation with OS (). The results above suggested that the immune and matrix components in the TME could reflect the prognosis of LUAD patients, indicating the potential value of each component of TME.

Figure 1

Correlation of different scores with the survival time of LUAD patients. (A) Kaplan-Meier survival analysis for LUAD patients grouped by median of ImmuneScore; (B) Kaplan-Meier survival analysis for LUAD patients grouped by median of StromalScore; (C) Kaplan-Meier survival analysis for LUAD patients grouped by median of ESTIMATEScore. LUAD, lung adenocarcinoma.

Correlation of scores with clinical characteristics of LUAD patients

Clinical characteristics were downloaded from the TCGA database, and the correlation analysis between three types of scores and clinical characteristics was carried out (). As shown in , males had lower scores relative to females. Meanwhile, three kinds of scores were negatively correlated with tumor stages (). Further analysis demonstrated that scores decreased significantly from T1 to T4 in T classification (), M classification also showed a negative correlation with scores. However, no statistically significant changes of scores were observed in N classification. In a word, results above demonstrated that the immune and matrix components in the TME correlated with tumor progression and metastasis.

Figure 2

Correlation analysis between scores and clinical characteristics. (A-C) Correlation analysis between scores and gender; (D-F) correlation analysis between scores and tumor stage; (G-I) correlation analysis between scores and T classification; (J-L) correlation analysis between scores and N classification; (M-O) correlation analysis between scores and M classification.

DEGs screening

In order to analyze the changes in the gene expression profile of immune and matrix components, we divided ImmuneScore and StromalScore into high and low score groups according to the median of scores, and the DEGs of the two groups were presented in heat maps (), thus identifying 1,394 genes from ImmuneScore (including 1,097 up-regulated and 297 down-regulated genes), and 1,623 genes from StromalScore (including 1,422 up-regulated and 201 down-regulated genes). The intersection analysis of the two scores obtained 585 up-regulated and 107 down-regulated shared genes (). These shared DEGs may exert dual effect on the immune and matrix component of TME.

Figure 3

Heatmaps, Venn plots, and functional enrichment analyses for DEGs. (A) Heatmap for DEGs generated in ImmuneScore; (B) Heatmap for DEGs generated in StromalScore; (C,D) Venn plots showing shared up-regulated and down-regulated DEGs in ImmuneScore and StromalScore; (E,F) GO and KEGG enrichment analyses for shared DEGs in ImmuneScore and StromalScore. DEG, differentially expressed gene; GO, gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Functional enrichment analysis

Subsequently, we performed functional enrichment analysis for the shared DEGs, the results of gene ontology (GO) analysis showed that the functions of 692 shared DEGs in two groups were mainly concentrated on immune-related terms, including immune response-activating cell surface receptor signaling pathway, immune response-activating signal transduction, lymphocyte mediated immunity, and humoral immune response (). Meanwhile, Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis indicated DEGs primarily enriched in cytokine-cytokine receptor interaction, chemokine signaling pathway, hematopoietic cell lineage and B cell receptor signaling pathway (). To sum up, DEGs obtained from the intersection analysis of ImmuneScore and StromalScore were enriched in immune function of TME, which may contain the key factors affecting immune state.

Intersection analysis of PPI network and univariate Cox regression analysis

We then constructed PPI network using STRING database to further screen the 692 shared genes. The top 30 genes sorted by the number of nodes were shown in bar plot (). Univariate Cox regression analysis of the 692 shared genes were performed showing that 22 genes were related to the OS of LUAD (). Having performed the intersection analysis of PPI network and univariate Cox regression analysis, we ultimately obtained two overlapping genes CCR2 and PTPRC ().

Figure 4

PPI network and univariate Cox regression analysis. (A) Construction of PPI network based on STRING database; (B) the top 30 genes ordered by the number of nodes; (C) forest plot of univariate Cox regression analysis showing the genes affecting the prognosis of LUAD; (D) Venn plot showing the shared genes between the top 30 genes in PPI and significant factors in univariate COX regression analysis. PPI, protein-protein interaction; LUAD, lung adenocarcinoma.

Integrated analysis of CCR2 and PTPRC

Next, we performed integrated analysis of these two genes. As shown in , the expression of CCR2 and PTPRC in tumors were significantly lower than that of normal samples. Survival analysis showed that patients in the high expression group of the two key genes had longer OS than those in the low expression group (). Next, we further analyzed the correlation between the two genes and clinical characteristics (), manifesting that CCR2 expression was negatively correlated with tumor stage of LUAD (), and its expression decreased as the T and N classification increased (). Similar results were presented in the analysis of PTPRC (). In conclusion, we can conclude that these two genes may affect the prognosis by affecting the immune state in TME based on the previous analysis, and they could be used as reliable biomarkers of LUAD prognosis.

Figure 5

Integrated analysis of CCR2 and PTPRC. (A,B) Different expression of CCR2 and PTPRC in the normal and cancer tissue; (C,D) Kaplan-Meier survival analysis of LUAD between different CCR2 and PTPRC expression level and prognosis. LUAD, lung adenocarcinoma.

Figure 6

The correlation analysis between CCR2 and PTPRC expression and clinic characteristics. (A-E) Correlation analysis between CCR2 expression level and different clinic characteristics; (F-J) correlation analysis between PTPRC expression level and different clinic characteristics.

Correlation analysis between key genes and TICs

The CIBERSORT algorithm was performed on 535 tumor samples to calculate the relative content of 22 types of TICs in each sample, the relative content of TICs in tumor samples and correlations analysis between TICs were shown in . The correlation analysis between the high and low expression groups of 2 key genes showed that 15 TICs were closely related to CCR2, including 10 positively related TIC subtypes and 5 negatively related TIC subtypes (). 13 types of TICs were correlated with PTPRC, including 7 positively related TIC subtypes and 6 negatively related TIC subtypes (). We subsequently divided LUAD patients into high and low groups based on median content of each TIC subtype and found that patients with a higher proportion of γδT cells and plasma cells had a better prognosis, while higher proportion of resting NK cells led to diametrically opposite results (). Taken together, these two genes were key factor affecting immune state of LUAD.

Figure 7

TICs analysis in tumor samples and correlation analysis of different TICs subtypes. (A) The relative content of 22 kinds of TICs in each sample based on CIBERSORT algorithm; (B) correlation analysis of 22 kinds of TICs. TIC, tumor-infiltrating immune cell.

Figure 8

Correlation analysis between the relative contents of TICs and CCR2 expression. (A) Violin plot showing the difference in the contents of 22 types of TICs in tumor tissues of high and low CCR2 expression groups; (B) scatter plot showing the correlation analysis between the proportions of 15 types of TICs and the expression of CCR2. TIC, tumor-infiltrating immune cell.

Figure 9

Correlation analysis between the relative contents of TICs and PTPRC expression. (A) Violin plot showing the difference in the contents of 22 types of TICs in tumor tissues of high and low PTPRC expression groups; (B) scatter plot showing the correlation analysis between the proportions of 13 types of TICs and the expression of PTPRC. TIC, tumor-infiltrating immune cell.

Figure 10

Three key immune cells affecting patients’ survival.

Discussion

Lung cancer is the most common cause of cancer-related deaths in the world’s population with a 5-year survival rate ranges from 4% to 17% depending on tumor stage and regional variation (26). According to statistics, there were 2.1 million lung cancer cases globally leading to 1.761 million deaths in 2018, ranking first in the incidence of malignant tumors (27). Although many advances have been made in understanding the mechanism of this disease, more than 1 million people die from lung cancer each year worldwide (28,29). Currently, the role of TME in tumorigenesis and cancer progression has become clearer (30), making it possible to predict treatment responsiveness and survival prognosis. Notably, the immune components in TME could mediate anti-tumor effect, and studies have indicated the correlation between tumor immune microenvironment and the efficacy of immunotherapy (31,32).

In our study, key genes, CCR2 and PTPRC, were obtained by analyzing the matrix and immune components of TME in LUAD. CCR2 was a chemokine receptor that regulated the immune response by inducing the recruitment of macrophages and monocytes to the site of inflammation (20). Intriguingly, numerous studies have elucidated that CCR2 was correlated with the occurrence and development of many diseases (20,33). High expression of CCR2 may lead to the poor prognosis of various cancers (21,34,35), and inhibition of CCR2 can enhance the therapeutic effect of PD-1 on tumor suppression (36). Paradoxically, the study by Zheng et al. proved that CCL2 was significantly down-regulated in SCLC and led to the proliferation of cancer cells (22). Our study indicated that the expression of CCR2 in LUAD was significantly lower compared with the corresponding normal tissue, and patients with high CCR2 expression had longer OS. Further analysis of the relationship between the expression of CCR2 and tumor stage shed bright light on the results that with the progress of T classification, the expression of CCR2 decreased. In a word, the above results demonstrated that CCR2 may play an anti-tumor effect in LUAD. In addition, the other key gene, PTPRC, encoded an evolutionarily highly conserved protein tyrosine phosphatase CD45, which was only expressed on nucleated cells of the hematopoietic system and was considered to be a regulator of B and T cell receptor signaling (24,37,38). In recent years, the role of CD45 in cancer and TME remained controversial. Chen et al. indicated that the accumulation of CD45+CD71+ erythroid cells in liver cancer may play an immunosuppressive effect in TME (39). Studies have confirmed that a higher proportion of CD45+ cells was closely related to the poor prognosis of NSCLC patients. However, in other types of tumors, high expression of CD45 may be related to the good prognosis of the tumor, indicating that CD45 may play a completely distinctive role in different tumors (24,25). Our study revealed that the expression of PTPRC in LUAD was significantly down-regulated compared with normal tissue, and patients with high PTPRC expression had better prognosis. Further analysis of the relationship between the expression of PTPRC and clinic characteristics indicated that with the progress of tumor stage, the expression of PTPRC decreased significantly.

We further calculated relative content of TICs of LUAD samples using CIBERSORT algorithm, the results found that patients with high proportion of γδT cells presented a better prognosis. γδT lymphocytes are a subset of T lymphocytes, which can directly inhibit tumor cells through cytotoxicity (40,41). It has been confirmed that γδT cells participated in the anti-tumor effect of lung cancer (42) and prostate cancer (43). Previous studies have confirmed that γδT cells have the property of suppressing and inhibiting a variety of tumor cell lines (44,45). However, some γδT lymphocyte subtypes were unexpected promoters for tumorigenesis and cancer development as their functions were affected by the immunosuppressive signal of TME (46). Our research found that the proportion of γδT cells was significantly correlated with the expression of CCR2 and PTPRC (P<0.05). This may suggest that CCR2 and PTPRC were the key factors that drive the anti-tumor effect of γδT cells in the TME. In summary, CCR2 and PTPRC may affect the function of immune cells, such as γδT cells, by participating in the regulation of the immune activity of TME and exert an impact on the prognosis of LUAD. CCR2 and PTPRC can be used as biomarkers to predict the immune response of TME and provide new therapeutic targets for LUAD. Notwithstanding, our study has some limitations. First, we still need to validate the findings of the study through in vivo and in vitro experiments. Secondly, we didn’t clarify the way in which key genes affect tumor immunity and exert anti-tumor effects.

Conclusions

Through bioinformatics to evaluate data in TCGA, our research showed that CCR2 and PTPRC may be potential prognostic markers in LUAD, which may affect the immune function of γδT cells and other immune cells by participating in the regulation of TME immune activity.

The article’s supplementary files as