Key words
ACTH-independent macronodular adrenal hyperplasia - AIMAH - computational biology
- genes
Introduction
Adrenocorticotrophic hormone (ACTH)-independent macronodular adrenal hyperplasia
(AIMAH) is a clinically rare disease, which is characterized by adrenal cortical
hyperplasia. It accounts for<1% of all cases developing
Cushing’s syndrome (CS) [1]. Over the
past few years, the diagnostic imaging techniques are gradually developing and have
continuously increased the detection rate of AIMAH. AIMAH can present with diverse
clinical symptoms, with changes in the cortisol secretion level, accordingly,
treatment for AIMAH also varies. The clinical presentations in patients with AIMAH
are sometimes mild or concealed without typical symptoms, and in some patients,
there is no significant abnormality in endocrine examination or computed tomography
(CT) scan, making it challenging to diagnose this disease clinically. Furthermore,
it has been suggested that in the context of chronic polyclonal hyperplasia, the
larger adrenal lesions can accumulate more genomes leading to transcriptional
abnormalities and ultimately show abnormal expression of oncogenic pathways. The
study by Almeida et al. [2] supported this
view by analyzing the integrated transcriptomic and genomic data of AIMAH. Under
this background, a well understanding towards the pathogenesis of AIMAH and the
expression profiles of relevant genes is conducive to formulating efficient
diagnostic and therapeutic strategies for AIMAH.
In the present study, microarray data generated by the GPL625 platform were
downloaded from the Gene Expression Omnibus (GEO) database, and differentially
expressed genes (DEGs) between AIMAH and normal samples were identified using R
analysis. In the meantime, a series of analyses, including Gene Ontology (GO)
functional annotation, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway
enrichment analysis, protein-protein interaction (PPI) analysis and Gene Set
Enrichment Analysis (GSEA), were performed to unravel the candidate core genes of
clinical significance in AIMAH. This study just computationally analyzes the
microarray data obtained from the GEO database, but findings in this study will
promote the progress of research on AIMAH and shed new lights on the clinical
diagnosis and treatment of this disease.
Materials and Methods
Ethics approval
All information from the GEO database was deidentified and no personal
identifying information was used in our analysis, so informed consent was not
required in our study.
Microarray data
The GSE25031 dataset generated by the microarray platform GPL6255 via Almeida et
al. [2] was downloaded from the GEO
database. This dataset involved 7 adrenal nodules collected in 2 AIMAH cases as
well as 3 total RNA available pools collected in human adrenal gland. Patient 1
was the 42-year-old man who experienced weight gain (60 lbs in more than 10
years) and mild hypertension, violaceous striae, and plethoric face, and was
hospitalized at the NIH Clinical Center to receive CS examination. Hormonal
assessment suggested no suppression on plasma cortisol following dexamethasone
(DEX, 1 mg) test overnight as well as following DEX-ovine
corticotropin-releasing hormone (CRH) test, suppression on ACTH content
(<5 pg/ml), abnormality in midnight serum cortisol
content (11.6 μg/dl), and great 24-hour urinary free
cortisol content (159.1 μg;
normal,<90 μg/24 h). Meanwhile, CT
images suggested enlargement of both adrenal glands and presence of
macronodules. The diagnosis of ACTH-independent CS was made in this patient, and
he received bilateral adrenalectomy. Patient 2 was the 42-year-old woman who
experienced bruising susceptibility, secondary amenorrhea for 4 years, weight
gain (25 lbs), hirsutism, behavioral alterations, and muscle weakness. According
to biochemical assessment, there was no suppression on plasma cortisol following
DEX (1 mg) test overnight as well as following DEX-ovine CRH test, ACTH
content was suppressed (<5 pg/ml), and 24-hour urinary
free cortisol content was high (270 μg;
normal,<90 μg/24 h). On CT images,
macronodules were observed in both adrenal glands. The diagnosis of
ACTH-independent CS was made, and the woman received laparoscopic bilateral
adrenalectomy. The sample data could be accessed in original study [2].
Data preprocessing and identification of DEGs
Raw data from the GSE25031 dataset and annotation data from the GPL6255 platform
were obtained. First of all, raw data were preprocessed using the RMA algorithm
of R package “affy”. Thereafter, DEGs between AIMAH and healthy
samples were identified using R package “limma”, then
corresponding adjusted-P and fold change (FC) values were calculated. Later, a
volcano plot was generated according to those identified DEGs using R package
“ggplot2”. Finally, the top 50 upregulated and downregulated
genes were mapped into a heatmap using the “Heml” software.
GO functional annotation and KEGG pathway enrichment analyses
GO annotation describes protein functions from three aspects, namely, Biological
Process (BP), Cellular Component (CC), and Molecular Function (MF). KEGG is a
database that integrates genomic, chemical, and systemic functional information.
Moreover, this work applied the Database for Annotation, Visualization, and
Integrated Discovery (DAVID) tool for GO functional annotation and KEGG pathway
enrichment analyses on DEGs.
Gene set enrichment analysis (GSEA)
GSEA was performed to identify the significantly different pathways between AIMAH
and healthy samples using R package “ClusterProfiler”. The
significance thresholds were set at adjusted p<0.05 and false discovery
rate (FDR)<0.25.
Construction of a protein-protein interaction (PPI) network and
identification of hub genes
STRING is a database that involves both known and predicted PPIs. Here, a PPI
network was constructed based on DEGs that were projected onto the STRING
database. The interaction score>0.4 was considered statistically
significant. Cytoscape was employed to visualize the intermolecular
interactions. Meanwhile, the CytoHubba plug-in was utilized to determine hub
genes in the PPI network, according to the Degree, Radiality Centrality,
Closeness Centrality and EPC.
Results
Identification of DEGs
Box-plot, principal component analysis (PCA) and uniform manifold approximation
and projection (UMAP) plots were generated based on sample data derived from the
GSE25031 dataset ([Fig. 1a–c]).
The Box-plot revealed good normalization between samples, as demonstrated by
distribution of median of each sample to a horizontal line basically. On the
contrary, the PCA and UMAP plots showed significant inter-group difference.
Besides, differential analysis revealed 295 DEGs between AIMAH and healthy
samples, including 164 upregulated and 131 downregulated genes. Accordingly, a
volcano plot and a heatmap of the top 50 upregulated and down-regulated DEGs
were drawn ([Fig. 2a, b]).
Fig. 1 Box-plot, PCA, and UMAP plots of the sample data derived
from the GSE25031 dataset.
Fig. 2 Volcano plot and Heat map of the top 50 up-/downregulated
genes of DEGs in AIMAH.
Enrichment analysis
GO functional annotation was performed using the DAVID tool. As a result, the
upregulated genes were mainly annotated into response to steroid hormones,
response to metal ion, response to ketone, response to corticosteroids, response
to cadmium (BP terms); hemoglobin complex and haptoglobin-hemoglobin complex (CC
terms); and heme binding, tetrapyrrole binding, oxygen binding, haptoglobin
binding and RNA polymerase II activating transcription factor binding (MF terms)
([Fig. 3a–d]). Typically, the
upregulated DEGs were mainly annotated into BP and MF terms. KEGG pathway
enrichment analysis demonstrated that the upregulated DEGs were mostly enriched
into the ADM, DUSP1, CYB5A, and CYP3A5 pathways ([Table 1]). GSEA was later conducted to
identify the signaling pathways significantly different between the two groups
upon the thresholds of adjusted p<0.05 and FDR<0.25 in MSigDB
Collection (c2.all.v7.0). There were 5 pathways of statistical significance,
including nerve growth factor (NGF) stimulated transcription, protein
localization, mitogen-activated protein kinase (MAPK) signaling pathway, cystic
fibrosis (CFTR), and metabolism of polyamines ([Fig. 4a–e]).
Fig. 3 GO functional annotation of the upregulated genes and the
results in BPs, CCs and MFs.
Fig. 4 GSEA enrichment analysis of two groups and 5 pathways of
statistical significance.
Table 1 KEGG pathway enrichment analysis of the most
activated pathways involved in the upregulated genes.
ID
|
Term
|
Count
|
p-Vaule
|
Genes
|
GO:0048545
|
Response to steroid hormone
|
11
|
6.94409E-08
|
ADM/DUSP1/FOS/FOSB/NR4A2/RORA/SPP1/ZFP36/AKR1C3/SLIT2/TXNIP
|
GO:1901654
|
Response to ketone
|
8
|
3.04287E-07
|
DUSP1/FOS/FOSB/CCL21/SPP1/AKR1C3/SLIT2/TXNIP
|
GO:0010038
|
Response to metal ion
|
10
|
4.16809E-07
|
CYB5A/CYP11B2/DUSP1/FOS/FOSB/MT1A/MT2A/AKR1C3/TXNIP/SLC40A1
|
GO:0031960
|
Response to corticosteroid
|
7
|
1.33375E-06
|
ADM/DUSP1/FOS/FOSB/ZFP36/AKR1C3/SLIT2
|
GO:0046686
|
Response to cadmium ion
|
5
|
2.43448E-06
|
CYB5A/FOS/MT1A/MT2A/AKR1C3
|
GO:0020037
|
Heme binding
|
5
|
0.00010949
|
CYB5A/CYP3A5/CYP11B2/HBA1/HBB
|
GO:0046906
|
Tetrapyrrole binding
|
5
|
0.000153175
|
CYB5A/CYP3A5/CYP11B2/HBA1/HBB
|
GO:0019825
|
Oxygen binding
|
3
|
0.000269234
|
CYP3A5/HBA1/HBB
|
Construction of the PPI network and identification of hub genes
A PPI network based on DEGs was constructed using the STRING online tool and
visualized using the Cytoscape software. According to the Degree, Radiality
Centrality, Closeness Centrality and EPC, those identified DEGs were ranked. The
top 15 genes selected according to each standard as hub genes are listed in
[Table 2]. Thereafter, interaction
analysis was performed, and 10 overlapped hub genes were obtained, including
FOS, NR4A2, DUSP1, FOSB, ZFP36, PPP1R15 A, GADD45B, KLF6, SPP, and
CYP11B2 ([Fig. 5a, b]).
Fig. 5 PPI network in the DEGs and interaction according to the
Degree, Radiality Centrality, Closeness Centrality, and EPC.
Table 2 The top 15 genes were screened out by each
standard as hub genes according to the Degree, Radiality Centrality,
Closeness Centrality, and EPC.
Closeness
|
Degree
|
EPC
|
Radiality
|
FOS
|
FOS
|
FOS
|
FOS
|
NR4A2
|
PPP1R15A
|
FOSB
|
NR4A2
|
DUSP1
|
NR4A2
|
DUSP1
|
DUSP1
|
FOSB
|
ZFP36
|
ZFP36
|
FOSB
|
ZFP36
|
DUSP1
|
PPP1R15A
|
ZFP36
|
PPP1R15A
|
FOSB
|
GADD45B
|
GADD45B
|
GADD45B
|
GADD45B
|
KLF6
|
SPP1
|
KLF6
|
KLF6
|
NR4A2
|
PPP1R15A
|
SPP1
|
CYB5A
|
RGS1
|
KLF6
|
CYP11B2
|
RGS1
|
CYP11B2
|
CYP11B2
|
MT2A
|
SPP1
|
SPP1
|
MT2A
|
RGS1
|
FGF9
|
ERN1
|
MYB
|
CYB5A
|
CYP3A5
|
MT2A
|
PRKCA
|
MYB
|
AKR1C3
|
RORA
|
RORA
|
PRKCA
|
CYP11B2
|
CYB5A
|
TUBB3
|
Identification of possible core genes
Combined with the enrichment analysis results, 3 of the 10 hub genes, namely,
FOS, FOSB and DUSP1, were identified as the potential core genes of clinical
significance in AIMAH.
Discussion
Currently, the etiology of AIMAH remains elusive, while it was previously thought
to
be a gradual transition from ACTH-dependent to ACTH-independent. However, more and
more scholars hold a different opinion. Apart from ACTH, research has revealed that
AIMAH can also be induced by the ectopic expression of receptors of arginine
vasopressin (AVP), gastric inhibitory peptide (GIP), and catecholamine (CA). All
these receptors belong to the G protein-coupled receptors, which will also stimulate
the secretion of adrenal hormone when the above hormones act on them, finally
leading to adrenal hyperplasia [1].
Additionally, ACTH secreted by adrenal medullary chromaffin tissues can be applied
over the cortex [3], and the adrenaline
secreted can directly excite β receptor, induce the hypersecretion of renin
and thereby advance the globular accretion [4]. Vezzosi et al. [5] revealed that
AIMAH was an autosomal dominant disorder. In terms of pathology, there is a trend
of
adrenal capsule vascular stiffening and decreasing blood supply with age, which
results in local cortical atrophy, thus stimulating the
hypothalamic-pituitary-adrenal (HPA) axis and inducing peripheral cell hyperplasia
in a feed-back fashion [6]. At present, the
diagnosis of AIMAH is also facing great challenges. In clinical practice, most
patients present with hypertension alone, or dizziness, headache, and diabetes,
without significant specific symptoms. This may be attributed to the fact that
cortical hyperplasia causes an increased number of cells, but the biosynthetic
capacity of a cell fails to be concurrently enhanced [7]. This phenomenon is also considered as a
subclinical syndrome. In a majority of patients with adrenal hyperplasia, the
clinical presentations are diverse and complicated, including either primary
aldosteronism or hypercortisolisme manifestations. In adrenal hyperplasia, in
addition to the elevation of lead hormone, various corticosteroids can increase
simultaneously to varying degrees. Additionally, there may be some atypical symptoms
due to comorbidities such as nervous system or endocrine system disease and oral
medicine. However, some patients have normal endocrine levels in clinical
examination, which may be associated with the periodical secretion fashion of
hormones [8]. This demonstrates that the
elevation of hormone levels in patients with adrenal hyperplasia can be
impersistent. Besides, the consistent biosynthetic capacity of a cell following the
increase in cell number can also be an explanation [9]. Moreover, the cells that undergo hyperplasia can be interstitial
cells without endocrine functions. Neuroendocrine, environmental, and emotional
factors can also affect the secretion levels of hormones. Under this background,
multiple blood examinations are required in some cases to determine the levels of
hormones. In patients with normal endocrine levels, adrenal hyperplasia cannot be
completely excluded, and further imaging examinations such as CT are required. The
final diagnosis of AIMAH is made depending on pathological results. Additionally,
there are some patients presenting with clinical presentations or abnormality in
endocrine manifestations alone, without any obvious evidence of abnormal images.
This may be because that the mild adrenal lesions are only shown in microscopy, and
the pathological changes are too mild to induce any morphological changes.
Pathological changes may precede morphological changes, given the detection of
significant cortical thickening and hyperplasia on CT scans.
At the microscopic level, Zhao et al. [9]
identified 12 miRNAs differentially expressed between AIMAH and normal adrenal
tissues, including 7 upregulated genes (hsa-miR-663, hsa-miR-498, hsa-miR-638,
hsa-miR-501-5p, hsa-miR-585, hsa-miR-557, and hsa-miR-144) and 5 downregulated ones
(hsa-miR-744, hsa-miR-143, hsa-miR-26a, hsa-miR-22, and hsa-miR-29a). They also
reported 4 miRNAs of clinical significance, including hsa-miR-663, hsa-miR-498,
hsa-miR-557, and hsa-miR-744. In the present study, a total of 295 DEGs were
identified between AIMAH and healthy samples, involving 164 upregulated and 131
downregulated genes. GO functional annotation based on the DAVID tool revealed that
the up-regulated genes were mainly enriched in BP and MF terms. As demonstrated by
KEGG pathway enrichment analysis, the most activated pathways were ADM, DUSP1,
CYB5A, and CYP3A5. Further GSEA revealed that 5 pathways were significantly
different between two groups, including NGF stimulated transcription, protein
localization, MAPK signaling pathway, cystic fibrosis (CFTR) and metabolism of
polyamines. Furthermore, a PPI network was established based on the identified DEGs.
Using the CytoHubba plug-in, 10 hub genes (FOS, NR4A2, DUSP1, FOSB, ZFP36, PPP1R15A,
GADD45B, KLF6, SPP, and CYP11B2) were identified according to the Degree, Radiality
Centrality, Closeness Centrality and EPC. Eventually, 3 possible core genes, FOS,
FOSB, and DUSP1, were obtained combined with the enrichment results. Damina et al.
[10] suggested that FOS and FOSB
upregulated the expression of 11β-hydroxylase and aldosterone synthase.
Takahashi et al. [11] reported that DUSP1
deactivated mitogen-activated protein kinase (MAPK) and played a role in cell
proliferation in mouse models. Nonetheless, the roles of the 3 core genes and the
5
significant signaling pathways in AIMAH should be further identified. Inevitably,
certain limitations should be noted in this work. For example, the sample size was
small, which might be ascribed to the relatively rare clinical AIMAH cases recorded
and the even limited expression profile data of AIMAH. AIMAH study is usually case
report and exploration of how some gene mutations are involved in the related
pathways during the formation of AIMAH. However, with the development of
bioinformatics and the expansion of sample size, the conclusions of this study
should be further verified in the future.
Conclusion
This study first identified 295 DEGs between AIMAH and normal samples and selected
5
key pathways and 10 hub genes. Three core genes, including FOS, FOSB, and DUSP1,
were eventually identified in this work, which were the potential biomarkers for
AIMAH. Our results are of guiding significance for the clinical diagnosis and
treatment of this disease.
Notice
This article was changed according to the following Erratum
on July 22nd 2022.
Erratum
The abbreviation for “Gene set enrichment analysis” has to be
“GSEA”, this has been corrected in the Abstract and in the
Materials and Methods section.