Single-cell dissection reveals the role of aggrephagy patterns in tumor microenvironment components aiding predicting prognosis and immunotherapy on lung adenocarcinoma

Results

Heterogeneity of aggrephagy genes across cells

The workflow chart illustrates the general design of this study (Figure 1A). We found that the aggrephagy score was significantly increased in normal samples compared LUAD samples analyzed via the ssGSEA algorithm (Figure 1B, P < 0.001). We then examined the landscape of aggrephagy on each cell type by using the scRNA-seq data of LUAD. After using the Seurat pipeline, a total of 12,554 cells were divided into 18 clusters (Figure 1C) and were annotated with major cell types based on classical marker genes, including mast cells, macrophages, plasma cells, T cells, stromal cells, epithelial cells, fibroblasts, and smooth muscle cells (Figure 1D). Marker genes for each cell population showed a clear separation between each cell type (Figure 1E). To further investigate whether aggrephagy activity was dynamic among TME in LUAD scRNA-seq level, we quantified aggrephagy score using ssGSEA, AUCell, Ucell, addmodulescore, and singsore algorithms (Figure 1F). Interestingly, the outcomes demonstrated relatively heightened aggrephagy score in smooth muscle cells, fibroblasts, epithelial, stromal, and macrophage, while lower in T, plasma, and mast cells. Moreover, a comparison between average aggrephagy score in LUAD and normal tissues unveiled intriguing observations: macrophage, plasma, T, epithelial, and fibroblasts within tumors exhibited notably low average aggrephagy score, while stromal cells exhibited notably high average aggrephagy score (Figure 1J). As a result of comparing the expression patterns of aggrephagy genes across cell types in LUAD, we found that aggrephagy genes also showed high heterogeneity among cell types (Figure 1G). For instance, VIM exhibited higher expression levels in fibroblasts, and plasma cells, while showed lower expression levels in T cells and macrophages (Figure 1H). TUBA1C is highly expressed in macrophages and epithelial cells, but shared lower expression levels in other cells (Figure 1I), indicating the necessity for further investigation of aggrephagy genes by targeting specific cell types.

Figure 1. Overview of aggrephagy gene expression in the scRNA-seq for LUAD. (A) The overall design of this study. (B) The GSVA scores based on the aggrephagy gene set were compared between normal and LUAD samples. (C) t-SNE plot colored by 18 cluster of cells. (D) Cell type annotations clustered by Seurat t-SNE in the scRNA-seq data. (E) Dot plot showing representative marker genes for each cell type. (F) Bubble plot showing aggrephagy scores for each cell type. (G) Dot plot showing the expression level of aggrephagy genes in all cell types. (H, I) Expression of example genes VIM (H), and TUBA1C (I) in different cell types. (J) The difference in average aggrephagy score in each cell type in normal and LUAD samples. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Novel aggrephagy-mediated fibroblasts contribute to the TME of LUAD

Recent research has recognized CAFs as a key factor in TME, and they are emerging as a research hotspot [42]. Advancements in single-cell analysis techniques have led to a more comprehensive comprehension of the diversity and role of CAFs. For instance, according to one of the most well-established categorizations of CAFs, CAFs can be categorized as myCAF (myfibroblastic CAF), dCAF (development CAF), iCAF (inflammatory CAF), as well as pCAF (Pdpn CAF). Distinctive characteristics differentiate various subgroups of CAFs. Specifically, myCAF subset is situated in proximity to cancer cell nests and is characterized by elevated expression levels of both FAP+ and ꬰ-SMA. Conversely, the iCAF subset is localized far from cancerous cells and is marked by low expression levels of ꬰ-SMA but high expression levels of IL-6; the dCAF subset is distinguished by the activation of genes associated with diverse types of stem cells [43]. Similarly, to determine whether aggrephagy expression has an effect on CAFs, we performed a dimension reduction analysis. Pseudotime analysis reveals that aggrephagy genes are expressed at various developmental stages (Figure 2A). For instance, the early development stages of CAFs were characterized by the significant expression of TUBA4A, TUBA1C, and ARL13B, whereas TUBB4B, VIM, and TUBA1B were the feature of late development. Based on the NMF algorithm, CAFs were further divided into four clusters, and we identified that DYNC1I2+CAF−C1, DYNLL1+CAF−C2, PARK7+CAF−C3, as well as Non−Aggrephagy−CAF−C4 subtypes (Figure 2B, Supplementary Table 2). The developing status of NMF-based CAF clusters varied greatly as shown in pseudotime analysis (Figure 2C, 2D). Interestingly, using the Cellchat analysis, DYNC1I2+CAF−C1, DYNLL1+CAF−C2, and PARK7+CAF−C3 presented more and tighter connections with other cell types than Non-Aggre-CAF-C4 (Figure 2E–2G). GSVA was employed to derive scores for established classical CAF phenotype markers, with the aim of uncovering the possible phenotypes and functions of the aggrephagy CAF subtypes [44, 45] (Figure 2I). DYNC1I2+CAF−C1 showed the most prominent scores of pan-myCAF, pan-dCAF, and pan-pCAF, while PARK7+CAF−C3 exhibited the prominent scores of pan-iCAF. In contrast, Non-Aggrephagy-CAF-C4 had the lowest score among all classical CAF subtypes. In addition, we also investigated important CAF phenotype markers, like pro-inflammatory genes, neo-angiogenic genes, and MMPs. From the pathway heatmap (Figure 2H), aggrephagy-CAF subtypes had a significantly different expression of common pathway genes. PARK7+CAF−C3 had significantly higher expression of MMPs, ECM, Neo-Angio, and proinflammatory pathways genes, DYNC1I2+CAF−C1 had significantly higher expression of TGFb and RAS pathways gens, while Non-Aggre-CAF-C4 had lower expression of above pathways genes. Ultimately, we compared TF regulation features and subtype-specific TFs between aggrephagy-related CAFs as well as non-aggrephagy-related CAFs. We found that the average activities of TFs including FOS and ATF3 were exclusively higher in DYNC1I2+CAF−C1, while downregulated in DYNLL1+CAF−C2 except for FOXO3, NFIA, and CEBPB (Figure 2J). In aggregate, these discoveries imply that aggrephagy could exert a significant influence on CAF remodeling within the TME of LUAD.

Figure 2. Aggrephagy genes modified the features of CAF. (A) Pseudotime trajectory analysis of aggrephagy genes in CAFs. (B) NMF clustering and annotation in CAFs classified by aggrephagy gene expression features. (C, D) The developing status of NMF-based CAF clusters obtained in pseudotime analysis. (E) A heat map summarizing the outgoing (secreting) and incoming (target) signal pathways of NMF-based aggrephagy-related CAF subtypes and other cell types. (F, G) The number (F, number of interactions) and weight (G, interaction weights/strength) of cell-cell interactions between Agg-related CAF subtypes and other cell types. (H) Heatmap showing the different average expression of common signaling pathway genes in the aggrephagy-related CAF subtypes, including collagens, ECM, MMPs, TGFb, Neo-Angio, Contractile, RAS and Proinflammatory. (I) Different aggrephagy-related CAF subtypes were correlated with the previous signatures. (J) Heatmap showing the significantly different activities of TFs among aggrephagy-related CAF subtypes by comparing the average AUC using pySCENIC in Python.

Aggrephagy participate in TAMs metabolism and polarization

Subsequently, we investigated the potential impact of aggrephagy on the phenotypes and functions of TAMs. Throughout the pseudotime analysis, aggrephagy genes were observed to exhibit diverse expression patterns at different stages of macrophage development, implying their complex and dynamic roles in regulating the differentiation and function of macrophages (Figure 3A). NMF with aggrephagy genes separated macrophages into VIM+Mac−C1 and PARK7+Mac−C2 (Figure 3B, Supplementary Table 3). The developing status of NMF-based TAM clusters varied greatly as shown in pseudotime analysis (Figure 3C, 3D). Further, VIM+Mac−C1 and PARK7+Mac−C2 had comparatively stronger interactions with other cellular components, especially with epithelial cells (Figure 3E, 3F). We inferred the specific pathways of intercellular communication and found that mainly outgoing (secreting) of VIM+Mac−C1 and PARK7+Mac−C2 were EGF and VISFATIN, and incoming (target) signal pathways were PTN and MK (Figure 3G). A tumor microenvironment induces metabolic reprogramming of macrophages, leading to protumor macrophages having an overactive metabolism [37], we conducted scMetabolism algorithm to explore the metabolic heterogeneity in different aggrephay TAM subtypes (Figure 3H). Interestingly, the result revealed that the two TAM subtypes exhibit unique metabolic activation pathways, indicating that aggrephagy could be a critical modulator of metabolic regulation. For instance, PARK7+Mac−C2 was distinguished by the upregulation of multiple metabolic pathways such as glycolysis, TCA cycle, fatty acid degradation, as well as oxidative phosphorylation, while VIM+Mac−C1 was distinguished by the upregulation of fatty acid elongation and ether lipid metabolism. In addition, we computed the M1-like/M2-like polarization scores of different aggrephagy TAMs subtypes [46]. Notably, PARK7+Mac−C2 and VIM+Mac−C1 exhibited a preference for expressing genes that were upregulated in M1 macrophages, suggesting their enhanced M1-like anti-tumor functions (Figure 3I).

Figure 3. Tumor-associated macrophages (TAMs) differed in metabolism and polarization during aggrephagy. (A) Pseudotime trajectory analysis of aggrephagy genes in TAMs. (B) NMF clustering and annotation in TAMs classified by aggrephagy gene expression features. (C, D) The developing status of NMF-based TAM subtypes obtained in pseudotime analysis. (E, F) The number and weight of cell-cell interactions between aggrephagy-related TAM subtypes and other cell types. (G) A heat map summarizing the outgoing (secreting) and incoming (target) signal pathways of NMF-based aggrephagy-related TAM subtypes and other cell types. (H) Metabolic status of aggrephagy-related TAM subtypes analyzed by scMetabolism. (I) M1-like/M2-like phenotype scoring among different aggrephagy-related TAM subtypes.

Aggrephagy contributed to the subgroup formation in CD8+T cells

As CD8+T cells are known to significantly impact tumor metastasis and treatment outcomes, we conducted a comprehensive analysis of these cells [47]. Likewise, Monocle2 revealed that aggrephagy genes occupy different developmental stages (Figure 4A). TUBA1A, concordant with CAFs, was highly expressed during late development, while PARK7 and UBA52 were expressed during early development. Furthermore, CD8+T cells were categorized as four clusters via the NMF algorithm, and the developing status of NMF-based CD8+T cells varied greatly (Figure 4C, 4D, Supplementary Table 4), as well as named them as VIM+CD8+T cells C1, UBA52+CD8+T cells C2, TUBA4A+CD8+T cells C3, and TUBA1A+CD8+T cells C4, respectively (Figure 4B). Similarly, the aforementioned aggrephagy CD8+ T cells subtypes were found to have more extensive and robust interactions with other components of the TME, indicating their potential role in orchestrating the tumor immune microenvironment (Figure 4E). Furthermore, we inferred the specific pathways of intercellular communication and found that VIM+CD8+T cells C1 and TUBA1A+CD8+T cells C4 presented stronger activity in the CD40 pathways compared to UBA52+CD8+T cells C2, TUBA4A+CD8+T cells C3 (Figure 4F). We computed the T exhaustion and T cytotoxic scores for these aggrephagy related CD8+ T cell subtypes, utilizing a previously published panel for calculating scores that evaluate overall functions and phenotypes. These scores play a crucial role in determining the effects of ICB and its significance in LUAD progression [48]. Interestingly, we found that UBA52+CD8+T cells C2 showed higher T cytotoxic scores, TUBA4A+CD8+T cells C3 demonstrated higher T exhaustion scores but lower T cytotoxic scores, while VIM+CD8+T cells C1 and TUBA1A+CD8+T cells C4 presented lower both scores (Figure 4G), indicating TUBA4A likely to be involved in CD8+T cell exhaustion, result in tumor immune escape. In addition, we performed a comparison of the mean expression levels of genes related to T cell function and immune checkpoint inhibitors across these subtypes of CD8+ T cells. Consequently, VIM+CD8+T cells C1, UBA52+CD8+T cells C2, TUBA4A+CD8+T cells C3, and TUBA1A+CD8+T cells C4 exhibited different immune function-related terms (Figure 4H). Ultimately, we investigated the regulatory mechanisms through network regulatory analysis. TUBA1A+CD8+T cells C4 demonstrated apparent activation of TFs including FOS, JUNB, FOSB, BCLAF1, and BATF (Figure 4H). In conclusion, our investigation has illustrated the role of aggrephagy in the restructuring of CD8+ T cells within the TME of LUAD.

Figure 4. CD8+T cells differed in metabolism and polarization during aggrephagy. (A) Pseudotime trajectory analysis of aggrephagy genes in CD8+T cells. (B) NMF clustering and annotation in CD8+T cells classified by aggrephagy gene expression features. (C, D) The developing status of NMF-based CD8+T cells subtypes. (E) The number and weight of cell-cell interactions between aggrephagy-related CD8+T cells subtypes and other cell types. (F) A heat map summarizing the outgoing (secreting) and incoming (target) signal pathways of NMF-based aggrephagy-related CD8+T cells subtypes and other cell types. (G) Heatmap showing the comparison of CD8+T cell function signatures (exhaustion score and T cytotoxic score) between aggrephagy-related CD8+T cells subtypes. (H) Heatmap showing significantly different features among aggrephagy-related CD8+T cells subtypes, including Co-stimulations (left), Co-inhibitors (middle), and TFs (right).

Aggrephagy-mediated TME remodeling contributes to prognosis and immunotherapy response in LUAD

We have set out to investigate whether the newly defined aggrephagy subtypes have an impact on the survival rates of patients diagnosed with LUAD. Initially, we depicted the intercellular communication among all the aggrephagy-related subtypes, which elucidated all the possible connections between these subtypes (Figure 5A). Then we utilized ssGSEA algorithm to compute the enrichment score of each aggrephagy subtype based on the corresponding differentially expressed genes (Supplementary Table 5) and investigated their prognostic significance in LUAD patients across multiple cohorts. Through the implementation of univariate Cox regression analysis, we acquired the hazard ratio of each aggrephagy cell subtype in TCGA, GSE68465, GSE50081, GSE37745, GSE3141, GSE31210, Meta-GEO, and Meta cohorts, which showed that PARK7-Mac-C2, UBA52+CD8+T_cells−C2, PARK7+CAF-C3, and VIM+Mac-C1 were poor prognostic factors (Figure 5B, Supplementary Table 6). In addition, certain aggrephagy-related cell subtypes could distinguish the OS in TCGA, GSE68465, GSE50081, GSE37745, GSE3141, GSE31210, and Meta-GEO cohorts (Supplementary Figure 1A–1AJ). Remarkably, all cell subtypes related to aggrephagy were able to differentiate overall survival (OS) among meta-cohorts, which included 1512 patients with LUAD (Figure 5D–5L). Having established the potential role of aggrephagy in shaping the TME, we proceeded to investigate whether the aggrephagy-based patterns of TME had an impact on the response to ICB therapy. We employed the TIDE algorithm to predict the response of each patient to ICB therapy across multiple cohorts. Additionally, we observed that in the TCGA cohort, DYNC1I2+CAF−C1, DYNLL1+CAF−C2, and PARK7+CAF−C3, TUBA4A+CD8+T_cells−C3 were downregulated in responders, suggesting that these cells may be linked to ICB resistance (Figure 5M), which was verified in other cohorts (Supplementary Figure 2). Utilizing logistic regression, we obtained the odds ratio value for response prediction, unveiling the detrimental impact of these cells (Figure 5C, Supplementary Table 7). To verify the presence of the above cell types, we used IHC to observe cellular localization of DYNC1I2, TUBA1A, TUBA4A, UBA52 and VIM. We discovered that some of the fibroblasts in stromal did express DYNC1I2. (Figure 6A, first panel). Similarly, TUBA1A+ T cells were observed, showing significant positivity, from which immune cells were seen to exude (Figure 6A, second panel). TUBA4A was expressed on the cell nuclear of some T cells, with stromal cell nearby (Figure 6A, third panel). As a highly conserved nuclear and cytoplasmic protein, UBA52 was expressed in some T cells. As an immune cells, we also observed a high abundance of fibroblasts around UBA52+ T cells in pathological sections (Figure 6A, fourth panel). As we known, VIM was an epithelial-mesenchymal transition biomarker, highly expressed in stromal cell and tumor cell. In pathological sections, we found VIM was expressed in some macrophages, with tumor cell nearby. (Figure 6A, fifth panel).

Figure 5. Multiple aggrephagy cell subtypes influenced the prognosis and immunotherapy response of LUAD patients. (A) The landscape of all cell-cell communication within all identified aggrephagy subtypes. (B) A hazard ration of different aggrephagy subtypes in 8 LUAD cohorts. (C) Odd ratio produced by logistic regression of different subtypes in 8 cohorts. (D–L) All aggrephagy cell subtypes could distinguish the survival of patients in Meta cohort. (M) Various aggrephagy-related subtypes demonstrated significantly different infiltration in responders and non-responders of ICB in TCGA cohorts, as predicted by TIDE. Mann Whitney-Wilcoxon test was applied between responders and non-responders. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Figure 6. Identification of the hub predictive genes to construct ADPS. (A) IHC staining images of five crucial aggrephagy genes. Panel1: DYNC1I2+ fibroblasts were marked with a red arrow. Protein was mainly expressed on cell membrane (source HPA. Patient ID: 448, female, age 76 years). Panel2: TUBA1A+ T cells were marked with a red arrow. Protein was mainly expressed on cytoplasmic/membranous (source HPA. Patient ID: 448, female, age 76 years). Panel3: TUBA4A+ T cells were marked with a red arrow. Protein was mainly expressed on cytoplasmic/membranous (source HPA. Patient ID: 2403, female, age 65 years). Panel4: UBA52+ T cells were marked with red arrows. Protein was mainly expressed in the nucleus (source HPA. Patient ID: 2393, female, age 54 years). Panel5: VIM + Macrophages were marked with red arrows. Protein was mainly expressed in the cytoplasmic/membranous (source HPA. Patient ID: 1421, female, age 76 years). (B–I) Kaplan-Meier survival curves of the ADPS regarding OS in the TCGA (B), GSE3141 (C), GSE31210 (D), GSE37745 (E), GSE50081 (F), GSE68465 (G), MetaGEO (H), and Meta cohorts (I). (J–Q) Time-dependent ROC curves of the ADPS regarding 1-, 3-, and 5-year OS in the TCGA (J), GSE3141 (K), GSE31210 (L), GSE37745 (M), GSE50081 (N), GSE68465 (O), MetaGEO (P), and Meta cohorts (Q).

Aggrephagy-derived prognostic score (ADPS)

Subsequently, considering the concordant performances of aggrephagy subtypes in anticipating the overall survival outcomes as well as immunotherapy response in multiple LUAD cohorts, we retired marker genes and a total of 571 genes from aggrephagy-related subtypes were obtained (Supplementary Table 8). A total of 114 genes were found to have prognostic value through the assessment of each gene using univariate Cox regression analysis (Supplementary Figure 3A). To streamline the number of genes, we utilized LASSO Cox regression analysis, which resulted in 32 remaining genes at a lambda value of 0.0252 (Supplementary Figure 3B, 3C). The ADPS was ultimately generated by performing multivariate Cox regression analysis (Supplementary Figure 3D, Supplementary Table 9). We calculated the ADPS for each patient according to the expression as well as weighted regression coefficients of the ADPS-related genes, and then partitioned them to high ADPS and low ADPS groups after z-mean normalization. K-M survival analyses revealed that the mortality rate in the high ADPS group was significantly higher than the low ADPS group in the training cohort (TCGA-LUAD, n = 500, P < 0.005), and other seven validation cohorts GSE3141 (n = 111, P < 0.0001), GSE31210 (n = 226, P < 0.0001), GSE37745 (n = 106, P < 0.0001), GSE50081 (n = 126, P < 0.0001), GSE68465 (n = 442, P < 0.0001), MetaGEO (n = 1012, P < 0.0001), and Meta cohort (n = 1512, P < 0.0001) (Figure 6B–6I). Moreover, in the TCGA cohort, the ADPS exhibited excellent performance, with time-dependent AUCs of 0.77, 0.75, and 0.71 at 1, 3, and 5 years ((Figure 6J). Comparable results were also acquired across the validation cohorts GSE3141 (0.75/0.86/0.78), GSE31210 (0.94/0.90/0.94), GSE37745 (0.76/0.75/0.77), GSE50081 (0.84/0.91/0.84), GSE68465 (0.73/0.67/0.64), MetaGEO (0.75/0.75/0.70), as well as Meta cohort (0.73/0.75/0.7), respectively (Figure 6K–6Q).

Predictive value of ADPS for immunotherapy

In addition, we evaluated the prognostic significance of ADPS about immunotherapy, utilizing the real-world cohort including IMvigor210 and GSE78220 datasets. Within the IMvigor210 cohort, which consisted of 298 patients, responses to anti-PD-L1 receptor blockers ranged from the complete response (CR) and partial response (PR) to stable disease (SD) and progressive disease (PD). A significant clinical benefit was observed in the low ADPS group of the IMvigor210 cohort compared with the high ADPS group (Figure 7A). CR/PR patients presented lower ADPS than SD/PD patients (Figure 7B, P < 0.0001). The low ADPS group displayed a higher percentage of CR/PR compared to the high ADPS group (Figure 7C). Survival differences between different ADPS groups were significant specifically in patients with Stage I+II (Figure 7D, P < 0.001), and in Stage III+IV patients (Figure 7E, P < 0.001). Our finding also indicated that patients with a low ADPS had a significantly better overall survival outcome compared to those with a high ADPS in the GSE78220 cohort (Figure 7F, P < 0.001). Besides, CR/PR patients also presented lower ADPS than SD/PD patients (Figure 7G, P < 0.01). Additionally, a higher percentage of SD/PD was observed in the high ADPS group as compared to the low ADPS group (Figure 7H). In addition to the abundance of immune checkpoints, we also noted a significantly higher expression of these checkpoints in the low ADPS group as compared to the high ADPS group (Figure 7I). Utilizing the TIDE web tool, we discovered that the low ADPS group had notably lower TIDE scores and higher rates of immunotherapy response compared to the high ADPS group (Figure 7J). There was a significant positive correlation observed between ADPS and TIDE score (Figure 7K). Responders to ICIs were found to be more prevalent among patients belonging to the high ADPS group, as predicted by the TIDE algorithm (Figure 7L, P < 0.001). Parallel to these findings, the subclass mapping analysis (Submap, modules in GenePattern, https://cloud.genepattern.org) also revealed similarities in expression patterns between patients with low ADPS and those with melanoma who responded to immunotherapy [49, 50] (A melanoma dataset that responded to immunotherapy was selected as the reference, while default settings were applied) responding to ICB (Figure 7M).

Figure 7. Immunotherapy response prediction of the ADPS. (A) Kaplan-Meier survival curves of the ADPS regarding OS in the IMvigor210 cohort. (B) Differences in ADPS among immunotherapy responses in the Imvigor210 cohort. (C) Distribution of immunotherapy responses among ADPS groups in the Imvigor210 cohort. (D) Prognostic differences between ADPS groups in early stage (stage I+II) patients in the Imvigor210 cohort. (E) Prognostic differences between ADPS groups in advanced stage (stage III+IV) patients in the Imvigor210 cohort. (F) Kaplan-Meier survival curves of the ADPS regarding OS in the GSE78220 cohort. (G) Differences in ADPS among immunotherapy responses in the GSE78220 cohort. (H) Distribution of immunotherapy responses among ADPS groups in the GSE78220 cohort. (I) Boxplot of relative expression levels at immune checkpoints between the high and low ADPS patients. (J) Boxplot of TIDE score between the high and low ADPS patients. (K) The relationship between the TIDE score and ADPS. (L) Contingency table between immunotherapy responses and ADPS groups based on TIDE algorithm. (M) Submap analysis of the two groups and 47 pretreated patients with comprehensive immunotherapy annotations. In submap analysis, a smaller p-value indicated a greater similarity of paired expression profiles. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Immune landscape for ADPS

To explore the immune characteristics reflected by the ADPS, we investigated the association between the ADPS and immune cell type, immune scores, stromal scores, as well as critical immune checkpoint scores. Figure 8A, 8B demonstrated that the group with low ADPS exhibited elevated levels of immune infiltrating cells as well as immune modulators. This suggests the presence of an inflamed yet relatively immune-supportive microenvironment, which may be more conducive to benefiting from immunotherapy [51]. The immune score, stromal score, and ESTIMATE score had a negative correlation to ADPS, respectively (Figure 8C–8E). In parallel, patients in the low ADPS group had higher stromal, immune, and ESTIMATE scores compared to the high ADPS group (Figure 8F). Besides, we compared the status of the TRS score, Cytolytic activity, and Th1/IFN Score, which was more related to a more immunoreactive microenvironment between the two ADPS score groups [52, 53]. The results revealed that all of these indicators were significantly upregulated in the low ADPS group (Figure 8G–8I). According to 28 immune cells infiltration [54] assessed by ssGSEA (Figure 8J) and CIBERSORT algorithm (Supplementary Figure 4A), we further confirmed that the low ADPS group had markedly higher overall infiltration abundance than the high ADPS group. Thus, we defined the high ADPS group as “immune-cold” tumors and the low ADPS group as “immune-hot” tumors. Considering the upregulation of immune-related characteristics observed within the group of patients with a low ADPS, we sought to delve deeper into the biological mechanisms responsible for this phenomenon. KEGG-based GSVA analysis indicated the ADPS had a strong positive correlation with many pathways facilitating tumor growth, such as the P53 signaling pathway, cell cycle, glycolysis, and WNT signaling, while the immunological pathways, such as B/T cell receptor signaling pathways, cytokine receptor interaction pathway showed a weaker correlation with ADPS (Supplementary Figure 4B). Additionally, we further explored the immune characteristics from gene level, as displayed in Figure 9A, 9B, most ADPS-related genes were positively correlated with immune score as well as immune cell types.

Figure 8. Immune characteristics of the ADPS in the TCGA dataset. (A) The correlation between the ADPS and immune modulators. (B) Heatmap exhibiting the immune score, stromal score, critical checkpoints, and cell types calculated through CIBERSORT analysis of the high and low ADPS groups. (C–E) Correlations between ADPS and immune score (C), stromal score (D), and ESTIMATE score (E). (F) Boxplot of relative stromal, immune, and ESTIMATE score between high and low ADPS groups. (G) Boxplot of relative TRS score between high and low ADPS groups. (H) Boxplot of relative cytolytic activity between high and low ADPS groups. (I) Boxplot of relative Th1/IFN score between high and low ADPS groups. (J) Boxplot of relative infiltrate abundance of 28 immune cell types between high and low ADPS groups. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Figure 9. Identification of candidate agents with higher drug sensitivity in high ADPS group. (A) The correlation of gene expression and immune score. (B) The correlation of gene expression and immune cell types. (C) A Venn diagram for summarizing included compounds from CTRP and PRISM datasets. (D) Barplot of gemcitabine, paclitaxel, and BI-2536 CMap scores in patients with high ADPS. (E) Schematic outlining the strategy to identify agents with higher drug sensitivity in high ADPS patients. (F, G) The result of Spearman’s correlation analysis and differential drug response analysis of CTRP-derived compounds (F) and profiling relative inhibition simultaneously in mixtures (PRISM)-derived compounds (G). (H, I) The results of differential drug response analysis of CTRP -derived compounds (H) and PRISM -derived compounds (I), The lower the value of the y-axis, the greater the drug sensitivity. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Searching for potential therapeutic agents for the high ADPS patients

The CTRP and PRISM databases comprise gene expression as well as drug sensitivity profiles from hundreds of cancer cell lines (CCLs), providing an opportunity to construct a model for predicting drug response. After removing duplicates, a total of 1770 unique compounds were obtained in the two datasets, with 160 compounds being shared between them (Figure 9C). Figure 9E demonstrates our strategic approach to developing potential agents for high ADPS patients. Subsequently, we implemented this methodology to pinpoint promising agent candidates for the high ADPS patients and generated three CTRP-derived agents (BI-2536, paclitaxel, and SB-743921), as well as five PRISM-derived agents (cabazitaxel, epothilone−b, gemcitabine, ispinesib, and SNS−314). The AUC values of these agents, as estimated, showed a statistically significant negative correlation with ADPS scores and were significantly lower in the high ADPS group (Figure 9F–9I). Furthermore, we leveraged the differential expression profiles of LUAD patients as well as normal samples to identify potential candidate compounds using the Cmap tool. Specifically, we aimed to identify agents whose gene expression patterns were opposite to the LUAD-specific expression patterns, indicating potential efficacy in suppressing LUAD tumor growth (i.e., gene expression increased in tumor tissues but decreased by treatment of certain compounds). Following the cross-referencing of the results from CTRP and PRISM, we were left with a trio of prospective agents: BI-2536, a PLK inhibitor, paclitaxel, and gemcitabine (Figure 9D). Among them, BI-2536, with a CMap score of −87.58, exhibited high sensitivity in LUAD patients, implying its potential as a therapeutic agent for those with high ADPS.

Multi-omics alteration characteristics targeting ADPS

The differences in frequently altered chromosomes were detected in two ADPS groups (Figure 10A, 10B). We conducted an integrated analysis of mutations and copy number alterations (CNA, Figure 10C) to examine the genomic heterogeneity of the high and low ADPS groups. In the high ADPS group, we observed a higher frequency of mutations in classical tumor suppressor genes TP53 and CSMD3, as well as the oncogene KRAS, compared to the low ADPS group (Figure 10C). Additionally, compared to the low ADPS group, the high ADPS group showed significantly higher levels of amplification or deletion at the focal and chromosome arm levels, such as the amplification of 8q24.21, 12p12.1, as well as the deletion of 11p15.5, and 9p21.3 (Figure 10C). This finding was further supported at the gene level by the clear amplification of the oncogene MYC located at 8q24.21 and the distinct deletion of the tumor suppressor genes CDKN2A located at 9p21.3 (Figure 10C). Furthermore, we explored the mutation frequency of 10 major oncogenic pathways [55] between the high and low ADPS groups, and result found that both groups had detectable mutations in most of the oncogenic pathways, such as RTK-RAS, PI3K, TP53, NOTCH, and Hippo pathways (Figure 10D, 10E). We further explored the relationship between ADPS-related genes and LUAD by analyzing their correlations with different molecular signatures (Figure 10F). Our findings showed that FHL1 and ADGRD1 were strongly and positively correlated with molecular signatures associated with genomic instability in LUAD, including Aneuploidy Score, Homologous Recombination Defects, Fraction Altered, Number of Segments, as well as Nonsilent Mutation Rate. Additionally, the high ADPS group exhibited a higher TMB score in comparison to the low ADPS group (Figure 10G, P < 0.01). Only a small fraction of samples exhibited copy number variation (CNV) in the ADPS-related genes (Figure 10H). Additionally, we explored the single-nucleotide variant (SNV) mutations of the ADPS-related genes. It showed that COLA4A1, COL1A2, and LUZP2 had SNV mutations in more samples (Supplementary Figure 4C). Finally, we investigated the co-occurrence probability between ADPS-related genes and the top 10 genes with the highest mutation frequency. Supplementary Figure 4D reveals significant co-occurrence probabilities between CACYBP and ENU1, as well as between RELA and AP2M.

Figure 10. Multi-omics alteration characteristics of the ADPS in the TCGA dataset. (A, B) GISTIC 2.0-based chromosome amplifications and deletions in high (A) and low (B) ADPS groups. (C) Genomic alteration landscape according to ADPS. Tumor mutation burden (TMB), relative contribution of four mutational signatures, top 20 mutated genes and broad-level CNA (>20%). The proportion of the high and low ADPS groups in each alteration is presented in the right bar charts. (D, E) The mutation frequencies of nine common oncogenic pathways in the high ADPS (D) and low (E) groups. (F) Correlation heatmap of ADPS-related genes with Aneuploidy Score, Homologous Recombination Defects, Fraction Altered, Number of Segments, and Nonsilent Mutation Rate. (G) Difference of TMB score between high and low ADPS group. (H) CNV mutations (gain, loss, none) of ADPS-related genes. ^*p < 0.05; ^**p < 0.01; ^***p < 0.001. Abbreviation: ns: not significant.

Nomogram based on ADPS and clinical features

To better apply ADPS in a clinical setting and optimize its predictive performance, we performed univariate and multivariate Cox regression analysis to integrate clinicopathological characteristics as well as ADPS. ADPS was identified as the most significant independent prognostic factor of LUAD according to the results of univariate Cox analysis (hazard ratio (HR) = 2.718, 95% confidence interval (CI): 2.25–3.284, p < 0.001), as well as multivariate Cox regression analysis (hazard ration (HR) = 2.687, 95% confidence interval (CI): 2.127–3.394, p < 0.001), respectively (Figure 11A, 11B). The circos plot illustrated a significant correlation between ADPS and survival status, tumor stage, as well as TNM staging system in the TCGA cohort (Figure 11C). Moreover, the ADPS showed superior predictive accuracy, as indicated by the higher C-index values, compared to traditional clinical factors including age, gender, tumor stage, and TNM staging system, in various independent LUAD cohorts including TCGA, GSE31210, GSE37745, GSE50081, as well as GSE68465 (Figure 11D). Thus, a nomogram was created by integrating both the stage and ADPS, as illustrated in Figure 11G. The calibration plot illustrated the effective ability of the nomogram in predicting the actual survival outcomes (Figure 11E). Furthermore, the decision curve analysis (DCA) indicated that the nomogram had a superior ability to identify high ADPS patients compared to both the ADPS as well as stage alone, as illustrated in Figure 11F. TimeROC analysis indicated that the ADPS and nomogram had higher AUC values than other indicators (Supplementary Figure 4E). Our research further revealed that ADPS was significantly increased in the T stage (Supplementary Figure 4F), N stage (Supplementary Figure 4G), and tumor stage (Supplementary Figure 4I), without significant differences in M stage (Supplementary Figure 4H).

Figure 11. Prognostic value of the ADPS and nomogram construction. (A, B) Univariate Cox (A) and multivariate Cox (B) analysis of the ADPS and clinicopathological characteristics. (C) Circos plot of different clinical factors in two ADPS groups. (D) The C-index of the ADPS and various clinical factors in the TCGA, GSE31210, GSE37745, GSE50081, and GSE68465 datasets. (E) Calibration curves for 1, 3, and 5 years of nomogram. (F) Decision curve for nomogram. (G) Nomogram model integrating the ADPS and stage was constructed.

Author Contributions

XTS and FM designed the study. XTS, MYN, and YW analyzed the data, participated in data collection, and prepared the manuscript. HF, CLL, and PZ helped the analysis with constructive discussions. All authors critically revised the manuscript.

Conflicts of Interest

The authors declare no conflicts of interest related to this study.

Funding

This work was supported by the Beijing-Tianjin-Hebei Basic Research Cooperation Project (20JCZXJC00190) and Tianjin Key Medical Discipline (Specialty) Construction Project (TJYXZDXK-061B).

References

• 1. Gridelli C, Rossi A, Carbone DP, Guarize J, Karachaliou N, Mok T, Petrella F, Spaggiari L, Rosell R. Non-small-cell lung cancer. Nat Rev Dis Primers. 2015; 1:15009. https://doi.org/10.1038/nrdp.2015.9 [PubMed]
• 2. Siegel RL, Miller KD, Fuchs HE, Jemal A. Cancer Statistics, 2021. CA Cancer J Clin. 2021; 71:7–33. https://doi.org/10.3322/caac.21654 [PubMed]
• 3. Dagogo-Jack I, Shaw AT. Tumour heterogeneity and resistance to cancer therapies. Nat Rev Clin Oncol. 2018; 15:81–94. https://doi.org/10.1038/nrclinonc.2017.166 [PubMed]
• 4. Parzych KR, Klionsky DJ. An overview of autophagy: morphology, mechanism, and regulation. Antioxid Redox Signal. 2014; 20:460–73. https://doi.org/10.1089/ars.2013.5371 [PubMed]
• 5. Kocaturk NM, Akkoc Y, Kig C, Bayraktar O, Gozuacik D, Kutlu O. Autophagy as a molecular target for cancer treatment. Eur J Pharm Sci. 2019; 134:116–37. https://doi.org/10.1016/j.ejps.2019.04.011 [PubMed]
• 6. Miller DR, Thorburn A. Autophagy and organelle homeostasis in cancer. Dev Cell. 2021; 56:906–18. https://doi.org/10.1016/j.devcel.2021.02.010 [PubMed]
• 7. Malampati S, Song JX, Chun-Kit Tong B, Nalluri A, Yang CB, Wang Z, Gopalkrishnashetty Sreenivasmurthy S, Zhu Z, Liu J, Su C, Krishnamoorthi S, Iyaswamy A, Cheung KH, et al. Targeting Aggrephagy for the Treatment of Alzheimer's Disease. Cells. 2020; 9:311. https://doi.org/10.3390/cells9020311 [PubMed]
• 8. Ma X, Lu C, Chen Y, Li S, Ma N, Tao X, Li Y, Wang J, Zhou M, Yan YB, Li P, Heydari K, Deng H, et al. CCT2 is an aggrephagy receptor for clearance of solid protein aggregates. Cell. 2022; 185:1325–45.e22. https://doi.org/10.1016/j.cell.2022.03.005 [PubMed]
• 9. Gong H, Yang X, Zhao Y, Petersen RB, Liu X, Liu Y, Huang K. Amyloidogenicity of p53: a hidden link between protein misfolding and cancer. Curr Protein Pept Sci. 2015; 16:135–46. [PubMed]
• 10. Ishimaru D, Andrade LR, Teixeira LS, Quesado PA, Maiolino LM, Lopez PM, Cordeiro Y, Costa LT, Heckl WM, Weissmüller G, Foguel D, Silva JL. Fibrillar aggregates of the tumor suppressor p53 core domain. Biochemistry. 2003; 42:9022–7. https://doi.org/10.1021/bi034218k [PubMed]
• 11. Kanapathipillai M. Treating p53 Mutant Aggregation-Associated Cancer. Cancers (Basel). 2018; 10:154. https://doi.org/10.3390/cancers10060154 [PubMed]
• 12. Obradovic A, Chowdhury N, Haake SM, Ager C, Wang V, Vlahos L, Guo XV, Aggen DH, Rathmell WK, Jonasch E, Johnson JE, Roth M, Beckermann KE, et al. Single-cell protein activity analysis identifies recurrence-associated renal tumor macrophages. Cell. 2021; 184:2988–3005.e16. https://doi.org/10.1016/j.cell.2021.04.038 [PubMed]
• 13. Devalaraja S, To TKJ, Folkert IW, Natesan R, Alam MZ, Li M, Tada Y, Budagyan K, Dang MT, Zhai L, Lobel GP, Ciotti GE, Eisinger-Mathason TSK, et al. Tumor-Derived Retinoic Acid Regulates Intratumoral Monocyte Differentiation to Promote Immune Suppression. Cell. 2020; 180:1098–114.e16. https://doi.org/10.1016/j.cell.2020.02.042 [PubMed]
• 14. Young TM, Reyes C, Pasnikowski E, Castanaro C, Wong C, Decker CE, Chiu J, Song H, Wei Y, Bai Y, Zambrowicz B, Thurston G, Daly C. Autophagy protects tumors from T cell-mediated cytotoxicity via inhibition of TNFα-induced apoptosis. Sci Immunol. 2020; 5:eabb9561. https://doi.org/10.1126/sciimmunol.abb9561 [PubMed]
• 15. Bai J, Liu T, Tu B, Yuan M, Shu Z, Fan M, Huo S, Guo Y, Wang L, Wang H, Zhao Y. Autophagy loss impedes cancer-associated fibroblast activation via downregulating proline biosynthesis. Autophagy. 2023; 19:632–43. https://doi.org/10.1080/15548627.2022.2093026 [PubMed]
• 16. Wang D, Jiao Z, Ji Y, Zhang S. Elevated TUBA1A Might Indicate the Clinical Outcomes of Patients with Gastric Cancer, Being Associated with the Infiltration of Macrophages in the Tumor Immune Microenvironment. J Gastrointestin Liver Dis. 2020; 29:509–22. https://doi.org/10.15403/jgld-2834 [PubMed]
• 17. Jin W. Novel Insights into PARK7 (DJ-1), a Potential Anti-Cancer Therapeutic Target, and Implications for Cancer Progression. J Clin Med. 2020; 9:1256. https://doi.org/10.3390/jcm9051256 [PubMed]
• 18. Chang YH, Lee SH, Chang HC, Tseng YL, Lai WW, Liao CC, Tsay YG, Liao PC. Comparative secretome analyses using a hollow fiber culture system with label-free quantitative proteomics indicates the influence of PARK7 on cell proliferation and migration/invasion in lung adenocarcinoma. J Proteome Res. 2012; 11:5167–85. https://doi.org/10.1021/pr300362g [PubMed]
• 19. Garcia-Mata R, Gao YS, Sztul E. Hassles with taking out the garbage: aggravating aggresomes. Traffic. 2002; 3:388–96. https://doi.org/10.1034/j.1600-0854.2002.30602.x [PubMed]
• 20. Dost AFM, Moye AL, Vedaie M, Tran LM, Fung E, Heinze D, Villacorta-Martin C, Huang J, Hekman R, Kwan JH, Blum BC, Louie SM, Rowbotham SP, et al. Organoids Model Transcriptional Hallmarks of Oncogenic KRAS Activation in Lung Epithelial Progenitor Cells. Cell Stem Cell. 2020; 27:663–78.e8. https://doi.org/10.1016/j.stem.2020.07.022 [PubMed]
• 21. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012; 28:882–3. https://doi.org/10.1093/bioinformatics/bts034 [PubMed]
• 22. Mariathasan S, Turley SJ, Nickles D, Castiglioni A, Yuen K, Wang Y, Kadel EE III, Koeppen H, Astarita JL, Cubas R, Jhunjhunwala S, Banchereau R, Yang Y, et al. TGFβ attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature. 2018; 554:544–8. https://doi.org/10.1038/nature25501 [PubMed]
• 23. Hugo W, Zaretsky JM, Sun L, Song C, Moreno BH, Hu-Lieskovan S, Berent-Maoz B, Pang J, Chmielowski B, Cherry G, Seja E, Lomeli S, Kong X, et al. Genomic and Transcriptomic Features of Response to Anti-PD-1 Therapy in Metastatic Melanoma. Cell. 2016; 165:35–44. https://doi.org/10.1016/j.cell.2016.02.065 [PubMed]
• 24. Satija R, Farrell JA, Gennert D, Schier AF, Regev A. Spatial reconstruction of single-cell gene expression data. Nat Biotechnol. 2015; 33:495–502. https://doi.org/10.1038/nbt.3192 [PubMed]
• 25. Mo S, Shen X, Wang Y, Liu Y, Sugasawa T, Yang Z, Gu W, Nakajima T. Systematic single-cell dissecting reveals heterogeneous oncofetal reprogramming in the tumor microenvironment of gastric cancer. Hum Cell. 2023; 36:689–701. https://doi.org/10.1007/s13577-023-00856-z [PubMed]
• 26. Gao Y, Wang H, Chen S, An R, Chu Y, Li G, Wang Y, Xie X, Zhang J. Single-cell N⁶-methyladenosine regulator patterns guide intercellular communication of tumor microenvironment that contribute to colorectal cancer progression and immunotherapy. J Transl Med. 2022; 20:197. https://doi.org/10.1186/s12967-022-03395-7 [PubMed]
• 27. Butler A, Hoffman P, Smibert P, Papalexi E, Satija R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nat Biotechnol. 2018; 36:411–20. https://doi.org/10.1038/nbt.4096 [PubMed]
• 28. Bray NL, Pimentel H, Melsted P, Pachter L. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016; 34:525–7. https://doi.org/10.1038/nbt.3519 [PubMed]
• 29. Qiu X, Mao Q, Tang Y, Wang L, Chawla R, Pliner HA, Trapnell C. Reversed graph embedding resolves complex single-cell trajectories. Nat Methods. 2017; 14:979–82. https://doi.org/10.1038/nmeth.4402 [PubMed]
• 30. Vento-Tormo R, Efremova M, Botting RA, Turco MY, Vento-Tormo M, Meyer KB, Park JE, Stephenson E, Polański K, Goncalves A, Gardner L, Holmqvist S, Henriksson J, et al. Single-cell reconstruction of the early maternal-fetal interface in humans. Nature. 2018; 563:347–53. https://doi.org/10.1038/s41586-018-0698-6 [PubMed]
• 31. Guan X, Polesso F, Wang C, Sehrawat A, Hawkins RM, Murray SE, Thomas GV, Caruso B, Thompson RF, Wood MA, Hipfinger C, Hammond SA, Graff JN, et al. Androgen receptor activity in T cells limits checkpoint blockade efficacy. Nature. 2022; 606:791–6. https://doi.org/10.1038/s41586-022-04522-6 [PubMed]
• 32. Van de Sande B, Flerin C, Davie K, De Waegeneer M, Hulselmans G, Aibar S, Seurinck R, Saelens W, Cannoodt R, Rouchon Q, Verbeiren T, De Maeyer D, Reumers J, et al. A scalable SCENIC workflow for single-cell gene regulatory network analysis. Nat Protoc. 2020; 15:2247–76. https://doi.org/10.1038/s41596-020-0336-2 [PubMed]
• 33. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics. 2013; 14:7. https://doi.org/10.1186/1471-2105-14-7 [PubMed]
• 34. Sun X, Li Z, Meng F, Huang X, Wang J, Song J, Sun L, Zhang P. Cuproptosis associated genes affect prognosis and tumor microenvironment infiltration characterization in lung adenocarcinoma. Am J Cancer Res. 2022; 12:4545–65. [PubMed]
• 35. Galbo PM Jr, Zang X, Zheng D. Molecular Features of Cancer-associated Fibroblast Subtypes and their Implication on Cancer Pathogenesis, Prognosis, and Immunotherapy Resistance. Clin Cancer Res. 2021; 27:2636–47. https://doi.org/10.1158/1078-0432.CCR-20-4226 [PubMed]
• 36. Bischoff P, Trinks A, Obermayer B, Pett JP, Wiederspahn J, Uhlitz F, Liang X, Lehmann A, Jurmeister P, Elsner A, Dziodzio T, Rückert JC, Neudecker J, et al. Single-cell RNA sequencing reveals distinct tumor microenvironmental patterns in lung adenocarcinoma. Oncogene. 2021; 40:6748–58. https://doi.org/10.1038/s41388-021-02054-3 [PubMed]
• 37. Wu Y, Yang S, Ma J, Chen Z, Song G, Rao D, Cheng Y, Huang S, Liu Y, Jiang S, Liu J, Huang X, Wang X, et al. Spatiotemporal Immune Landscape of Colorectal Cancer Liver Metastasis at Single-Cell Level. Cancer Discov. 2022; 12:134–53. https://doi.org/10.1158/2159-8290.CD-21-0316 [PubMed]
• 38. Yu L, Shen N, Shi Y, Shi X, Fu X, Li S, Zhu B, Yu W, Zhang Y. Characterization of cancer-related fibroblasts (CAF) in hepatocellular carcinoma and construction of CAF-based risk signature based on single-cell RNA-seq and bulk RNA-seq data. Front Immunol. 2022; 13:1009789. https://doi.org/10.3389/fimmu.2022.1009789 [PubMed]
• 39. Yang C, Huang X, Li Y, Chen J, Lv Y, Dai S. Prognosis and personalized treatment prediction in TP53-mutant hepatocellular carcinoma: an in silico strategy towards precision oncology. Brief Bioinform. 2021; 22:bbaa164. https://doi.org/10.1093/bib/bbaa164 [PubMed]
• 40. Subramanian A, Narayan R, Corsello SM, Peck DD, Natoli TE, Lu X, Gould J, Davis JF, Tubelli AA, Asiedu JK, Lahr DL, Hirschman JE, Liu Z, et al. A Next Generation Connectivity Map: L1000 Platform and the First 1,000,000 Profiles. Cell. 2017; 171:1437–52.e17. https://doi.org/10.1016/j.cell.2017.10.049 [PubMed]
• 41. Mayakonda A, Lin DC, Assenov Y, Plass C, Koeffler HP. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 2018; 28:1747–56. https://doi.org/10.1101/gr.239244.118 [PubMed]
• 42. Mao X, Xu J, Wang W, Liang C, Hua J, Liu J, Zhang B, Meng Q, Yu X, Shi S. Crosstalk between cancer-associated fibroblasts and immune cells in the tumor microenvironment: new findings and future perspectives. Mol Cancer. 2021; 20:131. https://doi.org/10.1186/s12943-021-01428-1 [PubMed]
• 43. Kanzaki R, Pietras K. Heterogeneity of cancer-associated fibroblasts: Opportunities for precision medicine. Cancer Sci. 2020; 111:2708–17. https://doi.org/10.1111/cas.14537 [PubMed]
• 44. Kasashima H, Duran A, Martinez-Ordoñez A, Nakanishi Y, Kinoshita H, Linares JF, Reina-Campos M, Kudo Y, L'Hermitte A, Yashiro M, Ohira M, Bao F, Tauriello DVF, et al. Stromal SOX2 Upregulation Promotes Tumorigenesis through the Generation of a SFRP1/2-Expressing Cancer-Associated Fibroblast Population. Dev Cell. 2021; 56:95–110.e10. https://doi.org/10.1016/j.devcel.2020.10.014 [PubMed]
• 45. Affo S, Nair A, Brundu F, Ravichandra A, Bhattacharjee S, Matsuda M, Chin L, Filliol A, Wen W, Song X, Decker A, Worley J, Caviglia JM, et al. Promotion of cholangiocarcinoma growth by diverse cancer-associated fibroblast subpopulations. Cancer Cell. 2021; 39:883. https://doi.org/10.1016/j.ccell.2021.05.010 [PubMed]
• 46. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst. 2015; 1:417–25. https://doi.org/10.1016/j.cels.2015.12.004 [PubMed]
• 47. Cui C, Wang J, Fagerberg E, Chen PM, Connolly KA, Damo M, Cheung JF, Mao T, Askari AS, Chen S, Fitzgerald B, Foster GG, Eisenbarth SC, et al. Neoantigen-driven B cell and CD4 T follicular helper cell collaboration promotes anti-tumor CD8 T cell responses. Cell. 2021; 184:6101–18.e13. https://doi.org/10.1016/j.cell.2021.11.007 [PubMed]
• 48. Guo X, Zhang Y, Zheng L, Zheng C, Song J, Zhang Q, Kang B, Liu Z, Jin L, Xing R, Gao R, Zhang L, Dong M, et al. Global characterization of T cells in non-small-cell lung cancer by single-cell sequencing. Nat Med. 2018; 24:978–85. https://doi.org/10.1038/s41591-018-0045-3 [PubMed]
• 49. Hubble J, Demeter J, Jin H, Mao M, Nitzberg M, Reddy TB, Wymore F, Zachariah ZK, Sherlock G, Ball CA. Implementation of GenePattern within the Stanford Microarray Database. Nucleic Acids Res. 2009; 37:D898–901. https://doi.org/10.1093/nar/gkn786 [PubMed]
• 50. Lu X, Jiang L, Zhang L, Zhu Y, Hu W, Wang J, Ruan X, Xu Z, Meng X, Gao J, Su X, Yan F. Immune Signature-Based Subtypes of Cervical Squamous Cell Carcinoma Tightly Associated with Human Papillomavirus Type 16 Expression, Molecular Features, and Clinical Outcome. Neoplasia. 2019; 21:591–601. https://doi.org/10.1016/j.neo.2019.04.003 [PubMed]
• 51. Wu W, Liu Y, Zeng S, Han Y, Shen H. Intratumor heterogeneity: the hidden barrier to immunotherapy against MSI tumors from the perspective of IFN-γ signaling and tumor-infiltrating lymphocytes. J Hematol Oncol. 2021; 14:160. https://doi.org/10.1186/s13045-021-01166-3 [PubMed]
• 52. Rooney MS, Shukla SA, Wu CJ, Getz G, Hacohen N. Molecular and genetic properties of tumors associated with local immune cytolytic activity. Cell. 2015; 160:48–61. https://doi.org/10.1016/j.cell.2014.12.033 [PubMed]
• 53. Danilova L, Ho WJ, Zhu Q, Vithayathil T, De Jesus-Acosta A, Azad NS, Laheru DA, Fertig EJ, Anders R, Jaffee EM, Yarchoan M. Programmed Cell Death Ligand-1 (PD-L1) and CD8 Expression Profiling Identify an Immunologic Subtype of Pancreatic Ductal Adenocarcinomas with Favorable Survival. Cancer Immunol Res. 2019; 7:886–95. https://doi.org/10.1158/2326-6066.CIR-18-0822 [PubMed]
• 54. Charoentong P, Finotello F, Angelova M, Mayer C, Efremova M, Rieder D, Hackl H, Trajanoski Z. Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep. 2017; 18:248–62. https://doi.org/10.1016/j.celrep.2016.12.019 [PubMed]
• 55. Sanchez-Vega F, Mina M, Armenia J, Chatila WK, Luna A, La KC, Dimitriadoy S, Liu DL, Kantheti HS, Saghafinia S, Chakravarty D, Daian F, Gao Q, et al, and Cancer Genome Atlas Research Network. Oncogenic Signaling Pathways in The Cancer Genome Atlas. Cell. 2018; 173:321–37.e10. https://doi.org/10.1016/j.cell.2018.03.035 [PubMed]
• 56. Marsh T, Debnath J. Autophagy suppresses breast cancer metastasis by degrading NBR1. Autophagy. 2020; 16:1164–5. https://doi.org/10.1080/15548627.2020.1753001 [PubMed]
• 57. Poillet-Perez L, White E. Role of tumor and host autophagy in cancer metabolism. Genes Dev. 2019; 33:610–9. https://doi.org/10.1101/gad.325514.119 [PubMed]
• 58. Janssens R, Struyf S, Proost P. The unique structural and functional features of CXCL12. Cell Mol Immunol. 2018; 15:299–311. https://doi.org/10.1038/cmi.2017.107 [PubMed]
• 59. Marković I, Savvides SN. Modulation of Signaling Mediated by TSLP and IL-7 in Inflammation, Autoimmune Diseases, and Cancer. Front Immunol. 2020; 11:1557. https://doi.org/10.3389/fimmu.2020.01557 [PubMed]
• 60. Mulder K, Patel AA, Kong WT, Piot C, Halitzki E, Dunsmore G, Khalilnezhad S, Irac SE, Dubuisson A, Chevrier M, Zhang XM, Tam JKC, Lim TKH, et al. Cross-tissue single-cell landscape of human monocytes and macrophages in health and disease. Immunity. 2021; 54:1883–900.e5. https://doi.org/10.1016/j.immuni.2021.07.007 [PubMed]
• 61. Vitale I, Manic G, Coussens LM, Kroemer G, Galluzzi L. Macrophages and Metabolism in the Tumor Microenvironment. Cell Metab. 2019; 30:36–50. https://doi.org/10.1016/j.cmet.2019.06.001 [PubMed]
• 62. Barrett CS, Millena AC, Khan SA. TGF-β Effects on Prostate Cancer Cell Migration and Invasion Require FosB. Prostate. 2017; 77:72–81. https://doi.org/10.1002/pros.23250 [PubMed]
• 63. Zhang X, Zhang C, Qiao M, Cheng C, Tang N, Lu S, Sun W, Xu B, Cao Y, Wei X, Wang Y, Han W, Wang H. Depletion of BATF in CAR-T cells enhances antitumor activity by inducing resistance against exhaustion and formation of central memory cells. Cancer Cell. 2022; 40:1407–22.e7. https://doi.org/10.1016/j.ccell.2022.09.013 [PubMed]
• 64. Borst J, Ahrends T, Bąbała N, Melief CJM, Kastenmüller W. CD4⁺ T cell help in cancer immunology and immunotherapy. Nat Rev Immunol. 2018; 18:635–47. https://doi.org/10.1038/s41577-018-0044-0 [PubMed]
• 65. Wu F, Zhang W, Shao H, Bo H, Shen H, Li J, Liu Y, Wang T, Ma W, Huang S. Human effector T cells derived from central memory cells rather than CD8(+)T cells modified by tumor-specific TCR gene transfer possess superior traits for adoptive immunotherapy. Cancer Lett. 2013; 339:195–207. https://doi.org/10.1016/j.canlet.2013.06.009 [PubMed]
• 66. Malta TM, Sokolov A, Gentles AJ, Burzykowski T, Poisson L, Weinstein JN, Kamińska B, Huelsken J, Omberg L, Gevaert O, Colaprico A, Czerwińska P, Mazurek S, et al, and Cancer Genome Atlas Research Network. Machine Learning Identifies Stemness Features Associated with Oncogenic Dedifferentiation. Cell. 2018; 173:338–54.e15. https://doi.org/10.1016/j.cell.2018.03.034 [PubMed]
• 67. Zhou J, Yang Q, Lu L, Tuo Z, Shou Z, Cheng J. PLK1 Inhibition Induces Immunogenic Cell Death and Enhances Immunity against NSCLC. Int J Med Sci. 2021; 18:3516–25. https://doi.org/10.7150/ijms.60135 [PubMed]
• 68. Uras IZ, Moll HP, Casanova E. Targeting KRAS Mutant Non-Small-Cell Lung Cancer: Past, Present and Future. Int J Mol Sci. 2020; 21:4325. https://doi.org/10.3390/ijms21124325 [PubMed]
• 69. Mogi A, Kuwano H. TP53 mutations in nonsmall cell lung cancer. J Biomed Biotechnol. 2011; 2011:583929. https://doi.org/10.1155/2011/583929 [PubMed]
• 70. Li B, Xiang Z, Xiong F, Yan B, Huang Q. C-terminal binding protein-2 is a prognostic marker for lung adenocarcinomas. Medicine (Baltimore). 2020; 99:e20076. https://doi.org/10.1097/MD.0000000000020076 [PubMed]
• 71. Wang H, Wang X, Xu L, Lin Y, Zhang J. CCT6A and CHCHD2 Are Coamplified with EGFR and Associated with the Unfavorable Clinical Outcomes of Lung Adenocarcinoma. Dis Markers. 2022; 2022:1560199. https://doi.org/10.1155/2022/1560199 [PubMed]
• 72. He T, Shen H, Wang S, Wang Y, He Z, Zhu L, Du X, Wang D, Li J, Zhong S, Huang W, Yang H. MicroRNA-3613-5p Promotes Lung Adenocarcinoma Cell Proliferation through a RELA and AKT/MAPK Positive Feedback Loop. Mol Ther Nucleic Acids. 2020; 22:572–83. https://doi.org/10.1016/j.omtn.2020.09.024 [PubMed]
• 73. Li H, Liu J, Shen S, Dai D, Cheng S, Dong X, Sun L, Guo X. Pan-cancer analysis of alternative splicing regulator heterogeneous nuclear ribonucleoproteins (hnRNPs) family and their prognostic potential. J Cell Mol Med. 2020; 24:11111–9. https://doi.org/10.1111/jcmm.15558 [PubMed]
• 74. Jin R, Wang X, Zang R, Liu C, Zheng S, Li H, Sun N, He J. Desmoglein-2 modulates tumor progression and osimertinib drug resistance through the EGFR/Src/PAK1 pathway in lung adenocarcinoma. Cancer Lett. 2020; 483:46–58. https://doi.org/10.1016/j.canlet.2020.04.001 [PubMed]
• 75. Zhang P, Pei S, Gong Z, Feng Y, Zhang X, Yang F, Wang W. By integrating single-cell RNA-seq and bulk RNA-seq in sphingolipid metabolism, CACYBP was identified as a potential therapeutic target in lung adenocarcinoma. Front Immunol. 2023; 14:1115272. https://doi.org/10.3389/fimmu.2023.1115272 [PubMed]