mRNA data acquisition and processingData on mRNA and clinical characteristics of the TCGA cohort were obtained from the TCGA database. (https://tcga-data.nci.nih.gov/), which includes 547 samples (501 primary tumor samples, 2 metastasis samples, and 44 normal samples). On the basis of the AJCC stage, all patients are categorized as early stage (stages I and II) or advanced stage (stages III and IV). The mRNA data and clinical characteristics of the validation cohort (GSE65858) are acquired from the GEO database (https://www.ncbi.nlm.nih.gov/geo/), which includes 270 tumor samples and downloaded using R package GEOquery. The human ERS-related genes are taken from the MSigDB database (http://www.gseamsigdb.org/gsea/msigdb), and 258 genes are used in following study after removing duplicates and two microRNA genes (MIR199A1 and MIR200C, Table S1). R package DESeq2 is used for differential expression analysis and count data from the TCGA-HNSC cohort15, and genes with a p value less than 0.05 were considered as differentially expressed genes.Weighted correlation network analysis (WGCNA)WGCNA is an algorithm to group genes with similar expression patterns, analyzing the relationship between gene modules and phenotypes by converting gene expression into gene modules16. R package WGCNA is performed to divide ERS-related genes into different gene modules and select the most relevant gene modules with phenotypes.Establishment and validation of the ERS-related signatureLASSO-Cox regression (performed by R package glmnet) is chosen to overcome the multicollinearity problem between variables and identify important predictors. In the regression, the independent variable is the normalized expression matrix of selected genes, the dependent variable is the patient’s survival status and survival time, and the penalty regularization parameter λ is determined via the ten-fold cross validation. Multiple Cox regression is utilized to further discriminate the genes with prognosis value by R package survival, and the results are visualized by forest plot (ggforest function in R package survival). Then the ERS-related prognostic gene signature is constructed based on corresponding regression coefficients and normalized mRNA expression of the genes in multiple Cox regression. R package survival and survminer are utilized to overall survival and visualization of the results. The nomogram and calibration curves were constructed by R package rms.Tumor mutation burden analysisSomatic mutation data of the HNSC cohort is acquired from the TCGA database (including 510 samples and clinical information in 495 patients was consistent with the TCGA-HNSC cohort). The mutational landscape, tumor mutation burden (TMB) calculation, and visualization of protein domains mutation sites are implemented in R package maftools17.The assessment of immune infiltrationThe CIBERSORT method in R18 is utilized to evaluate the relative abundance of 22 types of human immune cells in mRNA samples in the HNSC cohort. An evaluation of immune infiltration, stromal content, and tumor purity can be performed with R package ESTIMATE19.scRNA-seq data processing, Gene Set Variation Analysis (GSVA), trajectory analysis and cell–cell communicationR package Seurat and harmony were performed for integrating and clustering three scRNA-seq datasets (GSE162025, GSE150430, GSE150321, which included 25 patients with nasopharyngeal cancer and 2 patients with laryngeal cancer.) from GEO database (https://www.ncbi.nlm.nih.gov/geo/). We utilized the Harmony package to integrate single-cell data and remove batch effects between datasets. R package COSG was used to identified marker genes in each group20. Risk score was computed by the average risk score of all cells in the sample, and then divided into risk groups by median. GSVA was performed to quantify pathway activation. Trajectory analysis was also performed using R package monocle2. Cellchat method was performed to visualize intercellular communications21.Cell cultureHuman normal and HNSC cell lines (NOK, SAS, UM1, CAL27, CAL33) were purchased from iCell Bioscience Inc. (Shanghai, China) and cultured in DMEM ( Invitrogen, CA, USA) with 10% FBS ( Sigma-Aldrich, MO, USA) and 1% penicillin–streptomycin (Gibco, NY, US). All cell lines were authenticated by Short Tandem Repeat (STR) profiling and kept in a humidified atmosphere of 5% CO2 at 37 °C.Total RNA isolation, reverse transcription and real-time quantitative PCRTrizol kit (Invitrogen, US) is chosen to isolate total cellular RNA, and the concentration and purity of total RNA are measured by absorbance at 260/280 nm. Total RNA is reverse transcribed using a reverse transcription kit (Yishan Science, Shanghai, China). Quantitative Real-Time PCR (qRT-PCR) is performed with 2 × RealStar Green Fast Mixture (Genstar, Beijing, China). The expression of target genes is normalized to the mRNA expression of GAPDH, and the 2 − ΔΔCT method is used to calculate fold changes. Primers sequences of target genes and GAPDH are presented as follow:GAPDH forwardTGACTTCAACAGCGACACCCAGAPDH reverseCACCCTGTTGCTGTAGCCAAAATF6 forwardTCCTCGGTCAGTGGACTCTTAATF6 reverseCTTGGGCTGAATTGAAGGTTTTGTRIB3 forwardGCCTTTTTCACTCGGACCCATTRIB3 reverseCAGCGAAGACAAAGCGACACUBXN6 forwardGGAGCGCATTAACTGCCTGUBXN6 reverseGCTCAGCACGTAGAACTCCTCWestern blotTotal protein is harvested using RIPA buffer (Beyotime, Jiangsu, China) and quantified by the bicinchoninic acid (BCA) protein Assay Kit (Beyotime, Jiangsu, China). Equal volume of samples is electrophoretically separated using 10% SDS-PAGE and then transferred to a PVDF membrane (Millipore, Bedford, MA, USA). Subsequently, the PVDF membrane is incubated with the primary antibody overnight at 4 °C and then incubated with the second antibody for 1.5h at room temperature. Finally, enhanced chemiluminescence assays are used to detect the signal. Primary antibodies are listed below: TRIB3(Abcam, Cambridge, UK), UBXN6(Proteintech, Wuhan, China), α-tubulin (Used as the loading control, Proteintech, Wuhan, China). The original western blot images in this study all were provided in Figure S1.SiRNA transfectionThe transfection was performed when cells reached 50% confluence in a 6-well culture dish. SiRNA was transfected in cells with jetPRIME transfection reagent (polyplus-transfection, New York, USA), and the medium is replaced with DMEM(+10%FBS) after 8h. After 48h, cells are collected for further assays. The sequences of siRNAs (Si#1 and Si#2) as follows:SiTRIB3#1GACAACUUAGAUACCGAGCGUTTACGCUCGGUAUCUAAGUUGUCTTSiTRIB3#2GAUCUCAAGCUGUGUCGCUUUTTAAAGCGACACAGCUUGAGAUCTTColony-forming assayCells are digested and counted in an automated cell counter (Nexcelom, Lawrence, USA). Then 800 cells were seeded into a 6-well plate and cultured for 7–10 days. 4% paraformaldehyde was performed to fix colonies (10 min), and the colonies were stained with 0.01% crystal violet solution for 10 min. Finally, we clean the plate with water and captured it with the camera.Wound healing assayWhen the cells are full of the 6-well plates, we straightly scratch three lines with a 200 μL pipette tip in the well. After scratching, we wash the cells with 1xPBS to remove any floating cells, observe them under a microscope, and photographs were taken in 0h and 24h.CCK-8 proliferation assay1000 cells per well were seeded in a 96-well plate and incubated with 100μL DMEM(+ 10%FBS). Every day 10 μl CCK-8 (Takara, Tokyo, Japan) was added in each well of the plate. Then the 96-well plate was incubated for 2h and measured the OD value at 450nm.Statistical analysis and graphs plottingStatistical graphs are plotted by R package ggpubr, R package ggplot2, R package corrplot, and GraphPad Prism (version 8.0). All statistical analyses were performed in R (venison 4.2.1) and GraphPad Prism. Student’s t-test(two-tail) is used for paired samples, and non-parametric Wilcoxon rank-sum test is used for unpaired samples. Pearson analysis was used to assess the correlation. Significance in survival analysis was assessed using the Log-Rank test, and the Gehan-Breslow method was employed when survival curves intersected. P values less than 0.05 (*p < 0.05) and 0.01 (**p < 0.01) were considered significant.
