We acquired the gene expression matrix of the GBM dataset from the Cancer Genome Atlas (TCGA) database (https://www.cancer.gov/) [23]. The dataset included 174 samples of GBM, with clinical information available for all samples. This comprised tumor tissues from 169 GBM patients (Tumor) and partially matched paracancerous tissues from 5 GBM patients (Normal). The count sequencing data of these samples were normalized to Fragments Per Kilobase per Million (FPKM) format. By integrating clinical data with sequencing data and excluding incomplete clinical information, this study obtained complete progression-free survival clinical data for 168 GBM patients.

In addition, the GBM dataset of TCGA combined with Genotype-Tissue Expression (GTEx) (TCGA-GTEx-GBM) was uniformly processed into the Transcripts Per Kilobase of exon model per Million mapped reads (TPM) format. This study used 1323 samples, including 1152 normal samples from the GTEx dataset, 5 paraneoplastic samples, and 166 tumor samples from the TCGA-GBM dataset. These datasets were utilized to conduct differential expression analysis of the HOX family genes and perform receiver operating characteristic (ROC) curve analysis.

To investigate how HOXD family genes work in GBM and their related biological features and pathways, we used the Mann-Whitney U test (Wilcoxon rank sum test) to compare expression differences between the Tumor and Normal groups in the GBM dataset. We visually presented the results using subgroup comparative plots. Additionally, we analyzed the prognosis of the HOXD gene family in GBM using the TCGA-GBM dataset. We presented the results using the Kaplan-Meier (KM) survival curves to show the relationship between HOXD gene expression and patient survival in GBM. To further confirm the correlation between the HOXD family gene expression in GBM and different groups within the GBM dataset, we used the Gene Expression Profiling Interactive Analysis (GEPIA) tool [24]. GEPIA is an online platform offering diverse analytical functionalities utilizing data from TCGA and GTEx. In this study, GBM data from TCGA and GTEx were utilized within the GEPIA website to investigate the distinct expression patterns of the HOX family genes in tumor and non-tumor samples, and to evaluate their association with prognosis. Overall survival and Disease specific survival were chosen as the prognostic endpoints, and the correlation between gene expression levels and prognosis was assessed through hazard ratio (HR) analysis.

We used the survival ROC package (version 1.80.0) in R software (version 4.1.2, R Foundation for Statistical Computing, Vienna, Austria) to create ROC curves that show the connection between the expressions of the HOXD family genes and the survival time and status of GBM patients. We then calculated the Area Under Curve (AUC) to evaluate how well the expression of HOXD family genes can predict the survival outcomes of GBM patients.

Initially, the genes in the TCGA-GBM dataset were divided into high and low phenotypic relevance groups based on their correlation. Following this, enrichment analysis of all differential genes from both high and low phenotypic relevance groups was carried out using the cluster Profiler package (version 4.12.6). The GSEA was performed using the following parameters (Table 1). For enrichment analysis, we used the c2.cp.v7.2.symbols gene set from the Molecular Signatures Database (MSigDB). Significant enrichment was determined based on the criteria of p 0.05 and a false discovery rate (FDR) value (q.value) 0.05.

We conducted univariate Cox regression analysis to evaluate the expression of key genes in the GBM dataset. Following this, we employed a multivariate Cox regression analysis that encompassed all the identified key genes to construct strong multivariate Cox regression models. Finally, we utilized the R package ggDCA (version 1.1) [25] to produce a decision curve analysis (DCA) plot, enabling us to assess the predictive performance of the nomogram model on survival outcomes for GBM patients.

The human microglia cell line HMC3 (AW-CNH003, Abiowell, Changsha, China), human glioma cell lines U87-MG (AW-CCH001, Abiowell, Changsha, China), U251-MG (AW-CCH034, Abiowell, Changsha, China) and LN229 (AW-CCH094, Abiowell, Changsha, China) were acquired from Changsha Abiowell Biotechnology Co., Ltd. HMC3 cells were cultured in Minimum Essential Medium (MEM) medium containing 10% fetal bovine serum (FBS; 10099141, Gibco, Grand Island, NY, USA) and 1% penicillin-streptomycin (SV30010, Beyotime, Shanghai, China). U87-MG, U251-MG, and LN229 cells were cultured in a Dulbecco’s modified Eagle’s medium (DMEM) medium containing 10% FBS and 1% penicillin-streptomycin. All cell lines were validated by short-tandem repeat (STR) profiling and tested negative for mycoplasma. Cells were all cultured in a humidified incubator at 37 °C and 5% CO₂.

Cells in logarithmic phase were transfected using Lipofectamine 2000 (11668019, Invitrogen, Carlsbad, CA, USA) as per the kit instructions after cell adhesion. The cells were transfected or co-transfected with the following: overexpressing negative control plasmid (oe-NC; LV-NC01, Honorgene, Changsha, China), HMGB1 overexpressing plasmid (oe-HMGB1; HG-HO002128, Honorgene, Changsha, China), oe-HOXD9 (HG-HO014213, Honorgene, Changsha, China), si-NC (LV-shNC, Honorgene, Changsha, China), si-HOXD9 (HG-HS014213, Honorgene, Changsha, China), oe-PFKFB3 (HG-HO004566, Honorgene, Changsha, China). The si-HOXD9 target sequence was 5^′-GACTCTGTATTTGCTCGTTTA-3^′.

The in vitro experiments in this study were divided into the following 8 groups, with U87-MG and U251-MG cells serving as controls for each other: Control, Hypoxia, Hypoxia+si-NC, Hypoxia+si-HOXD9, NC, HMGB1, HMGB1+3PO, HMGB1+si-HOXD9. In the Control group, cells were cultured normally. In the Hypoxia group, cells were cultured under hypoxic conditions (1% O₂ and 5% CO₂) at 37 °C for 48 h [21]. In the Hypoxia+si-NC group, cells were transfected with si-NC under hypoxic conditions for 48 h. In the Hypoxia+si-HOXD9 group, cells were transfected with si-HOXD9 under hypoxic conditions for 48 h. In the NC group, cells were transfected with si-NC under hypoxic conditions for 48 h. In the HMGB1 group, cells were transfected with oe-HMGB1 under hypoxic conditions for 48 h. In the HMGB1+3PO group, cells were transfected with oe-HMGB1 and treated with 10 mM PFKFB3 inhibitor under hypoxic conditions for 48 h [22]. In the HMGB1+si-HOXD9 group, cells were transfected with oe-HMGB1 and si-HOXD9 under hypoxic conditions for 48 h.

In accordance with the protocol of Total RNA Extraction (Trizol; 15596026, Thermo, Waltham, MA, USA), total RNA was isolated from tissues and cells. Subsequently, the absorbance readings of the total RNA were assessed using an ultraviolet spectrophotometer. The following procedures proceeded once the OD260/OD280 ratio fell within the range of 1.8–2.0. RT-qPCR analysis was conducted utilizing the SYBR method on a Quantstudio1 (ABI, Foster, CA, USA). The primer sequences are shown in Table 2.

The cells were washed with chilled PBS buffer and subsequently treated with radioimmunoprecipitation assay (RIPA) lysate (AWB0136, Abiowell, Changsha, China) for lysis. Following the determination of protein concentration using a bicinchoninic acid (BCA) kit (AWB0104, Abiowell, Changsha, China), the appropriate amount of protein was combined with a 5 $\times$ loading buffer (AWB0055, Abiowell, Changsha, China) for thorough mixing. The membrane transfer current was 300 mA, with transfer time being adjusted based on specific antibodies utilized. Upon completion of the membrane transfer, it was washed once in 1 $\times$ phosphate-buffered saline with Tween-20 (PBST), followed by complete immersion in a blocking solution and agitation on a shaker for 90 min. The primary antibodies were diluted in 1 $\times$ PBST at the following ratios: HOXD9 (1:1000; AWA44737, Abiowell, Changsha, China), HMGB1 (1:1000; AWA43119, Abiowell, Changsha, China), PFKFB3 (1:1000; 13763-1-AP, Proteintech, Chicago, IL, USA), and $\beta$-actin (1:5000; AWA80002, Abiowell, Changsha, China). The membrane was incubated with the primary antibodies overnight at 4 °C. Subsequently, the ECL western blotting substrate (AWB0005, Abiowell, Changsha, China) was used to treat the membrane for 1 min. Finally, the membrane was visualized using a chemiluminescence imaging system (ChemiScope6100, Clinx Science Instruments, Shanghai, China).

We used an enzyme-linked immunosorbent assay (ELISA) kit (CSB-E08223h, Cusabio, Wuhan, China) to measure HMGB1 secretion levels. In short, cells were cultured at a density of 2 $\times$ 10⁶/mL, and then the cell supernatant (100 µL) was collected. Then, 100 µL of standards and samples were added to each well and incubated at 37 °C for 2 h. After removing the liquid in the wells, 100 µL of biotin-labeled antibody working solution was added to each well, followed by incubation at 37 °C for 1 h. Subsequently, the liquid in the wells was discarded, and 100 µL of horseradish peroxidase-labeled affinity protein working solution was added and incubated at 37 °C for 1 h. After discarding the liquid, 90 µL of substrate solution was added, and the color was allowed to develop at 37 °C for 15–30 min while protecting it from light. To stop the reaction, 50 µL of termination solution was added. Within 5 min of terminating the reaction, the optical density (OD) of each well was measured at 450 nm using an enzyme marker.

Cell viability was assessed using a CCK-8 kit (NU679, Dojindo, Kumamoto, Japan). Cells in the logarithmic growth phase were seeded at 1 $\times$ 10⁴ per well in a 24-well plate with 300 µL per well. After the culture adhered and was treated as described above, 30 µL CCK-8 was added to each well for the corresponding time. Then, 300 µL CCK-8 containing medium was added. After continued incubation, the absorbance value at 450 nm was analyzed with a microplate reader (MB-530, HEALES, Shenzhen, China), and the mean value was used to create a bar graph.

First, remove the fixed overnight sample. To remove ethanol, add 1 mL of PBS to resuspend the cells. Centrifuge the suspension at 800 rpm for 5 min to collect the cells. Then, add 150 µL of propidium iodide (PI) working solution (MB2920, Meilunbio, Dalian, China) to the cells and stain them for 30 min at 4 °C while shielding them from light. Finally, transfer the cells to a flow-through detector tube and analyze them using the FACSCanto Flow Cytometer (A00-1-1102, Beckman Coulter, Indianapolis, IN, USA).

In the Transwell chamber (3428, Corning Inc., Corning, NY, USA), fill the lower layer with 500 µL of complete medium. After the cells were digested with trypsin to form individual cells, the serum-free medium was resuspended to 2 $\times$ 10⁶ cells/mL, and 100 µL of cells were added to each well. The cells were placed in the incubator at 37 °C for 48 h. After the cells were fixed with 4% paraformaldehyde, cells were stained with 0.1% crystal violet (AWC0333, Abiowell, Changsha, China) for 5 min. After staining, wash the cells and count them using flow cytometry (A00-1-1102, Beckman Coulter, Indianapolis, IN, USA).

We used the Glucose kit (A154-1-1, Njjcbio, Nanjing, China) to measure the glucose content in tissue and cell supernatants. After sample processing, 2.5 µL the sample and 250 µL of working solution were added to the 96-well plate. The absorbance value was then detected at a wavelength of 505 nm.

Additionally, we utilized the Lactic Acid assay kit (A019-2-1, Njjcbio, Nanjing, China) to detect the lactate levels in tissue and cell supernatants. The test sample of 0.02 mL, enzyme working solution of 1 mL, and 0.2 mL chromogen were mixed, followed by the addition of 2 mL termination solution. The absorbance was detected at a wavelength of 530 nm.

After fixing the slides, we permeated them with 0.3% Triton X-100 for 30 min. Subsequently, we added 5% BSA-PBS to the block for 60 min. We then added the primary antibodies HOXD9 (1:50; sc-137134, Santa, Dallas, TEX, USA) and PFKFB3 (1:50; ab181861, Abcam, Cambridge, MA, USA), and incubated them overnight at 4 °C. 4^′,6-diamidino-2-phenylindole (DAPI) staining solution (AWI0429, Abiowell, Changsha, China) was applied, and the samples were observed using fluorescence microscopy. HOXD9 exhibited red fluorescence, PFKFB3 showed green fluorescence, and DAPI displayed blue fluorescence.

In this study, we followed the ChIP Assay Kit (ab500, Abcam, Cambridge, MA, USA) protocol. We used a final concentration of 1.1% formaldehyde to crosslink the DNA and proteins in the cells. Then, we fragmented the crosslinked DNA into smaller pieces using ultrasonication. We added antibodies of interest or positive control (ab1791, Abcam, Cambridge, MA, USA) were added to the supernatant for immunoprecipitation. Finally, we used the eluted and purified DNA as a template to amplify the fragment region of the PFKFB3 promoter region prediction site by RT-qPCR with PFKFB3-specific primer set, forward primer (5^′-GTGGGGAGGCAGAGCTGCAGG-3^′) and the reverse primer (5^′-TGCTGGGTACGAGGAGG-3^′), yielding a product of 273 bp.

We obtained 40 healthy male nude mice (4 weeks of age) from Hunan Slaike Jingda Experimental Animal Co., LTD (Changsha, China). After one week of adaptive feeding, we injected U87-MG and U251-MG cells (1 $\times$ 10⁶) [26, 27] subcutaneously into the nude mice. The mice were grouped as follows: Control, si-HOXD9, si-HOXD9+oe-NC, si-HOXD9+oe-PFKPB3. The mice were intratumorally injected once every 3 days with 100 µL interfering adenovirus [28] to interfere with HOXD9 expression or overexpression of PFKFB3 protein. We measured tumor size twice a week after tumor seeding. After 28 days, we photographed and harvested the tumors. Mice were sacrificed by intraperitoneal injection of sodium pentobarbital (150 mg/kg).

For data processing and analysis, we used R software (version 4.1.2) and conducted statistical analysis using GraphPad Prism 8.0 software (GraphPad Software, Inc., San Diego, CA, USA). Wilcoxon rank-sum tests were used to determine the differential expression of nine HOXD family genes between the Tumor group and the Normal group. If the data followed a normal distribution, we performed pairwise comparisons between the two groups using the T-test. Comparisons between multiple groups were performed using the one-way or two-way ANOVA test, and then Tukey’s multiple comparisons test was used. We considered results with a p 0.05 as significantly different.

In our investigation of the possible mechanism and associated biological features and pathways of the nine HOX family genes in GBM, we analyzed their expression differences in GBM using the Mann-Whitney U test. As depicted in Fig. 1A, except for HOXD1, the expressions of the other eight HOX family genes were highly statistically significant.

Fig. 1.

Differential expression analysis of HOX family genes in the TCGA-GBM dataset. (A) The comparative grouping figure displayed the expression difference analysis results for the HOXD gene family in the TCGA-GBM dataset. (B–J) The ROC curve represented the gene expressions of HOXD1, HOXD3, HOXD4, HOXD8, HOXD9, HOXD10, HOXD11, HOXD12, and HOXD13 in the TCGA-GBM dataset. *p 0.05, ***p 0.001. HOX, homeobox; TCGA, the Cancer Genome Atlas; GBM, Glioblastoma; ROC, Receiver operating characteristic; TPR, true positive rate; FPR, false positive rate; AUC, area under the curve; CI, confidence interval.

Subsequently, the ROC curve revealed that the expression level of HOXD1 (AUC = 0.554) had a low correlation between the Tumor and Normal groups. On the other hand, the expression level of HOXD12 (AUC = 0.891) demonstrated a certain correlation in the GBM dataset. Additionally, the expression of HOXD3 (AUC = 0.906), HOXD4 (AUC = 0.925), HOXD8 (AUC = 0.961), HOXD9 (AUC = 0.978), HOXD10 (AUC = 0.995), HOXD11 (AUC = 0.971) and HOXD13 (AUC = 0.976) in TCGA-GTEx-GBM exhibited a high correlation between the Tumor and Normal groups (Fig. 1B–J). Consequently, we identified seven candidate genes by plotting the ROC curves.

We analyzed the clinical information from GBM samples and created survival curves for nine HOX family genes in the TCGA-GBM dataset (FPKM format). These curves depict overall survival (OS) and disease-specific survival (DSS) outcomes. Our findings revealed that only HOXD9 (p = 0.028), HOXD10 (p = 0.042), and HOXD12 (p = 0.048) exhibited statistically significant correlations with the overall survival of GBM patients, based on their expression levels (p 0.05). Specifically, high expression levels of HOXD9, HOXD10, and HOXD12 implied lower overall survival in GBM patients. On the other hand, HOXD1, HOXD3, HOXD4, HOXD8, HOXD11, and HOXD13 did not show statistically significant associations with overall survival (p $\ge$ 0.05) (Fig. 2A–I).

Fig. 2.

Prognostic KM curve analysis of HOX family genes. (A–I) The OS analysis curves in the TCGA-GBM dataset. (J–R) DSS analysis curves from the TCGA-GBM dataset. Each point on the curve represents the survival rate of patients at that specific time point. The “+” symbol on the curve represented censored data. HOX, homeobox; TCGA, the Cancer Genome Atlas; GBM, Glioblastoma; KM, Kaplan-Meier; OS, overall survival; DSS, disease specific survival; HR, hazard ratio.

The DSS KM curves of the HOXD gene family in the FPKM format TCGA-GBM dataset indicated that the expression levels of HOXD9 (p = 0.022) and HOXD10 (p = 0.039) were significantly associated with the disease-specific survival of GBM patients (p 0.05). However, the expression levels of HOXD1, HOXD3, HOXD4, HOXD8, HOXD11, HOXD12, and HOXD13 did not show statistically significant associations with the disease-specific survival of GBM patients (p $\ge$ 0.05) (Fig. 2J–R). These results suggest that the expression levels of HOXD9 and HOXD10 were statistically significant in both prognostic OS and DSS of GBM patients.

To delve deeper into the correlation between the expression levels of HOX family genes in GBM and various subgroups within the GBM dataset, we conducted an online analysis using the GEPIA. The TCGA-GTEx-GBM dataset on the GEPIA website includes 163 samples from the tumor group (T) and 207 samples from the normal group (N). The analysis revealed that all nine gene expressions in the tumor group were higher in the GBM dataset. However, the expression differences of HOXD1, HOXD4, and HOXD12 in GBM did not show statistical significance (p $\ge$ 0.01). On the other hand, the HOXD3, HOXD8, HOXD9, HOXD10, HOXD11, and HOXD13 exhibited statistically significant increase in expression levels in GBM (Supplementary Fig. 1A–I).

Furthermore, we performed a prognostic KM curve analysis to assess the relationship between the expression levels of the 9 HOX family genes in GBM and the OS of GBM patients. The analysis revealed that the expressions of HOXD1, HOXD3, HOXD4, HOXD8, HOXD11, and HOXD13 did not exhibit statistically significant correlations with the prognosis of GBM patients (p (HR) $\ge$ 0.05). However, the expression levels of HOXD9 (p (HR) = 0.044), HOXD10 (p (HR) = 0.03), and HOXD12 (p (HR) = 0.012) exhibited a statistically significant association with the OS of GBM patients (Supplementary Fig. 1J–R). Specifically, higher HOXD9 expression implied poorer overall survival in GBM patients.

We selected nine genes from the HOXD family based on the results of the above differential analysis and prognostic KM curve analysis. We then summarized the results of our analyses, including differential analysis (Fig. 1A) and ROC curve analysis (Fig. 1B–J) in the GBM dataset of TCGA combined with GTEx, as well as the results of differential analysis of the expression of the HOXD family genes in the GBM dataset on the GEPIA website (Supplementary Fig. 1A–I) and then plotting a Venn diagram (Fig. 3A). As depicted in Fig. 3B, the outcomes of the prognostic survival KM curves for OS (Fig. 2A–I) and DSS KM curves (Fig. 2J–R) of the nine genes of the HOXD family in the TCGA-GBM dataset, as well as the results of the prognostic KM curves analysis of the HOXD family genes in the GBM dataset from the GEPIA website (Supplementary Fig. 1J–R) were also summarized in the form of a Venn diagram for presentation. The results showed that six of the nine HOXD family genes showed statistically significant (p 0.05) differences in expression between different subgroups in GBM, and these six genes were: HOXD3, HOXD8, HOXD9, HOXD10, HOXD11, and HOXD13 (Fig. 3A). Additionally, the expression of only two of the nine HOXD family genes, HOXD9, and HOXD10, was statistically significant (p 0.05) associated with the prognostic clinical of GBM patients (Fig. 3B). We further performed Venn analysis on the six differentially expressed genes and two differentially expressed genes with clinical prognostic relevance to GBM obtained, and obtained HOXD9 and HOXD10 (Fig. 3C). That is, HOXD9 and HOXD10 are HOXD candidate genes that satisfy both differential expression and clinical prognostic relevance. That is, both HOXD9 and HOXD10 had elevated expression levels in GBM and were positively correlated with poor prognosis.

Fig. 3.

Target gene screening and single gene co-expression analysis. (A) This Venn diagram depicted the results of the differential analysis performed on the TCGA_GTEx-GBM dataset, as well as the results obtained from the ROC curve analysis and GEPIA analysis of the HOX family genes in the GBM dataset. (B) A Venn diagram was used to summarize and illustrate the prognostic survival KM curves for both OS and DSS of the 9 HOX family genes in the TCGA-GBM dataset (FPKM format). Additionally, the Venn diagram included the prognostic KM curve analysis results for OS of the HOX family genes in the GBM dataset from the GEPIA website. (C) The Venn diagram results of the difference analysis and prognostic KM curve were summarized. (D) The forest plots showed the results of both the univariate and multivariate Cox regression analyses. (E) Nomogram analysis was used to evaluate the predictive ability of the Cox regression model. (F) Prognostic Calibration analyses were performed at 1, 2, and 3 years for variables in both univariate and multivariate Cox regressions. (G) Plotting the risk factors of the Cox regression prognostic model. (H) This study utilized DCA to assess how well the nomogram model predicts survival for patients with GBM. TCGA, the Cancer Genome Atlas; GTEx, Genotype-Tissue Expression; GBM, Glioblastoma; ROC, Receiver operating characteristic; HOX, homeobox; KM, Kaplan-Meier; OS, overall survival; DSS, disease specific survival; FPKM, Fragments Per Kilobase per Million; DCA, decision curve analysis; HR, hazard ratio; CI, confidence interval.

To further validate the association between HOXD gene expression levels and GBM occurrence, univariate and multivariate Cox regression analyses were conducted using the TCGA-GBM dataset. The results indicated a significant correlation between the expression levels of HOXD9 and HOXD10 and the clinical prognosis of GBM patients (Fig. 3D). The clinical information of GBM patients is displayed in Supplementary Table 1. To assess the prognostic potential of the Cox regression model, we constructed a prognostic nomogram to predict the 1-, 2-, and 3-year overall survival of GBM patients (Fig. 3E). HOXD9 was the largest contributing factor and the most significant predictor of GBM. Furthermore, prognostic calibration analyses were performed for 1, 2, and 3-year time points, and calibration curves were plotted for the variables in both univariate and multivariate Cox regressions (Fig. 3F). The calibration curve demonstrated that the prognostic nomogram was reliable and accurate. At the same time, a risk factor diagram was used to visualize the grouping of risk factors in the Cox regression model (Fig. 3G,H).

Taken together, the significance of HOXD9 and HOXD10 was relatively similar in the expression difference analysis. However, in the prognostic survival analysis, HOXD9 (p = 0.028, Fig. 2E) HOXD10 (p = 0.042, Fig. 2F) and HOXD9 (p = 0.022, Fig. 2N) HOXD10 (p = 0.039, Fig. 2O) were more prognostically relevant for HOXD9 than HOXD10. Furthermore, according to Fig. 3D,E, HOXD9 was the largest contributing factor and the most significant predictor of GBM. Based on the findings, HOXD9 was identified as a key gene for analyzing the prognosis of GBM.

PFKFB3, an isoform of fructose-6-phosphate-2 kinase (PFK2) and a regulatory glycolytic enzyme, is overexpressed in GBM, and selective inhibition of these enzymes has emerged as a new avenue for tumor therapy [17]. HOXD9 is known to bind directly to the promoter region of PFKFB3, thereby enhancing its transcription [18]. Hypoxia also promotes HMGB1 expression and secretion in glioma stem cells [21]. The lactate dehydrogenase inhibitor oxamic acid (OA) and the PFKFB3 inhibitor 3PO block HMGB1-induced fibroblast proliferation and glycolytic processes [22]. Importantly, the HOXD9-PFKFB3-glycolysis axis was found in non-small cell lung cancer (NSCLC) [18]. It is believed that HOXD9-mediated enhancement of PFKFB3 transcription may be involved in HMGB1 secretion and glycolytic process in hypoxia-induced GBM. To confirm the functional expression of HOXD9 in GBM cells, we assessed the expression of HOXD9, HMGB1, and PFKFB3 in human normal glial cells and glioma cell lines using RT-qPCR and western blotting analysis. We found a significant upregulation of HOXD9, HMGB1, and PFKFB3 in U87-MG and U251-MG cells compared to HMC3 cells, but no significant changes were observed in LN229 cells (Fig. 4A,B). Based on these findings, we selected the glioma cell lines U87-MG and U251-MG for further experiments. Subsequently, we exposed the cells to hypoxic conditions to investigate the role of HOXD9 in hypoxia-induced glycolytic changes in GBM. The results from RT-qPCR indicated a significantly elevated expression level of HOXD9 in the Hypoxia group compared with the Control group (Fig. 4D). To evaluate the regulatory role of HOXD9 on hypoxia-induced HMGB1 secretion and glycolysis, we transfected HOXD9 siRNA to GBM cells to silence HOXD9 expression and screened the targets with the highest silencing efficiency (Fig. 4C). As expected, HOXD9 expression in the Hypoxia+si-HOXD9 group was significantly lower than that in the Hypoxia+si-NC group, confirming the successful transfection (Fig. 4D). Notably, the knockdown of HOXD9 counteracted the hypoxia-induced upregulation of HOXD9, HMGB1, PFKFB3, and HIF-1$\alpha$ protein expression (Fig. 4E).

Fig. 4.

Silencing HOXD9 inhibited HMGB1 secretion and glycolysis in GBM cells under hypoxic conditions. (A,B) RT-qPCR and western blotting were used to detect the expression of HOXD9, HMGB1, and PFKFB3 in HMC3, U87-MG, U251-MG, and LN229 cells. (C) The HOXD9 expression was examined using western blotting. After hypoxia induction, (D) HOXD9 mRNA expression level in U87-MG and U251-MG cells was detected by RT-qPCR. (E) Western blotting was used to analyze HOXD9, HMGB1, PFKFB3, and HIF-1$\alpha$ protein expression. (F) HMGB1 contents in the cells were detected using ELISA. (G) A CCK-8 assay was used to detect the viability of U87-MG and U251-MG cells. (H) The levels of GLU intake and LD production in the cells were measured using appropriate kits. (I) Flow cytometry was performed to detect the distribution of two types of GBM cells in cell cycle stages after different treatments. (J) Transwell assay was used to determine cell migration ability. Scale bar = 100 µm. *p 0.05, vs. HMC3 or Control; #p 0.05, vs. Hypoxia+si-NC. n = 3. HOXD9, homeobox D9; HMGB1, high mobility group box 1; GBM, Glioblastoma; PFKB3, 6-phosphofructo-2-kinase/fructose-2,6-bisphosphatase 3; GLU, glucose; LD, lactic acid; RT-qPCR, real-time quantitative polymerase chain reaction; ELISA, enzyme-linked immunosorbent assay; si-NC, silencing negative control.

Subsequently, to explore the effect of HOXD9 on hypoxia-induced HMGB1 secretion, we quantified the levels of HMGB1 in the cells through ELISA. The results revealed a significant increase in HMGB1 levels in the Hypoxia group r compared to the Control group. Furthermore, reducing the expression of HOXD9 led to a significant decrease in HMGB1 levels (Fig. 4F). Based on these findings, we chose U87-MG and U251-MG cells for further experiments, maintaining them under hypoxic conditions for 48 h. In the Hypoxia group, cell viability, as measured by the CCK-8 assay, exhibited a significant increase compared to the Control group. However, knocking down HOXD9 significantly reduced cell viability (Fig. 4G). Additionally, we assessed glucose uptake and lactate production levels in U87-MG and U251-MG cells. This finding revealed a significant increase in both glucose uptake and lactate production levels in the Hypoxia group compared to the Glucose uptake and Lactate production groups. Conversely, knocking down HOXD9 significantly decreased these levels (Fig. 4H). We also examined the proliferative and migratory abilities of the cells. Flow cytometry analysis revealed a notable decrease in the number of cells in the G1 phase in the Hypoxia group compared to the Control group. Interestingly, silencing HOXD9 led to a significant increase in cell count. In the G2 and S phases, the cell number significantly increased after hypoxia induction compared with the Control group, while the cell number significantly decreased after silencing HOXD9 (Fig. 4I). Transwell assays demonstrated a significant elevation in migration levels in the Hypoxia group compared to the Control group. However, knocking down HOXD9 significantly reduced these levels (Fig. 4J). In summary, these findings suggest that HOXD9 may enhance the proliferative activity and migratory capability of glioma cells under hypoxic conditions while promoting HMGB1 secretion and glycolysis in these cells.

We conducted the following verification to investigate further the potential regulation of the functional phenotype of GBM cells under hypoxic conditions through PFKFB3. Our RT-qPCR and western blotting results showed a significant increase in the expression of HOXD9, PFKFB3, and HMGB1 in the HMGB1 group compared to the NC group. After the intervention of PFKFB3 inhibitor 3PO, we did not observe significant changes in the levels of HOXD9 and HMGB1. Silencing of HOXD9 led to a significant reduction in the expression of both HOXD9 and PFKFB3, while the expression of HMGB1 remained unchanged in the HMGB1 group (Fig. 5A,B). In order to detect the regulatory role of HOXD9/PFKFB3 on HMGB1 secretion, ELISA was used to detect the secretion level of HMGB1 in U87-MG and U251 cells. As shown in Fig. 5C, HMGB1 overexpression significantly upregulated the secretion level of HMGB1 in GBM cells. However, there was no significant change in the secretion level of HMGB1 upon inhibition of HXOD9 or inhibition of PFKFB3. This may therefore be due to the fact that HMGB1 enhances the sequence-specific DNA-binding activity of the HOXD9 protein in vitro and its transcriptional activation in vivo in a dose-dependent manner [29]. The downregulation of HOXD9 resulted in a decrease in the uptake of HMGB1, and at the same time, the promotion of HMGB1 secretion by HOXD9/PFKFB3 was inhibited. The two phases resulted in no significant change in HMGB1 levels in GBM cell supernatants. The HMGB1 group exhibited increased cell viability, changes in cell cycle distribution, and enhanced migration levels compared to the NC group. However, following intervention with the 3PO inhibitor and si-HOXD9, the cell viability and migration levels decreased significantly, and there was an increased percentage of cells in the G1 phase of the cell cycle (Fig. 5D–F). These results suggested an interaction between HMGB1 and HOXD9. We verified the expression of HOXD9 and PFKFB3 using IF double staining, which showed a significant increase in the expression levels after HMGB1 overexpression. Notably, PFKFB3 inhibitors did not significantly affect HOXD9 expression, but silencing HOXD9 inhibited the expression of PFKFB3 (Fig. 5G and Supplementary Fig. 2). Furthermore, in U87-MG and U251-MG cells under hypoxic conditions, the HMGB1 group exhibited significantly elevated levels of glucose uptake and lactic acid production. However, in the HMGB1+3PO and HMGB1+si-HOXD9 groups, the levels were significantly lower than those observed in the HMGB1 group (Fig. 5H). ChIP results showed that HOXD9 was bound to the PFKFB3 promoter region (Fig. 5I). As displayed in Fig. 5J, transfection of HOXD9 overexpression plasmid resulted in a significant increase of HOXD9 mRNA and protein levels in U87-MG and U251 cells. Overexpression of HOXD9 was observed by ChIP-qPCR data to significantly enhance the ability of HOXD9 to bind to the promoter region of PFKFB3 (Fig. 5K). These findings suggest that HOXD9 could promote cell proliferation and migration by activating PFKFB3 transcription and participating in HMGB1 secretion and glycolysis in glioma cells under hypoxic conditions. Conversely, HMGB1 could promote the expression of HOXD9.

Fig. 5.

HOXD9 promoted HMGB1 secretion and glycolysis in GBM cells under hypoxic conditions by transcriptional activation of PFKFB3. (A,B) The expression levels of HMGB1, HOXD9, and PFKFB3 were determined using RT-qPCR and western blotting. (C) The levels of HMGB1 secretion were evaluated by ELISA. (D) The viability of the cells was measured by CCK-8 assay. (E) Flow cytometry was performed to detect the distribution of two types of GBM cells in cell cycle stages after different treatments. (F) The migration abilities of the cells were evaluated using a transwell assay. Scale bar = 100 µm. (G) The expression of HOXD9 and PFKFB3 was detected using IF. (H) Relevant kits were used to measure the cells’ GLU intake and LD production levels. (I) ChIP was used to analyze HOXD9 binding to the PFKFB3 promoter region. (J) The mRNA and protein levels of HOXD9 were evaluated in U87-MG and U251 cells. (K) ChIP was utilized to analyze HOXD9 binding to the PFKFB3 promoter region. *p 0.05, vs. oe-NC; #p 0.05, vs. oe-HMGB1. n = 3. HOXD9, homeobox D9; HMGB1, high mobility group box 1; GBM, Glioblastoma; PFKB3, 6-phosphofructo-2-kinase/fructose-2,6-bisphosphatase 3; GLU, glucose; LD, lactic acid; ChIP, chromatin immunoprecipitation, 3PO, 3-(3-pyridinyl)-1-(4-pyridinyl)-2-propen-1-one.

To further investigate the role of HOXD9 in GBM formation, we conducted tumor formation experiments in nude mice. The results showed a significant reduction in both tumor volume and mass in the si-HOXD9 group compared to the Control group. Conversely, overexpression of PFKFB3 led to a significant increase in tumor volume and mass (Fig. 6A,B). Subsequently, we detected the expression of HOXD9, PFKFB3, and HMGB1 in tumor tissues. Compared to the Control group, the expressions of both HOXD9 and PFKFB3 were significantly decreased in the si-HOXD9+oe-PFKFB3 group, while the expressions in the si-HOXD9+oe-NC group were higher. Overexpression of PFKFB3 partially counteracted the inhibitory effect of silencing HOXD9 on HMGB1 protein expression (Fig. 6C,D). In addition, the silencing of HOXD9 in tumor tissues led to a considerable increase in residual glucose levels and a significant decrease in lactic acid production levels. In comparison to the si-HOXD9+oe-NC group, the residual glucose level was significantly reduced in the si-HOXD9+oe-PFKFB3 group, while the level of lactic acid production was significantly increased (Fig. 6E). Our findings suggest that the HOXD9/PFKFB3 pathway promoted HMGB1 secretion, glycolysis process, and tumor formation in GBM under hypoxic conditions.

Fig. 6.

HOXD9/PFKFB3 promoted tumor formation by promoting HMGB1 secretion and glycolysis in GBM under hypoxic conditions. (A,B) After the tumor formation experiment in nude mice, the changes in tumor volume and mass were photographed and recorded. (C) RT-qPCR was used to determine the mRNA levels of HOXD9 and PFKFB3 in tumor tissues. (D) The protein expression of HOXD9, HMGB1, and PFKFB3 in tumor tissues was detected using western blotting. (E) The levels of GLU and LD in tumor tissue were measured. *p 0.05, vs. Control, #p 0.05, vs. si-HOXD9+oe-NC. n = 5. HOXD9, homeobox D9; PFKB3, 6-phosphofructo-2-kinase/fructose-2,6-bisphosphatase 3; HMGB1, high mobility group box 1; GBM, Glioblastoma; GLU, glucose; LD, lactic acid.