MEDICAL SCIENCE MONITOR
| Received: | 2020.09.13 |
| Accepted: | 2020.12.15 |
| Available online: | 2020.12.31 |
| Published: | 2021.03.05 |
Identification of NDRG Family Member 4 (NDRG4) and CDC28 Protein Kinase Regulatory Subunit 2 (CKS2) as Key Prognostic Genes in Adrenocortical Carcinoma by Transcriptomic Analysis
Department of Urology, Shanghai Yangsi Hospital, Shanghai, P.R. China
Authors’ Contribution: Study Design A Data Collection B
Statistical Analysis C Data Interpretation D Manuscript Preparation E Literature Search F Funds Collection G
BE Zhengqing Yang
CD Hui Cheng
F
Yazhou Zhang
AG Yan Zhou
Corresponding Author: Source of support:
Yan Zhou, e-mail: zhou.yanys@hotmail.com Departmental sources
Background: Adrenocortical carcinoma (ACC) is an aggressive cancer with heterogeneous outcomes. In this study, we aimed to investigate genomic and prognostic features of ACC.
Material/Methods: Clinical, pathologic, and transcriptomic data from 2 independent datasets derived from ACC samples (TCGA- ACC dataset, GEO-GSE76021 dataset) were collected. Weighted gene co-expression network analysis (WGCNA) and survival analyses were performed to identify prognostic genes. Pathway analysis was performed for mech- anistic analysis. xCell deconvolution was performed for tumor microenvironment analysis.
Results: In the TCGA-ACC cohort, WGCNA identified a prognostic module of 5408 genes. Differential expression analysis identified 1969 genes that differed in expression level between long-term and short-term survivors. Univariate Cox regression model analysis identified 8393 genes with prognostic value. The intersection of these gene sets included 820 prognostic genes. Similar protocols were performed for the GSE76021 dataset, and 5 candidate genes were identified. Further intersection of these genes finally identified NDRG4 and CKS2 as key prognos- tic genes. Multivariate Cox regression model analysis validated the prognostic value of NDRG4 (HR=0.61, 95% CI 0.46-0.80) and CKS2 (HR=2.52, 95% CI 1.38-4.60). Moreover, NDRG4 and CKS2 expression predicted surviv- al in patients treated with mitotane (P<0.001). Further mechanism exploration found an association between CKS2 and DNA mismatch repair pathways. Moreover, NDRG4 positively correlated with CD8+ T cell infiltration, while CKS2 negatively correlated with it.
Conclusions:
We identified NDRG4 and CKS2 expression as key prognostic genes in ACC, which may help in risk stratifica- tion of ACC. Moreover, a close relationship was found between CKS2 and mismatch repair pathways. Moreover, immune cell infiltration differed according to NDRG4 and CKS2 expression.
Keywords: Adrenocortical Carcinoma . Biological Markers . Prognosis . Transcriptome
Full-text PDF:
https://www.medscimonit.com/abstract/index/idArt/928523
3256 1 7
2] 41
Background
Adrenocortical carcinoma (ACC) is an aggressive cancer orig- inating in the cortex of the adrenal gland [1]. Almost half of ACC patients will experience recurrence or metastasis [2], re- sulting in a relatively low 5-year overall survival (OS) rate, rang- ing from 16% to 60% [3,4]. Traditional prognostic factors for ACC include tumor stage [5], resection status [6], and Ki67 in- dex of tumor proliferation [7]. However, these traditional fac- tors still have limitations. Recent studies have identified nov- el molecular biomarkers for predicting ACC survival, such as a 5-gene (TOP2A, NDC80, CEP55, CDKN3, and CDK1) signature, as well as DNA damage and cell cycle pathways [8-12], indi- cating that the molecular pathogenesis of ACC may also be involved in disease progression and survival.
Advances in molecular biology have uncovered part of the pathogenesis mechanism of ACC [13-17], but the precise mech- anisms are still far from known. Recently, genomic approach- es have provided a more comprehensive landscape for ACC pathogenesis [18,19]. Assié et al were the first to apply exome sequencing and SNP array analysis in a relatively large num- ber of ACC patients [19]. They identified 2 distinct molecular ACC subgroups with opposite outcome [19]. Another study by Zheng et al applied a comprehensive pan-genomic analysis of
ACC patients, which identified PRKAR1A, RPL22, TERF2, CCNE1, and NF1 as new driver genes of ACC [18]. They also divided ACC into 3 subtypes according to a DNA-methylation signature, which may be helpful in the risk stratification of ACC [18]. These high-throughput genomic studies helped to expand knowledge boundaries with regard to ACC development.
Mitotane (1,1-(dichlorodiphenyl)-2,2-dichloroethane) is a ste- roidogenesis inhibitor and cytostatic antineoplastic medication for ACC treatment, which is recommended as systemic therapy for metastatic or unresectable ACC, as well as adjuvant ther- apy for high-risk postoperative patients [20]. However, not all patients respond to mitotane. Recently, the pharmacokinetics of mitotane has been found to be a predictive factor in ACC prognosis [21,22]. A recent study also identified the correla- tion between SOAT1 expression and response to mitotane in ACC patients [23]. However, no robust predictive biomarker has been applied in the clinic.
In the present study, by re-evaluation of 2 open-access tran- scriptomics datasets, we aimed to identify key prognostic genes in ACC, which may be predictive of the efficacy of mi- totane therapy as well. Moreover, we attempted to discover the molecular mechanisms underlying the functions of these prognostic genes.
12.5
Down Stable
-
10.0
Up
TCGA-ACC
-log10 (p-value)
7.5
GSE76021
5.0
.
25
0.0
WGCNA
Cox regression
DEGS
-2.5
log2 (fold change)
0.0
2.5
1.00
NDRGA level +high+low
1.00
CK52 level
-high +-low
100
Survival probabil ity
Survival probabil ity
CKS2 NDRG4
50
50
80
1.
0.25
p -< 0.001
0.25
p<0.001
60
0.00
0 1 2 3 4 5 6 7 8 9101112 Tomat (ytum)
0.00
Sensivity
01234567231711
40
1.00
NDAG4 +high+lov
1.00
CK52 +high +low
Survival probability
Survival probabil ity
D.TS
Multivariate cox
Kaplan-Meier analysis
Receiver operating characteristic analysis
20
Stage
AUC
Pvalue
0.2001
50
0.50
NORG4 Risk score
0.906
0.974
0.1182
0.0016
0.75
0
p<0.001
0.25
p<0.001
0
20
40
60
80
100
0.00
1123456713101112
0.00
8123456713112
100-specificity
Pathway analysis
4
2
0
-
Immunologic analysis
8
2
2
2
0
2
6
0
6
8
0
6
8
0
6
8
4
4
4
4
BRCA1
2
8
2
2
2
0
6
0
6
8
0
6
8
0
6
8
4
4
4
4
| Male (%) | Female (%) | Total (%) | |
|---|---|---|---|
| Sex | 31 (39.2%) | 48 (60.8%) | 79 (100%) |
| Age (years) | 48.74±15.67 | 45.38+15.69 | 46.70±15.67 |
| Status | |||
| Dead | 11 (13.9%) | 17 (21.5%) | 28 (35.4%) |
| Alive | 20 (25.3%) | 31 (39.2%) | 51 (64.6%) |
| T stage | |||
| 1 | 3 (3.8%) | 6 (7.6%) | 9 (11.4%) |
| 2 | 15 (19.0%) | 27 (34.2%) | 42 (53.2%) |
| 3 | 4 (5.0%) | 4 (5.0%) | 8 (10.1%) |
| 4 | 7 (8.9%) | 11 (13.9%) | 18 (22.8%) |
| N stage | |||
| 0 | 27 (34.2%) | 41 (51.9%) | 68 (86.1%) |
| 1 | 2 (2.5%) | 7 (8.9%) | 9 (11.4%) |
| M stage | |||
| 0 | 24 (30.4%) | 38 (48.1%) | 62 (78.5%) |
| 1 | 5 (6.3%) | 10 (12.7%) | 15 (19.0%) |
| Stage | |||
| I | 3 (3.8%) | 6 (7.6%) | 9 (11.4%) |
| II | 15 (19.0%) | 22 (27.8%) | 37 (46.8%) |
| III | 6 (7.6%) | 10 (12.7%) | 16 (20.3%) |
| 6IV | 5 (6.3%) | 10 (12.7%) | 15 (19.0%) |
| Events | |||
| Local recurrence | 5 (6.3%) | 5 (6.3%) | 10 (12.7%) |
| Distant metastasis | 6 (7.6%) | 20 (25.3%) | 26 (32.9%) |
| New primary tumor | 1 (1.3%) | 0 (0.0%) | 1 (1.3%) |
Material and Methods
Study Populations and Data Acquisition
The study design is illustrated in Figure 1. Two independent datasets, the TCGA-ACC dataset and the GEO-GSE76021 da- taset, were downloaded for the study. The TCGA-ACC cohort included 79 patients (Table 1) and the GEO-GSE76021 co- hort included 29 patients. For the TCGA-ACC cohort, clinical, pathologic, and transcriptomic data were downloaded from the UCSC Xena (https://xena.ucsc.edu/) website. For the GEO- GSE76021 cohort, clinical, pathologic, and transcriptomic data were downloaded from the Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo/) website. R software (version 4.0.2; https://www.r-project.org/) and Bioconductor packages (http://www.bioconductor.org/) were used for data analysis.
The “sva” package [24] was used for the normalization of the GSE76021 dataset. For duplicates, RNA expression was aver- aged, and genes with relatively low abundance were discard- ed. We followed the access policies of the TCGA and GEO da- tabases strictly during the process of the study.
Weighted Gene Co-expression Network Analysis
Weighted gene co-expression network analysis (WGCNA) was used as a systematic method to describe the pattern of gene correlation, as previously described by Langfelder et al [25]. We set OS as the phenotype for WGCNA analysis. Potential biomarker genes in ACC were identified through association between gene sets and OS. In the current study, gene expres- sion matrices for the TCGA-ACC cohort and the GEO-GSE76021 cohort were used as input data. The “WGCNA” package in R
software was used for WGCNA processing. An adjacency ma- trix and a topological overlap matrix were constructed. After that, we constructed a gene dendrogram, and performed mod- ule identification with a dynamic tree cut. The correlation be- tween OS and genomic data was then calculated. The most predominant module associated with OS was selected as the key module, and genes in the key module were used in the following analysis.
Survival Analysis
Univariate Cox analysis was used to screen survival-related genes. Hazard ratio (HR) was used to classify genes as either protective or deleterious. These genes were selected as po- tential prognostic genes. A multivariate Cox regression mod- el including gender, age, tumor stage, NDRG Family Member 4 (NDRG4), and CDC28 Protein Kinase Regulatory Subunit 2 (CKS2) was used to identify independent prognostic factors for OS in ACC. Kaplan-Meier survival analysis and log-rank test were used to validate the prognostic value of target genes.
Differential Expression Analysis
For differential expression analysis, long-term survivors were defined as those with OS >5 years, and short-term survivors were defined as those with OS <3 years. Differentially expressed genes (DEGs) were determined by the “Limma” package mod- erated t test, and genes with P value <0.001 and fold change >2 were regarded as DEGs [26].
Establishment of regression Model and Construction Of Risk Score
The risk score for predicting postoperative survival was calcu- lated based on a linear combination of tumor stage and gene expression in the TCGA-ACC cohort by logistic regression. The following calculation formula was used for the analysis:
Risk score= Σ Ni=1 (Expi*Coei)
N, Coei, and Expi represented gene number, coefficient value, and level of gene expression, respectively. Receiver operating characteristic (ROC) analysis was used to evaluate the accu- racy of the prognostic model.
Gene Set Enrichment Analysis
Gene set enrichment analysis (GSEA) was performed us- ing the molecular signatures database (MSigDB) (http:// www.broadinstitute.org/gsea/msigdb). The “clusterProfiler” R package was used during the process. Significantly differ- ent pathways were defined as P< 0.05 in the gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG)
databases. Enrichment scores and P values were based on 1000 permutations.
Tumor Microenvironment Analysis
We applied xCell (http://xCell.ucsf.edu/), a gene signature-based deconvolution method, to infer 64 immune and stromal cell types for ACC microenvironment analysis [27].
Statistical Analysis
Wilcoxon signed ranks test or t test was used for analysis of continuous variables, as appropriate. Chi-square test, Fisher’s exact test, or Cochran-Mantel-Haenszel test was applied for categorical variables. Correlational analyses were performed by Spearman’s correlation test. Two-tailed P value <0.05 was considered as statistically significant for all statistical analyses.
Results
Identification of Potential Prognostic Genes in the TCGA-ACC Cohort
Figure 1 shows the study design. We first performed weight- ed correlation network analysis to identify prognostic genes in the TCGA-ACC cohort. OS was set as the phenotype. Expression matrices for 79 samples were used for WGCNA construction. Figure 2A shows the cluster dendrogram of the samples in the TCGA-ACC cohort. After quality assessment for each ex- pression matrix, R2=0.9 was selected to ensure a scale-free network (Figure 2B). We also set a threshold of 0.4 to merge similar modules. After that, a total of 15 modules were iden- tified (Figure 2C). Genes in the grey module were removed in the subsequent analysis. Next, we calculated eigengenes for all of the modules and clustered them based on their correla- tion (Figure 2D, 2E). A module eigengene heatmap showed that the 7 modules were mainly divided into 3 clusters (Figure 2D). A network heatmap was also performed to analyze the inter- actions between the modules (Figure 2E). The results showed that the skyblue2 module was relatively independent from the other modules, with a high-scale independence of transcrip- tomic expression. Figure 2F shows a heatmap of the correla- tion between module eigengenes and patients’ survival in the TCGA-ACC cohort. Close relationships are represented by red color. We found that the skyblue2 module was significantly correlated with OS (Figure 2F, P=9×10-8). Taken together, the skyblue2 module was selected as the main prognostic module.
We further performed differential expression analysis to de- fine 3769 DEGs between long-term survivors (OS >5 years) and short-term survivors (OS <3 years) in the TCGA-ACC co- hort (Figure 2G). Moreover, we performed a univariate Cox
A
B
D
Sample denrogram and trait heatmap
Scale independence
121314151617181820
Scale free topology model fit, signed R2
1.0
1.0
160
0.8
0.8
0.6
3
0.4
0.6
Height
140
MEp lum
MEbrovin
MEcoral2
120
TOGA-OR-A5J8-01A
TOGA-OR-ASKT-01A TOGA-OR-ASKED
ICGA-OR-ASJG-01AOHASLO-01A
TOGA-OR- ASLO-OTA
KOGA- OR-A 515-01A
ICGA-OR-ASL3-01A
0.2
0.4
2
MEantiquevihite4
MEmediumorchid
MElightsteelblue1
MEdarkred
MEskyblue1
MEI avenderblush
MEbrown4
MEbarkseagreen
MEcoral
MEpurple.
MEskyblue?
0.0
1
100
KOGA-OR-SJ8-01A KGA-P6-ASOG-01A
TOGA-OR-ASLC-01A TOG A-OR-ASK7-01A
TOGA-OR-ASK9-01A
TOGA-OR-ASLE-01A
TOGA-08-15JM-01A
5
10
15
20
80
TOGA-OR-ASLK- 01A
TEGA -OR- ASAN-01A
TOG A-OR-A511-01.
ICIGA -OR- ASJV-01A
KOGA-OH-A SJ0-01A
1.0
TOGI-OH-ASJR-01A
HUN OR ASLH- 01A
KOGA-OR-ASAL-01A
TOCA-PR- ASHB-01A
LOGAN -OR -ASKS-01A TOGA-OR-ASLMOTA
ICGA-OR-A SKY-01
GAOR ASTE OL
IEGA -OR- ASJF-01
TOGA-OR ASI&LOLA
KOGA-OR-ASKW-01
TOGA-OR-ASI7-01A
TOGA-OR-ASK2-01A
TOGA-OR -ASJE-01A
Soft threshold (power) Scale independence
3000
1
0.8
Mean connectivity
0.6
2000
0.4
OS time
1000
2
0.2
500
3
0
0.0
5 6 7 8 9 1011121314151617181820
5
10
15
20
C
Cluster dendogram
Soft threshold (power)
E
Network heatmap plot, selected genes
1.0-
0.8
Height
0.6-
0.4
0.2
Dynamic tree cut
Merged dynamic
TCGA-COX
TCGA-DEGS
F
G
H
2787
2949
0
MEplum
0.17
12.5
MEbrown
(0.3)
0.16
Down
.
820
MEcoral2
(03)
1.0
Stable
Up
1837
0
MEantiquewhite4
-0.066
10.0
MEmediumorchid
(0.7)
MElightsteelblue 1
0.19 (0.2)
MEdarkred
0.18 (0.3)
0.5
2751
-log10 (p-value)
7.5
MEskyblue1
0.15
MEbrown4
(0.3)
0.0
WGCNA-skyblue2
MEdarkseagreen4
5.0
Size of each list
0.44
MElavenderblush3
10.004)
8393
(0.2)
4196.5
8393
MEcoral1
0.47
-0.5
(0.002)
2.5
0
3769
5408
MEpurple
0.63
MEskyblue2
(9e-06]
TCGA-COX TCGA-DEXs WGCNA-skyblue2 Number of lements: specific (1) or shared by 2,3, … lists
0.72
MEgrey
(ge-08)
-0.048
(0.8)
-1.0
0.0
4786
5538
OS time
-2.5
0.0
2.5
3
2
1
log2 (fold change)
DATABASE ANALYSIS
A
B
D
Sample denrogram and trait heatmap
Scale independence
1.1
Scale free topology model fit, signed R2
0.8
7 8 91011121914151617181820
80
0.9
75
GSM 1972921
0.6
Height
GSM1972924
0.7
70
GSM1972939
GSM197233
0.4
0.5
MEnavajowhite1
65
GSM11972927
GSM11972931
GSM1972937
GSM19/2935
0.2
MEbrovinz
MEbrown
MEgreen
MEdarlorange
MEdarkturqu ob
MEdarkgry
MBightyanl
ROPA
MEblack
MEplum3
MElightsteelbluet
ME avenderblush2
MEorangered4
MEantiquewhite4
MEbluež
1
60
GSM1972926
GSM1972928
GSM 1972944
GSM1972919
GSM1972929
GSM1972938
GSM1972942
5
10
15
55
Soft threshold (power)
20
50
GSM1972920
GSM1972922
GSM1972925
GSM1972923
GSM1972940
GSM1972934
GSM1972930
GSM1972947
3000
Scale independence
1.0
45
1
0.8
Mean connectivity
2000
0.6
0.4
OS time
1000
0.2
500-
0-
0.0
5 6 7 8 9 1011121314151617181820
5
Soft threshold (power)
10
15
20
C
Cluster dendogram
E
Network heatmap plot, selected genes
1.00-
0.95-
0.90
Height
0.85
0.80
0.75
0.70-
0.65-
Dynamic tree cut
Merged dynamic
GEO-COX
GEO-DEGS
F
G
H
774
255
144
MEplum
0.21
(0.3)
Down
MEbrown
5
-0.29
MEcoral2
(0-2)
1.0
Stable
-0.17
Up
25
6
MEantiquewhite4
(944)
0_27
00.2)
MEmediumorchid
0-067
4
MElightsteelblue1
02
302
00.4)
0.5
MEdarkred
0.44
[0.02)
MEskyblue1
0.36
-log 10 (p-value)
[0.03)
MEbrown4
0.13 00.6)
WGCNA-plum3
0,34
0.0
MEdarkseagreen4
(0.1)
2
Size of each list
MElavenderblush3
-0.29 (0.2)
1059
-0.26
529.5
1095
MEcoral1
00.2)
0.076
410
338
MEpurple
(0.7)
0.5
0
-0.28
TCGA-COX
MEskyblue2
TCGA-DEXs WGCNA-skyblue2
00.2)
-0.47
0
Number of lements: specific (1) or shared by 2,3, … lists
MEgrey
[0.02)
[Be-040)
286
1220
OS time
1.0
-2
-1
0
1
2
3 (5) 2
1
log2 (fold change)
A
WGCNA-GEO
B
WGCNA-TCGA
Hazard ratio
31
199
14
TCGA-differ
1.32
11
373
Gender
(n=79)
(0.56-3.1)
H
4
0.53
2435
169
11
13
Age
(n=79)
1.01
(0.98-1.0)
0.549
19
1666
65
2
1.80
72
8
Stage
(n=79)
78
(1.14-2.8)
0.012*
11
14
GEO-differ
393
NDRG4
(n=79)
0.61 (0.46-0.8)
<0.001 ***
65
5090
485
GEO-COX
CKS2
(n=79)
2.52
(1.38-4.6)
0.003 **
Events: 27; Global p-value (Log-Rank): 3.0303e-11 aIC: 159.51; Concordance index: 0.88
TCGA-COX
0.5
1
2
5
C
NDRG4
CKS2
Size of each list
0.94
0.0092
0.0087
0.00011
8393
8393
0.23
12.5
0.71
5408
0.0043
1.9e-05
4196.5
Gene expression
10
0.093
10.0
0.24
1969
0.31
0.038
338
454
1059
0
:
WGCNA-TCGA
7.5
WGCNA-gGEO
TCGA-differ
GEO-differ
TCGA-COX
GEO-COX
5
:
Number of lements: specific (1) or shared by 2,3, … lists
5.0
2654
8568
2.5
6 (2)
5 (132)
2
1
0
T1
T2
T3
T4
T1
T2
T3
T4
4 (154)
3 (819)
regression model to define 8393 potential prognostic genes (Figure 2H). Finally, 820 intersection genes were defined be- tween the DEGs, Cox regression prognostic genes, and the WGCNA-skyblue2 module genes.
Identification of Potential Prognostic Genes in the GSE76021 Cohort
We downloaded clinicopathologic and transcriptomic data from the GEO dataset GSE76021, which consisted of data from 29 ACC patients. We also set OS as the phenotype for WGCNA. A cluster dendrogram of the samples in the GSE76021 cohort is illustrated in Figure 3A. After quality assessment for the ex- pression matrix, R2=0.9 was selected to ensure a scale-free net- work (Figure 3B). We also set a threshold of 0.4 to merge sim- ilar modules. After that, a total of 16 modules were identified (Figure 3C). Next, we calculated eigengenes for all modules and clustered them based on their correlation (Figure 3D, 3E). The modules were mainly divided into 2 clusters (Figure 3D), and the interactions among all of the modules were detect- ed by network heatmap (Figure 3E). In the end, the plum3 module was defined as relatively independent from the other
modules, with a high-scale independence of transcriptomic ex- pression. Moreover, a close relationship was found between the plum3 module and the patients’ OS in the GSE76021 co- hort (Figure 3F). Taken together, the plum3 module was se- lected as the main prognostic module in the GSE76021 cohort.
Differential expression analysis was further performed in the GSE76021 cohort, and 410 DEGs were found between long- term survivors (OS >5 years) and short-term survivors (OS <3 years) (Figure 3G). Moreover, univariate Cox analysis was per- formed, and 1059 potential prognostic genes were identified (Figure 3H). Finally, 5 intersection genes were defined between the DEGs, Cox regression prognostic genes, and WGCNA-plum3 module genes in the GSE76021 cohort.
Identification of NDRG4 and CKS2 as Prognostic Genes for ACC
Figure 4A shows a Venn diagram of the intersection of prog- nostic genes in the TCGA-ACC cohort and the GSE76021 co- hort. We next constructed a multivariate Cox regression mod- el for OS in ACC, including gender, age, tumor stage, NDRG4,
DATABASE ANALYSIS
A
TCGA-ACC cohort
B
GSE76021 dataset
NDRG4 level + High+Low
CKS2 level
+ High+Low
NDRG4 level + High+Low
CKS2 level
+High+Low
1.00
1.00
1.00-
1.00
Survival probability
0.75
Survival probability
0.75
Survival probability
0.75
Survival probability
0.75
0.50
0.50
0.50
0.50
0.25
p<0.001
0.25
p<0.001
0.25
p=0.002
0.25
p
0.001
0.00
0.00
0.00
0.00
0
T
2
3
4
5
6
7
8
9101112
0
1
2
3
4
5
6
7
8
9
10 11 12
0
2.
3
1
6
8
9 1011 12 131415 16
0
2
4
5 6
8
9 9 10
01112131415 16
Time (years)
Time (years)
Time (years)
Time (years)
C
Disease-free survival
D
Mitotane therapy
NDRG4 +High+Low
CKS2 +High+Low
NDRG4
High+Low
CKS2 +High+Low
1.00
1.00
1
1.00-
1.00-
Survival probability
0.75
Survival probability
0.75
Survival probability
0.75
Survival probability
0.75
0.50
0.50
0.50
0.50
0.25
p<0.001
0.25
p<0.001
0.25
p<0.001
0.25
p<0.001
0.00
1
2
7
8
9101112
0.00
0.00
0.00
0
3
4
5
6
0
1
2
3
4
5
6 7
8
9 1011 12
0
1
2
3
4
5
6
7
8
9101112
0
1
2
3
4 5
6
7
8
9 101112
Time (years)
Time (years)
Time (years)
Time (years)
E
3-year survival
F
5-year survival
100
100
80
80
Sensivity
60
Sensivity
60
40
40
20
AUC
P value
Stage
0.804
reference
20
AUC
P value
0.891
CKS2
0.882
0.2091
Stage
reference
0.1182
CKS2
0.892
0.9763
NDRG4
0.906
NDRG4
0.866
0.6718
0
Risk score
0.974
0.0016
0
Risk score
0.969
0.0313
0
20
40
60
80
100
0
20
40
60
80
100
100-specificity
100-specificity
A
NDRG4
B
NDRG4 low
NDRG4 high
cell_list
- KEGG_CELL_CYCLE
CHEK2
4
BADSNC
0.6
- KEGG_DNA_REPLICATION
BARD1
2
CHECKT
0
Enrichment score
- KEGG_HEDGEHOG_SIGNALING_PATHWAY
- KEGG_HOMOLOGOUS_RECOMBINATION
-2
0.4
-4
.
FANCDZ
FANCA
0.2
ELM
FANCE
FANCF
FANCL
BAD52 ATM
0.0
MRE11A
ATR
FANCM
BAD518 COK12
PALB2
BADSO
high expression
low expression
C
0.29, p = 0.0097
R =- 0.41, p= 2e-04
R = - 0.56, p= 1.2e-07
R= · 0.44, p = 5.4e-05
BRCA1
BRCA2
EXO1
MSH2
NDRG4
NDRG4
NDRG4
NDRG4
D
E
CKS2
CKS2 low
CKS2 high
0.8
cel_list
4
KEGG_DNA_REPLICATION
MLH3
PMSZ
2
Enrichment score
- KEGG_MISMATCH_REPAIR
MSH3
0.6
- KEGG_P53_SIGNALING_PATHWAY
RFCT
POLD4
0
- KEGG_SPLICEOSOME
MLH1
RPA1
-2
0.4
SSEPT
RPA3
4
RPA4
POMA
RFC4
0.2
POLD3
REC2
RPAZ
MSH2
0.0
MISHG
RFC3
RFCS
LIGT
POLD1
YCOLA CREASKK01A
TOGA OR ASI
TOO.A. OR ASKEDIN
QUIASPIORA TOGA OFLASJBL.DIA
TOGA OR ASJYCIA
CIEASJFORA
TOO.A. DR ALICIA
TOGA GRASNO DIA
TODA CILASJ7.01A
TOGA OR ASKEDIA
TOUA OR ASKCZOŁA
TOG.A. ORLASKYDIA
high expression
low expression
F
R=0.46, p = 2.3e-05
R= 0.55, p = 1.9e-07
R=0.68, p < 2.2e-16
R= 0.46, p = 2.48-05
BRCA1
BRCA2
EXO1
MSH2
CKS2
CKS2
CKS2
CKS2
DATABASE ANALYSIS
A
p=0.002
0.5
p<0.001
0.4
p<0.001
p=0.008
Fraction
0.3
p<0.001
p=0.003
0.2
P<0.001
p<0.001
p<0.001
0.1
p=0.089
p=0.006
p=0.036
p=0.001
p<0.001
p=0.633
p=0.284
p=0.577
0.0
CD4+ T-cells
CD8+ T-cells
CD8+ Tcm
CD8+ Tem
Fibroplast
MSC
Macrophages
Macrophages M1
Macrophages M2
NK cells
NKT
Plateles
Th1 cells
Th2 cells
Tregs
CDC
DC
B
R=0.47,p=1e-05
R=0.29,p=0.0092
0.20
R=0.29,p=0.0083
R=0.59,p=8.4e-09
0.15
0.10
0.06
CD8+T cells
0.15
0.10
DC
0.05
NKT
0.10
Tregs
0.04
0.05
0.05
0.02
0.00
0.00
0.00
0.00
0
2
4
6
8
0
2
4
6
8
0
2
4
6
8
0
2
4
6
8
NDRG4
NDRG4
NDRG4
NDRG4
0.15
R=0.56, p=9.2e-08
0.3
R=0.39, p=0.00044
0.075
‘R=0.23,p=0.042
0.3
R=0.42,p=0.00011
0.10
Th1 cells
Th2 cells
0.2
Macrophages
0.050
Fibroblasts
0.2
0.05
0.1
0.025
0.1
0.00
0.0
0.000
0.0
0
2
4
6
8
0
2
4
6
8
0
2
4
6
8
0
2
4
6
8
NDRG4
NDRG4
NDRG4
NDRG4
C
P<0.001
0.5
p<0.001
0.4
p=0.002
P<0.001
Fraction
0.3
0.2
p<0.001
P<0.001
p=0.002
P<0.001
p<0.001
0.1
p=0.120
p=0.016
p=0.228
p<0.001
p<0.001
p=0.301
p=0.308
p=0.002
0.0
CD4+T-cells
CD8+T-cells
CD8+ Tcm
CD8+ Tem
Fibroplast
MSC
Macrophages
Macrophages M1
Macrophages M2
NK cells
NKT
Plateles
Th1 cells
Th2 cells
Tregs
DC
CDC
4
6
8
6 8
4
6
8
4
6
8
Indexed in: [Current Contents/Clinical Medicine] [SCI Expanded] [ISI Alerting System] [ISI Journals Master List] [Index Medicus/MEDLINE] [EMBASE/Excerpta Medica] [Chemical Abstracts/CAS]
D
R=0.51,p=2e-06
0.075
0.15
R =- 0.59, p=7.7e-09
0.20
R =- 0.37,p=0.00077
R =- 0.65, p=9.3e-11
0.10
0.15
CD8+T cells
0.050
0.10
DC
0.05
NKT
0.10
Tregs
0.025
0.05
0.05
0.00
0.00
0.000
0.00
4
6
8
4
6
8
4
6
8
4
6
8
CKS2
CKS2
CKS2
CKS2
0.15
R=0.61,p=2.9e-09
0.3
R=0.55,p=1.6e-07
0.075
R=0.38,p=0.00056
0.3
R =- 0.41,p=0.00014
0.10
0.2
Macrophages
0.05
0.050
0.2
Th 1 cells
Th2 cells
Fibroblasts
0.1
0.00
0.0
0.025
0.1
0.0
-0.1
0.000
4
6
8
4
6
8
4
6
8
4
6
8
CKS2
CKS2
CKS2
CKS2
and CKS2 (Figure 4B). In this multivariate Cox regression model, NDRG4 and CKS2 expression were both independent- ly associated with the patients’ OS (NDRG4: HR 0.61, 95% CI 0.46-0.80, P<0.001; CKS2: HR 2.52, 95% CI 1.38-4.60, P=0.003; Figure 4B). Moreover, NDRG4 expression decreased when tu- mor stage increased (Figure 4C). In contrast, CKS2 expression increased when tumor stage increased. The results indicat- ed that NDRG4 and CKS2 expression are potential prognos- tic genes in ACC, and are associated with ACC progression.
We also performed Kaplan-Meier survival analysis for the ACC patients. Kaplan-Meier survival curves and log-rank analysis showed that high expression of NDRG4 was associated with better survival in both the TCGA-ACC cohort (P<0.001; Figure 5A) and the GSE76021 cohort (P=0.002; Figure 5B). On the other hand, CKS2 high expression was associated with worse sur- vival in both the TCGA-ACC cohort (P<0.001; Figure 5A) and the GSE76021 cohort (P<0.001; Figure 5B). For disease-free survival, Kaplan-Meier survival analysis also identified NDRG4 (P<0.001) and CKS2 (P<0.001) expression as prognostic factors (Figure 5C). Moreover, NDRG4 (P<0.001) and CKS2 (P<0.001) expression could also be prognostic factors for patients treat- ed with mitotane therapy (Figure 5D).
We further performed ROC analysis to identify whether CKS2 and NDRG4 could be used as predictors for 3-year survival and 5-year survival of ACC patients. We constructed a combined risk score by logistic regression modelling for ACC survival
according to tumor stage, NDRG4 expression, and CKS2 ex- pression. We found that this combined risk score could be used for prediction of 3-year survival in ACC, and this worked better than tumor stage alone (AUC=0.974 vs AUC=0.804, P=0.0016, Figure 5E). The combined risk score could also be used for pre- diction of 5-year survival, a stronger prediction than for tumor stage alone (AUC=0.969 vs AUC=0.891, P=0.0313; Figure 5F).
Biologic Pathways Related to NDRG4 and CKS2 Expression
In order to clarify biologic pathways related to NDRG4 and CKS2 expression, we further performed GSEA analysis in the TCGA- ACC cohort. We found several pathways enriched in NDRG4 low-expression samples including the cell cycle pathway, the DNA replication pathway, the hedgehog signaling pathway, and the homologous recombination pathway (Figure 6A). Moreover, expression of genes related to the homologous recombina- tion pathway was mostly elevated in the NDRG4 low-expres- sion samples (Figure 6B). In addition, a negative correlation was found between NDRG4 expression and homologous re- combination genes including BRCA1, BRCA2, EXO1, and MSH2 (Figure 6C).
We also performed GSEA analysis between low and high ex- pression levels of CKS2. We found DNA replication mismatch repair, p53 signaling pathway, and spliceosome pathway genes were enriched in CKS2 high-expression samples (Figure 6D). Moreover, expression of genes related to DNA mismatch repair
pathways was mostly elevated in the CKS2 high-expression samples (Figure 6E). In addition, a positive correlation was found between NDRG4 expression and mismatch repair genes such as BRCA1, BRCA2, EXO1, and MSH2 (Figure 6F).
Immune Cell Infiltration Associated with NDRG4 and CKS2 Expression
Immunotherapy has been recently applied in multiple cancer types. We aimed to discover whether NDRG4 and CKS2 ex- pression was associated with immune cell infiltration in ACC. By using xCell cell type enrichment analysis, we identified im- mune cell infiltration in the TCGA-ACC samples. The infiltra- tion of CD8+ T cells, dendritic cells, natural killer T cells, regu- latory T cells, macrophages, and fibroblasts was significantly more abundant in NDRG4 high-expression samples, while the infiltration of Th1 and Th2 cells was more abundant in NDRG4 low-expression samples (Figure 7A). Correlation analysis also identified a positive correlation between NDRG4 expression and infiltration of CD8+ T cells, dendritic cells, natural killer T cells, regulatory T cells, macrophages, and fibroblasts; and a negative correlation between NDRG4 expression and infiltra- tion of Th1 and Th2 cells (Figure 7B).
We also analyzed the relationship between CKS2 expression and immune cell infiltration. The infiltration of CD8+ T cells, dendritic cells, natural killer T cells, regulatory T cells, macro- phages, and fibroblasts was significantly more abundant in CKS2 low-expression samples, while the infiltration of Th1 and Th2 cells was more abundant in CKS2 high-expression sam- ples (Figure 7C). Correlation analysis also identified a nega- tive correlation between CKS2 expression and infiltration of CD8+ T cells, dendritic cells, natural killer T cells, regulatory T cells, macrophages, and fibroblasts; and a positive correlation between CKS2 expression and infiltration of Th1 and Th2 cells (Figure 7D). Our results suggested immune cell infiltration dif- fers according to NDRG4 and CKS2 expression.
Discussion
ACC is an aggressive cancer with heterogeneous outcomes, but its molecular mechanisms have not yet been clarified. Recently, genomic technologies have been applied for the discovery of mechanisms of ACC. In the present study, we applied transcrip- tional analysis of available transcriptomic data from 2 data- sets of ACC. We identified NDRG4 and CKS2 as key prognos- tic genes for ACC, with potential predictive value for mitotane therapy as well. We further identified mismatch repair path- ways as dominant biologic pathways associated with CKS2. Moreover, immune cell infiltration in ACC differed according to NDRG4 and CKS2 expression.
NDRG4 is a member of the N-myc downregulated gene family, which belongs to the alpha/beta hydrolase superfamily [28]. The protein encoded by this gene is a cytoplasmic protein that is required for cell cycle progression and survival [29], and may be involved in the regulation of mitogenic signaling in vascu- lar smooth muscles cells [30]. Research into the connection between NDRG4 and cancer is gaining more and more atten- tion, although conflicting results have been reported. NDRG4 can show either tumor-suppressive or oncogenic function, de- pending on the tumor types [31-34]. Currently, the function of NDRG4 in ACC has not been determined. Our study found a relationship between high expression of NDRG4 and favor- able survival in ACC (Figure 5A-5C), suggesting that it may be a potential tumor suppressor. Since NDRG4 is required for tu- mor cell cycle regulation [29], and cell cycle regulation plays a crucial role in ACC progression [12], we believe this could be a potential mechanism underlying the correlation between NDRG4 and ACC development. Moreover, Agosta et al (2018) found that NDRG4 may be a potential target of miR-139-5p, and that the miR-139-5p/NDRG4 pathway promotes ACC ag- gressiveness by mediating epithelial-to-mesenchymal transi- tion, which results in tumor cell invasiveness and anchorage- independent growth [35]. However, this finding needs further in vitro and in vivo validation.
CKS2 is a protein-coding gene that binds to the catalytic sub- unit of the cyclin-dependent kinases and is essential for their biological function [36]. CKS2 is involved in the process of tu- morigenesis of gastric cancer [37], prostate cancer [38], blad- der cancer [39] and others. Our study also found an associ- ation between high expression of CKS2 and poor survival in ACC (Figure 5A-5C). Moreover, we found CKS2 expression to be correlated with elevated expression of mismatch repair genes such as BRCA1, BRCA2, EXO1, and MSH2 (Figure 6D-6F). Since CKS2 overexpression has been found to be able to override the DNA damage response barrier triggered by activated oncopro- teins [40], we believe this could be a potential mechanism un- derlying the correlation between CKS2 and survival of ACC pa- tients. Moreover, through in vitro experiments, Chehade et al found that miR-7, which participates in cell cycle arrest in ACC, acts as a key regulator of CKS2 [41]. According to these stud- ies, we speculate that miRNAs might play an important role in the regulation of NDRG4 and CKS2, resulting in tumor cell proliferation and aggressiveness, as well as ACC progression. However, further in vitro and in vivo functional studies are still needed. In vitro studies needed include gene overexpression and CRISPR knock-out in tumor cells, and evaluation of tumor cell growth, invasion, and other related biological processes. In vivo studies could be performed on tumor-bearing mice, includ- ing evaluation of tumor growth and mouse survival. Moreover, the functional status of immune cells in tumor-bearing mice should also be evaluated, since our study also suggested a re- lationship between CKS2/NDRG4 and intratumoral immunity.
Mitotane is recommended as systemic therapy for metastat- ic or unresectable ACC, and also as adjuvant therapy in high- risk postoperative patients [20], but not all patients respond to this drug. The development of robust predictive biomarkers for mitotane therapy efficacy may have potential clinical sig- nificance. Our study identified NDRG4 and CKS2 as predictive biomarkers for mitotane therapy success in ACC (Figure 5D), which may be helpful in clinical management. However, further validation studies in larger, prospective cohorts are still needed.
Sequence-based analyses, including genomic, transcriptomic, and epigenomic profiling, have been applied recently in ACC studies [18]. Such studies have identified novel driver genes and molecular classifications for ACC, which may guide fur- ther functional studies [18]. Our study integrated multiple da- tasets from previous studies, which provided a more compre- hensive view of the pathogenesis of ACC.
The present study still has some limitations. The limited sam- ple size is the major limitation. Further validation studies in larger, prospective cohorts are needed. Moreover, the meth- odology of the current study could not clarify the underlying molecular mechanism of the 2 prognostic genes. Future func- tional and mechanical studies are needed to fully understand the underlying mechanisms. In addition, the study indicated
References:
1. Crona J, Beuschlein F. Adrenocortical carcinoma - towards genomics guid- ed clinical care. Nat Rev Endocrinol, 2019;15(9):548-60
2. Datta J, Roses RE. Surgical management of adrenocortical carcinoma: An evidence-based approach. Surg Oncol Clin N Am, 2016;25(1):153-70
3. Else T, Kim AC, Sabolch A, et al. Adrenocortical carcinoma. Endocr Rev, 2014;35(2):282-326
4. Erickson LA, Rivera M, Zhang J. Adrenocortical carcinoma: Review and up- date. Adv Anat Pathol, 2014;21(3):151-59
5. Fassnacht M, Johanssen S, Quinkler M, et al. Limited prognostic value of the 2004 International Union Against Cancer staging classification for ad- renocortical carcinoma: Proposal for a Revised TNM Classification. Cancer, 2009;115(2):243-50
6. Margonis GA, Kim Y, Prescott JD, et al. Adrenocortical carcinoma: Impact of surgical margin status on long-term outcomes. Ann Surg Oncol, 2016;23(1):134-41
7. Beuschlein F, Weigel J, Saeger W, et al. Major prognostic role of Ki67 in lo- calized adrenocortical carcinoma after complete resection. J Clin Endocrinol Metab, 2015;100(3):841-49
8. Xiao H, Xu D, Chen P, et al. Identification of five genes as a potential bio- marker for predicting progress and prognosis in adrenocortical carcinoma. J Cancer, 2018;9(23):4484-95
9. Lippert J, Appenzeller S, Liang R, et al. Targeted molecular analysis in adre- nocortical carcinomas: A strategy toward improved personalized prognos- tication. J Clin Endocrinol Metab, 2018;103(12):4511-23
10. Sbiera S, Sbiera I, Ruggiero C, et al. Assessment of VAV2 expression refines prognostic prediction in adrenocortical carcinoma. J Clin Endocrinol Metab, 2017;102(9):3491-98
11. Subramanian C, Cohen MS. Over expression of DNA damage and cell cycle dependent proteins are associated with poor survival in patients with ad- renocortical carcinoma. Surgery, 2019;165(1):202-10
12. Pereira SS, Monteiro MP, Bourdeau I, et al. Mechanisms of endocrinol- ogy: Cell cycle regulation in adrenocortical carcinoma. Eur J Endocrinol, 2018;179(2):R95-110
a relationship between immunogenic status and prognosis in ACC. Further studies are needed to clarify whether or not im- mune cell infiltration is involved in ACC pathogenesis and de- velopment, and whether immune cell infiltration could be a potential therapeutic target in ACC. The prognostic genes iden- tified in this preliminary study require validation by further in- terventional and functional studies.
Conclusions
This study applied a comprehensive transcriptomic approach to analyze available transcriptomic data to identify key prognostic genes in ACC. We identified NDRG4 and CKS2 as key prognos- tic genes in ACC that also may have potential predictive val- ue for mitotane therapy. We further identified DNA mismatch repair pathways as dominant biologic pathways associated with CKS2. Moreover, immune cell infiltration in ACC differed according to NDRG4 and CKS2 expression. These findings re- quire validation by further interventional and functional studies.
Conflicts of Interest
None.
13. Gara SK, Lack J, Zhang L, et al. Metastatic adrenocortical carcinoma dis- plays higher mutation rate and tumor heterogeneity than primary tumors. Nat Commun, 2018;9(1):4172
14. Doghman M, Cazareth J, Lalli E. The T cell factor/beta-catenin antagonist PKF115-584 inhibits proliferation of adrenocortical carcinoma cells. J Clin Endocrinol Metab, 2008;93(8):3222-25
15. Bernini GP, Moretti A, Viacava P, et al. Apoptosis control and proliferation marker in human normal and neoplastic adrenocortical tissues. Br J Cancer, 2002;86(10):1561-65
16. Fay AP, Signoretti S, Callea M, et al. Programmed death ligand-1 expression in adrenocortical carcinoma: An exploratory biomarker study. J Immunother Cancer, 2015;3:3
17. Caramuta S, Lee L, Ozata DM, et al. Clinical and functional impact of TARBP2 over-expression in adrenocortical carcinoma. Endocr Relat Cancer, 2013;20(4):551-64
18. Zheng S, Cherniack AD, Dewal N, et al. Comprehensive pan-genomic char- acterization of adrenocortical carcinoma [published correction appears in Cancer Cell, 2016;30(2):363]. Cancer Cell, 2016;29(5):723-36
19. Assié G, Letouzé E, Fassnacht M, et al. Integrated genomic characterization of adrenocortical carcinoma. Nat Genet, 2014;46(6):607-12
20. Lubitz JA, Freeman L, Okun R. Mitotane use in inoperable adrenal cortical carcinoma. JAMA, 1973;223(10):1109-12
21. Terzolo M, Baudin AE, Ardito A, et al. Mitotane levels predict the outcome of patients with adrenocortical carcinoma treated adjuvantly following rad- ical resection. Eur J Endocrinol, 2013;169(3):263-70
22. Hermsen IG, Fassnacht M, Terzolo M, et al. Plasma concentrations of o,p’DDD, o,p’DDA, and o,p’DDE as predictors of tumor response to mitotane in adre- nocortical carcinoma: Results of a retrospective ENS@T multicenter study. J Clin Endocrinol Metab, 2011;96(6):1844-51
23. Sbiera S, Leich E, Liebisch G, et al. Mitotane inhibits Sterol-O-Acyl Transferase 1 triggering lipid-mediated endoplasmic reticulum stress and apoptosis in adrenocortical carcinoma cells. Endocrinology, 2015;156(11):3895-908
DATABASE ANALYSIS
24. Leek JT, Johnson WE, Parker HS, et al. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics, 2012;28(6):882-83
25. Langfelder P, Horvath S. WGCNA: An R package for weighted correlation network analysis. BMC Bioinformatics, 2008;9:559
26. Ritchie ME, Phipson B, Wu D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res, 2015;43(7):e47
27. Aran D, Hu Z, Butte AJ. xCell: Digitally portraying the tissue cellular hetero- geneity landscape. Genome Biol, 2017;18(1):220
28. Zhou RH, Kokame K, Tsukamoto Y, et al. Characterization of the human NDRG gene family: A newly identified member, NDRG4, is specifically ex- pressed in brain and heart. Genomics, 2001;73(1):86-97
29. Schilling SH, Hjelmeland AB, Radiloff DR, et al. NDRG4 is required for cell cycle progression and survival in glioblastoma cells. J Biol Chem, 2009;284(37):25160-69
30. Yang X, An L, Li X. NDRG3 and NDRG4, two novel tumor-related genes. Biomed Pharmacother, 2013;67(7):681-84
31. Ding W, Zhang J, Yoon JG, et al. NDRG4 is downregulated in glioblastoma and inhibits cell proliferation. OMICS, 2012;16(5):263-67
32. Xiao W, Zhao H, Dong W, et al. Quantitative detection of methylated NDRG4 gene as a candidate biomarker for diagnosis of colorectal cancer. Oncol Lett, 2015;9(3):1383-87
33. Chu D, Zhang Z, Zhou Y, et al. NDRG4, a novel candidate tumor suppressor, is a predictor of overall survival of colorectal cancer patients. Oncotarget, 2015;6(10):7584-96
34. Kotipatruni RP, Ren X, Thotala D, Jaboin JJ. NDRG4 is a novel oncogenic pro- tein and p53 associated regulator of apoptosis in malignant meningioma cells. Oncotarget, 2015;6(19):17594-604
35. Agosta C, Laugier J, Guyon L, et al. MiR-483-5p and miR-139-5p promote ag- gressiveness by targeting N-myc downstream-regulated gene family mem- bers in adrenocortical cancer. Int J Cancer, 2018;143(4):944-57
36. Spruck CH, de Miguel MP, Smith AP, et al. Requirement of Cks2 for the first metaphase/anaphase transition of mammalian meiosis. Science, 2003;300(5619):647-50
37. Kang MA, Kim JT, Kim JH, et al. Upregulation of the cycline kinase subunit CKS2 increases cell proliferation rate in gastric cancer. J Cancer Res Clin Oncol, 2009;135(6):761-69
38. Lan Y, Zhang Y, Wang J, et al. Aberrant expression of Cks1 and Cks2 con- tributes to prostate tumorigenesis by promoting proliferation and inhibit- ing programmed cell death. Int J Cancer, 2008;123(3):543-51
39. Chen R, Feng C, Xu Y. Cyclin-dependent kinase-associated protein Cks2 is as- sociated with bladder cancer progression. J Int Med Res, 2011;39(2):533-40
40. Liberal V, Martinsson-Ahlzén HS, Liberal J, et al. Cyclin-dependent kinase subunit (Cks) 1 or Cks2 overexpression overrides the DNA damage re- sponse barrier triggered by activated oncoproteins. Proc Natl Acad Sci USA, 2012;109(8):2754-59
41. Chehade M, Bullock M, Glover A, et al. Key MicroRNA’s and their targetome in adrenocortical cancer. Cancers (Basel), 2020;12(8):2198