International Journal of Molecular Sciences

MDPI

Article Identification of Molecular Subtypes and Prognostic Characteristics of Adrenocortical Carcinoma Based on Unsupervised Clustering

Yuan Zhang, Cong Zhang, Kangjie Li, Jielian Deng, Hui Liu ®, Guichuan Lai OD, Biao Xie *DD and Xiaoni Zhong *

Department of Epidemiology and Health Statistics, School of Public Health, Chongqing Medical University, Yixue Road, Chongqing 400016, China; 2022110598@stu.cqmu.edu.cn (Y.Z.); 2022120781@stu.cqmu.edu.cn (C.Z.); 2021160020@stu.cqmu.edu.cn (K.L.); 2021120746@stu.cqmu.edu.cn (J.D.); 2020121415@stu.cqmu.edu.cn (H.L.); 2020111425@stu.cqmu.edu.cn (G.L.)

* Correspondence: kybiao@cqmu.edu.cn (B.X.); zhongxiaoni@cqmu.edu.cn (X.Z.)

check for updates

Citation: Zhang, Y .; Zhang, C .; Li, K .; Deng, J .; Liu, H .; Lai, G .; Xie, B .; Zhong, X. Identification of Molecular Subtypes and Prognostic Characteristics of Adrenocortical Carcinoma Based on Unsupervised Clustering. Int. J. Mol. Sci. 2023, 24, 15465. https://doi.org/10.3390/ ijms242015465

Academic Editor: Nam Deuk Kim

Received: 21 September 2023

Revised: 11 October 2023

Accepted: 18 October 2023

Published: 23 October 2023

CC

İ BY

Copyright: @ 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Abstract: Adrenocortical carcinoma (ACC) is a rare endocrine malignancy with a poor prognosis. Increasing evidence highlights the significant role of immune-related genes (IRGs) in ACC progression and immunotherapy, but the research is still limited. Based on the Cancer Genome Atlas (TCGA) database, immune-related molecular subtypes were identified by unsupervised consensus clustering. Univariate Cox analysis and Least Absolute Shrinkage and Selection Operator (LASSO) regression were employed to further establish immune-related gene signatures (IRGS). An evaluation of immune cell infiltration, biological function, tumor mutation burden (TMB), predicted immunotherapy response, and drug sensitivity in ACC patients was conducted to elucidate the applicative efficacy of IRGS in precision therapy. ACC patients were divided into two molecular subtypes through consistent clustering. Furthermore, the 3-gene signature (including PRKCA, LTBP1, and BIRC5) based on two molecular subtypes demonstrated consistent prognostic efficacy across the TCGA and GEO datasets and emerged as an independent prognostic factor. The low-risk group exhibited heightened immune cell infiltration, TMB, and immune checkpoint inhibitors (ICIs), associated with a favorable prognosis. Pathways associated with drug metabolism, hormone regulation, and metabolism were activated in the low-risk group. In conclusion, our findings suggest IRGS can be used as an independent prognostic biomarker, providing a foundation for shaping future ACC immunotherapy strategies.

Keywords: adrenocortical carcinoma; unsupervised clustering; immune-related genes; subtype; biomarkers; prognosis

1. Introduction

Adrenocortical carcinoma (ACC) is a rare and aggressive tumor that starts in the adrenal cortex. It is worth mentioning that approximately 40% to 60% of patients show symptoms and signs that suggest an excessive production of adrenal steroids [1,2]. Accord- ing to reports, the annual incidence of ACC is approximately 0.7-2 cases/million people, and the 5-year average survival rate is 20-25%. Nevertheless, the percentage of patients who survive for five years in stage IV is only 6-13% [3]. Currently, the sole treatment choice for individuals with non-metastatic ACC is complete removal of the tumor [4]. Continued exploration of ACC’s pathogenesis remains crucial. The existing comprehension of the disease’s mechanisms is rather unclear. The lack of clarity poses challenges for newly diagnosed ACC patients when assessing their risk levels in their initial appointments [5]. Hence, conducting thorough investigations to identify novel predictive biomarkers can not only provide deeper insights into the development of ACC metastasis but also facilitate the identification of fresh treatment targets [6,7].

Immune checkpoint inhibitor (ICIs) treatment has been an effective therapeutic tool that has dramatically changed the treatment pattern of solid tumors, including ACC,

head and neck squamous cell carcinoma (HNSCC), and colorectal cancer (CRC) [8-10]. Previous studies have found that patients’ age, tumor grade, and tumor range may have an impact on the prognosis of ACC [11-13]. However, the range of predictive markers has broadened in recent studies. The response of patients treated with ICIs may be influenced by immune-related genes (IRGs), such as programmed death-1 (PD-1) and PD-L1 in the tumor microenvironment (TME) [14,15], tumor-infiltrating immune cells (TICs) [16,17], tumor mutation burden (TMB) [18], and the expression of immune-related signaling pathways [19] in the TME, which may have an impact on the clinical outcome of patients treated with ICIs [20]. Zhang et al. pointed out that the presence of a V domain Ig T cell activation inhibitor (VISTA) in tumor cells is linked to the advancement of ACC and an elevated overall mortality rate. The findings provide new evidence supporting the potential of VISTA as a viable target for immunotherapy in the treatment of ACC [21]. Interfering with the binding between PD-1 and PD-L1 has the potential to enhance T cell proliferation, release cytokines, eliminate infected cells, lower the viral load, and potentially boost the immune response against tumors [22]. According to certain research, the capacity of PD-L1 expression in both tumor cells and immune cells infiltrating tumors to serve as indicators of the response to immunomodulatory agents has been well documented [23]. Furthermore, the activation of Common y-receptor-dependent cytokines and their JAK/STAT pathways is commonly observed in the majority of T-cell malignancies, potentially leading to a complete transformation in the approach to treating such malignancies [24]. In addition, Peng et al. discovered a pair of IRGs that exert a significant influence on tumor progression and the body’s immune response. These findings provide valuable information about the tumor microenvironment in ACC and indicate the possible usefulness of these genes as markers for immunotherapy in ACC [25]. The changing environment implies that the clinical results of patients receiving immunotherapy may be influenced not only by conventional demographic and histopathological factors but also by complex immune factors [26]. Hence, it is crucial to examine and assess the connection between the aforementioned factors and ACC immunotherapy, aiming to offer fresh perspectives for enhancing clinical diagnosis, prevention, and the advancement of novel treatment approaches [27].

In this work, using the unsupervised clustering method, we first classified ACC patients into distinct molecular subtypes according to IRGs. Additionally, we developed a three-gene signature using IRGs that were differentially expressed, enabling accurate evaluation of the prognostic risk in ACC patients. To investigate the potential significance of the three-gene signature and gain fresh perspectives on the disease’s pathogenesis and potential immunotherapy, we delved deeper into the correlation between gene signature and TICs, TMB, functional enrichment, drug sensitivity, and predictive immunotherapy.

2. Results

2.1. Identification of Characteristics and Biological Functions of Molecular Subtypes

A grand total of 2483 IRGs were acquired from the Immunology Database and Analysis Portal (ImmPort) database. Based on the expression patterns of 2483 IRGs in TCGA, we conducted consensus unsupervised clustering on 79 samples from TCGA-ACC in order to investigate novel molecular classifications for ACC patients. Assess the clustering values (k) for a range of 2 to 7 outcomes. The CDF Delta plot shows that the clustering outcomes exhibit greater stability when k = 2, and the PAC algorithm and a consistent heatmap determine that the optimal number of subtypes is 2 (Figure 1A-C; Table S1). In the comparison between patients of subtype 2 and subtype 1, subtype 2 patients exhibited significantly longer overall survival (OS), as indicated by the Kaplan-Meier (KM) survival curves (p = 0.0025, Figure 1D). To identify important IRGs in ACC, we utilized TCGA-ACC data to screen differentially expressed genes (DEGs) between subtypes 1 and 2. Additionally, 340 DEGs were identified based on a significance threshold of p-values < 0.05. Out of these, a total of 202 genes exhibited up-regulation while 138 genes showed down-regulation (Figure 1E and Table S2). Through the examination of the Gene Ontology (GO) function of these immune-related DEGs, it was discovered that, in terms of biological processes (BP), they were primarily concentrated

in cytokine-mediated signaling pathways and the regulation of leukocyte proliferation. Regarding cellular components (CC), the DEGs were predominantly engaged in the biological processes occurring on the external side of the plasma membrane and the MHC protein complex. Furthermore, molecular function (MF) exhibited activation in cytokine-mediated signaling pathways and receptor ligand activity (Figure 1F; Table S3). Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis showed that a total of 325 differential pathways were identified, of which the ten most significant pathways were Neuroactive ligand-receptor interaction, Hematopoietic cell lineage, Cytokine-cytokine receptor interaction, Cell adhesion molecules, Calcium signaling pathway, Graft-versus-host disease, Type I diabetes mellitus, Allograft rejection, Viral protein interaction with cytokine and cytokine receptor, PI3K-Akt signaling pathway (Figure 1G; Table S4). The results indicate that DEGs are linked to functions and pathways related to the immune system in ACC.

Figure 1. Characterization of immune properties and biological roles of molecular subtypes. (A-C) Identification of two molecular subtypes. (D) Kaplan-Meier survival analysis of OS between

A

B

C

consensus matrix k=2

consensus CDF

1.0

0.5

Delta area

relative change in area

n

Cluster 1

0.8-

0.4

Cluster 2

CDF

0.6-

0.3

0.4

0.2

0.2

0.0

0

0.0-

0.2-

0.4-

0.6

0.8

1.0-

2

3

4

5

6

7

consensus index

k

D

E

1.00

35

Survival probablity

Subtype 1

0.75

Subtype 2

30-

.0.50

-log10(pvalue)

25-

CDHR5

TNR

UGT2A3

p=0.0025

20

GLBIL3

0.25

· PRSS8

Down

0.00

15

CST2

None

Up

0

1000

2000

Days

3000

4000

5000

DACT2

10

Number at risk

SOX1

TECRL

Subtype 1

49

27

13

4

2

0

Subtype 2

30

19

9

4

0

0

·DDX53

0

0

1000

2000

3000

4000

5000

-10

-5

Days

0

5

10

log2FC

F

G

Top 20 of Pathway Enrichment

BP

Neuroactive ligand-receptor interaction-

Cytokine-cytokine receptor interaction

PI3K-Akt signaling pathway

Calcium signaling pathway

Cell adhesion molecules

Hematopoietic cell lineage

Count

CC

Phagosome

Viral protein interaction with cytokine and

30

cytokine receptor

Amoebiasis

60

Rheumatoid arthritis

90

Staphylococcus aureus infection-

Aldosterone synthesis and secretion-

p.adjust

Bile secretion

Antigen processing and presentation

0.0020

Inflammatory bowel disease

0.0015

MF

Graft-versus-host disease

0.0010

Type I diabetes mellitus

0.0005

Allograft rejection

logFC

Nicotine addiction

GO Terms

-3

3

Asthma

cytokine-mediated signaling pathway

regulation of leukocyte proliferation

0.025

0.050

0.075

0.100

external side of plasma membrane

MHC protein complex

receptor ligand activity

signaling receptor activator activity

GeneRatio

molecular subtypes. (E) The volcano plot of DEGs between two subtypes. (F) GO analysis showing the pathways enriched in DEGs. (G) KEGG analysis showing the pathways enriched in DEGs.

2.2. Construction and Verification of IRGs Gene Signature

At present, 340 DEGs have been identified in TCGA-ACC samples. To further decrease the feature dimension, we acquired 120 genes associated with OS (p < 0.05) that were identi- fied as prognostic genes in ACC patients by univariate Cox regression analysis. Afterwards, the 120 genes were incorporated into the Least Absolute Shrinkage and Selection Operator (LASSO) regression analysis, which was subsequently narrowed down to only three genes. In conclusion, a prognostic model (immune-related gene signature, IRGS) for ACC patients was established using three genes related to the immune system (PRKCA, LTBP1, and BIRC5). We then computed risk scores for each patient with ACC in both the training set and validation cohorts, using the subsequent mathematical equation: Riskscore = - 0.00487 x PRKCA + 0.18795 x LTBP1 + 0.30888 × BIRC5. Based on the median risk score, the patients were categorized into high-risk and low-risk groups. The prognostic effectiveness of IRGS was evaluated using time-dependent receiver operating characteristic (ROC) curves and KM curves. The survival rate in the training set showed a notable disparity between the high-risk and low-risk groups (p < 0.001). In addition, patients in the high-risk group demonstrated a significant decrease in their overall survival duration (Figure 2A). In the TCGA-ACC, the area under the curve (AUC) values for 1-year, 3-year, and 5-year OS were 0.851, 0.954, and 0.898, respectively, indicating that IRGS has excellent precision in predicting survival duration for both high-risk and low-risk groups (Figure 2E). The correlation between the risk score and the duration and status of survival is shown in Figure 2D.

Figure 2. The prognostic characteristics of IRGS. (A) The Kaplan-Meier analysis of TCGA-ACC patients. (B) The Kaplan-Meier analysis of the GSE33371 dataset. (C) The Kaplan-Meier analysis of

A

D

Survival probability

1.00

0.4

0.75

p < 0.0001

0.50

- High

Risk Score

0.3

Riskscore

+ Low

0.25

0.2

· High

· Low

0.00

0.1

0

1000

2000

3000

4000

5000

Riskscore

Number at risk

Days

0.0

High

39

18

6

0

0

0

20

40

60

80

Low

40

28

16

8

2

0

0

0

1000

2000

Days

3000

4000

5000

id

B

4000-

3000

Status

Survival probability

Survival Time

1.00 -

2000

Alive

0.75

p=0.0028

+ High

+ Low

1000

Dead

0.50

0.25

0

0

20

40

60

80

0.00

id

0

5

Years

10

15

Riskscore

Number at risk

High

11

1

1

E

Low

12

5

1

1

0

0

5

Years

10

15

1.00

C

1.00

True positive rate

0.75

Survival probablity

0.75

p=0.0048

High

Low

0.50

0.50

0.25

0.00

0.25

years = 0.851

0

2.5

5

7.5

Number at risk

Years

10

12.5

3 years = 0.954

Riskscore

High

0.00

5 years = 0.898

12

1

1

1

0

0

Low

12

6

2

1

1

0

5

7.5

10

0.25

0.50

0.75

1.00

0

2.5

12.5

0.00

Years

False positive rate

the GSE10927 dataset. (D) Risk score and survival time distribution in high- and low- risk groups. (E) The receiver operating characteristic curve in the TCGA-ACC.

To confirm the effectiveness of the IRGS, it was validated in the GSE33371 and GSE10927 datasets. Consistent with the results observed in the TCGA-ACC training cohort, the high-risk group demonstrated a notably reduced survival rate in the validation cohorts compared to the low-risk group (GSE33371: p = 0.0028, GSE10927: p = 0.0048; Figure 2B,C). Moreover, we confirmed the robust performance of the ROC curve during our validation process (Figure S1A,B). To thoroughly assess the model’s ability to differen- tiate, we generated a risk score and survival time distribution map for the GSE33371 and GSE10927 cohorts (Figure S1C,D). To summarize, these findings indicate the strength of IRGS in predicting outcomes for ACC patients.

2.3. IRGS and Prognostic Analysis of Clinical Characteristics

Next, we examined the association between the risk scores and different clinical attributes of the individuals. Figure 3A illustrates the utilization of a Sankey diagram to display the allocation of individuals across various risk categories. The Wilcoxon test results indicated significant variations in risk scores among various T stages, M stages, and clinical stages (Figure S2A-F). Nevertheless, there was no variation in the average risk scores across gender, age, and N stage. Next, we investigated the significance of IRGS in biological processes. The high-risk and low-risk groups underwent GO enrichment and KEGG pathway analysis. The results of GO enrichment included the cytokine-mediated signaling pathway, cell chemotaxis, receptor ligand activity, signaling receptor activator activity, the external side of the plasma membrane, and the secretory granule lumen (Figure 3B; Table S5). Furthermore, a grand total of 218 distinct KEGG pathways were discovered when comparing the high-risk and low-risk groups. Among these, the top ten pathways of utmost importance included cytokine-cytokine receptor interaction, viral protein interaction with cytokine and cytokine receptor, neuroactive ligand-receptor interaction, hematopoietic cell lineage, chemokine signaling pathway, IL-17 signaling pathway, rheumatoid arthritis, natural killer cell mediated cytotoxicity, NF-kappa B signaling pathway, and graft-versus- host disease (Figure 3C; Table S6).

2.4. Independent Prognostic Indicator of IRGS

In order to determine if the IRGS can function as an independent prognostic indicator, we incorporated various clinical factors, including age, gender, T stage, M stage, N stage, and clinical stage, into the evaluation of the IRGS in TCGA-ACC patients. Finally, prognosis was found to be associated with T stage, M stage, clinical stage, and risk score in the univariate Cox analysis (Figure 4A). Moreover, in the multivariate Cox analysis, the risk score was identified as a separate prognostic indicator for ACC patients (hazard ratio: 0.011-0.261, p < 0.001, Figure 4B). These findings were consistently validated in the GEO cohort (Figure S3A,B). The aforementioned findings suggest that the IRGS holds promise as a reliable and independent prognostic marker for individuals with ACC.

Figure 3. Functional analysis of IRGS. (A) The distribution of clinical features of TCGA-ACC samples. (B) The circular cluster diagram showing the results of GO pathway enrichment analysis. (C) The bar chart showing the top 20 most significant pathways in the KEGG pathway enrichment analysis.

A

C

80

Old

Cytokine-cytokine receptor interaction

Dead

Viral protein interaction with cytokine and cytokine receptor

60

Female

Neuroactive ligand-receptor interaction

count

Riskscore

40

Young

High

Hematopoietic cell lineage

Alive

Low

Chemokine signaling pathway

20

Male

IL-17 signaling pathway

Rheumatoid arthritis

0

Natural killer cell mediated cytotoxicity

-log10(pvalue)

Age

Gender

Status

B

Pathway

NF-kappa B signaling pathway

Graft-versus-host disease Allograft rejection

-

50

40

MAPK signaling pathway

30

Autoimmune thyroid disease

20

Type I diabetes mellitus

10

Amoebiasis

logFC

Viral myocarditis

-3

3

Intestinal immune network for IgA production

GO Terms

Malaria

· cytokine-mediated signaling pathway

JAK-STAT signaling pathway

external side of plasma membrane

receptor ligand activity

Phagosome

a cell chemotaxis

secretory granule lumen

0

20

40

60

· signaling receptor activator activity

Number of Gene

Figure 4. Univariate and multivariate Cox analysis results in TCGA-ACC patients. (A) The forest plot illustrating the results of the univariate Cox regression analysis in the TCGA-ACC patients. Green color blocks show the results of univariate Cox analysis. (B) The forest plot illustrating the findings of multivariate Cox regression analysis in the TCGA-ACC patients. Red color blocks show the results of multivariate Cox analysis. (p < 0.05 was considered significant).

A

B

pvalue

Hazard ratio

pvalue

Hazard ratio

Age

0.357

1.011(0.987-1.036)

Age

0.698

1.006(0.976-1.036)

Gender

0.999

1.001(0.469-2.137)

Gender

0.389

1.495(0.599-3.735)

T stage

<0.001

10.286(3.976-26.608)

T stage

0.077

9.358(0.784-111.718)

M stage

<0.001

6.150(2.710-13.959)

M stage

0.655

1.268(0.448-3.589)

N stage

0.152

2.038(0.769-5.400)

N stage

0.221

2.327(0.601-9.005)

H

Clinical stage

<0.001

6.476(2.706-15.498)

Clinical stage

0.342

0.296(0.024-3.652)

Riskscore

<0.001

(High vs. Low)

0.061(0.018-0.204)

Riskscore (High vs. Low)

<0.001

0.054(0.011-0.261)

0

5

10

15

20

25

0

20

40

60

80

100

Hazard ratio

Hazard ratio

2.5. Analysis of Tumor Immune Infiltration

We applied multiple algorithms to investigate variations in the level of immune infiltration in TCGA-ACC cohorts. According to our findings, the group with low risk showed decreased levels of tumor infiltration but displayed elevated immune and stromal scores, along with composite scores. In contrast, the high-risk group displayed a contrasting trend, characterized by increased levels of tumor infiltration and reduced immune and

stromal scores, along with composite scores. We utilized the CIBERSOR algorithm to determine the infiltration percentage of 22 immune cells. In the high-risk group, there was a notable increase in the presence of activated dendritic cells, and macrophages M0. Conversely, the low-risk group exhibited elevated levels of M2 macrophages, activated NK cells, and CD8+T cells (Figure 5A). The proportion of infiltration by 64 immune cells was calculated using the xCell algorithm. The low-risk group had significantly higher immunity scores, stromal scores, and microenvironmental scores compared to the high-risk group. The risk score (Figure 5B-D) showed a negative correlation with the immune score, stromal score, and microenvironmental score. Afterwards, the ESTIMATE algorithm was employed to assess the tumor’s purity. Typically, when there is a higher presence of immune cells and stromal cells, the tumor’s purity tends to decrease. During this investigation, a comparison between the high-risk and low-risk groups revealed that patients in the latter group displayed improved immune scores and stromal scores (Figure S3C-E). Furthermore, a significant inverse correlation was found between tumor purity and risk score (r = - 0.31, p = 0.0061; Figure 5E-F). To summarize, these results indicate that the tumor’s ability to suppress the immune system could be a contributing factor to the unfavorable prognosis of ACC patients in the high-risk group.

Figure 5. Exploring the landscape of immune cell infiltration in IRGS (A) CIBERSORT algorithm showing 22 immune cell infiltration levels. (B-D) xCell algorithm evaluating the infiltration levels of 64 immune cell types. (E,F) ESTIMATE algorithm showing the relationship between tumor purity and risk score. Data in (A-E) were analyzed by the Wilcoxon test; (ns, no significance, * p < 0.05, ** p < 0.01, *** p < 0.001, and **** p < 0.0001).

A

B

**

0.6

ns *

*

*

*


*

E

ns ns ns ns


ns

ns

ns

ns

ns

ns ns

ns

0.2

0.1

Riskscore

Proportion

ImmuneScore

0.4

₿ High

₿ Low

0.0

0.2

High

Low

High

C

Low


0.0

0.20

B cells memory

B cells naive

Dendritic cells activated

Dendritic cells resting

Eosinophils Macrophages MO

Macrophages M1

Macrophages M2

Mast cells activated

Mast cells resting

Monocytes

Neutrophils NK cells activated

NK cells resting

Plasma cells

T cells CD4 memory activated

T cells CD4 memory resting

T cells CD4 naive

T cells CD8

T cells follicular helper

T cells gamma delta T cells regulatory (Tregs)

StromalScore

0.15

0.10

Riskscore

0.05

₿ High

₿ Low

0.00

High

Low

D

E

F



MicroenvironmentScore

1.00

0.4

0.4

0.3

TumorPurity

0.75

0.3

r =- 0.31

0.2

0.50

TumorPurity

p =0.0061

0.2

0.1

0.25

0.1

0.0

High

0.00

Low

High

Low

0.0

Riskscore = High

Riskscore = High

0.00

0.05

0.10

0.15

0.20

0.25

@ Low

= Low

Riskscore

2.6. The Correlation between IRGS and TMB

To assess the correlation between IRGS and TMB, we utilized the R package “maftools”. The genes MUC16, TP53, CTNNB1, TTN, ASXL3, CNTNAP5, PCDH15, PKHD1, DST, and HMCN1 were found to have the highest mutation frequencies in patients with ACC (Figure 6A). Upon computation of the TMB for patients with ACC, it was observed that the TCGA-ACC cohort had a median TMB of 0.42/MB (Figure 6B), leading to the categorization of patients into high-TMB and low-TMB groups. To establish the inherent relationship between TMB and IRGS, we conducted a comparison of the TMB variation in different risk groups. The results of our study showed a markedly elevated TMB in the high-risk cohort, which demonstrated a direct association with the risk score of patients with ACC (Figure 6C,F). Moreover, the low-TMB group showed a tendency towards improved OS (Figure 6D). The KM curve indicated that TMB status had no impact on the prediction of IRGS, and consistently, the low-risk group demonstrated a superior survival benefit (Figure 6E).

Figure 6. The correlation between IRGS and TMB. (A) Waterfall plot illustrating the top 20 highly mu- tated genes of TCGA-ACC. (B) Median value of TMB distribution in TCGA-ACC patients. (C) Com- parison of TMB scores in different risk groups. (D) Kaplan-Meier curves of OS for high- and low- TMB groups. (E) KM curves of OS for ACC stratified by both risk score and TMB. (F) Correlation analysis between risk scores and TMB. Data in (C) were analyzed by the Wilcoxon test; ( **** p < 0.0001).

A

B

872

Mutation Burden

TMB

Altered in 51 (56.67%) of 90 samples

2.0

0

0

No. of samples

12

Riskscore

TMB/MB (log10)

1.5

MUC16

13%

1.0

TP53

CTNNB1

13%

High

13%

0.5

TTN

9%

Low

ASXL3

6%

0.0

CNTNAP5

7%

PCDH15

4%

-0.5

7%

Age

PKHD1

DST

4%

Old

-1.0-

HMCN1

7%

LRP1

4%

Young

-1.5-

NF1

6%

Median: 0.42/MB

SVEP1

6%

APOB

6%

CMYA5

4%

Gender

C

CSMD1

4%

FAT4

4%

female

FBN2

4%


FRAS1

4%

male

KMT2B

6%

Riskscore

Age

Status

15

Gender

Alive

Status

Dead

Missense_Mutation

Frame_Shift_Del

TMB

10

Splice_Site

Frame_Shift_Ins

5

Riskscore

In_Frame_Del

Multi_Hit

= High

Nonsense_Mutation

Low

0

High

Low

D

E

F

1.00 -

1.00-

Survival probablity

0.75

Survival probability

0.75

0.50

+ High

0.50

- H-H-TMB

15

+ Low

+ H-L-TMB

+ L-H-TMB

0.25

p<0.0001

0.25

+ L-L-TMB

total_perMB

r =0.59

p<0.0001

10

p =1.2×10-8

0.00

0.00

0

1000

2000

3000

4000

5000

0

1000

2000

3000

4000

Days

5000

Number at risk

5

Number at risk

Days

TMB

High- 39

21

7

1

0

0

H-H-TMB 27

11

0

0

0

H-L-TMB

11

2

L-H-TMB

12

7

0

0

O

Low

38

24

15

7

2

0

10

L-L-TMB

5

1

0

0

0

27

17

11

7

2

0

0

1000

2000

3000

4000

5000

0

1000

2000

3000

4000

5000

0.0

0.1

0.2

0.3

0.4

Days

Days

Riskscore

2.7. GSVA Analysis of KEGG and Immunologic Pathways between Different Risk Groups

By analyzing the enrichment score (ES) of 186 KEGG and 4872 immune gene sets, we discovered a combined count of 69 distinct KEGG pathways and 2385 immune path- ways among the two risk groups. The low-risk group showed significant enrichment in tyrosine metabolism, steroid hormone biosynthesis, drug metabolism by cytochrome P450,

metabolism of xenobiotics by cytochrome P450, aldosterone-regulated sodium reabsorption, primary bile acid bilsynthesis, the calcium signaling pathway and so on among the top 20 pathways of KEGG (Figure 7A; Table S7). Meanwhile, we computed the enrichment score (ES) for 2385 gene sets associated with the immune system using Gene Set Variation Analy- sis (GSVA). The gene sets belonging to the high-risk group were markedly up-regulated in the top 20 gene sets, whereas the gene sets of the low-risk group exhibited significant down-regulation (Figure 7B, Table S8).

A

B

Status

Status

Clinical stage

Riskscore

Clinical stage

Gender

Riskscore

Gender

Age

High

High

Riskscore

Low

2

Age

Riskscore

Low

2

KEGG_TYROSINE_METABOLISM

KEGG_STEROID_HORMONE_BIOSYNTHESIS

Age

1

Age

KEGG_DRUG_METABOLISM_CYTOCHROME_P450

Old

Old

1

KEGG_METABOLISM_OF_XENOBIOTICS_BY_CYTOCHROME_P450

Young

0

Young

KEGG_ALDOSTERONE_REGULATED_SODIUM_REABSORPTION

Gender

0

KEGG_PRIMARY_BILE_ACID_BIOSYNTHESIS

Gender

female

-1

KEGG_CALCIUM_SIGNALING_PATHWAY

female

KEGG_OOCYTE_MEIOSIS

male

-2

male

-1

KEGG_PROGESTERONE_MEDIATED_OOCYTE_MATURATION

KEGG_MISMATCH_REPAIR

Clinical stage

Clinical stage

KEGG_HOMOLOGOUS_RECOMBINATION

not reported

not reported

-2

KEGG_CELL_CYCLE

stage i

stage i

KEGG_DNA_REPLICATION

stage ii

stage ii

KEGG_BASE_EXCISION_REPAIR

stage iii

stage iv

stage iii

KEGG_NUCLEOTIDE_EXCISION_REPAIR

stage iv

KEGG_N_GLYCAN_BIOSYNTHESIS

Status

Status

KEGG_PROTEIN_EXPORT

Alive

Alive

KEGG_ONE_CARBON_POOL_BY_FOLATE

Dead

KEGG_PYRIMIDINE_METABOLISM

Dead

KEGG_PURINE_METABOLISM

Figure 7. KEGG and immunologic pathways in different risk groups and clinical characteris- tics. (A) The distribution of the top 20 most significant KEGG pathways in different risk groups. (B) Heatmap showing the immunologic pathways in different risk groups. (C) The box plot illus- trating the expression levels of ten immune checkpoints in different risk groups. Data in (C) were analyzed by the Wilcoxon test; (* p < 0.05, ** p < 0.01).

C

*

**

*

**

*

*

**

*

1.00

**

**

Expression level

0.75

Riskscore

0.50

High

Low

:

0.25

0.00

ji

CD27

CD274

CD40LG

CD8A

CTLA4

ICOS

LAIR1

PDCD1LG2

PRF1

TNFSF14

2.8. Correlation of IRGS with Immunotherapy Indicators

According to the Wilcox test analysis, several immune checkpoints, such as CD8A, exhibit high expression in the low-risk group in immunotherapy (Figure 7C). To evaluate the predictive capability of IRGS to predict the response to ACC immunotherapy, the Immunophenoscore (IPS) algorithm was employed. The findings indicated a significant

increase in the IPS score of the low-risk group, and a negative correlation was found between the risk score and the IPS score (Figure 8A,D). This implies that within the low-risk category, there could be a greater variety of immunophenotypes. Furthermore, patients classified as low-risk exhibited an elevated T cell-inflamed gene expression profile (GEP) score (Figure 8B,E) and MHC I association immunoscore (MAIS) score compared to those in the high-risk group (Figure 8C,F). The observation indicated that individuals in the low-risk category might exhibit a more positive reaction to immunotherapy.

Figure 8. The correlation between IRGS and immunotherapy indicators. (A-C) Violin plot of three ICI-related indicators between the high-risk group and the low-risk group in the TCGA-ACC dataset. (D-F) The relationship between three ICI-related indicators (IPS, GEP, and MAIS) and risk scores, respectively. (G-I) Comparison of three chemotherapeutic drugs between high- and low-risk groups in the TCGA-ACC dataset. Data in (A-C,G-I) were analyzed by the Wilcoxon test; ( ** p < 0.01, *** p < 0.001).

A

B

C

Riskscore

High

Riskscore

High

Riskscore

Low

High

Low

Low

**

**


10.0

0.8

1.0

E

IPS Score

GEP Score

7.5

0.4

MIAS Score

0.5

5.0

0.0

0.0

High

Low

High

Low

High

Low

D

E

F

10

1.00

IPS Score

9

0.75

GEP Score

MIAS Score

r =- 0.36

r =- 0.35

0.75

r =- 0.37

8

0.50

p = 0.00077

7

p =0.0012

0.50

0.25

p =0.0017

6

0.25

0.00

5

0.00

0.0

0.2

0.4

0.6

0.8

0.0

0.2

0.4

0.6

0.8

0.0 0.2 0.4 0.6 0.8

Riskscore

Riskscore

Riskscore

G

H

I

Cisplatin

Docetaxel

Gemcitabine

100

p=0.00024

0.03

p=2.2×10-7

p=0.0047

75

0.02

4

IC50

50

IC50

0.01

IC50

25

2

0

0.00

0

High

Low

High

Low

High

Low

Riskscore

High

Riskscore

High

Riskscore

Low

High

Low

Low

2.9. Drug Sensitivity Analysis

Given the fact that certain individuals display resistance towards traditional medi- cations, we investigated how patients with ACC respond to three frequently prescribed chemotherapy drugs. Notably, patients classified as high-risk may derive greater benefits

from these commonly prescribed medications, as evidenced by the higher IC50 found in patients from the low-risk group. This finding provides valuable insights for guiding clinical treatment strategies (Figure 8G-I).

3. Discussion

The prognosis for ACC varies greatly. Some tumors can be cured by complete surgery, while others cannot be removed, and they grow rapidly. The possibility of metastasis and diffusion is great, resulting in a poor prognosis [28]. Despite the utilization of various clinical biomarkers in the assessment of ACC prognosis [29-31], the adequacy of risk stratification for ACC is still lacking. Hence, it is crucial to discover additional potent biomarkers for precise treatment and prognosis estimation in ACC. This study employed information from the TCGA and GEO databases to investigate the molecular classification, prognostic factors, immune infiltration, gene signature construction, and drug sensitivity of ACC from the perspective of IRGs and provided strong support for treatment decision- making and prognostic assessment of ACC patients.

Unsupervised clustering analysis categorized the patients into two subtypes by 2483 IRGs. Our findings revealed that the OS in patients with subtype 2 exhibited a notably greater value compared to subtype 1, suggesting that these IRGs may potentially impact the prognosis of ACC. To enhance the investigation of the involvement of these IRGs in ACC risk classification, we employed two subtypes as the primary characteristics. Through univariate Cox and LASSO regression analysis, we identified three genes (PRKCA, LTBP1, BIRC5). Subsequently, we developed an ACC IRG signature to shed light on the biological understanding of IRGS in ACC prognosis. The IRGS demonstrated its effectiveness in risk stratification across different cohorts. Furthermore, it demonstrated the ability to autonomously forecast the outcome of ACC by merging with additional clinical characteristics. It is noteworthy that individuals classified into the high-risk group experience a less favorable outcome compared to the low-risk category. To summarize, the signature is anticipated to improve our comprehension of the origin and progression of ACC, ultimately aiding in the advancement of precise clinical management approaches.

Based on the research, it was revealed that there are three IRGS genes that are impli- cated in signaling pathways associated with both the immune system and the formation of tumors. Prior research has shown that the initiation of PRKCA prompts subsequent tumor- promoting consequences in various forms of cancer. This happens via the MAPK/ERK [32] and PI3K/AKT signaling pathways [33]. These impacts include the suppression of cell death and the stimulation of cell growth, movement, infiltration, and blood vessel forma- tion. Consequently, these processes collectively contribute to facilitating tumor progression. Regarding the prognosis, the presence of PRKCA in oral tongue squamous cell carcino- mas (OTSCC) was linked to a lack of smoking and alcohol consumption history, leading to reduced OS [34]. PRKCA overexpression has been associated with reduced survival outcomes not only in ACC patients but also in lung adenocarcinoma [35]. The primary function of LTBP1 is to have a significant impact on the extracellular matrix and attach to the extracellular trap substances of the transforming growth factor-ß (TGF-ß) group, including TGF-ß and BMPs (bone morphogenetic proteins), thereby controlling their acti- vation and release [36]. Elevated LTBP1 levels in NK/T cell lymphoma (NKTCL) result in heightened TGFß-1 secretion and release within the TME. Blocking the interaction between TGF-ß and LTBP1 can impede the activity of natural killer/T cells by deactivating the TGF- ß/Smad and p38MAPK signaling pathways, consequently impacting the advancement of lymphoma [37]. A previous investigation has presented proof that LTBP1 could function as a predictive risk element in instances of esophageal squamous cell carcinoma, and the elevated LTBP1 expression has been observed to display a favorable association with lym- phatic metastasis in individuals suffering from esophageal squamous cell carcinoma [38]. In addition, Fu et al. reported that glioblastoma (GBM) cells showing increased LTBP1 expression demonstrate heightened abilities in terms of proliferation and migration when compared to cells with lower LTBP1 levels. Furthermore, this occurrence is associated

with a more negative prognosis [39]. BIRC5 is widely acknowledged as a molecule that inhibits cell death and promotes cell proliferation and the advancement of tumors. Hence, it is widely recognized as a hopeful contender for therapeutic treatment [40]. Increased expression of BIRC5 in cancerous tissues has been demonstrated to enhance the forma- tion of new blood vessels, stimulate cellular proliferation, and impede programmed cell death [41]. Multiple research studies have indicated that a significant number of cancer types frequently exhibit elevated levels of BIRC5, which is linked to resistance against chemicals and an unfavorable prognosis in individuals with cancer, such as lung adenocar- cinoma [42], neuroblastoma [43], and glioma [44]. Furthermore, BIRC5 holds promise as a prospective candidate for identifying and predicting biomarkers for the identification, assessment, or prediction of breast cancer in patients [45].

In terms of signaling pathways, we noted that DEGs in the high-risk and low-risk groups were primarily enriched in immune regulation, inflammation, pathology, and bi- ological processes associated with cell signal transduction. The GO enrichment results comprised pathways involving cytokine-mediated signaling, cell chemotaxis, receptor ligand activity, and signaling receptor activator activity, among others. The KEGG analysis comprises pathways such as cytokine-cytokine receptor interaction, viral protein interaction with cytokine and cytokine receptor, the chemokine signaling pathway, and neuroactive ligand-receptor interaction, among others. Our findings are also supported by molecular typing using IRGs in head and neck cancer [46]. To summarize, the findings of these enrichment analyses emphasize the contrasting functional characteristics of high-risk and low-risk groups, specifically in relation to immune responses, cell signaling, and potentially viral interactions. The Janus kinase-signal transducer and activator of transcription (JAK- STAT) pathway is activated by cytokine-mediated signaling, leading to changes in gene expression within cells. Consequently, this pathway regulates the activity, differentiation, and function of immune cells [47]. Various diseases, including AML/MDS, have been thor- oughly investigated for immune system activation by extensively studying cytokines and their corresponding receptors. The potential to provide valuable insights into individual prognosis can be found in IL-6 and IL-10, as well as their corresponding receptors [48]. The presence of chemokines in tumors was linked to inflammation caused by T-cells [49]. In related studies, the neuroactive ligand-receptor interaction signaling pathway plays a direct and key role in the function of the nervous system [50].

Based on the attributes of IRGS, ACC patients were classified into different risk groups. The results showed that the group with low risk exhibited a higher levels of M2 macrophages, activated NK cells, and CD8+T cells. The results strongly indicate that there are different immune profiles between the high-risk and low-risk groups. Furthermore, our findings indicate a negative correlation between elevated immune score, stromal score, and microenvironment score and a low-risk score in patients. Similar findings were noted in patients with TME-associated subtypes of ACC [5]. In the same way, the risk score showed a strong correlation with elevated immune scores and the presence of tumor-infiltrating immune cells in melanoma. Additionally, the low-risk group exhibited the presence of CD8+T cells and NK cells [51]. Furthermore, based on the findings, we employed various techniques for immune infiltration to guarantee result stability. In the past few years, the rise of ICIs and targeted immunotherapy has become a crucial method of personalized treatment for patients with ACC. Recent studies indicate that the therapy of ICIs boosts the function of T cells by obstructing CTLA-4, PD-1, or PD-L1. There exists a direct relationship between the increased concentrations of ICIs and the efficacy of immunotherapy. This finding is of considerable importance in the context of prognosticating treatment responses [52]. Through an examination of the correlation between IRGS and ICIs manifestation, it was discovered that ten prevalent ICIs exhibited significant expression in the low-risk category, indicating heightened immune system activity within this group. In a highly inflammatory tumor setting, the increased expression of these ICIs could trigger a distinct feedback mechanism, fostering a dynamic immune environment that has the potential to enhance patient outcomes [53]. The high expression of ICIs may be linked to a better survival

outcome in the low-risk group, considering our findings. Furthermore, apart from utilizing these typical indicators, we have also chosen several potential immunotherapy predictors to evaluate the effectiveness of our signature in forecasting the reaction of ACC patients to ICIs. We discovered that individuals in the low-risk group exhibited elevated IPS, GEP, and MAIS scores, with all three scores displaying an inverse relationship with the risk scores. The results indicate that there are notable variations in the immune microenvironment between high-risk and low-risk groups, with the low-risk group having a higher probability of gaining advantages from immunotherapy [54].

The efficacy of ICIs in different types of tumors, such as non-small-cell lung cancer (NSCLC), metastatic gastrointestinal cancer, and colorectal cancer, has been assessed using TMB as a potential biomarker [55-58]. A higher TMB results in an increased number of novel antigens, enhancing the likelihood of T cell detection and improving the clinical response to ICIs [59]. In line with prior research on ACC, it was observed that individuals with high TMB experienced poorer survival results and exhibited a negative association with their risk score. Conversely, those in the low TMB category demonstrated favorable survival outcomes. Notably, these findings align with the outcomes of diffuse glioma and neuroblastoma, yet contrast with the findings of NSCLC, suggesting that TMB displays distinct prognostic traits across various cancer types [55,60]. Hence, a comprehensive examination of pan-cancer is necessary to completely understand its function. Furthermore, following risk stratification, the TMB status had no impact on the prediction of IRGS. The prognosis of the low-risk group has been better consistently, suggesting that IRGS possessed superior prognostic prediction capability, aligning with the findings of low- grade glioma [61]. We additionally investigated the association between IRGS and the effectiveness of antineoplastic medications for the created IRGS. Three frequently employed chemotherapeutic medications were chosen for the treatment of ACC. It was observed that the IC50 of cisplatin, gemcitabin, and docetaxel in the low-risk group was higher than that in the high-risk group, suggesting that high-risk patients may derive greater benefits from these commonly prescribed medications. Based on the current treatment guidelines, it is recommended to use etoposide, docetaxel, and cisplatin in combination with oral mitotane (EDP-M) as the primary treatment for metastatic ACC [62]. However, research has indicated that the combination of cisplatin and gemcitabin is also effective in treating ACC [63].

The current investigation extensively examined the possible influence of immune- associated genes on the prognosis of ACC individuals, unveiling their potential in terms of prognosis, immune attributes, and immunotherapy from various angles. Our main emphasis was on three IRGs, and we projected the OS of ACC patients by creating feature labels using these genes. Furthermore, the validation of the signature’s predictive capability is reaffirmed simultaneously. Of course, some limitations should be noted in this study. To begin with, ACC is an uncommon and diverse form of cancer, making it challenging to acquire extensive patient samples [64]. Nevertheless, the gene signature we have built can also be validated in two separate external validation sets. Secondly, this is a retrospective study, and it is inevitable that there may be selection bias and low statistical ability.

4. Materials and Methods

4.1. Data Acquisition and Pre-Processing

4.1.1. Data Download

In this work, we obtained RNA-seq data and the most recent clinical data of ACC patients from The Cancer Genome Atlas (TCGA) (https://xenabrowser.net/, accessed on 1 June 2023), and a grand total of 79 RNA-seq data samples were acquired. Key clinical features comprise age, gender, clinical stage, T stage, M stage, and N stage (Table S9). Furthermore, the validation datasets GSE33371 and GSE10927 were obtained from the GEO (https://www.ncbi.nlm.nih.gov/geo/, accessed on 2 June 2023) repository. According to the screen criteria, GSE33371 contains 23 samples with clinical features, and GSE10927 contains 24 samples with clinical features. Table S7 displays the detailed clinicopathological

characteristics of patients with ACC. Subsequently, we downloaded the data of 2483 IRGs from the ImmPort (https://immport.niaid.nih.gov, accessed on 23 July 2023) for analysis.

4.1.2. Data Pre-Processing

For this investigation, gene expression profile data were measured on the Illumina platform. We marked these genes utilizing the annotation document named “gencode.v22. annotation.gene.probeMap”. The measurement of gene expression profiles was done through the estimation of transcripts per million (TPM) and fpkm values converted to a log2 scale. Differential expression analysis utilized data in count format.

4.1.3. Data Screen Criteria

First, excluded were samples that did not have clinical follow-up information and those with a survival period shorter than 30 days in TCGA-ACC patients. Second, the gene expression data for all samples were complete. Finally, normal tissue samples were not considered, and only tumor samples were included. Figure 9 displays the workflow chart.

Figure 9. The workflow chart in this study.

TCGA-ACC RNA-seq (n=79)

Intersection with immune- related genes

Consensus clustering analysis

Identification of two molecular subtypes

Subtype 1

Subtype 2

1. KM survival analysis

2. ROC curve analysis

Identification immune-related DEGs (n=340)

3. GO enrichment and KEGG pathway analysis

Cox regression analysis and LASSO

Candidate genes

1. KM survival analysis 2. ROC curve analysis 3. Comparing the clinical features by univariate and multivariate Cox analysis

GSE33371

Construct immune-related gene signature (IRGS)

Validation

GSE10927

High-risk group

Low-risk group

GO, KEGG analysis

Immune analysis

TMB score

GSVA

Drug sensisitivity

1. Comparison of immune checkpoint expression 2. “CIBERSORT” algorithm 3. “ESTIMATE” algorithm

1. Enrichment scores of 186 KEGG pathways by GSVA-KEGG algorithm

4. xCell algorithm 5. IPS, GEP, MAIS score

2. GSVA analysis of 4872 immunologic pathways between high- and low-risk groups

4.2. Molecular Subtype Analysis

After intersecting the IRGs with the gene expression profile in TCGA-ACC, we used the “ConsensusClusterPlus” package [65] for unsupervised clustering and identified molecular subtypes. Then, we applied the “pam” algorithm and Pearson’s method for similarity measurement in hierarchical clustering. Sampling 80% of the data was done in each repetition of the clustering process, which was repeated 1000 times. The ideal value for the clustering number k was determined by analyzing the CDF, Delta area plot, and PAC algorithm [66], and only the stable clustering results were chosen.

4.3. Differential Expression Analysis

In this study, we employed the “DESeq2” package [67] to examine the variation in expression among various subtypes. The significant DEGs were identified by applying screening criteria of |log2FC| > 1 and p < 0.05. These DEGs were then compared with the IRGs. Afterward, the volcano map is employed for visualization.

4.4. Development and Verification of IRGs Signature

Initially, we conducted univariate Cox regression analysis on the potential IRGs us- ing the “survival” R package, taking into account genes with p-values below 0.05 as statistically noteworthy. Afterwards, we employed LASSO regression analysis to fur- ther identify prognostic-associated genes and obtained noteworthy outcomes from uni- variate Cox regression [68]. The lambda values were determined using a 10-fold cross- validation approach. Here, lambda = 0.20 with an appropriate partial likelihood residual was selected. The risk score of each sample was obtained by the following formula: Riskscore = ExpressionIRGS1 x CoefficientIRGS1 + ExpressionIRGS2 x CoefficientIRGS2 + … + ExpressionIRGSn × CoefficientIRGSn [69]. The risk assessment of TCGA-ACC cancer pa- tients resulted in the categorization of individuals into high- and low-risk groups based on the median risk score. Ultimately, the identical approach is employed to authenticate the outcomes in the validation set. To ascertain the independence of IRGS as a predictive marker, univariate Cox analysis was conducted to assess the risk score in conjunction with relevant clinical variables.

4.5. Analysis of Immune Microenvironment Characteristics

In our study, the “CIBERSORT” algorithm was used to calculate 22 representative immune cells based on the LM22 profile and normalized expression matrix [70]. Simultane- ously, we employed the “ESTIMATE” algorithm to quantify the immune microenvironment of ACC patients. To ensure the robustness and consistency of our findings, we also in- tegrated the xCell algorithm into our analysis, enabling the evaluation of the infiltration ratios of 64 specific immune cell types and facilitating the comparative assessment of their levels of infiltration among distinct risk groups.

4.6. Analysis of Somatic Mutations

To determine the mutational burden of ACC, we employed the R package “TCGAbi- olinks” [71] to fetch ACC mutation data. Then, we used the “maftools” package [72] to analyze the mutation data and obtain TMB for each patient. Visualizing the results with a waterfall diagram. Based on the median, TMB was categorized into two groups: high-TMB and low-TMB, and the mutation characteristics were subsequently compared between the different risk groups. Additionally, the KM curve was utilized to examine the relationship between TMB, different risk levels of TMB, and the prognosis of ACC.

4.7. Biological Functional Analysis

To gain insight into the possible molecular mechanisms underlying the DEGs in IRG subtypes and risk groups, the “clusterprofiler” package was used to conduct “GO” enrichment and “KEGG” pathway analyses [73]. The differential functions and pathways were selected via the “BH” method with an adjusted p value < 0.05. A bubble diagram

was used to depict the major functions and pathways. Furthermore, we used GSVA to directly access the downloaded C2 reference gene sets, which consist of KEGG data and C7 immunologic signature gene sets. The pathway difference was analyzed using the “limma” package, and the analysis results were visualized through a heatmap.

4.8. Predictive Analysis of Immunotherapy

In order to enhance our comprehension of the correlation between IRGS and specific prognostic immunotherapy indicators, we performed a comparative examination between different risk groups. The immunogenicity of tumors can be evaluated, and the effectiveness of different types of ICIs therapies can be predicted using IPS. The evaluation of four different immunophenotypes, namely antigen presentation, effector cells, inhibitory cells, and checkpoints, is encompassed by IPS. The immunogenicity of the IPS was assessed using a z-score ranging from 0 to 10, where a higher value suggests an increased immunogenic response. Simultaneously, we employed GEP [74] and MIAS (Microenvironment, Immune Response, and Anti-tumor Sensitivity) scores as indicators for forecasting patients’ reactions to immunotherapy.

4.9. Drug Sensitivity Analysis

The presence of diversity in ACC leads to significant differences in how patients respond to the same pharmacological treatments. We utilized the Genomics of Drug Sensitivity in Cancer (GDSC) database as our main training dataset to examine the variation in the responsiveness of anti-cancer medications among patients belonging to diverse risk categories. Afterwards, we utilized the R “oncoPredict” package to determine the IC50 values of three traditional chemotherapy drugs used in the clinical treatment of ACC patients [75]. It should be emphasized that cisplatin is an ACC treatment drug that has been extensively researched. In addition, in vitro studies have also confirmed the effectiveness of taxanes [76]. Currently, several investigations have conducted research on the utilization of cisplatin in conjunction with gemcitabine and docetaxel for treating individuals diagnosed with advanced adrenocortical carcinoma.

5. Conclusions

We identified two new subtypes of ACC based on IRGs, constructed a prognostic gene signature for IRG-related risks, including PRKCA, LTBP1, and BIRC5, and verified its predictive ability. The IRGS could potentially be a dependable predictive biomarker for individualized care of ACC patients, offering a theoretical foundation for accurate immunotherapy in ACC patients.

Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ijms242015465/s1.

Author Contributions: Conceptualization, Y.Z .; methodology, Y.Z. and C.Z .; software, Y.Z. and K.L .; validation, Y.Z. and C.Z .; formal analysis, Y.Z. and J.D .; resources, Y.Z., C.Z. and H.L .; data curation, G.L .; writing-original draft preparation, Y.Z .; writing, review and editing, X.Z. and B.X .; visualization, Y.Z .; supervision, B.X .; project administration, X.Z .; funding acquisition, B.X. All authors have read and agreed to the published version of the manuscript.

Funding: This research was funded by the National Youth Science Foundation Project, grant number (82204159), and Science and Technology Research Program of Chongqing Municipal Education Commission (Grant No. KJQN202300423).

Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.

Data Availability Statement: The corresponding data and results were generated by the TCGA (https://xenabrowser.net/, accessed on 1 June 2023), GEO (http://www.ncbi.nlm.nih.gov/geo/, accessed on 2 June 2023), and ImmPort databases (https://immport.niaid.nih.gov, accessed on 23 July 2023).

Acknowledgments: Here, we thank the TCGA, GEO, and ImmPort databases for providing relevant data to support our studies.

Conflicts of Interest: The authors declare no conflict of interest.

Abbreviations

ACCAdrenocortical carcinoma
ICIsImmune checkpoint inhibitors
HNSCCHead and neck squamous cell carcinoma
CRCColorectal cancer
IRGsImmune-related genes
PD-1Programmed death-1
TICsTumor infiltrating immune cells
TMBTumor mutation burden
TMETumor microenvironment
VISTAV domain Ig T cell activation inhibitor
NSCLCNon-small cell lung cancer
OSOverall survival
DEGsDifferentially expressed genes
KMKaplan-Meier
GOGene ontology
KEGGKyoto Encyclopedia of Genes and Genomes
BPBiological process
MFMolecular function
CCCellular component
LASSOLeast Absolute Shrinkage and Selection Operator
ROCReceiver operating characteristic
AUCArea under the curve
GEOGene Expression Omnibus
ESEnrichment score
GSVAGene Set Variation Analysis
OTSCCOral tongue squamous cell carcinomas
NKTCLNK/T cell lymphoma
GBMGlioblastoma
JAK-STATJanus kinase-signal transducer and activator of transcription
NSCLCnon-small-cell lung cancer
TCGAThe Cancer Genome Atlas
EDP-MEtoposide, Docetaxel, and cisplatin plus oral mitotane
ImmPortImmunology database and analysis portal
TPMPer million transcript
CDFCumulative Distribution Function
PACProportion of ambiguous clustering
LASSOLeast abolute shrinkage and selection operator
IRGSImmune-related gene signature,
BHBenjamini-Hochberg
IPSImmunophenoscore
GEPGene expression profile
MIASMHC I association immunoscore
GDSCGenomics of Drug Sensitivity in Cancer

References

1. Landwehr, L .- S .; Altieri, B .; Schreiner, J .; Sbiera, I .; Weigand, I .; Kroiss, M .; Fassnacht, M .; Sbiera, S. Interplay between glucocorti- coids and tumor-infiltrating lymphocytes on the prognosis of adrenocortical carcinoma. J. Immunother. Cancer 2020, 8, e000469. [CrossRef] [PubMed]

2. Allolio, B .; Fassnacht, M. Clinical review: Adrenocortical carcinoma: Clinical update. J. Clin. Endocrinol. Metab. 2006, 91, 2027-2037. [CrossRef] [PubMed]

3. Shetty, I .; Venzon, D.J .; Mauda-Havakuk, M .; Thomas, B .; Bernstein, D .; Flowers, C .; Anderson, V .; Raygada, M .; Levy, E.B .; Wood, B.J .; et al. Effect of local therapies on survival in patients with metastatic adrenocortical carcinoma. J. Clin. Oncol. 2020, 38, e17114. [CrossRef]

4. Kerkhofs, T.M .; Verhoeven, R.H .; Van der Zwan, J.M .; Dieleman, J .; Kerstens, M.N .; Links, T.P .; Van de Poll-Franse, L.V .; Haak, H.R. Adrenocortical carcinoma: A population-based study on incidence and survival in the Netherlands since 1993. Eur. J. Cancer 2013, 49, 2579-2586. [CrossRef]

5. Lai, G .; Liu, H .; Deng, J .; Li, K .; Zhang, C .; Zhong, X .; Xie, B. The characteristics of tumor microenvironment predict survival and response to immunotherapy in adrenocortical carcinomas. Cells 2023, 12, 755. [CrossRef]

6. Crona, J .; Beuschlein, F. Adrenocortical carcinoma-Towards genomics guided clinical care. Nat. Rev. Endocrinol. 2019, 15, 548-560. [CrossRef]

7. Zheng, S .; Cherniack, A.D .; Dewal, N .; Moffitt, R.A .; Danilova, L .; Murray, B.A .; Lerario, A.M .; Else, T .; Knijnenburg, T.A .; Ciriello, G .; et al. Comprehensive Pan-Genomic Characterization of Adrenocortical Carcinoma. Cancer Cell 2016, 29, 723-736. [CrossRef]

8. Brabo, E.P .; Moraes, A.B .; Neto, L.V. The role of immune checkpoint inhibitor therapy in advanced adrenocortical carcinoma revisited: Review of literature. J. Endocrinol. Investig. 2020, 43, 1531-1542. [CrossRef]

9. Burtness, B .; Harrington, K.J .; Greil, R .; Soulières, D .; Tahara, M .; de Castro, G., Jr .; Psyrri, A .; Basté, N .; Neupane, P .; Bratland, A .; et al. Pembrolizumab alone or with chemotherapy versus cetuximab with chemotherapy for recurrent or metastatic squamous cell carcinoma of the head and neck (KEYNOTE-048): A randomised, open-label, phase 3 study. Lancet 2019, 394, 1915-1928. [CrossRef]

10. Oliva, M .; Spreafico, A .; Taberna, M .; Alemany, L .; Coburn, B .; Mesia, R .; Siu, L. Immune biomarkers of response to immune- checkpoint inhibitors in head and neck squamous cell carcinoma. Ann. Oncol. 2019, 30, 57-67. [CrossRef]

11. Assié, G .; Antoni, G .; Tissier, F .; Caillou, B .; Abiven, G .; Gicquel, C .; Leboulleux, S .; Travagli, J .- P .; Dromain, C .; Bertagna, X .; et al. Prognostic Parameters of Metastatic Adrenocortical Carcinoma. J. Clin. Endocrinol. Metab. 2007, 92, 148-154. [CrossRef] [PubMed]

12. Miller, B.S .; Gauger, P.G .; Hammer, G.D .; Giordano, T.J .; Doherty, G.M. Proposal for modification of the ENSAT staging system for adrenocortical carcinoma using tumor grade. Langenbeck’s Arch. Surg. 2010, 395, 955-961. [CrossRef]

13. Gonzalez, R.J .; Tamm, E.P .; Ng, C .; Phan, A.T .; Vassilopoulou-Sellin, R .; Perrier, N.D .; Evans, D.B .; Lee, J.E. Response to mitotane predicts outcome in patients with recurrent adrenal cortical carcinoma. Surgery 2007, 142, 867-875. [CrossRef]

14. Thompson, R.H .; Dong, H .; Kwon, E.D. Implications of B7-H1 Expression in Clear Cell Carcinoma of the Kidney for Prognostica- tion and Therapy. Clin. Cancer Res. 2007, 13, 709s-715s. [CrossRef]

15. Grimm, M .; Gasser, M .; Koenigshausen, M .; Stein, C .; Lutz, J .; Thiede, S.K .; Heemann, U .; Waaga-Gasser, A. Clinical significance and therapeutic potential of programmed death-1 ligand-1 and programmed death-1 ligand-2 expression in human colorectal cancer. J. Clin. Oncol. 2008, 26, 15005. [CrossRef]

16. Tian, X .; Xu, W .; Wang, Y .; Anwaier, A .; Wang, H .; Wan, F .; Zhu, Y .; Cao, D .; Shi, G .; Zhu, Y .; et al. Identification of tumor-infiltrating immune cells and prognostic validation of tumor-infiltrating mast cells in adrenocortical carcinoma: Results from bioinformatics and real-world data. Oncoimmunology 2020, 9, 1784529. [CrossRef]

17. Pereira, S.S .; Costa, M.M .; Guerreiro, S.G .; Monteiro, M.P .; Pignatelli, D. Angiogenesis and Lymphangiogenesis in the Adrenocor- tical Tumors. Pathol. Oncol. Res. 2017, 24, 689-693. [CrossRef]

18. Lin, G .; Hari, P .; Chhabra, S .; Hamadani, M .; D’Souza, A .; Szabo, A .; Dhakal, B. The significance of beta-II microglobulin (B2M) and International Staging System (ISS) in multiple myeloma (MM) patients (pts.) with renal impairment (RI). J. Clin. Oncol. 2020, 38 (Suppl. S15), 8544. [CrossRef]

19. Zhou, Y .; Xu, J .; Luo, H .; Meng, X .; Chen, M .; Zhu, D. Wnt signaling pathway in cancer immunotherapy. Cancer Lett. 2022, 525, 84-96. [CrossRef]

20. Bai, R .; Lv, Z .; Xu, D .; Cui, J. Predictive biomarkers for cancer immunotherapy with immune checkpoint inhibitors. Biomark. Res. 2020, 8, 34. [CrossRef]

21. Zhang, Z .; Li, M .; Wang, J .; Liu, M.S .; Chen, H .; Lou, Y .; Wang, Y .; Sun, Q .; Zhu, D.L .; Li, P .; et al. Expression and clinical significance of VISTA and PD-L1 in adrenocortical carcinoma. Endocr. Relat. Cancer 2022, 29, 403-413. [CrossRef] [PubMed]

22. Vilain, R.E .; Menzies, A.M .; Wilmott, J.S .; Kakavand, H .; Madore, J .; Guminski, A .; Liniker, E .; Kong, B.Y .; Cooper, A.J .; Howle, J.R .; et al. Dynamic Changes in PD-L1 Expression and Immune Infiltrates Early During Treatment Predict Response to PD-1 Blockade in Melanoma. Clin. Cancer Res. 2017, 23, 5024-5033. [CrossRef] [PubMed]

23. Taube, J.M .; Klein, A .; Brahmer, J.R .; Xu, H .; Pan, X .; Kim, J.H .; Chen, L .; Pardoll, D.M .; Topalian, S.L .; Anders, R.A. Association of PD-1, PD-1 ligands, and other features of the tumor immune microenvironment with response to anti-PD-1 therapy. Clin. Cancer Res. 2014, 20, 5064-5074. [CrossRef] [PubMed]

24. Waldmann, T.A .; Chen, J. Disorders of the JAK/STAT Pathway in T Cell Lymphoma Pathogenesis: Implications for Immunother- apy. Annu. Rev. Immunol. 2017, 35, 533-550. [CrossRef] [PubMed]

25. Peng, Y .; Song, Y .; Ding, J .; Li, N .; Zhang, Z .; Wang, H. Identification of immune-related biomarkers in adrenocortical carcinoma: Immune-related biomarkers for ACC. Int. Immunopharmacol. 2020, 88, 106930. [CrossRef]

26. Assié, G .; Jouinot, A .; Fassnacht, M .; Libé, R .; Garinet, S .; Jacob, L .; Hamzaoui, N .; Neou, M .; Sakat, J .; de La Villéon, B .; et al. Value of Molecular Classification for Prognostic Assessment of Adrenocortical Carcinoma. JAMA Oncol. 2019, 5, 1440-1447. [CrossRef]

27. Liu, W .; Xie, X .; Qi, Y .; Wu, J. Exploration of Immune-Related Gene Expression in Osteosarcoma and Association with Outcomes. JAMA Netw. Open 2021, 4, e2119132. [CrossRef]

28. Lippert, J .; Dischinger, U .; Appenzeller, S .; Prete, A .; Kircher, S .; Skordilis, K .; Elhassan, Y.S .; Altieri, B .; Fassnacht, M .; Ronchi, C.L. Performance of DNA-based biomarkers for classification of adrenocortical carcinoma: A prognostic study. Eur. J. Endocrinol. 2023, 189, 262-270. [CrossRef]

29. Mete, O .; Gucer, H .; Kefeli, M .; Asa, S.L. Diagnostic and Prognostic Biomarkers of Adrenal Cortical Carcinoma. Am. J. Surg. Pathol. 2018, 42, 201-213. [CrossRef]

30. Salama, M.F .; Liu, M .; Clarke, C.J .; Espaillat, M.P .; Haley, J.D .; Jin, T .; Wang, D .; Obeid, L.M .; Hannun, Y.A. PKCa is required for Akt-mTORC1 activation in non-small cell lung carcinoma (NSCLC) with EGFR mutation. Oncogene 2019, 38, 7311-7328. [CrossRef]

31. Martini, M .; De Santis, M.C .; Braccini, L .; Gulluni, F .; Hirsch, E. PI3K/AKT signaling pathway and cancer: An updated review. Ann. Med. 2014, 46, 372-383. [CrossRef] [PubMed]

32. Moon, H .; Ro, S.W. MAPK/ERK Signaling Pathway in Hepatocellular Carcinoma. Cancers 2021, 13, 3026. [CrossRef] [PubMed]

33. Xiao, C.L .; Yin, W.C .; Zhong, Y.C .; Luo, J.Q .; Liu, L.L .; Liu, W.Y .; Zhao, K. The role of PI3K/ Akt signalling pathway in spinal cord injury. Biomed. Pharmacother. 2022, 156, 113881. [CrossRef] [PubMed]

34. Parzefall, T .; Schnoell, J .; Monschein, L .; Foki, E .; Liu, D.T .; Frohne, A .; Grasl, S .; Pammer, J .; Lucas, T .; Kadletz, L .; et al. PRKCA Overexpression Is Frequent in Young Oral Tongue Squamous Cell Carcinoma Patients and Is Associated with Poor Prognosis. Cancers 2021, 13, 2082. [CrossRef]

35. Jiang, H .; Fu, Q .; Song, X .; Ge, C .; Li, R .; Li, Z .; Zeng, B .; Li, C .; Wang, Y .; Xue, Y .; et al. HDGF and PRKCA upregulation is associated with a poor prognosis in patients with lung adenocarcinoma. Oncol. Lett. 2019, 18, 4936-4946. [CrossRef]

36. Wang, Y .; Liu, S .; Yan, Y .; Li, S .; Tong, H. SPARCL1 promotes C2C12 cell differentiation via BMP7-mediated BMP/TGF-ß cell signaling pathway. Cell Death Dis. 2019, 10, 852. [CrossRef]

37. Lin, R .; Li, X .; Wu, S .; Qian, S .; Hou, H .; Dong, M .; Zhang, X .; Zhang, M. Suppression of latent transforming growth factor-ß (TGF-ß)-binding protein 1 (LTBP1) inhibits natural killer/T cell lymphoma progression by inactivating the TGF-ß/Smad and p38MAPK pathways. Exp. Cell Res. 2021, 407, 112790. [CrossRef]

38. Cai, R .; Wang, P .; Zhao, X .; Lu, X .; Deng, R .; Wang, X .; Su, Z .; Hong, C .; Lin, J. LTBP1 promotes esophageal squamous cell carcinoma progression through epithelial-mesenchymal transition and cancer-associated fibroblasts transformation. J. Transl. Med. 2020, 18, 139. [CrossRef]

39. Fu, X .; Zhang, P .; Song, H .; Wu, C .; Li, S .; Li, S .; Yan, C. LTBP1 plays a potential bridge between depressive disorder and glioblastoma. J. Transl. Med. 2020, 18, 391. [CrossRef]

40. Frazzi, R. BIRC3 and BIRC5: Multi-faceted inhibitors in cancer. Cell Biosci. 2021, 11, 8. [CrossRef]

41. Margulis, V .; Lotan, Y .; Shariat, S.F. Survivin: A promising biomarker for detection and prognosis of bladder cancer. World J. Urol. 2008, 26, 59-65. [CrossRef] [PubMed]

42. 2. He, D .; Huang, K .; Liang, Z. Prognostic value of baculoviral IAP repeat containing 5 expression as a new biomarker in lung adenocarcinoma: A meta-analysis. Expert Rev. Mol. Diagn. 2021, 21, 973-981. [CrossRef]

43. Hagenbuchner, J .; Kiechl-Kohlendorfer, U .; Obexer, P .; Ausserlechner, M.J. BIRC5/Survivin as a target for glycolysis inhibition in high-stage neuroblastoma. Oncogene 2016, 35, 2052-2061. [CrossRef] [PubMed]

44. Kahm, Y.J .; Kim, R.K. BIRC5: A novel therapeutic target for lung cancer stem cells and glioma stem cells. Biochem. Biophys. Res. Commun. 2023, 682, 141-147. [CrossRef] [PubMed]

45. Adinew, G.M .; Messeha, S .; Taka, E .; Soliman, K.F.A. The Prognostic and Therapeutic Implications of the Chemoresistance Gene BIRC5 in Triple-Negative Breast Cancer. Cancers 2022, 14, 5180. [CrossRef]

46. Wang, P .; Wang, Y .; Jiang, Y .; Li, M .; Li, G .; Qiao, Q. Immune Cluster and PPI Network Analyses Identified CXCR3 as a Key Node of Immunoregulation in Head and Neck Cancer. Front. Oncol. 2020, 10, 564306. [CrossRef]

47. . Inci, N .; Akyildiz, E.O .; Bulbul, A.A .; Turanli, E.T .; Akgun, E .; Baykal, A.T .; Colak, F .; Bozaykut, P. Transcriptomics and Proteomics Analyses Reveal JAK Signaling and Inflammatory Phenotypes during Cellular Senescence in Blind Mole Rats: The Reflections of Superior Biology. Biology 2022, 11, 1253. [CrossRef]

18. Kupsa, T .; Vanek, J .; Zak, P .; Jebavy, L .; Horacek, J.M. Serum Levels of Cytokines and Cytokine Receptors Are Associated with Outcome in Newly Diagnosed AML Patients. Blood 2016, 128, 5246. [CrossRef]

49. Garris, C.S .; Luke, J.J. Dendritic Cells, the T-cell-inflamed Tumor Microenvironment, and Immunotherapy Treatment Response. Clin. Cancer Res. 2020, 26, 3901-3907. [CrossRef]

50. Duan, J .; Yu, Y .; Li, Y .; Li, Y .; Liu, H .; Jing, L .; Yang, M .; Wang, J .; Li, C .; Sun, Z. Low-dose exposure of silica nanoparticles induces cardiac dysfunction via neutrophil-mediated inflammation and cardiac contraction in zebrafish embryos. Nanotoxicology 2016, 10, 575-585. [CrossRef]

51. Zhang, C .; Dang, D .; Wang, Y .; Cong, X. A Nomogram Combining a Four-Gene Biomarker and Clinical Factors for Predicting Survival of Melanoma. Front. Oncol. 2021, 11, 593587. [CrossRef] [PubMed]

52. Qu, J .; Jiang, M .; Wang, L .; Zhao, D .; Qin, K .; Wang, Y .; Tao, J .; Zhang, X. Mechanism and potential predictive biomarkers of immune checkpoint inhibitors in NSCLC. Biomed. Pharmacother. 2020, 127, 109996. [CrossRef] [PubMed]

53. Saleh, R.R .; Peinado, P .; Fuentes-Antrás, J .; Pérez-Segura, P .; Pandiella, A .; Amir, E .; Ocaña, A. Prognostic Value of Lymphocyte- Activation Gene 3 (LAG3) in Cancer: A Meta-Analysis. Front. Oncol. 2019, 9, 1040. [CrossRef] [PubMed]

54. Yang, L .; Wei, S .; Zhang, J .; Hu, Q .; Hu, W .; Cao, M .; Zhang, L .; Wang, Y .; Wang, P .; Wang, K. Construction of a predictive model for immunotherapy efficacy in lung squamous cell carcinoma based on the degree of tumor-infiltrating immune cells and molecular typing. J. Transl. Med. 2022, 20, 364. [CrossRef] [PubMed]

55. Kao, C .; Powers, E .; Datto, M.B .; Green, M .; Harrison, M.R .; Armstrong, A.J .; George, D.J .; Clarke, J.M .; Strickler, J.H .; Zhang, T. Tumor mutational burden (TMB) as a predictive biomarker of immune checkpoint blockade (ICB) in metastatic solid tumors. J. Clin. Oncol. 2020, 38, 80. [CrossRef]

56. Nandimandalam, S .; Sharma, N. Correlation between tissue PD-L1, TMB, and blood PD-L1, MSI biomarkers in patients with advanced-stage non-small cell lung cancer (NSCLC). J. Clin. Oncol. 2020, 38, e21564. [CrossRef]

57. Lentz, R.W .; Friedrich, T .; Hu, J .; Leal, A.D .; Kim, S.S .; Davis, S.L .; Purcell, T .; Messersmith, W.A .; Lieu, C.H. Tissue tumor mutational burden (TMB) as a biomarker of efficacy with immune checkpoint inhibitors (ICI) in metastatic gastrointestinal (mGI) cancers. J. Clin. Oncol. 2021, 39, e14559. [CrossRef]

58. Huang, X .; Tse, J.Y .; Protopopov, A .; Russell, M .; Weeraratne, D .; Ring, J.E .; Bjonnes, A .; Pei, S .; Sun, R .; Lvova, M .; et al. Character- ization of tumor mutation burden (TMB) and microsatellite instability (MSI) interplay for gastroesophageal adenocarcinoma (GA) and colorectal carcinoma (CRC). J. Clin. Oncol. 2018, 36, 22. [CrossRef]

59. Jardim, D.L .; Goodman, A .; de Melo Gagliato, D .; Kurzrock, R. The Challenges of Tumor Mutational Burden as an Immunotherapy Biomarker. Cancer Cell 2021, 39, 154-173. [CrossRef]

60. Wang, L .; Ge, J .; Lan, Y .; Shi, Y .; Luo, Y .; Tan, Y .; Liang, M .; Deng, S .; Zhang, X .; Wang, W .; et al. Tumor mutational burden is associated with poor outcomes in diffuse glioma. BMC Cancer 2020, 20, 213. [CrossRef]

61. Lai, G .; Zhong, X .; Liu, H .; Deng, J .; Li, K .; Xie, B. Development of a Hallmark Pathway-Related Gene Signature Associated with Immune Response for Lower Grade Gliomas. Int. J. Mol. Sci. 2022, 23, 11971. [CrossRef] [PubMed]

62. Shuayb, M .; Das, A .; Uddin, M.N. Long-term survival in recurrent adrenocortical cancer. BMJ Support. Palliat. Care 2019, 9, 47-50. [CrossRef] [PubMed]

63. Menefee, M .; Edgerly, M .; Velarde, M .; Herbert, K .; Fojo, A.T. The efficacy of combination chemotherapy with cisplatin and gemcitabine in patients with advanced adrenal cortical carcinoma (ACC). J. Clin. Oncol. 2006, 24, 12033. [CrossRef]

64. Subramanian, C .; Cohen, M.S. Identification of novel lipid metabolic biomarkers associated with poor adrenocortical carcinoma prognosis using integrated bioinformatics. Surgery 2022, 171, 119-129. [CrossRef] [PubMed]

65. Wilkerson, M.D .; Hayes, D.N. ConsensusClusterPlus: A class discovery tool with confidence assessments and item tracking. Bioinformatics 2010, 26, 1572-1573. [CrossRef]

66. Şenbabaoğlu, Y .; Michailidis, G .; Li, J.Z. Critical limitations of consensus clustering in class discovery. Sci. Rep. 2014, 4, 6207. [CrossRef] [PubMed]

67. Love, M.I .; Huber, W .; Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014, 15, 550. [CrossRef]

68. Engebretsen, S .; Bohlin, J. Statistical predictions with glmnet. Clin. Epigenet. 2019, 11, 123. [CrossRef]

69. Zhong, Z .; Wang, Y .; Yin, J .; Ni, S .; Liu, W .; Geng, R .; Liu, J .; Bai, J .; Yu, H. Identification of Specific Cervical Cancer Subtypes and Prognostic Gene Sets in Tumor and Nontumor Tissues Based on GSVA Analysis. J. Oncol. 2022, 2022, e6951885. [CrossRef]

70. Newman, A.M .; Liu, C.L .; Green, M.R .; Gentles, A.J .; Feng, W .; Xu, Y .; Hoang, C.D .; Diehn, M .; Alizadeh, A.A. Robust enumeration of cell subsets from tissue expression profiles. Nat. Methods 2015, 12, 453-457. [CrossRef]

71. Colaprico, A .; Silva, T.C .; Olsen, C .; Garofano, L .; Cava, C .; Garolini, D .; Sabedot, T.S .; Malta, T.M .; Pagnotta, S.M .; Castiglioni, I .; et al. TCGAbiolinks: An R/Bioconductor package for integrative analysis of TCGA data. Nucleic Acids Res. 2016, 44, e71. [CrossRef] [PubMed]

72. Mayakonda, A .; Lin, D .- C .; Assenov, Y .; Plass, C .; Koeffler, H.P. Maftools: Efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 2018, 28, 1747-1756. [CrossRef] [PubMed]

73. Yu, G .; Wang, L.G .; Han, Y .; He, Q.Y. clusterProfiler: An R package for comparing biological themes among gene clusters. OMICS 2012, 16, 284-287. [CrossRef] [PubMed]

74. Trujillo, J.A .; Sweis, R.F .; Bao, R .; Luke, J.J. T Cell-Inflamed versus Non-T Cell-Inflamed Tumors: A Conceptual Framework for Cancer Immunotherapy Drug Development and Combination Therapy Selection. Cancer Immunol. Res. 2018, 6, 990-1000. [CrossRef]

75. Maeser, D .; Gruener, R.F .; Huang, R.S. oncoPredict: An R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief. Bioinform. 2021, 22, bbab260. [CrossRef]

76. Bonacci, R .; Inserm, R.C .; Gigliotti, A .; Baudin, E .; Wion-Barbot, N .; Emy, P .; Bonnay, M .; Cailleux, A .; Nakib, I .; Schlumberger, M. Cytotoxic therapy with etoposide and cisplatin in advanced adrenocortical carcinoma. Br. J. Cancer 1998, 78, 546-549. [CrossRef]

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.