Long non-coding RNA LINC02613 is a prognostic biomarker for breast cancer and correlates with the cell cycle and immune infiltration based on TCGA data
Original Article

Long non-coding RNA LINC02613 is a prognostic biomarker for breast cancer and correlates with the cell cycle and immune infiltration based on TCGA data

Mengyao Cui1#, Yanbiao Liu1#, Lei Cui2

1Department of Breast Surgery, the 1st Affiliated Hospital, China Medical University, Shenyang, China; 2College of Health Management, China Medical University, Shenyang, China

Contributions: (I) Conception and design: M Cui, L Cui; (II) Administrative support: L Cui; (III) Provision of study materials or patients: M Cui, Y Liu; (IV) Collection and assembly of data: M Cui, Y Liu; (V) Data analysis and interpretation: M Cui, Y Liu; (VI) Manuscript writing: All authors; (VII) Final approval of manuscript: All authors.

#These authors contributed equally to this work.

Correspondence to: Cui Lei, MD. College of Health Management, China Medical University, 77 Puhe Road, Shenbei Xin District, Shenyang, China. Email: lcui@cmu.edu.cn.

Background: The development of transcriptomics has provided a new perspective to understand cancers. However, the role of long non-coding RNAs (lncRNAs) has not been fully elucidated. In the field of breast cancer, exploring the biological role and impact of lncRNAs on patient prognosis is of great importance.

Methods: Expression data of the lncRNA, LINC02613, were downloaded from The Cancer Genome Atlas (TCGA) dataset. The Kruskal-Wallis rank sum test, Wilcoxon rank sum test and logistic regression analysis were used to investigate the relationship between LINC02613 and clinicopathological characteristics. Cox regression and Kaplan-Meier analyses were used to assess the prognostic value of LINC02613. Gene ontology (GO) analysis, gene set enrichment analysis (GSEA) and immune infiltration analysis were used to explore the potential biological mechanisms of LINC02613.

Results: The expression of LINC02613 was downregulated in tumor samples and was related to poor overall survival (OS) in breast cancer patients. Statistical analysis revealed that low expression of LINC02613 was associated with Asian ethnicity, older age, higher histological grade and higher expression of estrogen receptor (ER), progesterone receptor (PR) and human epidermal growth factor receptor 2 (HER2). Positive expression of immunohistochemical markers (ER, PR and HER2), high pathological stage and N stage were independent predictors for low LINC02613 expression. Decreased expression of LINC02613 caused an increase in Wnt2 and Wnt3 protein expression, which activated the Wnt signaling pathway, thereby promoting cell cycle progression, DNA replication and glycolysis. At the same time, the antitumor immunity of the body was also weakened when the expression of LINC02613 was decreased.

Conclusions: Low expression of LINC02613 predicts poor OS in breast cancer patients. LINC02613 is involved in the regulation of the cell cycle, DNA replication and glycolysis via the Wnt signaling pathway, and it is related to antitumor immunity.

Keywords: LINC02613; breast cancer; cell cycle; immune infiltration and signaling pathway


Submitted Nov 07, 2021. Accepted for publication Feb 20, 2022.

doi: 10.21037/tcr-21-2479


Introduction

Breast cancer is the most diagnosed malignant tumor and one of the leading causes of cancer-related mortality among females worldwide. In 2020, the number of newly diagnosed breast cancer cases was 2.3 million, representing 11.7% of all cancer cases. With 690,000 deaths caused by breast cancer, it is the fifth leading cause of cancer-related mortality (1). From 2017 to 2019, China was one of the top 3 countries in the world with the highest breast cancer incidence and deaths (2,3).

In clinical settings, immunohistochemical markers [estrogen receptor (ER), progesterone receptor (PR), human epidermal growth factor receptor 2 (HER2) and Ki-67] are used to classify patients into 4 different categories (luminal A, luminal B, triple negative breast cancer (TNBC) and Her2-enriched), which serve as a basis for decision making for surgery and adjuvant therapies. Although this strategy has improved patient prognosis to some extent in recent decades, it is still far from benefiting every patient (4). According to a study in 2019, the global mortality rate of breast cancer continued to show a slow upward trend from 1990 to 2015 with recurrence and metastasis remaining the leading cause of death among patients (5). Advancements in second-generation sequencing technology and bioinformatics technology have allowed continual exploration of the biological behavior of tumor cells and their relationship with immune cells in the tumor microenvironment (TME) as well as searches for additional potential regulatory targets to achieve more precise and individualized medical care.

Long intergenic non-coding RNAs (lncRNAs) are defined as RNA transcripts composed of more than 200 nucleotides. As nomenclature has evolved, they are simply labeled lncRNAs (6). The omitted term, ‘intergenic’, indicates their origin as follows: transcribed from regions of the genome that do not encode any proteins (7,8). LncRNAs were once considered useless because they lack open-reading frames and their expression is much lower than that of protein-encoding genes (9). However, the importance of lncRNAs is gradually being recognized as they play important roles in a series of biological processes, including gene expression control, scaffold formation and epigenetic control. LncRNAs can encode some short polypeptides (10,11). In addition, their dysregulation has been shown to be related to tumorigenesis (12).

LINC02613 is transcribed from the 2p22.1 region of the human genome. In this study, we downloaded the expression data of LINC02613 from The Cancer Genome Atlas (TCGA) dataset and explored its impact on the prognosis of breast cancer patients. On the basis of differentially expressed genes (DEGs), we demonstrated that LINC02613 affected the Wnt signaling pathway by influencing the expression of Wnt2 and Wnt3, which in turn regulated the cell cycle and DNA replication as well as promoted glycolysis, and we also demonstrated that LINC02613 was associated with antitumor immunity. We present the following article in accordance with the REMARK reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-21-2479/rc).


Methods

Acquisition and preprocessing of RNA-sequencing (RNA-seq) data

We downloaded RNA-seq data in HTSeq-fragments per kilobase million (FPKM) format from UCSC XENA (https://xenabrowser.net/datapages/). Data without clinical information were omitted. Data in FPKM format were transformed then into TPM format to perform the analysis. Our study was in accordance with the publication guidelines provided by TCGA (https://cancergenome.nih.gov/publications/publicationguidelines) and did not contain any experiments involving human participants or animals performed by any of the authors. The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).

DEGs and functional enrichment analysis

All tumor samples were divided into two groups (high expression group and low expression group) according to the median expression level of LINC02613. The high expression group was labeled the test group, and the low expression group was labeled the control group. The DESeq2 package was used to compare expression data between the two groups (13). |Log fold change (logFC)| >2 and adjusted P<0.05 were regarded as the thresholds for DEGs. Functional enrichment analysis was performed using the Metascape database and its online analysis tools (http://metascape.org) (14).

Gene set enrichment analysis (GSEA)

GSEA is a computational method for verifying the extent to which a previously defined gene set is consistent with the expressed gene set. In the present study, signaling pathways with significant differences between the high and low LINC02613 expression groups were identified. The c2.all.v7.0.symbols.gmt curated from the Molecular Signatures Database (MSigDB) collection was used as the reference gene set. Permutations of gene sets were performed 1,000 times for each analysis. A false discovery rate (FDR) <0.25 and adjusted P<0.05 were regarded as statistically significant enrichments.

Immune infiltration analysis

Marker genes of immune cells, including dendritic cells (DCs), T cells, T follicular helper cells, cytotoxic cells, effector memory T cells, CD8+ T cells, Th1 cells, plasmacytoid DCs (pDCs), macrophages, activated DCs (aDCs), neutrophils, immature DCs (iDCs), NK CD56dim cells, NK cells, central memory T cells, gamma delta T cells, Tregs, mast cells, T helper cells, NK CD56bright cells, Th17 cells, eosinophils, Th2 cells and B cells, were obtained from a previously published paper (15). The single-sample Gene Set Enrichment Analysis (ssGSEA) method in the GSVA package was used to perform immune infiltration analysis by checking the expression level of these genes in published gene lists. Spearman correlation was then used to investigate the relationship between LINC02613 and immune cells.

Statistical analysis

The Wilcoxon rank sum test and Wilcoxon signed-rank test were used to compare LINC02613 expression between tumor samples and paratumor samples. The Kruskal-Wallis rank sum test and Wilcoxon rank sum test were used to analyze the relationship between multiple clinical features and the expression of LINC02613. Survival curves were constructed by the Kaplan-Meier method to assess the prognostic value of LINC02613 in breast cancer patients with overall survival (OS) (period from the date of randomization to the date of death or the last follow-up) as the primary endpoint as well as disease-specific survival (DSS) (period from the date of randomization to the date of death caused by breast cancer) and progression-free survival (PFS) (period from the date of randomization to the date of tumor progression or death from any cause) as the secondary end points. Cox regression models were used to determine independent prognostic factors. Variables with P<0.1 in univariate analysis were included in subsequent multivariate analysis. Hazard ratios (HRs) with 95% confidence intervals (CIs) were measured in all tests. All hypothesis tests were two-tailed, and a P<0.05 was considered statistically significant. Receiver operating characteristic (ROC) analysis was used to evaluate the discrimination ability of LINC02613 (16). All statistical analyses were performed in R (v.3.6.2).

Construction and assessment of nomogram

Based on the independent prognostic factors established by Cox regression analysis, a nomogram was constructed to show the predicted survival probabilities for 1-, 3- and 5-year OS using the rms package (https://cran.rproject.org/web/packages/rms/index.html). A calibration plot was then constructed to assess the predictive efficiency of the nomogram.


Results

Expression of LINC02613 is downregulated in tumors

After preliminary analysis, LINC02613 was found to be expressed at low levels in 24 tumor subtypes, including bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma (CESC), endocervical adenocarcinoma, colon adenocarcinoma (COAD), diffuse large B-cell lymphoma (DLBC), esophageal carcinoma (ESCA), kidney chromophobe (KICH), kidney renal clear cell carcinoma (KIRC), acute myeloid leukemia (LAML), low-grade glioma (LGG) of the brain, liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), ovarian serous cystadenocarcinoma (OV), pancreatic adenocarcinoma (PAAD), pheochromocytoma and paraganglioma (PCPG), prostate adenocarcinoma (PRAD), rectum adenocarcinoma (READ), skin cutaneous melanoma (SKCM), stomach adenocarcinoma (STAD), testicular germ cell tumors (TGCT), thyroid carcinoma (THCA), thymoma (THYM), uterine corpus endometrial carcinoma (UCEC) and uterine carcinosarcoma (UCS) (Figure 1). We then explored the expression of LINC02613 in breast cancers. The Wilcoxon rank sum test and Wilcoxon signed rank test confirmed that LINC02613 expression was downregulated in breast cancer tissues compared to normal tissues and adjacent noncancerous tissues (Figure 2). In addition, ROC analysis was performed to measure the discriminative power of LINC02613 in breast cancer. The area under the curve (AUC) was 0.925 (95% CI: 0.908–0.941), representing an efficient power of LINC012613 to discriminate tumors from normal tissues (Figure 3).

Figure 1 Pancancer analysis of LINC02613. ns, P≥0.05; *, P<0.05; **, P<0.01; ***, P<0.001. TCGA, The Cancer Genome Atlas.
Figure 2 LINC02613 is downregulated in breast cancer tissues compared to normal and adjacent noncancerous tissues. (A) The expression of LINC02613 was significantly lower in breast cancer tissues than in normal tissues according to the Wilcoxon rank sum test (P<0.001); (B) LINC02613 showed significantly lower expression in breast cancer tissues than in adjacent noncancerous tissues according to the Wilcoxon signed-rank test (P<0.001).
Figure 3 Receiver operating characteristic (ROC) curve showing the discrimination power of LINC02613. The AUC was 0.925, suggesting an efficient power to discriminate tumors from normal tissues. AUC, area under the curve; TPR, true positive rate; FPR, false positive rate.

Relationship between LINC02613 and clinicopathological characteristics

Clinicopathological characteristics were compared between patients in the high expression and low expression groups (Table 1). A total of 1,065 patients were enrolled in our study, and 976 of whom had racial information, including 60 Asian patients, 179 African American patients and 737 Caucasian patients. Statistical results revealed that a greater proportion of Asian patients expressed lower LINC02613 than others (P=0.006). The mean age of the low expression group was 61 years (ranging from 50 to 70 years), while the mean age of the high expression group was 55 years (ranging from 47 to 64 years) with a statistically significant difference (P<0.001). Among all patients, those with left-sided neoplasms accounted for 51.9%, and those with right-sided neoplasms accounted for 48.1%. For pathological stage, lower expression of LINC02613 was related to higher pathological stage (stages III and IV vs. stages I and II, P=0.016). In addition, there were more patients with low expression of LINC02613 in the T4 stage and N2 stage (71.4% for T4, P=0.04; and 63.8% for N2, P=0.01). Low expression of LINC02613 was related to patients with positive expression of ER, PR and HER2 (ER positive vs. ER negative: 54.9% vs. 31.6%, P<0.001; PR positive vs. PR negative: 58.2% vs. 38.8%, P<0.001; HER2 positive vs. HER2 negative: 66.2% vs. 44.5%, P<0.001). The above results are shown in Figure S1. Further logistic analysis revealed that positive expression of ER [odds ratio (OR) =0.38; 95% CI: 0.28–0.52; P<0.001], PR (OR =0.52; 95% CI: 0.40–0.68; P<0.001) and HER2 (OR =0.41; 95% CI: 0.28–0.59; P<0.001), higher pathological stage (OR =0.71; 95% CI: 0.53–0.94; P=0.016) and N stage (OR =0.78; 95% CI: 0.61–1.00; P=0.047) were independent predictors for low expression of LINC02613 (Table 2).

Table 1

Characteristics of patients

Characteristics Level Low expression group (n=533) High expression group (n=532) P
T stage, n (%) T1 130 (47.3) 145 (52.7) 0.04
T2 313 (50.9) 302 (49.1)
T3 63 (46.0) 74 (54.0)
T4 25 (71.4) 10 (28.6)
N stage, n (%) N0 236 (46.5) 271 (53.5) 0.01
N1 175 (50.9) 174 (49.1)
N2 74 (63.8) 42 (36.2)
N3 35 (47.3) 39 (52.7)
M stage, n (%) M0 449 (50.5) 440 (49.5) 1.000
M1 10 (50.0) 10 (50.0)
Pathological stage, n (%) Stage I & II 374 (47.6) 412 (52.4) 0.016
Stage III & IV 144 (56.3) 112 (43.7)
PR status, n (%) Negative 131 (38.8) 207 (61.2) <0.001
Positive 370 (58.2) 304 (41.8)
ER status, n (%) Negative 75 (31.6) 162 (68.4) <0.001
Positive 427 (54.9) 351 (45.1)
HER2 status, n (%) Negative 244 (44.5) 304 (55.5) <0.001
Positive 104 (66.2) 53 (43.8)
Race, n (%) Asian 41 (68.3) 19 (31.7.) 0.006
Black or African American 82 (17.4) 97 (19.2)
White 349 (73.9) 388 (77.0)
Anatomic neoplasm subdivisions, n (%) Left 275 (49.7) 278 (50.3) 0.877
Right 258 (50.4) 254 (49.6)
Age, median (IQR), years 61.00 (50.00, 70.00) 55.50 (47.00, 64.00) <0.001

P<0.05 indicates statistical significance. ER, estrogen receptor; PR, progesterone hormone receptor; HER2, human epidermal growth factor receptor 2; T, topography distribution; N, lymph node metastasis; M, distant metastasis; IQR, interquartile range.

Table 2

Results of logistic regression analysis

Characteristics OR in LINC02613 expression OR P value
T stage (T3 & T4 vs. T1 & T2) 1,062 0.95 (0.68–1.31) 0.739
N stage (non-N0 vs. N0) 1,046 0.78 (0.61–1.00) 0.047
M stage (M1 vs. M0) 909 1.02 (0.41–2.51) 0.964
Pathological stage (stage III & stage IV vs. stage I & stage II) 1,042 0.71 (0.53–0.94) 0.016
PR status (positive vs. negative) 1,012 0.52 (0.40–0.68) <0.001
ER status (positive vs. negative) 1,015 0.38 (0.28–0.52) <0.001
HER2 status (positive vs. negative) 705 0.41 (0.28–0.59) <0.001

P<0.05 indicates statistical significance. OR, odds ratio; ER, estrogen receptor; PR, progesterone hormone receptor; HER2, human epidermal growth factor receptor 2; T, topography distribution; N, lymph node metastasis; M, distant metastasis.

Low expression of LINC02613 indicates poor OS in breast cancer patients

Kaplan-Meier analysis revealed that low expression of LINC02613, although having no effect on DSS and PFS, predicted worse OS in patients [HR =0.65 (0.49–0.90), P=0.01] (Figure 4). Subgroup analysis of the TNM staging system revealed that low LINC02613 expression was significantly associated with poor OS in patients with stage T2 disease [HR =0.581 (0.364–0.928), P=0.023], stage N0 disease [HR =0.540 (0.297–0.982), P=0.043] and stage M0 disease [HR =0.660 (0.459–0.949), P=0.025] (Figure 5). T2 indicates that the maximum diameter of the primary tumor was >2.0 cm and <5.0 cm. N0 indicates that there was no infiltration of regional lymph nodes, and M0 indicates that no distant metastasis was discovered.

Figure 4 Kaplan-Meier analysis. (A) The Kaplan-Meier plot of relationship between LINC02613 and OS; (B) The Kaplan-Meier plot of relationship between LINC02613 and DSS; (C) The Kaplan-Meier plot of relationship between LINC02613 and PFS. OS, overall survival; DSS, disease-specific survival; PFS, progression-free survival.
Figure 5 Subgroup analysis revealing that low LINC02613 expression is significantly associated with poor OS of patients in the stage T2, stage N0 and stage M0 groups. If the HR of a group was infinite, the result was removed from the plot. OS, overall survival.

Prognostic factors identified by Cox regression analysis

Six independent prognostic factors were identified by univariate Cox regression, including T stage [HR =1.673 (1.152–2.429), P=0.007), N stage [HR =2.145 (1.497–3.073), P<0.001], M stage [HR =4.327 (2.508–7.465), P<0.001], pathological stage [HR =2.519 (1.787–3.549), P<0.001], age [HR =2.036 (1.468–2.822), P<0.001] and LINC02613 [HR =0.650 (0.469–0.901), P=0.010]. ER status (P=0.062) and HER2 status (P=0.059) were marginally significant. Further multivariate analysis confirmed that five of these factors, including T stage [HR =2.508 (1.293–4.864), P=0.007], M stage [HR =2.939 (1.141–7.545), P=0.026], ER status [HR =2.738 (1.597–4.694), P<0.001], age [HR =3.326 (2.001–5.501), P<0.001] and LINC02613 [HR =0.572 (0.332–0.984), P=0.544], were independent prognostic factors (Table 3).

Table 3

Results of Cox regression analysis

Characteristics Total (N) Univariate analysis HR (95% CI) P value Multivariate analysis HR (95% CI) P value
T stage (T3 & T4 vs. T1 & T2) 1,061 1.673 (1.152–2.429) 0.007 2.508 (1.293–4.864) 0.007
N stage (N1 & N2 & N3 vs. N0) 1,045 2.145 (1.497–3.073) <0.001 1.304 (0.699–2.435) 0.404
M stage (M1 vs. M0) 909 4.327 (2.508–7.465) <0.001 2.934 (1.141–7.545) 0.026
Pathological stage (stage III & stage IV vs. stage I & stage II) 1,041 2.519 (1.787–3.549) <0.001 1.953 (0.932–4.095) 0.076
PR status (negative vs. positive) 1,011 1.312 (0.931–1.849) 0.120
ER status (negative vs. positive) 1,014 1.420 (0.983–2.052) 0.062 2.738 (1.597–4.694) <0.001
HER2 status (negative vs. positive) 705 0.621 (0.378–1.019) 0.059 1.210 (0.657–2.227) 0.541
Age (>60 vs. ≤60 years) 1,064 2.036 (1.468–2.822) <0.001 3.326 (2.011–5.501) <0.001
Race (Caucasian vs. Asian and African American) 975 0.880 (0.593–1.306) 0.526
Anatomic neoplasm subdivisions (right vs. left) 1,064 0.776 (0.559–1.077) 0.130
LINC02613 (high vs. low) 1,064 0.650 (0.469–0.901) 0.010 0.572 (0.332–0.984) 0.044

P<0.05 indicates statistical significance. Variables with P<0.1 in univariate analysis were included in multivariate analysis. ER, estrogen receptor; PR, progesterone hormone receptor; HER2, human epidermal growth factor receptor 2; T, topography distribution; N, lymph node metastasis; M, distant metastasis; HR, Hazard ratio.

Construction and evaluation of the nomogram

A nomogram was constructed to provide surgeons with a quantitative approach to predict the prognosis of breast cancer patients (Figure 6). In this nomogram, each variable corresponds to a scaled line segment, in which length represents the rescaled range and impact on prognosis. The individual score of each variable and their sum can be easily obtained, and the probability of survival at different time nodes can be obtained. A calibration plot was then constructed to intuitively evaluate the prediction performance of the nomogram (Figure 7). The abscissa represents the probability of prognosis obtained from the nomogram, and the ordinate represents the actual prognosis (0 to 1 indicates that the probability of events is 0% to 100%). The colored lines are fitted lines, and the diagonal represents the ideal situation. Increased overlapping indicated better prediction performance. The C-index of the nomogram was 0.713 (0.687–0.738), indicating a medium level of predictive accuracy.

Figure 6 Nomogram for predicting 1-, 3-, and 5-year OS for breast cancer patients. The top row shows the point value for each variable. Rows 2–6 indicate each variable included in this nomogram. Row 7 shows total points of all 5 variables. Row 9–11 shows corresponding survival probability.
Figure 7 Calibration plot of the nomogram for predicting 1-, 3-, and 5-year OS for breast cancer patients. OS, overall survival.

Identification of DEGs and gene ontology (GO) analysis

Because the importance of LINC02613 for breast cancer patients has been well established, we were interested in further exploring its functions in organisms. Differential expression analysis revealed that the expression of 382 genes met the threshold of DEGs, of which 126 were upregulated and 256 were downregulated. Ten of these genes were selected for heatmapping (Figure 8). GO analysis revealed that genes regulated by LINC02613 were mainly enriched in biological processes related to tissue development and energy metabolism, including cornification, regulation of neuroblast proliferation, pentose and glucuronate interconversions and galactose metabolism. The detailed results and a network of their interactions are provided in Figure S2.

Figure 8 Identification of DEGs. (A) Volcano plot of LINC02613-related DEGs; (B) heatmap of LINC02613-related DEGs. DEGs, differentially expressed genes.

GSEA

Traditional functional enrichment analysis only focuses on a few genes with significant differential expression, omitting the genes with great biological importance but inconspicuously different expression. Importantly, GSEA fills this gap. The GSEA results showed that six datasets were significantly enriched in the LINC02613 low expression group, including DNA replication, cell cycle, glycolysis gluconeogenesis, antigen processing and presentation, chemokine signaling pathway and WNT signaling pathway (Figure S3). Tumor cells are characterized by their unrestricted ability to divide and proliferate; DNA replication and smooth cell cycle are necessary processes for cell proliferation. In addition, adequate nutrition is a prerequisite, but the energy provided by normal intracellular aerobic oxidation is insufficient to support this abnormal process. According to Warburg's theory, in such a condition, tumor cells will generate energy by activating glycolysis in the presence of oxygen (17). Antigen processing and presentation as well as chemokine secretion are important components of antitumor immunity (18). Combined with previous studies, we hypothesized that low expression of LINC02613 regulates tumor cell cycle progression, DNA replication and glycolysis through the Wnt signaling pathway, and we also hypothesized that low expression of LINC02613 weakens antitumor immunity.

Relationship between LINC02613 and cell cycle

The cell cycle refers to the activity of a cell from the beginning of one mitosis to the end of the next mitosis. The cell cycle is divided into the following five phases: G1 phase (pre-DNA synthesis), S phase (DNA synthesis), G2 phase (post-DNA synthesis), M phase (mitosis) and G0 phase (pre-DNA synthesis). DNA replication occurs in S phase. When growth factors bind to cell surface receptors, cells in the G0 phase are activated and express cyclins and cyclin-dependent kinases (CDKs). Cyclin D1, the first discovered and most studied cyclin molecule, is a cell cycle initiator and growth-sensing factor. In the early G1 phase, cyclin D1 binds to CDK4 or CDK6 to drive cell cycle progression and DNA replication (19). In this study, we examined the relationship of CCND1 (the gene encoding cyclin D1), CDK4 and CDK6 with LINC02613. The results showed that the expression of CCND1 was negatively correlated with LINC02613 (Figure S4).

Relationship between LINC02613 and the Wnt signaling pathway

The Wnt pathway is an evolutionarily conserved signaling pathway, and it can be divided into β-catenin-dependent and β-catenin-independent subtypes. Aberrant activation of the β-catenin-dependent Wnt pathway has been shown to result in tumorigenesis and increased cellular glycolysis (20,21). When Wnt proteins bind to Frizzled receptors, repressed β-catenin is released, activating downstream TEF/LEF and affecting the transcription of numerous genes, including CCND1. There are many members of the Wnt protein family. Wnt2, Wnt3 and Wnt3a are all considered to be initiators of the β-catenin-dependent Wnt signaling pathway. In this study, we investigated the relationship between the Wnt2, Wnt3 and Wnt3a coding genes and LINC02613. The results revealed that the expression of Wnt2 and Wnt3 was negatively correlated with LINC02613 (Figure S5). Considering the relationship between LINC02613 and CCND1, the downregulation of LINC02613 may promote the expression of Wnt2 and Wnt3, leading to activation of the Wnt signaling pathway, which promotes the expression of downstream cyclin D1, leading to cell cycle progression and DNA replication of tumor cells.

Relationship between LINC02613 and immune infiltration

The balance between the tumoricidal effect of immune cells and the immune escape of tumor cells in the microenvironment plays an important role in the occurrence and metastasis of cancer. We analyzed the relationship between the expression of LINC02613 and infiltration of immune cells in the TME by ssGSEA. The results suggested that the expression of LINC02613 was positively related to the infiltration of 20 out of 24 immune cells, including CD8+ T cells, NK cells, macrophages and DCs, which remain the workhorses in eliminating tumors. In addition, other innate immune cells, such as mast cells and neutrophils, as well as functional T cells, including Tfh, Tem, Th1, Tcm, Tgd, Treg and T helper cells, were also shown to be highly expressed in the test group (Figure 9). These results suggested that the expression of LINC02613 is positively related to immune infiltration in the TME and that the downregulation of LINC02613 may weaken the antitumor effect in breast cancer patients.

Figure 9 Relationship between LINC02613 and immune infiltration. (A) Immune infiltration in test group; (B) LINC02613 was positively related to CD8+ T cells; (C) LINC02613 was positively related to macrophages; (D) LINC02613 was positively related to NK cells; (E) LINC02613 was positively related to DCs. DCs, dendritic cells; pDCs, plasmacytoid DCs; aDCs, activated DCs; iDCs, immature DCs

Discussion

Treatment for breast cancer is entering the era of individual-based precision medicine. LncRNAs have great potential to stratify the risk of recurrence and predict the prognosis of patients (22).

Dysregulation of lncRNA expression has been shown to be related to the occurrence and metastasis of tumors. HOTAIR has been shown to be highly expressed in breast cancer tissues and to promote tumor cell metastasis by reprogramming the chromatin state (23). LUNAR has been shown to regulate T-acute lymphoblastic leukemia (T-ALL) through the Notch signaling pathway (24). Another non-coding RNA, MALAT1, has been shown to be associated with tumor cell metastasis in a variety of cancer types, including lung, breast and colon cancers (25). In the present study, we identified a potential breast cancer suppressor gene, LINC02613.

Analysis of TCGA database revealed that LINC02613 expression was significantly lower in tumor samples compared to normal samples. Low expression of LINC02613 was related to ethnicity (Asian vs. non-Asian, P<0.001), mean age (61.0 vs. 55.5 years, P<0.001), pathological stage (stages I and IV vs. stages I and II, P=0.012), T and N stage (T4 vs. non-T4, P=0.012; non-N0 vs. N0, P=0.038) and expression of immunohistochemical markers (ER, PR and HER2). Logistic regression analysis revealed that positive expression of immunohistochemical markers (ER, PR and HER2), higher pathological stage (stages III and IV vs. stages I and II, P=0.016) and N stage (non-N0 vs. N0, P=0.047) were independent predictors for low expression of LINC02613. These results indicated that LINC02613 may play a tumor suppressive role and that low expression of LINC02613 implies a more aggressive and metastatic ability of breast cancer cells.

Based on TCGA database, Kaplan-Meier survival analysis and Cox regression analysis indicated that low expression of LINC02613 was associated with worse OS in breast cancer patients. We identified the following five independent prognostic factors: age, T stage, LINC02613, metastasis and ER status. These results indicated that LINC02613 may serve as a candidate marker to assist survival prediction in clinical settings. Therefore, a nomogram including these prognostic factors was subsequently built and showed good predictive value.

Because the importance of LINC02613 for breast cancer patients has been well established, we further explored its functional mechanisms. In conventional GO analysis, LINC02613 was mainly associated with energy metabolism and tissue development. In GSEA, six gene sets were significantly enriched in the LINC02613 low expression group, including DNA replication, cell cycle, glycolysis gluconeogenesis, Wnt signaling pathway, antigen presentation and chemokine signaling pathway. Further gene-to-gene correlation analysis showed that LINC02613 expression was inversely correlated with Wnt2, Wnt3 and CCND1.

The eukaryotic cell cycle is a highly organized and strictly regulated process that ensures the correct distribution of genetic material in the nucleus during cell division (26). Cyclin D1, as one of the key regulators of the cell cycle, binds to CDK4 or CDK6, enters the nucleus and phosphorylates downstream molecules to promote the transcription of genes required for cell proliferation (27). CCND1, as the encoding gene of cyclin D1, has been shown to be overexpressed in tumor cells and is regarded as a negative prognostic biomarker (28). Increasing evidence has revealed that CCND1 is regulated by the Wnt signaling pathway (29,30). When the Wnt pathway is activated, the downstream TEF/LEF complex promotes the transcription of CCND1, which induces cell cycle progression and DNA replication (31). In the present study, the expression of Wnt2, Wnt3 and CCND1 was negatively correlated with LINC02613, which implied that low expression of LINC02613 contributed to the unrestricted proliferation of tumor cells by upregulating Wnt protein, thus activating the Wnt pathway and promoting transcription of downstream CCND1.

Metabolic reprogramming is a prominent feature of tumors, and abnormalities in glucose metabolism were the first identified examples of metabolic reprogramming (32,33). Under aerobic conditions, tumor cells still rely on glycolysis for energy, which is also known as the Warburg effect (17). Recent studies have revealed that the canonical Wnt signaling pathway promotes tumor cell glycolysis in a variety of ways, including activating the expression of MCT-1, CYC1 and ATP synthase through the downstream transcription factor, TCF/LEF, which promotes the overexpression of c-Myc and regulates the activity of rate-limiting enzymes in the energy metabolism pathway (21). In the present study, both glycolysis gluconeogenesis and the Wnt signaling pathway were enriched in the low LINC02613 expression group, implying that low expression of LINC02613 may promote cellular glycolysis by activating the Wnt pathway.

In addition, the present study showed that LINC02613 expression was positively correlated with several types of antitumor immune cells, including CD8+ T cells, macrophages, DCs and Treg cells, which are recognized to play vital roles in antitumor immune responses (34). It is conceivable that the number of immune cells in the TME decreases dramatically when LINC02613 is expressed at low levels, further exacerbating tumor growth and metastasis.

In the present study, we first explored the relationship between LINC02613 and breast cancer. Despite the important findings in the present study, there were several limitations. First, although advances in sequencing and database technologies have greatly broadened information access compared to the past, the sample size in this study was still limited due to the single source of sample information. Second, due to objective constraints, we only explored the expression of LINC02613 on the bioinformatics platform, which is insufficient for rigorous scientific study. Thus, there is a clear need for additional multicenter and large-scale clinical trials before applying our nomogram in clinical practice. In addition, further cytological and in vivo experiments are required to validate our findings in terms of molecular mechanisms.


Conclusions

LINC02613 is a potential tumor suppressor gene. Low expression of LINC02613 predicts poor OS in breast cancer patients. LINC02613 not only affects the Wnt signaling pathway by influencing the expression of Wnt2 and Wnt3 proteins, which in turn regulates the cell cycle, regulates DNA replication and promotes glycolysis, but is also associated with antitumor immunity.


Acknowledgments

Funding: None.


Footnote

Reporting Checklist: The authors have completed the REMARK reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-21-2479/rc

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-21-2479/prf

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-21-2479/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.


References

  1. Sung H, Ferlay J, Siegel RL, et al. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin 2021;71:209-49. [Crossref] [PubMed]
  2. Yi M, Li T, Niu M, et al. Epidemiological trends of women's cancers from 1990 to 2019 at the global, regional, and national levels: a population-based study. Biomark Res 2021;9:55. [Crossref] [PubMed]
  3. Li N, Deng Y, Zhou L, et al. Global burden of breast cancer and attributable risk factors in 195 countries and territories, from 1990 to 2017: results from the Global Burden of Disease Study 2017. J Hematol Oncol 2019;12:140. [Crossref] [PubMed]
  4. Gannon LM, Cotter MB, Quinn CM. The classification of invasive carcinoma of the breast. Expert Rev Anticancer Ther 2013;13:941-54. [Crossref] [PubMed]
  5. Azamjah N, Soltan-Zadeh Y, Zayeri F. Global Trend of Breast Cancer Mortality Rate: A 25-Year Study. Asian Pac J Cancer Prev 2019;20:2015-20. [Crossref] [PubMed]
  6. St Laurent G, Wahlestedt C, Kapranov P. The Landscape of long noncoding RNA classification. Trends Genet 2015;31:239-51. [Crossref] [PubMed]
  7. Boon RA, Jaé N, Holdt L, et al. Long Noncoding RNAs: From Clinical Genetics to Therapeutic Targets? J Am Coll Cardiol 2016;67:1214-26. [Crossref] [PubMed]
  8. Deniz E, Erman B. Long noncoding RNA (lincRNA), a new paradigm in gene expression control. Funct Integr Genomics 2017;17:135-43. [Crossref] [PubMed]
  9. Cabili MN, Trapnell C, Goff L, et al. Integrative annotation of human large intergenic noncoding RNAs reveals global properties and specific subclasses. Genes Dev 2011;25:1915-27. [Crossref] [PubMed]
  10. Ruiz-Orera J, Messeguer X, Subirana JA, et al. Long non-coding RNAs as a source of new peptides. Elife 2014;3:e03523. [Crossref] [PubMed]
  11. Rinn JL, Chang HY. Genome regulation by long noncoding RNAs. Annu Rev Biochem 2012;81:145-66. [Crossref] [PubMed]
  12. Younger ST, Kenzelmann-Broz D, Jung H, et al. Integrative genomic analysis reveals widespread enhancer regulation by p53 in response to DNA damage. Nucleic Acids Res 2015;43:4447-62. [Crossref] [PubMed]
  13. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014;15:550. [Crossref] [PubMed]
  14. Zhou Y, Zhou B, Pache L, et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun 2019;10:1523. [Crossref] [PubMed]
  15. Bindea G, Mlecnik B, Tosolini M, et al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity 2013;39:782-95. [Crossref] [PubMed]
  16. Robin X, Turck N, Hainard A, et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinformatics 2011;12:77. [Crossref] [PubMed]
  17. WARBURG O. On the origin of cancer cells. Science 1956;123:309-14. [Crossref] [PubMed]
  18. Jhunjhunwala S, Hammer C, Delamarre L. Antigen presentation in cancer: insights into tumour immunogenicity and immune evasion. Nat Rev Cancer 2021;21:298-312. [Crossref] [PubMed]
  19. Musgrove EA, Caldon CE, Barraclough J, et al. Cyclin D as a therapeutic target in cancer. Nat Rev Cancer 2011;11:558-72. [Crossref] [PubMed]
  20. Bugter JM, Fenderico N, Maurice MM. Mutations and mechanisms of WNT pathway tumour suppressors in cancer. Nat Rev Cancer 2021;21:5-21. [Crossref] [PubMed]
  21. Mo Y, Wang Y, Zhang L, et al. The role of Wnt signaling pathway in tumor metabolic reprogramming. J Cancer 2019;10:3789-97. [Crossref] [PubMed]
  22. Malih S, Saidijam M, Malih N. A brief review on long noncoding RNAs: a new paradigm in breast cancer pathogenesis, diagnosis and therapy. Tumour Biol 2016;37:1479-85. [Crossref] [PubMed]
  23. Gupta RA, Shah N, Wang KC, et al. Long non-coding RNA HOTAIR reprograms chromatin state to promote cancer metastasis. Nature 2010;464:1071-6. [Crossref] [PubMed]
  24. Trimarchi T, Bilal E, Ntziachristos P, et al. Genome-wide mapping and characterization of Notch-regulated long noncoding RNAs in acute leukemia. Cell 2014;158:593-606. [Crossref] [PubMed]
  25. Zhang A, Xu M, Mo YY. Role of the lncRNA-p53 regulatory network in cancer. J Mol Cell Biol 2014;6:181-91. [Crossref] [PubMed]
  26. Otto T, Sicinski P. Cell cycle proteins as promising targets in cancer therapy. Nat Rev Cancer 2017;17:93-115. [Crossref] [PubMed]
  27. Montalto FI, De Amicis F. Cyclin D1 in Cancer: A Molecular Connection for Cell Cycle Control, Adhesion and Invasion in Tumor and Stroma. Cells 2020;9:2648. [Crossref] [PubMed]
  28. Zhang LQ, Jiang F, Xu L, et al. The role of cyclin D1 expression and patient's survival in non-small-cell lung cancer: a systematic review with meta-analysis. Clin Lung Cancer 2012;13:188-95. [Crossref] [PubMed]
  29. Vlad-Fiegen A, Langerak A, Eberth S, et al. The Wnt pathway destabilizes adherens junctions and promotes cell migration via β-catenin and its target gene cyclin D1. FEBS Open Bio 2012;2:26-31. [Crossref] [PubMed]
  30. Zhang J, Gill AJ, Issacs JD, et al. The Wnt/β-catenin pathway drives increased cyclin D1 levels in lymph node metastasis in papillary thyroid cancer. Hum Pathol 2012;43:1044-50. [Crossref] [PubMed]
  31. Stamatakos M, Palla V, Karaiskos I, et al. Cell cyclins: triggering elements of cancer or not? World J Surg Oncol 2010;8:111. [Crossref] [PubMed]
  32. Cantor JR, Sabatini DM. Cancer cell metabolism: one hallmark, many faces. Cancer Discov 2012;2:881-98. [Crossref] [PubMed]
  33. Granger A, Mott R, Emambokus N. Hacking Cancer Metabolism. Cell Metab 2016;24:643-4. [Crossref] [PubMed]
  34. Kaymak I, Williams KS, Cantor JR, et al. Immunometabolic Interplay in the Tumor Microenvironment. Cancer Cell 2021;39:28-37. [Crossref] [PubMed]
Cite this article as: Cui M, Liu Y, Cui L. Long non-coding RNA LINC02613 is a prognostic biomarker for breast cancer and correlates with the cell cycle and immune infiltration based on TCGA data. Transl Cancer Res 2022;11(4):615-628. doi: 10.21037/tcr-21-2479

Download Citation