Bioinformatics analysis and identification of cuproptosis-related long non-coding RNAs in colorectal cancer

Abstract

Objective

Identifying precise biomarkers for colorectal cancer (CRC) detection and management remains challenging. Here, we developed an innovative prognostic model for CRC using cuproptosis-related long non-coding RNAs (lncRNAs).

Methods

In this retrospective study, CRC patient transcriptomic and clinical data were sourced from The Cancer Genome Atlas database. Cuproptosis-related lncRNAs were identified and used to develop a prognostic model, which helped categorize patients into high- and low-risk groups. The model was validated through survival analysis, risk curves, independent prognostic analysis, receiver operating characteristic curve analysis, decision curves, and nomograms. In addition, we performed various immune-related analyses. LncRNA expression levels were examined in normal human colorectal epithelial cells (FHC) and CRC cells (HCT-116) using quantitative polymerase chain reaction (qPCR).

Results

Six cuproptosis-related lncRNAs were identified: ZKSCAN2-DT, AL161729.4, AC016394.1, AC007128.2, AL137782.1, and AC099850.3. The prognostic model distinguished between high-/low-risk populations, demonstrating excellent predictive ability for survival outcomes. Immunocorrelation analysis showed significant differences in immune cell infiltration and functions, immune checkpoint expression, and m⁶A methylation-related genes. The qPCR results showed significant upregulation of ZKSCAN2-DT, AL161729.4, AC016394.1, AC007128.2 in HCT-116 cells, while AL137782.1 and AC099850.3 expression patterns were significantly downregulated.

Conclusion

Cuproptosis-related lncRNAs can potentially serve as reliable diagnostic and prognostic biomarkers for CRC.

Introduction

Colorectal cancer (CRC) has been identified as one of the most prevalent malignant tumor types, ranking as the third most frequent cancer, second most severe risk to human life, and fourth deadliest cancer worldwide. Notably, new favorable subsets of cancers of undefined origin (CUP) have emerged, including CRC CUP. This new clinical entity is treated as CRC, which contributes to the current increased incidence of CRC. ¹ Both the CRC incidence and mortality rates have shown a consistent upward trend annually, resulting in 881,000 and 940,000 recorded deaths globally in 2018 and 2020, respectively.^2–6 Current CRC treatment methods include surgical treatment, targeted therapy, and other approaches.^7,8 Because bloody stools and anemia do not typically appear until the late stages of CRC, patients often miss the optimal treatment opportunity during their initial diagnosis when the tumor is most susceptible to further deterioration and metastasis. Patients who receive an early diagnosis of CRC could potentially reduce their mortality rate by 90%.^9–11

The goals of CRC screening programs around the world are to increase early detection, reduce patient mortality, and improve prognosis of this disease.^12–15 However, no consensus has been reached if age affects CRC survival outcomes. The prognosis of older patients may be confounded by factors such as differences in stage at presentation, tumor site, preexisting comorbidities, and type of treatment received. ¹⁶ Classical biomarkers, like carcinoembryonic antigen and carbohydrate antigen 19-9, have been integrated into medical practice following technological advancements. However, these markers still do not meet all the necessary requirements for widespread clinical implementation. ¹⁷ Among several major genetic mutations in CRC, RAS mutations correlate with the oncological aggressiveness and response to chemotherapy. There is also growing evidence that inflammation drives disease development. As a result, many studies have investigated the predictive and prognostic roles of various blood-based inflammatory markers, including the neutrophil–lymphocyte ratio, lymphocyte–monocyte ratio, and platelet–lymphocyte ratio. Finally, microRNAs (miRNAs) can have roles as tumor suppressor genes and oncogenes, with their diagnostic, prognostic, and predictive implications now being explored. ¹⁸ To ensure tailored care for individuals with CRC, it is essential to identify specific biomarkers that can aid in the early detection and assessment of treatment outcomes for advanced stages of the disease.

The regulation of cell death by copper is a unique and emerging concept in the cell biology field. This form of cell death sets itself apart from the established pathways of cellular demise. ¹⁹ Called cuproptosis, this phenomenon occurs when copper ions directly attach to lipid acylated sections of the tricarboxylic acid cycle during mitochondrial respiration. This process triggers the clumping together of lipid acylated proteins and decreased iron-sulfur cluster protein levels, which induces proteotoxic stress and eventual cell demise.^20,21 Cuproptosis therefore offers a distinct viewpoint for exploring cancer treatment methods with significant promise.

Long non-coding RNAs (lncRNAs) are a distinct group of RNA molecules categorized by their extensive length (more than 200 bases) and lack of protein-coding ability. They primarily accumulate in the cell nuclei and cytoplasm. LncRNAs are transcribed from spacer or gene regions as main products of RNA polymerase II. ²² They can compete with other endogenous RNA molecules, such as miRNAs, to regulate gene expression and thereby influence a wide range of biological processes. ²³ LncRNAs are associated with numerous methods of tumor cell death. For example, lncRNAs can influence the miR-513a-3p/SLC7A11 axis, which in turn can induce ferroptosis and inhibit esophageal squamous cell cancer progression. ²⁴ Previous studies have demonstrated that lncRNAs can enhance NLRP3-mediated pyroptosis in hepatocellular carcinoma cells through its interaction with the miR-34a/SIRT1 axis. ²⁵ Research findings have demonstrated that lncRNAs can also impede the apoptosis and autophagy abilities of bladder cancer cells by affecting the PI3K/AKT/mTOR signaling pathway. ²⁶ Specific research on the connection between cuproptosis and lncRNAs is still in its early stages. The number of lncRNA-related studies has notably expanded in recent years, largely driven by its promise as a biomarker for the clinical detection, management, and prognosis of CRC.^27–29 Data have suggested that lncRNAs may control the glioblastoma (GBM) microenvironment, specifically by influencing the activity of cytokines and growth factors. The expression profiles of different lncRNAs can also be used to determine the grade of glioma and its subtype. Importantly, the dynamics of lncRNAs circulating in the blood can help determine the GBM patient prognosis and response to treatment. ³⁰

In this study, we analyzed the transcriptomic and clinical survival data of CRC patients using The Cancer Genome Atlas (TCGA) database. We aimed to identify any cuproptosis-associated lncRNAs through gene co-expression analysis. A one-way Cox analysis was used to screen for prognostically relevant data, after which we constructed prognostic models and assessed if they could function as independent prognostic factors. Subsequently, Gene Set Enrichment Analysis (GSEA) was conducted on the prognostic model, with the results then validated using receiver operating characteristic (ROC) curve analysis and decision curve analysis (DCA). We compared the variations in immunological checkpoint genes, immune infiltration, immune function, and m⁶A methylation-related genes, then validated the results through in vitro experiments. Ultimately, this study introduces a new direction for the development of gene therapy, targeted therapy, and individually tailored therapeutic approaches for CRC patients.

Methods

Data acquisition and processing

This study was a retrospective study. The reporting of this study conforms to STROBE guidelines. ³¹ The purpose of TCGA is to support research by allowing the comparison of data from 34 prevalent human cancer types with normal biological tissue samples that are associated with the tumors. Researchers can access various types of data from this database, including information on genomic copies, outlier genetics, and gene expression patterns. Both genome and transcriptome sequencing data, which include lncRNAs and miRNAs, as well as clinical data such as age, sex, ethnicity, pathological histology categorization, tumor stage, recurrence, survival time, and clinical outcome, are available for free download. Users can access the complete TCGA dataset through the Genomic Data Commons Data Portal (https://portal.gdc.cancer.gov/), which is provided by the National Cancer Institute. This research used several datasets, including the RNA sequencing (RNA-seq) dataset, clinical dataset, and metadata dataset, containing 473 tumor samples and 41 normal samples. We have de-identified all patient details. This study was examined by the Anxi County Hospital Ethics Review Committee and was granted an ethical exemption.

The RNA-seq dataset and metadata information were merged using the Perl programming language to construct the gene expression matrix file. The gene expression matrix obtained from the TCGA database was completed using the FPKM tool. Genes with the same gene name were combined into a single row by calculating the mean value, while genes with zero expression in any sample were excluded. The following clinical data were retrieved from the dataset using Perl programming language: age, sex, TMN stage, survival status, and other variables.

Differential analysis and screening of cuproptosis-related lncRNAs

Cuproptosis-related genes were identified from the literature. To determine the differential expression of cuproptosis-related genes between tumor samples and normal samples, we used the “limma” package in R software (www.r-project.org). Our analysis employed the criteria of |log2 (fold change)| > 1 and P < 0.05 to identify cuproptosis-related genes. In addition, we used the “ggplot2” package in R software to perform gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses.

An association analysis was conducted to identify co-expression associations between cuproptosis-related genes and lncRNAs. The “coxph” function from the “survival” package in R software was used to visualize the expression files of single-factor significant lncRNAs. From the screening results, two categories of lncRNA variables were identified: risk and protective lncRNAs. Using the expression data of these lncRNAs, forest plots were generated to show the risk ratio and identify lncRNAs associated with prognosis. The aforementioned significant single-factor lncRNA expression data underwent regression analysis using the least absolute shrinkage and selection operator (Lasso) function, which led to the identification of significant lncRNAs. Subsequent research focused on these lncRNAs, with network maps created using Cytoscape software.

Prognostic model construction

To calculate the risk score for each sample, the lncRNA expression levels were multiplied by the outcomes derived from the Lasso Cox regression analysis. Subsequently, the samples were divided into two groups based on the median value of the risk score. By merging the clinical data of the high- and low-risk groups and examining their differences in prognosis and survival, we could evaluate the discriminative ability of this model for distinguishing between individuals.

Single-factor and multi-factor independent prognostic analyses were conducted using the risk scores and all clinical traits collectively to determine if this model could serve as an independent prognostic factor for predicting patient survival time and status, regardless of other clinical factors.

Risk score DCA was performed using the “ggDCA” software. The prediction accuracy of the prognostic model was evaluated using the ROC function. A larger area under the ROC curve (AUC) value indicates an improved forecast accuracy and efficiency for the model.

Using the screened clinical data, a scatter plot diagram was generated using the rms function. The variables were categorized by sex, age, risk score, total risk score, linear predictive value, and the 1-, 3-, and 5-year survival rates. Line segments with scales were used to depict the relationships between the variables based on a specific ratio plotted in the same plane. This nomogram can be employed for predicting patient survival.

GSEA enrichment analysis

When processing the standardized expression matrix and conducting the GSEA analysis using the GSEA 4.0.3 software, the following criteria were used for filtering: a normalized enrichment score >1, a nominal P < 0.05, and a false discovery rate q < 0.25.

Immune cell infiltration and immune function analyses

Immune infiltration data from CRC patients in the TCGA database were obtained from the TIMER database (http://timer.cistrome.org/). To analyze immune cell infiltration, both the TIMER and CIBERSORT methods were used, focusing specifically on B cells, CD4⁺ T cells, and CD8⁺ T cells. R software was employed to visualize the immune infiltration. The “heatmap” package in R was used to generate a heat map.

Thirteen distinct immune function scores were aggregated, then the immune function score variations were calculated. Subsequently, significant differences in immune function scores were identified. Then, single sample GSEA (ssGSEA) analysis was conducted using the gene set variation analysis command, with the ssGSEA scores adjusted to obtain the immunological scores for the 13 immune functions. The differential analysis of immune function expression was carried out using the sigVec command, with all immune functions examined for variability. A box line graph was generated to assess the significance of the observed differences.

Immune checkpoint and m⁶A-related gene analyses

The m⁶A-related genes and immunological checkpoints were identified by reviewing the literature. To assess their expression patterns, we used the “limma” package in the R software. In addition, the R program “ggplot2” was employed to compare the immunological checkpoints that exhibited significant differences. Box plots were generated to compare the outcomes of the m⁶A-related genes.

In vitro experimental validation

Human normal colorectal epithelial cells (FHC) and CRC cells (HCT-116) were purchased from American Type Culture Collection (ATCC; Manassas, VA, USA) and cultured in Dulbecco’s Modified Eagle Medium containing 10% fetal bovine serum (Life Technologies, Beijing, China) at 37°C in a 5% CO₂ incubator. Total RNA was extracted from the cells in each group using TRIzol reagent (Thermo Fisher Scientific, Waltham, MA, USA), then the RNA purity and concentration were determined using the Nano Drop 2000 system (Thermo Fisher Scientific). The total RNA was reverse transcribed into cDNA using Prime ScriptRT Master Mix reagent (Thermo Fisher Scientific), then quantitative polymerase chain reaction (qPCR) was performed using the miScript SYBR Green PCR kit (Qiagen, Hilden, Germany). The relative expression levels of the target genes were calculated by the 2^−ΔΔCt method with β-actin as the internal reference gene. The primer sequences used in this study are listed in Table 1.

Table 1.

Primer sequences used for qPCR.

RNA molecule	DNA sequence (5′ to 3′)
β-actin mRNA	Forward Primer TCCGGCACTACCGAGTTATC
	Reverse Primer GATCCGGTGTAGCAGATCGC
AL161729.4	Forward Primer CCTTTCTGGTCTGTGGGAGG
	Reverse Primer TCAGTGTTCTTGCCAACTGC
ZKSCAN2-DT	Forward Primer GGATTCTAGGGAACCTGCCG
	Reverse Primer AGGTGACACTGTTGCCTCTG
AL137782.1	Forward Primer CAGGGAGTACAGCCTTGTCG
	Reverse Primer TTTCTCCTGCACTCTCCACG
AC016394.1	Forward Primer CGAAGGCAGCAGTGTTTGTT
	Reverse Primer GGGAGCATTGATGAAACCTATAC
AC099850.3	Forward Primer CCCAGGTTCAAGTAACTGGGAC
	Reverse Primer GACAATGCCTTGCCAAGGAATC
AC007128.2	Forward Primer ACCAGCCGGAATAGCAACC
	Reverse Primer GTCTGCTGGAGGGTAGAAGC

Statistical analysis

The data were analyzed using R software. Group differences were analyzed using one-way ANOVA. Two-by-two comparisons were conducted using the LSD test, unpaired t-test, and Wilcoxon test for non-normally distributed data. Differences were considered statistically significant when P < 0.05.

Results

Identification of cuproptosis-related lncRNAs

After reviewing the literature, we identified ten cuproptosis-related genes: FPX1, LIAS, LIPI1, DLD, DLAT, PDHA1, PDHB, MTF1, GLS, and CDKNIA. The GO enrichment data are shown in Figure 1(a). For biological processes, the cuproptosis-related genes were primarily involved in the acetyl-coA biosynthetic process from pyruvate, acetyl-coA biosynthetic process, and acetyl-coA metabolic process. For cellular components, these genes were mainly involved with the mitochondrial matrix, oxidoreductase complex, and mitochondrial protein-containing complex. For molecular function, these genes were primarily involved in oxidoreductase activity and sulfurtransferase activity. Figure 1(b) shows the KEGG enrichment results. The citrate cycle, pyruvate metabolism, glycolysis/gluconeogenesis, and other processes were primarily impacted by the cuproptosis related-genes. Gene co-expression analysis was then conducted on the genes and lncRNAs that were differentially expressed in CRC, resulting in 2578 cuproptosis-related lncRNAs. The clinical survival data were integrated with these findings and assessed using one-way Cox analysis (Figure 1(c)). Overall, 23 lncRNAs were selected for this research, which were the key participants in the Lasso regression analysis. Subsequently, the lncRNAs underwent further filtering and visualization using Cytoscape software (Figure 1(d)). This led to the identification of six significant lncRNAs: AC016394.1, AC099850.3, ZKSCAN2-DT, AL137782.1, AC007128.2, and AL161729.4.

Figure 1.

Cuproptosis-related long non-coding RNAs (lncRNAs) were screened and subjected to enrichment analysis. (a) Gene ontology (GO) enrichment analysis. (b) Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis. (c) Forest plot depicting the one-way Cox analysis results and (d) network map illustrating the interactions between lncRNAs (yellow) and genes (green).

Constructing the prognostic model

The samples were classified into high- and low-risk groups according to their risk score and lncRNA expression levels. Kaplan–Meier survival curves were then generated, as shown in Figure 2(a). The model exhibited increased precision in forecasting the survival outlook of CRC patients, as evident from the reduced survival forecast in the high-risk classification group (P < 0.001). The scatter plots in Figure 2(b–d) depict the correlations between the risk score and model-generated survival status prognosis. The results showed a clear correlation between the increased patient mortality rates and higher risk scores, validating the model's predictions. A visual representation as a risk heat map was then produced using the levels of lncRNAs included in model construction for every individual sample. The high-risk lncRNAs identified were AL161729.4, ZKSCAN2-DT, AC016394.1, and AC007128.2, while the low-risk lncRNAs were AL137782.1 and AC099850.3.

Figure 2.

Survival and risk curves. (a) Survival curves. (b) Risk group curves. (c) Patient survival status scatter plot and (d) risk heat map.

Evaluation of the prognostic model performance

The model’s position as an independent prognostic factor was validated through single-factor and multi-factor independent prognostic analyses, as displayed in Figure 3(a–b), which indicated that the model's predictive ability had increased. Risk score ROC curve analysis was conducted using the “timeROC” package in the R software, as depicted in Figure 3(c–d). In this analysis, a higher AUC value indicates greater accuracy. The AUC values were determined to be 0.648, 0.713, and 0.717 at 1, 3, and 5 years, respectively. The ROC curve analysis of this model, in combination with other clinical traits, revealed that the model was more accurate for predicting CRC patient survival. The DCA data, shown in Figure 3(e), illustrates that the risk score decision curve was significantly different from the All curve, indicating that the model was a stronger predictor of survival compared with other clinical characteristics. The risk score decision curve was the farthest away from the All curve (Figure 3(e)), demonstrating that the model predicted survival better than other clinical traits examined. Figure 3(f) displays a nomogram illustrating various clinical characteristics, such as gender, stage, age, and risk score, plotted against the 3-year survival rate (0.764%), 5-year survival rate (92.5%), and 1-year survival rate (43.0%).

Figure 3.

Evaluation of prognostic model performance. (a) Single-factor independent prognostic analyses. (b) Multi-factor independent prognostic analyses. (c) Receiver operating characteristic (ROC) curve. (d) Examination of the model's ROC curves at 1, 3, and 5 years. (e) Decision curves and (f) Nomogram.

GSEA enrichment analysis

The GSEA analysis outcomes are displayed in Figure 4. Oocyte meiosis, cell cycle, nod-like receptor signaling pathway, ubiquitin-mediated proteolysis, and progesterone-mediated oocyte maturation were significantly active in the low-risk group, while basal cell carcinoma, tyrosine metabolism, maturity-onset diabetes of the young, taste transduction, and linoleic acid metabolism were significantly active in the high-risk group.

Figure 4.

Gene set enrichment analysis (GSEA).

Immune-related analysis

The TIMER and CIBERSORT methods were used to analyze immune cell infiltration, including B cells, CD4⁺ T cells, CD8⁺ T cells, and other immune cells (Figure 5(a)). APC co-inhibition, APC co-stimulation, CCR, Checkpoint, Cytolytic activity, Inflammation-promoting, MHC class I, Parainflammation, T cell co-inhibition, T cell co-stimulation, Type I IFN Response, and Type II IFN Response were found to indicate significantly different immune functions, as shown in Figure 5(b). Thirteen m⁶A-related genes and 47 immunological checkpoints were identified by integrating data from the literature. The “limma” package of R software was used to compare the expression patterns of m⁶A-related genes, with the results revealing that YTHDF2, WTAP, and HNRNPC were significantly different (Figure 5(c)). ICOS, CTAL4, CD80, CD28, TNFRSF9, TIGIT, CD274, PDCD1LG2, TNFRSF14, TNFRSF25, HHLA2, IDO1, HAVCR2, CD48, BTLA, CD244, NRP1, and CD86 were also found to be significantly different (Figure 5(d)).

Figure 5.

Immune-related analysis. (a) Immune cell infiltration. (b) Immune function analysis. (c) M⁶A methylation-related gene expression analysis and (d) immune checkpoint expression analysis.

qPCR results

Figure 6 demonstrates that the expression levels of ZKSCAN2-DT, AL161729.4, AC016394.1, and AC007128.2 were significantly upregulated in the HCT-116 CRC cell line compared with the control FHC cells (P < 0.001). Conversely, the expression levels of AL137782.1 and AC099850.3 were significantly downregulated compared with the control cells (P < 0.001).

Figure 6.

qPCR results. Compared with the control group, *P < 0.05; **P < 0.01; ***P < 0.001.

Discussion

Cuproptosis is distinct from oxidative stress and does not depend on reactive oxygen species (ROS) pathways like other forms of cell death do, such as ferroptosis, apoptosis, and necroptosis. ³² Although CRC maintains low ROS levels, which can promote tumor development and boost tumor cell resistance, it is crucial to investigate cuproptosis-related genes and their possible use as CRC biomarkers.^33,34 Numerous studies have demonstrated how lncRNAs can serve as tumor molecular markers, playing significant roles in CRC invasion, growth, and apoptosis.^35–37 A new gene therapy-based therapeutic approach for CRC patients may be provided through examining cuproptosis-related lncRNAs. This will offer fresh perspectives on the diagnosis, care, and prognosis of this disease.

Cuproptosis-related lncRNAs, including AL161729.4, ZKSCAN2-DT, AC016394.1, AC007128.2, AL137782.1, and AC099850.3, were screened in this study. We determined the high-risk lncRNAs to be AL161729.4, ZKSCAN2-DT, AC016394.1, and AC007128.2, while the low-risk lncRNAs were AL137782.1 and AC099850.3. AC099850.3 has shown potential as a tumor prognosis biomarker.^38–40 The present research's findings support those of earlier studies, emphasizing the dependability of the screened biomarkers.

Here, Lasso regression analysis was used to construct a prognostic model, which supported sample classification into high- and low-risk groups using the median risk score. The resulting Kaplan–Meier survival curves clearly indicated that the low-risk group had notably higher survival rates compared with the high-risk group, which validated the efficacy of our model. The risk curves also revealed higher patient mortality rates with increased risk scores, aligning with the Kaplan–Meier survival curve data. Moreover, our assessment of the model’s predictive value, both through individual factors and multiple independent prognostic analyses, demonstrated that the risk score can potentially serve as a standalone prognostic indicator. The DCA and ROC curve analysis, which display CRC patient 1-, 3-, and 5-year survival rates, were used to confirm the predictive ability of the model. The effectiveness of the model was significantly higher than that of the conventional approach. In addition, enhanced precision was noted in identifying individuals at both elevated and reduced risk levels. GSEA was also used to explore the variations in biological functions and associated signaling pathways. Interestingly, basal cell carcinoma, tyrosine metabolism, and maturity-onset diabetes of the young emerged as notably enriched in the high-risk cohort, while oocyte meiosis, cell cycle, and nod-like receptor signaling pathway were significantly active in the low-risk group. The nod-like receptor signaling pathway is the primary cytosolic pattern recognition receptor directly implicated in innate immunity. It can identify various microbiological components and promptly stimulate immune cells, including dendritic cells, CD8⁺ T cells, and natural killer (NK) cells, to infiltrate the tumor microenvironment (TME) and eradicate tumor cells. ⁴¹ This research also examined the immune-related differences in deeper quantity.

Components of the TME are potential therapeutic targets for immunotherapy in different cancer types, including CRC, mainly because of their significant contributions to treatment resistance and disease advancement. ⁴² The TME is composed of a diverse array of cell types, including stromal cells, such as fibroblasts and adipocytes, as well as immune cells, like lymphocytes and macrophages. Within the TME, tumor infiltrating immune cells (TIICs) have a significant impact on the development and progression of cancer. ⁴³ Research on TIICs has advanced considerably in recent years. Multiple studies have shown that CRC patients with high TIIC infiltration levels tend to have a more favorable prognosis. ⁴⁴ Our data in this study are consistent with these previous investigations, demonstrating that the presence of immune cells, such as CD8⁺ T cells, CD4⁺ T cells, and NK cells, was notably higher in the low-risk group.

Immune function differences, as determined by the ssGSEA analysis, revealed that the high-risk group exhibited significantly higher levels of CCR and checkpoint expression. We investigated the role of various immune checkpoints, including both stimulatory and suppressive checkpoints, in regulating T cell recognition of antigens during the immune response. Immune checkpoints have been extensively explored in recent years for tumor therapy applications. The suppressive checkpoints act as brakes to prevent immune overactivation, which may harm normal cells and lead to autoimmune disorders. ⁴⁵ Tumor cells have developed an immunological escape mechanism by exploiting the binding capabilities of immune cells and proteins to evade T cell-mediated elimination. ⁴⁶ More potent immune checkpoint inhibitors have been developed, including anti-CTLA-4 and anti-PD-1 antibodies. These inhibitors work by disrupting tumor cell defense mechanisms against immune checkpoints, thereby boosting the immune system's ability to effectively eliminate cancer cells. ⁴⁷ Immune cell PD-L1 expression levels are significantly higher in mismatch repair (MMR)-deficient (MSI-H) CRC compared with MMR-proficient (MSI-L) tumors, with no differences among the different MSI-H molecular subtypes. The recommended screening methods for defective DNA mismatch repair includes immunohistochemistry (IHC) assays and/or MSI tests. However, there are challenges in obtaining usable data from the biological and technical heterogeneity of MSI testing. IHC testing of the mismatch repair machinery can reportedly yield different results for a given germline mutation, which has been suggested to be from somatic mutations. ⁴⁸ We evaluated immunological checkpoint variations and observed significant differences in ICOS and CTLA-4. These findings may have implications for individuals with CRC undergoing immunotherapy.

Upon its discovery in 1974, the prevailing modification observed in eukaryotic mRNA molecules was the introduction of a methyl group at the N6 position of adenine, commonly referred to as m⁶A. ⁴⁹ This process can impact mRNA transcription, translation, degradation, and other associated physiological and pathological processes. ⁵⁰ Downregulating genes associated with m⁶A may accelerate the growth of certain cancers. ⁵¹ In the current study, significantly diminished YTHDF2, WTAP, and HNRNPC expression levels were noted in the high-risk group. Previous research has demonstrated that YTHDF2, belonging to the YTH structural domain family 2, can potentially regulate the m⁶A methylation levels of the lncRNA XIST. Moreover, YTHDF2 may also play a role in promoting tumor metastasis, while also increasing CRC cell growth and invasion capabilities. ⁵² WTAP, also known as WT1 related protein, induces PDK4 m⁶A modifications, leading to enhanced CRC development, higher tumor cell proliferation, migration, and invasion rates, and accelerated xenograft growth. ⁵³ HNRNPC, also referred to as heterogeneous nuclear ribonucleoprotein-C, is an essential factor involved in regulating selective cleavage and polyadenylation mechanisms. These processes have a significant impact on cell proliferation and metastasis by influencing the expression of genes associated with carcinogenic properties. ⁵⁴ YTHDF2, WTAP, and HNRNPC might be useful resources for CRC gene therapy.

Finally, using the CRC cell line HCT-116, qPCR results from in vitro experiments showed that the ZKSCAN2-DT, AL161729.4, AC016394.1, and AC007128.2 transcript expression levels were significantly upregulated compared with the control group. In addition, the AL137782.1 and AC099850.3 transcript expression levels were significantly downregulated compared with the control.

The current study has some limitations. First, bioinformatics analyses rely heavily on data from public databases, such as TCGA. These datasets may have a limited sample size, population diversity, and data quality, which may affect the accuracy and generalizability of the results. Second, while in vitro experiments can simulate biological processes to a certain extent, they cannot fully replicate the complex in vivo environment. Therefore, there may be discrepancies between our in vitro experimental results and actual clinical situations. The results of this study are mainly based on bioinformatics analyses and in vitro experiments, lacking multi-center, large-sample clinical validation. Further validation is necessary to confirm the reliability and validity of these results.

Conclusion

We identified six cuproptosis-related lncRNAs: AL161729.4, ZKSCAN2-DT, AC016394.1, AC007128.2, AL137782.1, and AC099850.3. These lncRNAs were used to construct a prognostic model for CRC patients. The accuracy of this model was validated and a dependable nomogram was developed to predict patient outcomes. This study offers recommendations for personalized treatment of CRC patients and provides valuable information for gene therapy and immunotherapy development.