mRNA expression data of six members of DLX gene family in The Cancer Genome Atlas (TCGA) dataset (HTSeq Counts) were acquired from University of California, Santa Cruz (UCSC) Xena website (xenabrowser.net) (23); the project number was TCGA-KIRC. The data were processed using RSEM normalization methods and expressed as log₂(x+1) transformed. Clinical information of 533 patients, including age, sex, pathological stage, pathological grade, TNM stage, overall survival and survival status were obtained from UCSC Xena. DFS of patients with ccRCC was obtained from the cBio Cancer Genomics Portal (http://cbioportal.org).

The differential expression of six members of the DLX gene family was compared between ccRCC and adjacent normal renal tissues in TCGA dataset. The differential expression of DLX4 was also compared between the different clinical subgroups (advanced vs. early stage, high vs. low grade, N0 vs. N1, M0 vs. M1). Kaplan-Meier curves and log-rank tests were performed to evaluate the prognostic significance, and receiver operating characteristic (ROC) curves were generated to assess the diagnostic efficiency (24). Univariate and multivariate Cox hazard analyses were performed to assess the overall survival and DFS of patients with ccRCC.

A total of 40 pairs of ccRCC tumors and matched normal renal tissues were retrieved from the Department of Urology at the Second People's Hospital of Wuhu (Wuhu, China) and The First Affiliated Hospital of Anhui Medical University (Hefei, China) between May 2019 and May 2020. Adjacent normal renal tissues, ≥2 cm away from the location of the tumor, were cut and collected. Samples were placed in RNAlater stabilization solution (Invitrogen; Thermo Fisher Scientific, Inc.) and frozen in liquid nitrogen. These samples were evaluated by two pathologists to confirm the histopathological results of ccRCC. Inclusion criteria were ccRCC, and exclusion criteria were other types of RCC and benign renal tumors. Written informed consent was provided by all patients. The present study was approved by the Ethics Committee of Human Research of The Second People's Hospital of Wuhu and The First Affiliated Hospital of Anhui Medical University (approval no. PJ2019-14-22). The study methodology conformed with the standards of the Declaration of Helsinki.

Total RNA was extracted from ccRCC tissue using the TRIzol^® (Invitrogen; Thermo Fisher Scientific, Inc.) method. The purity and concentration of total RNA solution were detected using a NanoDrop spectrophotometer (NanoDrop Technologies; Thermo Fisher Scientific, Inc.). Extracted RNA was reverse transcribed into cDNA using a PrimeScript™ RT with gDNA Eraser kit (Takara Bio, Inc.) according to the manufacturer's protocol. The cDNA was subjected to qPCR by a SYBR Green Mix (Takara Bio, Inc.). Thermocycling conditions were as follows: Initial denaturation at 95 for 30 sec, followed by 40 cycles of 95°C for 5 sec and 60°C for 34 sec, then final extension at 95°C for 15 sec, 60°C for 1 min and 95°C for 15 sec. A final reaction volume of 20 µl was used, and the ABI 7500 Real-Time PCR System (Thermo Fisher Scientific, Inc.) was used to measure the reactions. Relative gene expression of DLX4 were analyzed using the 2^−ΔΔCq method (25). GAPDH was used as endogenous control to normalize the relative expression of DLX4. The primers were purchased from Sangon Biotech Co., Ltd., and the sequences were as follows: DLX4 forward, 5′-CAGCACCTAAACCAGCGTTTC-3′ and reverse, 5′-GAGCTTCTTATACTTGGAGCGTT-3′; and GAPDH forward, 5′-GGGAGCCAAAAGGGTCAT-3′ and reverse, 5′-GAGTCCTTCCACGATACCAA-3′.

Paraffin-embedded tissue specimens were sliced into 4-µm-thick sections. The slides were placed in a 60°C thermostat and baked for 20 min to dewax. Then, the slides were hydrated in xylene I (30 min), xylene II (30 min), anhydrous ethanol (5 min) twice and 95, 90, 80 and 70% ethanol (5 min each). Then, the slides were heated at 100°C for 10 min in citric acid buffer (0.01 M, pH 6.0) for antigen retrieval. The slides were incubated in 3% hydrogen peroxide solution (cat. no. SP 9000; Beijing Zhongshan Jinqiao Biotechnology Co., Ltd.) for 15 min at room temperature and washed three times in phosphate-buffered saline (PBS; pH 7.4). The slides were blocked with 3% bovine serum albumin (cat. no. G5001; Servicebio) for 1 h at room temperature, then incubated with anti-DLX4 antibody (1:100; cat. no. DF3387; Affinity Biosciences, Ltd.) overnight at 4°C. After three washes with PBS, the slides were incubated with biotinylated goat anti-rabbit IgG (1:200; cat. no. G23303; Servicebio) for 2 h at room temperature. Finally, DAB (Beyotime Institute of Biotechnology) was used to detect the immune complexes and the sections were counterstained with hematoxylin for 3 min at room temperature. The images were captured using a light microscope (Olympus Corporation). DLX4 expression was measured based on the intensity of immune staining (intensity score) with ImageJ version 6.0 software (National Institutes of Health), as previously described (26).

Patients in TCGA data set were divided into high and low expression groups according to the median values of DLX4 gene expression. GSEA was conducted to explore the possible functional mechanisms of DLX4 in ccRCC with GSEA 3.0 software (Broad Institute; broad.mit.edu/gsea) (27). The reference gene set file was h.all.v7.0.symbols.gmt; normalized P-value <0.05 and FDR <0.25 were set as threshold values.

Pearson correlation coefficients was calculated to determine the association between DLX4 gene expression and genes in TCGA dataset. Genes that strongly correlated with the expression levels of DLX4 (correlation coefficient R>0.4) were obtained as previously described (28). A total of 559 positively correlated genes and 226 negatively correlated genes were identified. Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Ontology (GO) enrichment analyses were performed using The Database for Annotation, Visualization and Integrated Discovery (https://david.ncifcrf.gov) to investigate the possible molecular mechanisms of these selected genes. The cutoff value was set at P<0.05.

To quantify the cellular composition of the immune infiltrates in each patient in TCGA dataset, a set of metagenes, including non-overlapping sets of genes that are representative of 28 specific immune cell subpopulations, was obtained from a previous study (29). Subsequently, single-sample GSEA was used to quantify the 28 types of immune cells based on the metagenes. In the tumor microenvironment, immune and stromal cells are the two main non-tumor components; these components have been proposed to be valuable for tumor treatment and prognosis (30). To assess the tumor microenvironment associated with DLX4 expression levels, the immune and the stromal scores (which reflect the infiltration levels of non-tumor cells) for the TCGA dataset were calculated using the ESTIMATE package (30). The differences in immune cell composition and immune and stromal scores were compared between the high and low DLX4 expression groups, which were divided based on the median expression values. In addition, the correlations between DLX4 expression levels and immunosuppressive immune cells were evaluated.

Tumor mutation burden (TMB) is defined as the total amount of coding errors of somatic genes, base substitutions, insertions or deletions detected per million bases (31). In the present study, 38 megabases was used as the exon length. TMB was calculated as (the number of variants)/(exon length) for each patient with ccRCC using Perl scripts (32). The somatic mutation status data of ccRCC samples (workflow type: VarScan2 Variant Aggregation and Masking) were downloaded from the TCGA data portal (https://portal.gdc.cancer.gov/repository) in May 2021. The cytolytic activity scores were calculated as the geometric mean of the granzyme A and perforin expression (13). The TMB scores and cytolytic activity scores were also compared between high and low DLX4 expression groups.

For equivalent variables with a normal distribution, an unpaired Student's t-test was performed for comparisons between the two groups, and one-way ANOVA followed by Tukey's post hoc test was used for comparison between multiple groups. The paired Student's t-test was used to compare paired samples. For non-normally distributed variables, the Wilcoxon rank-sum test was performed for comparisons between two groups. Pearson's χ² test was used to evaluate the distributions of clinical factors between the high and low-DLX4 expression groups. Log-rank tests were used for Kaplan-Meier curves to assess survival differences. The correlations between DLX4 expression and tumor immunosuppressive cells were evaluated using Pearson's correlation analysis, and ‘general linear model’ method was used to fit curves. P<0.05 was used to indicate a statistically significant difference.

To investigate the roles of DLX family genes in ccRCC, the mRNA expression levels of six members of the DLX family in the TCGA dataset were investigated. A heat map revealed the overall expression levels of the six members of the DLX family between normal and tumor samples in TCGA dataset (Fig. 1A). Furthermore, the expression levels of six members were compared between normal and tumor samples. Downregulated levels of DLX3 and upregulated levels of DLX1, DLX2, DLX4, DLX5 and DLX6 were observed in tumor samples compared with normal samples in TCGA dataset (Fig. 1B). The expression levels of paired normal and tumor samples in TCGA dataset were also compared, which revealed the same expression patterns (Fig. 1C).

Patients in the TCGA data set were divided into high and low expression group according to the median values of the mRNA expression levels for each gene of the DLX family. A total of 511 patients with ccRCC with overall survival and 422 with DFS information were included. The differences in overall and DFS status for six members of DLX family were compared using the χ² test. The distribution for overall survival status was found to be significantly different for DLX2 and DLX4, whereas the distribution for DFS status was significantly different for DLX2, DLX4 and DLX5 (Table I). Kaplan-Meier curves were used to evaluate the prognostic significances of the DLX family in ccRCC. Kaplan-Meier curves for overall survival revealed poor prognosis in high DLX2 expression group (Fig. 2B), DLX4 (Fig. 2D) and DLX5 (Fig. 2E) for patients with ccRCC. However, the overall survival rates between the high and low expression group of DLX1 (Fig. 2A), DLX3 (Fig. 2C) and DLX6 (Fig. 2F) were not significantly different. DFS was also compared between the high and low expression groups. A significantly lower DFS rate was identified in the high expression group of DLX2 (Fig. 2H), DLX4 (Fig. 2J) and DLX5 (Fig. 2K) relative to the low expression group. However, the DFS rate was not significantly different between DLX1 (Fig. 2G), DLX3 (Fig. 2I) and DLX6 (Fig. 2L). The expression of DLX2 and DLX4 was revealed to be associated with overall and DFS in ccRCC patients in both survival distribution difference assessment and Kaplan-Meier curves. Difference in either overall survival or DFS were not identified for several genes, including DLX3 and DLX5. Therefore, DLX2 and DLX4 were selected for subsequent analyses.

To investigate the association of DLX2 and DLX4 with the distributions of clinical phenotypes, patients were divided into two subgroups according to the median mRNA expression levels. A high proportion of males, advanced pathological stage, higher grade, higher T stage, N1 stage, and M1 stage were identified in the high expression group of DLX4 (Table II); however, no significant difference was identified for age distribution (Table II). In addition, the expression of DLX2 was not shown to be associated with the clinicopathological variables examined for patients with ccRCC (Table II). Therefore, DLX4 was selected for further studies. Univariate Cox regression analysis showed that the expression level of DLX4 was a risk factor for overall survival (Table III) and DFS (Table IV) in ccRCC patients. Based on multivariate Cox hazard analysis, the expression level of DLX4 was an independent risk variable for the overall survival of ccRCC patients (Table III) and DFS (Table IV) after integration with multiple clinicopathological characteristics. These findings indicated that the expression level of DLX4 was as an independent prognostic predictor for ccRCC patients, which may provide a supplement for clinical factors.

The relationships between clinical variables and DLX4 expression levels were investigated in the TCGA dataset. The expression levels of DLX4 were significantly higher in stage 3&4 vs. 1&2 (Fig. 3A and D), grade 3&4 vs. 1&2 (Fig. 3B and E) and T 3&4 vs. T 1&2 (Fig. 3C and F). In addition, high expression levels of DLX4 were observed in N1 vs. N0 stage (Fig. 3G), M1 vs. M0 stage (Fig. 3H) and dead vs. alive (Fig. 3I). These results suggested that the expression levels of DLX4 were highly associated with tumor malignancy.

Based on the prognostic values of DLX4 for the overall and DFS of patients with ccRCC in TCGA dataset, it was sought to further confirm the prognostic significances of DLX4 in different subgroups of patients. A worse overall survival rate was revealed in the high DLX4 expression group compared with the low expression group for male, female, stage I and II, stage III and IV, grade 1 and grade 2, grade 3 and grade 4, low age (≤60 years) and elderly (>60 years) subgroups (Fig. 4A). In addition, a worse DFS rate was identified in high compared with low expression group for male, female, stage I and II, grade 1 and grade 2, high age (>60 years) and low age (≤60 years) subgroups (Fig. 4B). However, DFS was not significantly different between stage III and IV and grade 3 and grade 4 subgroups. These results suggested that DLX4 may be a potential prognostic biomarker for ccRCC. To determine whether the expression of DLX4 has diagnostic values in ccRCC, ROC curves were generated and AUCs were calculated to determine the diagnostic efficiency. DLX4 could sufficiently differentiate ccRCC from normal tissues with an AUC of 0.863 in the entire TCGA dataset (Fig. 4C). In addition, the diagnostic values of DLX4 expression levels were evaluated for paired normal and adjacent tumor samples in TCGA dataset, with an AUC of 0.837 (Fig. 4D). These results suggested that DLX4 may be a potential biomarker for ccRCC patients with favorable diagnostic performance.

To further confirm the expression levels of DLX4 in ccRCC vs. normal renal tissue, RT-qPCR and immunohistochemistry were used to verify the results of the TCGA databases at the mRNA and protein levels. The results revealed higher DLX4 mRNA expression levels in ccRCC samples compared with matched normal renal samples (Fig. 5A and B). Furthermore, the expression levels of DLX4 in ccRCC samples and the matched normal renal samples were detected using immunohistochemistry, which revealed higher protein levels of DLX4 in ccRCC tissues (Fig. 5C and D). These results indicated that DLX4 was highly expressed in ccRCC tissues, which aligns with the TCGA database analysis results aforementioned.

To investigate the functional mechanism of DLX4 in ccRCC, genes that are highly associated with the expression levels of DLX4 in ccRCC were first selected; of which 328 positively related genes and 210 negatively associated genes were selected (Fig. 6A). Functional enrichment analysis of these genes was carried out to determine the potential mechanisms. The GO terms, including ‘positive regulation of ubiquitin protein ligase activity’, ‘positive regulation of protein phosphorylation’, ‘mitotic nuclear division’, ‘microtubule cytoskeleton organization’, and ‘cell division’ were enriched for biological process (Fig. 6B); ‘extracellular region’, ‘endoplasmic reticulum lumen’, ‘extracellular exosome’, and ‘cytoplasmic microtubule’ were enriched for cellular component (Fig. 6C) and ‘protein domain specific binding’, ‘identical protein binding’, ‘RNA polymerase II regulatory region sequence-specific DNA binding’, ‘protein kinase binding’, and ‘enzyme binding’ were enriched for molecular function (Fig. 6D). The terms ‘cell cycle’, ‘PPAR signaling pathway’ and ‘fatty acid degradation’ were enriched for KEGG pathways (Fig. 6E). In addition, the terms that are more closely related to the cell cycle and complement and coagulation cascades were selected and highlighted in red in Fig. 6.

Patients with ccRCC in the TCGA dataset were divided into high and low expression groups based on median expression level of DLX4. GSEA was then performed to determine the statistical significance of the biological pathways associated with the expression levels of DLX4 (Fig. 7A, C and E). The results revealed that the expression of DLX4 was associated with biological pathways related to G2M checkpoint, epithelial-mesenchymal transition (EMT), glycolysis, inflammatory response and TNFα signaling via NF-κB. The distributions of cell cycle-related genes for high and low DLX4 expression groups are presented in Fig. 7B. Based on the heat map, patients in the high DLX4 expression group exhibited a trend of higher expression levels of cell cycle-related genes compared with the low DLX4 expression group. Genes involved in EMT were expressed at notably higher expression levels in the high DLX4 expression group (Fig. 7D). The distributions of the inflammatory molecules, including chemokines, cytokines and their receptors, between the high and low expression groups of DLX4 were also compared (Fig. 7F). A trend towards higher expression of inflammatory molecules was identified in patients in the high DLX4expression group. These findings suggested that cell cycle-related pathways, EMT and inflammatory response may be a potential mechanism of DLX4 in ccRCC tumorigenesis.

Based on the associations of DLX4 expression with inflammatory and immune responses, the associations between the expression levels of DLX4 and immune cell infiltration were investigated by examining the differences in the expression profiles of 28 immune cell types separated into high and low DLX4 expression groups of DLX4 (Fig. 8A and B). Patients with high DLX4 expression in TCGA dataset had a high percentage of immune cells, including immature dendritic cells, activated dendritic cells, gamma delta T cells, macrophages, myeloid-derived suppressor cells (MDSCs), plasmacytoid dendritic cells, regulatory T cells and T follicular helper cells. In addition, a positive correlation was determined between the expression levels of DLX4 and tumor immunosuppressive cells, including immature dendritic cells, activated dendritic cells, gamma delta T cells, macrophages, MDSCs, plasmacytoid dendritic cells, regulatory T cells and T follicular helper cells (Fig. 8C). The immune and stromal scores for TCGA cohort were calculated using the ESTIMATE algorithm. Patients in the high DLX4 expression group had significantly higher immune and stromal scores compared with the low DLX4 expression group (Fig. 8D and E, respectively). Furthermore, patients in the high DLX4 expression group had significantly higher TMB and cytolytic activity scores compared with those in the DLX4 low expression group (Fig. 8F and G, respectively).