Skip to main content

An immune scoring system predicts prognosis and immune characteristics in lung adenocarcinoma brain metastases by RNA sequencing

Abstract

Background

Previous studies have reported that the tumor immune microenvironment (TIME) was associated with the prognosis of lung cancer patients and the efficacy of immunotherapy. However, given the significant challenges in obtaining specimens of brain metastases (BrMs), few studies explored the correlation between the TIME and the prognosis in patients with BrMs from lung adenocarcinoma (LUAD).

Methods

Transcript profiling of archival formalin-fixed and paraffin-embedded specimens of BrMs from 70 LUAD patients with surgically resected BrMs was carried out using RNA sequencing. An immune scoring system, the green-yellow module score (GYMS), was developed to predict prognosis and immune characteristics in both BrMs and primary LUAD using Weighted Correlation Network analysis (WGCNA) and GSVA analysis. We comprehensively evaluated the immunological role of GYMS based on gene expression profile of LUAD BrMs by systematically correlating GYMS with immunological characteristics and immunotherapy responsiveness in the BrMs. Immunohistochemistry was applied for validation.

Results

We found that the high-GYMS group had better clinical prognosis and inflamed immune landscape including high infiltrations of various immune cells, increased immunomodulatory expression, and enriched immune-related pathways by using RNA-seq and immunohistochemical analysis. Low-GYMS group presented a lacked immune infiltration characteristic. Besides, the high-GYMS group had lower TIDE score and higher T-cell inflamed score than low-GYMS group. The GYMS has been validated in independent BrMs cohorts and primary NSCLC cohort treated with anti-PD-1/PD-L1, showing strong reproducibility and stability in both primary LUAD and BrMs. In addition, we construct a GYMS-related risk signature for patients with LUAD BrMs to predict prognosis.

Conclusions

We identified two immune-related subtypes which used to estimate prognosis and immune characteristics and developed a reliable GYMS-related risk signature in LUAD BrMs. These results will enhance the understanding of the immune microenvironment in LUAD BrMs and lay the theoretical foundation for the development of personalized therapies for LUAD patients with BrMs.

Introduction

Brain metastases (BrMs) are the most prevalent malignant brain tumors, even tenfold more common than primary malignant brain tumors [1, 2]. Among all types of cancer, lung adenocarcinoma (LUAD) is the most common cancer that metastasizes to the brain [3]. Unfortunately, patients with BrMs have significant mortality with poor median overall survival ranging from 7 to 13 months[4]. The conventional treatment pattern, including a combination of surgery, chemotherapy, and radiotherapy, has shown limited efficacy [5].

Recently, immune checkpoint inhibitors (ICIs), single-agent or combination strategies, is generally become the treatment option for patients with advanced NSCLC [6, 7]. Subgroup analyses in clinical trials and retrospective analyses provided evidence for the effective treatment of ICIs in NSCLC patients with BrMs, although there is currently a lack of large clinical trials specifically focusing on BrMs [8]. The tumor immune microenvironment is a key factor in immunotherapy. BrMs exhibit a distinct tumor immune microenvironment compared to primary lung cancer, as demonstrated by previous studies [9, 10]. Meanwhile, there was heterogeneity in the tumor immune microenvironment among different individuals. Therefore, it is crucial to characterize the specific immune microenvironment of each patient with BrMs to tailor individualized therapeutic approaches. New therapeutic strategies guided by biomarkers for the stratification of LUAD patients with BrMs are urgently needed to improve overall survival of BrMs patients.

Few prognostic and immune biomarkers in lung cancer with BrMs have been previously reported. It has been observed that BrMs patients with high expression of TIM-3 and LAG-3 in CD3+ T cell was also associated with longer survival (P < 0.05) [10]. Additionally, a tumor microenvironment classification, based on PD-L1 status and tumor-infiltrating lymphocytes (TILs), has been proposed [11]. A retrospective study showed that both PD‐L1 positivity on tumor cells and high infiltration of CD8 + T cell might have better post‐surgery outcomes in lung cancer patients after intracranial resection of BrMs [12]. However, there is still a lack of systematic studies to explore immune-related prognostic models in BrMs. Therefore, it is imperative to investigate the immune-related prognostic biomarkers in both BrMs and primary lung cancer. In this study, we aimed to develop a biomarker using RNA-seq analysis for 70 LUAD patients with BrMs to predict prognosis and characterize the immune landscape in both BrMs and primary LUAD.

Materials and methods

Figure 1 shows the workflow of this study.

Fig. 1
figure 1

The workflow of this study

Study cohort

Seventy Advanced LUAD patients with BrMs were diagnosed and enrolled at the Central South University Xiangya Hospital between 2017 and 2021 (Table S1). Resected BrMs samples and clinical data were collected. Overall survival (OS) was measured as the interval between the date of the surgical resection of brain metastases (BrM) and the date of death or the end of follow-up period. The study was approved by the ethics committee of Xiangya hospital, Central South University (No.202207391). We also obtained two independent NSCLC cohorts, including TCGA-LUAD and GSE135222. The Cancer Genome Atlas (TCGA) data of LUAD cohort: the RNA-seq data (COUNT) and clinical data (https://portal.gdc.cancer.gov/) was downloaded by utilizing "TCGAbiolinks" R package [13]. Then, the COUNT values were transformed into transcripts per kilobase million (TPM) values. An immunotherapy cohort of NSCLC from Samsung Medical Center (SMC cohort), GSE135222, with the gene expression profile and survival data were downloaded from the Gene Expression Omnibus (GEO). In the SMC cohort, DCB/responder were derived from patients who achieved partial response (PR) or stable disease (SD) that lasted more than 6 months. NDB/non-responder were derived from patients who had progressive disease (PD) or SD that lasted less than 6 months [14]. In addition, we obtained three independent BrMs cohorts for external validation of the immune scoring system. A NSCLC BrMs cohort (n = 43) with gene expression data from Sun Yat-Sen University Cancer Center was downloaded from OMIX (https://ngdc.cncb.ac.cn/omix/preview/Spe7IDiX, ID: OMIX575). The gene expression data and clinical data from the GSE50493 (Melanoma BrMs) and GSE173661(Breast cencer BrMs) also were downloaded from GEO.

Next-generation sequencing

BrMs tissues were collected and prepared into formalin-fixed paraffin-embedded (FFPE) samples. A total of 15 FFPE slides/sample were utilized, with a minimum tumor content of 20% confirmed by histological examination by pathologists. RNA was extracted from FFPE brain metastasis samples. All RNA with a percentage of RNA fragments > 200 nucleotides (DV200) ≤ 50% skipped fragmentation and proceeded to library preparation. Then, cDNA synthesis and NGS library preparation were performed using NEBNext® Ultra™ II Directional RNA Library Prep Kit (NEB#E7760L). The library was quantitated using Qubit 3.0 (life Invitrogen, USA), and quality was assessed with LabChip GX Touch (PerkinElmer, USA). After library quantification, the samples were subjected to high-throughput sequencing using DNBSEQ-T7. After removal of terminal adaptor sequences and low-quality data by using fastp (version: 0.19.5) [15] and removal of rRNA reads through aligning clean reads to the rRNA database (download from NCBI) by using bowtie2 (version:2.2.8), clean reads without known rRNA were aligned to the reference human genome (hg19) through STAR (version 020201) [16]. A series of quality control metrics were computed by using RNA-SeQC assessment [17]. A threshold of ≥ 80 million mapped reads and ≥ 10 million junction reads per sample was set. The procedure details can be found in Supplementary Method.

WGCNA construction and functional annotation

The gene expression data of our Xiangya BrMs cohort with 70 BrMs samples was used to construct a co-expression network by applying Weighted Gene Co-expression Network Analysis (WGCNA). Data processing of WGCNA based on the protocol of Horvath Lab UCLA (https://horvath.genetics.ucla.edu/html/CoexpressionNetwork/Rpackages/WGCNA/). 149 genes from the green-yellow module (Table S2) were chosen for gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses by using the “clusterprofiler” R package with a cutoff value of adjust P < 0.05. Then, the Cytoscape (version 3.9.0) [18] was used to map the network of the green-yellow module and screened and visualized the top 10 genes through the cytoHubba method [19].

Construction of immune scoring system and GSEA

To generate robust biomarkers to predict the prognosis and immune phenotype in LUAD BrMs, we constructed an immune scoring system by utilizing the GSVA method based on the expression of green-yellow module genes, which named the green-yellow module score (GYMS).

To explore the difference in biological process and pathway terms between the high-GYMS and low-GYMS groups, the GSEA was used via the R package "clusterProfiler". The gene sets of "c5.go.bp.v7.4.entrez.gmt", "c2.cp.kegg.v7.4.entrez.gmt", and "h.all.v7.4.1.entrez.gmt" were downloaded from the MSigDB (https://www.gsea-msigdb.org/gsea/msigdb) for performing the GSEA.

Immunological characteristics based on RNA-seq data in LUAD BrMs

The "ESTIMATE" R package was used to evaluate the general immune infiltration in BrMs samples, which is a tool that can predicts tumor purity by using gene signatures [20]. Then, to estimate the detailed immune cell infiltrations, we used six independent algorithms through the TIMER2.0 tool, including the Timer, CIBERSORT, MCP-COUNTER, QUANTISEQ, xCell, and ssGSEA [21,22,23]. CIBERSORT evaluate the proportions of 22 immune cells based on RNA-seq data by applying deconvolution algorithm [24]. MCP-COUNTER calculate 8 immune-cell lineage scores [25]. Meanwhile, the single sample gene set enrichment analysis (ssGSEA) algorithm was utilized to estimate the infiltration of 28 immune cells based on previously published gene sets from Charoentong et al. study [26]. Moreover, we collected gene sets of astrocyte, mature astrocyte, reactive astrocyte, and microglia signature from V. K. Han et al. study [27], and then using ssGSEA to estimate the infiltration abundance of four brain resident cells.

The prediction for therapeutic response

The Tumor Immune Dysfunction and Exclusion (TIDE) algorithm and subclass mapping (SubMap) method were utilized to predict clinical responses to ICIs. TIDE (http://tide.dfci.harvard.edu) is a computational method used to assess tumor immune escape potential and predict response to immunotherapy [28]. Responders were derived from patients who achieved complete responses, partial responses or stable disease control in response to anti-PD-1 therapy. Non-responders were derived from patients who had progressive disease [29]. The SubMap (https://www.genepattern.org/) was used to compare the similarity of the gene expression data of our group to previous melanoma cohort treated PD-1 inhibitor or CTLA-4 inhibitor [30]. In addition, a pan-cancer T cell-inflamed score was developed by Ayers et al., including 18 inflammatory genes, which could be an inflammatory biomarker and predict the clinical response of ICIs [31]. Here, T cell-inflamed score was calculated as described previously [32].

Gene expression profile data of human cancer cell lines were obtained from the Cancer Cell Line Encyclopedia (CCLE) project (https://depmap.org/portal/). Drug sensitivity data of human cancer cell lines were obtained from the Cancer Therapeutics Response Portal (CTRP) (https://portals.broadinstitute.org/ctrp.v2.1/) and  Profiling Relative Inhibition simultaneously in Mixtures (PRISM) Repurposing dataset (19Q4) (https://depmap.org/portal/prism/). To estimate the drugs response in LUAD BrMs, the R package "pRRophetic" was employed to calculate the area under the dose–response curve of each BrMs sample [33, 34]. The area under the curve (AUC) values were a measure indicator of drug sensitivity, and the lower AUC values mean the higher sensitivity to drugs in the two datasets. K-nearest neighbor (k-NN) imputation was used to impute the missing AUC values, as the previous study described [33]. Compounds with more than 20 percent of missing data were excluded before imputation.

Establishment of GYMS-related risk signature

The differentially expressed analysis was conducted through "DESeq2" R package. Based on the differentially expressed genes (DEGs) between high- and low-GYMS BrMs, Univariate Cox analysis was used to screen for prognostic genes in the Xiangya cohort. Then, prognostic genes were put into the least absolute shrinkage and selection operator (LASSO) Cox regression analysis utilizing the "glmnet" R package [35, 36]. Finally, a GYMS-related risk signature of LUAD BrMs was developed by weighting the Cox regression coefficients to estimate a risk score for each BrMs patient, classifying patients as low- and high-risk groups based on the median values of 70 BrMs samples. Kaplan–Meier analysis with the log-rank test was used to appraise the prognostic performance of the GYMS-related risk signature. To assess the sensitivity and specificity of the risk signature, the Time-dependent receiver operating characteristic (ROC) curve was displayed using the "timeROC" R package.

Immunohistochemistry

All tumor slides were stained with antibodies against CD8 (Cat. ab209775, 1:2000, Abcam), CD68(ZM-0060 Ready-to-Use, ZSGB-BIO), and CD163(ZA-0428 Ready-to-Use, ZSGB-BIO). The average from five independent areas containing the greatest abundance was calculated as the density of CD8+ T cells, CD68+ cells, and CD163+ cells under 200 × magnification (Table S3).

Statistical analysis

All statistical analyses were performed by using R software (Version 4.1.0). Non-parametric Mann–Whitney U test was used to compare gene expression level, ESTIMATE score, immune cells abundance, TIDE score, or sensitivity score. The correlation was analyzed by using the Spearman correlation test. The fisher exact test was used to compare the proportion of responders to immunotherapy. Survival analyses were conducted using the Kaplan–Meier method based on the log-rank test. The multiple hypothesis correction was performed for the comparison of multiple categories based on the false discovery rate (FDR) method between the two groups in Xiangya cohort (Table S4).

Results

Identification of the green-yellow module involved in CD8 + T cells and CD274 expression in LUAD BrMs via WGCNA

The CD8+ T cell level and PD-L1 expression are essential factors in determining the tumor microenvironment classification and the responsiveness of ICIs in NSCLC [37]. Thus, we performed weighted gene co-expression network analysis (WGCNA) by integrating the infiltration level of CD8+ T cells, cytotoxic lymphocytes, and CD274 expression. In our study, the top-5000 variation genes in RNA-seq data and power of β = 4 (scale free R2 = 0.89) as the soft threshold were selected to build the co-expression network. The analysis indicated that 13 co-expression modules were identified by applying the dynamic tree cut method (Fig. 2A, Figure S1A-F). Furthermore, we explored the relationships between the module gene and molecular traits, including CD8+ T cell infiltrations, cytotoxicity lymphocytes infiltrations, and CD274 expression. The heatmap of module trait correlations showed that the green-yellow module was highly positively correlated with CD8+ T cells and CD274 levels in LUAD BrMs, which meant that the green-yellow module (149 genes) may play an essential role in estimating the tumor immune microenvironment classification and the clinical response of ICIs (Fig. 2B). The scatter plots in Fig. 2C, D showed similar results. Thus, the green-yellow module was recognized as the key module.

Fig. 2
figure 2

Module-trait relationships and green-yellow-related enrichment analysis. A Cluster dendrogram of DEGs based on different metrics. Each color depicts a module that possesses weighted co-expressed genes. B Correlation of module eigengenes to CD8+ T cells and PD-L1 in LUAD BrMs. The values in the cells of each column are presented as p-value. Red represents a positive correlation, and blue represents a negative correlation. *: < 0.05, **: P < 0.01, ***: < 0.001, ****: < 0.0001. C, D Scatter plot show the correlation between Gene significance for CD8+ T cells or CD274 and module eigengenes in greenyellow module. E, F Gene Ontology and KEGG enrichment analysis for 149 genes in greenyellow module

Furthermore, the GO and KEGG analysis demonstrated that the genes of the green-yellow module were enriched in various immune processes and pathways. GO analysis indicated that lymphocyte activation, regulation of T cell activation, T cell activation and positive regulation of leukocyte activation were enriched in the green-yellow gene set (Fig. 2E). KEGG analysis showed that hematopoietic cell lineage, natural killer cell mediated cytotoxicity, cytokine-cytokine receptor interaction, T cell receptor signaling pathway were enriched in green-yellow module (Fig. 2F). These pathways had been demonstrated to be linked to immune process and immunotherapy [37, 38]. This shown the importance of the green-yellow module genes for the immune microenvironment and immunotherapy in LUAD BrMs, and it is required to further explore its characteristics and roles.

Construction and immune prediction of GYMS

Considering the individual differences in the immune microenvironment and fully exploiting the predictive immune value of green-yellow module genes, we constructed an immune scoring system (green-yellow module score, GYMS) by applying the GSVA algorithm based on 149 genes of green-yellow module in Xiangya BrMs cohort (n = 70). These patients were divided into low-GYMS (n = 35) and high-GYMS (n = 35) group based on the median value 0.0051 as a constant threshold across all cohort. To explore the performance of GYMS to predict the outcome in LUAD BrMs patients, we conducted a survival analysis for patients with survival information (n = 55). The survival analysis indicated that the high-GYMS group (n = 29) had a better OS benefit than low-GYMS group (n = 26) in our Xiangya BrMs cohort (P = 0.025, Fig. 3A). Moreover, the result of survival analysis had been verified in an immunotherapy cohort in lung cancer from Samsung Medical Center (high-GYMS vs low-GYMS; P = 0.06, Figure S2A) and the TCGA-LUAD cohort (high-GYMS vs low-GYMS; P = 0.031, Figure S3A). Furthermore, we investigated the correlation between the GYMS and clinical molecular information. Our analysis revealed no significant difference in the proportion of sex, age, smoking status, and molecular mutations between the high- and low-GYMS groups (Figure S4). Nevertheless, we observed a higher proportion of patients with KRAS mutations in the high-GYMS group (high-GYMS vs low-GYMS, 23.5% vs 6.1%), although statistical significance was not reached (Figure S4).

Fig. 3
figure 3

Immune characteristics of the green-yellow module score (GYMS). A Kaplan–Meier curve shows the effect of GYMS on OS with the constant threshold of 0.0051. B Heatmap shows the spearman’s correlation between GYMS and co-inhibitors. C Heatmap shows the expression in immunomodulators (MHC I, MHC II, co-stimulators, chemokines, and receptors) between low- and high-GYMS. D Correlation between GYMS and immune cells. E Correlation between GYMS and the infiltration levels of seven types of immune cells (CD8+ T cells, NK cells, CTL, macrophages, dendritic cells, MDSC, and Treg), which were estimated using six independent algorithms. F Heatmap shows the spearman's correlation between GYMS and ESTIMATE score. G Comparison of ESTIMATE score, immune score, and stromal score between low- and high-GYMS

Apart from the prognostic value, the GSEA was performed to explore the transcriptomic features produced by high-GYMS BrMs. By using GO, KEGG, and Hallmark gene sets as references, GSEA was performed between the high-GYMS and low-GYMS groups in patients from the Xiangya cohort. The level of GYMS in LUAD BrMs patients was positively correlated with the immune-related biological process and pathways, such as adaptive immune response, antigen processing and presentation, chemokine signaling pathway, natural killer cell mediated cytotoxicity, IFNG response, and inflammatory response (Figure S5; Table S5). Therefore, high-GYMS BrMs were correlated with high antigen presentation capacity and immune cell infiltrations.

Then, we evaluated the correlation between immunological characteristics and GYMS in LUAD BrMs. The GYMS positively correlated with the multiple immunomodulators, including MHC, Co-stimulators, and chemokines, which was consistent with the GSEA results (Fig. 3C). The high-GYMS group had higher expression of MHC I/II molecules than the low-GYMS group, showing a solid antigen presentation capacity (Figure S6D). A majority of chemokines and acceptor, including CXCL9, CXCL10 and CXCR3, were significantly upregulated in the high-GYMS group (P < 0.01; Figure S6A), promoting the recruitment of effector immune cells. Generally, the high expression of immune checkpoint molecules such as PD-L1/PD-1/CTLA4 was reported in inflamed TIME [39]. We found that the GYMS was positively correlated with many immune checkpoint molecules, including PD-L1, PD-1, CTLA-4, HAVCR2, LAG3, IDO1, PDCD1LG2, TIGIT, and BTLA (P < 0.01; Fig. 3B).

The immune score, stromal score and ESTIMATE score were also positively correlated with the GYMS (P < 0.001, Fig. 3F). High-GYMS had higher immune score, stromal score, and ESTIMATE score than the low-GYMS (P < 0.0001, Fig. 3G). Next, the correlation between GYMS and immune cells infiltration was evaluated using six different algorithms. We found that GYMS was positively correlated with overall level of immune cells (Fig. 3D, E). In addition to showing superior immune-activated cells (e.g. cytotoxic lymphocytes, CD8+ T cells, NK cells), high-GYMS BrMs also indicated a more increased abundance of immunosuppressive cells (e.g. M2-type TAMs, MDSC, Treg, and fibroblast). Meanwhile, the consistent results were observed in primary lung cancer in two external cohorts (the SMC immunotherapy cohort and the TCGA-LUAD cohort; Figure S2D-G, Figure S3B-E). Then, resident cells (astrocyte and microglia) in brain also were explored in this study. High-GYMS BrMs had higher infiltrations of the microglia and reactive astrocyte than low-GYMS BrMs (P < 0.01, Figure S7A,B). Therefore, these results demonstrated that high-GYMS BrMs were correlated with the inflamed tumor microenvironment that may well be susceptible to immunotherapy. However, it also indicates significant immune escape in high-GYMS BrMs.

The above results were validated in protein level. 66 FFPE samples from BrMs of LUAD was obtained to conduct immunohistochemical staining for immune cell markers (Table S3). High-GYMS BrMs had a higher infiltrating density of CD8+ T cells, CD68+ macrophages, and CD163+ M2 macrophages compared to low-GYMS BrMs (P < 0.05, Fig. 4B, D, F). Representative images are shown in Fig. 4A, C, E.

Fig. 4
figure 4

Immunohistochemical analysis shows the immune infiltration characteristics. Comparison of infiltrating density of CD8+ T cells (A, B), CD68+ macrophages (C, D), and CD163+ M2 macrophages (E, F) between high- and low-GYMS group. All tumor slides were stained with antibodies against CD8 (Cat. ab209775, 1:2000, Abcam), CD68(ZM-0060 Ready-to-Use, ZSGB-BIO), and CD163(ZA-0428 Ready-to-Use, ZSGB-BIO). Scale bars: 100 μm in 200x and 50 μm in 400x

To further explore new immune-related targets, the Cytohubba functions of the Cytoscape were used to discover ten most essential hub genes in the green-yellow module, including CD8A, LAG3, TMEM176B, TNFRSF13B, SH2D2A, KLHL6, PARP15, CXCL10, IGLL5, and SLAMF6 (Figure S1G). The CD8A is the marker gene of CD8+ T cells and the LAG3 is the important immune checkpoint molecule, which further demonstrates the accuracy of the above results. Given the clinical utility of 10 genes would be easier than 149 genes, utilizing the expression profiles of these 10 genes, we developed a 10-gene scoring system using the GSVA method. Nonetheless, our analysis revealed no significant difference in prognosis between the high and low score groups, although there was a trend toward longer median survival in the high score group (P > 0.05; Figure S1H).

The potential role of GYMS in ICIs responsiveness

To investigate whether GYMS can predict the clinical response to ICIs, we calculated the TIDE score in our Xiangya BrMs cohort. We found that the GYMS was positively correlated with CD274 (r = 0.34; P = 0.0042), IFNG (r = 0.71; P < 0.0001), CD8 (r = 0.64; P < 0.0001), and Dysfunction (r = 0.59; P < 0.0001) and negatively correlated with MDSC (r = − 0.32; P = 0.0061), TAM M2 (r = − 0.54; P < 0.0001), and Exclusion (r = − 0.37; P = 0.0018) (Fig. 5B). Besides, a strong negatively correlation (r = − 0.37; P = 0.0015) between GYMS and TIDE score had been found in the Xiangya BrMs cohort (Fig. 5B). The high-GYMS group had a borderline significant lower TIDE score than the low-GYMS group (P = 0.076, Fig. 5D). It has been reported that the low-TIDE score is a surrogate biomarker to predict a good response to cancer immunotherapy [28]. As expected, higher proportion of responder was found in the high-GYMS group compared to low-GYMS group (Fig. 5C). Meanwhile, we also used the SubMap algorithm to evaluate the clinical response to ICIs (PD-1 and CTLA-4 inhibitors) in high- and low-GYMS patients with LUAD BrMs. In this study, the high-GYMS BrMs displayed more promise in response to anti-PD-1 inhibitors (Bonferroni-corrected P = 0.008, Fig. 5E). In addition, we found that GYMS was positively correlated with T-cell inflamed score (r = 0.85; P < 0.0001) that can predict the clinical response of ICIs (Fig. 5A). Meanwhile, the expression of several therapeutic immune checkpoints, including PD-L1, CTLA-4, TIM-3, and LAG-3, was significantly higher in the high-GYMS group (P < 0.05, Figure S6B). Consistent with above results, in the anti-PD-1/PD-L1 NSCLC cohort, the high-GYMS group had a higher proportion of responders with durable clinical benefits (DCB) than the low-GYMS group (Figure S2B,C). Therefore, these results indicated that the GYMS is a potential biomarker with applications in predicting immunotherapy response in both LUAD patients with and without BrMs.

Fig. 5
figure 5

Immunotherapeutic responsiveness of the green-yellow module score (GYMS). A Scatter plot shows the correlation between GYMS and T-cell inflamed score. B Bar plot shows the correlation between GYMS and TIDE score. C Comparison of the proportion of the responder between low- and high-GYMS. D Comparison of the TIDE score between low- and high-GYMS. E Immunotherapeutic responses to anti-CTLA-4 and -PD-1 treatments in low- and high-GYMS patients

GYMS predicts immune phenotypes in external validation cohort of BrMs

In a NSCLC BrMs cohort (n = 43) from Sun Yat-Sen University Cancer Center, we observed that high-GYMS (n = 25) had higher expression of Co-stimulators and chemokine than low-GYMS (n = 18), such as CXCR3, CXCL9, CXCL10, CCL4, and CCL3 (Fig. 6C). The GYMS also was positively correlated with immune score (P < 0.001; Fig. 6A). Meanwhile, GYMS was uncovered to be positively correlated with T cells, CD8+ T cells, Cytotoxic lymphocytes, NK cells, Myeloid dendritic cells, monocytes, Treg, and MDSC (P < 0.05, Fig. 6D, E). GYMS was positively correlated with the T-cell inflamed score, showing high-GYMS BrMs had an inflamed tumor microenvironment (r = 0.58, P < 0.0001, Fig. 6F). Moreover, GYMS was negatively correlated with TIDE score; Reduced TIDE score was found in the high-GYMS group (P < 0.05, Fig. 6G). As expected, borderline significant higher proportion of responder was found in the high-GYMS group compared to low-GYMS group (P = 0.064; Fig. 6H). Then, the SubMap algorithm was used to evaluate the clinical response to ICIs, showing the high-GYMS BrMs might have better response to anti-PD-1 inhibitors (Bonferroni-corrected P = 0.024, Fig. 6I). Furthermore, GYMS was also positively correlated with expression of several immune checkpoints, including CD274, PDCD1, CTLA4, TIM3, LAG3, TIGIT (P < 0.05, Fig. 6B). Therefore, the GYMS can distinguish immune phenotypes and immunotherapy responsiveness in the external validation cohort of NSCLC BrMs.

Fig. 6
figure 6

Immune characteristics and immunotherapeutic responsiveness of GYMS in OMIX575 cohort (43 BrMs). A Heatmap shows the spearman’s correlation between GYMS and ESTIMATE score. B Heatmap shows the spearman’s correlation between GYMS and co-inhibitors. C Heatmap shows the expression in immunomodulators (MHC I, MHC II, co-stimulators, chemokines, and receptors) between low- and high-GYMS. D Correlation between GYMS and immune cells. E Correlation between GYMS and the infiltration levels of seven types of immune cells (CD8 + T cells, NK cells, CTL, macrophages, dendritic cells, MDSC, and Treg), which were estimated using six independent algorithms. F Scatter plot shows the correlation between GYMS and T-cell inflamed score. G Bar plot shows the correlation between GYMS and TIDE score. H Comparison of the proportion of the responder between low- and high-GYMS. I Immunotherapeutic responses to anti-CTLA-4 and -PD-1 treatments in low- and high-GYMS patients

To broaden the application of GYMS, we explored the performance in Pan-cancer BrMs cohorts. In melanoma BrMs and breast cancer BrMs, the T-cell inflamed score was highly positively correlated with the GYMS (P < 0.0001, Figure S8A,E). Besides, the higher CD274 expression and lower TIDE score was found in the high-GYMS group compared to low-GYMS group (Figure S8 B,C,F,G). As expected, the high-GYMS BrMs may benefit more from anti-PD1 therapy in Pan-cancer BrMs by using SubMap algorithm (P < 0.01, Figure S8D,H).

Identification of potential therapeutic drugs for low-GYMS group

To improve the poor prognosis of low-GYMS BrMs patients, we employed two methods to identify potential drugs based on the drug response data from the CTRP and PRISM database. On the one hand, differential drug response between low-GYMS (bottom quintile) and high-GYMS (top quintile) groups was performed to determine drugs of lower estimated AUC values in the low-GYMS group with log2FoldChange > 0.10 and Wilcox rank sum test P < 0.05 (Fig. 7A). On the other hand, correlation analysis between the GYMS and estimated AUC values was utilized to search drugs with positive correlation coefficient with Spearman’s correlation coefficient > 0.25 (Fig. 7A). Then, 2 CTRP-derived compounds (ABT-737 and navitoclax) and 1 PRISM-derived compound (romidepsin) had lower estimated AUC values in low-GYMS group and positively correlated with GYMS (Fig. 7B,C). Therefore, the three compounds were identified as potential candidate drugs for low-GYMS BrMs patients.

Fig. 7
figure 7

Identification of potential therapeutic drugs for low-GYMS BrMs patients. A The workflow for identifying drugs with higher drug sensitivity in low-GYMS BrMs patients. B Spearman’s correlation analysis and differential drug response analysis of 2 CTRP-derived compounds. C Spearman's correlation analysis and differential drug response analysis of 1 PRISM-derived compounds

Development and validation of a GYMS-related risk signature

To identify a biomarker to predict prognostic value, we constructed a GYMS-related risk signature in the Xiangya BrMs cohort. Differential expression analysis was used between the high-GYMS (n = 35) and low-GYMS (n = 35). A total of 354 differential expression genes (DEGs) with log2FC > 1 and adjust P < 0.05, including 283 up-regulated genes and 71 down-regulated genes (Fig. 8A; Table S6). Then, in 55 BrMs samples with survival information, the GYMS-related risk signature composed of 4 genes (MAB21L1, C14orf132, IL12RB1, and CCR2) was developed based on the 354 DEGs by applying Univariate Cox and LASSO Cox analysis (Fig. 8B, C; Table S7). Risk scores were calculated for each sample (risk score = 0.084* MAB21L1 + 0.029* C14orf132 + (− 0.005)* IL12RB1 + (− 0.077)*CCR2). BrMs patients of LUAD were divided into low-risk group (n = 28) and high-risk group (n = 27) using the median value as threshold (Fig. 8D). As shown in Fig. 8E, BrMs patients with high-risk had significantly shorter overall survival than those patients with low-risk (log-rank P < 0.0001), and consistent results were acquired in two external validation sets (the immunotherapy cohort and the TCGA-LUAD cohort; P < 0.05, Figure S9A, S10A). The time-dependent AUC of the risk signature was 0.66, 0.74, and 0.77 at 1-, 2-, and 3-year, respectively, which showed robust predictive value of the GYMS-related risk signature (Fig. 8F). Simultaneously, the AUC values at 1, 2, and 3 years for the GYMS-related risk signature outperformed those markers derived from immune score, CD8+ T cells, and CTL (Fig. 8F). The prognostic value of the risk signature was further validated in two external cohorts (Figure S9B,S10B). Then, the effect of the risk signature on the prognosis of BrMs patients was further validated. We found that the higher risk score correlated with poor prognosis in the univariate Cox regression analysis, which was statistical significance (P < 0.001; HR = 67.604, 95% CI 8.819–518.264, Fig. 8G). Next, three factors (Postoperative radiotherapy, Age, and GYMS-related risk signature) with statistical significance in univariate Cox analysis were imported into Multivariate Cox regression analysis (Fig. 8GH), and indicated that the GYMS-related risk signature was an independent prognostic biomarker in LUAD BrMs and primary LUAD by using the Xiangya cohort, the SMC immunotherapy cohort in NSCLC (Figure S9C,D), and the TCGA-LUAD cohort (Figure S10C,D) (P < 0.05). In addition to GYMS-related risk signature, Age and Postoperative radiotherapy are also independent prognostic factors in LUAD patients with BrMs (P < 0.05; Fig. 8H). Therefore, our findings confirmed that the GYMS-related risk signature could be used as a valuable prognostic indicator in LUAD patients with or without BrMs.

Fig. 8
figure 8

Developing a GYMS-related prognostic signature using LASSO Cox regression. A Volcano plot of differently expression genes between the high- and low-GYMS group. The x-axis range is from − 6 to 10. B, C LASSO Cox analysis identified genes most correlated to overall survival in the Xiangya cohort. D Visualization of the association of the risk scores with survival status and gene expression profiles in LUAD BrMs. E Kaplan–Meier survival analysis of low- and high-risk group in LUAD BrMs. F Comparison of the 1-,2-,3-year AUC value between the GYMS-related signature and other similar score, including Stromal score, Immune score, ESTIMATE score, CD8+ T cell infiltration, and CTL infiltration in LUAD BrMs. G, H Forest plot of the GYMS-related risk signature in univariate and multivariate cox analysis. CTL: Cytotoxic T Lymphocyte. KPS: Karnofsky Performance Status

Discussion

Overall, brain metastases (BrMs) exhibit a distinct tumor immune microenvironment [9]. Meanwhile, there exists heterogeneity in the tumor immune microenvironment among different individuals with BrMs. We developed an immune scoring system (GYMS) to characterize the immune landscape of each LUAD patients with BrMs and to predict prognosis. High-GYMS BrMs presented elevated infiltration of both innate and adaptive immune cells, indicative of an inflamed immune landscape, thus defining them as the high-immunity subtype. Low-GYMS BrMs show a scarcity of infiltrating immune cells, defined as the low-immunity subtype. The GYMS demonstrates its prediction power for immunotherapy responsiveness in both BrMs and primary LUAD. We further validated the reproducibility and stability of the GYMS in independent external cohorts. Moreover, based on the DEGs between the two GYMS subtypes, we established a GYMS-related risk signature that served as an independent prognosis factor for LUAD patients with BrMs. To our knowledge, this study represents the current largest cohort based on bulk RNA-seq data to introduce an immune scoring system for predicting prognosis and characterizing immune landscape in BrMs patients with LUAD.

The GYMS was developed and validated to assess the immune status and predict the potential responsiveness to immunotherapy in LUAD BrMs. We identified ten hub genes in the green-yellow module, which consists of CD8A, TMEM176B, TNFRSF13B, LAG3, SH2D2A, KLHL6, PARP15, CXCL10, IGLL5, and SLAMF6. LAG3 (lymphocyte-activation gene 3) is a known immune checkpoint gene, and its expression product serving as an inhibitory receptor on activated T cells [40]. Ligand binding to LAG3 expressed on T-cell surface leads to the inhibition of T-cell proliferation and loss of T-cell effector function. Thus, high expression of LAG3 ligands contributes to tumor immune escape [41]. Blockade of LAG3 enhances the function of CD8+ T cell, which becomes a target for immunotherapy. CD8A serves as a biomarker of intratumoral prevalence of CD8+ T cells and has been identified as a predictive marker for response to immunotherapy [42]. CXCL10 is a crucial molecule in mediating CD8+ T cells chemotaxis and activation of T cells by binding to the CXCR3 receptor [43]. A previous study suggested that promoting CXCL10 secretion by cancer cells and recruiting CD8+ T cells could sensitize the tumor to ICIs [44]. Accordingly, CD8+ T cells may play a crucial role in determining the prognosis of LUAD patients with BrMs, as well as the response to immunotherapy. As expected, our study revealed higher levels of CD8+ T cell infiltration in high-GYMS tumors.

The two GYMS subtypes exhibited distinct immune profiles. The high-GYMS subtype showed increased accumulation of adaptive and innate immune cells, such as CD8+ T cells, cytotoxic T lymphocytes (CTLs), and NK cells, as well as enhanced expression of co-stimulatory molecules. These findings can be corroborated through the application of various R algorithms in training and validation cohorts. Thus, the high-GYMS tumors depict an inflamed TIME. In addition, the overexpressed immune checkpoints (eg. PD-L1) is another essential marker of an inflamed TIME, which could be induced by TILs [39]. Avoiding excessive immune response is the primary role of these immune checkpoint molecules by inhibiting pre-existing cancer immunity [45]. However, we observed that high-GYMS BrMs also displayed immune-suppressive features. High infiltration of immunosuppressive cells, such as Tregs, MDSC, and M2-type TAMs, along with pro-tumor cytokines, were also detected in high-GYMS tumors, potentially facilitating tumor immune escape. Meanwhile, high-GYMS tumors exhibited increased infiltrations of reactive astrocytes compared to low-GYMS tumors. Previous studies have reported that tumor-associated reactive astrocytes promote the formation of a tumor immunosuppressive environment in CNS cancer and accelerate the metastatic potential of cancer cells to the brain [46, 47]. Moreover, the significantly increased expression of immune checkpoints in high-GYMS tumors revealed functionally exhausted T cells enriched, causing immune evasion of malignant tumors after activation. Then, we further explored the specific immune characteristics of the low-GYMS BrMs, which characterized by depleted cytotoxic immune cells, including CD8+ T cells and NK cells, and deficient phagocytosis resulting from decreased macrophages and antigen presentation capacity. Thus, the low-GYMS group belongs to the "cold" tumor-immune phenotypes.

Given the limited efficacy of current therapies for LUAD BrMs, it is crucial to identify novel therapeutic strategies and targets for improving treatment outcomes. Clinical trials have revealed that immunotherapy based on immune checkpoint blockade may offer potential benefit for LUAD patients with BrMs [6, 48]. Nevertheless, not all patients exhibit a positive response to immunotherapy, and there exists a discordant immune profile between primary tumors and BrMs [9]; thus, exploring novel and effective biomarkers is required to recognize LUAD patients with BrMs who may benefit from immunotherapy. In this study, the high-GYMS BrMs had higher expression of multiple immune checkpoints than the low-GYMS, potentially attributed to the increased infiltration of TILs. Besides, low TIDE score and high T-cell inflamed score were indicated in the high-GYMS BrMs. Therefore, high-GYMS BrMs may benefit from anti-PD-(L)1 therapy. Considering immunotherapy was required to play a role in both primary tumors and BrMs, the predicted value of the GYMS has also been verified in primary LUAD. As expected, a better immunotherapeutic response was observed in high-GYMS group of primary LUAD compared to low-GYMS group in an immunotherapy cohort. The pivotal role of TIME in immunotherapy underscores the potential to improve the efficacy of immunotherapy for LUAD BrMs through targeted interventions directed at specific cells and regulatory molecules within the TIME. Hence, reversing the immunosuppressive effects potentially mediated by M2-type macrophages and reactive astrocytes could serve as an effective strategy [27]. For example, targeting STAT3 can turn the pro-tumor M2 into a tumor suppressive M1 phenotype and inhibit astrocyte activation [47]. Collectively, for high-GYMS LUAD BrMs, co-inhibiting multiple immune checkpoints and combined regulator of M2-type TAM and astrocyte activation may yield a more potent therapeutic outcome compared to single-agent immunotherapy.

This study has several limitations. The study is retrospective, which restricts the application range of our results. Further clinical trials are needed to determine the predictive value for survival and the response to immunotherapy. Moreover, in the future, we will continue to follow up and reassess the survival outcomes of the Xiangya BrMs cohort given the majority of patients in the cohort don't reach the clinical endpoints.

Conclusions

In summary, we identified the most relevant gene modules for immunotherapy utilizing the WGCNA algorithm and subsequently established an immune scoring system for LUAD BrMs. The system could serve as a reliable biomarker for prognosis and predicting immune status in LUAD patients with BrMs.

Availability of data and materials

The raw sequence data reported in this paper have been deposited in the Genome Sequence Archive (Genomics, Proteomics & Bioinformatics 2021) [49] in National Genomics Data Center (Nucleic Acids Res 2022) [50], China National Center for Bioinformation / Beijing Institute of Genomics, Chinese Academy of Sciences (GSA-Human: HRA003286) that are publicly accessible at https://ngdc.cncb.ac.cn/gsa-human.

Abbreviations

AUC:

The area under the curve

BrMs:

Brain metastases

CCLE:

The cancer cell line encyclopedia

CTRP:

The cancer therapeutics response portal

CTL:

Cytotoxic T lymphocyte

DEGs:

Differentially expressed genes

FFPE:

Formalin-fixed paraffin-embedded

GYMS:

Green-yellow module score

GEO:

The gene expression omnibus

GO:

Gene ontology

ICIs:

Immune checkpoint inhibitors

KEGG:

Kyoto encyclopedia of genes and genomes

KPS:

Karnofsky performance status

LASSO:

The least absolute shrinkage and selection operator

LUAD:

Lung adenocarcinoma

ROC:

Receiver operating characteristic

RNA-seq:

RNA-sequencing

SubMap:

Subclass mapping

TCGA:

The cancer genome atlas

TIDE:

The tumor immune dysfunction and exclusion

TIME:

Tumor immune microenvironment

WGCNA:

Weighted correlation network analysis

References

  1. Lamba N, Wen PY, Aizer AA (2021) Epidemiology of brain metastases and leptomeningeal disease. Neuro Oncol 23(9):1447–1456

    Article  PubMed  PubMed Central  Google Scholar 

  2. Ostrom QT, Wright CH, Barnholtz-Sloan JS (2018) Brain metastases: epidemiology. Handb Clin Neurol 149:27–42

    Article  PubMed  Google Scholar 

  3. Pocha K et al (2020) Surfactant expression defines an inflamed subtype of lung adenocarcinoma brain metastases that correlates with prolonged survival. Clin Cancer Res 26(9):2231–2243

    Article  CAS  PubMed  Google Scholar 

  4. Kondziolka D et al (2005) Long-term survivors after gamma knife radiosurgery for brain metastases. Cancer 104(12):2784–2791

    Article  PubMed  Google Scholar 

  5. Achrol AS et al (2019) Brain metastases. Nat Rev Dis Primers 5(1):5

    Article  PubMed  Google Scholar 

  6. Xiao G et al (2021) Immune checkpoint inhibitors for brain metastases in non-small-cell lung cancer: from rationale to clinical application. Immunotherapy 13(12):1031–1051

    Article  CAS  PubMed  Google Scholar 

  7. Grant MJ, Herbst RS, Goldberg SB (2021) Selecting the optimal immunotherapy regimen in driver-negative metastatic NSCLC. Nat Rev Clin Oncol 18(10):625–644

    Article  PubMed  Google Scholar 

  8. Eguren-Santamaria I et al (2020) PD-1/PD-L1 blockers in NSCLC brain metastases: challenging paradigms and clinical practice. Clin Cancer Res 26(16):4186–4197

    Article  CAS  PubMed  Google Scholar 

  9. Kudo Y et al (2019) Suppressed immune microenvironment and repertoire in brain metastases from patients with resected non-small-cell lung cancer. Ann Oncol 30(9):1521–1530

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Lu BY et al (2021) Spatially resolved analysis of the T cell immune contexture in lung cancer-associated brain metastases. J Immunother Cancer 9(10):e002684

    Article  PubMed  PubMed Central  Google Scholar 

  11. Teng MW et al (2015) Classifying cancers based on T-cell infiltration and PD-L1. Cancer Res 75(11):2139–2145

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Li LL et al (2022) An integrated biomarker of PD-L1 expression and intraepithelial CD8(+) T cell infiltration was associated with the prognosis of lung cancer patients after intracranial resection of brain metastases. Thorac Cancer 13(13):1948–1960

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Colaprico A et al (2016) TCGAbiolinks: an R/Bioconductor package for integrative analysis of TCGA data. Nucleic Acids Res 44(8):e71

    Article  PubMed  Google Scholar 

  14. Jung H et al (2019) DNA methylation loss promotes immune evasion of tumours with high mutation and copy number load. Nat Commun 10(1):4278

    Article  PubMed  PubMed Central  Google Scholar 

  15. Chen S et al (2018) fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 34(17):i884–i890

    Article  PubMed  PubMed Central  Google Scholar 

  16. Dobin A et al (2013) STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29(1):15–21

    Article  CAS  PubMed  Google Scholar 

  17. DeLuca DS et al (2012) RNA-SeQC: RNA-seq metrics for quality control and process optimization. Bioinformatics 28(11):1530–1532

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Shannon P et al (2003) Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res 13(11):2498–2504

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Chin CH et al (2014) cytoHubba: identifying hub objects and sub-networks from complex interactome. BMC Syst Biol 8(Suppl 4):S11

    Article  PubMed  PubMed Central  Google Scholar 

  20. Yoshihara K et al (2013) Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun 4:2612

    Article  PubMed  Google Scholar 

  21. Li T et al (2020) TIMER2.0 for analysis of tumor-infiltrating immune cells. Nucleic Acids Res 48(1):509–514

    Article  Google Scholar 

  22. Finotello F et al (2019) Molecular and pharmacological modulators of the tumor immune contexture revealed by deconvolution of RNA-seq data. Genome Med 11(1):34

    Article  PubMed  PubMed Central  Google Scholar 

  23. Li B et al (2016) Comprehensive analyses of tumor immunity: implications for cancer immunotherapy. Genome Biol 17(1):174

    Article  PubMed  PubMed Central  Google Scholar 

  24. Chen B et al (2018) Profiling tumor infiltrating immune cells with CIBERSORT. Methods Mol Biol 1711:243–259

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Becht E et al (2016) Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol 17(1):218

    Article  PubMed  PubMed Central  Google Scholar 

  26. Charoentong P et al (2017) Pan-cancer immunogenomic analyses reveal genotype-immunophenotype relationships and predictors of response to checkpoint blockade. Cell Rep 18(1):248–262

    Article  CAS  PubMed  Google Scholar 

  27. Zhang Q et al (2022) The spatial transcriptomic landscape of non-small cell lung cancer brain metastasis. Nat Commun 13(1):5983

    Article  PubMed  PubMed Central  Google Scholar 

  28. Jiang P et al (2018) Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med 24(10):1550–1558

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Hugo W et al (2016) Genomic and transcriptomic features of response to Anti-PD-1 therapy in metastatic melanoma. Cell 165(1):35–44

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Hoshida Y et al (2007) Subclass mapping: identifying common subtypes in independent disease data sets. PLoS ONE 2(11):e1195

    Article  PubMed  PubMed Central  Google Scholar 

  31. Ayers M et al (2017) IFN-γ-related mRNA profile predicts clinical response to PD-1 blockade. J Clin Invest 127(8):2930–2940

    Article  PubMed  PubMed Central  Google Scholar 

  32. Cristescu R et al (2018) Pan-tumor genomic biomarkers for PD-1 checkpoint blockade-based immunotherapy. Science 362(6411):eaar3593

    Article  PubMed  PubMed Central  Google Scholar 

  33. Yang C et al (2021) Prognosis and personalized treatment prediction in TP53-mutant hepatocellular carcinoma: an in silico strategy towards precision oncology. Brief Bioinform 22(3):bbaa164

    Article  PubMed  Google Scholar 

  34. Geeleher P, Cox N, Huang RS (2014) pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS ONE 9(9):e107468

    Article  PubMed  PubMed Central  Google Scholar 

  35. Tibshirani R (1997) The lasso method for variable selection in the Cox model. Stat Med 16(4):385–395

    Article  CAS  PubMed  Google Scholar 

  36. Friedman J, Hastie T, Tibshirani R (2010) Regularization paths for generalized linear models via coordinate descent. J Stat Softw 33(1):1–22

    Article  PubMed  PubMed Central  Google Scholar 

  37. Chen X et al (2021) CD8(+) T effector and immune checkpoint signatures predict prognosis and responsiveness to immunotherapy in bladder cancer. Oncogene 40(43):6223–6234

    Article  CAS  PubMed  Google Scholar 

  38. Spangler JB et al (2015) Insights into cytokine-receptor interactions from cytokine engineering. Annu Rev Immunol 33:139–167

    Article  CAS  PubMed  Google Scholar 

  39. Spranger S et al (2013) Up-regulation of PD-L1, IDO, and T(regs) in the melanoma tumor microenvironment is driven by CD8(+) T cells. Sci Transl Med 5(200):200ra116

    Article  PubMed  PubMed Central  Google Scholar 

  40. Huard B et al (1994) Lymphocyte-activation gene 3/major histocompatibility complex class II interaction modulates the antigenic response of CD4+ T lymphocytes. Eur J Immunol 24(12):3216–3221

    Article  CAS  PubMed  Google Scholar 

  41. Lino AC et al (2018) LAG-3 inhibitory receptor expression identifies immunosuppressive natural regulatory plasma cells. Immunity 49(1):120-133.e9

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Tumeh PC et al (2014) PD-1 blockade induces responses by inhibiting adaptive immune resistance. Nature 515(7528):568–571

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Luo R et al (2019) Cisplatin facilitates radiation-induced abscopal effects in conjunction with PD-1 checkpoint blockade through CXCR3/CXCL10-mediated T-cell recruitment. Clin Cancer Res 25(23):7243–7255

    Article  CAS  PubMed  Google Scholar 

  44. Limagne E et al (2022) MEK inhibition overcomes chemoimmunotherapy resistance by inducing CXCL10 in cancer cells. Cancer Cell 40(2):136-152.e12

    Article  CAS  PubMed  Google Scholar 

  45. Hu J et al (2021) Siglec15 shapes a non-inflamed tumor microenvironment and predicts the molecular subtype in bladder cancer. Theranostics 11(7):3089–3108

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Henrik Heiland D et al (2019) Tumor-associated reactive astrocytes aid the evolution of immunosuppressive environment in glioblastoma. Nat Commun 10(1):2541

    Article  PubMed  PubMed Central  Google Scholar 

  47. Priego N et al (2018) STAT3 labels a subpopulation of reactive astrocytes required for brain metastasis. Nat Med 24(7):1024–1035

    Article  CAS  PubMed  Google Scholar 

  48. Goldberg SB et al (2020) Pembrolizumab for management of patients with NSCLC and brain metastases: long-term results and biomarker analysis from a non-randomised, open-label, phase 2 trial. Lancet Oncol 21(5):655–663

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Chen T et al (2021) The genome sequence archive family: toward explosive data growth and diverse data types. Genom Proteom Bioinform 19(4):578–583

    Article  Google Scholar 

  50. Database Resources of the National Genomics Data Center (2022) China National Center for Bioinformation in 2022. Nucleic Acids Res 50(D1):D27-d38

    Article  Google Scholar 

Download references

Acknowledgements

We thank the investigators and research staff involved in this study.

Funding

Research reported in this publication was supported by Advanced Lung Cancer Research Fund for targeted therapy in China (CTONG-YC20210303), Chen Xiao-Ping Foundation for the Development of Science and Technology of Hubei Province (CXPJJH121005-01), National Multidisciplinary Cooperative Diagnosis and Treatment Capacity (lung cancer z027002), Health Research Project of Hunan Provincial Health Commission (W20242005), and Beijing Xisike Clinical Oncology and Research Foundation (Y-HR2019-0185).

Author information

Authors and Affiliations

Authors

Contributions

Conception and design of the study: Rongrong Zhou, Gang Xiao; Data acquisition: Zhiyuan Liu; Data analysis and interpretation: Gang Xiao, Lifeng Li; IHC Staining: Guilong Tanzhu; Gang Xiao; Drafting of the manuscript: Gang Xiao; Critical comments and suggestions for revision of the manuscript: Lifeng Li, Xuan Gao, Rongrong Zhou, Xuefeng Xia. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Rongrong Zhou.

Ethics declarations

Ethics approval and consent to participate

All tissues were collected with the approval of the ethics committee of Xiangya Hospital, Central South University (No.202207391).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have 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

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, 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 you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. 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-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Xiao, G., Tanzhu, G., Gao, X. et al. An immune scoring system predicts prognosis and immune characteristics in lung adenocarcinoma brain metastases by RNA sequencing. acta neuropathol commun 12, 181 (2024). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40478-024-01895-9

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40478-024-01895-9

Keywords