Identification of key target genes and pathway analysis in nonalcoholic fatty liver disease via integrated bioinformatics analysis

doi:10.21203/rs.3.rs-49190/v1

Download PDF

Research article

Identification of key target genes and pathway analysis in nonalcoholic fatty liver disease via integrated bioinformatics analysis

https://doi.org/10.21203/rs.3.rs-49190/v1

This work is licensed under a CC BY 4.0 License

Journal Publication

published 01 Jun, 2022

Read the published version in Balkan Journal of Medical Genetics →

Version 1

posted

You are reading this latest preprint version

Background

This study aimed to explore the mechanisms underlying nonalcoholic fatty liver disease (NAFLD) and develop new diagnostic biomarkers for nonalcoholic steatohepatitis (NASH).

Methods

The microarray dataset GES83452 was downloaded from the NCBI-GEO database, and the differentially expressed RNAs (DERs) were screened between the NAFLD and non-NAFLD samples of the baseline and 1-year follow-up time point group based on the Limma package. Subsequently, the lncRNA–miRNA–mRNA regulation network was constructed, and Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses on the regulated target genes in the ceRNA regulatory network were performed based on DAVID. Finally, PharmGKB database was used to search for the gene-related drug molecules, and the gene–drug connection network was constructed.

Results

A total of 561 DERs (268 downregulated and 293 upregulated) were screened in the baseline time point group, and 1163 DERs (522 downregulated and 641 upregulated) were screened in the 1-year follow-up time point group. A total of 74 lncRNA–miRNA pairs and 523 miRNA–mRNA pairs were obtained to construct lncRNA–miRNA–mRNA regulatory network. Subsequently, functional enrichment analysis revealed 28 GO and 9 KEGG pathways in the ceRNA regulatory network. LEPR and CXCL10 are involved in the Cytokine–cytokine receptor interaction (P = 1.86E-02), and the FOXO1 is involved in both the Insulin signaling pathway (P = 1.79E-02) and the pathways in cancer (P = 2.87E-02).

Conclusion

LEPR, CXCL10, and FOXO1 were the characteristic target genes for NAFLD.

Medical Genetics

NAFLD

ceRNA regulation network

PharmGKB

Nonalcoholic fatty liver disease (NAFLD) is the most common type of chronic liver disease with a prevalence rate of 25% worldwide[1]. Several lifestyle-related factors are associated with incident fatty liver such as alcohol intake, lower physical activity, smoking, and shift work. Bad lifestyles is the important cause of fatty liver, such as smoking, drinking, lack of physical activity, and shift work, etc. In addition, high triglycerides, type 2 diabetes mellitus, obesity and hypertension are associated with incident fatty liver. Therefore, lifestyle modification is strongly recommended to prevent fatty liver [2, 3]. It is difficult to detect in the earlier stages of the disease, and may further develop into advanced liver diseases, such as cirrhosis and hepatocellular carcinoma, which brings clinical challenges to the treatment of NAFLD [4]. In the literature, NAFLD patients with type 2 diabetes and obesity at the same time, the severity of NAFLD will be significantly affected, increasing the degree of deterioration of liver fibrosis and the possibility of further development of end-stage liver disease [5–7]. Likewise, studies have shown that when NAFLD patients suffer from cardiovascular diseases and dyslipidemia, they have a strong negative impact on the natural progression of NAFLD [8–10]. Nearly 40% of patients with NAFLD die of complications, as previously reported [1]. However, detailed mechanisms under which NAFLD develops remain largely unknown. Diet adjustment and weight loss can improve NAFLD, but it is difficult to maintain. Moreover, the theory of insulin resistance has been widely accepted clinically. Insulin sensitizers have a certain therapeutic effect, but they can cause adverse reactions such as increased body weight and its therapeutic target is too limited. Therefore, this study aimed to find new molecular targets to provide a theoretical basis for new and effective treatment methods of NAFLD.

Long non-coding RNA (lncRNA) is the main component of the human transcriptome. It plays an important role in regulating cell migration, proliferation, invasion and metastasis, and can be used as a diagnostic marker or therapeutic target for malignant tumors and other diseases. Competitive endogenous RNA (ceRNA) is a transcript with the same miRNA response element, which binds to miRNA to compete and regulate its target gene, thereby affecting the biological behavior of the disease. Studies have confirmed that the mutual regulation between lncRNA and miRNA and their downstream target genes plays an important role in the occurrence and development of diseases [11].

As we all know, the inflammatory component of nonalcoholic steatohepatitis (NASH) is more difficult to capture with ultrasound-assisted techniques. Although more and more technologies are applied in clinical practice such as quantitative and contrast-enhanced ultrasound, there are still many technical barriers to be broken through, and not all technologies have been successful in clinical and research practice [12]. Due to the limitations of liver biopsy, searching for non-invasive and reliable diagnostic biomarkers for NAFLD is the research hotspot and priority of current research. Bioinformatics has been widely used to explore biomarkers of different diseases, but NAFLD-related biomarkers need to be further explored to help early diagnosis and prognosis evaluation of NAFLD [13].

In this study, the human samples from the Gene Expression Omnibus (GEO) database was used to identify key genes related to NAFLD and non-NASH samples during baseline and 1-year follow-up time point, and to explore the underlying mechanism of NAFLD and develop new NAFLD diagnostic biomarkers. Then, the lncRNA–miRNA–mRNA network related to NAFLD was constructed by mapping the differentially expressed RNAs (DERs) into a global triple network by starBase and miRcode databases to identify which RNAs can be used as a sensitive and specific marker for NAFLD. Furthermore, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses were performed to explore the potential regulatory functions of RNAs. Finally, PharmGKB database was used to search and obtain gene-related drug molecules in the ceRNA regulatory network and then build a gene–drug connection network to screen out important gene molecules and KEGG signaling pathways involved in genes.

Microarray data and data preprocessing

GES83452 [14] in the NCBI-GEO (https://www.ncbi.nlm.nih.gov/) [15] database were downloaded on April 10, 2020, which included a total of 231 samples, including 159 patients at baseline (44 no NASH, 104 NASH, and 4 undefined) and 79 patients at 1-year follow-up (54 no NASH, 22 NASH, and 3 undefined) based on platforms GPL16686 [HuGene-2_0-st] Affymetrix Human Gene 2.0 ST Array [transcript (gene) version].

Screening significantly differentially expressed RNAs and functional enrichment analyses

The mRNA and lncRNA in GES83452 datasets were reannotated using the HUGO Gene Nomenclature Committee (http://www.genenames.org/) [16] based on information of Transcript ID, RefSeq ID, etc., which contain 4600 lncRNAs and 19195 protein coding gene. The Limma package (version 3.34.0, https://bioconductor.org/packages/release/bioc/html/limma.html) [17] in R was used to identify DERs between the NAFLD and non-NAFLD samples of the baseline and 1-year follow-up time point group. False discovery rate (FDR) < 0.05 and |log₂ fold change (FC)| > 0.5 were used as the cutoff criteria to define DERs, and the ggplot2packages in R was used to visualization the volcano plots. The heatmap was plotted using the pheatmap package (version 1.0.8, https://cran.r-project.org/package=pheatmap)[18] in R and was presented by two-way hierarchical clustering heatmaps [19] based on Euclidean distance [20]. P < .05 was considered statistically significant. The Venn software online (http://bioinformatics.psb.ugent.be/webtools/Venn/) was used to detect overlapping DERs among the baseline and 1-year follow-up time point groups. Then, GO and KEGG enrichment analyses were performed on intersection mRNAs that commonly contained DERs using the online tool DAVID (version 6.8, https://david.ncifcrf.gov/) [21, 22] P < .05 was considered as significant enrichment.

Construction of ceRNA network

The miRNAs related to NAFLD included in the Human MicroRNA Disease Database (HMDD) database (http://www.cuilab.cn/hmdd) were downloaded [23]. We constructed a ceRNA network based on NAFLD directly related to lncRNAs and miRNAs, as well as the miRNAs with significantly consistent expression. Firstly, we downloaded the connection relationship pairs of lncRNA-miRNA in the DIANA–LncBase (version 2, http://carolina.imis.athena-innovation.gr/diana_tools/web/index.php) [24]. The regulatory relationship between significantly DElncRNA and NAFLD-related differentially expressed miRNA (DEmiRNA) was retained, with retention connection score (miRNA target gene score (miTG–score): the target gene score of DEmiRNA; the higher the value, the greater the probability of targeting) higher than 0.6, thereby the lncRNA–miRNA connection network was constructed. Then, the starBase database (version 2.0, http://starbase.sysu.edu.cn/) [25] was used to predict target genes regulated by miRNA linked to lncRNA, and the comprehensive target gene prediction information from five databases (targetScan, picTar, RNA22, PITA and miRanda) was provided in the StarBase database. The target miRNA regulatory target gene relationship pair was selected in at least one of the databases, and the miRNA–mRNA pairs with the opposite significant differential expression direction was retained to construct the miRNA–mRNA connection network. Finally, a ceRNA regulation network composed of lncRNA–miRNA–mRNA was constructed by combining lncRNA–miRNA and miRNA–mRNA, and the ceRNA network was visualizatied by using the cytoscape (version 3.6.1, http://www.cytoscape.org/).

The screened target genes in the ceRNA regulatory network were submitted to DAVID 6.8 online tool to perform functional annotation based on GO biological processes and KEGG pathway enrichment analysis, the P value < 0.05 as the significance threshold.

Construction of drug–gene regulation network

The pharmacogenetics and pharmacogenomics knowledge base (PharmGKB) (https://www.pharmgkb.org/) [26] collected the most complete genotype and phenotype information related to the drug genome and was classified systematically, which contained 27,007 genes related to 3579 drugs and 3410 diseases. In this study, the PharmGKB database was used to search for and obtain the gene-related drug molecules in the regulated ceRNA network; then the gene–drug connection network was constructed, the important gene molecules were screened out, and the KEGG signaling pathway of those genes participated in in-depth analysis.

Data preprocessing and DERs screening

A total of 9698 mRNAs and 1116 lncRNAs were detected after data preprocessing, and a total of 561 DERs (48 lncRNA and 513 mRNA; 268 downregulated and 293 upregulated) were screened in the baseline time point group, and 1163 DERs (114 lncRNA and 1049 mRNA; 522 downregulated and 641 upregulated) were screened in the 1-year follow-up time point group, with FDR < 0.05 and |log₂FC| > 0.5 as the cutoff criteria. We identified all DERs shown in the volcano map according to the value of |log₂FC| and displayed DERs on a heatmap (Fig. 1A, B). The expression values of DERs were 2-way hierarchically clustered, and the color contrast indicated that there was a significant difference in expression levels between the NAFLD and non-NAFLD samples (Fig. 1C, D). Subsequently, a total of 220 overlapping DERs were identified between the baseline and 1-year follow-up time points, which were used by the Venn diagram software (Fig. 2).

In addition, the functional enrichment analysis of the overlapping DERs based on online DAVID analyses revealed 22 significantly related GO biological processes and 9 KEGG pathways, with P < .05 as the cutoff criteria (Table 1). We found that chemotaxis (GO, 0006935; P = 3.110E-04), unsaturated fatty acid biosynthetic process (GO, 0006636; P = 4.770E-04), and cell-cell signaling (GO, 0007267; P = 1.513E-03) were the three most significant pathways in GO biological processes. Meanwhile, fatty acid metabolism (hsa01212, P = 2.300E-04), PPAR signaling pathway (hsa03320, P = 1.090E-03), and Toll-like receptor signaling pathway (hsa04620, P = 1.479E-03) were the three most significant pathways in KEGG signaling pathways.

Table 1

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis for the differentially expressed RNAs (DERs).
Category	Term	Count	P-Value	Genes
Biology Process	GO:0006935 ~ chemotaxis	8	3.110E-04	MAP2K1, CXCL9, AMOT, etc
	GO:0006636 ~ unsaturated fatty acid biosynthetic process	4	4.770E-04	FADS1, SCD, FADS2, etc
	GO:0007267 ~ cell-cell signaling	10	1.513E-03	SH2D1A, ADM, FADS1, etc
	GO:0007568 ~ aging	8	1.848E-03	TYMS, ADM, APOD, etc
	GO:0007166 ~ cell surface receptor signaling pathway	10	2.538E-03	MARCO, PRLR, TSPAN3, etc
	GO:0055114 ~ oxidation-reduction process	15	4.070E-03	KDM6A, HSD17B2, PYROXD2, etc
	GO:0016337 ~ single organismal cell-cell adhesion	6	4.312E-03	MPZL2, PKHD1, FAT1, etc
	GO:0006629 ~ lipid metabolic process	7	6.425E-03	APOD, PLIN1, APOF,etc
	GO:0070098 ~ chemokine-mediated signaling pathway	5	6.698E-03	TFF2, CXCL9, ACKR3,etc
	GO:0032496 ~ response to lipopolysaccharide	7	7.895E-03	ADM, DUSP10, CXCL9, etc
	GO:0010508 ~ positive regulation of autophagy	4	8.546E-03	RNF152, FOXO1, TRIM22, etc
	GO:0035338 ~ long-chain fatty-acyl-CoA biosynthetic process	4	9.779E-03	SCD, FASN, ELOVL6, ACSL5
	GO:0006915 ~ apoptotic process	13	1.716E-02	PEG10, RASSF6, LITAF, etc
	GO:0002250 ~ adaptive immune response	6	2.033E-02	SH2D1A, EOMES, CD1C, etc
	GO:0006959 ~ humoral immune response	4	2.223E-02	SH2D1A, BST1, LTF, CD28
	GO:0006968 ~ cellular defense response	4	2.767E-02	SH2D1A, CXCL9, LBP, etc
	GO:0060326 ~ cell chemotaxis	4	3.124E-02	CXCL9, CCL5, DOCK4, etc
	GO:0045766 ~ positive regulation of angiogenesis	5	3.346E-02	ADM, LRG1, RHOB, etc
	GO:0032868 ~ response to insulin	4	3.374E-02	ADM, INSIG2, FADS1, PCK1
	GO:0008203 ~ cholesterol metabolic process	4	3.504E-02	INSIG2, APOF, LEPR, ERLIN1
	GO:0001889 ~ liver development	4	4.331E-02	COBL, ONECUT1, DBP, RPGRIP1L
	GO:0006954 ~ inflammatory response	9	4.788E-02	LXN, AOX1, LYZ, etc
KEGG Pathway	hsa01212:Fatty acid metabolism	6	2.300E-04	FADS1, SCD, FASN, etc
	hsa03320:PPAR signaling pathway	6	1.090E-03	PLIN1, SCD, FADS2, etc
	hsa04620:Toll-like receptor signaling pathway	7	1.479E-03	CTSK, MAP2K1, CXCL9, etc
	hsa01040:Biosynthesis of unsaturated fatty acids	4	2.350E-03	FADS1, SCD, FADS2, ELOVL6
	hsa04640:Hematopoietic cell lineage	6	3.475E-03	CR1, CD2, CD1C, etc
	hsa04152:AMPK signaling pathway	6	1.465E-02	LEPR, SCD, FASN, etc
	hsa04668:TNF signaling pathway	5	3.695E-02	MAP2K1, IL15, CREB3L3, etc
	hsa00760:Nicotinate and nicotinamide metabolism	3	4.529E-02	BST1, ENPP3, AOX1
	hsa04060:Cytokine-cytokine receptor interaction	7	4.659E-02	PRLR, LEPR, CXCL9, etc

Construction of ceRNA regulation network

A total of 77 miRNAs directly related to NAFLD were downloaded from the HMDD database. After the lncRNA and miRNA connection relationship pairs were downloaded and the regulation relationship between significantly DElncRNA and NAFLD-related DEmiRNAs were selected, a total of 74 connection pairs were retained to construct an lncRNA–miRNA connection network with a connection coefficient higher than 0.6. In addition, after the target genes of the miRNA were screened, then comparing the regulated target genes with the significant DEmRNAs in target modules and retaining the opposite relationship pairs of the expressions of significant differential direction. The miRNA–mRNA regulation network was constructed by using a total of 523 pairs of regulation relationships. Finally, as shown Fig. 3, a ceRNA regulation network was constructed.

A total of 28 GO biological processes and 9 KEGG pathways of the mRNAs in the ceRNA regulatory network were obtained, with P < .05 as the significance threshold (Table 2). We found that lipid biosynthetic process (GO, 0008610; P = 1.15E-06), steroid metabolic process (GO, 0008202; P = 1.26E-05), and steroid biosynthetic process (GO, 0006694; P = 9.34E-05) were the three most significant pathways in GO biological processes. Meanwhile, biosynthesis of unsaturated fatty acids (hsa01040, P = 2.39E-05), terpenoid backbone biosynthesis (hsa00900, P = 1.090E-03), and heparan sulfate biosynthesis (hsa00534, P = 1.38E-02) were the three most significant pathways in KEGG signaling pathways.

Table 2

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis for the mRNA in the ceRNA regulatory network.
Category	Term	Count	P-Value	Genes
Biology Process	GO:0008610 ~ lipid biosynthetic process	12	1.15E-06	PRLR, ADM, ISPD, etc
	GO:0008202 ~ steroid metabolic process	9	1.26E-05	PRLR, INSIG2, ADM, etc
	GO:0006694 ~ steroid biosynthetic process	6	9.34E-05	PRLR, ADM, FDPS, etc
	GO:0006633 ~ fatty acid biosynthetic process	5	8.64E-04	FADS1, SCD, FASN, etc
	GO:0046649 ~ lymphocyte activation	6	4.38E-03	PRLR, EOMES, IL15, etc
	GO:0016053 ~ organic acid biosynthetic process	5	9.79E-03	FADS1, SCD, FASN, etc
	GO:0046394 ~ carboxylic acid biosynthetic process	5	9.79E-03	FADS1, SCD, FASN, etc
	GO:0045321 ~ leukocyte activation	6	9.84E-03	PRLR, EOMES, IL15, etc
	GO:0006355 ~ regulation of transcription, DNA-dependent	18	1.19E-02	CAMTA2, MAP2K1, LITAF, etc
	GO:0030334 ~ regulation of cell migration	5	1.31E-02	MAP2K1, CLIC4, AMOT, etc
	GO:0051252 ~ regulation of RNA metabolic process	18	1.47E-02	CAMTA2, MAP2K1, etc
	GO:0007267 ~ cell-cell signaling	9	1.52E-02	ADM, FADS1, FAT1, etc
	GO:0051094 ~ positive regulation of developmental process	6	1.71E-02	MAP2K1, EOMES, AMOT, etc
	GO:0001775 ~ cell activation	6	1.93E-02	PRLR, EOMES, IL15, etc
	GO:0040012 ~ regulation of locomotion	5	2.00E-02	MAP2K1, CLIC4, AMOT, etc
	GO:0051270 ~ regulation of cell motion	5	2.04E-02	MAP2K1, CLIC4, AMOT, etc
	GO:0009719 ~ response to endogenous stimulus	7	2.19E-02	MAP2K1, ADM, FADS1, etc
	GO:0006631 ~ fatty acid metabolic process	5	2.21E-02	FADS1, SCD, FASN, etc
	GO:0045860 ~ positive regulation of protein kinase activity	5	3.23E-02	PRLR, MAP2K1, CD24, etc
	GO:0060429 ~ epithelium development	5	3.42E-02	STX2, MAP2K1, RPGRIP1L, etc
	GO:0033674 ~ positive regulation of kinase activity	5	3.61E-02	PRLR, MAP2K1, CD24, etc
	GO:0051347 ~ positive regulation of transferase activity	5	4.06E-02	PRLR, MAP2K1, CD24, etc
	GO:0001568 ~ blood vessel development	5	4.33E-02	LEPR, FOXO1, AMOT, etc
	GO:0051254 ~ positive regulation of RNA metabolic process	7	4.51E-02	CAMTA2, MAP2K1, CSRNP1, etc
	GO:0001944 ~ vasculature development	5	4.66E-02	LEPR, FOXO1, AMOT, etc
	GO:0009725 ~ response to hormone stimulus	6	4.80E-02	MAP2K1, ADM, FADS1, etc
	GO:0030030 ~ cell projection organization	6	4.84E-02	STX2, MAP2K1, ADM, etc
	GO:0007243 ~ protein kinase cascade	6	4.94E-02	PRLR, MAP2K1, DUSP10, etc
KEGG Pathway	hsa01040:Biosynthesis of unsaturated fatty acids	4	2.39E-05	FADS1, SCD, FADS2, ELOVL6
	hsa00900:Terpenoid backbone biosynthesis	2	8.23E-03	FDPS, MVK
	hsa00534:Heparan sulfate biosynthesis	2	1.38E-02	EXT1, HS3ST3B1
	hsa04910:Insulin signaling pathway	3	1.79E-02	MAP2K1, FASN, FOXO1
	hsa04060:Cytokine-cytokine receptor interaction	4	1.86E-02	PRLR, LEPR, IL15, CXCL10
	hsa04010:MAPK signaling pathway	4	1.93E-02	MAP2K1, DUSP10, MAP3K13, FGF2
	hsa04630:Jak-STAT signaling pathway	3	2.21E-02	PRLR, LEPR, IL15
	hsa05200:Pathways in cancer	4	2.87E-02	MAP2K1, FZD1, FOXO1, FGF2
	hsa03320:PPAR signaling pathway	2	3.28E-02	SCD, FADS2

Construction of drug regulation gene network

The gene-related drug molecule connection pairs were downloaded from the PharmGKB database. A total of 154 connection pairs were obtained by selecting the parts related to the genes in the ceRNA regulatory network, and a gene–drug connection network was constructed (Fig. 4). Compared with the pathway in which the RNAs significantly participated in the ceRNA regulatory network constructed in the previous step, the leptin receptor (LEPR) and CXCL10 are involved in the Cytokine–cytokine receptor interaction (P = 1.86E-02), and FOXO1 is involved in both the Insulin signaling pathway (P = 1.79E-02) and the pathways in cancer (P = 2.87E-02).

The characteristics of NAFLD are necrotizing inflammation and lipid accumulation in the liver, and the continuous improvement of living standards lead to overnutrition. In addition, bad living habits lead to the incidence of NAFLD on a global scale [27]. The specific mechanism of the transition from benign steatosis to steatohepatitis in NAFLD is not fully understood, and there are currently no pharmacological options for the treatment of NAFLD. Therefore, the treatment of NASH mainly depends on changing lifestyles, such as strengthening exercise, reducing weight, and light diet [28]. Although current studies have shown that weight loss improves the histological characteristics of NAFLD, but most patients have not achieved the goal of curing NAFLD. However, some potentially valuable molecules are currently being clinically evaluated [29]. For example, PNPLA3 [30], TM6SF2 [31], MBOAT7 [32], and HSD17B13 [33] that predispose an individual to the spectrum of NAFLD-related disease have been found to play a role in macrophage phagocytosis, immune response, oxidative stress, and inflammation, insulin signaling, and lipid metabolism in NAFLD susceptibility and progression [34]. But there is no unmet clinical need for drug discovery and development for patients with NAFLD. Increased levels of toxic lipids (free fatty acids or free cholesterol) can lead to liver cell damage and trigger inflammation is the pathogenesis of NAFLD is currently understood. In addition, oxidative stress, pro-inflammatory chemokines and cytokines have been proven to be lead to liver inflammation, which in turn leads to the damage and fibrosis of liver. Therefore, the identification of pro-inflammatory cytokines related to lipotoxicity may improve our understanding of the pathogenesis of NAFLD, thereby helping to develop new pharmacological methods.

In this study, a total of 220 overlapping DERs were identified between the baseline and 1-year follow-up time points. In addition, functional enrichment analysis of overlapping DERs based on online DAVID analyses revealed 22 significantly related GO biological processes and 9 KEGG pathways, with P < .05 as the cutoff criteria. We found that chemotaxis (P = 3.110E-04), unsaturated fatty acid biosynthetic process (P = 4.770E-04), and cell-cell signaling (P = 1.513E-03) were the three most significant pathways in GO biological processes. Meanwhile, fatty acid metabolism (P = 2.300E-04), PPAR signaling pathway (P = 1.090E-03), and Toll-like receptor signaling pathway (P = 1.479E-03) were the three most significant pathways in KEGG signaling pathways. Then a ceRNA regulatory network was constructed. The GO and pathway enrichment analysis indicated that mRNAs of ceRNA regulatory network were involved in various important biological functions and metabolic pathways associated with NAFLD, including lipid biosynthetic process, steroid metabolic process, steroid biosynthetic process, biosynthesis of unsaturated fatty acids, terpenoid backbone biosynthesis, heparan sulfate biosynthesis, Cytokine–cytokine receptor interaction, Insulin signaling pathway, and the pathways in cancer. To further understand the functional mechanism of the ceRNA network, the drug regulation gene network was constructed, which included 154 gene–drug connection pairs. Subsequently, LEPR, CXCL10, and FOXO1 were investigated using the PharmGKB database. It was revealed that the therapeutic effect of antipsychotics, atorvastatin, valproic acid, risperidone, clozapine, olanzapine, simvastatin, and quetiapine were produced may by targeting to LEPR through the Cytokine–cytokine receptor interaction pathway. The therapeutic effect of Peginterferon alfa-2a and peginterferon alfa-2b were produced may by targeting the CXCL10 through the Cytokine–cytokine receptor interaction pathway. The therapeutic effect of Epirubicin, cyclophosphamide, and fluorouracil were produced may by targeting the FOXO1 through the Insulin signaling pathway or the pathways in cancer.

LEPR is responsible for encoding the leptin receptor that binds to leptin in target tissues. Due to its role in regulating lipid metabolism and insulin resistance, it is considered to be a candidate gene for NAFLD and coronary atherosclerosis [35]. Simultaneously, An et al. [36] found that LEPR Q223R polymorphism may lead to a significant risk of NAFLD and coronary atherosclerosis, which is consistent with the results of this study. The CXC motif chemokine ligand 10 (CXCL10) is a particularly important pro-inflammatory cytokine related to lipotoxicity, which can recruit inflammatory cells to the site of tissue injury [37, 38]. Studies have shown that CXCL10 is upregulated in NAFLD patients [39], and revealed that CXCL10 may be a key molecule that contributes to the transition from benign steatosis to steatohepatitis, promoting liver cell damage and inflammation [40].

Our study revealed that peginterferon alfa-2a and peginterferon alfa-2b can downregulate the expression of CXCL10, suggesting a potential role of CXCL10 in the development of intrahepatic inflammation through the Cytokine–cytokine receptor interaction pathway, and demonstrated that CXCL10 is an independent risk factor for patients with NAFLD. FOXO1 is an important transcriptional effector. It is widely expressed in various types of tissues and plays an important role in the signaling pathway of insulin and insulin-like growth factor 1 [41]. In addition, the expression levels of most genes related to adipocyte differentiation are affected by the coordination of FOXO1 [42]. Yue Li et al. [43] conducted a comprehensive analysis of the relevant information about the activity of FOXO1 in lipid metabolism, and found that FOXO1 has a significant inhibitory effect on the production of fibrotic effector cells, and pointed out that FOXO1 has the potential to become a target for the treatment of NAFLD, but the related mechanism needs to be further verified by experiments [44]. L. Valenti et al.[45] found that FOXO1 may affect the susceptibility of NAFLD, and regulating the level of FOXO1 mRNA to regulate the relevant cytokines in insulin signaling to promote the progression of liver injury. This study found that epirubicin, cyclophosphamide, and fluorouracil can downregulate the expression of FOXO1, suggesting that these drugs may produce therapeutic effect by targeting FOXO1 through the insulin signaling pathway or the pathways in cancer. However, further study is necessary to validate this hypothesis.

In conclusion, this study constructed and analyzed an ceRNA network, which may provide some evidence to future studies focusing on the molecular mechanisms of NALFD. LEPR, CXCL10, and FOXO1 may function as ceRNAs to serve critical roles in NALFD. In addition, antipsychotics, atorvastatin, valproic acid, risperidone, clozapine, olanzapine, simvastatin, quetiapine, peginterferon alfa-2a, peginterferon alfa-2b, epirubicin, cyclophosphamide, and fluorouracil may produce therapeutic effect for patients with NALFD.

Abbreviations

Nonalcoholic fatty liver disease (NAFLD)

Long noncoding RNAs (lncRNAs)

Competing endogenous RNAs (ceRNAs)

Nonalcoholic steatohepatitis (NASH)

Gene Expression Omnibus (GEO)

Differentially expressed RNAs (DERs)

Gene Ontology (GO)

Kyoto Encyclopedia of Genes and Genomes (KEGG)

False discovery rate (FDR)

Human MicroRNA Disease Database (HMDD)

Leptin receptor (LEPR)

Declarations

Ethics approval and consent to participate:

Not applicable

Consent for publication:

Not applicable

Availability of data and material

The raw data were collected and analyzed by the Authors, and are not ready to share their data because the data have not been published.

Competing interests:

The authors declare no conflict of interest.

Funding:

none

Authors' contributions

XC and LZ participated in the design of this study, and they both performed the statistical analysis. YYW and REL carried out the study and collected important background information. MY and LG drafted the manuscript. All authors read and approved the final manuscript.

Acknowledgements:

none

Erman H, Beydogan E, Cetin SI et al. (2020) Endocan: A Biomarker for Hepatosteatosis in Patients with Metabolic Syndrome. Mediators Inflamm 2020:3534042. doi:10.1155/2020/3534042
Aizawa M, Inagaki S, Moriyama M et al. (2019) Modeling the natural history of fatty liver using lifestyle-related risk factors: Effects of body mass index (BMI) on the life-course of fatty liver. PLoS One 14 (10):e0223683. doi:10.1371/journal.pone.0223683
Zelber-Sagi S, Lotan R, Shlomai A et al. (2012) Predictors for incidence and remission of NAFLD in the general population during a seven-year prospective follow-up. Journal of hepatology 56 (5):1145-1151. doi:10.1016/j.jhep.2011.12.011
Sookoian S, Pirola CJ (2019) Review article: shared disease mechanisms between non-alcoholic fatty liver disease and metabolic syndrome - translating knowledge from systems biology to the bedside. Aliment Pharmacol Ther 49 (5):516-527. doi:10.1111/apt.15163
Friedman SL, Neuschwander-Tetri BA, Mary R et al. (2018) Mechanisms of NAFLD development and therapeutic strategies.
Rinella MEJJtJotAMA (2015) Nonalcoholic Fatty Liver Disease: A Systematic Review. 313 (22):2263-2273
Noureddin M, Rinella MEJCiLD (2015) Nonalcoholic Fatty Liver Disease, Diabetes, Obesity, and Hepatocellular Carcinoma. 19 (2):361-379
Fabbrini E, Magkos F, Mohammed BS et al. (2009) Intrahepatic fat, not visceral fat, is linked with metabolic complications of obesity. 106 (36):15430-15435
Konerman MA, Jones JC, Harrison SAJJoH (2017) Pharmacotherapy for NASH: Current and Emerging.S0168827817323942
Pirola CJJH (2010) Epigenetic regulation of insulin resistance in nonalcoholic fatty liver disease: Impact of liver methylation of the peroxisome proliferator–activated receptor γ coactivator 1α promoter. 52
Liu NK, Xu X-M MicroRNA in central nervous system trauma and degenerative disorders. 43 (10):571-580
Piazzolla VA, Mangia A (2020) Noninvasive Diagnosis of NAFLD and NASH. Cells 9 (4). doi:10.3390/cells9041005
Zou B, Yeo YH, Nguyen VH et al. (2020) Prevalence, characteristics and mortality outcomes of obese, nonobese and lean NAFLD in the United States, 1999-2016. J Intern Med. doi:10.1111/joim.13069
Philippe L, Fanny L, Eric B et al. (2017) Interspecies NASH disease activity whole-genome profiling identifies a fibrogenic role of PPARα-regulated dermatopontin. 2 (13):e92264-
Barrett T, Troup DB, Wilhite SE et al. (2007) NCBI GEO: mining tens of millions of expression profiles--database and tools update. 35 (Database):D760-D765
Bryony, Braschi, Paul et al. Genenames.org: the HGNC and VGNC resources in 2019.
Ritchie ME, Belinda P, Di W et al. (2015) \n limma powers differential expression analyses for RNA-sequencing and microarray studies. (7):7
Wang L, Cao C, Ma Q et al. RNA-seq analyses of multiple meristems of soybean: novel and alternative transcripts, evolutionary and functional implications. 14 (1)
Bonow RO Braunwald`s Heart Disease: A Textbook of Cardiovascular Medicine, Single Volume, 9th Edition.
Szekely GJ, Classification MLRJJo Hierarchical Clustering via Joint Between-Within Distances: Extending Ward\"s Minimum Variance Method. 22 (2):151-183
Da WH, Sherman BT, Protocols RALJN (2008) Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. 4 (1):44-57
Huang DW, Sherman BT, Lempicki RAJNAR (2008) Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists. (1):1
Zhou, Huang, Jiangcheng et al. HMDD v3.0: a database for experimentally supported human microRNA-disease associations.
Paraskevopoulou MD, Vlachos IS, Dimitra K et al. (2015) DIANA-LncBase v2: indexing microRNA targets on non-coding transcripts. (D1):D1
Li JH, Liu S, Zhou H et al. (2013) Starbase v2.0: Decoding miRNA-ceRNA, miRNA-ncRNA and protein-RNA interaction networks from large-scale clip-seq data. 42 (D1)
Thorn CF, Klein TE, Altman RBJMiMB (2005) PharmGKB: The Pharmacogenomics Knowledgebase. 311:179-191
Chen HJ, Liu J (2018) Actein ameliorates hepatic steatosis and fibrosis in high fat diet-induced NAFLD by regulation of insulin and leptin resistant. Biomed Pharmacother 97:1386-1396. doi:10.1016/j.biopha.2017.09.093
Pan X, Zheng M, Zou T et al. (2018) The LEPR K109R and Q223R Might Contribute to the Risk of NAFLD: A Meta-Analysis. Curr Mol Med 18 (2):91-99. doi:10.2174/1566524018666180705110412
Zhang Y, Xiang D, Hu X et al. (2019) Identification and study of differentially expressed miRNAs in aged NAFLD rats based on high-throughput sequencing. Ann Hepatol. doi:10.1016/j.aohep.2019.12.003
Romeo S, Kozlitina J, Xing C et al. (2008) Genetic variation in PNPLA3 confers susceptibility to nonalcoholic fatty liver disease. 40 (12):1461-1465
Kozlitina J, Smagris E, Stender S et al. (2014) Exome-wide association study identifies a TM6SF2 variant that confers susceptibility to nonalcoholic fatty liver disease. 46 (4):352-356
Mancina RM, Dongiovanni P, Petta S et al. (2016) The MBOAT7-TMC4 Variant rs641738 Increases Risk of Nonalcoholic Fatty Liver Disease in Individuals of European Descent. 150 (5):1219-1230.e1216
Interferon-λ rs12979860 genotype and liver fibrosis in viral and non-viral chronic liver disease %J Nature Communications (2015). 6:6422
Petta S, Valenti L, Marra F et al. (2016) MERTK rs4374383 polymorphism affects the severity of fibrosis in non-alcoholic fatty liver disease. 64 (3):682-690
Aller R, Luis DAD, Izaola O et al. (2012) Lys656Asn polymorphism of leptin receptor, leptin levels and insulin resistance in patients with non alcoholic fatty liver disease. 16 (3):335-341
An BQ, Lu LL, Yuan C et al. (2016) Leptin Receptor Gene Polymorphisms and the Risk of Non-Alcoholic Fatty Liver Disease and Coronary Atherosclerosis in the Chinese Han Population. Hepat Mon 16 (4):e35055. doi:10.5812/hepatmon.35055
Neville LF, Mathiak G, ., Bagasra O, . %J Cytokine et al. (1997) The immunobiology of interferon-gamma inducible protein 10 kD (IP-10): a novel, pleiotropic member of the C-X-C chemokine superfamily. 8 (3):207
Luster AD, Unkeless JC, Ravetch JVJN (1985) γ-Interferon transcriptionally regulates an early-response gene containing homology to platelet proteins. 315 (6021):672-676
Bertola A, Bonnafous S, Anty R et al. (2010) Hepatic Expression Patterns of Inflammatory and Immune Response Genes Associated with Obesity and NASH in Morbidly Obese Patients. 5 (10):e13577
Zhang X, Shen J, Man K et al. (2014) CXCL10 plays a key role as an inflammatory mediator and a non-invasive biomarker of non-alcoholic steatohepatitis. Journal of hepatology 61 (6):1365-1375. doi:10.1016/j.jhep.2014.07.006
Dong XC, Copps KD, Guo S et al. (2008) Inactivation of Hepatic Foxo1 by Insulin Signaling Is Required for Adaptive Nutrient Homeostasis and Endocrine Growth Regulation. 8 (1):0-76
Munekata K, Cellular KSJV, Animal DB (2009) Forkhead Transcription Factor Foxo1 Is Essential for Adipocyte Differentiation. 45 (10):642-651
Li Y, Ma Z, Jiang S et al. (2017) A global perspective on FOXO1 in lipid metabolism and lipid-related diseases. Prog Lipid Res 66:42-49. doi:10.1016/j.plipres.2017.04.002
Xin Z, Ma Z, Hu W et al. (2018) FOXO1/3: Potential suppressors of fibrosis. Ageing Res Rev 41:42-52. doi:10.1016/j.arr.2017.11.002
Valenti L, Dongiovanni P, Rametta R et al. (2009) FOXO1 genotype influences the susceptibility to and severity of NAFLD by modulating FOXO1 expression. 41 (3):A2–A3

Download PDF

Journal Publication

published 01 Jun, 2022

Read the published version in Balkan Journal of Medical Genetics →

Version 1

posted

You are reading this latest preprint version

Identification of key target genes and pathway analysis in nonalcoholic fatty liver disease via integrated bioinformatics analysis

Status:

Journal Publication

Version 1

Abstract

Background

Methods

Results

Conclusion

Figures

Background

Materials And Methods

Microarray data and data preprocessing

Screening significantly differentially expressed RNAs and functional enrichment analyses

Construction of ceRNA network

Results

Data preprocessing and DERs screening

Construction of ceRNA regulation network

Construction of drug regulation gene network

Discussion

Conclusion

Declarations

Ethics approval and consent to participate:

Funding:

Acknowledgements:

References

Status:

Journal Publication

Version 1