PeerJ
A transcription factor signature predicts the survival of patients with adrenocortical carcinoma
Jianyu Zhao1, Bo Liu1 and Xiaoping Li2
1 Department of Endocrinology, China-Japan Union Hospital of Jilin University, Changchun, Jilin, China
2 Department of Pediatrics Endocrinology, The First Hospital of Jilin University, Changchun, Jilin, China, Jilin, China
ABSTRACT
Background: Adrenocortical carcinoma (ACC) is a rare endocrine cancer that manifests as abdominal masses and excessive steroid hormone levels and is associated with poor clinical outcomes. Transcription factors (TFs) deregulation is found to be involved in adrenocortical tumorigenesis and cancer progression. This study aimed to construct a TF-based prognostic signature for the prediction of survival of ACC patients.
Methods: The gene expression profile and clinical information for ACC patients were downloaded from The Cancer Genome Atlas (TCGA, training set) and Gene Expression Omnibus (GEO, validation set) datasets after obtained 1,639 human TFs from a previously published study. The univariate Cox regression analysis was applied to identify the survival-related TFs and the LASSO Cox regression was conducted to construct the TF signature based on these survival-associated TFs candidates. Then, multivariate analysis was used to reveal the independent prognostic factors. Furthermore, Gene Set Enrichment Analysis (GSEA) was performed to analyze the significance of the TFs constituting the prognostic signature.
Results: LASSO Cox regression and multivariate Cox regression identified a 13-TF prognostic signature comprised of CREB3L3, NR0B1, CENPA, FOXM1, E2F2, MYBL2, HOXC11, ZIC2, ZNF282, DNMT1, TCF3, ELK4, and KLF6. The risk score based on the TF signature could classify patients into low- and high-risk groups. Kaplan-Meier analyses showed that patients in the high-risk group had significantly shorter overall survival (OS) compared to the low-risk patients. Receiver operating characteristic (ROC) curves showed that the prognostic signature predicted the OS of ACC patients with good sensitivity and specificity both in the training set (AUC > 0.9) and the validation set (AUC > 0.7). Furthermore, the TF-risk score was an independent prognostic factor.
Conclusions: Taken together, we identified a 13-TF prognostic marker to predict OS in ACC patients.
Submitted 10 June 2021 Accepted 13 October 2021 Published 9 December 2021 Corresponding author Xiaoping Li, xiaopingli@jlu.edu.cn Academic editor Katherine Mitsouras Additional Information and Declarations can be found on
page 14
DOI 10.7717/peerj.12433
cc Copyright 2021 Zhao et al.
Distributed under Creative Commons CC-BY 4.0
Subjects Bioinformatics, Oncology, Medical Genetics Keywords Adrenocortical carcinoma, TCGA, GEO, Transcription factor, Prognosis
OPEN ACCESS
INTRODUCTION
Adrenocortical carcinoma (ACC) is a rare endocrine cancer with an annual incidence of 0.7-2.0 cases per million (Bilimoria et al., 2008; Else et al., 2014). It usually affects adults aged around 40-50 years and children younger than 10 years (Koschker et al., 2006; Rodgers et al., 2006). Clinical manifestations of ACC include abdominal masses and elevated steroid hormones and result in overall poor outcomes with 5-year survival ranging from 32% to 45% (Icard et al., 2001). Therefore, it is essential to identify prognostic markers of ACC to screen for patients at high risk.
Transcription factors (TFs) are regulatory proteins that bind to the promoter sequences of genes and decrease or increase their transcription (Latchman, 1997), and thus control cell differentiation (Arendt et al., 2016), proliferation (Bouchard, Staller & Eilers, 1998), and death (Shaulian & Karin, 2002). Given the important role of TFs in determining cell fate, the contribution to tumor development and progression is easily made. Not surprisingly, the genes encoding TFs are often aberrantly expressed in ACC (Brown et al., 2018), and the relationship between the TFs and survival has been implicated in ACC patients (Parviainen et al., 2013b). For instance, Snail is overexpressed in numerous ACC patients and associated with decreased survival (Waldmann et al., 2008). In addition, TGF-ß pathway components including GATA-6 and SF-1 are also correlated with poor outcomes in ACC patients (Parviainen et al., 2013a). A recent study has shown that high expression of HOXB9 is associated with a poorer prognosis of ACC patients (Francis et al., 2021). Although these reports demonstrated that the TFs are closely correlated with the prognosis of ACC patients, fewer studies have investigated the prognostic value of the TF signature in ACC patients. Therefore, it is necessary to investigate the correlation between the TFs and the survival of ACC patients and construct a prognostic TF signature for predicting the survival of ACC patients.
The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) datasets have enabled researchers to correlate clinical outcomes with transcriptomic profiles. To this end, we systematically analyzed the gene expression data of ACC datasets using univariate and multivariate Cox regression models. Based on the survival analysis, we then developed a 13-TF prognostic signature and validated this model in an independent microarray data set from GEO.
MATERIALS AND METHODS
Data collection and preprocessing
A total of 1,639 human TFs were identified from a previously published study (Lambert et al., 2018). TCGA ACC level 3 RNAseqv2 data (RSEM_genes_normalized file) and corresponding clinical information were downloaded from the TCGA database (https://tcga-data.nci.nih.gov/). A total of 79 patients with ACC from TCGA were included after excluding those lacking complete clinical and survival data, which served as training sets. The gene expression profiling of 22 patients with ACC was selected via the accession number GSE19776 from GEO (http://www.ncbi.nlm.nih.gov/geo/) using the GPL570 platform [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array,
which served as validation set. The clinical data included age, gender, and tumor stage. Firstly, we obtained 1,554 overlapped TFs between TFs in literature and TCGA dataset, and 1,508 overlapped TFs between TFs in literature and GEO dataset. A total of 118 TFs were excluded from the TCGA dataset by removing lower expression genes. Next, we acquired the common TFs between the TCGA and GEO datasets.
Survival analysis
After obtained common TFs between TCGA and GEO datasets, univariate Cox survival analysis was performed using Cox proportional hazards regression model. Only TFs with Cox P < 0.001 and Log Rank P < 0.001 on univariate analysis were incorporated into the Lasso Cox regression analysis. Kaplan-Meier method was used to analyze the correlation between overall survival (OS) and TF expression, and the OS of different patient groups were compared using the Log-Rank test. The “survival” package (Therneau & Lumley, 2015) in R software was used for survival analysis and the time-dependent ROC (Receiver Operating Characteristic) curve analysis was performed using the “survival ROC” package (Robin et al., 2011).
Construction of ACC prognostic signature
LASSO Cox regression model was widely used for high-dimensional predictor identification. In the present study, OS-associated TFs were used to select the significant TFs associated with the OS of ACC patients according to the coefficient value. These factors were incorporated in the multivariate Cox regression model to construct the ACC prognostic signature. The risk score for each TF-encoding gene was calculated as follows: Risk score = ""1 expi * Bi, where n is the number of selected genes, expi is the expression level of gene i, and ß; is the coefficient of gene i.
mRNA expression of 13 TFs in ACC cell line
Gene expression data of all 13 TFs in the prognostic signature from the SW13 cell line were downloaded directly from the Cancer cell line encyclopedia (CCLE) database (https://portals.broadinstitute.org/ccle/home), which is an open-access database covering large-scale deep sequencing raw data of more than 900 human cancer cell lines.
Gene set enrichment analysis (GSEA)
GSEA was performed to analyze the significance of the 13 TFs constituting the prognostic signature using GSEA v2.0.12 (http://www.broadinstitute.org/gsea) by computing the enrichment score for each gene set (Subramanian et al., 2005). The distributions of these TFs against the rank-ordered gene ontology (GO) hallmarks were characterized using GSEA with the default settings. Positive and negative normalized scores indicated enrichment in the high-risk and low-risk groups respectively.
RESULTS
Construction of the TF prognostic signature
The scheme for developing the TF signature is outlined in Fig. 1. After initially identifying 1,639 TFs by literature search, the down-regulated genes were removed and 1,304 common
Transcription factors in literature (1639)
Overlap with TGCA gene expression dataset
Overlap with GEO gene expression dataset
Transcription factors (1554)
Transcription factors (1508)
Remove lower expression gene
Transcription factors (1366)
Transcription factors (1508)
Overlap
Transcription factors (1304)
Univariate analysis of OS Cox P< 0.001 & Log-Rank P < 0.001
Survival correlated transcription factors (23)
Lasso regression
Genes for constructing the prognostic signature (13)
Risk score (Ei-1 exp; * B;) High risk threshold determination
Training set TCGA (n=79)
Validating set GSE19776 (n=22)
Univariate & Multivariate Cox analysis
ROC analysis
Clinical associated analysis
TFs were screened from TCGA and GEO datasets. Using the TCGA training set (n = 79), univariate regression analysis showed that 23 TFs were correlated with OS (Cox-P < 0.001 and Log-Rank P < 0.001), of which 13 were identified by the Lasso regression analysis, and the risk score was calculated by multivariate cox analysis.
The X value was selected in the LASSO Cox regression analysis when the median of the sum of squared residuals was the smallest (Fig. 2A). The following survival-related TFs with non-zero coefficients were then screened: CREB3L3 (cAMP-responsive element-binding protein 3-like 3), NR0B1 (nuclear receptor subfamily 0, group B, member 1), CENPA (centromere protein-A), FOXM1 (Forkhead Box M1), E2F2 (E2F transcription factor 2), MYBL2 (v-myb myeloblastosis viral oncogene homolog (avian)- like 2), HOXC11 (homeobox C11), ZIC2 (Zic family member 2), ZNF282 (zinc finger protein 282), DNMT1 (DNA methyltransferase 1), TCF3 (transcription factor 3), ELK4 (ETS transcription factor ELK4) and KLF6 (Krüppel-like factor 6) (Fig. 2B). Only CREB3L3 and NR0B1 were negatively correlated with the remaining TFs (Fig. 2C), while CENPA, FOXM1, and E2F2 displayed a strong correlation (Table 1). As shown in Fig. 3A, a total of 79 patients from the TCGA database were divided into a high-expression
A
B
0 20 19 17 16 14 13 14 9 8 1
ZNF282
ZIC2
ZBTB34
TET1
Partial Likelihood Deviance
8
TCF3
STAT5B
OTX1
NROB1
2
MYBL2
LHX2
KLF6
HOXC11
2
HIC2
FOXM1
ELK4
E2F2
Y
DPF1
DOT1L
DNMT1
10
CREB3L3
CENPA
-7
-6
-5
-4
-3
-2
-1
-0.1
0.0
0.1
0.2
0.3
Log(lambda)
Coefficients
C
CREB3L3
NROB1
CENPA
FOXM1
E2F2
MYBL2
HOXC11
ZIC2
ZNF282
DNMT1
TCF3
ELK4
KLF6
1
CREB3L3
1.00
NROB1
0.29
1.00
X
0.75
CENPA
-0.39
-0.44
1.00
FOXM1
0.5
-0.30
-0.38
0.86
1.00
E2F2
-0.25
-0.32
0.84
0.87
1.00
0.25
MYBL2
-0.35
-0.40
0.84
0.86
0.87
1.00
HOXC11
-0.34
0.48
0.47
0.51
0.44
1.00
0
ZIC2
-0.45
0.43
0.38
0.44
0.42
0.62
1.00
ZNF282
0.25
-0.43
-0.30
0.37
0.32
0.29
0.32
0.43
0.42
1.00
DNMT1
-0.47
-0.35
0.76
0.70
0.61
0.66
0.38
0.39
0.47
1.00
-0.5
TCF3
-0.52
-0.43
0.57
0.45
0.48
0.57
0.47
0.39
0.55
0.64
1.00
ELK4
-0.50
-0.47
0.39
0.33
0.27
0.33
0.32
0.41
0.30
0.49
0.32
1.00
0.75
KLF6
-0.37
-0.37
0.40
0.31
0.27
0.46
0.56
0.37
0.42
1.00
-1
DOI: 10.7717/peerj.12433/fig-2
| Table 1 Information of the 13 transcription factors for constructing the prognostic signature. | |||
|---|---|---|---|
| Gene symbol | Stable ID | Gene type | Chr position (start-end) |
| CENPA | ENSG00000115163 | Protein coding | 2 (26764289-26801067) |
| CREB3L3 | ENSG00000060566 | Protein coding | 19 (4153631-4173054) |
| DNMT1 | ENSG00000130816 | Protein coding | 19 (10133345-10231286) |
| E2F2 | ENSG00000007968 | Protein coding | 1 (23506438-23531233) |
| ELK4 | ENSG00000158711 | Protein coding | 1 (205597556-205631962) |
| FOXM1 | ENSG00000111206 | Protein coding | 12 (2857681-2877155) |
| HOXC11 | ENSG00000123388 | Protein coding | 12 (53973126-53977643) |
| KLF6 | ENSG00000067082 | Protein coding | 10 (3775996-3785281) |
| MYBL2 | ENSG00000101057 | Protein coding | 20 (43667019-43716495) |
| NR0B1 | ENSG00000169297 | Protein coding | X (30304206-30309598) |
| TCF3 | ENSG00000071564 | Protein coding | 19 (1609290-1652615) |
| ZIC2 | ENSG00000043355 | Protein coding | 13 (99981784-99986765) |
| ZNF282 | ENSG00000170265 | Protein coding | 7 (149195546-149226238) |
Note: Chr, chromosome.
group (n = 40) and a low-expression group (n = 39). Significantly, high expression of CREB3L3 (HR = 0.663, Cox P= 5e-05, Log-Rank P = 1.97e-07) and NR0B1 (HR = 0.799, Cox P = 6.93e-05, Log-Rank P = 3.18e-06) were associated with good prognosis, and that of other TFs with poor prognosis (HR > 1, P < 0.001).
TF signature risk score can predict prognosis of ACC patients
The clinical relevance of these TFs was further assessed by multivariate Cox regression analysis based on the TCGA training set, and the risk score based on their expression levels and coefficients was calculated. The 13-TF risk score classified the patients from TCGA training set into the high- (n = 39) and low-risk (n = 40) groups (Fig. 4A). Except for CREB3L3 and NR0B1, all TFs were overexpressed in the high-risk group (Fig. 4B). Furthermore, ACC patients in the high-risk group had significantly shorter survival compared to the low-risk patients (HR = 16.95 (5.02-57.2); Cox P = 5.11e-06; Log Rank P = 2.09e-09) (Fig. 4C). The sensitivity and specificity of the 13-TF signature were determined using time-dependent receiver operating characteristic (ROC) analysis, and the area under curves (AUCs) at all follow-up time points were greater than 0.9 (Fig. 4D). The predictive model was then validated in a GEO dataset, and as shown in Fig. 5A, the high-risk group (n = 11) had worse survival compared to the low-risk group (n = 11). In addition, the AUC values of the signature were greater than 0.75 (Fig. 5B). Next, we investigated the mRNA levels of all 13 TFs in the ACC cell line through the CCLE database (Fig. 5C). The results showed that the mRNA levels of CREB3L3, NR0B1, and HOXC11 were almost undetectable in SW13 cells. In contrast, the other 10 TFs were highly expressed in ACC cells. Taken together, the 13-TF signature can predict the survival of ACC patients with high sensitivity and specificity.
CENPA
CREB3L3
DNMT1
E2F2
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
0.8
0.8
0.8
0.8
OS
0.6
OS
0.6
OS
0.6
OS
0.6
0.4
0.4
0.4
0.4
0.2
0.2
0.2
0.2
HR=2.07
HR=0.663
HR=4.651
HR=2.063
0.0
Cox P=1.78e-07
Log-Rank P=6.23e-06
0.0
Cox P=5e-05
Log-Rank P=1.97e-07
0.0
Cox P=2.8e-07
Log-Rank P=3.51e-05
0.0
Cox P=1.65e-06
Log-Rank P=1.69e-06
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
ELK4
FOXM1
HOXC11
KLF6
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
0.8
0.8
0.8
0.8
#
OS
0.6
OS
0.6
OS
0.6
OS
0.6
0.4
0.4
0.4
0.4
0.2
0.2
0.2
0.2
HR=2.625
HR=2.445
HR=1.478
HR=2.301
0.0
Log-
Cox P=3.61e-05
Rank P=3.07e-06
0.0
Cox P=3.07e-06 Log-Rank P=4.32e-05
0.0
Cox P=7.55e-08
Log-Rank P=8.74e-08
0.0
Cox P=3.51e-05
Log-Rank P=4.18e-05
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
NROB1
MYBL2
TCF3
ZIC2
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
1.0
Low (n=39)
High (n=40)
0.8
0.8
0.8
0.8
OS
0.6
OS
0.6
OS
0.6
OS
0.6
0.4
0.4
0.4
0.4
0.2
0.2
0.2
0.2
HR=0.799
Cox P=6.93e-05
HR=1.72
Cox P=5.11e-07
HR=3.583
HR=1.568
0.0
Log-Rank P=3.18e-06
0.0
Log-Rank P=4.54e-07
0.0
Cox P=6.09e-07
Log-Rank P=1.11e-06
0.0
Cox P=1.76e-06
Log-Rank P=6.96e-05
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
0
2
Time (years)
4
6
8
10
12
ZNF282
1.0
Low (n=39)
High (n=40)
0.8
OS
0.6
0.4
0.2
HR=7.799
0.0
Cox P=9e-05
Log-Rank P=1.57e-05
0
2
Time (years)
4
6
8
10
12
Figure 3 Survival curves of ACC patients demarcated on the basis of TF risk score. The median risk score was used as the cutoff value. ACC, Adrenocortical carcinoma; OS, overall survival. Full-size
DOI: 10.7717/peerj.12433/fig-3
A
5
Low
High
Risk score
0
Alive
Dead
“?
0
20
40
60
80
Patients
B
Risk
CREB3L3
NROB1
Low
CENPA
High
FOXM1
E2F2
4
MYBL2
2
HOXC11
ZIC2
0
ZNF282
DNMT 1
-2
TCF3
4
ELK4
KLF6
C
Low risk (n=40)
D
1.0
1.0
High risk (n=39)
0.8
Overall survival
0.8
HR=16.95(5.02-57.2)
0.6
Cox P=5.11e-06
Sensitivity
0.6
Log Rank P=2.09e-09
0.4
0.4
0.2
10-year AUC = 0.929
0.2
5-year AUC = 0.927
4-year AUC = 0.941
3-year AUC = 0.959
0.0
0.0
2-year AUC = 0.96
1-year AUC = 0.922
0
2
4
6
8
10
12
14
0.0
0.2
0.4
0.6
0.8
1.0
Years
1-Specificity
Full-size DOI: 10.7717/peerj.12433/fig-4
The risk score is an independent prognostic factor of ACC
The prognostic relevance of the 13-TF signature was further validated by multivariate Cox regression analysis after normalizing for variables including age, gender, and pathological stage. In both the training and validation ACC cohorts, the 13-TF risk score was an independent prognostic factor (Table 2). However, no significant correlation was seen between OS and age, gender, or pathological stage. The high- and low-risk groups of
A
B
Low risk (n=11)
1.0
1.0
High risk (n=11)
HR=16.57(2.02-135.76)
Cox P=8.89e-03
0.8
Overall survival
0.8
Log Rank P=6.8e-04
0.6
Sensitivity
0.6
0.4
0.4
0.2
0.2
10-year AUC = 0.934
5-year AUC = 0.849
4-year AUC = 0.849
3-year AUC = 0.843
0.0
0.0
2-year AUC = 0.763
1-year AUC = 0.754
0
5
10
15
0.0
0.2
0.4
0.6
0.8
1.0
Years
1-Specificity
C
8
6
Log2(TPM+1)
4
2
0
CENPA
CREB3L3
DNMT1
E2F2
ELK4
FOXM1
HOXC11
KLF6
MYBL2
NROB1
TCF3
ZIC2
ZNF282
both training and validation datasets were further divided into subgroups based on age (≤50 vs >50 years), gender (male vs female), and pathological stage (I-II vs III-IV). As shown in the Figs. 6A and 6B, patients in the high-risk group had poor prognosis and significantly shorter survival compared to those in the low-risk group regardless of other variables. Thus, the 13-TF signature is an independent prognostic predictor of ACC.
The 13-TF signature is associated with cancer-related functions
GSEA results showed that four hallmarks including G2M_CHECKPOINT (P = 0.021), E2F_TARGETS (P=0.023), SPERMATOGENESIS (P= 0.046), and MITOTIC_SPINDLE (P= 0.048) were significantly enriched in high-risk patients, suggesting a mechanistic basis of the prognostic role of the 13-TF signature in ACC (Fig. 7).
| Table 2 Cox regression analysis of the correlation between clinicopathological factors and overall survival of ACC patients. | |||||||
|---|---|---|---|---|---|---|---|
| Variables | Group | Patients | Univariate analysis | Patients (N) | Multivariate analysis | ||
| (N) | HR (95% CI) | P | HR (95% CI) | P | |||
| TCGA | |||||||
| Risk score | Low/High | 40/39 | 16.95 [5.02-57.2] | 5.11E-06 | 39/38 | 22.59 [4.55-112.22] | 1.38E-04 |
| Age | <= 50/>50 | 41/38 | 1.8 [0.85-3.82] | 1.27E-01 | 40/37 | 1.47 [0.63-3.45] | 3.73E-01 |
| Gender | F/M | 48/31 | 1.00 [0.47-2.14] | 9.99E-01 | 48/29 | 1.96 [0.79-4.84] | 1.46E-01 |
| Stage | I-II/III-IV | 46/31 | 6.48 [2.71-15.5] | 2.72E-05 | 46/31 | 1.78 [0.72-4.39] | 2.13E-01 |
| GSE19776 | |||||||
| Risk score | Low/High | 11/11 | 16.57 [2.02-135.76] | 8.89E-03 | 11/11 | 18.85 [2.09-170.08] | 8.89E-03 |
| Age | <= 50/>50 | 15/7 | 1.72 [0.60-4.88] | 3.11E-01 | 15/7 | 1.23 [0.38-3.97] | 7.33E-01 |
| Gender | F/M | 11/11 | 1.26 [0.46-3.44] | 6.50E-01 | 11/11 | 0.77 [0.26-2.26] | 6.32E-01 |
| Stage | I-II/III-IV | 14/8 | 1.21 [0.44-3.29] | 7.12E-01 | 14/8 | 1.46 [0.51-4.14] | 4.81E-01 |
Note:
HR, Hazard ratio; CI, Confidence interval.
DISCUSSION
ACC is a rare endocrine cancer with limited therapeutic options and poor clinical outcomes. Studies have increasingly identified molecular diagnostic and prognostic signatures of various cancers by screening multiple databases via high-throughput technologies such as microarrays and next-generation sequencing. Transcription factors (TFs) are often aberrantly expressed in tumors, and correlated with cancer prognosis (Lutz, Leon & Eilers, 2002; Ozaki & Nakagawara, 2011). Recent studies have identified specific TFs that are independent prognostic factors in various cancers (Jiang et al., 2010; Span et al., 2002; Wolf et al., 2005). The zinc-finger transcription factor Snail is associated with decreased survival of ACC patients and a higher risk of distant metastasis (Waldmann et al., 2008). However, a TF-related prognostic signature has not yet been identified for ACC.
We analyzed the gene expression profiles of ACC patients deposited in GEO and TCGA databases and constructed a prognostic signature of 13 TFs, including CREB3L3, NR0B1, CENPA, FOXM1, E2F2, MYBL2, HOXC11, ZIC2, ZNF282, DNMT1, TCF3, ELK4, and KLF6. NR0B1, also known as DAX1, is an atypical orphan nuclear hormone receptor that is expressed in the adrenal glands and at all levels of the hypothalamusry- gonadal axis (Guo et al., 1996). Previous studies have shown that NR0B1 was associated with a variety of cancers, although its role in promoting or suppressing tumors is not consistent (Ranhotra, 2013). NR0B1 silencing decreased the in vitro invasiveness of the lung cancer cell line A549 and inhibited xenograft growth without affecting cell proliferation (Oda et al., 2009). In addition, NR0B1 is overexpressed in cervical cancer and promotes cancer cell proliferation via the Wnt/ß-catenin pathway (Liu et al., 2018). It is also overexpressed in pediatric adrenocortical tumors (de Sousa et al., 2015), ovarian cancer (Abd-Elaziz et al., 2003), breast cancer (Conde et al., 2004), endometrial cancer (Saito et al., 2005), and prostate cancer (Tong et al., 2014). Interestingly, NR0B1 was down-regulated in hepatocellular carcinoma tissues and cell lines and overexpression of
A
Age ⇐ 50 & Low risk (n=22)
Female & Low risk (n=21)
Stage I-II & Low risk (n=33)
1.0
Age ⇐ 50 & High risk (n=19)
1.0
Female & High risk (n=27)
1.0
Stage I-II & High risk (n=13)
Log-Rank P=1.27e-04
Log-Rank P=2.5e-04
Log-Rank P=3.37e-06
Overall survival
0.8
Overall survival
0.8
Overall survival
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.0
0
2
4
6
8
10
12
14
0
2
4
6
8
10
12
14
0
2
4
6
8
10
12
14
Years
Years
Years
Age>50 & Low risk (n=18)
Age>50 & High risk (n=20)
Male & Low risk (n=19)
Male & High risk (n=12)
Stage III-IV & Low risk (n=6)
1.0
1.0
1.0
Stage III-IV & High risk (n=25)
Log-Rank P=9.97e-06
Log-Rank P=0.02
Overall survival
0.8
Overall survival
0.8
Log-Rank P=6.89e-08
Overall survival
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.0
0
2
4
6
8
10
12
14
0
2
4
6
8
10
0
2
4
6
8
10
12
Years
Years
Years
B
Age ⇐ 50 & Low risk (n=5)
Female & Low risk (n=7)
Stage I-II & Low risk (n=4) Stage I-II & High risk (n=4)
1.0
Age ⇐ 50 & High risk (n=2)
1.0
Female & High risk (n=4)
1.0
Log-Rank P=8.23e-03
Log-Rank P=0.04
Overall survival
0.8
Overall survival
0.8
Log-Rank P=0.012
Overall survival
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.0
0
5
10
15
0
5
10
15
0
5
10
15
Years
Years
Years
Age>50 & Low risk (n=6)
Age>50 & High risk (n=9)
Male & Low risk (n=4)
Stage III-IV & Low risk (n=7)
1.0
1.0
Male & High risk (n=7)
1.0
Stage III-IV & High risk (n=7)
Log-Rank P=0.038
Log-Rank P=9.75e-03
Overall survival
0.8
Overall survival
0.8
Log-Rank P=0.025
Overall survival
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.0
0
5
10
15
0
5
10
15
0
5
10
15
Years
Years
Years
Full-size DOI: 10.7717/peerj.12433/fig-6
Enrichment plot: HALLMARK_G2M_CHECKPOINT
0.7
Enrichment score (ES)
P=0.021
0.6
0.5
0.4
0.3
0.2
0.1
0.0
Ranked list metric (Signal2 Noise)
1.0
‘High’ (positively correlated)
0.5
0.0
Zero cross at 10200
-0.5
-1.0
‘Low’ (negatively correlated)
0
2,500
5,000
7,500
10,000
12,500
15,000
17,500
20,000
Rank in Ordered Dataset
Enrichment profile - Hits
Ranking metric scores
Enrichment plot: HALLMARK_E2F_TARGETS
0.7
Enrichment score (ES)
P=0.023
0.6
0.5
0.4
0.3
0.2
0.1
0.0
Ranked list metric (Signal2Noise)
1.0
‘High’ (positively correlated)
0.5
0.0
Zero cross at 10200
-0.5
-1.0
‘Low’ (negatively correlated)
0
2,500
5,000
7,500
10,000
12,500
15,000
17,500
20,000
Rank in Ordered Dataset
Enrichment profile - Hits
Ranking metric scores
Enrichment plot: HALLMARK_SPERMATOGENESIS
0.40
Enrichment score (ES)
P=0.046
0.35
0.30
0.25
0.20
0.15
0.10
0.05
0.00
-0.05
Ranked list metric (Signal2Noise)
1.0
‘High’ (positively correlated)
0.5
0.0
Zero cross at 10200
-0.5
-1.0
‘Low’ (negatively correlated)
0
2,500
5,000
7,500
10,000
12,500
15,000
17,500
20,000
Rank in Ordered Dataset
Enrichment profile - Hits
Ranking metric scores
Enrichment plot: HALLMARK_MITOTIC_SPINDLE
0.45
Enrichment score (ES)
0.40
P=0.048
0.35
0.30
0.25
0.20
0.15
0.10
0.05
0.00
Ranked list metric (Signal2 Noise)
1.0
‘High’ (positively correlated)
0.5
0.0
Zero cross at 10200
-0.5
-1.0
‘Low’ (negatively correlated)
0
2,500
5,000
7,500
10,000
12,500
15,000
17,500
20,000
Rank in Ordered Dataset
Enrichment profile - Hits
Ranking metric scores
Full-size ☒ DOI: 10.7717/peerj.12433/fig-7
NR0B1 could inhibit cell proliferation (Jiang et al., 2014), indicating that it may be a tumor suppressor in hepatocellular carcinoma. CREB3L3, also called CREB-H, was originally identified as a TF specifically expressed in the liver (Omori et al., 2001). The role of CREB3L3 in the liver is mainly related to triglyceride metabolism (Lee et al., 2011).
A previous study has shown that loss of CREB3L3 function in hepatocellular carcinoma might contribute to the occurrence and/or progression of cancer (Chin et al., 2005). These results were consistent with our study that low expression of NR0B1 and CREB3L3 was associated with a worse prognosis of ACC patients. However, their potential roles in ACC have not been investigated.
CENPA is a histone-H3 variant that regulates cell division by establishing kinetochore assembly and ensuring proper centromere segregation and is associated with cancer progression (Athwal et al., 2015; Vardabasso et al., 2014). CENPA expression level is a potential biomarker of poor prognosis in cancer patients (Sun et al., 2016). FOXM1 and E2F2 are the upstream regulators of CENPA and play critical roles in cell cycle progression and tumorigenesis (Nakahata et al., 2014; Wang et al., 2011). Both TFs can potentially bind to the CENPA promoter sequence indicating that they regulate CENPA transcription (Grant et al., 2013). Previous studies have correlated CENPA with poor prognosis in ACC patients and identified E2F2 as an ACC-related TF (Xia et al., 2019; Zhang et al., 2013). Although few studies have associated CREB3L3, MYBL2, HOXC11, ZIC2, ZNF282, DNMT1, TCF2, FLK4, and KLF6 with ACC progression (DiFeo, Martignetti & Narla, 2009; Jain, Rechache & Kebebew, 2012), there is no evidence linking FOXM1 to ACC. We found that high levels of CREB3L3 and NR0B1 were correlated with good prognosis, while that of other TFs were correlated with poor prognosis in the ACC patients. The tumor stage was correlated with the OS of patients in the training set but not in the validation set. Furthermore, the 13-TF signature was an independent prognostic factor in both TCGA and GEO datasets.
The GSEA results showed that G2M-CHECKPOINT and E2F-TARGET were significantly enriched in the high-risk group. The G2/M checkpoint is frequently impaired in cancer cells, which promotes genomic instability and tumorigenesis (Lobrich & Jeggo, 2007). Since the E2F transcription factors regulate DNA replication and are aberrantly expressed in almost all human cancers (Johnson & Schneider-Broussard, 1998), targeting E2Fs could be a generic approach in anti-cancer treatment.
To the best of our knowledge, this report is the first to investigate the cancer-specific TFs and their association with clinical outcomes in ACC patients. The 13-TF signature showed accurate predictive ability and is a promising prognostic biomarker for clinical applications. However, the in-silico results were not validated by PCR or Western blotting. Future studies should focus on validating these survival-related TFs through molecular and functional assays, and determine the mechanistic basis.
ABBREVIATIONS
| ACC | adrenocortical carcinoma |
| TF | transcription factor |
| TCGA | The Cancer Genome Atlas |
| GEO | Gene Expression Omnibus |
| ROC | Receiver Operating Characteristic |
| OS | overall survival |
| CCLE | Cancer cell line encyclopedia |
| CREB3L3 | cAMP-responsive element-binding protein 3-like 3 |
| NR0B1 | nuclear receptor subfamily 0, group B, member |
| CENPA | centromere protein-A |
| FOXM1 | Forkhead Box M1 |
| E2F2 | E2F transcription factor 2 |
| MYBL2 | v-myb myeloblastosis viral oncogene homolog (avian)-like 2 |
| HOXC11 | homeobox C11 |
| ZIC2 | Zic family member 2 |
| ZNF282 | zinc finger protein 282 |
| DNMT1 | DNA methyltransferase |
| TCF3 | transcription factor 3 |
| ELK4 | ETS transcription factor ELK4 |
| KLF6 | Krüppel-like factor 6 |
| GSEA | Gene Set Enrichment Analysis |
ADDITIONAL INFORMATION AND DECLARATIONS
Funding
The authors received no funding for this work.
Competing Interests
The authors declare that they have no competing interests.
Author Contributions
· Jianyu Zhao conceived and designed the experiments, performed the experiments, analyzed the data, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft.
· Bo Liu performed the experiments, analyzed the data, prepared figures and/or tables, and approved the final draft.
· Xiaoping Li conceived and designed the experiments, authored or reviewed drafts of the paper, and approved the final draft.
Data Availability
The following information was supplied regarding data availability: All the raw data is available at TCGA (search term: ACC) and GEO: GSE19776.
REFERENCES
Abd-Elaziz M, Akahira J, Moriya T, Suzuki T, Yaegashi N, Sasano H. 2003. Nuclear receptor DAX-1 in human common epithelial ovarian carcinoma: an independent prognostic factor of clinical outcome. Cancer Science 94(11):980-985 DOI 10.1111/j.1349-7006.2003.tb01388.x.
Arendt D, Musser JM, Baker CVH, Bergman A, Cepko C, Erwin DH, Pavlicev M, Schlosser G, Widder S, Laubichler MD, Wagner GP. 2016. The origin and evolution of cell types. Nature Reviews Genetics 17(12):744-757 DOI 10.1038/nrg.2016.127.
Athwal RK, Walkiewicz MP, Baek S, Fu S, Bui M, Camps J, Ried T, Sung MH, Dalal Y. 2015. CENP-A nucleosomes localize to transcription factor hotspots and subtelomeric sites in human cancer cells. Epigenetics & Chromatin 8(1):2 DOI 10.1186/1756-8935-8-2.
Bilimoria KY, Shen WT, Elaraj D, Bentrem DJ, Winchester DJ, Kebebew E, Sturgeon C. 2008. Adrenocortical carcinoma in the United States: treatment utilization and prognostic factors. Cancer 113(11):3130-3136 DOI 10.1002/cncr.23886.
Bouchard C, Staller P, Eilers M. 1998. Control of cell proliferation by Myc. Trends in Cell Biology 8(5):202-206 DOI 10.1016/S0962-8924(98)01251-3.
Brown TC, Murtha TD, Rubinstein JC, Korah R, Carling T. 2018. SLC12A7 alters adrenocortical carcinoma cell adhesion properties to promote an aggressive invasive behavior. Cell Communication and Signaling 16(1):27 DOI 10.1186/s12964-018-0243-0.
Chin KT, Zhou HJ, Wong CM, Lee JM, Chan CP, Qiang BQ, Yuan JG, Ng IO, Jin DY. 2005. The liver-enriched transcription factor CREB-H is a growth suppressor protein underexpressed in hepatocellular carcinoma. Nucleic Acids Research 33(6):1859-1873 DOI 10.1093/nar/gki332.
Conde I, Alfaro JM, Fraile B, Ruiz A, Paniagua R, Arenas MI. 2004. DAX-1 expression in human breast cancer: comparison with estrogen receptors ER-alpha, ER-beta and androgen receptor status. Breast Cancer Research 6(3):R140-R148 DOI 10.1186/bcr766.
de Sousa GR, Soares IC, Faria AM, Domingues VB, Wakamatsu A, Lerario AM, Alves VA, Zerbini MC, Mendonca BB, Fragoso MC, Latronico AC, Almeida MQ. 2015. DAX1 Overexpression in pediatric adrenocortical tumors: a synergic role with SF1 in tumorigenesis. Hormone and Metabolic Research 47:656-661 DOI 10.1055/s-00000025.
DiFeo A, Martignetti JA, Narla G. 2009. The role of KLF6 and its splice variants in cancer therapy. Drug Resistance Updates 12(1-2):1-7 DOI 10.1016/j.drup.2008.11.001.
Else T, Kim AC, Sabolch A, Raymond VM, Kandathil A, Caoili EM, Jolly S, Miller BS, Giordano TJ, Hammer GD. 2014. Adrenocortical carcinoma. Endocrine Reviews 35(2):282-326 DOI 10.1210/er.2013-1029.
Francis JC, Gardiner JR, Renaud Y, Chauhan R, Weinstein Y, Gomez-Sanchez C, Lefrancois-Martinez AM, Bertherat J, Val P, Swain A. 2021. HOX genes promote cell proliferation and are potential therapeutic targets in adrenocortical tumours. British Journal of Cancer 124(4):805-816 DOI 10.1038/s41416-020-01166-z.
Grant GD, Brooks L 3rd, Zhang X, Mahoney JM, Martyanov V, Wood TA, Sherlock G, Cheng C, Whitfield ML. 2013. Identification of cell cycle-regulated genes periodically expressed in U2OS cells and their regulation by FOXM1 and E2F transcription factors. Molecular Biology of the Cell 24(23):3634-3650 DOI 10.1091/mbc.e13-05-0264.
Guo W, Burris TP, Zhang YH, Huang BL, Mason J, Copeland KC, Kupfer SR, Pagon RA, McCabe ER. 1996. Genomic sequence of the DAX1 gene: an orphan nuclear receptor responsible for X-linked adrenal hypoplasia congenita and hypogonadotropic hypogonadism. The Journal of Clinical Endocrinology & Metabolism 81(7):2481-2486 DOI 10.1210/jcem.81.7.8675564.
Icard P, Goudet P, Charpenay C, Andreassian B, Carnaille B, Chapuis Y, Cougard P, Henry JF, Proye C. 2001. Adrenocortical carcinomas: surgical trends and results of a 253-patient series from the French Association of Endocrine Surgeons study group. World Journal of Surgery 25(7):891-897 DOI 10.1007/s00268-001-0047-y.
Jain M, Rechache N, Kebebew E. 2012. Molecular markers of adrenocortical tumors. Journal of Surgical Oncology 106(5):549-556 DOI 10.1002/jso.23119.
Jiang J, Tang Y, Zhu G, Zheng M, Yang J, Liang X. 2010. Correlation between transcription factor Snail1 expression and prognosis in adenoid cystic carcinoma of salivary gland. Oral Surgery,
Oral Medicine, Oral Pathology, Oral Radiology, and Endodontology 110(6):764-769 DOI 10.1016/j.tripleo.2010.06.015.
Jiang HL, Xu D, Yu H, Ma X, Lin GF, Ma DY, Jin JZ. 2014. DAX-1 inhibits hepatocellular carcinoma proliferation by inhibiting beta-catenin transcriptional activity. Cellular Physiology and Biochemistry 34(3):734-742 DOI 10.1159/000363038.
Johnson DG, Schneider-Broussard R. 1998. Role of E2F in cell cycle control and cancer. Frontiers in Bioscience 3(4):d447-d448 DOI 10.2741/A291.
Koschker AC, Fassnacht M, Hahner S, Weismann D, Allolio B. 2006. Adrenocortical carcinoma- improving patient care by establishing new structures. Experimental and Clinical Endocrinology & Diabetes 114(2):45-51 DOI 10.1055/s-2006-923808.
Lambert SA, Jolma A, Campitelli LF, Das PK, Yin Y, Albu M, Chen X, Taipale J, Hughes TR, Weirauch MT. 2018. The human transcription factors. Cell 172(4):650-665 DOI 10.1016/j.cell.2018.01.029.
Latchman DS. 1997. Transcription factors: an overview. The International Journal of Biochemistry & Cell Biology 29(12):1305-1312 DOI 10.1016/S1357-2725(97)00085-X.
Lee JH, Giannikopoulos P, Duncan SA, Wang J, Johansen CT, Brown JD, Plutzky J, Hegele RA, Glimcher LH, Lee AH. 2011. The transcription factor cyclic AMP-responsive element-binding protein H regulates triglyceride metabolism. Nature Medicine 17(7):812-815 DOI 10.1038/nm.2347.
Liu XF, Li XY, Zheng PS, Yang WT. 2018. DAX1 promotes cervical cancer cell growth and tumorigenicity through activation of Wnt/beta-catenin pathway via GSK3beta. Cell Death & Disease 9(3):339 DOI 10.1038/s41419-018-0359-6.
Lobrich M, Jeggo PA. 2007. The impact of a negligent G2/M checkpoint on genomic instability and cancer induction. Nature Reviews Cancer 7(11):861-869 DOI 10.1038/nrc2248.
Lutz W, Leon J, Eilers M. 2002. Contributions of Myc to tumorigenesis. Biochimica et Biophysica Acta 1602(1):61-71 DOI 10.1016/s0304-419x(02)00036-7.
Nakahata AM, Suzuki DE, Rodini CO, Fiuza ML, Okamoto OK. 2014. RNAi-mediated knockdown of E2F2 inhibits tumorigenicity of human glioblastoma cells. Oncology Letters 8(4):1487-1491 DOI 10.3892/ol.2014.2369.
Oda T, Tian T, Inoue M, Ikeda J, Qiu Y, Okumura M, Aozasa K, Morii E. 2009. Tumorigenic role of orphan nuclear receptor NR0B1 in lung adenocarcinoma. American Journal of Pathology 175(3):1235-1245 DOI 10.2353/ajpath.2009.090010.
Omori Y, Imai J, Watanabe M, Komatsu T, Suzuki Y, Kataoka K, Watanabe S, Tanigami A, Sugano S. 2001. CREB-H: a novel mammalian transcription factor belonging to the CREB/ATF family and functioning via the box-B element with a liver-specific expression. Nucleic Acids Research 29(10):2154-2162 DOI 10.1093/nar/29.10.2154.
Ozaki T, Nakagawara A. 2011. p53: the attractive tumor suppressor in the cancer research field. Journal of Biomedicine and Biotechnology 2011(5):603925 DOI 10.1155/2011/603925.
Parviainen H, Schrade A, Kiiveri S, Prunskaite-Hyyrylainen R, Haglund C, Vainio S, Wilson DB, Arola J, Heikinheimo M. 2013a. Expression of Wnt and TGF-beta pathway components and key adrenal transcription factors in adrenocortical tumors: association to carcinoma aggressiveness. Pathology Research and Practice 209(8):503-509 DOI 10.1016/j.prp.2013.06.002.
Parviainen H, Schrade A, Kiiveri S, Prunskaite-Hyyryläinen R, Heikinheimo M. 2013b. Expression of Wnt and TGF-ß pathway components and key adrenal transcription factors in adrenocortical tumors: association to carcinoma aggressiveness. Pathology-Research and Practice 209(8):503-509 DOI 10.1016/j.prp.2013.06.002.
Ranhotra HS. 2013. The orphan nuclear receptors in cancer and diabetes. Journal of Receptors and Signal Transduction 33(4):207-212 DOI 10.3109/10799893.2013.781624.
Robin X, Turck N, Hainard A, Tiberti N, Lisacek F, Sanchez JC, Muller M. 2011. pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinformatics 12(1):77 DOI 10.1186/1471-2105-12-77.
Rodgers SE, Evans DB, Lee JE, Perrier ND. 2006. Adrenocortical carcinoma. Surgical Oncology Clinics of North America 15(3):535-553 DOI 10.1016/j.soc.2006.05.005.
Saito S, Ito K, Suzuki T, Utsunomiya H, Akahira J, Sugihashi Y, Niikura H, Okamura K, Yaegashi N, Sasano H. 2005. Orphan nuclear receptor DAX-1 in human endometrium and its disorders. Cancer Science 96(10):645-652 DOI 10.1111/j.1349-7006.2005.00101.x.
Shaulian E, Karin M. 2002. AP-1 as a regulator of cell life and death. Nature Cell Biology 4(5):E131-E136 DOI 10.1038/ncb0502-e131.
Span PN, Manders P, Heuvel JJ, Thomas CM, Bosch RR, Beex LV, Sweep CG. 2002. Expression of the transcription factor Ets-1 is an independent prognostic marker for relapse-free survival in breast cancer. Oncogene 21(55):8506-8509 DOI 10.1038/sj.onc.1206040.
Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, Paulovich A, Pomeroy SL, Golub TR, Lander ES, Mesirov JP. 2005. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proceedings of the National Academy of Sciences of the United States of America 102(43):15545-15550 DOI 10.1073/pnas.0506580102.
Sun X, Clermont PL, Jiao W, Helgason CD, Gout PW, Wang Y, Qu S. 2016. Elevated expression of the centromere protein-A(CENP-A)-encoding gene as a prognostic and predictive biomarker in human cancers. International Journal of Cancer 139(4):899-907 DOI 10.1002/ijc.30133.
Therneau TM, Lumley T. 2015. Package ‘survival’[J]. R Top Doc 128(10):28-33.
Tong SJ, Liu J, Wang X, Qu LX. 2014. microRNA-181 promotes prostate cancer cell proliferation by regulating DAX-1 expression. Experimental and Therapeutic Medicine 8(4):1296-1300 DOI 10.3892/etm.2014.1846.
Vardabasso C, Hasson D, Ratnakumar K, Chung CY, Duarte LF, Bernstein E. 2014. Histone variants: emerging players in cancer biology. Cellular and Molecular Life Sciences 71(3):379-404 DOI 10.1007/s00018-013-1343-z.
Waldmann J, Feldmann G, Slater EP, Langer P, Buchholz M, Ramaswamy A, Saeger W, Rothmund M, Fendrich V. 2008. Expression of the zinc-finger transcription factor Snail in adrenocortical carcinoma is associated with decreased survival. Br J Cancer 99(11):1900-1907 DOI 10.1038/sj.bjc.6604755.
Wang Z, Park HJ, Carr JR, Chen YJ, Zheng Y, Li J, Tyner AL, Costa RH, Bagchi S, Raychaudhuri P. 2011. FoxM1 in tumorigenicity of the neuroblastoma cells and renewal of the neural progenitors. Cancer Research 71(12):4292-4302 DOI 10.1158/0008-5472.CAN-10-4087.
Wolf D, Wolf AM, Rumpold H, Fiegl H, Zeimet AG, Muller-Holzner E, Deibl M, Gastl G, Gunsilius E, Marth C. 2005. The expression of the regulatory T cell-specific forkhead box transcription factor FoxP3 is associated with poor prognosis in ovarian cancer. Clinical Cancer Research 11(23):8326-8331 DOI 10.1158/1078-0432.CCR-05-1244.
Xia WX, Yu Q, Li GH, Liu YW, Xiao FH, Yang LQ, Rahman ZU, Wang HT, Kong QP. 2019. Identification of four hub genes associated with adrenocortical carcinoma progression by WGCNA. PeerJ 7(9):e6555 DOI 10.7717/peerj.6555.
Zhang B, Xu ZW, Wang KH, Lu TC, Du Y. 2013. Complex regulatory network of microRNAs, transcription factors, gene alterations in adrenocortical cancer. Asian Pacific Journal of Cancer Prevention 14(4):2265-2268 DOI 10.7314/APJCP.2013.14.4.2265.