Abstract
Hox genes, a highly conserved homolog in most animals, play vital functions in cell development and organ formation. In recent years, researchers have discovered that it can act as a tumor regulator, and its members can participate in tumorigenesis by regulating receptor signaling, cell differentiation, apoptosis, migration, EMT, and angiogenesis. Hox genes and which major members play a vital role in the progress of head and neck squamous cell carcinoma (HNSCC) is still unclear. After analyzing the expression differences and prognostic value of all Hox genes through the TCGA-HNSC database, we use histochemistry stains in 52 pairs of HNSCC slices to verify the expression level of the key member-HOXA1. In correlation analysis, we found that high HOXA1 expression is related to poor pathological grade (p = 0.0077), advanced T stage (p = 0.021) and perineural invasion (PNI) (p = 0.0019). Furthermore, we used Cox univariate and multivariate regression analysis to confirm the independent predictive power of HOXA1 expression. To explore the underlying mechanisms behind HOXA1, we ran GSVA and GSEA and found fourteen mutual signaling pathways, including neuroprotein secretion and transport, tumor-associated signaling pathways, cell adhere junction and metabolic reprogramming. Finally, we found that the high expression of HOXA1 is significantly related to the decrease of CD8+ T cell infiltration and the decline of DNA methylation level. Our findings demonstrated that HOXA1, as a notable member of the HOX family, maybe an independent prognostic indicator in HNSCC.
Introduction
Head and neck cancer was one of the most common malignant upper respiratory tract tumors globally, of which squamous cell carcinoma (HNSCC) accounts for about 90% (Bray et al., 2018). Although the continuous improvement of treatment such as surgery, radiotherapy, chemotherapy, targeted therapy and immunotherapy, the 5-year survival rate of HNSCC is still stagnant (Cancer Genome Atlas Network, 2015). Senile patients are mainly those who are habitually ingesting alcohol and tobacco, while HPV infection has led to an increase in young patients recently, especially oropharyngeal squamous cell carcinoma (Chaturvedi et al., 2011). Among them, patients with HPV-negative suffer from irreversible damage due to long-term smoking and drinking, who are more tolerant of the treatments with a worse prognosis.
Hox genes is a subgroup of the homeobox transcription factor family, which aggregates in 4 types of chromosomal loci (A-D), with 39 genes in humans. Hox gens plays its main function by encoding transcription factors or directly targeting downstream genes as monomers or homodimers (Alexander et al., 2009). Members of Hox genes can participate in cell development, organ formation and even carcinogenesis by regulating receptor signaling, cell differentiation, apoptosis, migration, epithelial-mesenchymal transition (EMT) and angiogenesis (Cho et al., 2012). Hox genes is tissue-specific and may play a key role in developmental regulation in one tissue type (Economides and Capecchi, 2003), while in another tissue, it is an accomplice of tumor malignant biological behavior (Wang et al., 2007). Previous studies have found that many Hox genes are abnormally expressed in hematological malignancies and various solid tumors. There are inadequate studies on the involvement of Hox genes in the occurrence and development of HNSCC. Therefore, studying the role of Hox genes in HNSCC is of great significance.
In this study, we used the method of bioinformatics mining database of TCGA-HNSC combined with immunohistochemical (IHC) verification to find a novel prognostic biomarker, HOXA1 (Homeobox A1) in Hox genes. HOXA1 is the first Hox gene to be discovered that can be upregulated during carcinogenesis by gastrin (Hofsli et al., 2005). In breast cancer, HOXA1 has been identified as a breast epithelial oncogene whose forced expression is sufficient to induce the transformation of immortalized human mammalian epithelial cells into aggressive cancer cells (Zhang et al., 2006; Mohankumar et al., 2008). In previous reports, the abnormal expression of HOXA1 is related to the prognosis of breast cancer, lung cancer, liver cancer, prostate cancer, oral cancer, gastric cancer and melanoma (Shah and Sukumar, 2010). We evaluated the abnormal expression of HOXA1 in HNSCC and confirmed its independent predictive power. At the same time, we analyzed the correlation between the expression level and clinicopathological characteristics. According to the expression level, we ran gene set variation analysis (GSVA) and gene set enrichment analysis (GSEA) to further analyze the biological pathways of HNSCC pathogenesis related to the expression of HOXA1. Finally, results of TIMER, CIBERSORT and MEXPRESS analysis indicated that HOXA1 was correlated with CD8+ T cell infiltration and DNA methylation level.
Materials and Methods
HOX Genes Expression Datasets and Relevant Clinicopathological Features
HNSCC dataset (528 cases) contain RNA-Seq (HTSeq) and clinicopathological characteristics from the cancer genome atlas (TCGA-HNSC) were downloaded (https://cancergenome.nih.gov/). Screened out 499 cases, including the expression of 39 Hox genes and corresponding clinicopathological features are used in the subsequent analysis. R Package “survminer” were used to calculate the cutoff value of Hox genes. The 499 cases were divided into a high-expression group (n = 396) and a low-expression group (n = 103) according to the cutoff value of HOXA1 expression (overall survival, OS: 1.0525, disease-free survival, DFS: 2.9749). Clinical features of HNSCC patients had been shown in Table 1.
Table 1
| Clinical characteristics | Total (n = 499) | Expression of HOXA1 | p-Value* | |
|---|---|---|---|---|
| High | Low | |||
| (n = 396) | (n = 103) | |||
| Age, n (%) | 0.8741 | |||
| <65 | 311 (62.3) | 248 (62.6) | 63 (61.2) | |
| ≥65 | 188 (37.7) | 148 (37.4) | 40 (38.8) | |
| Gender, n (%) | 0.0086 | |||
| Female | 133 (26.7) | 110 (27.8) | 23 (22.3) | |
| Male | 366 (73.3) | 286 (72.2) | 80 (77.7) | |
| Histologic grade, n (%) | 0.0023 | |||
| G1-2 | 362 (72.5) | 300 (75.8) | 62 (60.2) | |
| G3-4 | 121 (24.2) | 84 (21.2) | 37 (35.9) | |
| NA | 16 (3.2) | 12 (3.0) | 4 (3.9) | |
| AJCC Stage, n (%) | 1.0000 | |||
| I-II | 114 (22.8) | 91 (23.0) | 23 (22.3) | |
| III-IV | 371 (74.3) | 296 (74.7) | 75 (72.8) | |
| NA | 14 (2.8) | 9 (2.3) | 5 (4.9) | |
| T stage, n (%) | 0.0080 | |||
| T1-2 | 176 (35.3) | 129 (32.6) | 47 (45.6) | |
| T3-4 | 308 (61.7) | 258 (65.2) | 50 (48.5) | |
| NA | 15 (3.0) | 9 (2.3) | 6 (5.8) | |
| N classification, n (%) | 0.1761 | |||
| N0 (negative) | 238 (47.7) | 197 (49.7) | 41 (39.8) | |
| N1-3(positive) | 239 (47.9) | 185 (46.7) | 54 (52.4) | |
| NA | 22 (4.4) | 14 (3.5) | 8 (7.8) | |
| M classification, n (%) | 1.0000 | |||
| M0 | 469 (37.3) | 374 (94.4) | 95 (92.2) | |
| M1 | 5 (0.2) | 4 (1.0) | 1 (1.0) | |
| NA | 25 (50.3) | 18 (4.5) | 7 (6.8) | |
| HPV_STATUS_P16, n (%) | 1.0000 | |||
| Positive | 34 (6.8) | 22 (5.6) | 12 (11.7) | |
| Negative | 78 (15.6) | 61 (15.4) | 17 (16.5) | |
| NA | 387 (75.56) | 313 (79.0) | 74 (71.8) | |
| ANGIOLYMPHATIC_ INVASION, n (%) | 0.5906 | |||
| YES | 218 (43.7) | 92 (23.2) | 28 (27.2) | |
| NO | 120 (24.0) | 174 (43.9) | 44 (42.7) | |
| NA | 161 (32.3) | 130 (32.8) | 31 (30.1) | |
| PERINEURAL_INVASION, n (%) | 0.00197 | |||
| YES | 185 (37.1) | 142 (35.9) | 23 (22.3) | |
| NO | 165 (33.1) | 133 (33.6) | 52 (50.5) | |
| NA | 149 (29.8) | 121 (30.6) | 28 (27.2) | |
| Smoking pack / years, n (%) | 0.5314 | |||
| ≥40 | 161 (32.3) | 127 (32.1) | 34 (33.0) | |
| <40 | 126 (25.3) | 104 (26.3) | 22 (21.4) | |
| NA | 212 (42.4) | 165 (41.7) | 47 (45.6) | |
| ALCOHOL_HISTORY_ DOCUMENTED, n (%) | 1.0000 | |||
| YES | 331 (66.3) | 262 (66.2) | 69 (67.0) | |
| NO | 157 (31.5) | 125 (31.5) | 32 (31.1) | |
| NA | 11 (2.2) | 9 (2.3) | 2 (1.9) | |
| Primary Tumor Site, n (%) | 0.4934 | |||
| Oropharynx | 378 (75.8) | 296 (74.7) | 82 (79.6) | |
| Larynx | 111 (22.2) | 91 (23.0) | 20 (19.4) | |
| Hypopharynx | 10 (2.0) | 9 (2.3) | 1 (1.0) | |
Clinical features of TCGA-HNSC patients (n = 499).
The chi-square test. Bold values indicate statistically significant, P < 0.05.
Differentially Expressed Genes and Clinical Pathological Parameters Analysis
First, the R package “limma” was used to analyze the differential expression of 39 Hox genes between tumor and normal samples in 499 patients and differentially expressed genes (DEGs) with |log2 fold change (FC)| > 1 and adjusted P-values < 0.05 were selected for subsequent analysis. The prognostic genes with the predictive ability of OS and DFS were screened out by the analysis of the Kaplan–Meier survival curve and the log-rank test (P < 0.05). The correlation analysis between the expression of HOXA1 and clinicopathological characteristics were evaluated by the R package “combn.” Furthermore, Cox univariate and multivariate regression analysis were performed on the individual clinical variables [age (years ≥65/ <65), gender (male/female), grade (G1&2/G3&4), AJCC Stage (I&II/III&IV), T stage (T1–2/T3–4), N classification (N0/N1–3), M classification (M0/M1), HPV_STATUS_P16 (positive/negative), angiolymphatic invasion (Yes/No), perineural invasion (Yes/No), smoking (packs/year ≥40/ <40), alcohol history documented (Yes/No) and primary tumor site (Oropharynx/ Larynx/Hypopharynx)]with the HOXA1 expression value and high or low-expression level. Hazard ratios (HRs) and 95% confidence intervals (CIs) were calculated using the R package “coxph.” The chi-square test was used to assess the statistical significance of expression value between the two groups.
Analysis of Immune Infiltration and DNA Methylation
The fraction of twenty-two immune cell types of each TCGA-HNSC sample was yielded through R package “CIBERSORT” by estimating relevant subsets of RNA profiles (https://cibersort.stanford.edu/) (Supplementary Material 1). Only samples with P < 0.05 were conceived to perform the further analysis of comparing differential the tumor-infiltrating immune cells (TICs) level between the high and low-HOXA1 level group. In addition, another immune infiltration data is derived from TIMER (Tumor Immune Estimation Resource, https://cistrome.shinyapps.io/timer/) and DNA methylation data came from MEXPRESS database (https://mexpress.be/).
Gene Set Variation Analysis (GSVA) and Gene Set Enrichment Analysis (GSEA)
The R package “GSVA” was used to evaluate the signaling pathways (5529 gene sets) and Gene Ontology (GO) terms (10,192 gene sets) in the TCGA-HNSC database (Hänzelmann et al., 2013). Then, the R package “limma” was used to analyze the significant differential gene sets between the two groups. Besides, 499 HNSCC cases were divided according to the median expression value of HOXA1. GSEA was then performed to detect the gene sets that were enriched in the gene rank for identifying potential the pathways most associated with HOXA1. The annotated gene sets of c2.cp.kegg. v7.1. symbols in the Molecular Signatures Database (MSigDB) was selected in GSEA version 4.1.0 (https://www.gsea-msigdb.org/gsea/msigdb/index.jsp). We performed 1,000 times of permutations and filtered the results by normalized enrichment scores (NES), nominal p-value and false discovery rate (FDR) q-value (Subramanian et al., 2005).
Tissue Specimen Collection and Immunohistochemistry (IHC)
The study was approved by the Clinical Research Ethics Committee of The First Affiliated Hospital of Bengbu Medical College. Written informed consent was obtained from all-volunteer patients. In all, 52 pairs of HNSCC paraffin sections were collected from The First Affiliated Hospital of Bengbu Medical College (Supplementary Material 2). In the progress of antigen retrieval, the sections were placed in 10 mM citrate buffer solution (pH = 6.0) and boiling in a pressure cooker for 15 min. After cooling to 37°C, 3% H2O2 was added and incubated for 10 min to deplete the activity of endogenous peroxidase. After blocking by Goat serum, the sections were stained with an anti-HOXA1 polyclonal antibody (1:100, Thermo Fisher Scientific, catalog # PA5-36164, RRID AB_2553401) for a night at 4°C-refrigerator. After soaking and rinsing thoroughly with PBS three times, the sections were stained with secondary antibodies. Finally, the positive expression of the HOXA1 protein was observed by a microscope after adding substrate and hematoxylin staining. The expression density was scored by two experienced pathologists, following the single-blind trial technique. HOXA1 protein staining was assessed by the sum of the percentage of positive nuclear cells: − (0–25%), + (26–50%), ++ (51–100%). At last, we evaluated HOXA1 expression into two levels, low (− to +) and high (++). Fifty-two pairs of IHC pictures were taken for each slide, and images (200×) were analyzed using software ImageJ (US National Institute of Health, Bethesda, Maryland, USA, https://imagej.nih.gov/).
Statistical Analysis
Statistical analyses were conducted with software R (v3.6.0) and GraphPad Prism 7 (GraphPad Software, California, USA). Differences among variables were analyzed with two-tailed Student's t-test and chi-square test. All significance tests were carried out at the 0.05 level.
Results
Differentially Expression and Prognostic Value of HOX Genes in HNSCC
The analysis flow chart of this article is shown in Figure 1. Most Hox genes (30/39, Supplementary Material 3) have significant differences in expression between 499 tumor samples and 43 tumor-adjacent normal samples (Figure 2A). Among DEGs, only HOXA1 have predictive ability both in OS and DFS (Supplementary Figure 1), and worse clinical outcomes in patients with high expression of HOXA1 (Figures 2B–E).
Figure 1
Figure 2
HOXA1 Is Upregulated in HNSCC
In TCGA-HNSC, HOXA1 expression in tumor tissues was significantly higher than adjacent tissues (p = 1.241e-10) (Figure 3A), the differential analysis of HOXA1 expression in 43 paired samples got the same result (p = 9.984e−09) (Figure 3B). We performed IHC in 52 pairs of HNSCC paraffin sections for further verification (Figures 3E–G). HOXA1 protein immunostaining was mainly found in the cell nucleus but widely distributed in cancer nest. Among them, HOXA1 protein was detected in 45 cases (86.5%, + to ++) of tumor samples, while HOXA1 protein was detected in only ten sections (19.2%, + to ++) of adjacent pieces. There were 38 cases with low positive staining (− to +) and 14 patients with high positive staining (++) among the tumor samples. Our results also showed that HNSCC samples have a higher mean percentage of HOXA1-positive area (Figure 3C) and HOXA1-positive cells (Figure 3D) than paracancerous samples (P < 0.0001). The numeric dimension of the areas (mean ± SD, median) used for the IHC evaluation were shown in Table 2.
Figure 3
Table 2
| The numeric dimension | Adjacent tissues (mean ± SD, median) | HNSCC tissues (mean ± SD, median) | P-value |
|---|---|---|---|
| HOXA1-positive area | 2.86 ± 1.29, 2.94 | 10.44 ± 3.24, 10.48 | <0.0001 |
| HOXA1-positive cells | 16.69 ± 7.60, 15.50 | 40.87 ± 11.62, 42.00 | <0.0001 |
The numeric dimension of the areas (mean ± SD, median) used for the IHC evaluation.
Association of HOXA1 Expression With Clinicopathological Features
The expression of HOXA1 was highly statistically significant correlated with pathology grade (G1&2 vs. G3&4, p = 0.0077), local invasion (T1–2 vs. 3–4, p = 0.023), perineural invasion (No vs. Yes, p = 0.0019) and HPV infection (positive vs. negative, p = 0.001) (Figures 4A–L). The expression of HOXA1 in oropharyngeal cancer was significantly lower than that in laryngeal cancer (p = 0.021), and there was no significant difference between the other two tumor locations (Figure 5A). We conducted further survival analysis in three tumor locations and found that the expression of HOXA1 only has prognostic value in oropharyngeal cancer (P < 0.05), while no significant predictive value in hypopharyngeal and laryngeal cancer (p > 0.05) in TCGA-HNSC patients (Figures 5B–D).
Figure 4
Figure 5
Independent Predictive Power of HOXA1 Expression Value and Expression Level
In Cox univariate regression analysis, the HRs of high-expression vs. low-expression of HOXA1 and expression value of HOXA1 for overall survival (OS) was 2.316 (P < 0.001, CI: 1.515–3.540) and 1.236 (P < 0.001, CI: 1.115–1.369), respectively (Figure 6A). In further multivariate Cox regression analysis, the expression level and expression value of HOXA1 still have a good predictive ability for OS [expression level (HR = 1.927, p = 0.010, CI: 1.173–3.166) and expression value (HR = 1.144, p = 0.047, CI: 1.002–1.307)] (Figure 6B). The results of the above analysis indicated HOXA1 was an independent prognostic factor, and the expression level of HOXA1 seem better.
Figure 6
GSVA and GSEA Reveal HOXA1-Related Pathways and Molecular Functions
To explore the potential carcinogenic function of HOXA1, we used both GSVA and GSEA both approaches. After scoring 499 patients by the GSVA, the significant differences in signaling pathways and GO functions associated with the expression level of HOXA1 were obtained, as shown in Figures 7A,B. Besides, we ran GSEA to figure out the significant differential signaling pathways related to the expression value of HOXA1 (Table 3), including neuroprotein secretion and transport, the tumor-associated signaling pathway, cell adhere junction and metabolic reprogramming (Figures 8A–D). After the intersection of the results of GSVA and GSEA, 14 mutual pathways mainly focusing on the following four aspects: neurotrophin, vesicular transport, protein export; Ubiquitin mediated proteolysis, Glycan biosynthesis; gap junction, adherens junction, focal adhesion, regulation of actin cytoskeleton; Pathways in cancer including WNT/ERBB/ NOTCH/TGF-beta.
Figure 7
Table 3
| MSigDB collection | Gene set name | NES | NOM p-val | FDR q-val |
|---|---|---|---|---|
| c2.cp.kegg.v7.1.symbols.gmt | KEGG_WNT_SIGNALING_PATHWAY | 1.961 | 0 | 0.075 |
| KEGG_GLYCOSYLPHOSPHATIDYLINOSITOL_GPI_ANCHOR_BIOSYNTHESIS | 1.952 | 0.002 | 0.054 | |
| KEGG_ERBB_SIGNALING_PATHWAY | 1.760 | 0.003 | 0.137 | |
| KEGG_UBIQUITIN_MEDIATED_PROTEOLYSIS | 1.829 | 0.004 | 0.100 | |
| KEGG_N_GLYCAN_BIOSYNTHESIS | 2.057 | 0.004 | 0.038 | |
| KEGG_ADHERENS_JUNCTION | 1.875 | 0.005 | 0.079 | |
| KEGG_PURINE_METABOLISM | 1.679 | 0.006 | 0.141 | |
| KEGG_OTHER_GLYCAN_DEGRADATION | 1.758 | 0.008 | 0.125 | |
| KEGG_REGULATION_OF_ACTIN_CYTOSKELETON | 1.712 | 0.015 | 0.135 | |
| KEGG_PROTEIN_EXPORT | 1.720 | 0.019 | 0.138 | |
| KEGG_NOTCH_SIGNALING_PATHWAY | 1.664 | 0.026 | 0.132 | |
| KEGG_GAP_JUNCTION | 1.614 | 0.028 | 0.155 | |
| KEGG_TGF_BETA_SIGNALING_PATHWAY | 1.674 | 0.029 | 0.137 | |
| KEGG_FOCAL_ADHESION | 1.696 | 0.031 | 0.132 | |
| KEGG_INOSITOL_PHOSPHATE_METABOLISM | 1.629 | 0.034 | 0.161 | |
| KEGG_PATHWAYS_IN_CANCER | 1.591 | 0.038 | 0.146 | |
| KEGG_NEUROTROPHIN_SIGNALING_PATHWAY | 1.623 | 0.039 | 0.160 | |
| KEGG_SNARE_INTERACTIONS_IN_VESICULAR_TRANSPORT | 1.588 | 0.045 | 0.142 |
Gene sets enriched related with expression of HOXA1.
Gene sets with NOM P-value < 0.05 and FDR q-value < 0.25 are considered as significant.
MsigDB, Molecular Signatures Database; NES, normalized enrichment score; NOM, nominal; FDR, false discovery rate.
Figure 8
Immune Infiltration and DNA Methylation Correlating With HOXA1 Expression
More and more studies have found that some oncogenes significantly affect the tumor immune microenvironment, including infiltration of B cell, CD8+ T cell, CD4+ T cell, Macrophage, Neutrophil, Dendritic cell, etc. The HOXA1 expression in the TIMER database is associated with a decrease of B cell / CD8+T cell infiltration and an increase of CD4+ T cell infiltration in total HNSCC patients, a reduction of B cell and CD8+T cell infiltration in HPV(+) HNSCC patients, and a decrease of CD8+T cell infiltration and an increase of CD4+ T cell infiltration in HPV(−) HNSCC patients (Figures 9A–C). From the results of immune infiltration in total and HPV subtype HNSCC patients, it is found that the high expression of HOXA1 significantly leads to the reduction of CD8+ T cell infiltration. Besides, we performed the CIBERSORT analysis to screen the differential TICs of the high and low-HOXA1 expression groups and found that the results were similar to the TIMER analysis results. The CIBERSORT analysis results showed that the proportion of M0 macrophages was positively correlated with HOXA1 expression level, and the ratio of CD8+T cells, naïve B cells, CD4 memory activated T cells and follicular helper T cells were negatively associated with the level of HOXA1 expression (Figure 9D).
Figure 9
The correlation analysis between DNA methylation level and HOXA1 expression showed that high HOXA1 expression with less promoter methylation in tumor samples (Pearson correlation coefficients ranging from −0.166 to −0.528 for promoter region probes) in the MEXPRESS database (Figure 10). We also found that the probes such as cg03116258, cg07450037 and cd12686016 of the seventh chromosome had the strongest correlation in DNA methylation of HOXA1, which was worthy of further exploration and research.
Figure 10
Discussion
There are only sporadic reports in the literature on the role of HOX genes in HNSCC. Jinyun Li et al. used bioinformatics to analyze the HOXA genes' role in laryngeal squamous cell carcinoma (LSCC), in which HOXA10, HOXA11 and HOXA13 expression and clinical characteristics (M stage and gender) are independent prognostic biomarkers for the prognosis of LSCC patients (Li et al., 2020). The expression of HOXA members in LSCC is negatively correlated with the level of methylation. HOXA7, HOXA10 and HOXA13, which are highly expressed in LSCC tissues, have a weak positive correlation with tumor purity. However, there is no significant correlation between the TICs and HOXA members. Studies in oral squamous cell carcinoma (OSCC) showed that HOXA genes are more enriched than normal mucosal tissues. High expression of HOXA1 can promote OSCC proliferation and maybe an independent prognostic indicator for patients (Bitu et al., 2012). Both AMOT-p80 and miR-146a, identified as HOXD10 targets in HNSCC, can represent therapeutic targets for specific tumor stages (Hakami et al., 2014). High expression of HOXD10 promotes tumor proliferation and migration, while low expression of HOXD10 can promote tumor invasion and metastasis by regulating AMOT-p80. In this paper, we showed the expression profile of 39 members of the HOX genes in HNSCC, as well as the most important independent prognostic indicator and potential biological role of HOXA1 in HNSCC, providing novel insights for further research on HOXA1 as a potential target of HNSCC.
The expression of HOXA1 is tissue-specific (Fagerberg et al., 2014) (https://www.proteinatlas.org), and as tumor modulators, it may cover a wide range of functionality in types of tumor dependent on the cellular context (Shah and Sukumar, 2010). HOXA1 played a role of an oncogene in breast carcinoma (Liu et al., 2019), lung carcinoma (Abe et al., 2006), oral squamous cell carcinoma (Bitu et al., 2012), hepatocellular carcinoma (Zha et al., 2012), melanoma (Wardwell-Ozgo et al., 2014), gastric cancer (Yuan et al., 2016), leukemia (Chen et al., 2019) and retinoblastoma (RB) (Lyv et al., 2020). Inversely, it is shown as anti-oncogenes in pancreatic carcinoma (Ohuchida et al., 2012) and small cell lung cancer (SCLC) (Xiao et al., 2014). Nilay et al. pointed out that the abnormal expression of Hox genes triggers oncogenesis or tumor suppression, which may be caused by her mechanism of action (Shah and Sukumar, 2010). On the level of expression, Hox genes are associated with oncogenesis when expressed at an increased level in tumor tissues not normally seen in that tissue type. In the results of our HNSCC IHC, the expression of HOXA1, which is mainly located in the cell nucleus, is extremely low in adjacent normal tissues while it is widely distributed in cancer nest. Thus, it can be seen that the abnormal up-regulation of HOXA1 in HNSCC, especially in the nucleus, indicates that HOXA1 mainly plays a role in oncogenesis.
E-cadherin activated Rac signaling pathway through the transcription target of HOXA1 (Zhang et al., 2006), which can promote breast cancer cell survival by upregulating cyclinD1, c-Myc, and BCL-2 (Mohankumar et al., 2007, 2008). HOXA1 also takes effect independent of transcriptional regulation (Taminiau et al., 2016). In proteomics research, it is found that HOXA1 also participates in cell signal transduction as a functional protein, such as pathways of focal adhesion, receptor tyrosine kinase function, TGF/BMP and TNF/NF-kB (Lambert et al., 2012). In breast cancer, HOXA1 transforms immortalized human breast epithelial cells stimulated by p44/42 MAPK pathway, signal transducer and activator of transcription 3 (STAT3) and STAT5B into an invasive malignant phenotype. WNT, TGF-beta and FGF pathways may synergistically affect HOXA1 expression, and WNT7A directly maintains expression of HOXA1, making for lung cancer recurrence (Calvo et al., 2000). HOXA1 dedifferentiates melanoma cells into a highly aggressive cell state accompanied by TGF activation (Wardwell-Ozgo et al., 2014). Down-regulation of HOXA1 can increase the expression of E-cadherin while reducing the expression of Snail and MMP-3, which can inhibit ERK1/2 and AKT to suppress the growth, migration, invasion and metastasis of prostate cancer (Wang et al., 2015).
The methylation of DNA is inversely proportional to transcription activity, which also affects the genome's stability. Yeon et al. evaluated the methylation status of 15 genes of promoter CpG island, including HOXA1, that are involved in breast cancer progression (Park et al., 2012). The expression level of HOXA1 in subtypes of breast cancers are related to the disordered epigenetic activity, in which the methylation frequency of Luminal B is the highest (< 0.001). Meanwhile, HOXA1 hypermethylation is related to the absence of HER2 neu expression and proliferation of breast cancer cells (P < 0.05) (Pilato et al., 2013). Through the MEXPRESS database, we found that more HOXA1 expression in HNSCC samples with less promoter methylation, which possibly foreshadowing the lower DNA methylation level of oncogenes, may lead to the occurrence of tumors.
In addition to DNA methylation, many non-coding RNAs (ncRNAs) as another way of apparent epigenesis have been found to participate in the control of gene expression and activation of pathways, which in turn participates in the development of cancer and other diseases. Many studies have shown that circRNAs, lncRNAs and miRNAs play a key role in tumorigenesis, invasion, metastasis and drug resistance through the function of competitive RNA (ceRNA). Here we summarize Non-coding genes act as regulators of HOXA1 in cancers (Table 4).
Table 4
| Non-coding form | Non-coding RNAs | Cancer types | Functions | Ref. (PMID) |
|---|---|---|---|---|
| microRNA | miR-193a-5p | Breast cancer | Tumor promotor | 32497022 |
| miR-100 | Lung cancer | Chemotherapy resistance | 24559685 | |
| Nasopharyngeal Carcinoma | Growth and proliferation | 32021301 | ||
| miR-218 | Lung cancer | Gefitinib Resistance | 32084702 | |
| Cervical cancer | Proliferation, colony formation, migration and invasion | 30896864 | ||
| miR-181b | SCLC | Regulatory of p53 signaling pathway | 29658571 | |
| miR-181a-5p | MM | Apoptosis resistance | 29228867 | |
| miR-216b-5p | Cervical cancer | Proliferation | 31114990 | |
| miR-210 | SCLC/melanoma | Immune escape | 22962263 | |
| miR-145 | Oral squamous cell carcinoma | Cancer inhibition | 31138758 | |
| miR-10a | PC | Tumor promotor | 22407312 | |
| miR-30c | Malignant giant-cell tumor of bone | Metastasis and growth | 29164581 | |
| miR-30b | Esophageal cancer | Growth, migration and invasion | 28189678 | |
| miR-99a | HCC | Migration and invasion | 31186723 | |
| miR-99 family | Epithelial cell | Proliferation and migration | 24312487 | |
| miR-515 | RB | Tumor promotor | 32561925 | |
| LncRNA | LncRNA SNHG1 | Breast cancer | Tumor promotor | 32497022 |
| LncRNA CCAT1 | Lung cancer | Gefitinib Resistance | 32084702 | |
| MM | Apoptosis resistance | 29228867 | ||
| LINC00152 | Cervical cancer | Proliferation | 31114990 | |
| lncRNA HOTAIR | SCLC | Multidrug resistance | 26707824 | |
| LncRNA HOTAIRM1 | Breast cancer | Avoiding silence of HOXA1 | 32284737 | |
| Endometrial Cancer | Proliferation, Migration, Invasion and EMT | 31853186 | ||
| GBM | Growth and invasion | 30376874 | ||
| LncRNA ZFPM2-AS1 | RB | Tumor promotor | 32561925 | |
| CircRNA | CircEIF4G2 | Cervical cancer | Proliferation, colony formation, migration and invasion | 30896864 |
Non-coding genes act as regulators of HOXA1 in cancers.
SCLC, small cell lung cancer; RB, retinoblastoma; GBM, glioblastoma multiforme; PC, Pancreatic cancer; MM, Multiple Myeloma; HCC, hepatocellular carcinoma cells; EMT, epithelial–mesenchymal transition.
Tumor immunity research is in full swing, and we are still full of unknowns on how the abnormal expression of HOXA1 affects immune function. In the study of endometrial cancer, high expression of HOXA1 in Myeloid-derived suppressor cells (MDSCs) promotes tumor progression and hinders anti-tumor immune response (Li et al., 2019). Tumors adapt to hypoxic stress and create an immunosuppressive microenvironment suitable for survival (Keith et al., 2011). Studies have found that in hypoxic areas of SCLC and melanoma, miR-210 regulates the expression of HOXA1 and inhibits the anti-tumor effect of cytotoxic T lymphocyte (Noman et al., 2012). Regarding the results of TIMER and CIBERSORT analysis, we suspected that the high expression of HOXA1 in HNSCC might reduce CD8+ T cell infiltration, leading to an immunosuppressive environment not conducive to patient's survival outcomes.
Previous studies of pancreatic cancer and glioma showed that the close interaction between tumors and neurons is not conducive to the clinical prognosis. The density of newborn nerve cells in the tumor microenvironment is often associated with tumor growth and highly invasive biological behaviors (Magnon et al., 2013; Venkataramani et al., 2019). The latest study by Moran Amit et al. found that the loss of TP53 in HNSCC can promote tumor growth by affecting neuronal remodeling (Amit et al., 2020). So far, we have found that the function of p53 in SCLC and osteosarcoma is related to the regulation of HOXA1 (Zhao et al., 2015; Zhang et al., 2018). In mouse development studies, it was found that at embryonic stage (D7-9), HOXA1 expression precisely regulates the development of neural crest precursor cells, which are essential for tumor microenvironment remodeling (Gavalas et al., 2001). In this paper, we found that high HOXA1 expression is related to PNI by mining TCGA-HNSC data. According to the analysis of GSVA and GESA, we speculated that HNSCC might act on neurons by secreting exosomes containing neurotrophic proteins. However, this result lacked evidence concerning that no cases that meet PNI were found after careful screening of the paraffin blocks by pathologists. Whether HNSCC hijacks surrounding neurons to promote tumor growth through the regulation of HOXs is still worth exploring.
In summary, HOXA1 may mainly exert tumor regulation functions through the following four mechanisms: spatiotemporal expression disorder (high expression in cancer tissues, especially in the cell nucleus), abnormal protein-signaling pathways (neuronal protein secretion), mediation of immunosuppression (inhibition of CD8+ T cell infiltration) and aberrant epigenetic activity (high expression of HOXA1 leads to a decrease in DNA methylation level). Our study demonstrated that HOXA1 is highly expressed in HNSCC, which may be an independent prognostic factor in TCGA-HNSC patients, especially for oropharyngeal cancer patients. At present, HOXA1 small interfering RNA (siRNA) nanoparticles have achieved initial success in mouse models of breast cancer and are effective (75%) to reduce the incidence of tumors (Brock et al., 2014). More random cohort studies and mechanism experiments will help confirm the regulatory effects of HOXA1 on HNSCC.
Statements
Data availability statement
Publicly available datasets were analyzed in this study. This data can be found here: https://portal.gdc.cancer.gov/repository.
Ethics statement
The studies involving human participants were reviewed and approved by Clinical Research Ethics Committee of The First Affiliated Hospital of Bengbu Medical College. The patients/participants provided their written informed consent to participate in this study.
Author contributions
The project was proposed by HL and SM. HL provided paraffin tissue sections and wrote the manuscript. MW assisted in editing the manuscript. JZ downloaded and mined the relevant data of TCGA. XW and MZ did immunohistochemistry and assisted in organizing data. All authors contributed to the article and approved the submitted version.
Funding
This work was funded by the Anhui Provincial Key Research and Development Project (202004j07020007), Anhui Provincial Natural Science Foundation (1608085MH236), Key Project of Natural Science of Bengbu Medical College (BYKY2019086ZD), and Science Research Project of Bengbu Medical College (BYKF1855).
Acknowledgments
Thanks to Zhaoyi Lu (Department of Otolaryngology-Head and Neck Surgery, Xiangya Hospital, Central South University) for suggestions on the manuscript.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2020.602068/full#supplementary-material
Supplementary Figure 1Kaplan-Meier survival analysis of (A) overall survival and (B) disease free survival) according to the expression of all HOX genes in TCGA-HNSC patients. The log-rank test was used to calculate the p-value.
Supplementary Material 1CIBERSORT results of the fraction of twenty-two immune cell types of each TCGA-HNSC sample.
Supplementary Material 2Clinicopathological characteristics of 52 HNSCC patients collected from The First Affiliated Hospital of Bengbu Medical College.
Supplementary Material 3Results of differential expression analysis in all HOX genes.
References
1
AbeM.HamadaJ.-I.TakahashiO.TakahashiY.TadaM.MiyamotoM.et al. (2006). Disordered expression of HOX genes in human non-small cell lung cancer. Oncol. Rep.15, 797–802. 10.3892/or.15.4.797
2
AlexanderT.NolteC.KrumlaufR. (2009). Hox genes and segmentation of the hindbrain and axial skeleton. Annu. Rev. Cell Dev. Biol.25, 431–456. 10.1146/annurev.cellbio.042308.113423
3
AmitM.TakahashiH.DragomirM. P.LindemannA.Gleber-NettoF. O.PickeringC. R.et al. (2020). Loss of p53 drives neuron reprogramming in head and neck cancer. Nature578, 449–454. 10.1038/s41586-020-1996-3
4
BituC. C.Destro MFdS. S.CarreraM.da SilvaS. D.GranerE.KowalskiL. P.et al. (2012). HOXA1 is overexpressed in oral squamous cell carcinomas and its expression is correlated with poor prognosis. BMC Cancer12:146. 10.1186/1471-2407-12-146
5
BrayF.FerlayJ.SoerjomataramI.SiegelR. L.TorreL. A.JemalA. (2018). Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin.68, 394–424. 10.3322/caac.21492
6
BrockA.KrauseS.LiH.KowalskiM.GoldbergM. S.CollinsJ. J.et al. (2014). Silencing HoxA1 by intraductal injection of siRNA lipidoid nanoparticles prevents mammary tumor progression in mice. Sci. Transl. Med.6:217ra2. 10.1126/scitranslmed.3007048
7
CalvoR.WestJ.FranklinW.EricksonP.BemisL.LiE.et al. (2000). Altered HOX and WNT7A expression in human lung cancer. Proc. Natl. Acad. Sci. U.S.A.97, 12776–12781. 10.1073/pnas.97.23.12776
8
Cancer Genome Atlas Network (2015). Comprehensive genomic characterization of head and neck squamous cell carcinomas. Nature517, 576–582. 10.1038/nature14129
9
ChaturvediA. K.EngelsE. A.PfeifferR. M.HernandezB. Y.XiaoW.KimE.et al. (2011). Human papillomavirus and rising oropharyngeal cancer incidence in the United States. J. Clin. Oncol.29, 4294–4301. 10.1200/JCO.2011.36.4596
10
ChenS.-L.QinZ.-Y.HuF.WangY.DaiY.-J.LiangY. (2019). The role of the gene family in acute myeloid leukemia. Genes10:621. 10.3390/genes10080621
11
ChoH.-S.ToyokawaG.DaigoY.HayamiS.MasudaK.IkawaN.et al. (2012). The JmjC domain-containing histone demethylase KDM3A is a positive regulator of the G1/S transition in cancer cells via transcriptional regulation of the HOXA1 gene. Int. J. Cancer131, E179–E189. 10.1002/ijc.26501
12
EconomidesK. D.CapecchiM. R. (2003). Hoxb13 is required for normal differentiation and secretory function of the ventral prostate. Development130, 2061–2069. 10.1242/dev.00432
13
FagerbergL.HallströmB. M.OksvoldP.KampfC.DjureinovicD.OdebergJ.et al. (2014). Analysis of the human tissue-specific expression by genome-wide integration of transcriptomics and antibody-based proteomics. Mol. Cell. Proteomics MCP13, 397–406. 10.1074/mcp.M113.035600
14
GavalasA.TrainorP.Ariza-McNaughtonL.KrumlaufR. (2001). Synergy between Hoxa1 and Hoxb1: the relationship between arch patterning and the generation of cranial neural crest. Development128, 3017–3027.
15
HakamiF.DardaL.StaffordP.WollP.LambertD. W.HunterK. D. (2014). The roles of HOXD10 in the development and progression of head and neck squamous cell carcinoma (HNSCC). Br. J. Cancer111, 807–816. 10.1038/bjc.2014.372
16
HänzelmannS.CasteloR.GuinneyJ. (2013). GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics14:7. 10.1186/1471-2105-14-7
17
HofsliE.ThommesenL.YadetieF.LangaasM.KusnierczykW.FalkmerU.et al. (2005). Identification of novel growth factor-responsive genes in neuroendocrine gastrointestinal tumour cells. Br. J. Cancer92, 1506–1516. 10.1038/sj.bjc.6602535
18
KeithB.JohnsonR. S.SimonM. C. (2011). HIF1α and HIF2α: sibling rivalry in hypoxic tumour growth and progression. Nat. Rev. Cancer12, 9–22. 10.1038/nrc3183
19
LambertB.VandeputteJ.RemacleS.BergiersI.SimonisN.TwizereJ.-C.et al. (2012). Protein interactions of the transcription factor Hoxa1. BMC Dev. Biol.12:29. 10.1186/1471-213X-12-29
20
LiJ.YeM.ZhouC. (2020). Expression profile and prognostic values of family members in laryngeal squamous cell cancer. Front. Oncol.10:368. 10.3389/fonc.2020.00368
21
LiX.PangL.YangZ.LiuJ.LiW.WangD. (2019). LncRNA HOTAIRM1/HOXA1 axis promotes cell proliferation, migration and invasion in endometrial cancer. Onco. Targets. Ther.12, 10997–11015. 10.2147/OTT.S222334
22
LiuJ.LiuJ.LuX. (2019). HOXA1 upregulation is associated with poor prognosis and tumor progression in breast cancer. Exp. Ther. Med.17, 1896–1902. 10.3892/etm.2018.7145
23
LyvX.WuF.ZhangH.LuJ.WangL.MaY. (2020). Long noncoding RNA ZFPM2-AS1 knockdown restrains the development of retinoblastoma by modulating the microRNA-515/HOXA1/Wnt/β-catenin axis. Invest. Ophthalmol. Vis. Sci.61:41. 10.1167/iovs.61.6.41
24
MagnonC.HallS. J.LinJ.XueX.GerberL.FreedlandS. J.et al. (2013). Autonomic nerve development contributes to prostate cancer progression. Science341:1236361. 10.1126/science.1236361
25
MohankumarK. M.PerryJ. K.KannanN.KohnoK.GluckmanP. D.EmeraldB. S.et al. (2008). Transcriptional activation of signal transducer and activator of transcription (STAT) 3 and STAT5B partially mediate homeobox A1-stimulated oncogenic transformation of the immortalized human mammary epithelial cell. Endocrinology149, 2219–2229. 10.1210/en.2007-1320
26
MohankumarK. M.XuX. Q.ZhuT.KannanN.MillerL. D.LiuE. T.et al. (2007). HOXA1-stimulated oncogenicity is mediated by selective upregulation of components of the p44/42 MAP kinase pathway in human mammary carcinoma cells. Oncogene26, 3998–4008. 10.1038/sj.onc.1210180
27
NomanM. Z.BuartS.RomeroP.KetariS.JanjiB.MariB.et al. (2012). Hypoxia-inducible miR-210 regulates the susceptibility of tumor cells to lysis by cytotoxic T cells. Cancer Res.72, 4629–4641. 10.1158/0008-5472.CAN-12-1383
28
OhuchidaK.MizumotoK.LinC.YamaguchiH.OhtsukaT.SatoN.et al. (2012). MicroRNA-10a is overexpressed in human pancreatic cancer and involved in its invasiveness partially via suppression of the HOXA1 gene. Ann. Surg. Oncol.19, 2394–2402. 10.1245/s10434-012-2252-3
29
ParkS. Y.KwonH. J.ChoiY.LeeH. E.KimS.-W.KimJ. H.et al. (2012). Distinct patterns of promoter CpG island methylation of breast cancer subtypes are associated with stem cell phenotypes. Modern Pathol.25, 185–196. 10.1038/modpathol.2011.160
30
PilatoB.PintoR.De SummaS.LamboR.ParadisoA.TommasiS. (2013). HOX gene methylation status analysis in patients with hereditary breast cancer. J. Hum. Genet.58, 51–53. 10.1038/jhg.2012.118
31
ShahN.SukumarS. (2010). The Hox genes and their roles in oncogenesis. Nat. Rev. Cancer10, 361–371. 10.1038/nrc2826
32
SubramanianA.TamayoP.MoothaV. K.MukherjeeS.EbertB. L.GilletteM. A.et al. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. U.S.A.102, 15545–15550. 10.1073/pnas.0506580102
33
TaminiauA.DraimeA.TysJ.LambertB.VandeputteJ.NguyenN.et al. (2016). HOXA1 binds RBCK1/HOIL-1 and TRAF2 and modulates the TNF/NF-κB pathway in a transcription-independent manner. Nucleic Acids Res.44, 7331–7349. 10.1093/nar/gkw606
34
VenkataramaniV.TanevD. I.StrahleC.Studier-FischerA.FankhauserL.KesslerT.et al. (2019). Glutamatergic synaptic input to glioma cells drives brain tumour progression. Nature573, 532–538. 10.1038/s41586-019-1564-x
35
WangH.LiuG.ShenD.YeH.HuangJ.JiaoL.et al. (2015). HOXA1 enhances the cell proliferation, invasion and metastasis of prostate cancer cells. Oncol. Rep.34, 1203–1210. 10.3892/or.2015.4085
36
WangZ.DahiyaS.ProvencherH.MuirB.CarneyE.CoserK.et al. (2007). The prognostic biomarkers HOXB13, IL17BR, and CHDH are regulated by estrogen in breast cancer. Clin. Cancer Res.13, 6327–6334. 10.1158/1078-0432.CCR-07-0310
37
Wardwell-OzgoJ.DogrulukT.GiffordA.ZhangY.HeffernanT. P.van DoornR.et al. (2014). HOXA1 drives melanoma tumor growth and metastasis and elicits an invasion gene expression signature that prognosticates clinical outcome. Oncogene33, 1017–1026. 10.1038/onc.2013.30
38
XiaoF.BaiY.ChenZ.LiY.LuoL.HuangJ.et al. (2014). Downregulation of HOXA1 gene affects small cell lung cancer cell survival and chemoresistance under the regulation of miR-100. Eur. J. Cancer50, 1541–1554. 10.1016/j.ejca.2014.01.024
39
YuanC.ZhuX.HanY.SongC.LiuC.LuS.et al. (2016). Elevated HOXA1 expression correlates with accelerated tumor cell proliferation and poor prognosis in gastric cancer partly via cyclin D1. J. Exp. Clin. Cancer Res.35:15. 10.1186/s13046-016-0294-2
40
ZhaT.-Z.HuB.-S.YuH.-F.TanY.-F.ZhangY.ZhangK. (2012). Overexpression of HOXA1 correlates with poor prognosis in patients with hepatocellular carcinoma. Tumour Biol.33, 2125–2134. 10.1007/s13277-012-0472-6
41
ZhangX.EmeraldB. S.MukhinaS.MohankumarK. M.KraemerA.YapA. S.et al. (2006). HOXA1 is required for E-cadherin-dependent anchorage-independent survival of human mammary carcinoma cells. J. Biol. Chem.281, 6471–6481. 10.1074/jbc.M512666200
42
ZhangY.LiX.-J.HeR.-Q.WangX.ZhangT.-T.QinY.et al. (2018). Upregulation of HOXA1 promotes tumorigenesis and development of non-small cell lung cancer: a comprehensive investigation based on reverse transcription-quantitative polymerase chain reaction and bioinformatics analysis. Int. J. Oncol.53, 73–86. 10.3892/ijo.2018.4372
43
ZhaoS.NiuF.XuC.-Y.YeL.BiG.-B.ChenL.et al. (2015). Microarray and ChIP-seq data analysis revealed changes in p53-mediated transcriptional regulation in Nutlin-3-treated U2OS cells. Mol. Med. Rep.12, 4284–4290. 10.3892/mmr.2015.3933
Summary
Keywords
hox genes, homeobox A1, head and neck squamous cell carcinoma, perineural invasion, prognosis, biomarker, tumor
Citation
Li H, Wang X, Zhang M, Wang M, Zhang J and Ma S (2021) Identification of HOXA1 as a Novel Biomarker in Prognosis of Head and Neck Squamous Cell Carcinoma. Front. Mol. Biosci. 7:602068. doi: 10.3389/fmolb.2020.602068
Received
02 September 2020
Accepted
12 November 2020
Published
08 March 2021
Volume
7 - 2020
Edited by
Satyendra Chandra Tripathi, All India Institute of Medical Sciences Nagpur, India
Reviewed by
José Manuel Lopes, Universidade Do Porto, Portugal; Ajay Dixit, University of Minnesota Twin Cities, United States; Bin Zhao, Xiamen University, China
Updates
Copyright
© 2021 Li, Wang, Zhang, Wang, Zhang and Ma.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Shiyin Ma mashiyinent@sina.com
This article was submitted to Molecular Diagnostics and Therapeutics, a section of the journal Frontiers in Molecular Biosciences
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.