ELSEVIER

Surgery

journal homepage: www.elsevier.com/locate/surg

SURGERY

MWEMBER 2018

$ -

===

Multi-genomic analysis of 260 adrenocortical cancer patient tumors identifies novel network BIRC5-hsa-miR-335-5p-PAX8-AS1 strongly associated with poor survival

Check for updates

Chitra Subramanian, PhD, MBAª, Reid McCallister, MSª, Mark S. Cohen, MD, FACSa,b,c,*

ª Department of Surgery, Michigan Medicine, Ann Arbor, MI

b Department of Pharmacology, University of Michigan, Ann Arbor, MI

” Department of Biomedical Engineering, University of Michigan, Ann Arbor, MI

ARTICLE INFO

Article history: Accepted 9 August 2022 Available online 3 October 2022

ABSTRACT

Background: Adrenocortical carcinoma is a rare endocrine cancer with poor overall survival. Linking survival outcomes to a common target across multiple genomic datasets incorporating microRNA-long non-coding RNA dysregulation have not been well described. We hypothesized that a multi-database analysis of microRNA-long noncoding RNA-messenger RNA regulatory networks associated with sur- vival will identify novel biomarkers.

Methods: Significantly dysregulated genes or microRNA in adrenocortical carcinoma compared to normal adrenal was identified from sequencing data for 260 human adrenocortical carcinomas using GEO2R. The miRnet identified hub microRNA and genes and long noncoding RNA and microRNA asso- ciated with survival genes. The R2 generated Kaplan-Meier curves. The database miRTarBase linked genes associated with poor survival and dysregulated microRNA.

Results: Analysis of genes and microRNAs differentially regulated in >50% of datasets revealed 75 genes and 12 microRNAs were upregulated, and 167 genes and 12 microRNAs were downregulated (bonf. P < .05). Kyoto Encyclopedia of Genes and Genomes pathway analysis revealed cell cycle, P53 signaling, arachidonic acid and innate immune response, and PI3/Akt are altered in adrenocortical carcinoma. A microRNA-target interaction network of differentially regulated microRNAs identified upregulated miRNA107, 103a-3p and 27a-3p, 16-5p, and downregulated 335-5p to have the highest degree of interaction with upregulated (ie, TPX2, CDK1, BIRC5, PRC1, CCNB1, GINS1) and downregulated (ie, RSPO3, NR2F1, TLR4, HOXA5, USP53, SLC16A9) hub genes as well as hub long noncoding RNAs XIST, NEAT1, KCNQ1OT1, and PAX8-AS1. Survival analysis revealed that the hub genes are associated with poor overall survival (P < . 05) of adrenocortical carcinoma in the Cancer Genome Atlas data.

Conclusion: A messenger RNA-microRNA-long noncoding RNA network analysis identified the BIRC5- miR335-5p-PAX8-AS1 network as one that was associated with poor overall survival in adrenocortical carcinoma, warranting further validation as a potential therapeutic target.

@ 2022 Published by Elsevier Inc.

Accepted for Presentation at the 42nd Annual Meeting of the American Asso- ciation of Endocrine Surgeons, April 25-27, 2022.

Chitra Subramanian and Mark S. Cohen contributed equally.

* Reprint requests: Mark S. Cohen MD, FSSO, FACS, Professor of Surgery, Phar- macology, and Biomedical Engineering, Vice Chair in Surgery for Clinical Opera- tions, Director, Medical School Pathway of Excellence in Innovation and Entrepreneurship, Director Center for Surgical Innovation, Department of Surgery, University of Michigan Hospital and Health Systems, 2920K Taubman Center, SPC 5331, 1500 East Medical Center Drive, Ann Arbor, MI 48109-5331.

E-mail address: cohenmar@med.umich.edu (M.S. Cohen); Twitter: @MarkCohenFACS

https://doi.org/10.1016/j.surg.2022.08.025

0039-6060/ 2022 Published by Elsevier Inc.

Introduction

Adrenocortical carcinoma (ACC) is a rare, clinically aggressive endocrine neoplasm of the adrenal cortex, with an annual inci- dence between 1 to 2 cases per million.1 Despite surgical inter- vention and existing adjuvant therapies, prognosis for patients diagnosed with ACC remains extremely poor, with an estimated 5- year overall survival of 30% to 47%.2,3 Complete surgical resection is the only potentially curative treatment available, though 70% of patients develop disease recurrence after initial neoplasm resec- tion.4 For patients with stage IV metastatic disease, treatment op- tions are severely limited, and median survival is <1 year.5

Figure 1. (A) Differential expression of genes in adrenocortical carcinomas compared to normal adrenal gland in 5 different datasets, namely GSE12368, GSE14922, GSE19750, GSE143383, and GSE90713 by Venn diagram. (B and C) Heatmaps of the top 25 upregulated and downregulated genes in adrenocortical carcinoma compared to normal adrenal gland.

A

GSE14922

GSE12368

401

54

2

206

3

3

97

21

1

8

24

516

33

GSE19750

11

24

3

16

0

2

0

4

35

14

4

37

23

38

GSE90713

18

76

19

108

GSE143383

B

Upregulated Genes

C

Downregulated Genes

KIAA0101 -

RSPO3

CDKI -

C2orf40-

CCNB1 -

ADH1B-

NCAPG-

PRC1 -

5

FMO2-

WFDC1-

TPX2-

HOXA5-

BIRC5-

NR2F1-

RACGAPI -

CYP11B2-

ZWINT-

DAPL1-

GINS1 -

SLITRK4-

DEPDCIB-

SLC16A9-

ASPM -

STEAP4-

TOPZA-

4

C7-

DTL-

GIPC2-

PBK-

CRHBP-

ANLN-

FAM65C-

BUBI -

PHYHIP-

CDKN3-

KLHL4-

MYBLI -

FNDC4-

CENPK-

PLCXD3-

CCNB2-

CAB39L-

CEP55-

3

ECHDC3-

MND1 -

IGF1-

ZNF367-

DLG2-

UBE2T-

DUSP26-

NDC80-

FNDC5-

CCNA2-

LRRC32-

CENPF-

ZNF185-

KIF11-

LMOD1-

TYMS-

2

SLC37A2-

MAD2L1 -

EMCN-

AURKA-

EPHX2-

ESMI-

USP53-

CDCA5-

TLR4

T

T

T

T

GSE143383

GSE9073

GSE12368

GSE14922

GSE19750

GSE143383

GSE9073

GSE12368

GSE14922

GSE19750

-2

-4

-6

Mitotane, either as a monotherapy or in combination with etopo- side, doxorubicin, and cisplatin (EDP-M, also known as the Italian Protocol), is the only approved therapy for patients at high risk of tumor recurrence and in cases not amenable to complete surgical resection. However, the efficacy of current chemotherapeutic

options, including EDP-M, remains limited, especially in advanced disease, and drug toxicities pose a significant barrier to therapy completion.7

Advances in genetic sequencing and immunotherapy have resulted in targeted therapies and ongoing clinical trials. However,

Figure 2. (A-D) David functional annotation to identify the pathways targeted by differentially expressed genes. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Ge- nomes pathway analysis of differentially regulated genes to identify GO biological process, GO cellular component (CC) and GO molecular function and associated pathways (P < .05). GO, gene ontology.

A

KEGG Pathway Analysis of Dysregulated Genes

B

GO Biosynthetic Pathway Analysis

Cell cycle

cell division

Oocyte meiosis

mitotic nuclear division

Progesterone-mediated oocyte maturation

sister chromatid cohesion

p53 signaling pathway

cell proliferation

Upregulated

MicroRNAs in cancer

G2/M transition of mitotic cell cycle

Downregulated

Folate biosynthesis

G1/S transition of mitotic cell cycle

protein phosphorylation

One carbon pool by folate

DNA replication

Focal adhesion

regulation of cell cycle

PI3K-Akt signaling pathway

mitotic cytokinesis

Complement and coagulation cascades

cell adhesion

Amoebiasis

innate immune response

inflammatory response

Phagosome

extracellular matrix organization

Mineral absorption

response to drug

Arachidonic acid metabolism

positive regulation of gene expression

Gastric acid secretion

aging

Pertussis

regulation of cell growth

response to hypoxia

Tyrosine metabolism

lipid metabolic process

0

2

4

6

8

10

12

0

5

10

15

20

25

Number of Genes

Number of Genes

C

D

GO Cellular Component Analysis

GO Molecular Function Analysis

nucleus

protein binding

nucleoplasm

ATP binding

cytoplasm

Upregulated

protein kinase binding

cytosol

Downregulated

chromatin binding

membrane

microtubule binding

midbody

protein serine/threonine kinase activity

centrosome

protein kinase activity

microtubule

ubiquitin-protein transferase activity

spindle

condensed chromosome kinetochore

drug binding

plasma membrane

cyclin-dependent protein serine/threonine kinase activity

extracellular exosome

protein binding

extracellular region

calcium ion binding

extracellular space

receptor binding

integral component of plasma membrane

actin binding

endoplasmic reticulum

catalytic activity

proteinaceous extracellular matrix

Upregulated

heparin binding

focal adhesion

cell surface

Downregulated

calmodulin binding

growth factor activity

extracellular matrix

integrin binding

0

10

20

30

40

50

60

extracellular matrix structural constituent

Number of Genes

0

10

20

30

40

50

Number of Genes

60

70

80

90

targeted therapeutics have yet to produce satisfactory results in ACC, demonstrating the need to identify novel targets for inter- vention.8 The noncoding RNAs (ncRNAs), namely microRNAs (miRNAs) and long ncRNAs (lncRNAs), do not encode proteins but regulate gene expression and protein function and represent a promising therapeutic target for treating ACC. MicroRNA are short ncRNAs that mediate gene expression by targeting sites in the untranslated region of mRNAs, leading to translation repression and mRNA degradation.9 Global repression of miRNA function is a hallmark of multiple cancers, and individual miRNA can function as either an oncogene or a tumor suppressor.10 Due to their role in dysregulation of cancer genetics, miRNA can be prognostic in- dicators or targets for novel drug therapies.11 Long ncRNAs can be precursors to miRNA or degrade miRNA through target RNA- directed miRNA degradation.9 In this way, lncRNA modulates miRNA, which in turn regulates mRNA translation and gene expression. Analysis of the lncRNA-miRNA-mRNA network pro- vides a more detailed picture of genetic dysregulation in cancer and identifies potential targets for novel therapeutics.

Independent studies have identified differentially expressed genes, miRNA, and lncRNA in ACC.12,13 However, the full picture of ACC pathogenesis and pathway dysregulation remains unclear. The integration of data from multiple studies targeting lncRNA, miRNA, and gene expression generates a more robust analysis of differen- tial expression in ACC and allows for the inspection of the larger regulatory network. Key genes associated with poor survival in ACC can be identified and cross-referenced with differentially expressed miRNA and lncRNA to construct a network of dysregulation

depicting potential drivers of ACC. Herein, we conducted an inte- grated analysis of genetic dysregulation data from multiple inde- pendent studies to identify pathogenic drivers and potential therapeutic targets in ACC.

Methods

Data acquisition and identification of differentially expressed genes and miRNAs

Genetic sequencing data for 260 ACC patients and 37 normal adrenal cortex samples (NAC) were acquired from 8 independent databases on the National Center for Biotechnology Information gene expression omnibus (GSE143385, GSE22816, GSE49279, GSE143383, GSE19750, GSE12368, GSE90713, GSE14922). Three databases contained miRNA sequencing data, and 5 contained gene expression data. Gene expression omnibus databases with ACC tumor samples and NAC were included in the analysis. The GEO2R analysis tool identified genes and miRNA with differential expres- sion between ACC and NAC. Genes were filtered using absolute value of log fold change in expression >1.5, P value < . 05 after Bonferroni correction, and dysregulation in >50% of datasets. The Bioinformatics and Evolutionary Genomics platform (VIB-UGent Center for Medical Biotechnology, Ghent, Belgium) generated a Venn diagram. Heat maps of the 25 most significantly up and downregulated genes were plotted in PRISM (GraphPad Software, Inc, San Diego, CA).

Figure 3. (A-B) miRNet identifying hub microRNAs (miRNAs) regulating differentially expressed genes. (C) miRNA-long ncRNA network. (D) Expression of hub miRNAs in adrenocortical carcinoma compared to normal adrenal.

A

B

C

CCNB1

CENPF

ZNF367

GINS1

NR2F1

CEP55

hsa-mir-16-5p

hsa-mir-16-5

hsa-mir-103a-3p

HOXA5

hsa-mir-107

CDK1

hsa-mir-147a

hsa-mir-27a

SLC

6A9

hsa-mir-107

hsa-mir-129-2-3p

RACGAP1

hsa-mir

KCNQ10 1

hsa mir 3

335-5p

hsa

-27

P

NE

hsa-mir-34a-5p

PAXB-

PRC1

hsa-mir-195-5p

hsa-mir-16

6:5p

USP53

hsa-mi

03a-3p

TLR4

hsa-mir-1-3p

hsa-mir-124-3p

CCNA2

nir- 103a

BIRC5

hsa-mir-335-5p

ASPM

RS

03

TPX2

TOP2A

DTL

D

hsa-miR-16-5p

hsa-miR-107

hsa-miR-103a-3p

1.5

2

2

1

1.5

1.5

0.5

1

0

1

-0.5

0.5

0.5

-1

0

0

-1.5

-0.5

-0.5

-2

ACC

Normal

-1

ACC

Normal

-1

ACC

Normal

hsa-miR-27a-3p

hsa-miR-335-5p

3

2

2

1

1

0

0

-1

-2

-1

-3

-2

-4

ACC

Normal

ACC

Normal

Functional enrichment analysis of the differentially expressed genes

The Database for Annotation Visualization and Integrated Dis- covery (DAVID) identified biochemical pathways enriched by dys- regulated genes common to >50% of datasets.14,15 The DAVID functional clustering annotation tool enabled enrichment analysis of genes in the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases. Pathways with significant enrichment and biological significance (P < . 001 and Fisher Exact test < . 001) were identified. Ten pathways containing the most

upregulated and downregulated genes were selected for further analysis.

Construction of miRNet network and the identification of hub miRNAs and lncRNAs

miRNet (https://www.mirnet.ca) used miRNA-target gene in- teractions from 4 well-annotated databases (ie, miRanda, miRTar Base, TarBase, and miRecords) to identify hub miRNAs involved in regulating the top differentially regulated genes. Deregulated

miRNAs from the miRNA databases with high degree interactions and betweenness centrality with significantly differentially regu- lated genes were identified. This network identified the hub miR- NAs critical to the regulation of top differentially regulated genes. These miRNAs were used as the input to determine which lncRNAs associated with hub miRNAs.

Survival analysis of hub genes and their validation

Hub gene expression in ACC compared to NAC were validated in The Cancer Genome Atlas (TCGA) dataset using University of Cali- fornia Santa Cruz Xena (https://xena.ucsc.edu). Genomic analysis and visualization platform R2 (https://hgserver1.amc.nl/cgi-bin/r2/ main.cgi) further validated differential expression and evaluated the prognostic effects of consistently dysregulated genes using the TCGA set of 79 ACC samples. Kaplan-Meier analysis determined the association of consistently dysregulated genes, including hub genes, with overall survival (P < . 05) in R2: the Genomics Analysis and Visualization Platform.16

Results

Comprehensive data analysis for the identification of the differentially expressed genes in ACC

Comprehensive analysis of the 5 ACC databases after standard- ization using GEO2R with a filtration threshold of absolute log-fold change >1.5 and risk-adjusted P value of ≤ .05 revealed that 773, 411, 807, 455, and 289 genes were differentially regulated in the databases GSE12368, GSE14922, GSE19750, GSE143383, and GSE90713, respectively. Comparison between the datasets revealed that 24 genes (10 upregulated and 14 downregulated) were present in all 5 datasets analyzed. Eighty-three (23 upregulated and 50 downregulated) and 147 (42 upregulated and 105 downregulated) genes were common between 4 and 3 datasets, respectively, out of the 5 as shown (Figure 1, A). Heat maps of the top upregulated and downregulated genes in the analyzed datasets are presented in Figure 1, B and C.

Functional annotation of differentially regulated genes in KEGG and GO enrichment pathways using DAVID

To evaluate the pathways modulated by the genes differentially regulated in >50% of the datasets, we performed KEGG pathway analysis and GO enrichment analysis. The top 10 biological path- ways that are enriched in ACC compared to normal adrenal gland (P < . 05) by upregulated genes in GO pathway mainly include cell cycle, mitotic cytokinesis, sister chromatid cohesion, G1/S, and G2/ M transition of mitotic cell cycle (Figure 2, A). Conversely, the top 10 biological pathways that are uniquely enriched in the down- regulated genes include the following: inflammatory response, innate immune response, response to hypoxia, and lipid metabolic processes (Figure 2, A). Cellular component analysis in GO revealed that upregulated genes in ACCs were enriched in the cytoplasm and nucleus, whereas the downregulated gene were enriched in the plasma membrane and endoplasmic reticulum (Figure 2, B). The GO molecular functions of differentially regulated genes were enriched in protein, adenosine triphosphate, microtubule, and kinase bind- ing, as well as catalytic activities (Figure 2, C). A KEGG pathway enrichment analysis of upregulated genes confirmed the role of cell cycle and P53 signaling, whereas the downregulated genes were enriched in PI3-Akt signaling, arachidonic acid metabolism, and tyrosine metabolism (Figure 2, D).

Figure 4. Box plot of differentially regulated hub genes in adrenocortical carcinoma compared to normal tissue (P < . 05). ACC, adrenocortical carcinoma.

20

ACC

Normal

Log, normalized count

15.

10

0

5-

0

TPX2

CDK1

BIRC5

PRCI

CCNB1

NS1

GIN

SPO3

NR2F

P53 HO

HOXA

LR4SLC1

16A9

Genes

Identification of hub miRNAs in ACC using miRnet

The miRNA databases GSE143385, GSE22816, and GSE49279 containing ACC tumor samples and normal adrenal tissue in the database of differentially expressed miRNAs in cancer were iden- tified using the query term adrenal cancer (dbDEMC [biosino.org]). These databases were analyzed to identify hub miRNAs involved in the differential regulation of key genes in ACC compared to normal adrenal cortex. The results from the analysis indicated that 146 (76 upregulated and 70 downregulated), 82 (30 upregulated and 52 downregulated), and 56 (23 upregulated and 36 downregulated) miRNAs, respectively, were differentially regulated in these data- sets with a significant risk-adjusted P value of ≤ .05. Analysis of miRNAs that were differentially regulated in 2 or more databases showed 12 upregulated, and 12 downregulated miRNAs were common to 2 of the 3 datasets. To identify the hub miRNAs that regulate the expression of the key differentially regulated genes, we constructed an miRnet (miRNet) as shown in Figure 3, A and B. The miRnet of upregulated genes that were expressed in >50% of the databases resulted in a network from which miRNAs and genes that have high degree and betweenness centrality were selected. This identified hsa-miR-103a-3p, hsa-miR-107, and hsa-miR-16-5p as the hub miRNAs that target the upregulated genes present in all 5 datasets (Figure 3, A). The miRnet of genes downregulated in > 50% of the databases had degree centrality of 60, 48, 32, 27, and 25 and betweenness centrality of 2,039.80, 1,490.02, 693.14, 295.13, 247.73 for hsa-miR-16-5p, hsa-miR-27a-3p, hsa-miR-335-5p, hsa-miR107, and hsa-miR-103a-3p, respectively (Figure 3, B). The expression levels of the hub miRNAs examined using dbDEMC showed upre- gulation of the hub miRNAs hsa-miR-16-5p, hsa-miR-27a-3p, hsa- miR107, and hsa-miR-103a-3p and downregulation of hsa-miR- 335-5p in ACC compared to normal adrenal gland (Figure 3, C). An miRNA-lncRNA network constructed in miRNet identified XIST, NEAT1, KCNQ1OT1, and PAX8-AS1 (Figure 3, D) as the hub miRNA interacting lncRNAs.

Identification, validation, and Kaplan-Meier survival analysis of miRNA regulated hub genes

To identify the hub genes that are regulated by miRNAs, we first considered differentially regulated genes that are expressed in ≥4 of the 5 datasets, then identified the genes that had high degree and betweenness centrality to the hub miRNAs. Finally, these hub genes were evaluated for significant associations with the overall survival of ACC. From this analysis, we identified a novel set of 12 hub genes regulated by miRNAs significantly modulated in ACCs. These include the following: TPX2, CDK1, BIRC5, PRC1, CCNB1, and GINS1 that are upregulated in all 5 datasets, and RSPO3, NR2F1, HOXA5,

Figure 5. (A and B) Overall Kaplan-Meyer survival analysis for upregulated and downregulated hub genes (P < . 05).

A

TPX2

CDK1

BIRC5

Overall survival probability

1.00

Overall survival probability

1.00

0.90

0.90

Overall survival probability

1.00

0.90

0.80

0.80

0.80

0.70

High

0.70

High

0.70

High

Low

0.60

Low

0.60

Low

0.60

0.50

0.50

0.50

0.40

0.40

0.40

0.30

0.30

0.30

0.20

0.20

0.20

0.10

raw p 7.4e-14

0.10

raw p 7.4e-14

0.10

0.00

bonf p 4.6e-12

0.00

bonf p 4.6e-12

raw p 5.1e-13

0.00

bonf p 3.2e-11

0

24

48

72

96

120

144

0 24 48 72 96 120 144

0

24

48

72

96

120

144

Follow up in months

Follow up in months

Follow up in months

PRC1

1.00

Overall survival probability

1.00

CCNB1

1.00

GINS1

Overall survival probability

0.90

0.90

Overall survival probability

High

0.90

0.80

Low

0.80

0.80

High

0.70

0.70

High

0.70

Low

0.60

0.60

Low

0.60

0.50

0.50

0.50

0.40

0.40

0.40

0.30

0.30

0.30

0.20

0.20

0.20

0.10

raw p 4.6e-10 bonf p 2.9e-08

0.10

raw p 1.1e-10

0.10

bonf p 7.2e-09

raw p 2.5e-07

0.00

0.00

bonf p 1.6e-05

0

24

48

72

96

120

144

0

24

48

72

96

120

144

0.00

0

24

48

72

96

120

144

Follow up in months

Follow up in months

Follow up in months

USP53, TLR4, and SLC16A9 that are downregulated in 4 or 5 data- sets. The differential expression levels of the miRNA-regulated hub genes were further validated in the TCGA dataset using the Uni- versity of California Santa Cruz Xena functional genome explorer. A box plot of expression of the 6 upregulated and 6 downregulated genes with significant differential expression in adrenal cancer compared to normal adrenal tissue (P < . 001) is shown in Figure 4. The overall survival of miRNA-regulated genes in TCGA dataset was further evaluated using the R2 database. Correlative analysis with Kaplan-Meier Survival data in ACC patients (Figure 5) using an automatic scan cutoff identified 6 upregulated and 6 down- regulated hub genes that significantly (P <. 01) correlated with poor overall survival. These genes were the differentially upregulated genes CDK1, CCNB1, PRC1, BIRC5, GINS1, and TPX2, as well as the differentially downregulated genes RSPO3, NR2F1, TLR4, USP53, HOXA5, and SLC16A9. This hub network is represented in Figure 6.

Discussion

Advancements in genetic sequencing have allowed for the characterization of molecular alterations that may occur in ACC, including TP53, APC, CTNNB1, IGF-2, CDK4, CDKN2A, RB1, MEN1, ZNRF3, DAXX, TERT, GNAS, RPL22, and PRKAR1A,13 leading to the proposal of therapeutic targets or prognostic indicators. Abnormal expression of such genes results in deregulation of WNT-6-catenin signaling, cell-cycle progression, telomere maintenance, chromatin remodeling, DNA repair, and steroid metabolism in ACC. Despite current understanding of ACC molecular mechanisms, it remains a

deadly disease, with a small population of patient samples available in individual microarray analysis. Integrated bioinformatic analysis of the microarray data can provide insights into the pathogenesis of ACC, leading to novel prognostic and diagnostic biomarkers.

In the present study, we used ACC databases to identify hub genes from a total of 242 differentially expressed genes using the following filter: present in >4 of the 5 total datasets, log fold change >1.5, bonf. P value ≤ .05, and associated with overall survival of ACC patients. Upregulated hub genes identified by this analysis, CDK1, CCNB1, PRC1, BIRC5, GINS1, and TPX2, are associated with cell-cycle pathways. Previous analyses of ACC genetics have identified most of the upregulated hub genes validating the findings of our analysis. However, upregulation of BIRC5 has not been previously reported in ACC. BIRC5, also known as survivin, is known to inhibit apoptosis and is associated with the immune pathway. Drugs targeting BIRC5, such as YM155, are already in clinical trials for use in cancer, and early trial reports do not indicate toxicity, making BIRC5 an exciting novel biomarker and therapeutic target for ACC patients.17 Our analysis also generated downregulated hub genes NR2F1 and USP53, which are novel identifications in ACC, in addition to HOXA5, TLR4, SLC16A9, and RSPO3, which were previously iden- tified as downregulated genes in ACC. Studies by Kanczkowski et al18 have shown downregulation of TLR4 in ACC tumor samples compared to normal adrenal tissue. They have also found that the introduction of TLR4 sensitized the ACC cells to apoptosis, vali- dating our findings. RSPO3 and TLR4 modulate WNT-6-catenin signaling and downregulation of NR2F1, a downstream component of TLR signaling, leading to activation of dormant cancer cells and

Figure 5. (continued).

B

RSPO3

NR2F1

USP53

1.00

1.00

Overall survival probability

0.90

Overall survival probability

0.90

Overall survival probability

1.00

0.90

0.80

0.80

0.80

0.70

High

0.70

High

0.70

High

0.60

Low

0.60

Low

0.60

Low

0.50

0.50

0.50

0.40

0.40

0.40

0.30

0.30-

0.30

0.20

0.20-

0.20

0.10

raw p 3.4e-04

0.10-

raw p 1.0e-03

0.10

raw p 5.0e-05

0.00

bonf p 0.021

0.00

bonf p 0.064

0.00

bonf p 3.2e-03

0

24

48

72

96

120

144

0

24

48

72

96

120

144

0

24 48 72 96 120 144

Follow up in months

Follow up in months

Follow up in months

HOXA5

Overall survival probability

1.00

Overall survival probability

1.00

TLR4

SLC16A9

Overall survival probability

1.00

0.90

0.90

0.90

0.80

0.80

High

0.80

0.70

High

0.70

Low

0.70

High

Low

0.60

0.60

0.60

Low

0.50

0.50

0.50

0.40

0.40

0.40

0.30

0.30

0.30

0.20

0.20

0.20

0.10

raw p 1.1e-03

0.10

raw p 6.5e-06

0.10

raw p 1.1e-06

0.00

bonf p 0.067

0.00

bonf p 4.1e-04

bonf p 6.7e-05

0

24

48

72

96

120

144

0 24 48 72 96 120 144

0.00

0

24 48 72 96 120 144

Follow up in months

Follow up in months

Follow up in months

metastasis in breast and head and neck squamous cell carcinoma.19 SLC16A9, a membrane transporter involved in adrenal steroid metabolism, has been identified as a biomarker for distinguishing ACC from adrenocortical adenoma, but its role in genetic dysre- gulation in ACC remains ill-defined and has yet to be investigated as a driver of ACC.20 HOXA5, found to be downregulated in our anal- ysis, may activate p53 signaling, thus inhibiting proliferation and promoting apoptosis of ACC cells.21 Hub gene USP53 is a ubiquitin- specific protease that has yet to be associated with ACC, but downregulation of USP53 has been reported in lung adenocarci- noma, esophageal carcinoma, and cervical cancer.22 In lung adenocarcinoma, increasing the expression levels of USP53 induced apoptosis and blocked glycolysis.23 Hence, better understanding of the role of USP53 in ACC will shed light on ACC pathogenesis.

Noncoding RNA, namely miRNA and lncRNA, are known to modulate genetic expression and can function as cancer bio- markers or novel therapeutic targets. Additionally, ncRNAs can easily be identified using liquid biopsy, allowing for noninvasive acquisition of diagnostic and prognostic data. In ACC, abnormal expression of miRNAs, including miR-503, miR-210, miR-195, miR- 335, miR-483-5p, and miR-483-3p,12 has been identified. However, these analyses identified abnormal miRNA expression independent of alterations in genetic deregulation. Therefore, integration of ACC sample data for both gene and miRNA expression enabled identi- fication and characterization of a complex network of key miRNAs associated with hub genes linked to overall survival in patients with ACC. Differentially regulated hub genes and miRNAs identified in patient tumors from 3 different datasets were used to construct an miRNet. Analysis of this network identified hsa-miR-16-5p, hsa-

miR-27a-3p, hsa-miR107, hsa-miR-103a-3p, and hsa-miR-335-5p as key miRNA interacting with hub genes in ACC samples. Down- regulation of hsa-miR-335-5p has been previously identified as a biomarker of ACC in comparison to adrenocortical adenoma, and its downregulation is associated with increased metastasis and inva- sion in multiple cancer types.24 However, our analysis demon- strated its association with deregulation of gene expression and overall survival of ACC, which were previously uncharacterized. The interaction between hub miRNA and the lncRNA XIST, NEAT1, KCNQ1OT1, and PAX8-AS1 completes the lncRNA-miRNA gene network of deregulation in ACC samples. Further, we identified PAX8-AS1, a potential therapeutic target in papillary thyroid car- cinoma and breast cancer,25 as a potential target for modulation of miR-335-5p and its downstream genes in ACC. Analysis of the overall network identified that our novel targets BIRC5, miRNA hsa- miR-335-5p, and lncRNA PAX8-AS1 are linked, and evaluation of this BIRC5-335-5p-PAX8-AS1 network axis may provide insights into targeted therapy for ACC.

Though this robust method of analysis has identified a novel network of hub genes and deregulated ncRNA, there were limita- tions to the study. First, this initial data analysis of gene expression and the interaction network effects of miRNA and lncRNA provided merely suggestive evidence of genetic network association and correlation with survival outcomes within the TCGA database but could not in itself establish a causative effect, which would require in vitro and in vivo validation of ACC cell lines and animal tumor models, which is beyond the scope of this manuscript. Second, 3 suitable databases containing miRNA expression data in both ACC and normal adrenal tissue were identified, limiting the level of

Figure 6. Hub messenger RNA-microRNA-long ncRNA network BIRC5-hsa-mir-335- 5p-PAX8-AS1 represented as star.

hsa

-5p

BIRC5

8-AS

hsa-mir- 16-5p

evidence for miRNA deregulation compared to the 5 databases used to analyze genetic regulation. Finally, deregulated miRNA and genes were filtered for significance and log fold change in expres- sion, which could result in additional genes/miRNA associated with ACC being excluded from the analysis.

In summary, the integrated analysis of 8 databases identified a novel set of 12 hub genes, 5 miRNA, and 4 lncRNA that form a network of dysregulation associated with ACC compared to normal adrenal tissue. Hub genes and miRNA were identified from 242 genes and 24 miRNA, with abnormal expression in >50% of data- bases based on interactions, significance, and association with overall survival in ACC. Of the hub genes identified, BIRC5 was identified as a potential therapeutic target for ACC. Four lncRNA were identified to be associated with all hub miRNA, completing the lncRNA-miRNA gene expression network. Hub miRNA are sig- nificant modulators of hub gene expression and represent novel therapeutic targets. Of the hub miRNAs, hsa-miR-335-5p was identified as a potential biomarker because it has previously been identified to be downregulated in ACC, but its role in ACC patho- genesis was not well studied. Long ncRNA associated with hub miRNA provides additional possibilities for novel ACC therapeutic targets, including PAX8-AS1, an exploratory therapeutic target in multiple other carcinomas. In total, mRNA-miRNA-lncRNA network analysis generated a more robust analysis of ACC tumor samples and identified the BIRC5-miR335-5p-PAX8-AS1 network as one that was associated with poor overall survival in ACC. Future in vitro validation is warranted to define its role as a potential therapeutic target.

Funding/Support

This research was funded in part by the National Institutes of Health, United States [R01 CA173292 and R01 CA216919 (M.S.C. and B.S.J.B.)] and the University of Michigan Department of Surgery, United States (M.S.C.).

Conflict of Interest/Disclosure

M.S.C. has equity in the following companies: Hylapharm LLC, WideawakeVR Inc, Pathware Inc, Ferroximend LLC, GigXR Inc, and Medguider LLC; he also has advisory roles in WideawakeVR Inc, Pathware Inc, and GigXR Inc. None of these disclosures have any relationship or conflict with the manuscript. The other authors have no conflicts of interests or disclosures to report.

References

1. Else T, Kim AC, Sabolch A, et al. Adrenocortical carcinoma. Endocr Rev. 2014;35: 282-326.

2. Allolio B, Fassnacht M. Clinical review: adrenocortical carcinoma: clinical up- date. J Clin Endocrinol Metab. 2006;91:2027-2037.

3. Pommier RF, Brennan MF. An eleven-year experience with adrenocortical carcinoma. Surgery. 1992;112:963-970;discussion 70-1.

4. Glenn JA, Else T, Hughes DT, et al. Longitudinal patterns of recurrence in pa- tients with adrenocortical carcinoma. Surgery. 2019;165:186-195.

5. Ayala-Ramirez M, Jasim S, Feng L, et al. Adrenocortical carcinoma: clinical outcomes and prognosis of 330 patients at a tertiary care center. Eur J Endo- crinol. 2013;169:891-899.

6. Fassnacht M, Dekkers OM, Else T, et al. European Society of Endocrinology Clinical Practice Guidelines on the management of adrenocortical carcinoma in adults, in collaboration with the European Network for the Study of Adrenal Tumors. Eur J Endocrinol. 2018;179:G1-G46.

7. Lagana M, Grisanti S, Cosentini D, et al. Efficacy of the EDP-M scheme plus adjunctive surgery in the management of patients with advanced adrenocor- tical carcinoma: the Brescia experience. Cancers (Basel). 2020;12:941.

8. De Filpo G, Mannelli M, Canu L. Adrenocortical carcinoma: current treatment options. Curr Opin Oncol. 2021;33:16-22.

9. Gebert LFR, MacRae IJ. Regulation of microRNA function in animals. Nat Rev Mol Cell Biol. 2019;20:21-37.

10. Lin S, Gregory RI. MicroRNA biogenesis pathways in cancer. Nat Rev Cancer. 2015;15:321-333.

11. Rupaimoole R, Slack FJ. MicroRNA therapeutics: towards a new era for the management of cancer and other diseases. Nat Rev Drug Discov. 2017;16: 203-222.

12. Hassan N, Zhao JT, Sidhu SB. The role of microRNAs in the pathophysiology of adrenal tumors. Mol Cell Endocrinol. 2017;456:36-43.

13. Crona J, Beuschlein F. Adrenocortical carcinoma - towards genomics guided clinical care. Nat Rev Endocrinol. 2019;15:548-560.

14. Huang da W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc. 2009;4: 44-57.

15. Huang da W, Sherman BT, Lempicki RA. Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists. Nucleic Acids Res. 2009;37:1-13.

16. Amsterdam University Medical Centers. R2: genomics analysis and visualiza- tion platform. http://r2.amc.nl. Accessed 5 June 2022.

17. Li F, Aljahdali I, Ling X. Cancer therapeutics using survivin BIRC5 as a target: what can we do after over two decades of study? J Exp Clin Cancer Res. 2019;38: 368.

18. Kanczkowski W, Tymoszuk P, Ehrhart-Bornstein M, Wirth MP, Zacharowski K, Bornstein SR. Abrogation of TLR4 and CD14 expression and signaling in human adrenocortical tumors. J Clin Endocrinol Metab. 2010;95:E421-E429.

19. Khalil BD, Sanchez R, Rahman T, et al. An NR2F1-specific agonist suppresses metastasis by inducing cancer cell dormancy. J Exp Med. 2022;219, e20210836.

20. Fernandez-Ranvier GG, Weng J, et al. Identification of biomarkers of adreno- cortical carcinoma using genomewide gene expression profiling. Arch Surg. 2008;143:841-846;discussion 846.

21. Chen D, Shen Z, Cheng X, et al. Homeobox A5 activates p53 pathway to inhibit proliferation and promote apoptosis of adrenocortical carcinoma cells by inducing Aldo-Keto reductase family 1 member B10 expression. Bioengineered. 2021;12:1964-1975.

22. Cheng W, Tang Y, Tong X, et al. USP53 activated by H3K27 acetylation regulates cell viability, apoptosis, and metabolism in esophageal carcinoma via the AMPK signaling pathway. Carcinogenesis. 2022;43:349-359.

23. Zhao X, Wu X, Wang H, Yu H, Wang J. USP53 promotes apoptosis and inhibits glycolysis in lung adenocarcinoma through FKBP51-AKT1 signaling. Mol Car- cinog. 2020;59:1000-1011.

24. Ye L, Wang F, Wu H, et al. Functions and targets of miR-335 in cancer. Onco Targets Ther. 2021;14:3335-3349.

25. Zhou P, Xu T, Hu H, Hua F. Overexpression of PAX8-AS1 inhibits malignant phenotypes of papillary thyroid carcinoma cells via miR-96-5p/PKN2 axis. Int J Endocrinol. 2021;2021, 5499963.

Discussion

Check for updates

Dr Martha Zeiger (Bethesda, MD): What are your true next steps? You say you are going to target this, but in what way? With CRISPR (Clustered Regularly Interspaced Short Palindromic Re- peats) technology? How are you going to target this pathway?

Dr Reid McCallister: There are some BIRC5-specific inhib- itors-one is YM155 and another is called Shepherdin; both would be of potential use for us in an in vitro study to see how they might modify the survival of adrenocortical carcinoma (ACC) in vitro. We also would like to create agents specific for miR-335-5p to see how modulation there affects ACC.

Dr Heather Wachtel (Philadelphia, PA): The clinical scenario we often see is an indeterminant adrenal mass. Do you see this as a potential diagnostic role for serum circulating micro RNAs (miR- NAs)? Do you have any plans to be studying this prospectively using patient serum? Did you have any consideration of using benign

adenoma tissue to see if we could identify targets on the spectrum of adenoma to carcinoma?

Dr Reid McCallister: That is an excellent question. miRNA 335-5p has been identified in ACC patients and in their serum, so it could represent a potential avenue for liquid biopsy and diagnostic tech- niques. At this time, that is not a future plan for our team, but I could see that happening as we move forward. I do think that the miRNA could be a differentiator between pathology and nonpathologic tissue.

Regarding the use of benign adenoma tissue, that was some- thing that we considered for this analysis. Ultimately, we decided to exclude adenoma from our normal group, but it is an analysis that could be done. We settled on this methodology because the data that we were able to access in the National Center for Biotech- nology Information database lent itself to comparing normal tissue to ACC, so we did not include adrenocortical adenoma at this time.