Integrative bulk and single-cell transcriptomics identify a recurrence-prognostic risk model with exploratory associations with treatment-response phenotypes in colorectal cancer

Abstract

Chemoresistance remains a major obstacle to improving outcomes in colorectal cancer (CRC), particularly among patients receiving 5-fluorouracil (5-FU)-based therapy. Here, we developed a READ-derived response-associated gene signature for recurrence stratification and performed exploratory cross-cohort evaluation in additional colorectal cancer cohorts, while further exploring its association with treatment-response phenotypes. In the GSE190826 cohort, pathological complete response (pCR) versus non-pCR status was used as a clinically annotated treatment-response phenotype for downstream screening. Using differential expression analysis in TCGA-READ, we identified 4277 differentially expressed genes and constructed a weighted gene co-expression network analysis (WGCNA) to define clinically relevant modules. Genes from READ-related modules were further screened by logistic regression in GSE190826 and univariate Cox regression in TCGA-READ, yielding six prognostic genes associated with pathological non-response. Least absolute shrinkage and selection operator (LASSO)–Cox regression subsequently refined these candidates to a four-gene signature, including ENHO, BLACAT1, LEMD1, and ASNS, which was used to establish a READ-derived response-associated gene (RDRG) risk score. High-risk patients showed poorer recurrence-free survival in the TCGA-READ training set and an exploratory TCGA-READ testing subset, and exploratory evaluations in additional colorectal cancer cohorts (TCGA-COAD and GSE17536) showed similar recurrence-stratification trends, with moderate time-dependent AUC values in TCGA-READ (3-year AUC = 0.766; 5-year AUC = 0.714). Multivariate Cox analyses showed that the risk score remained associated with RFS after adjustment for clinical covariates, and a nomogram incorporating disease stage provided additional exploratory prognostic stratification. The RDRG score also showed only modest and exploratory associations with treatment-response phenotypes across multiple datasets. Single-cell RNA-seq analysis revealed that malignant cells contributed the highest risk scores, with enrichment of epithelial-mesenchymal transition (EMT), KRAS signaling, coagulation, and hypoxia-related programs in high-risk tumor cells. Among the four genes, LEMD1 showed the strongest malignant-cell enrichment and was further examined in functional assays, where LEMD1 knockdown reduced CRC cell migration and increased 5-FU sensitivity, supporting its role as a candidate functional contributor rather than a validated mechanistic driver. Collectively, the RDRG risk score may serve as a READ-derived model for recurrence stratification, with exploratory cross-cohort applicability in broader colorectal cancer cohorts and exploratory associations with treatment-response phenotypes. The same risk-score formula was applied across cohorts, whereas high- and low-risk thresholds were cohort-specific rather than fixed from the training cohort.