Characterization of cellular subtype in COPD and LCCOPD
3 freshly resected lung specimens from non-small cell lung cancer patients within COPD were collected, and adjacent COPD tissue from a distal region within the same lobe was served as controls (Fig.1A). Following resection, tissues were rapidly digested to a single-cell suspension and unbiasedly subjected to scRNA-seq. These 6 samples (3 COPD tissues and 3 COPD-based lung cancer) were undertaken single-cell RNA sequencing and divided into COPD group and LCCOPD group. After quality filtering, a total of 27246 cells were detected and 61.8% of cells originated from malignant tissue. (Supplementary information, Tab.S1, Tab.S2). Principle component analysis (PCA) using variably expressed genes was used to generate t-SNE at different resolutions (Fig.1B, C). Marker genes were used to identify cell types (Supplementary information, Fig. S1). At low resolution, cells clustered based on different sources whereas at the high resolution they clustered based on patient identity and draw a heatmap showing the top enrichment gene of all identified cell groups (Fig.1F).
Gene-expression profiles of 27246 cells were retained after quality control and data filtering. The resulting cell clusters were annotated by established marker genes based on cell types., We classified our clusters with 12 broad cell types and roughly divided all cell subpopulations into structural cells and immune cell populations, depending on the role of cell subtypes in our body29. The structural cell population includes endothelial cells, epithelial cells, ciliated cells, and fibroblasts. The immune cell populations can be subdivided into natural immune cell populations, such as NK cells, monocytes, and DC cells, and acquired immune cell populations, such as macrophages, B cells, T cells, and plasma cells.
Overall, the main cell populations maintain a high degree of overlap in COPD and LCCOPD. Due to data discrepancies, we prioritize the normalization of data and visualized cellular composition. After normalization, monocyte, fibroblast cells, endothelial cells, T cells, and plasma cells account for a similar proportion in the cellular microenvironment under different disease conditions (COPD or LCCOPD). However, there is also some difference in cell subgroups between COPD and LCCOPD groups. In the LCCOPD groups, dendritic cells, epithelial cells, and B cells account for a larger proportion than COPD groups, while NK cells, ciliated cells, and macrophages account for more in COPD groups. Overall, in the LCCOPD environment, mainly structural cells and primary immune populations are significantly enriched, while in the COPD environment, there was a significant difference in the enrichment of powerful immune cell populations. The LCCOPD microenvironment showed a certain degree of tumor-killing ability deficiency (Fig.1D, E).
Tissue-specific T-reg cells showed significant changes in metabolic profiles and different FOXP3 expression
7939 transcriptomes (4505 for COPD, 3434 for LCCOPD) of T cells were identified and categorized as 10 clusters, each with corresponding marker genes in both COPD groups and LCCOPD groups (Fig.2A, Supplementary information, Fig.S2E). We draw the distribution graphs according to the type of organization (Fig.2B) and more than 45000 genes were detected and under-analyzed. T cells were usually divided into CD8+ T cell, CD4+ T cell, and T-regs, so we categorized all 10 clusters as specific T cell subgroups (Fig.2C), and gave each cell a risk score to estimate the contribution (Supplementary information, Fig. S2A, B). From the perspective of composition (Fig.2D, Supplementary information, Fig. S2C, D), we identified each cluster to be a tumor-relative or COPD-relative group based on the proportion of tissue origin.
We conducted a heatmap of top differential enrichment genes for respective biological functions. (Fig.2E). Compared to COPD tissues, tumor-related clusters expressed a large number of mitochondrial and ribosomal genes, such as MT-CO2, MT-CO1, MT-ND4, RPS27, etc. Highly expressed mitochondrial genes and ribosomal genes are often related to protein synthesis activities30. COPD-related cluster8 was marked by highly expressed HLA molecules, such as HLA-A, HLA-B, and HLA-C, which showed a high degree of tissue identification. Other clusters were not specifically enriched in COPD, nor enriched in cancer tissues. This result matched the conclusions of previous studies that CD4+ T cell has both anti-tumor effects and tumor-promoting effects31.
Compared to COPD tissues, T-reg cells of LCCOPD are mostly concentrated in cancer foci. Given the key role of T-reg cells in cancer, we re-clustered the T-reg cells and subdivided them into 10 clusters, namely c0-c9. (Fig.2F, Supplementary information Fig. S2G).c9, c8, and c0 were identified as the tumor-related T-regs, and c5 was identified as the COPD-related T-regs, because of the difference in composition. By visualizing the most differentially expressed non-ribosomal or mitochondrial genes in each cluster (Supplementary information, Fig. S2I), tumor-related clusters are marked by a large number of signal transduction factors, such as IL6ST, RGCC, IKZF3, etc., which are involved in the regulation of immune metabolism and cell-cycle process, while COPD-related clusters were mainly marked by genes such as CCL5, AHNAK, MALAT1, etc., which are involved in the process of immune regulation.
Comparing metabolic programs in different environments, we noticed that several metabolic processes were differentially expressed in COPD and LCCOPD groups (Fig.2G). In the LCCOPD groups, pathway enrichment analysis for T-reg cell function showed that metabolic pathways such as oxidative phosphorylation, glycosylation process, and glycosylation, inhibitory receptors such as CTLA4, PD-1, functional pathways such as IL2-stat, TGF-β, epithelial-mesenchymal transition, and cell cycle signaling were significantly up-regulated. As for the COPD groups, T-regs showed a highly up-regulated lipid metabolic process and were more active in the response to IL10, IFN-α, and VEGF.
After visualizing different metabolic and functional pathways for each T-reg clusters
(Fig.2H). It was found that the metabolic characteristics were not completely consistent among different clusters. Overall, Glucose-related metabolic pathways are significantly correlated with the acquisition of T-regs immunosuppressive ability and proliferation activity, while lipid metabolism is related to T-regs hypoxia sensitivity and pro-angiogenic activity.
To figure out the mechanism of metabolic differences, we constructed an upstream transcription factor network by analyzing targeted key functional genes. Through the Chea3 database, we identified a total of 1600 transcription factors that may be involved in functional regulation and mapped the top 10 functional transcription factors (Fig. 2I, 2G). It can be found that FOXP3 is the central and important transcription factor among them. We also mapped the expression levels of all 10 transcription factors and found significant differences in the expression of each transcription factor between clusters (Fig. 2K, Supplementary information, Fig. S2F, H). Under these two conditions, the IRF transcription factor family exhibits differential enhancement patterns. COPD-related T-regs showed a high enrichment on IRF8. Highly IRF8 expression in T-regs controlled the type Ⅰ response and IRF8 is an identity-keeper for suppressive Th1-like T-regs32. And LCCOPD-elated T-regs showed a high enrichment on IRF4 and IRF5. Highly IRF4 expression in T-regs controlled the type Ⅱ response, and IRF4 and FOXP3 synergistically participate in the differentiation of effector T-regs33. IRF4 is necessary for the inhibition of Th2-driven autoimmunity by T-regs34.
By comparing T-regs associated with COPD tissues and LCCOPD tissues, we found that LCCOPD-associated T-regs expressed stronger proliferation activities, such as up-regulated oxidative phosphorylation pathways and higher expression of FOXP3. At the same time, in the highly proliferative T-reg infiltration environment, the expression of inhibitory ability such as CTLA-4 and PD-1 was also up-regulated. It implied that FOXP3 in T-regs was progressively activated in the tumor environment and linked to dramatic changes in several metabolic processes and cytokine pathways.
Tissue-specific macrophages exhibit distinct polarization
A total of 6529 transcriptomes (5062 for COPD, 1467 for LCCOPD) in macrophages were identified and divided into 11clusters, each with corresponding marker genes in both COPD groups and LC groups (Fig.3A, Supplementary information, Fig. S3B). We draw distribution graphs according to the type of organization and the risk radio (Fig.3B, Supplementary information, Fig. S3A). Among all 6529 transcriptomes, more than 70% of them were from COPD tissue, significantly higher than that in cancer tissue. We deduced that this is related to the anti-tumor effect of macrophages. We visualized the number distribution for 11 clusters to find representative subgroups that are abundant in COPD or cancer tissues (Fig. 3C, D). For all 11 clusters, only c5 and c4 had a higher composition in LCCOPD groups, and identified as the LCCOPD-related macrophage. Other clusters, such as c0, c3, c9, c2, c1, and c8 almost comes from the COPD tissue environment, and had a specific relationships of COPD tissue. The differentially enriched gene analysis for each cluster revealed the HLA molecular symbol was significantly different. Among COPD-related clusters, macrophage HLA molecules were mainly HLA-DOB1, HLA-DRB1, HLA-B, HLA-C, etc. For LCCOPD-related macrophages, HLA-DOA2 and HLA-DRB6 were the main high-expressed HLA molecules. And LCCOPD-related macrophages had a significant up-regulation of LY86, CYBA, CCL13, IFITM3, RNASE6 and LYZ, which participated in immune response, including humoral immune response and innate immune response. COPD-related macrophages had a stronger antigen processing and presentation, with a significant up-regulation of HSP90AA1, HSP90AA1, B2M, HLA-A, and HLA-C. (Supplementary information, Fig. S3C)
We further deciphered the expression identify-markers of different types of macrophages (Fig. 3E). we found that COPD-specific macrophage populations showed activation of CD86, CD68, and MARCO, which is highly overlapped with M1 macrophages. The LCCOPD-specific macrophage populations were characterized by high expression of PDGFB, CCL2, CCR2, CD200R1, CSF1R, and other M2 or TAM phenotypes markers. In the COPD environment, macrophages had a conventional pro-inflammatory phenotype, namely M1 macrophages, as the main population, while in the tumor-bearing state, macrophages had an anti-inflammatory phenotype that promotes cancer as the main population, such as M2 type and tumor-associated macrophage population. This differential distribution was consistent with the specific pathophysiological environment.
Pseudo-time trajectory analysis revealed differentiation route of tissue-specific
macrophages
To further investigate the potential differentiation direction of macrophage populations in the two environments, we performed a pseudo-time trajectory analysis of macrophages by different conditions. We subset the initial dataset and analyzed differentially expressed genes (DEGs) across clusters to obtain a list of differentially expressed genes. Subsets of data from different tissues were then imported into the Monocle R package for separate pseudo-time trajectories, following a standard workflow using default parameters.
In COPD groups, macrophages showed five differentiation states and four differentiation outcomes (Fig. 3F). In the LCCOPD environment, may due to the existence of tumor cells, the differentiation path of macrophages was significantly increased, and there were seven different differentiation states and five differentiation outcomes (Fig. 3G). By downloading and marking the lung macrophage markers from the CellMarker database, we analyzed several cell markers on the trajectory (Supplementary information, Fig. S3D, E). Compared to COPD groups, tumor-related macrophages exhibit more diverse changes in the marker expression on the trajectory.
Combined with the clustering results, we found that COPD-related macrophages were dominated by M1 macrophages with stronger immune killing ability, and the killing factors such as TNF, IL1B, and IL6R gradually increased with the trajectory. LCCOPD-related macrophages were dominated by M2 macrophages involved in immune regulation, and the secretion of TNF, IL1B, and CXCL16 showed a downward trend with the trajectory (Fig. 3H, 3I).
For the functional differences of macrophages on the trajectories, we performed a time-sequential differential gene heatmap for DEGs. All DEGs were clustered into four clusters according to the expression changes on the trajectory. In the COPD group (Fig.4A), progressively increased expression happened in the most temporally differential genes. This indicated that most of the differential pathways were gradually enhanced in COPD. In the LCCOPD groups (Fig. 4B), cluster 3 whose expression was gradually down-regulated accounted for the most differential genes, indicating most of the pathway differences were gradually weakened.
To study the main biological processes of the genes in these gene clusters, we carried out GO and KEGG pathway enrichment analysis for each cluster to obtain the pathways with significant differences. In the COPD groups (Fig.4C), COPD-related macrophages exhibit progressively increased proliferative activities, proliferation-related pathways, such as cell cycle, DNA replication, mitosis, were enhanced along trajectories. As for metabolic profiles, the lipid metabolism-related processes gradually decreased from the early partial expression to changes in lipid catabolism, such as lipid catabolism, fatty acid metabolism, and cholesterol metabolism. For the LCCOPD group (Fig.4D), pathway enrichment on the trajectories was significantly different from the COPD group. Compared with the significantly up-regulated cell proliferation-related signals in the COPD group, in the LCCOPD group, the cell proliferation signal showed a short-term mid-term increase and remained in a low-proliferation state in the later stage. The most significant changes happened on the chemotactic ability of the macrophage along the trajectory, such as T cells, lymphocytes, and other leukocytes. The increased chemotactic ability shaped the special tumor immune microenvironment successfully. Different from COPD-related macrophages, the metabolic changes of lipid metabolism in LCCOPD-related macrophages were mainly fatty acid anabolism.
To investigate the clinical role of the different trajectories of macrophages identified in the present study, we selected different gene clusters for subsequent analysis. In COPD groups, gene cluster2 mainly represented a set of genes that were gradually down-regulated along the trajectory, while gene cluster 4 in the LCCOPD group was a set of genes that were gradually up-regulated along the trajectory. Therefore, we first performed the analysis for these two clusters of genes with opposite expression trends (Fig. 4E). It can be seen that there are eight crossover genes in the two clusters, namely FN1, PERP, SFTPA2, FGFBP2, KLRB1, SPON2, GZMK, and TUBA4A. Therefore, we mapped the protein interaction network for these eight repetitive genes and found that the strength of the interaction between them was not significant (Fig. 4F). TCGA survival analysis indicated that FN1, KLRB1, SPON2, and TUBA4A had significant prognostic differences in non-small cell lung cancer (Fig. 4G).
COPD status controlled the reprogramming of the immunosuppressive environment
The cancer microenvironment is characterized by a high degree of immunosuppression, and the main effector cell populations are highly infiltrating immunosuppressive cell populations, such as T-reg cells and M2 macrophages. In our data, we found that the number of T-regs and M2 macrophages in LCCOPD was significantly increased compared with COPD, and the related immunosuppressive pathways were also upregulated. It suggested that in LCCOPD groups, the microenvironment had a more powerful immunosuppressive ability. However, the association of the reshaping of this immunosuppressive microenvironment with COPD overlapping states need to be further complemented. To this end, we externally incorporated single-cell sequencing data from four sets of different settings to verify the contribution of COPD status to the immunosuppressive environment. We analyzed the external data from the normal lung tissue, COPD lung tissue, pure lung cancer tissue, and COPD combined lung cancer tissue and marked them with the Other-Nor group, Other-LC group, Other-COPD group, and Other-LCCOPD group.
First, we performed batch correction for all single-cell sequencing data from different tissue sources to eliminate the between-group differences and systematic errors. Dimensionality reduction analysis was performed on the T-reg population of the four groups as well as the macrophage population to validate the corrected results. A high degree of integration of cells in the four populations was observed, demonstrating that data correction can be used for subsequent analysis (Fig. 5A, B).
In terms of cell composition, we normalized the proportion of T-regs and macrophages in each group (Supplementary information, Fig. S4B). There was indeed a significant change in the number of cells, specifically, the number of T-regs and macrophages in the simple COPD group was higher than the normal lung tissue, and the number of T-regs and macrophages in the LCCOPD group was higher than the LC tissue. The number of T-regs in tumor-bearing tissues was significantly higher than that in non-tumor tissues. Compared with the population of macrophages, there was no obvious trend on the influence of COPD, which may be related to the abnormal number of macrophages.
The distribution difference study indicated that a large number of immune populations aggregated in the pure COPD state, and then we compared the intensity of pathways involved in the shaping of the inhibitory environment for T-regs and M2 macrophages (Fig. 5C, 5D). Compared to Other-LC groups, we found out that, in the presence of COPD, T-regs showed a dramatic decline in functional pathways and metabolic profiles. while all metabolic pathways were up-regulated in the other-LC group alone, only DNA demethylation, inflammatory response, and KRAS pathways were up-regulated in Other-LCCOPD groups, indicating that the existence of COPD environment leads to a certain decline in the function of T-regs in cancer tissues, including the expression of CTLA-4, PD-1, and the secretion of cytokines.
And in macrophages, the overlapping state of COPD also led to a large-scale decline in macrophage function in the lung cancer microenvironment. Overall, T-regs, as well as macrophages, exhibited significantly functional suppression in the lung cancer microenvironment with COPD overlap compared to LC. The secretion of cytokines and many metabolic pathways were significantly down-regulated.
In the tumor-bearing state, overlapping COPD results in a decline in cellular function. Compared to Other-Nor groups (Fig. 5C), only a small number of pathways were up-regulated in Foxp3+ T-reg cells in the other-COPD group, such as Notch pathway, KARS pathway, β fatty acid oxidation pathway, and DNA demethylation, the only major cytokine pathway was an increase in the intensity of the anti-inflammatory factor-IL10. The notch pathway was related to cancer. The relationship between K-Ras and lung cancer had been discovered. IL10 as an immunosuppressive factor participated in the survival, proliferation, and anti-apoptotic activity of non-small cell lung cancer. Based on these results, we concluded that T-reg cells exhibit significant immunosuppressive, cancer-promoting, and anti-cell death activities in the COPD state before the onset of cancer. For macrophages (Fig. 5D), mitochondrial respiratory dysfunction is one of the factors that induce the inflammatory differentiation of macrophages and change their anti-cancer phenotype35. By comparison, we found that a few macrophage pathways, such as the bile acid pathway, hypoxia pathway, and angiogenesis pathway, were up-regulated in other-COPD groups. Oxygen-dependent and mitochondrial-dependent metabolic pathways, such as the oxidative phosphorylation pathway, were upregulated in other-Nor groups. Mitochondrial metabolic program changes are often related to mitochondrial respiratory dysfunction36.
T-regs and macrophages had extensive metabolic profile changes in COPD tissue and LCCOPD tissue, and there was a correlation with function. For further verification, we selected lung sections from lung cancer patients with COPD and simple lung cancer patients, and performed immunofluorescence localization on T-regs and macrophages. The immunofluorescence co-localization (Figure.6) showed that compared with pure COPD and pure lung cancer tissue, T-regs with high expression of SREBP1(Figure6.A) and M2 macrophages with high expression of GLS1(Figure6.B), a key rate-limiting enzyme in glutamine metabolism, could be found in lung cancer tissue samples with COPD overlapped. The SREBP metabolic pathway in T-regs can regulate the immunosuppressive ability of T-regs . When the SREBP metabolic pathway in T-regs was blocked, the intensity of immunosuppression in the lung cancer microenvironment was significantly reduced, suggesting that the enhancement of the SREBP metabolic pathway was related to the enhanced immunosuppression of T-regs. Glutamine metabolism was considered to be one of the important metabolic markers of M2 polarization in macrophages. The high enrichment of T-regs and M2 macrophages confirmed the differences in the number of microenvironments in previous results, and also revealed some differences in the immunosuppressive capabilities of the three lung disease settings.