Single-Cell RNA Sequencing Reveals the Role of Epithelial Cell Marker Genes in Predicting the Prognosis of Colorectal Cancer Patients

Abstract

Single-cell RNA sequencing (scRNA-seq) is increasingly used in studies on gastrointestinal cancers. This study investigated the prognostic value of epithelial cell-associated biomarkers in colorectal cancer (CRC) using scRNA-seq data. We downloaded and analysed scRNA-seq data from four CRC samples from the Gene Expression Omnibus (GEO), and we identified marker genes of malignant epithelial cells (MECs) using CRC transcriptome and clinical data downloaded from The Cancer Genome Atlas (TCGA) and GEO as training and validation cohorts, respectively. In the TCGA training cohort, weighted gene correlation network analysis, univariate Cox proportional hazard model (Cox) analysis, and least absolute shrinkage and selection operator regression analysis were performed on the marker genes of MEC subsets to identify a signature of nine prognostic MEC-related genes (MECRGs) and calculate a risk score based on the signature. CRC patients were divided into high- and low-risk groups according to the median risk score. We found that the MECRG risk score was significantly correlated with the clinical features and overall survival of CRC patients, and that CRC patients in the high-risk group showed a significantly shorter survival time. The univariate and multivariate Cox regression analyses showed that the MECRG risk score can serve as an independent prognostic factor for CRC patients. Gene set enrichment analysis revealed that the MECRG signature genes are involved in fatty acid metabolism, p53 signalling, and other pathways. To increase the clinical application value, we constructed a MECRG nomogram by combining the MECRG risk score with other independent prognostic factors. The validity of the nomogram is based on receiver operating characteristics and calibration curves. The MECRG signature and nomogram models were well validated in the GEO dataset. In conclusion, we established an epithelial cell marker gene-based risk assessment model based on scRNA-seq analysis of CRC samples for predicting the prognosis of CRC patients.

Figure 1

Characterisation of scRNA-seq from 15,465 cells. (a) The cells were classified into 18 subsets using the t-SNE algorithm. (b) Distribution ratios of the cell subsets in the four CRC samples.

Figure 2

Characteristics of the cell subsets. (a) Annotated cell subsets. (b) Proportions of the various types of cells in the four CRC samples. (c) Heat map of the five most differentially expressed marker genes in each cell subset.

Figure 3

CNV analysis of epithelial cells from CRC patients. (a) Heat map of CNV in epithelial cells from the four CRC samples. (b) CNV score distribution among different epithelial cell subtypes.

Figure 4

Functional enrichment analysis of marker genes in the MEC subsets. Results of GO function (a) and KEGG pathway (b) enrichment analyses.

Figure 5

MEC marker genes associated with CRC. (a) Dendrogram of MEC subset marker genes obtained by WGCNA according to colour. (b) Correlations between characteristic genes of different modules and CRC.

Figure 6

Identification and validation of MECRGs in the TCGA-CRC training and GSE17538 validation cohorts. (a) Univariate Cox analysis was used to screen MECRGs with prognostic significance. (b) Distribution of MECRGs in cell subsets. (c and d) Kaplan–Meier survival curve showing the prognostic value of the MECRG signature in the training cohort (c) and the validation cohort (d). (e and f) Distribution of the MECRG risk scores and survival status of CRC patients in the training cohort (e) and validation cohort (f). ROC curve representing the efficiency of the MECRG signature in predicting 1-, 3-, and 5-year OS in CRC patients in the training cohort (g) and validation cohort (h).

Figure 7

Independent prognostic value of the MECRG signature. Univariate and multivariate Cox analyses of the MECRG risk score in the TCGA-CRC training cohort (a and b) and GSE17538 validation cohort (c and d).

Figure 8

Analysis of the relationships between the MECRG signature and clinical features using the TCGA-CRC training and GSE17538 validation cohorts. (a) Survival analysis of the MECRG signature in clinical features based on the training cohort. (b) Survival analysis of the MECRG signature in clinical features based on the validation cohort. (c and d) Heat maps showing the correlation between MECRG risk grouping and TNM stage in the training cohort (c) and validation cohort (d). ^∗∗∗P < 0.001, ^∗∗P < 0.01, and ^∗P < 0.05.

Figure 9

Functional enrichment analysis of the MECRG signature and distribution of enriched pathways in the cell subsets. (a) GSEA of the MECRG high- and low-risk groups. (b) Pathways enriched in the cell subsets of the MECRG risk groups.

Figure 10

Construction of a nomogram based on the MECRG signature. ROC curve analysis of the nomogram (a) using the TCGA-CRC training cohort (b) and the GSE17538 validation cohort (c). (d and e) Calibration curves for the nomogram using the training cohort (d) and the validation cohort (e).