基于生物信息学分析肝硬化与miR-125a-5p和miR-125b-5p的相关研究

Bioinformatics-Based Analysis of Liver Cirrhosis in Association with miR-125a-5p and miR-125b-5p

摘要: 目的：基于生物信息学分析肝硬化与miR-125a-5p和miR-125b-5p的功能富集通路以及关键基因，以期为肝硬化诊断与治疗提供新思路。方法：通过miRDB、miRWalk和TargetScan数据库获得miR-125a-5p和miR-125b-5p的靶基因，从GEO数据库获得肝硬化相关数据集GSE14323和GSE25097，筛选正常组织样本和肝硬化组织样本差异表达基因，将GSE25097和GSE14323所得差异表达基因取交集，获得肝硬化数据集的共表达基因，将miR-125a-5p和miR-125b-5p的靶基因与肝硬化数据集GSE14323和GSE25097的共表达基因取交集得到核心基因，对核心基因进行GO和KEGG富集以及构建PPI网络互作图。结果：GO富集分析中主要在糖氧代谢过程和线粒体外膜的活性中富集，KEGG富集分析中主要在脂肪酸生物合成和类固醇生物代谢等中富集，通过构建PPI网络互作图，得到miR-125a-5p和miR-125b-5p与肝硬化相关的4个关键基因NPL、H6PD、KIAA0319L和RFX5。结论：得到的相关通路以及4个关键基因NPL、H6PD、KIAA0319L和RFX5以期为肝硬化发病机制以及肝硬化发生发展过程中相关分子靶点和早期筛查的手段提供新思路，从而指导疾病诊疗计划和管理。

Abstract: Objective: To analyze the functional enrichment pathway and key genes of miR-125a-5p and miR-125b-5p in cirrhotic patients based on bioinformatics, so as to provide new ideas for the diagnosis and treatment of cirrhosis. Methods: The target genes of miR-125a-5p and miR-125b-5p were obtained from miRDB, miRWalk and TargetScan databases, the differentially expressed genes were screened from normal and cirrhotic tissues, and the differentially expressed genes from GSE25097 and GSE14323 were crossed to obtain the co-expressed genes from the cirrhotic data set, the target genes of miR-125a-5p and miR-125b-5p were intersected with the co-expressed genes of liver cirrhosis data sets GSE14323 and GSE25097 to obtain core genes, which were enriched by GO and KEGG, and a PPI network interaction map was constructed. Results: GO was enriched in the processes of glucose and oxygen metabolism and the activities of Outer mitochondrial membrane. KEGG was enriched in the biosynthesis of fatty acids and steroids, four key genes of miR-125a-5p and miR-125b-5p related to liver cirrhosis, NPL, H6PD, KIAA0319L and RFX5, were obtained by constructing PPI network interaction maps. Conclusion: The related pathways and four key genes, NPL, H6PD, KIAA0319L and RFX5, are expected to provide new ideas for the pathogenesis of liver cirrhosis, molecular targets and early screening of liver cirrhosis, to guide the disease diagnosis and treatment planning and management.

1. 引言

肝细胞脂质过载、胆固醇代谢失调以及慢性病毒感染等因素引起肝细胞死亡，在免疫细胞作用下，处于静息状态的肝星状细胞(hepatic stellate cell, HSC)活化，持续活化的HSCs破坏细胞外基质(extracellular matrix, ECM)沉积与溶解之间的平衡，引发渐进性肝硬化 [1] 。肝硬化多是由病毒性肝炎、自身免疫性肝病(autoimmune liver disease, AIH)、酒精性肝病(alcoholic liver disease, ALD)、非酒精性脂肪性肝病(Nonalcoholic fatty liver disease, NAFLD)等慢性肝脏疾病发展而来 [2] 。2019年肝硬化死亡人数占全球死亡人数的2.4%，成为全世界的主要死因 [3] 。代偿期肝硬化患者大多数没有症状或者症状不明显，临床上肝硬化患者多处于失代偿期。约4%~12%的肝硬化患者每年至少发生一次肝硬化失代偿事件，最常见的肝硬化失代偿事件是腹水、食管胃底静脉曲张出血和肝性脑病(hepatic encephalopathy, HE) [4] 。

微小RNA (miRNA)是长度为20~25个核苷酸非编码的小RNA分子，其通过调节细胞生长和组织分化等过程，在生长发育和疾病等发生发展过程中发挥重要作用。miR-125家族由miR-125a和miR-125b组成，miR-125a-5p与miR-125b-5p的成熟体核心序列具有高度一致性。肝硬化患者药物治疗前后分别检测miR-125a-5p的表达水平，药物治疗后miR-125a-5p表达水平降低，且随着治疗时间的延长，其表达水平持续下降 [5] 。miR-125a-5p在肝细胞癌(hepatocellular carcinoma, HCC)中表达量降低，其通过靶向作用于TER2P途径，抑制蛋白表达、HCC细胞的侵入能力以及HCC的增殖和迁移 [6] 。miR-125a-5p在乳腺癌、胃癌、宫颈癌、前列腺癌和膀胱癌中参与调节细胞增殖、凋亡、迁移、浸润和转移等过程 [7] 。一项关于姜黄素的研究中，miR-125b-5p过表达通过降低正常肝细胞抗氧化能力，降低过氧化氢酶(CAT)和超氧化物歧化酶(SOD)的活性，起到促进氧化的作用 [8] 。在构建的肝硬化小鼠和细胞模型中，miR-125b-5p的表达水平降低，其通过调控ITGA8过表达，激活RhoA信号通路，促进肝纤维化的发生 [9] 。肝纤维化是肝硬化发展的重要阶段，HSC是肝纤维化发生发展过程中起重要作用的细胞，连翘提取物(PHI)通过抑制JAK1/JAK2-STAT1和Notch1信号通路抑制M1巨噬细胞极化，以及巨噬细胞来源的外泌体中miR-125b-5p靶基因Stard13表达，抑制NF-KB激活，下调HSC中纤维化因子的表达，减少HSC的活化，抑制肝纤维的进展 [10] 。在慢性乙型肝炎、乙型肝炎相关肝硬化的乙型肝炎阳性HCC中miR-125b-5p的表达水平分别上调2.85倍、2.46倍和1.89倍 [11] 。miR-125b-5p在肺癌、乳腺癌、膀胱癌和胃癌中的表达量增高；在乳腺癌中，miR-125b-5p抑制乳腺癌细胞的增殖、迁移和入侵；在膀胱癌中，通过抑制PI3K/AKT通路，降低膀胱癌中miR-125b-5p的表达，抑制膀胱癌的进展 [12] 。

随着对表观遗传学的不断深入研究，越来越多研究开始尝试从基因水平上探究疾病发生发展过程中的发病机制。生物信息学是一门研究生物学数据的收集、存储、管理、分析和解释的学科。它结合了生物学、计算机科学和统计学的原理和方法，旨在从大规模的生物学数据中提取有意义的信息和知识。生物信息学的研究内容包括基因组学、转录组学、蛋白质组学、代谢组学等领域，通过使用计算机算法和工具，可以对生物学数据进行序列比对、基因预测、蛋白质结构预测、基因表达分析等各种分析和解释。生物信息学预测和计算机技术成为在检测遗传改变和发现许多疾病生物标记物的有用且有前景的工具 [13] 。基因表达综合(gene expression omnibus, GEO)是一个国际公共信息库，负责记录和分发研究界提交的微阵列、下一代测序和其他形式的高通量功能基因组学数据。本研究通过生物信息学分析miR-125a-5p和miR-125b-5p关键基因，为探索肝硬化的诊断和治疗提供新思路。

2. 材料与方法

2.1. 数据获取

选定所有样本均为人类，通过数据库miRDB、miRWalk和TargetScan (见表1)分别获取miR-125a-5p和miR-125b-5p的靶基因。基于GEO数据库检索肝硬化基因表达谱芯片数据集，该数据集满足以下筛选条件：① 它们使用人类肝组织样本；② 每个数据集至少包含3个样本。得到数据集GSE14323 (正常肝组织样本19个，肝硬化样本41个)和GSE25097 (正常肝组织样本19个，肝硬化样本40个)。

2.2. 数据处理

将miRDB、miRWalk和TargetScan数据库预测的miR-125a-5p和miR-125b-5p的靶基因通过Bioinformatics绘制韦恩图。选取校正后P < 0.05和|logFC| > 0.585的基因作为差异表达基因(FC：fold change，差异倍数)，通过Bioinformatics对GSE825097和GSE14323所得差异表达基因取交集，得到共表达基因(differentially expressed genes, DEGs)，对GSE825097和GSE14323所得差异基因使用R语言中的相关数据包绘制热图和火山图。通过Bioinformatics对miR-125a-5p和miR-125b-5p靶基因和DEGs取交集得到核心基因。

2.3. 京都基因的基因组数据库(Kyoto Encyclopedia of Genes and Genomes, KEGG) 和基因本体论(Gene Ontology, GO)富集分析

对核心基因通过Metascape和WebGestalt数据库进行GO富集分析，通过微生信进行GO和KEGG富集分析。

2.4. 蛋白相互作用(Protein-Protein Interaction Network, PPI)网络的构建

通过String对所得交集基因构建PPI网络。

Table 1. The bioinformatics websites and analysis tools

表1. 所用的生物信息学网站及分析工具如表

3. 结果

3.1. miR-125a-5p和miR-125b-5p靶基因预测

通过数据库miRDB、miRWalk和TargetScan分别预测miR-125a-5p和miR-125b-5p的靶基因，将3个数据库所得miR-125a-5p和miR-125b-5p的靶基因分别取交集，最终获得miRDB 918、miRWalk 14,419、TargetScan 931个靶基因。将3个数据库所得靶基因通过Bioinformatics绘制韦恩图(见图1)，3个数据库重合的靶基因有515个。

Figure 1. Venn diagram of miRDB, miRWalk and TargetScan target genes

图1. miRDB、miRWalk和TargetScan靶基因Venn图

3.2. GEO数据库获取肝硬化数据集和筛选DEGs

共纳入了2个数据集，分别为GSE14323 (正常肝组织样本19个，肝硬化样本41个)，GSE25097 (正常肝组织样本19个，肝硬化样本40个)；选取校正后 P < 0.05 和|logFC| > 0.585的基因作为差异表达基因(FC：fold change，差异倍数)。以log FC > 0.585为上调基因，以log FC < −0.585为下调基因，获得GSE14323差异基因1789个，上调基因有1030个，下调的基因有759个，获得GSE25097差异基因2528个，上调基因有1236个，下调的基因有1292个。用Bioinformatics将GSE25097和GSE14323所得差异基因取交集，得到DEGs。使用R语言，绘制GSE25097和GSE14323热图和火山图(见图2)。

(a) (b)(c) (d)注：其中红色代表上调基因，绿色代表下调基因

Figure 2. Figures A and B are heat maps and volcano maps of GSE14323, respectively; Figures C and D are heat maps and volcano maps of GSE25097, respectively

图2. 图A、B分别为GSE14323热图和火山图；图C、D分别为GSE25097的热图和火山图

3.3. miR-125a-5p、miR-125b-5p靶基因与DEGs获取核心基因

将数据库miRDB、miRWalk和TargetScan所得miR-125a-5p和miR-125b-5p的515个靶基因和GEO数据库筛选出的肝硬化649个DEGs通过Bioinformatics制作Venn图(见图3)，获得核心基因17个，去除2个不符合筛选标准的基因，最后获得15个基因，其中包括9个下调基因，6个上调基因。用chiplot在线网站绘制15个核心基因的热图(见图4)。

Figure 3. Venn diagram of miRNAs and DEGs

图3. miRNA和DEGs的Venn图

注：其中红色代表上调基因，蓝色代表下调基因，颜色越深代表log FC绝对值越大

Figure 4. Heat map of intersecting genes

图4. 交集基因热图

3.4. 核心基因进行GO和KEGG功能富集分析

使用微生信、Metascape和WebGestalt数据库对核心基因进行GO分析(见图5)，该分析表明交集基因在生物过程(biological process)中主要在己糖代谢和单糖代谢等中富集，在细胞组成(cellular component)中主要线粒体外膜、血小板致密颗粒腔和细胞器外膜等中富集，在分子功能(molecular function)中主要在羟酸酯水解酶的活性和甘露糖寡糖甘露糖苷酶活性等中富集。在KEGG富集分析中(见图6)主要在脂肪酸生物合成和类固醇生物代谢等中富集。

Figure 5. Core gene GO enrichment analysis map

图5. 核心基因GO富集分析图

Figure 6. Map of KEGG enrichment analysis of core genes

图6. 核心基因KEGG富集分析图

3.5. 蛋白互作(PPI)网络构建

通过String对核心基因构建蛋白互作网络(见图7)，最终得到4个节点，4条边，得到4个关键基因，分别为NPL、H6PD、KIAA0319L和RFX5，其中NPL和RFX5为上调基因，H6PD和KIAA0319L为下调基因。

注：其中红色代表上调基因，蓝色代表下调基因

Figure 7. Interaction diagram of PPI network

图7. PPI网络互作图

4. 讨论

本研究基于数据库miRDB、miRWalk和TargetScan将miR-125a-5p和miR-125b-5p的515个靶基因与GEO数据库中肝硬化数据集GSE14323和GSE25097的647个靶基因取交集，得到17个核心基因，根据筛选标准，去除2个不符合标准的基因，最终获得15个核心基因，对所得核心基因进行GO和KEGG富集分析以及PPI网络构建，最终得到4个关键基因。

GO富集分析中主要在己糖代谢过程和线粒体外膜等中富集。人体摄入的70%的果糖由肝脏中肝细胞、HSC和库普弗细胞代谢，合成脂肪酸和甘油三酯，果糖超过肠道清除能力时，其被诱导溢出到肝脏 [14] 。果糖摄入通过增加肠道细胞色素P450-2E1 (CYP2E1)和硝基氧化应激，肠道微生物组发生变化、肠上皮细胞凋亡以及肠道TJ和AJ蛋白的酪氨酸硝化，导致肠道屏障功能障碍和内毒素血症升高，最终导致大鼠肝脏TLR4的激活以及NAFLD/非酒精性脂肪性肝炎(nonalcoholic steatohepatitis，NASH)和肝纤维化的发生 [15] 。CCl4诱导的小鼠肝脏中，肝葡萄糖及其代谢物己糖和果糖显著减少 [16] [17] 。谷胱甘肽过氧化物酶(GPX4)通过线粒体蛋白输入系统–外膜转位酶/内膜转位酶(TOM/TIM)复合体进入线粒体，GPX4主要通过线粒体自噬以及ROS诱导的线粒体受损而被降解，导致肝细胞铁死亡 [18] 。KEGG富集分析中主要在脂肪酸生物合成和类固醇生物合成中富集。在NAFLD的发病过程中，在白细胞介素6 (IL-6)的驱动下，向肝脏输送的脂肪酸增加，促进肝脏脂肪变性、胰岛素抵抗和血脂异常 [19] 。酒精摄入会改变肝脏中脂肪酸的数量和组成，在酒精喂养的啮齿动物中，饱和脂肪酸转变为不饱和脂肪酸，以及18碳脂肪酸(C18:0、C18:1、C18:2和C18:3)和22碳六烯酸(C22:6)在肝脏组织中升高 [20] 。类固醇激素分泌过多或不足与肝脏代谢或功能性疾病之间存在明显的双向联系，肝脏提供类固醇激素，通过调节合成和生物活性，类固醇激素也有助于维持正常的肝功能 [21] 。

通过PPI网络互作图的构建获得NPL、H6PD、KIAA0319L和RFX5 4个关键基因。N-乙酰神经氨酸醛缩酶(NPL)能够催化N-乙酰神经氨酸(Neu5Ac)裂解成N-乙酰–甘露糖胺(ManNAc)和丙酮酸，是唾液酸循环的关键酶之一，一般出现在需要分解唾液酸的生物中；炎性或恶性肿瘤的发生过程中血清和其它体液中唾液酸含量通常发生一定的改变 [22] 。细胞信息传递是指细胞内外信息的传递和处理过程。它涉及到多种信号分子、信号通路和细胞器的相互作用。其一般过程有六部分，分别为信号分子的释放、受体的识别、信号传导、信号放大、信号响应、负反馈调节。该过程是一个复杂而精密的调控网络，在维持细胞功能和适应环境变化中起着重要的作用 [23] 。唾液酸作为糖蛋白和糖脂的重要成分，在细胞信息传递、感染和转移过程至关重要，有研究表明唾液酸分解代谢途径在哺乳动物肌肉生理学中有重要作用 [24] 。己糖-6-磷酸脱氢酶(H6PD)，11β羟基类固醇脱氢酶1型(HSD11B1)是调节糖皮质激素活性的关键酶，主要表现为还原型烟酰胺腺嘌呤二核苷酸磷酸(NADPH)依赖型的还原酶活性，H6PD可以催化烟酰胺腺嘌呤二核苷酸(NAD)还原为NADPH，两者协同作用影响糖皮质激素的活性，因此HSD11B1与H6PD是糖、脂代谢的关键基因 [25] 。H6PD通过影响11β-羟基类固醇脱氢酶活性，调节脂肪组织胰岛素敏感性的皮质醇外周代谢，还可能通过影响肥胖、胰岛素抵抗和高雄激素血症而促成多囊卵巢综合征 [26] 。11-β羟基类固醇脱氢酶1型(HSD11B1)在H6PD介导的过程中将无活性可的松转化为活性皮质醇，两者的多态性可能在预防2型糖尿病(diabetes mellitus type 2, T2DM)个体肥胖方面相互作用 [27] 。KIAA0319L在系统性硬化症(systemic Sclerosis, SSC)和系统性红斑狼疮(systemic lupus erythematosus, SLE)患者的外周血细胞中过表达，KIAA0319L可作为SSC和SLE的易感位点 [28] 。对高级别鳞状上皮内病变(HSIL)和低级别鳞状上皮内病变(LSIL)患者的外周血样本进行基因表达谱分析，与LSIL相比，HSIL的KIAA0319L基因被上调，由此说明，KIAA0319L可能参与宫颈癌前病变，可用作宫颈鳞状上皮内病变(SIL)的诊断生物标志物和靶向治疗 [29] 。调节因子X-5 (RFX5)是免疫系统中的关键转录调节因子，与非肿瘤组织相比，KDM4A在HCC肿瘤组织中的表达水平显著升高 [30] 。RFX5在HCC和细胞系中的mRNA表达水平显著上调，RFX5可以与三肽基肽酶1 (TPP1)启动子区域结合，增加转录活性，TPP1在HCC组织中过表达，与HCC患者的不良预后显着相关 [31] 。

5. 结论

综上所述，本研究基于生物信息学对miR-125a-5p和miR-125b-5p以及肝硬化相关数据集的处理得到的相关通路以及筛选出的核心基因和最终的4个关键基因，初步探索肝硬化的发生机制，此结果有助于为肝硬化发病机制以及诊疗过程中的研究提供新思路。然而，本研究基于大量数据库对miR-125a-5p和miR-125b-5p以及肝硬化进行了一定的研究，缺乏肝硬化相关临床数据和实验的进一步验证，后续的研究应该进一步结合相关临床数据和数据对上述研究结果进行验证。