Data acquisition and analysis

The transcriptomic data of 473 COAD tumors and 41 normal samples were downloaded from the TCGA database [30]. Then, we separated the expression of 14,056 lncRNAs and 19,573 mRNAs in COAD samples by Strawberry Perl. The clinical information of 421 COAD patients was also obtained from TCGA database after excluding samples with short-term survival (less than 30 days) or missing follow-up days (Table 1). By data merging, 417 COAD patients were finally included in the present analysis. We collected 18 cuproptosis-related mRNAs from previous literature [31,32,33,34,35,36], including FDX1, DLD, PDHA1, PDHB, MTF1, GLS, CDKN2A, DLAT, LIAS, LIPT1, LIPT2, ATP7A, ATP7B, SLC31A1, SLC31A2, DLST, NFE2L2, NLRP3 and extracted the expression of those CRGs from COAD samples accordingly.

Table 1 The clinical characteristics of COAD patients

Serum of 150 COAD patients and 135 healthy controls was collected. The diagnosis of COAD patients was confirmed by histopathology or biopsy and recruited from the Department of General Surgery, Qilu Hospital of Shandong University, from April 2018 to October 2020. The healthy controls were enrolled from the Department of Physical Examination Center, Qilu Hospital of Shandong University. Serum samples were separated by centrifugation at 6000 g for 10 min followed by another centrifugation at 12,000 g for 10 min and then stored at − 80 °C for further analysis.

Identification of differentially expressed cuproptosis-related LncRNAs

After the Pearson correlation algorithm with the filter of |coefficient|> 0.3 and P < 0.001, we selected lncRNAs that were related with cuproptosis genes and considered as CRLs. Next, we identified differentially expressed lncRNAs (|Log2 fold change (FC)|> 1, false discovery rate (FDR) < 0.05) in COAD tumor tissues comparing with normal tissue using differential analysis by R package “limma” [37].

Establishment and evaluation of CRLs prognostic signature

After the COAD patients were randomly divided into training and test groups, we first performed univariate Cox analysis (P < 0.01) to screen CRLs associated with prognosis. Next, we established a prognostic signature by multivariate Cox regression analysis in the training group. Hence, the risk score of each COAD patient could be calculated according to the coefficient and CRLs expression in our prognostic signature. And the COAD patients were divided into the high-risk and low-risk groups by the median value of the risk score.

We used the Kaplan–Meier (K–M) and receiver operating characteristic (ROC) curves to evaluate the value of the prognostic signature in both training and test groups by R packages “survival,” “survminer” and “timeROC” [38]. Besides, the relationship between the risk score and prognosis of COAD patients was also displayed by heatmap jointly, risk score curve, and survival status diagram. Moreover, whether the risk score was related with clinical parameters was also examined.

Independent prognostic analysis and development of nomogram model

The univariate and multivariate Cox regression analyses were performed to identify the independent risk factors of COAD, including risk score and clinical characteristics. Subsequently, the nomogram model was constructed based on independent risk factors using the R package “rms.” Then, we used calibration curves to estimate the prediction power of the model.

Gene enrichment analysis by GSEA

To identify pathway enrichment in two risk groups, we used GSEA software (4.2.2) to perform the enrichment analysis of the Kyoto Encyclopedia of Genes and Genomes (KEGG) (c2.cp.kegg.v7.5.1.symbols.gmt) [39]. The random sample permutation number was set as 1,000, and the significance threshold was P < 0.05.

Immune infiltration analysis by single sample gene set enrichment analysis (ssGSEA)

The enrichment score of infiltration estimation and immune function of different immune cells between two risk groups was compared using ssGSEA analysis [40] by R packages “GSVA,” “GSEABase” and “Limma.” So, we could easily explore the association between risk score, immune infiltration and immune function. And the significance threshold was FDR < 0.05.

The value of risk score in predicting response of patients to immunotherapy and chemotherapy

We also analyzed the differential expression of 47 immune checkpoint genes, including IDO1, LAG3, CTLA4, TNFRSF9, ICOS, CD80, PDCD1LG2, TIGIT, CD70, TNFSF9, ICOSLG, KIR3DL1, CD86, PDCD1 (PD1), LAIR1, TNFRSF8, TNFSF15, TNFRSF14, IDO2, CD276, CD40, TNFRSF4, TNFSF14, HHLA2, CD244, CD274 (PD-L1), HAVCR2, CD27, BTLA, LGALS9, TMIGD2, CD28, CD48, TNFRSF25, CD40LG, ADORA2A, VTCN1, CD160, CD44, TNFSF18, TNFRSF18, BTNL2, C10orf54, CD200R1, TNFSF4, CD200, NRP1 between two risk groups. Besides, we collected 15 commonly used drugs for the clinical treatment of gastrointestinal tumors, including Epothilone B, Sorafenib, Cisplatin, Doxorubicin, Etoposide, Imatinib, Lapatinib, OSI.906, PHA.665752, ABT.888, Camptothecin, Docetaxel, Mitomycin C, Paclitaxel, and Sunitinib. The half-maximal inhibitory concentration (IC50) of drugs was used to evaluate the therapy response of patients in two risk groups by R package “pRRophetic.” The significance threshold of all the above analyses results was P < 0.05 except for the multiple hypothesis test which used FDR to adjust.

RNA extraction and RT-qPCR

The total RNA was extracted from serum samples using TRIzol LS Reagent (Invitrogen, Eugene, OR, USA). The concentration of RNA was measured using a NanoDrop spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). RNA was reverse transcription into cDNA using SureScript RTase Mix and RT Reaction Buffer, and qPCR was performed using Blaze Taq qPCR Mix (GeneCopoeia, Guangzhou, China). The relative expression of target lncRNAs was normalized to the glyceraldehyde-3-phosphate dehydrogenase (GAPDH) and calculated 2 − ΔΔCt. The primer sequences are shown in Additional file 1: Table S1.