Identification of key miRNAs and their targets in peripheral blood mononuclear cells of IgA nephropathy using bioinformatics analysis.

BACKGROUND: Currently, renal biopsy is the gold standard for clinical diagnosis and evaluation the degrees of IgA nephropathy. However, renal biopsy is an invasive examination and not suitable for long-term follow-up IgA nephropathy. The activation of peripheral blood mononuclear cells (PBMCs) are related to IgA nephropathy, but the key molecular marker and target of PBMCs for evaluating the progression and prognosis of IgA nephropathy is still unclear.
METHODS: We downloaded gene expression omnibus series 25590 (GSE25590) datasets, of which PBMCs from IgA nephrology (IgAN) and healthy patients, from the gene expression omnibus (GEO) database. Differentially expressed miRNAs (DEMs) between IgAN and healthy patients were identified. The Funrich software was used to predict the differentially expressed genes (DEGs). Gene ontology (GO) and kyoto encyclopedia of genes and genomes (KEGG) analyzes of overlapping genes were analyzed at the function level on DAVID 6.8. We used search Tool for the retrieval of interacting genes (STRING) online database constructed the protein-protein interaction (PPI) network. Then we further analyzed the hub genes by Cytoscape software and the hub miRNA by TargetScan.
RESULTS: We identified 418 DEMs from the GSE25590 datasets. The upstream transcription factors SP1 regulates most DEMs. According to the GO and KEGG results, the DEGs were enriched in the MAPK signaling pathway and small GTPase mediated signal transduction. SYN1, SYT4, RBFOX1, KCNC1, VAMP2, FBXO11, ASB9, SYT9, KLHL5, and KRAS were identified as hub genes. Hsa-miR-532-5p, hsa-miR-92a, hsa-miR-328, hsa-miR-137, hsa-miR-153, hsa-miR-9-5p, hsa-miR-140-5p, hsa-miR-217, hsa-miR-155, and hsa-miR-212 were predicted as hub miRNAs.
CONCLUSIONS: The DEMs and DEGs re-analysis provided potential key genes and hub miRNA of PMBCs, which may help to monitor the happening and prognosis of IgAN.

Introduction

IgA nephropathy is the most common primary glomerular disease in Asia, accounting for 1/3 of all renal biopsy cases.[ Epidemiological data show that more than 30% of IgA nephropathy will eventually develop into end stage renal disease.[ Currently, renal biopsy is the gold standard for clinical diagnosis and evaluation the degrees of IgA nephropathy.[ However, renal biopsy is an invasive examination and is not suitable to biopsy patients without an indication with the sole purpose to observe the renal tissue changes in IgA nephropathy. So finding the key molecules during the occurrence and development of IgA nephropathy can provide new markers for the IgA nephropathy assessment. Also, it is still worldwide scarcity of centers having the capacity to perform routine renal biopsies.

To date, the pathogenesis of IgA nephropathy is still not fully understood. Generally, the pathogenesis of IgA nephropathy includes 3 links, the formation of pathogenic IgA1 molecules and the deposition of IgA1 in the mesangial region of the glomerulus, which triggers local inflammation in the glomerulus, leading to the formation and development of IgA nephropathy.[ Most researchers accepted the “Four-Hit Hypothesis,” which including abnormal glycosylation process of IgA1, autoantibody formation, immune complex formation, immune complexes and polymeric IgA1 glomerular deposition. Observation after kidney transplantation found that patients with IgA nephropathy can still relapse after kidney transplantation.[ Studies have found that after the donor kidney of IgA nephropathy is implanted into the body of non-IgA nephropathy, the IgA deposited in the glomeruli will disappear, suggesting that IgA nephropathy may be related to abnormal IgA immune system and peripheral blood leucocytes.[ In particular, studies have found excessive activation of peripheral blood mononuclear cells (PBMC) from patients with IgA nephropathy.[ This indicates that peripheral blood mononuclear cells are related to IgA nephropathy, and may be used as a marker for evaluating the progression and prognosis of IgA nephropathy.

MicroRNAs are involved in the pathological process of many diseases including renal diseases.[ However, the role of microRNAs in the pathological process of IgA nephropathy is still unclear, especially in peripheral blood mononuclear cells. Therefore, exploring the relationship between the changes of microRNAs in peripheral blood and IgA nephropathy may bring new markers for the diagnosis, treatment and prognosis of IgA.

In recent years, the microarray technique combined with bioinformatics analysis has provided a powerful tool for the research between genes and renal diseases. In current study, we reanalyzed the miRNA (PBMC from IgAN patients) datasets from the gene expression omnibus (GEO) database to identify the differentially expressed miRNAs (DEMs) and differentially expressed genes (DEGs) between IgAN and normal human peripheral blood mononuclear cells. The enrichment of GO terms, pathways and protein-protein interaction (PPI) screen were also analyzed to find potential hub genes in the development of IgAN. We hypothesize that the miRNA of PMBCs may help to monitor the happening and prognosis of IgAN.

Materials and methods

Identification of differentially expressed miRNAs

The Gene Expression Omnibus (free on www.ncbi.nlm.nih.gov/geo) provides original submitter-offered records and curated datasets. We searched the datasets by the term “IgAN” and “peripheral blood mononuclear cells” from GEO database. And we selected the dataset with human samples and peripheral blood mononuclear cells: GSE25590. GSE25590 was submitted by Sharon Natasha Cox. GSE25590 was obtained from the GPL7731 platform (Agilent-019118 Human miRNA Microarray 2.0 G4470B). There were 7 normal and 7 IgAN patients were enrolled in this study and a total of 14 peripheral blood mononuclear cells samples were applied for this research.

We used the TXT files of raw data files for the analysis through GEO2R. The |logFC| > 1.0 and P < .05 were used as the cut-off criteria.

Prediction potential upstream transcription factors

We used FunRich software (Version 3.1.3) to predict the upstream transcription factors of DEMs. The screened up-regulation and down-regulation of DEMs were inputted FunRich, the result showed the top-10 predicted transcription factors.

Gene ontology and kyoto encyclopedia of genes and genomes enrichment analysis

We used FunRich software to predict the up-regulation and down-regulation of DEGs by up-regulation and down-regulation of DEMs respectively. Then we used DAVID 6.8 (https://david.ncifcrf.gov/) to analyze the GO enrichment and KEGG pathway of DEGs. After acquiring the results of GO enrichment and KEGG pathway, we used the online tool of Image GP (http://www.ehbio.com/ImageGP) to create the bubble images. We set the P < .05 and κ=1.0 as the statistically significant difference.

Protein-protein interaction network construction and analysis

The PPI network construction is a useful tool for exploring the potential inter-relationship among genes. We used search tool for the retrieval of interacting genes (STRING) (http://www.string-db.org/) database for the PPI network construction. We mapped all the 5276 DEMs into STRING and acquired the nodes of the network. Then we analyzed the data by Cytoscape software. The CytoHubba (version 0.1) plugin Cytoscape is used for identifying key genes from PPI network. We got the top 10 genes by the method of degree calculation. In order to further screen the miRNA-gene pairs of the hub genes, we searched the hub genes in TargetScan website (www.targetscan.org).

Ethical approval

All data in this study were obtained from open, public databases, so ethical approval was not necessary.

Results

Identification of differentially expressed miRNAs in peripheral blood mononuclear cells of IgA nephrology

GSE25590 contains 7 samples of peripheral blood mononuclear cells from patients with IgAN and 7 samples from healthy people. The miRNAs were detected by microRNA Array. We conducted a GEO2R analysis and there are 418 DEMs, of which 103 DEMs were up-regulated and 315 DEMs were down regulated. The volcano plot and heatmap of these DEMs were showed in Figure 1.

Figure 1

The expression of differentially expressed miRNA (DEMs) in GSE25590. (A) Volcano plot showing DEMs. The red dots represent up-regulated and the blue dots represent down-regulated miRNAs. (B) Heatmap showing the relative expression level of all the DEMs in each group. Red represents relatively high and blue represents relatively low.

Prediction upstream transcription factors and downstream target genes of differentially expressed miRNAs

In our research, we used Funrich software to predict upstream transcription factors of DEMs. The Figure 2 A and B showed the top10 up-regulated and 10 down-regulated predicted transcription factors respectively. According to the results, SP1 regulates most DEMs.

Figure 2

Predicted transcription factors of DEMs. (A) Top 10 transcription factors of up-regulated DEMs; (B) Top 10 transcription factors of down-regulated DEMs.

We also used Funrich software to predict downstream target genes of DEMs. There were 3754 targeted genes predicted by down-regulated DEMs and 1522 targeted genes predicted by up-regulated DEMs. Table 1 showed the top 10 up and down regulated targeted genes.

Table 1

Top 10 targeted genes predicted by up- regulated and down-regulated DEMs.

Predicted by up-regulated DEMs	Predicted by down-regulated DEMs
CDC42EP3	PCNP
PFN2	KCNT2
IGIP	KRTAP19–6
SHISA7	MSR1
RASL11B	NAT8L
BTC	ERBB4
MRPL34	NUP107
GSG1L	NOA1
RGN	CCDC126
GAD1	ELOVL5

Gene ontology term enrichment and kyoto encyclopedia of genes and genomes pathways analysis of target genes of differentially expressed miRNAs

In order to further explore the biological characteristics of these candidate miRNAs, we conducted GO annotation analysis and KEGG analysis of target genes of DEMs. The GO results showed that the related bioprocesses of up-regulated DEGs are considerably enriched in positive regulation of peptidyl-serine phosphorylation, negative regulation of transcription from RNA polymerase II promoter, focal adhesion assembly, behavioral fear response and small GTPase mediated signal transduction. As for cellular components, the upregulated DEGs are considerably enriched in cytoskeleton, cytoplasm, voltage-gated potassium channel complex, intracellular membrane-bounded organelle and axon. As for molecular functions, the up-regulated DEGs are considerably enriched in GDP binding, protein serine/threonine kinase activator activity, transcriptional activator activity, RNA polymerase II core promoter proximal region sequence-specific binding, protein binding and GTP binding (Fig. 3A).

Figure 3

Gene ontology function enrichment bubble diagram. (A) Top 5 of up-regulate differentially expressed genes. (B) Top 5 of down-regulate differentially expressed genes (P < .05).

For down-regulated DEGs, the bioprocesses of DEGs are considerably enriched in negative regulation of cell migration, neuron migration, circadian regulation of gene expression, small GTPase mediated signal transduction and regulation of exit from mitosis. As for cellular components analysis, the down-regulated DEGs are considerably enriched in cytoplasm, membrane, melanosome, axon and nucleolus. As for molecular functions, the up-regulated DEGs are considerably enriched in sequence-specific DNA binding, protein binding, poly(A) RNA binding, protein C-terminus binding and GTP binding (Fig. 3B).

In order to further exploring the pathways related with the DEGs, we performed the KEGG pathways analysis of all DEGs. KEGG analysis of overlapping genes showed they were mainly enriched in the Ras signaling pathway, Axon guidance, ErbB signaling pathway, Rap1 signaling pathway, T cell receptor signaling pathway, MAPK signaling pathway, RNA degradation, Circadian rhythm, Regulation of actin cytoskeleton, Sphingolipid signaling pathway, Vasopressin-regulated water reabsorption, Protein processing in endoplasmic reticulum, Histidine metabolism and RNA transport. (Fig. 4)

Figure 4

Kyoto encyclopedia of genes and genomes (KEGG) enrichment analysis results of all differentially expressed genes (P < .05).

Protein–protein interaction network construction and searching for hub genes

In order to further explore the functional interactions of the DEGs, we performed the PPI network construction through the STRING database and the combined score was set as >0.5. Then the PPI network built and visualization was performed through Cytoscape. There were 352 nodes and 725 edges. We screened the top 10 nodes with higher degrees by CytoHubba plug-in of Cytoscape software. The top 10 hub genes were SYN1, SYT4, RBFOX1, KCNC1, VAMP2, FBXO11, ASB9, SYT9, KLHL5, and KRAS (Fig. 5). We further constructed all the nodes and edges through Cytoscape, and the results showed in Figure 6. Then the predicted target miRNAs of hub genes were searched through TargetScan Website. The screened hub miRNAs showed in Table 2.

Figure 5

The top 10 hub genes. The color depth represents the absolute value of fold change.

Figure 6

Protein–protein interaction network. The larger the absolute value of fold change is, the darker the color is, and the smaller the absolute value of fold change is, the lighter the color is.

Table 2

A summary of miRNAs that regulate hub genes.

	Gene	Predicted miRNA		Gene	Predicted miRNA
1	KRAS	hsa-miR-532-5p hsa-miR-142-3p hsa-miR-155	2	RBFOX1	hsa-miR-92a hsa-miR-92b hsa-miR-25
3	VAMP2	hsa-miR-328 hsa-miR-185 hsa-miR-520a-3p	4	KCNC1	hsa-miR-137 hsa-let-7b hsa-let-7a
5	SYN1	hsa-miR-153	6	SYN4	hsa-miR-9-5p
7	KLHL5	hsa-miR-140-5p	8	ASB9	hsa-miR-217
9	FBXO11	hsa-miR-155	10	SYT9	hsa-miR-212 hsa-miR-133a hsa-miR-133b

Discussion

IgA nephropathy is the most common glomerulonephritis, which is mainly characteristic as immunoglobulin A (IgA) or IgA-based immunoglobulin deposition in the mesangial region of the glomerulus. However, the occurrence and progression of IgA nephropathy remains poorly understood at present. Studies have found that the activation of inflammation and immune factors in peripheral blood mononuclear cells is associated with IgA nephropathy. So exploring the key molecular of PBMCs in IgAN patients may provide the possibility of future molecularly targeted monitor disease progression and therapy, which may avoid the invasive renal biopsy and even prevent the occurrence of IgAN.

In this study, we identified 418 differentially expressed miRNAs (DEMs), including 103 up-regulated DEMs and 315 down regulated DEMs. SP1 and EGR1 were identified as upstream transcription factors that regulated the most of DEMs. Xu et al reported that SP1 could bind DEFA, which is associated with IgAN.[ EGR1 is a zinc-finger transcription factor expressed in various renal cells. Several study have shown that up-regulation of EGR1 was associated with renal fibrosis and inflammation.[

GO analysis showed that both up and down regulated DEGs are mainly enriched in biological processes of small GTPase mediated signal transduction. Kocher et al[ found that small GTPase is closely related with IgAN. So the PBMCs may affect the renal through small GTPase. KEGG enrichment analysis indicated MAPK pathway enriched by the DEGs. Milillo et al[ found that down regulation of the MAPK pathway represents a common mechanism leading to IgAN. Another study showed that IgA1 immune complex-mediated activation of the MAPK pathway in mesangial cells was associated with glomerular damage in IgA nephropathy.[ So the PBMCs shared the similar key pathways with IgAN progression.

Through PPI network analysis, we found top 10 hub gens of PBMCs from IgAN, including SYN1, SYT4, RBFOX1, KCNC1, VAMP2, FBXO11, ASB9, SYT9, KLHL5, and KRAS. KRAS gene is like a “switch” in vivo, which plays an important role in the signal transduction pathway of tumor cell growth and angiogenesis.[ Recent years, some studies showed that KRAS gene was associated with IgA.[ KRAS gene can active MAPK pathway and then participate in a variety of physiological processes.[ MAPK pathways are likely form the basis of IgAN progression.[ So the KRAS in PBMCs may be a prediction of potential targets for IgAN. We also screened hub miRNAs through TargetScan Website (Table 2). However, further experimental studies are needed to confirm these hub miRNAs.

There are still several limitations of this study. Firstly, the samples of PBMCs miRNA expression profiling is relative small. So larger sample sizes studies are still needed. Secondly, our result and conclusion are only based by the computational analysis. So future experimental studies are needed to confirm our results. Thirdly, analysis of stratify expression profiles according to KDIGO chronic kidney disease staging may give us more information about the miRNAs in the PBMCs of IgAN patients. However, there was no available datasets of stratify expression profiles in GEO. Fourthly, asymptomatic IgA deposition could be found in “normal” subjects, so this may bias the results.

In conclusion, our research screened out some key miRNAs and their targets that may be related with IgAN. These key genes and hub miRNA may help to monitor the happening and prognosis of IgAN. However, these predictions need further experimental studies to verify.

Author contributions

Conceptualization: Ling Liu.

Data curation: Ling Liu.

Formal analysis: YaZhen Yang.

Funding acquisition: Ling Liu, Dongrong Yu.

Investigation: Ling Liu.

Methodology: Ling Liu, YaZhen Yang.

Software: Ling Liu.

Writing – original draft: Ling Liu, YaZhen Yang.

Writing – review & editing: Dongrong Yu.