Construction of a new immune-related signature based on three lncRNAs as the factor for prognosis prediction of hepatocellular carcinoma.

Dear Editor,

Hepatocellular carcinoma (HCC) ranks the second place with regard to the cancer‐related death in the world, and there are approximately 800 000 new HCC cases reported every year. However, effective treatment is lacking at present. Recently, the treatment strategies to offset the mechanism of immunosuppression may possibly alter the HCC clinical outcomes. This reveals that more studies should be conducted to examine the association between abnormality in local immune status and the occurrence and progression of HCC. In addition, the long noncoding RNAs (lncRNAs) have exhibited diverse effects, whereas dysregulation in them, including abnormality in methylation modification, can result in HCC cell invasion, migration, proliferation, as well as epithelial‐mesenchymal transition. This letter was written aiming to shed more lights on the possible clinical significance of our immune‐related signature established using three lncRNAs in the stratification of patient prognosis, and whether they might be used as the biomarkers in the targeted therapy for HCC.

Altogether 331 immune genes and 14 142 lncRNAs were acquired based on the Molecular Signatures Database (Immune system process, M13664 and Immune response, M19817; http://www.broadinstitute.org/gsea/msigdb/index.jsp) and The Cancer Genome Atlas (TCGA) (https://cancergenome.nih.gov) for constructing our immune‐related signature. All TCGA‐derived tumor samples were randomly divided into training or test set, whereas the test and entire sets were adopted in signature verification. The relationship between the expression of immune‐related lncRNAs and the overall survival (OS) of patients was determined by univariate Cox regression. In addition, for preventing overfitting and deleting genes with close relation, the least absolute shrinkage and selection operator (LASSO) Cox regression was applied. Moreover, key genes were selected from those identified through univariate Cox regression. Then, for the genes of interest, their contributions to the prediction of patient prognosis were assessed by multivariate Cox regression.

Data in training set were analyzed and optimized, and a signature was established based on three new immune‐related lncRNAs, including AL365203.2, AC015908.3, and AC068987.4, and it was an accurate and convenient approach for predicting the prognosis for HCC patients (Table 1). The expression level was determined by the formula below, [expression quantity of AC015908.3 × (–0.4054)] + [expression quantity of AC068987.4 × (0.1099)] + [expression quantity of AL365203.2 × (0.0875)]. Figure 1A compares the differences in patient survival between both training set groups. In addition, the obtained results were additionally validated using the test and the entire sets, respectively (Figures 1B and 1C). The values of area under the curve (AUC) of training, test, and entire sets were 0.788, 0.739, and 0.796, respectively, indicating that our lncRNA‐based immune‐related prognostic model showed stable and moderate ability in monitoring patient survival (Figure 1D‐F). In addition, our as‐constructed signature showed the best AUC value relative to other conventional clinicopathological characteristics among TCGA cases, which also reflected that it had superb predictive ability. Figure 1G‐I presents risk score distribution, lncRNAs expression, along with survival status of training, test, and entire sets, separately. Gene set enrichment analysis was adopted for functional annotation, which indicated that those differentially expressed genes (DEGs) were upregulated within immune‐associated gene sets (Figure S1). Additionally, this study also analyzed the association between the as‐constructed lncRNA‐based signature and the infiltrating level of immunocytes, which revealed that the levels of CD8+T cells (Correlation coefficient = .249; P = 2.978 × 10–6), CD4+T cells (Correlation coefficient = .112; P = .039), B cells (Correlation coefficient = .143; P = .008), dendritic cells (Correlation coefficient = .274; P = 2.507 × 10–7), and neutrophils (Correlation coefficient = .312; P = 3.702 × 10–9) together with macrophages (Correlation coefficient = .470; P = 2.802 × 10–20) evidently increased in tumor microenvironment of high‐risk patients (P < .05). This observation suggested that there was difference in the immune status between both groups (Figure S2). The ESTIMATE method was used to obtain the ESTIMATE, stromal, and immune scores of high‐risk patients, which markedly increased compared with low‐risk patients (P < .05) (Figure 2A‐C). Nevertheless, difference in tumor purity was not statistically significant between both groups (Figure 2D). Noteworthily, expression of several human leukocyte antigen (HLA) genes remarkably increased among high‐risk patients relative to low‐risk counterparts (P < .05) (Figure 2E). We further analyzed the common T‐cell exhaustion‐associated genes and immune checkpoints in both groups to examine their expression. Remarkably, the results indicated that the PD‐1, CTLA‐4, TIGIT, PD‐L1, LAG3, and TIM‐3 contents of high‐risk patients significantly increased relative to those in low‐risk patients (P < .05) (Figure 2A‐F).

TABLE 1

lncRNAs in the risk assessment model

lncRNA	Coefficient	HR	Lower. 95	Upper. 95	P‐value
AC015908.3	–0.4054	0.67	0.5	0.89	.006
AC068987.4	0.1099	1.12	1.0	1.22	.019
AL365203.2	0.0875	1.09	1.0	1.18	.035

Abbreviations: HR, hazard ratio; lncRNA, long noncoding RNA.

FIGURE 1

Survival curves of patients in high‐risk group and low‐risk group of training set (A), testing set (B), and entire set (C). Patients in high‐risk group suffered shorter overall survival. D‐F, Survival‐dependent receiver operating characteristic (ROC) curves validation at 1 year of prognostic value of the prognostic index in the three sets (the training set, the testing set, and the entire set). Distribution of risk score, overall survival (OS), and gene expression in (A) training set, (B) testing set, and (C) entire set were also exhibited. Distribution of risk score, OS, and heat map of the expression of three signature lncRNAs in low‐risk and high‐risk groups are listed in the picture from top to bottom

FIGURE 2

Analysis of different immune status in high‐ and low‐risk groups of The Cancer Genome Atlas (TCGA) hepatocellular carcinoma (HCC) cohort. Comparison of (A) immune score, (B) stromal score, (C) ESTIMATE score, and (D) tumor purity between high‐ and low‐risk groups is shown. E, Comparison of the expression levels of human leukocyte antigen (HLA) genes between high‐ and low‐risk groups. F‐K, Box plots visualizing significantly different immune checkpoints between high‐ and low‐risk cases. Abbreviations: CTLA‐4, cytotoxic T‐lymphocyte associated protein 4; PD‐1 (PDCD1), programmed cell death 1; PD‐L1 (CD274), programmed death ligand 1; LAG3, lymphocyte activation gene‐3; TIGIT, T cell immunoreceptor with Ig and ITIM domains; TIM‐3 (HAVCR2), T‐cell immunoglobulin mucin receptor 3

In conclusion, the present work selects three lncRNAs (namely, AL365203.2, AC068987.4, and AC015908.3), which can be used to be the new lncRNA biomarkers for predicting the prognosis of HCC. In addition, difference in the immune status is statistically significant between both groups classified based on the median value of risk score. Nonetheless, more studies are warranted to examine the abovementioned lncRNAs together with related genes, so as to explore the potential mechanism underlying HCC occurrence.

CONFLICT OF INTEREST

The authors declare no conflict of interest.

AUTHOR CONTRIBUTIONS

Xin‐Ting Sang and Bo Hu created the idea for the paper. Bo Hu performed the collection and assembly of data, conducted the analysis, drafted the manuscript, and prepared the figures. Xiao‐Bo Yang and Xin‐Ting Sang revised the manuscript. All authors read and approved the final manuscript.

Supplementary Figure 1. Enrichment plots from gene set enrichment analysis (GSEA). GSEA results showing (A) GSE1460_intrathymic T progenitor vs naïve CD4+Tcell adult blood_UP, (B) GSE14308_Th2 vs Th17 UP, (C) GSE17974_CTRL vs act IL4 and anti IL12 12H CD4+Tcell_DN, (D) GSE20727_CTRL vs DNFB allergen treated DC_UP, (E) GSE24634_Treg vs Tconv post day5 IL4 conversion_UP, (F) GSE39110_untreated vs IL2 treated CD8+Tcell day3 post immunization_DN are positively differentially enriched in high risk group that screened out by the lncRNA‐related signature. (G) summarizes the above six gene sets.

Click here for additional data file.

Supplementary Figure 2. Relationships between the immune‐related prognostic index and infiltration abundances of six types of immune cells. The correlation was performed by using Pearson correlation analysis. (A) B cells; (B) CD4+Tcells; (C) CD8+T cells; (D) neutrophils; (E) macrophages; and (F) dendritic cells.

Click here for additional data file.