Screening key lncRNAs with diagnostic and prognostic value for head and neck squamous cell carcinoma based on machine learning and mRNA-lncRNA co-expression network analysis

Abstract

BACKGROUND:

Head and neck squamous cell carcinoma (HNSCC) is the seventh most common type of cancer around the world. The aim of this study was to seek the long non-coding RNAs (lncRNAs) acting as diagnostic and prognostic biomarker of HNSCC.

METHODS:

Base on TCGA dataset, the differentially expressed mRNAs (DEmRNAs) and lncRNAs (DElncRNAs) were identified between HNSCC and normal tissue. The machine learning and survival analysis were performed to estimate the potential diagnostic and prognostic value of lncRNAs for HNSCC. We also build the co-expression network and functional annotation. The expression of selected candidate mRNAs and lncRNAs were validated by Quantitative real time polymerase chain reaction (qRT-PCR).

RESULTS:

A total of 3363 DEmRNAs (1822 down-regulated and 1541 up-regulated mRNAs) and 32 DElncRNAs (13 down-regulated and 19 up-regulated lncRNAs) between HNSCC and normal tissue were obtained. A total of 13 lncRNAs (IL12A.AS1, RP11.159F24.6, RP11.863P13.3, LINC00941, FOXCUT, RNF144A.AS1, RP11.218E20.3, HCG22, HAGLROS, LINC01615, RP11.351J23.1, AC024592.9 and MIR9.3HG) were defined as optimal diagnostic lncRNAs biomarkers for HNSCC. The area under curve (AUC) of the support vector machine (SVM) model, decision tree model and random forests model and were 0.983, 0.842 and 0.983, and the specificity and sensitivity of the three model were 95.5% and 96.2%, 77.3% and 97.6% and 93.2% and 97.8%, respectively. Among them, AC024592.9, LINC00941, LINC01615 and MIR9-3HG was not only an optimal diagnostic lncRNAs biomarkers, but also related to survival time. The focal adhesion, ECM-receptor interaction, pathways in cancer and cytokine-cytokine receptor interaction were four significantly enriched pathways in DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs. But for most of the selected DEmRNAs and DElncRNAs, the expression was consistent with our integrated analysis results, including LINC00941, LINC01615, FOXCUT, TGA6 and MMP13.

CONCLUSION:

AC024592.9, LINC00941, LINC01615 and MIR9-3HG was not only an optimal diagnostic lncRNAs biomarkers, but also were a prognostic lncRNAs biomarkers.

Keywords

Head and neck squamous cell carcinoma diagnostic prognostic machine learning

1. Introduction

Head and neck squamous cell carcinoma (HNSCC) is one of the most malignant tumors, ranks among the seventh most common type of cancer worldwide [1]. In spite of considerable advances in surgery, radiotherapy, chemotherapy and multimodality therapy, prognosis of HNSCC remains poor, and the 5-year overall survival is $\sim$ 50% [2]. Thence, it is urgently needed to reveal molecular pathogenesis of HNSCC and seek novel diagnostic and prognostic biomarkers for HNSCC.

With advances in gene expression profiles, lncRNAs are reported to be as pivotal regulators involved in various biological processes in cancer, including tumor occurrence, development and metastasis [3, 4]. Recently, many lncRNAs are considered to be novel candidate biomarkers for the diagnosis and prediction of various cancers [5, 6, 7]. However, there are few study on lncRNA diagnostic and prognostic biomarkers for HNSCC is rare. Therefore, the identification of key tumor markers of HNSCC is required.

The aim of this study is to obtain biomarkers of lncRNAs for diagnostic and prognostic HNSCC patients, we used the TCGA datasets to obtain the differentially expressed mRNAs (DEmRNAs) and lncRNAs (DElncRNAs). We performed the machine learning and survival analysis to evaluate the potential diagnostic and prognostic value of lncRNAs for HNSCC. The DElncRNA-DEmRNA co-expression network and the functional annotation of the DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs in HNSCC was performed. The expression levels of candidate genes were verified by qRT-PCR. To our knowledge, this is first time to find diagnostic and prognostic lncRNAs biomarkers in HNSCC by using machine learning.

2. Materials and methods

2.1 TCGA data

The lncRNA and mRNA gene expression profiles and clinical data of HNSCC were download by the Cancer Genome Atlas (TCGA) (http://tcgadata.nci.nih. gov/). In this study, 500 HNSCC tissues and 44 normal adjacent samples from patients with HNSCC were included.

2.2 Identification of DEmRNAs and DElncRNAs

The undetectable lncRNAs and mRNAs (with read count value $=$ 0 in more than 20% HNSCC case) were filtered and deleted. The DElncRNAs and DEmRNAs in HNSCC compared to adjacent normal tissues were analyzed using the R package DESeq2. Benjamini and Hochberg multiple testing method was applied to acquire the false discovery rate (FDR). FDR $<$ 0.05 and $|$ Log2fold change $|$ $>$ 1 was set as the cut-off criteria to define DElncRNAs and DEmRNAs. Hierarchical clustering analysis of DElncRNAs and DEmRNAs were generated by using R package v3.3.3 (https://www.r- project.org/).

2.3 Identification of the optimal diagnostic lncRNAs biomarkers for HNSCC

To further identify the optimal diagnostic lncRNA biomarkers for HNSCC, we performed feature selection procedures as follows. (1) The importance value of each lncRNA was ranked according to mean decrease in accuracy from large to small by random forest algorithm. (2) The optimum number of features was found by adding a DElncRNA at a time in the top down forward-wrapper packaging method. (3) By using support vector machine (SVM) at each increment and the optimal diagnostic lncRNA biomarkers were identified for HNSCC.

The ‘random Forests’ packet (https://cran.r-project. org/web/packages/randomForest/) was used to establish the random forest model. The ‘rpart’ packet (https: //cran.r-project.org/web/packages/rpart/) was used to build the decision tree model. The e1071 package (https://cran.r-project.org/web/packages/e1071/index. html) in R was used to establish the SVM model. Diagnostic ability of these three models, and each lncRNA biomarker was evaluated by acquiring the area under the receiver operating characteristic (ROC) curve (AUC), sensitivity and specificity. By using pROC package in R language we performed the ROC analyses to assess the diagnostic value of lncRNA biomarker. The AUC under binomial exact confidence interval was calculated and ROC curve was produced.

2.4 Survival analysis of optimal diagnostic lncRNAs biomarkers for HNSCC

To determine the potential association between the identified DElncRNAs and survival in HNSCC patients, survival analysis (https://cran.r-project.org/web/ packages/survival/index.html) in R was performed.

2.5 DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

The correlation between the optimal diagnostic lncRNAs and DEmRNAs were analyzed by the Pearson correlation coefficient. The threshold for DElncRNA-DEmRNA co-expression pairs was $P<$ 0.05 and $R>$ 0.5. We used the cytoscape (http://www. cytoscape.org/) to build the DElncRNA-DEmRNA co-expression network.

2.6 Functional annotation of DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

Gene Ontology (GO) classification and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis were performed using GeneCoDis (http://genecodis.cnb.csic.es/analysis). FDR $<$ 0.05 was considered to be significantly enriched terms.

Table 1
Primer sequences used for qRT-PCR

Name	Sequence (5’ to 3’)
GAPDH-F	GGAGCGAGATCCCTCCAAAAT
GAPDH-R	GGCTGTTGTCATACTTCTCATGG
ACTB-F	CATGTACGTTGCTATCCAGGC
ACTB-R	CTCCTTAATGTCACGCACGAT
LINC00941-F	CATGCACCACTACACTCAGC
LINC00941-R	GCCTCCTTGCTGACTGATTC
LINC01615-F	CGGAGCAGCAGGAAATCTCT
LINC01615-R	GCAGGTTTTGTTGCGAAGTG
FOXCUT-F	CAGATGTTCGGTCCCCACTC
FOXCUT-R	TCGGCCTCCTTTGCTTAGAC
ITGA6-F	GGCGGTGTTATGTCCTGAGTC
ITGA6-R	AATCGCCCATCACAAAAGCTC
MMP13-F	TACACCTACACCGGCAAAAGC
MMP13-R	AAAACAGCTCCGCATCAACCT
FOXC1-F	TGTTCGAGTCACAGAGGATCG
FOXC1-R	ACAGTCGTAGACGAAAGCTCC

2.7 Confirmation by qRT-PCR

Base on the results of TCGA integration analysis, ITGA6, MMP13 and FOXC1, LINC00941, LINC01615 and FOXCUT were screened as candidate mRNAs and lncRNAs. A total of 6 HNSCC patients were enrolled in this study. Twelve tissues samples of HNSCC patients ( $n=$ 6) and normal adjacent ( $n=$ 6) were obtained from surgery. Tissue samples were frozen in liquid nitrogen and stored at $-$ 80 ${}^{\circ}$ C before RNA isolation. All the participants submitted the signed informed consents and the protocol were approved by the ethical committee of our Hospital. The study protocol conformed to the ethical guidelines of the Helsinki of Declaration.

Total RNA was extracted from samples with TRIzol reagent (Tiangen, China) and then subjected to reverse transcription and PCR reactions using Super Real PreMix Plus SYBR Green (Tiangen, China). The 2- $\Delta\Delta$ Ct method was used to analyze the relative quantification of mRNA and lncRNA levels. Each sample was analyzed in triplicate. The PCR primers used are listed in Table 1. The human ACTB and GAPDH were used as endogenous controls for mRNA and lncRNA expression in analysis.

Table 2
Top 40 DEmRNAs

Symbol	log2FoldChange		$P$ value	FDR	Up/down
CA9	5.	640938	2.63E-86	1.40E-82	Up
MMP11	4.	915466	2.20E-78	5.02E-75	Up
IL11	4.	212103	1.63E-73	2.37E-70	Up
COLGALT1	1.	440801	1.27E-70	1.56E-67	Up
COL10A1	5.	129459	2.26E-68	2.12E-65	Up
GPRIN1	2.	08661	1.17E-61	8.12E-59	Up
CDCA5	1.	915155	2.66E-61	1.66E-58	Up
GRIN2D	3.	749138	4.01E-60	2.21E-57	Up
MFAP2	2.	870943	2.08E-59	1.07E-56	Up
MYBL2	2.	074464	4.17E-57	1.80E-54	Up
EME1	1.	908775	6.67E-57	2.81E-54	Up
INHBA	3.	875959	8.86E-57	3.63E-54	Up
HOXA10	3.	481532	2.10E-56	8.18E-54	Up
LAMC2	3.	446091	8.44E-56	3.21E-53	Up
KIF2C	1.	785265	1.25E-55	4.66E-53	Up
GJC1	2.	721089	1.16E-54	4.12E-52	Up
HOXD10	3.	562357	1.22E-54	4.26E-52	Up
NETO2	1.	895504	1.50E-54	5.10E-52	Up
MMP13	5.	691391	1.81E-54	6.04E-52	Up
COL13A1	2.	358811	2.71E-54	8.85E-52	Up
PRH1	$-$ 7.	11762	2.89E-118	4.62E-114	Down
AGFG2	$-$ 2.	36505	1.91E-86	1.40E-82	Down
NUCB2	$-$ 2.	13622	8.68E-85	3.47E-81	Down
CAB39L	$-$ 2.	19268	6.04E-82	1.93E-78	Down
CCL28	$-$ 4.	09831	1.12E-79	2.98E-76	Down
NOSTRIN	$-$ 2.	47725	5.82E-78	1.16E-74	Down
GPD1L	$-$ 2.	72233	8.82E-76	1.57E-72	Down
EMP1	$-$ 3.	14287	1.19E-75	1.91E-72	Down
CLEC3B	$-$ 3.	51315	6.24E-71	8.31E-68	Down
MPO	$-$ 4.	50315	1.48E-69	1.70E-66	Down
PLIN1	$-$ 5.	32039	9.95E-69	1.06E-65	Down
GPD1	$-$ 5.	00532	1.23E-68	1.23E-65	Down
ETFDH	$-$ 1.	57206	9.04E-68	8.03E-65	Down
LPO	$-$ 7.	38367	2.02E-64	1.70E-61	Down
UBL3	$-$ 1.	6461	2.33E-63	1.86E-60	Down
SLC25A4	$-$ 3.	20134	3.90E-62	2.97E-59	Down
PRB3	$-$ 6.	97523	9.74E-62	7.08E-59	Down
PFKFB1	$-$ 2.	94389	1.50E-61	1.00E-58	Down
ZG16B	$-$ 6.	91906	2.71E-61	1.66E-58	Down
FAM166B	$-$ 3.	38307	6.07E-61	3.60E-58	Down

Table 3

All of DElncRNAs

Symbol	log2FoldChange		$P$ value	FDR	Up/down
AC005532.5	$-$ 3.	97031	7.59E-50	2.79E-47	Down
AC024592.9	$-$ 3.	01971	1.35E-26	7.10E-25	Down
AC108676.1	3.	602306	1.22E-28	7.99E-27	Up
AFAP1.AS1	5.	804364	1.49E-52	7.82E-50	Up
CASC9	3.	315859	1.36E-25	6.32E-24	Up
CTB.51J22.1	$-$ 3.	18292	2.15E-37	3.16E-35	Down
DUXAP8	3.	384062	1.07E-41	2.63E-39	Up
FALEC	$-$ 3.	27553	2.82E-34	3.34E-32	Down
FOXCUT	$-$ 3.	84581	8.77E-49	2.93E-46	Down
HAGLROS	3.	163708	3.67E-40	7.93E-38	Up
HCG22	$-$ 4.	73685	4.61E-41	1.06E-38	Down
IL12A.AS1	$-$ 4.	03334	4.44E-74	8.16E-71	Down
LINC00520	3.	025544	2.39E-25	1.08E-23	Up
LINC00668	4.	12077	2.28E-29	1.68E-27	Up
LINC00941	3.	85692	6.90E-50	2.79E-47	Up
LINC01405	$-$ 3.	80476	3.59E-25	1.59E-23	Down
LINC01615	3.	25345	1.78E-36	2.51E-34	Up
LINC01694	3.	1717	3.96E-20	9.46E-19	Up
MIR9.3HG	3.	031614	1.06E-25	5.00E-24	Up
MYHAS	$-$ 3.	64474	6.07E-30	4.85E-28	Down
RNF144A.AS1	3.	079565	1.43E-47	4.39E-45	Up
RP11.159F24.6	3.	252705	5.26E-63	4.83E-60	Up
RP11.218E20.3	3.	027075	6.64E-44	1.74E-41	Up
RP11.284F21.9	3.	37323	8.94E-25	3.82E-23	Up
RP11.351J23.1	$-$ 4.	41511	1.65E-31	1.59E-29	Down
RP11.54H7.4	3.	033871	1.72E-25	7.91E-24	Up
RP11.594N15.3	$-$ 3.	36	1.01E-28	6.90E-27	Down
RP11.68I18.10	$-$ 3.	03083	2.67E-26	1.38E-24	Down
RP11.863P13.3	4.	025531	2.02E-59	1.49E-56	Up
RP3.323A16.1	3.	641611	3.24E-35	4.25E-33	Up
SLC8A1.AS1	$-$ 4.	95298	1.81E-83	6.66E-80	Down
TSPEAR.AS2	3.	481242	4.39E-22	1.35E-20	Up

Figure 1.

Hierarchical clustering analysis of DElncRNAs and top 100 DEmRNAs between HNSCC and normal tissues. A. DEmRNAs. B. DElncRNAs. Row and column represented DElncRNAs/DEmRNAs and tissue samples, respectively. The color scale represented the expression levels.

Figure 2.

Identification of optimal lncRNA biomarkers for HNSCC. A. The importance value of each DElncRNA ranked according to the mean decrease in accuracy by using the random forest analysis. B. The variance rate of classification performance when increasing numbers of the predictive DElncRNAs. C. Hierarchical clustering analysis of 13 lncRNAs biomarkers.

3. Results

3.1 DEmRNAs and DElncRNAs in HNSCC

Base on the DEmRNAs analysis, we found 3363 DEmRNAs in HNSCC patients compared with normal, including 1822 down-regulated mRNAs and 1541 up-regulated mRNAs. Base on the DElncRNAs analysis, we defined 32 DElncRNAs in HNSCC patients compared with normal, including 13 down-regulated lncRNAs and 19 up-regulated lncRNAs. Hierarchical clustering analysis of top 100 DEmRNAs and all of DElncRNAs were illuminated in Fig. 1A and Fig. 1B, respectively. Top 40 DEmRNAs and all of DElncRNAs were listed in Tables 2 and 3, respectively.

3.2 Identification of the optimal diagnostic lncRNAs biomarkers for HNSCC

The random forest analysis was used to rank the 32 DElncRNAs according to the mean decrease in accuracy (Fig. 2A). Ten-fold cross-validation result showed that the average accuracy rate of 13 DElncRNAs (IL12A.AS1, RP11.159F24.6, RP11.863P13.3, LINC 00941, FOXCUT, RNF144A.AS1, RP11.218E20.3, HCG22, HAGLROS, LINC01615, RP11.351J23.1, AC024592.9 and MIR9.3HG) reached the higher score for the first time (Fig. 2B). Hierarchical clustering analysis of the 13 DElncRNAs were displayed in Fig. 2C. Thence, these 13 DElncRNAs were identified as the optimal diagnostic lncRNAs biomarkers for HNSCC which were used to build the random forests, decision tree and SVM models. Violin Plot showed the expression levels of these 13 DElncRNAs between HNSCC and normal tissues (Fig. 3). Among which, 5 DElncRNAs (IL12A.AS1, FOXCUT, HCG22, RP11.351J23.1 and AC024592.9) were down-regulated in HNSCC tissues compared with normal, and 8 DElncRNAs (RP11.159F24.6, RP11.863P13.3, LINC00941, RNF144A.AS1, RP11. 218E20.3, HAGLROS, LINC01615 and MIR9.3HG) were up-regulated in HNSCC tissues compared with normal. The AUC of the SVM model was 0.983 and the specificity and sensitivity of this model were 95.5% and 96.2% (Fig. 4A). The AUC of the decision tree model was 0.824 and the specificity and sensitivity of this model were 77.3% and 97.6%, respectively (Fig. 4B). The AUC of the random forests model was 0.983 and the specificity and sensitivity of this model were 93.2% and 97.8%, respectively (Fig. 4C). The ROC of all these 13 lncRNAs (IL12A.AS1, RP11. 159F24.6, RP11.863P13.3, LINC00941, FOXCUT, RNF144A.AS1, RP11.218E20.3, HCG22, HAGLROS, LINC01615, RP11.351J23.1, AC024592.9 and MIR9.3HG) were displayed in Fig. 4D–P. Taken together, the AUC of all these three lncRNAs and their combination were all greater than 0.567 which indicated these 13 DElncRNAs and their combination were related to HNSCC and could predict the development of HNSCC.

Figure 3.

Violin Plot displayed the expression levels of 13 lncRNAs biomarkers between HNSCC and normal tissues. The X-axis represented normal and HNSCC groups. The Y-axis represented gene expression levels.

Figure 4.

ROC analysis of three HNSCC-specific lncRNAs biomarkers. The x-axis shows 1-specificity and y-axis shows sensitivity.

Figure 5.

Survival analysis of lncRNAs biomarkers between HNSCC and normal tissues.

3.3 Survival analysis of optimal diagnostic lncRNAs biomarkers for HNSCC

According to the survival analysis, we found AC024592.9, LINC00941, LINC01615 and MIR9-3HG was significantly associated with the prognosis of patients with HNSCC (Fig. 5A–D).

3.4 DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

A total of 13 optimal DElncRNA biomarkers for HNSCC were co-expressed with 610 DEmRNAs, accounting for 1096 DElncRNA-DEmRNA co- expression pairs (Fig. 6). IL12A.AS1, RP11.159F24.6, RP11.863P13.3, LINC00941, FOXCUT, RNF144A. AS1, RP11.218E20.3, HCG22, HAGLROS, LINC 01615, RP11.351J23.1, AC024592.9 and MIR9.3HG were co-expressed with 9, 79, 184, 139, 4, 123, 105, 47, 3, 109, 74, 137 and 83 DEmRNAs, respectively.

Figure 6.

lncRNAs-mRNAs co-expression network. The ellipses and rhombuses were represented the mRNAs and lncRNAs, respectively. Green and pink color represented up-regulated and down-regulated, respectively. The black border indicates DEmRNAs and DElncRNAs.

Figure 7.

The function annotation of DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs. The x-axis shows -log P and y-axis shows GO terms or KEGG pathways. A. GO terms. B. KEGG pathways.

3.5 Functional annotation of DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

Base on GO enrichment analysis, we found that cell adhesion (FDR $=$ 2.83E-10), signal transduction (FDR $=$ 5.81E-05), plasma membrane (FDR $=$ 9.49E-22) and protein binding (FDR $=$ 1.59E-14) were significantly enriched GO terms (Fig. 7A–C). The KEGG pathways analysis suggested that focal adhesion (FDR $=$ 1.92E-21), ECM-receptor interaction (FDR $=$ 7.59E-24), pathways in cancer (FDR $=$ 3.94E-08) and cytokine-cytokine receptor interaction (FDR $=$ 0.006683) were four significantly enriched pathways.

3.6 Verification by qRT-PCR

Three DElncRNAs and three DEmRNAs were chose for confirmation. As Fig. 8 displayed, FOXCUT was down-regulated and LINC00941, LINC01615, ITGA6, MMP13 and FOXC1 were up-regulated in HNSCC compared with adjacent tissues. Base on the TCGA results, FOXCUT and FOXC1 were down-regulated while the other lncRNAs and DEmRNAs were up-regulated in HNSCC compared to adjacent tissues. In generally, the confirmation results of qRT-PCR were consistent with our integrated analysis.

Figure 8.

Validation optimal lncRNA biomarkers in HNSCC tissue by qRT-PCR. ${}^{*}P<$ 0.05.

4. Disscussion

HNSCC, a frequent cancer of the head and neck region, has the high incidence and mortality rate in developing countries nowadays [8]. Although the diagnosis and treatment have been advanced during recent years, HNSCC still ranks among the seventh most common type of cancer worldwide [1]. Thence, exploring for diagnosis and prognosis biomarkers of HNSCC is an urgent problem. Through the present study, we aimed to obtain DEmRNA and DElncRNA in HNSCC patients compared with normal base on TCGA dataset. The machine learning and survival analysis were used to uncover the potential diagnostic and prognostic value of lncRNAs for HNSCC. We obtained 13 optimal diagnostic lncRNA biomarkers including IL12A.AS1, RP11.159F24.6, RP11.863P13.3, LINC00941, FOXCUT, RNF144A.AS1, RP11.218E 20.3, HCG22, HAGLROS, LINC01615, RP11.351J 23.1, AC024592.9 and MIR9.3HG. Among which, AC024592.9, LINC00941, LINC01615 and MIR9-3HG were significantly associated with the prognosis of patients with HNSCC.

To our knowledge, except of IL12A.AS1, RP11.159 F24.6, LINC00941, FOXCUT, HCG22, HAGLROS, LINC01615 and AC024592.9, the present study was the first to identify the RP11.863P13.3, RNF144A.AS1, RP11.218E20.3, RP11.351J23.1 and MIR9.3HG in HNSCC. Li et al. reported that IL12A.AS1 and HCG22 were down-regulated and RP11.159F24.6 was up-regulated between laryngeal squamous-cell carcinoma and adjacent tissues [9]. Zou et al. found that HCG22 and AC024592.9 were down-regulated in HNSCC carcinoma compared with adjacent tissues [10]. HAGLROS was elevated in nasopharyngeal carcinoma, and HAGLROS contributes to nasopharyngeal carcinoma progress via regulating miR-100/ATG14 axis-mediated PI3K/AKT/mTOR signals [11]. Sassenberg et al. displayed that LINC01615 was up-regulated in HNSCC carcinoma compared with adjacent tissues [12]. In this study, IL12A.AS1, HCG22 and AC024592.9 were down-regulated and RP11.159F24.6, LINC01615 and HAGLROS was up-regulated in HNSCC patients compared with normal, which were consistent with reports of other researchers, indicating the our integration analysis data were reliable.

Long intergenic non-protein coding RNA 941 (LINC00941) is an lncRNAs considered as a candidate biomarkers of multiple cancer. LINC00941 is overexpressed in liver cancer and correlated with poor prognosis, and LINC00941 overexpression promoted epithelial to mesenchymal transition of liver cancer [13]. In lung cancer, LINC00941 showed prognostic values and regulated focal adhesion and PI3K-AKT signaling pathway [14]. LINC00941 is overexpressed in gastric cancer tissues and may be a potential novel biomarker for therapeutic of gastric cancer [15]. Liu et al. uncovered that LINC00941 plays a key oncogenic function in gastric cancer and may be used as a potential biomarker for diagnosis and prognosis of gastric cancer [16]. It is reported that LINC00941 was up-regulated in HNSCC carcinoma compared with adjacent tissues [17]. In this study, LINC00941 was up-regulated in both TCGA integration analysis and qRT-PCR validation. The machine learning showed that LINC00941 was a one of optimal diagnostic lncRNAs biomarkers. The survival analysis results displayed that LINC00941 was significantly associated with the prognosis of patients with HNSCC. Therefore, LINC00941 was considered as a potential biomarker for diagnosis and prognosis of HNSCC.

Here, integrin subunit alpha 6 (ITGA6), belongs to the integrin family, was enriched in focal adhesion and ECM-receptor interaction pathway. Silencing of ITGA6 expression was shown to inhibit cell migration and invasion in head and neck cancer cells [18]. Yang et al. reported that ITGA6 in ECM-receptor interaction signaling pathway may play an crucial role in HNSCC [19]. ITGA6 was identified to be significantly correlated with HNSCC [20]. In this study, ITGA6 co-expressed with LINC00941 and RP11.159F24.6 was enriched in focal adhesion and ECM-receptor interaction pathway. Thence, we speculated that LINC00941 and RP11.159F24.6 might be involved in the occurrence of HNSCC by regulating focal adhesion and ECM-receptor interaction pathway.

Jin et al. showed that matrix metallopeptidase 13 (MMP13) was DEGs in HNSCC, suggesting that MMP13 may play a key roles in HNSCC progression [21]. MMP13 was up-regulated in HNSCC, and it has potential value as a marker predictive of the decreased efficacy of anti-PD-1 therapy [22]. In the study, MMP13 was up-regulated in HNSCC and co-expressed with RP11.863P13.3, RP11.218E20.3 and LINC01615. Therefore, we concluded that RP11. 863P13.3, RP11.218E20.3 and LINC01615 might be involved in the progression of HNSCC by regulating MMP13.

Collagen type IV alpha 1 chain (COL4A1), belongs to collagen family, is associated with cell adhesion. Chen et al. found that COL4A1 might be as a biomarker to distinguish oral squamous cell carcinoma and controls [23]. COL4A1 was high expressed in HNSCC which indicates that COL4A1 may participate in HNSCC cell adhesion [24]. In this study, COL4A1 was co-expressed with RP11.159F24.6, RP11.863P13.3, LINC00941 and RNF144A.AS1, and COL4A1 was enriched in pathway of focal adhesion, ECM-receptor interaction and pathways in cancer. These findings demonstrated that COL4A1 might play a pivotal role in the progression of HNSCC by regulating pathway of focal adhesion, ECM-receptor interaction and pathways in cancer.

It has been reported that the fork box C1 (FOXC1) is involved in various cancer cell proliferation, migration and invasion [25, 26, 27]. Gao et al. displayed that FOXC1 is involved in the invasion and metastasis of laryngeal squamous cell carcinoma, suggesting that FOXC1 might be a risk factor for laryngeal squamous cell carcinoma [28]. In this study, FOXC1 co-expressed with FOXCUT in HNSCC. FOXCUT may play a key role in the tumor progression and survival of esophageal squamous s cell carcinoma patients by modulating FOXC1 [29]. Functional of FOXC1 and FOXCUT is correlated with tumor progression of oral squamous cell carcinoma, and FOXC1 and FOXCUT may serve as novel biomarkers in oral squamous cell carcinoma patients [30]. One study found that FOXCUT expression contributed to the proliferation and migration of nasopharyngeal carcinoma by targeting FOXC1 and that FOXCUT may be as a novel nasopharyngeal carcinoma biomarker [31]. We provide further evidence that FOXCUT was a novel diagnosis biomarkers of HNSCC.

In conclusion, 13 DElncRNAs were identified as potential biomarkers of HNSCC. Among which, 4 DElncRNAs were significantly associated with the prognosis of patients with HNSCC. Overall, the confirmation results of qRT-PCR were generally consistent with our integrated analysis. Several limitations still detected in our study, the sample size for qRT-PCR confirmation was small, and large numbers of HNSCC samples are needed for further research. Identifying diagnostic and prognostic biomarkers for HNSCC is a pilot study, and further experiments are needed to understand the biological of key lncRNAs in HNSCC.

Footnotes

Acknowledgments

This study was supported by the Natural science foundation of hunan province (No: 2017JJ2178).

Conflict of interest

The authors declare that they have no conflict of interest.

References

Dhull

A.K.

Atri

Dhankhar

Chauhan

A.K.

and Kaushal

, Major risk factors in head and neck cancer: a retrospective analysis of 12-year experiences, World J Oncol 9(3) (2018), 80–84.

Ferlay

Soerjomataram

Dikshit

Eser

Mathers

Rebelo

Parkin

D.M.

Forman

and Bray

, Cancer incidence and mortality worldwide: sources, methods and major patterns in GLOBOCAN 2012, Int J Cancer 136(5) (2015), E359–386.

Fan

Zhao

Tan

Hua

and Wu

, Identification of aberrantly expressed long non-coding RNAs in stomach adenocarcinoma, Oncotarget 8(30) (2017), 49201–49216.

Dey

B.K.

Mueller

A.C.

and Dutta

, Long non-coding RNAs as emerging regulators of differentiation, development, and disease, Transcription 5(4) (2014), e944014.

Zeng

J.H.

Liang

R.Q.

Tang

R.X.

Cai

X.Y.

Chen

J.Q.

Luo

D.Z.

and Chen

, Comprehensive investigation of a novel differentially expressed lncRNA expression profile signature to assess the survival of patients with colorectal adenocarcinoma, Oncotarget 8(10) (2017), 16811–16828.

Pan

Liu

Wang

and Li

, Analysis of lncRNA-Mediated ceRNA crosstalk and identification of prognostic signature in head and neck squamous cell carcinoma, Front Pharmacol 10 (2019), 150.

Wei

Liang

Liu

Huang

Chen

Guo

and Li

, Bioinformatics profiling utilized a nine immune-related long noncoding RNA signature as a prognostic target for pancreatic cancer, J Cell Biochem (2019).

Siegel

Naishadham

and Jemal

, Cancer statistics, 2012, CA Cancer J Clin 62(1) (2012), 10–29.

Wang

F.L.

Fei

Y.H.

Wang

Y.T.

Zhang

J.E.

H.Y.

and Zhang

, Aberrant expressed long non-coding RNAs in laryngeal squamous-cell carcinoma, Am J Otolaryngol 40(5) (2019), 615–625.

10.

Zou

Guan

Zhang

Wang

Yuan

D.Y.

Meng

and Zhang

, Integrated analysis of lncRNA CTD-2357A8. 3 expression and its potential roles in head and neck squamous cell carcinoma, Oncology Letters.

11.

Zhang

and Xi

, Upregulation of lncRNA HAGLROS enhances the development of nasopharyngeal carcinoma via modulating miR-100/ATG14 axis-mediated PI3K/AKT/mTOR signals, Artif Cells Nanomed Biotechnol 47(1) (2019), 3043–3052.

12.

Sassenberg

Droop

Schulz

W.A.

Dietrich

Loick

S.M.

Wiek

Scheckenbach

Gaisa

N.T.

and Hoffmann

M.J.

, Upregulation of the long non-coding RNA CASC9 as a biomarker for squamous cell carcinoma, BMC Cancer 19(1) (2019), 806.

13.

Yan

Zhang

Qian

Hao

Yan

Zhu

Huang

Luo

Liu

Chen

and Fan

, Mesenchymal stem cells promote hepatocarcinogenesis via lncRNA-MUF interaction with ANXA2 and miR-34a, Cancer Res 77(23) (2017), 6704–6716.

14.

Wang

Zhao

Deng

Bai

Xiao

and Zhang

, Systematic identification of lincRNA-based prognostic biomarkers by integrating lincRNA expression and copy number variation in lung adenocarcinoma, Int J Cancer 144(7) (2019), 1723–1734.

15.

Luo

Tao

Zhang

Zhu

Minyao

D.N.

Haleem

Dong

Zhang

Zhao

and Liao

, Regulatory network analysis of high expressed long non-coding RNA LINC00941 in gastric cancer, Gene 662 (2018), 103–109.

16.

Liu

Zhang

Zhong

Zhang

Guo

Nie

and Liu

, Long non-coding RNA LINC00941 as a potential biomarker promotes the proliferation and metastasis of gastric cancer, Front Genet 10 (2019), 5.

17.

Zou

A.E.

Honda

T.K.

Kuo

S.Z.

Zheng

Xuan

Saad

M.A.

Hinton

Brumund

K.T.

Lin

J.H.

Wang-Rodriguez

and Ongkeko

W.M.

, Transcriptome sequencing uncovers novel long noncoding and small nucleolar RNAs dysregulated in head and neck squamous cell carcinoma, RNA 21(6) (2015), 1122–1134.

18.

Kinoshita

Nohata

Hanazawa

Kikkawa

Yamamoto

Yoshino

Itesako

Enokida

Nakagawa

Okamoto

and Seki

, Tumour-suppressive microRNA-29s inhibit cancer cell migration and invasion by targeting laminin-integrin signalling in head and neck squamous cell carcinoma, Br J Cancer 109(10) (2013), 2636–2645.

19.

Yang

Chen

Huang

Han

and Li

, Identification of potential biomarkers and analysis of prognostic values in head and neck squamous cell carcinoma by bioinformatics analysis, Onco Targets Ther 10 (2017), 2315–2321.

20.

Zhao

Chi

Cao

Cui

Meng

Guo

and Wang

, Screening and clinical significance of tumor markers in head and neck squamous cell carcinoma through bioinformatics analysis, Mol Med Rep 19(1) (2019), 143–154.

21.

Jin

and Yang

, Identification and analysis of genes associated with head and neck squamous cell carcinoma by integrated bioinformatics methods, Mol Genet Genomic Med (2019), e857.

22.

Hira-Miyazawa

Nakamura

Hirai

Kobayashi

Kitahara

Bou-Gharios

and Kawashiri

, Regulation of programmed-death ligand in the human head and neck squamous cell carcinoma microenvironment is mediated through matrix metalloproteinase-mediated proteolytic cleavage, Int J Oncol 52(2) (2018), 379–388.

23.

Chen

Mendez

Houck

Fan

Lohavanichbutr

Doody

Yueh

Futran

N.D.

Upton

Farwell

D.G.

Schwartz

S.M.

and Zhao

L.P.

, Gene expression profiling identifies genes predictive of oral squamous cell carcinoma, Cancer Epidemiol Biomarkers Prev 17(8) (2008), 2152–2162.

24.

Kuang

Zhao

Dang

and Li

, Identification of potential therapeutic target genes and mechanisms in head and neck squamous cell carcinoma by bioinformatics analysis, Oncol Lett 11(5) (2016), 3009–3014.

25.

Ray

P.S.

Wang

Sim

M.S.

Shamonki

Bagaria

S.P.

Liu

Elashoff

Hoon

D.S.

Walter

M.A.

Martens

J.W.

Richardson

A.L.

Giuliano

A.E.

and Cui

, FOXC1 is a potential prognostic biomarker with functional significance in basal-like breast cancer, Cancer Res 70(10) (2010), 3870–3876.

26.

Wei

L.X.

Zhou

R.S.

H.F.

Wang

J.Y.

and Yuan

M.H.

, High expression of FOXC1 is associated with poor clinical outcome in non-small cell lung cancer patients, Tumour Biol 34(2) (2013), 941–946.

27.

Z.Y.

Ding

S.M.

Zhou

Xie

H.Y.

Chen

K.J.

Zhang

Xing

C.Y.

Guo

H.J.

and Zheng

S.S.

, FOXC1 contributes to microvascular invasion in primary hepatocellular carcinoma via regulating epithelial-mesenchymal transition, Int J Biol Sci 8(8) (2012), 1130–1141.

28.

Gao

Zhang

Zhu

Chen

Liu

Wen

and Wang

, MicroRNA-204-5p inhibits invasion and metastasis of laryngeal squamous cell carcinoma by suppressing forkhead box C1, J Cancer 8(12) (2017), 2356–2368.

29.

Pan

Yao

Chen

Zhou

Geng

Mao

and Fang

, A novel long non-coding RNA FOXCUT and mRNA FOXC1 pair promote progression and predict poor prognosis in esophageal squamous cell carcinoma, Int J Clin Exp Pathol 7(6) (2014), 2838–2849.

30.

Kong

X.P.

Yao

Luo

Feng

F.K.

J.T.

Ren

Y.P.

Wang

D.L.

and Bu

R.F.

, The expression and functional role of a FOXC1 related mRNA-lncRNA pair in oral squamous cell carcinoma, Mol Cell Biochem 394(1-2) (2014), 177–186.

31.

Y.Z.

Chen

F.F.

Zhang

Zhao

Q.F.

Guan

X.L.

Wang

H.Y.

Song

S.S.

Zhou

and Li

X.J.

, The long noncoding RNA FOXCUT promotes proliferation and migration by targeting FOXC1 in nasopharyngeal carcinoma, Tumour Biol 39(6) (2017), 1010428317706054.

Screening key lncRNAs with diagnostic and prognostic value for head and neck squamous cell carcinoma based on machine learning and mRNA-lncRNA co-expression network analysis

Abstract

BACKGROUND:

METHODS:

RESULTS:

CONCLUSION:

Keywords

1. Introduction

2. Materials and methods

2.1 TCGA data

2.2 Identification of DEmRNAs and DElncRNAs

2.3 Identification of the optimal diagnostic lncRNAs biomarkers for HNSCC

2.4 Survival analysis of optimal diagnostic lncRNAs biomarkers for HNSCC

2.5 DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

2.6 Functional annotation of DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

Table 1 Primer sequences used for qRT-PCR

Table 2 Top 40 DEmRNAs

3.1 DEmRNAs and DElncRNAs in HNSCC

3.2 Identification of the optimal diagnostic lncRNAs biomarkers for HNSCC

3.4 DEmRNAs co-expressed with the identified optimal diagnostic lncRNAs

3.6 Verification by qRT-PCR

Footnotes

Acknowledgments

Conflict of interest

References

Table 1
Primer sequences used for qRT-PCR

Table 2
Top 40 DEmRNAs