Exploration of Core Genes in ACTH-Independent Macronodular Adrenal Hyperplasia

Abstract

This study explores the core genes involved in the pathogenesis of ACTH-independent macronodular adrenal hyperplasia (AIMAH), so as to provide robust biomarkers for the clinical diagnosis and treatment of this disease. Gene Expression Omnibus (GEO) database was used to obtain GSE25031 microarray dataset. R package “limma” was applied to identify differentially expressed genes (DEGs) between AIMAH and normal samples. The Database for Annotation, Visualization and Integrated Discovery (DAVID) was employed to perform Gene Ontology (GO) annotation for the DEGs, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was conducted. A protein-protein interaction network (PPI) was constructed using the STRING online website and visualized using the Cytoscape software. The key modules and hub genes were then identified. Finally, Gene Set Enrichment Analysis (GSEA) enrichment analysis was carried out to find the signaling pathways of significant clinical value in AIMAH. A total of 295 DEGs between AIMAH and healthy samples were screened out, including 164 upregulated genes and 131 downregulated genes. Combining enrichment analysis and PPI network construction, there were 5 signifiant pathways and 10 hub genes, among which 3 genes (FOS, FOSB, and DUSP1) were identified as potential core genes of clinical significance in AIMAH. In conclusion, the 3 core genes, FOS, FOSB, and DUSP1, identified here might be potential biomarkers for AIMAH, and the current study is of guiding significance for clinical diagnosis and treatment of this disease.

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.

Conflict of Interest

The authors declare that they have no conflict of interest.

• References

• 1 Hsiao HP, Kirschner LS, Bourdeau I. et al. Clinical and genetic heterogeneity, overlap with other tumor syndromes, and atypical glucocorticoid hormone secretion in adrenocorticotropin-independent macronodular adrenal hyperplasia compared with other adrenocortical tumors. J Clin Endocrinol Metab 2009; 94: 2930-2937
• 2 Almeida MQ, Harran M, Bimpaki EI. et al. Integrated genomic analysis of nodular tissue in macronodular adrenocortical hyperplasia: progression of tumorigenesis in a disorder associated with multiple benign lesions. J Clin Endocrinol Metab 2011; 96: E728-E738
• 3 Oh H, Koh JM, Kim MS. et al. A case of ACTH-producing pheochromocytoma associated with pregnancy. Endocr J 2003; 50: 739-744
• 4 Hu WL. Adrenal cortex-medulla mixed lesions. Chin J Laparoscop Urol (Electronic Edition) 2012; 6: 84-86
• 5 Vezzosi D, Cartier D, Régnier C. et al. Familial adrenocorticotropin-independent macronodular adrenal hyperplasia with aberrant serotonin and vasopressin adrenal receptors. Eur J Endocrinol 2007; 156: 21-31
• 6 Reinhard C, Saeger W, Schubert B. Adrenocortical nodules in post-mortem series. Development, functional significance, and differentiation from adenomas. Gen Diagn Pathol 1996; 141: 203-208
• 7 Chen WX, Lin Y, Yang CH. et al. Diagnosis and treatment of 12 cases of ACTH-independent adrenal nodular hyperplasia. Chin J Endocrinol Metab 2011; 5-8
• 8 Rogoff D, Bergada I, Venara M. et al. Intermittent hyperaldosteronism in a child due to an adrenal adenoma. Eur J Pediatr 2001; 160: 114-116
• 9 Zhao X, Zhang XB, Li HZ. Screening and preliminary validation of microRNAs associated with corticotropin-independent adrenal nodular hyperplasia. Chin J Urol 2015; 36: 921-924
• 10 Romero DG, Rilli S, Plonczynski MW. et al. Adrenal transcription regulatory genes modulated by angiotensin II and their role in steroidogenesis. Physiol Genomics 2007; 30: 26-34
• 11 Takahashi H, Ikematsu K, Tsuda R. et al. Increase in dual specificity phosphatase 1, TGF-beta stimulated gene 22, domain family protein 3 and Luc7 homolog ( S. cerevisiae )-like messenger RNA after mechanical asphyxiation in the mouse lung. Leg Med (Tokyo) 2009; 11: 181-185

Publication History

Received: 15 January 2022

Accepted after revision: 17 March 2022

Article published online:
09 May 2022