Prognostic value of lipid metabolism‐related genes in head and neck squamous cell carcinoma

Altered lipid metabolism is involved in the development of many tumors. However, the role of dissimilar lipid metabolism in head and neck squamous cell carcinoma (HNSCC) is not fully established.

mediate intercellular communication in tumor immune microenvironment (TIME). 4 Prostaglandin E 2 , produced by cyclooxy-genase-2, has been shown to promote TIME and block production of Type I interferon, preventing tumor elimination. 5 Head and neck squamous cell carcinoma (HNSCC) is a multisite malignancy that affects about 600,000 patients annually, worldwide. 6 Despite advances in treatment approaches, the 5-year survival rate remains about 50%. 7 In recent years, immunotherapy has gained a lot of interest for cancer treatment. However, it is only effective in 20%-30% of patients. 8 Thus, it is necessary to study multicomponent antitumor responses and develop new anti HNSCC therapeutics. Current studies have confirmed that lipid metabolism promotes HNSCC invasion. 9 However, other links between lipid metabolism and HNSCC are unclear.
Here, we sought to determine the survival role of lipid metabolism-related genes and their therapeutic potential against HNSCC. We also systematically assessed the relationship between lipid metabolism-related genes and overall HNSCC survival, as well as the correlation between prognostic index (PI) and six types of immune cells using bioinformatics. Findings from this study may improve personalized HNSCC treatment.

| Prognostic analysis
The R survival package was used to assess relationships between lipid metabolism-related genes and overall survival (p < .05) and to plot survival curves. Prognosis related lipid metabolism genes were visualized on forest maps. Multivariate Cox analysis was used to select signature genes and determine the PI for each patient (p < .05). Risk curves were drawn using the pheatmap package. ROC maps were created using the survival ROC package. Univariate and multivariate independent prognostic analysis by the survival package were used to determine if PI could be used as an independent prognostic indicator.

| Correlation analyses
Correlation maps between clinical features and signature genes were drawn using beeswarm package. Tumor infiltration immune cells data were obtained from tumor immune estimation resource (TIMER; https://cistrome. shinyapps.io/timer/) and included data on B-cells, CD4 + T-cells, CD8 + T-cells, neutrophils, macrophages, and dendritic cells. 14 These raw data were used to construct a relationship landscape between PI and six immune cell types (p < .05).

| Statistical analysis
R software was used for results visualization. Differences between two groups were compared using independent t test. p < .05 indicated statistically significant differences.

| Summary of results
Our analysis identified 4783 DEGs. Of these, 136 are associated with lipid metabolism and 63 are TFs. Differentially  lipid metabolism-related genes were subjected to GO term and KEGG pathway analysis. We then assessed the mutations of the prognosis-related lipid metabolism genes and constructed a regulatory network of TFs and lipid metabolism-related genes. We further established a prognosis model of lipid metabolism-related genes and divided all patients into two groups based on PI. The gene set was validated on HPA database and their relationship to clinical traits and immune cells determined.

| Identification of differentially expressed lipid metabolism-related genes
There were 4783 HNSCC DEGs, of which 136 are lipid metabolism-related genes. A total of 3602 DEGs were highly expressed and 1181 lowly expressed ( Figure 1A-C). 64 differentially expressed lipid metabolism-related genes were upregulated and 72 downregulated ( Figure 1D). GO term analysis of these lipid metabolism-related DEGs revealed that "fatty acid metabolic process," "lipid droplet," and "cofactor binding" were commonest in biological progress, cellular components, and molecular functions (MF), respectively ( Figure 1E). The KEGG pathway analysis confirmed that glycerol phospholipid metabolism was the most abundant pathway ( Figure 1F).

| Mutation of prognosis-associated lipid metabolism-related DEGs
Univariate analysis identified 23 prognosis-related lipid metabolism genes (Table 1). Forest map analysis found that 14 prognosis-related lipid metabolism genes were high-risk genes ( Figure 2). Mutation analysis revealed that many prognosis-related lipid metabolism genes had missense mutations ( Figure 3).

| Transcriptional regulation of prognosis-associated lipid metabolism-related genes
A total of 63 DEGs were TFs ( Figure 4A), of which 46 were highly expressed and 17 lowly expressed ( Figure 4B). Transcription regulation results showed that eigth highrisk lipid metabolism-related genes (GLA, ACAT1, DGKG, TXNRD1, PTDSS1, SMS, PIP5K1B, and CAV1) and eight TFs (PPAPG, MYH11, H2AFX, HEY1, CBX2, SNAI2, SPDEF, and FOXA2) were strongly associated with poor HNSCC prognosis. Moreover, five low-risk lipid metabolism-related genes, including LIPE, CYP2D6, ACACB, FLA2G2D, and ACSM3, were regulated by 13 TFs (FOXP3, HOXB13, LIN9, EZH2, E2F7, E2F1, DNMT1, CDK2, BRCA1, POU5F1, FBX1, MEF2C, and LMNB1; Figure 4C).  . Next, patients were divided into high-and lowrisk group based on median PI and a risk curve constructed ( Figure 5A-C). Figure 6A shows that low-risk patients had higher survival rates relative to high-risk patients. The PI based on lipid metabolism-related genes can distinguish patients with HNSCC based on potential discrete clinical outcomes. Thus, we established that PI effectively and accurately stratifies patients with HNSCC. Area under the curve values of the 1-, 3-, and 5-year were 0.664, 0.724, and 0.623, respectively ( Figure 6B-D), indicating a high prognostic value. Univariate and multivariate analyses were used to estimate the value of the acquired PI (Figure 7). Univariate risk analysis identified clinical stage, N stage and PI, as independent prognosis predictors. However, comprehensive analysis of all clinical information indicated that only PI was an independent predictor (p < .001). HPA database IHC data validated the risk model, except for CYP27B1 and PLA2G2D ( Figure 8).

| Correlation analysis of clinical and immune cell infiltration
Multivariate Cox regression analysis of the relationship between lipid metabolism-related genes and clinical characteristics like age, sex, pathological stage, clinical stage, T stage, and N stage ( Table 2) revealed statistically significant association with different genes (Figure 9). Analysis of the relationship between risk genes and clinical features revealed that CYP2D6, LPIN1, and PLA2G2D were highly expressed in G3&4 relative to G1&2 (p = .014, p = .005, and p = .022). PIP5K1B was highly expressed in females (p = .048). DGKG, PHYH, TRIB3, and risk score positively correlated with clinical stage (p = .008, p < .001, p < .001, and p = .007). Risk score and TRIB3 were positively correlated with T stage (p = .006 and p = .015), while RDH16 negatively correlated with T stage (p = .038). DGKG, LPIN1, and PHYH positively correlated with N stage (p = .017, p = .022, and p = .008). Examination of the relationship between PI and immune cell infiltration to determine if these risk genes precisely reflect HNSCC immune environment revealed that lipid metabolism-related genes PI negatively correlated with CD4 + T-cells, CD8 + T-cells, neutrophils cell, B-cells, and dendritic cells (p = .002, p = .004, p = .025, p < .001, and p = .005; Figure 10).

| DISCUSSION
HNSCC is characterized by a low cure rate and easy recurrence. Multiple studies have shown that aberrant lipid metabolism influences HNSCC tumorigenesis. In HNSCC models, cetuximab resistance has been linked to lipid metabolism alterations. 15 Reduced expression of metabolism-related genes may improve survival of HPVpositive patients with head and neck cancer (HNC). 16 FADS1 has emerged as mediator of lipid metabolism gene and drives laryngeal squamous cell carcinoma progression by activating AKT/mTOR signaling pathway. 17 In oral squamous cell carcinoma (OSCC), the Kennedy signaling pathway is upregulated by cholesterol and glycerophospholipid (GPL) metabolic changes. 18 While lipid metabolism is vital for HNC, its role is not completely understood in HNSCC.
Here, we used bioinformatics to comprehensively study lipid metabolism-related genes in HNSCC. Our data highlight the value of these genes in anti-HNSCC therapy. We confirmed that some lipid metabolismrelated genes are involved in HNSCC progression and RNA-seq analysis of lipid metabolism-related genes in HNSCC tumor tissue versus adjacent noncancer tissue found that lipid metabolism-related genes are altered in HNSCC. These DEGs exhibit enrichment for fatty acid metabolic processes. Elevated fatty acid synthesis has been observed in various cancers. Multiple studies have confirmed that lipogenesis is important for tumor growth. 19 The MFs analysis revealed cofactor binding enrichment. An active RNA editing complex is formed by APOBEC1 cytosine deaminase and cofactor A1CF. This complex acts on APOB RNA to regulate lipid metabolism. 20 The KEGG pathway analysis revealed GPL metabolism as the most enriched pathway. In OSCC, GPL metabolism is the most abundant form of lipid metabolism. 18 These data indicated the molecular mechanism underlying lipid metabolism in HNSCC.
Cox regression analysis revealed 11 lipid metabolismrelated genes associated with HNSCC prognosis. ROC indicated that our results are highly accurate. Lipid metabolism-related genes PI had a high prognosis value and were related to age, pathological tumor stage, TN stage, and gender. It is reported that CYP27B1, CYP2D6, and TRIB3 are directly or indirectly involved in HNSCC development and progression. [36][37][38][39] DGKG, DHCR7, LPIN1, PIP5K1B, PLA2G2D, and RDH16 are involved in the development of other tumors. 26,28,[40][41][42][43] However, ARSI and PHYH have not been previously studied in cancers. In addition, the gene set was validated using HPA data on HNSCC and other tumors. Further analysis of the relationship between these genes will offer theoretical basis for HNSCC treatment.
27-hydroxyl cholesterol indirectly affects breast cancer metastasis by increasing the number of metastasizing immune cells and suppressing CD8 + T-cells. 44 Increased TIME cholesterol may cause CD8 + T-cell exhaustion. 45 Together, these findings indicate that abnormal lipid metabolism may affect the number and function of tumor immune cells. Thus, we designed a simple protocol for monitoring immune status in patients with HNSCC based on lipid metabolism-associated genes. The results indicated that CD4 + T-cells, CD8 + T-cells, neutrophils, B-cells, and dendritic cells had a higher degree of infiltration in low-risk patients with HNSCC. Consistent with past studies, [46][47][48][49] we found that CD4 + T-cells and CD8 + T-cells were necessary for antitumor immunity and improved prognosis outcome. Neutrophils, B-cells, and dendritic cells were comprised of different subpopulations with pro-and antitumor function. Here, we found that these three immune cell types have antitumor roles in HNSCC. [50][51][52][53][54][55][56][57][58] In addition, tumor associated macrophages may be polarized into M1-like macrophages with anticancer activity, or M2 macrophages with pro-cancer roles. 59 We found that macrophages (predominantly M1) have a higher degree of infiltration in low-risk HNSCC. Our findings validate and enhance the recognition of immune cells functions in HNSCC. In summary, we conducted a comprehensive analysis of HNSCC mutations, the regulatory mechanisms of lipid metabolism-related genes, the relationship between these genes and clinical prognosis and immune infiltration, and determined the genes' prognostic value. The risk score of lipid metabolism-related genes may reflect the prognosis and immune status of patients. However, the prognostic value of lipid metabolism-related genes in HNSCC needs further experimental verification.

CONFLICTS OF INTEREST
All the authors declare that there are no conflicts of interest.

AUTHOR CONTRIBUTIONS
Conceptualization, software, and methodology: Ying Xiong and Yu Si. Validation, investigation, and data curation: Ying Xiong, Yu Si, and YiSi Feng. Writing-original draft preparation: Ying Xiong, Yu Si, Shipei Zhuo, and Bozhen Cui. Supervision and funding acquisition: Ying Xiong, Yu Si, and Zhigang Zhang. All authors have read and agreed to the published version of the manuscript.

DATA AVAILABILITY STATEMENT
The authors confirm that the data supporting the findings of this study are available within the article.