DNAJB11 predicts a poor prognosis and is associated with immune infiltration in thyroid carcinoma: a bioinformatics analysis.

OBJECTIVE: To investigate the prognostic value of the co-chaperone protein DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) in thyroid carcinoma (THCA).
METHODS: This bioinformatics analysis study evaluated the prognostic value of DNAJB11 mRNA levels in THCA based on data from The Cancer Genome Atlas (TCGA). The levels of DNAJB11 mRNA in THCA and normal tissues were compared with Wilcoxon signed rank test. Kaplan-Meier survival curve analysis and Cox regression analysis were performed to evaluate the correlation between DNAJB11 mRNA levels and survival. Gene Ontology (GO) enrichment analysis was used to elucidate the functional enrichment difference.
RESULTS: Data from the 502 patients with THCA from the TCGA database were analysed. DNAJB11 mRNA was downregulated in THCA tissues compared with normal tissues. Decreased levels of DNAJB11 mRNA were significantly correlated with T stage, N stage, pathological stage, histological type, extrathyroidal extension and BRAF gene status. The low levels of DNAJB11 mRNA were associated with a shorter progression-free interval. GO enrichment analysis showed that DNAJB11 was involved in immune-related biological processes.
CONCLUSION: Low levels of DNAJB11 mRNA were associated with poor prognosis in THCA.

Introduction

Thyroid carcinoma (THCA) is the most common endocrine malignant tumour, with an increasing global incidence each year. THCA is classified into four subtypes based on histology: papillary thyroid carcinoma, follicular thyroid carcinoma, medullary thyroid carcinoma and anaplastic thyroid carcinoma. Papillary thyroid carcinoma is responsible for 80% of thyroid carcinomas, which are more common in young women. Despite advances in detection and surgical management, such as surgical resection, radiotherapy and levothyroxine treatment, metastasis and recurrence are common, with a 5-year survival rate of approximately 59% in patients at an advanced stage. As a result, it is critical to identify reliable predictors of prognostic stratification and provide new therapeutic targets. Many different types of molecular signals (such as DNA mutations and abnormal mRNA and microRNA expression) have been identified independently as favourable or unfavourable prognostic biomarkers.

The DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) gene is a protein coding gene that encodes a co-chaperone protein involved in protein processing in the endoplasmic reticulum (ER). The protein performs multifaceted functions involved in coordinating ER and extracellular proteostasis.[7,8] DNAJB11 binds to misfolded proteins in the ER and transports them to the ER HSP70 Bip for ATP-dependent chaperoning.[9-11] Recent research has revealed that the endoplasmic reticulum co-chaperone ERdj3/DNAJB11 promotes hepatocellular carcinoma progression through suppressing AATZ degradation. Moreover, DNAJB11 is positively correlated with cirrhosis, advanced hepatocellular carcinoma status and poor survival outcomes. DNAJB11 is found in high concentrations in the cytoplasm of oral cavity squamous cell carcinoma cells. The ER stress response has been implicated in malignant cell survival and development. DNAJB11 protects against ER stress by stabilizing proteins. DNAJB11 is associated with a variety of signalling pathways, including aberrant signalling pathways associated with cancer.

This current study used bioinformatics analysis of publicly available sequencing data to gain a better understanding of the potential actions of DNAJB11 in THCA and its regulatory network. RNA sequence data from The Cancer Genome Atlas (TCGA) were used in this current study to assess the prognostic value of DNAJB11 gene expression in human THCA. The differences in the levels DNAJB11 mRNA between tumours from patients with THCA and normal tissues were analysed so that any correlation between DNAJB11 mRNA levels and clinicopathological parameters could be determined.

Materials and methods

Data acquisition and preprocessing

The RNA sequence data and the corresponding clinicopathological information of patients with THCA were obtained from the TCGA database (https://www.cancer.gov/about-nci/organization/ccg/research/structural-genomics/tcga) to forecast DNAJB11 mRNA levels; and the association between the levels and the overall survival (OS) of patients with THCA. The number of fragments per kilobase of non-overlapped exons per million fragments mapped (FPKM) was calculated and then the FPKM was subsequently transferred into transcripts per kilobase million (TPM) values. There were 502 patients with THCAs included in the further analysis, while patients without RNA sequence data (n = 10) and those with OS <30 days (n = 20) were excluded. As this study analysed bioinformatics data from a public database and did not involve human subjects, informed consent was not required. In addition, no approval from an institutional review board was required.

Differentially expressed gene (DEG) analysis

The DEGs between THCA and normal tissues were identified using HTSeq-Counts by the DESeq2 software by using R. The Benjamini–Hochberg procedure was introduced to reduce the false positive rate (FDR) in multiple comparisons. The cut-off criteria were set as |log2 Fold Change (FC)| > 1.5 and adjusted P-value <0.05.

GO enrichment analysis

Metascape online software (http://metascape.org) was used as a tool to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis of DNAJB11-related DEGs. The significantly enriched (P-value <0.01, minimum count of 3, and enrichment factor of >1.5) terms of biological process were visualized.

GSEA enrichment analysis

Gene Set Enrichment Analysis (GSEA) was performed using the gseGO and gseKEGG functions of R package clusterProfiler (3.8.0). In this study, GSEA was used to elucidate the significant function and pathway differences between high- and low-DNAJB11 groups. The analysis process was repeated 1000 times using the default weighted enrichment statistics methods. Genes with an adjusted P-value <0.05 and FDR q-value <0.25 were considered significantly related genes.

Immune infiltration analysis by ssGSEA

The immune infiltration analysis of DNAJB11 was done by applying a single-sample Gene Set Enrichment Analysis (ssGSEA) method from the Gene Set Variation Analysis (GSVA) package in R (version 3.6.3) for 24 types of immune cells in tumour samples. The following 24 types of immune cells were obtained: CD4+ T cells, CD8+ T cells, gamma delta T cells, B cells, CD56bright natural killer cell, CD56dim natural killer cells, natural killer cells, dendritic cells (DCs), activated dendritic cells (aDCs), immature dendritic cells (iDCs), plasmacytoid dendritic cells, T helper cells, regulatory T (Treg) cells, follicular helper T cells, type 1 helper T (Th1) cells, type 2 helper T (Th2) cells, type 17 helper T cells, central memory T (Tcm) cells, effector memory T cells (Tem), macrophages, mast cells, neutrophils, eosinophils and cytotoxic cells.

Statistical analyses

All statistical analysis and plots were analysed using the R statistical package (R version 3.6.3; R Foundation for Statistical Computing, Vienna, Austria). The relationship between clinicopathological features and the levels of DNAJB11 were analysed using Wilcoxon signed rank test. The Kaplan–Meier method and Cox regression analysis were used to evaluate prognostic factors. Multivariate Cox analysis was used to compare the influence of DNAJB11 levels on survival along with other clinical factors, including stage, grade, histological type, extrathyroidal extension, v-raf murine sarcoma viral oncogene homolog B1 (BRAF) gene status, neuroblastoma RAS viral oncogene homolog (NRAS) gene status and Harvey rat sarcoma viral oncogene homolog (HRAS) gene status. The median DNAJB11 level was regarded as the cut-off value. The hazard ratio (HR) with 95% confidence interval (CI) was measured to estimate the hazard risk of individual prognostic factors. A P-value <0.05 was considered statistically significant in all tests.

A nomogram was constructed based on the results of the multivariate analysis to predict the prognosis of patients with THCA by R packages rms. Multivariate Cox regression analysis and Akaike’s information criterion method were used to evaluate each variable. The TCGA THCA cohort was stratified into a high-risk group and low-risk group based on the median risk score. The receiver operating characteristic (ROC) curve was used to evaluate the predictive accuracy of the prognostic model.

Results

The data collected from the TCGA included 502 patients with THCA that had both clinical and gene expression data. The characteristics of the patients are shown in Table 1. The cancer status included 502 tumour tissues and 58 adjacent tissues. A total of 367 female patients and 135 male patients were analysed in this study. The median age at diagnosis was 46 years. The number of patients with pathological stage I, II, III and IV were 280 (55.8%), 55 (11.0%), 112 (22.4%) and 55 (11.0%). A total of 223 of 502 (44.4%) patients had lymph node metastases and nine of 502 (1.8%) had distant metastases.

Table 1.

Clinical and demographic characteristics of the patients (n = 502) with thyroid carcinoma for whom data were collected from The Cancer Genome Atlas (TCGA) database stratified according to the levels of DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) mRNA.

Characteristic	Level	n #	Low levels of DNAJB11 mRNA	High levels of DNAJB11 mRNA	Statistical analysisa	Test
T stage	T1	143	62 (43.4%)	81 (56.6%)	P < 0.001	Exact
	T2	164	63 (38.4%)	101 (61.6%)
	T3	170	105 (61.8%)	65 (38.2%)
	T4	25	19 (76.0%)	6 (24.0%)
N stage	N0	229	97 (42.4%)	132 (57.6%)	P < 0.001
	N1	223	132 (59.2%)	91 (40.8%)
M stage	M0	282	153 (54.3%)	129 (45.7%)	P = 0.01	Exact
	M1	9	5 (55.6%)	4 (44.4%)
Pathological stage	Stage I	280	120 (42.9%)	160 (57.1%)	P < 0.001
	Stage II	55	21 (38.2%)	34 (61.8%)
	Stage III	112	73 (65.2%)	39 (34.8%)
	Stage IV	55	35 (63.6%)	20 (36.4%)
Residual tumour	R0	384	185 (48.2%)	199 (51.8%)	P = 0.03	Exact
	R1	52	34 (65.4%)	18 (34.6%)
	R2	4	3 (75.0%)	1 (25.0%)
Histological type	Classical	356	182 (51.1%)	174 (48.9%)	P < 0.001	Exact
	Follicular	101	29 (28.7%)	72 (71.3%)
	Other	9	6 (66.7%)	3 (33.3%)
	Tall cell	36	34 (94.4%)	2 (5.6%)
Sex	Female	367	191 (52.0%)	176 (48.0%)	NS
	Male	135	60 (44.4%)	75 (55.6%)
Race	Asian	51	21 (41.2%)	30 (58.8%)	NS
	Black or African American	27	17 (63.0%)	10 (37.0%)
	White	332	179 (53.9%)	153 (46.1%)
Neoplasm location	Bilateral	86	42 (48.8%)	44 (51.2%)	NS
	Isthmus	22	9 (40.9%)	13 (59.1%)
	Left lobe	175	82 (46.9%)	93 (53.1%)
	Right lobe	213	114 (53.5%)	99 (46.5%)
Primary neoplasm focus type	Multifocal	226	107 (47.3%)	119 (52.7%)	NS
	Unifocal	266	141 (53.0%)	125 (47.0%)
Thyroid gland disorder history	Lymphocytic thyroiditis	71	33 (46.5%)	38 (53.5%)	NS
	Nodular hyperplasia	68	30 (44.1%)	38 (55.9%)
	Normal	280	143 (51.1%)	137 (48.9%)
	Other, specify	25	15 (60.0%)	10 (40.0%)
Extrathyroidal extension	No	331	138 (41.7%)	193 (58.3%)	P < 0.001
	Yes	153	105 (68.6%)	48 (31.4%)
BRAF status	Mut	285	191 (67.0%)	94 (33.0%)	P < 0.001
	WT	200	53 (26.5%)	147 (73.5%)
NRAS status	Mut	39	8 (20.5%)	31 (79.5%)	P < 0.001
	WT	446	236 (52.9%)	210 (47.1%)
HRAS status	Mut	17	1 (5.9%)	16 (94.1%)	P < 0.001	Exact
	WT	468	243 (51.9%)	225 (48.1%)
Age, years			50.00 (38.00–60.00)	43.00 (33.00–56.00)	P = 0.001	Non-norm

Data presented as n of patients (%) or median (interquartile range).

aWilcoxon signed rank test.

#Data are missing for some patients.

T, tumour; N, node; M, metastasis; BRAF, v-raf murine sarcoma viral oncogene homolog B1; MUT, mutant; WT, wild type; NRAS, neuroblastoma RAS viral oncogene homolog; HRAS, Harvey rat sarcoma viral oncogene homolog; NS, no significant between-group difference (P ≥ 0.05).

There were significantly lower levels of DNAJB11 mRNA in THCA tissues compared with non-paired normal tissues and paired adjacent samples (P < 0.001 for both comparisons) (Figures 1A and 1B). The ROC curve analysis showed that DNAJB11 mRNA levels can effectively discriminate THCA from normal tissue in the TCGA database. The area under the curve of DNAJB11 mRNA level was 0.739 (Figure 1C). Compared with the normal tissues, differential expression was observed in bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC), cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), oesophageal carcinoma (ESCA), glioblastoma multiforme (GBM), head and neck squamous cell carcinoma (HNSC), kidney chromophobe (KICH), kidney renal clear cell carcinoma (KIRC), kidney renal papillary cell carcinoma (KIRP), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), prostate adenocarcinoma (PRAD), rectum adenocarcinoma (READ), sarcoma (SARC), stomach adenocarcinoma (STAD), THCA and uterine corpus endometrial carcinoma (UCEC) (P < 0.05 for all comparisons) (Figure 1D).

Figure 1.

Identification of the differentially expressed DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) gene. (a) Levels of DNAJB11 mRNA from The Cancer Genome Atlas (TCGA) database in non-paired thyroid carcinoma (THCA) and normal samples were analysed using Wilcoxon signed rank test. Lower levels of DNAJB11 mRNA were observed in THCA compared with normal samples. (b) Levels of DNAJB11 mRNA from the TCGA database in paired THCA and adjacent samples were analysed using Wilcoxon signed rank test. Lower levels of DNAJB11 mRNA were observed in THCA compared with adjacent samples. (c) Receiver operating characteristic curve analysis of DNAJB11 mRNA levels showing promising discriminatory power between THCA and non-tumour tissues. (D) Human DNAJB11 levels in different tumour types from the TCGA database. Wilcoxon signed rank test, ns, P ≥ 0.05; *P < 0.05; **P < 0.01; ***P < 0.001. ACC, adrenocortical carcinoma; BLCA, bladder urothelial carcinoma; BRCA, breast invasive carcinoma; CESC, cervical squamous cell carcinoma and endocervical adenocarcinoma; CHOL, cholangiocarcinoma; COAD, colon adenocarcinoma; DLBC, diffuse Large B-cell; ESCA, oesophageal carcinoma; GBM, glioblastoma multiforme; HNSC, head and neck squamous cell carcinoma; KICH, kidney chromophobe; KIRC, kidney renal clear cell carcinoma; KIRP, kidney renal papillary cell carcinoma; LAML, acute myeloid leukaemia; LGG, low grade glioma; LIHC, liver hepatocellular carcinoma; LUAD, lung adenocarcinoma; LUSC, lung squamous cell carcinoma; MESO, mesothelioma; OV, ovarian cancer; PAAD, pancreatic adenocarcinoma; PCPG, pheochromocytoma and paraganglioma; PRAD, prostate adenocarcinoma; READ, rectum adenocarcinoma; SARC, sarcoma; SKCM, ▪; STAD, stomach adenocarcinoma; TGCT, testicular germ cell tumour; THYM, thymoma; UCEC, uterine corpus endometrial carcinoma; UCS, uterine carcinosarcoma; UVM, uveal melanoma. The colour version of this figure is available at: http://imr.sagepub.com.

The Kruskal-Wallis test and Wilcoxon signed rank test showed that lower levels of DNAJB11 mRNA were significantly correlated with T stage (P < 0.001), N stage (P < 0.001), pathological stage (P < 0.001), histological type (P < 0.001), extrathyroidal extension (P < 0.001) and BRAF status (P < 0.001) (Figure 2).

Figure 2.

Correlation between the levels of DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) mRNA and clinicopathological characteristics: (a) tumour (T) stage; (b) node (N) stage; (c) pathological stage; (d) histological type; (e) extrathyroidal extension; (f) v-raf murine sarcoma viral oncogene homolog B1 (BRAF) gene status. Wilcoxon signed rank test. The heavy central black horizontal lines for each group are the means. The extremities of the box are the 25th and 75th percentiles. The error bars represent the minimum and maximum outliers. The colour version of this figure is available at: http://imr.sagepub.com.

To explore the prognostic value of DNAJB11 mRNA levels in THCA outcome, the median value of DNAJB11 mRNA was used to divide patients into DNAJB11-high and DNAJB11-low level groups. The survminer R package was used to perform a Kaplan–Meier analysis. Compared with the DNAJB11-high level group, the DNAJB11-low level group was significantly associated with shorter progression-free survival (PFS) (HR 0.41; 95% CI 0.22, 0.74; P = 0.003) (Figure 3A). To provide clinicians with a quantitative method for predicting the prognosis of THCA patients, a nomogram was constructed to integrate DNAJB11 mRNA level and other independent prognostic variables including pathological stage and age (Figure 3B). The C-index for the nomogram was 0.678 with 1000 bootstrap replicates (95% CI 0.644, 0.713). The bias-corrected line in the calibration plot was found to be close to the ideal curve (the 45-degree line), which indicated good agreement between the prediction and the observation (Figure 3C). The nomogram showed that DNAJB11 levels and pathological stage were independent predictive factors of progression-free interval in THCA patients, while age was not (Figure 3D).

Figure 3.

The prognostic value of the levels of DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) mRNA in thyroid carcinoma (THCA). (a) Kaplan–Meier survival curve analysis of patients with THCA with respect to the DNAJB11 mRNA levels. The progression-free interval in The Cancer Genome Atlas (TCGA) database cohort was performed for all patients. (b) Low levels of DNAJB11 mRNA were associated with poor progression-free interval in patients with THCA. (c) Calibration plot of the nomogram for predicting the probability of progression-free interval (PFI) at 1, 5 and 10 years. (d) Nomogram for predicting the probability of 1-, 5- and 10-year PFI for patients with THCA. The colour version of this figure is available at: http://imr.sagepub.com.

The Cox regression model was used to analyse the risk factors for the prognosis of patients with THCA. Cox regression analysis showed that DNAJB11 mRNA levels were an independent prognostic factor (HR 0.410; 95% CI: 0.237, 0.766; P = 0.009) in THCA (Table 2). In addition, T stage (HR 2.977; 95% CI 1.032, 8.586; P = 0.044), pathological stage (HR 7.197; 95% CI 2.316, 22.366; P < 0.001) and residual tumour (HR 3.356; 95% CI 1.008, 11.182; P = 0.049) were also independent prognostic factors of PFS in THCA.

Table 2.

Univariate and multivariate Cox regression analysis of risk factors for progression-free survival in patients with thyroid carcinoma.

Characteristic	n	Univariate analysis		Multivariate analysis
		HR (95% CI)	P-value	HR (95% CI)	P-value
T stage, T3&T4 versus T1&T2	502	2.977 (1.032, 8.586)	P = 0.044	0.288 (0.038, 2.161)	NS
N stage, N1 versus N0	452	1.438 (0.469, 4.407)	NS
M stage, M1 versus M0	291	4.252 (0.908, 19.922)	NS	0.000 (0.000, Inf)	NS
Pathological stage, III&IV versus I&II	502	7.197 (2.316, 22.366)	P < 0.001	11.876 (1.554, 90.748)	P = 0.017
Residual tumour, R1&R2 versus R0	440	3.356 (1.008, 11.182)	P = 0.049	4.226 (0.803, 22.226)	NS
Histological type, classical versus follicular & other & tall cell	502	4.930 (0.649, 37.425)	NS
Age, >45 versus ≤45	502	773760021.346 (0.000, Inf)	NS
Sex, male versus female	502	1.969 (0.712, 5.444)	NS
Race, white versus Asian/black or African-American	410	1.328 (0.301, 5.855)	NS
Neoplasm location, bilateral versus isthmus & left lobe & right lobe	496	1.156 (0.255, 5.235)	NS
Primary neoplasm focus type, multifocal versus unifocal)	492	0.255 (0.058, 1.130)	NS	0.237 (0.028, 2.003)	NS
Thyroid gland disorder history, lymphocytic thyroiditis & nodular hyperplasia & other, specify versus normal)	444	0.182 (0.024, 1.404)	NS
Extrathyroidal extension, yes versus no	484	2.241 (0.831, 6.040)	NS
BRAF status, Mut versus WT	485	0.608 (0.213, 1.738)	NS
NRAS status, Mut versus WT	485	0.000 (0.000, Inf)	NS
HRAS status, Mut versus WT	485	1.911 (0.247, 14.760)	NS
DNAJB11 mRNA levels, high versus low	502	0.410 (0.237, 0.766)	P = 0.009

HR, hazard ratio; CI, confidence interval; T, tumour; N, node; M, metastasis; BRAF, v-raf murine sarcoma viral oncogene homolog B1; MUT, mutant; WT, wild type; NRAS, neuroblastoma RAS viral oncogene homolog; HRAS, Harvey rat sarcoma viral oncogene homolog; DNAJB11, DnaJ Heat Shock Protein Family (Hsp40) Member B11; NS, no significant association (P ≥ 0.05).

To elucidate whether DNAJB11 could play a role in promoting THCA occurrence, RNA sequence gene expression analysis was performed to compare the gene expression profiles of the DNAJB11-high and -low groups. Using the cut-off criteria of adjust P-value <0.05, |logFC| > 1.5, the data from the TCGA database were analysed using the DESeq2 package in R. A total of 177 upregulated mRNA and 88 downregulated mRNAs were detected between the DNAJB11-high and -low groups (Figure 4A). The co-occurrence differential mRNA of the DNAJB11-high and -low groups are illustrated by a Heat map (Figure 4B). To predict the functional enrichment information of interactive mRNA of DNAJB11, clusterProfiler package was used for gene associated DEG in THCA. The top 20 GO enrichment items were classified into three functional groups: biological process group (13 items) (Figure 4C), cellular components (two items) (Figure 4D) and molecular functions analysis (five items) (Figure 4E). The biological processes of DNAJB11 included extracellular structure organization, peptide hormone secretion, extracellular matrix organization, signal release, hormone secretion, hormone transport, cartilage development, multi-multicellular organism process, positive regulation of synapse assembly, female gonad development, regulation of hormone secretion, connective tissue development and development of primary female sexual characteristics. The cellular components for these genes were involved in the collagen-containing extracellular matrix and ER lumen. The molecular functions for DNAJB11 were mainly transcription regulation by receptor ligand activity, cytokine receptor binding, cytokine activity, metal ion transmembrane transporter activity and extracellular matrix structural constituent. The analysis also found that DNAJB11 was involved in immune-related biological process, including neutrophil migration and neutrophil chemotaxis. The DNAJB11 co-occurrence derived protein–protein interaction network was constructed using the Search Tool for the Retrieval of Interacting Genes (Figure 4F).

Figure 4.

Identifying genes with differential expression between DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11)-high and -low patients with thyroid carcinoma. (a) Volcano plot of the differentially expressed mRNA. (b) Heat map of the co-occurrence differential mRNA of the DNAJB11-high and -low groups. Bubble graphs showing significantly enriched biological process (c), cellular component (d) and molecular function (e) of gene ontology (GO) annotations of DNAJB11-related mRNA. (f) The protein–protein interaction network of significant DNAJB11 co-occurrence genes. The colour version of this figure is available at: http://imr.sagepub.com.

To identify signalling pathways that were differentially activated in THCA, a GSEA was undertaken between the DNAJB11-low and -high data sets. GSEA revealed significant differences (FDR <0.25, adjust P-value <0.05) in enrichment of Molecular Signature Database Collection (h.all.v7.0.symbols). The results showed that the DNAJB11 low mRNA level phenotype consisted of many key pathways and correlated with tumorigenesis. The analysis showed that coagulation, apoptosis, complement, epithelial mesenchymal transition, mitotic spindle, p53 pathway, androgen response and oestrogen response were differentially enriched in the DNAJB11 low level phenotype.

The association between the level of DNAJB11mRNA and immune cell infiltration level quantified by ssGSEA in the THCA microenvironment was analysed using Spearman's rank correlation coefficient. The results showed that the DNAJB11 mRNA level was significantly inversely correlated with the abundance of innate immunocytes (DCs and Treg) (Figures 5A and 5B; P < 0.001). Wilcoxon signed rank test also showed that the enrichment scores of DCs, Treg cells, macrophages, aDCs, iDCs, neutrophils, Th1, Th2, mast cells, Tem, cytotoxic cells, T cells, B cells and Tcm cells were significantly higher in DNAJB11-low level samples than in DNAJB11-high level samples (Figure 5C; P < 0.001).

Figure 5.

The correlation between the level of DnaJ Heat Shock Protein Family (Hsp40) Member B11 (DNAJB11) mRNA and immune infiltration in the tumour microenvironment. (a) Correlation between the relative abundances of 24 immune cell types and DNAJB11 mRNA levels. The size of dots shows the absolute value of Spearman's rank correlation coefficient. (b) DNAJB11 mRNA had a significant inverse correlation with dendritic cells (DCs), regulatory T (Treg) cells and macrophage infiltrating levels. (c) The Wilcoxon signed rank test revealed that the enrichment scores of DCs, Treg cells and macrophages were significantly higher in low DNAJB11 mRNA samples than in the high DNAJB11 mRNA samples. Data presented as mean ± SD. aDCs, activated dendritic cells; iDCs, immature dendritic cells; Th1, type 1 helper T cell; Th2, type 2 helper T cell; Tem, effector memory T cell; Tcm, central memory T cell; NK, natura killer; Th17, type 17 helper T cell; Tfh, follicular helper T cell; Tgd, gamma delta T; pDCs, plasmacytoid dendritic cell. The colour version of this figure is available at: http://imr.sagepub.com.

Discussion

The ER stress response has been identified as a key driver in the survival of malignant cells and the progression of cancer. DNAJB11 protects against ER stress by stabilizing proteins. The function and mechanisms of DNAJB11 are only now being revealed. For example, DNAJB11 has been implicated in autosomal-dominant polycystic kidney disease. DNAJB11 is known to be closely related to a variety of signalling pathways, including abnormal signalling pathways that are prevalent in tumours. However, little is known about DNAJB11 mRNA levels and their prognostic significance in THCA. To gain a better understanding of the possible functions of DNAJB11 in THCA and its regulatory network, bioinformatics analysis of publicly available sequencing data was used to guide future research in THCA.

According to research, the role of DNAJB11 in tumorigenesis can be either inhibitory or proliferative depending on the specific tissues. DNAJB11 promotes the progression of hepatocellular carcinoma by inhibiting alpha-1-antitrypsin mutant Z degradation. The levels of DNAJB11 are increased in oral cavity squamous cell carcinoma tissues when compared with adjacent noncancerous epithelia. K1 is at the far-left end of the Kaposi sarcoma-associated herpesvirus genome that encodes a unique transmembrane glycoprotein; and DNAJB11 is required for its anti-apoptotic function. Overexpression of the GRP78 gene partially suppressed the apoptosis and induction of C/EBP homologous protein in the presence of celecoxib and this suppression was stimulated by co-expression of DNAJB11. Based on the TCGA database, this current study found that DNAJB11, compared with normal tissues, was present in significantly higher levels in BLCA, BRCA, CESC, CHOL, COAD, ESCA, GBM, HNSC, KIRC, KIRP, LIHC, LUAD, LUSC, PRAD, READ, SARC, STAD and UCEC, while it was present at significantly lower levels in KICH and THCA. These current findings showed that a lower DNAJB11 level predicted poor THCA outcomes. The prognostic value of all 502 patients with THCA from the TCGA database was evaluated. Previously, the link between DNAJB11 levels and the stage of THCA tissues has not been evaluated. However, in this current study, THCA had lower DNAJB11 levels than normal tissues and this lower DNAJB11 level was associated with an advanced clinical stage. This suggests that DNAJB11 in THCA could act as a tumour suppressor gene. In addition, DNAJB11 level was found to be correlated with poor THCA prognosis in pathological stages I to IV, T2 to T4 and N0 to N1, with the highest HR for poor PFS being associated with a low DNAJB11 level in THCA. Cox proportional hazard analysis also confirmed this result, demonstrating that DNAJB11 is an independent prognostic marker for THCA.

Disruption of ER function leads to abnormal accumulation of unfolded proteins and promotes apoptosis-induced cell death in a process known as ‘ER stress’. The ER stress response has been identified as a key driver in the survival of malignant cells and the progression of cancer, and the increased induction of the ER stress pathway has been identified as a hallmark of cancer cells. ER stress inhibition has been shown to slow the progression of cancers such as pancreatic carcinoma, breast carcinoma and multiple myeloma. DNAJB11 protects against ER stress by stabilizing proteins. To investigate the function of DNAJB11 in THCA, this current study undertook GSEA using data from the TCGA database. The current results showed that coagulation, apoptosis, complement, epithelial mesenchymal transition, mitotic spindle, p53 pathway, androgen response and oestrogen response are differentially enriched in the DNAJB11 low level phenotype. This suggests that DNAJB11 could be used as a prognostic indicator and therapeutic target in THCA. In addition, this current study evaluated the association between DNAJB11 and immune cell populations. DNAJB11 was found to have a significant correlation with DCs and Treg cells. DNAJB11 may be required for effective anti-tumour immune responses.

Furthermore, this current study found that DNAJB11 was involved in immune-related biological processes, including neutrophil migration and neutrophil chemotaxis. These analyses showed that the levels of immune infiltration and diverse immune marker sets were correlated with DNAJB11 level in THCA. These current findings demonstrated that there was a moderate-to-strong inverse correlation between DNAJB11 level and DC, Treg cell and macrophage infiltration. Treg cells are immunosuppressive immune cells that are involved in the development of tumours. The presence of Treg cells in tumour tissues has been shown to impair anti-tumour immunity, implying a poor prognosis for affected patients. Suppression of Treg cell function is considered to be an effective strategy for cancer therapy. In addition, the immune cells that are effective are suppressed in response to an increase in Treg cells. These immune effective cells are crucial in tumour escape. ER stress plays a critical role in the induction of inflammation in a variety of pathological conditions.

Recent research indicates that macrophages are primarily responsible for tumour promotion in the immune microenvironment. M2-type tumour associated macrophages have been linked to the aggressiveness of various cancers. Thus, this current research contributes to a better understanding the potential role of DNAJB11 in THCA immunology and its application as a cancer biomarker.

This current study had several limitations. First, the correlation between DNAJB11 mRNA and protein levels should be validated using cell experiments and clinical samples, as predicting protein levels solely on the basis of mRNA levels is far from perfect. Secondly, larger studies of clinical data should be used to assess the relationship between DNAJB11 levels and clinical features, THCA stage and prognosis.

In conclusion, DNAJB11 may be associated with a poor prognosis and may play important roles in the regulation of the ER stress pathway in THCA. This current study revealed some of the roles of DNAJB11 in THCA and identified a potential biomarker for THCA diagnosis and prognosis.