基于生物信息学方法分析心肌梗死大鼠miRNA芯片

姚道敏,谢 亮,宫剑滨,朱 静,刘 晶

(1．南京中医药大学药学院，江苏 南京 210023；2．南京大学医学院附属金陵医院心脏内科，江苏 南京 210000)

心肌梗死是导致心脏重塑和慢性心衰发生的常见不良心血管事件，发病率和致死率居高不下，患者病理性心肌改变的分子机制仍未得到充分阐释[1]。微小RNA(microRNA，miRNA)作为非编码RNA的一种，长度在20个核苷酸左右，广泛存在于机体的各个组织细胞内，在基因的表达调控中具有重要的作用[2]。miRNA通常与信使RNA(messenger RNA,mRNA)的3′非翻译区( 3′untranslated region，3′UTRs)或编码序列(coding sequence,CDS)区结合抑制翻译或诱导mRNA降解。miRNA与靶基因的5′非翻译区(5′untranslated regions，5′UTRs)结合能激活转录[3]。有关MicroRNA在心肌梗死发生发展中的研究日益增多。MiR-15、miR-195、miR-122 、miR-181a等miRNA在心肌梗死后表达增加，而MiR-24、miR-138、miR-214、miR-499等miRNA在心肌梗死后表达减少。实验表明根据这一变化规律，相应得使miRNA过表达或抑制能有效改善心肌梗死的症状。特殊的是，MiR-320在心肌梗死后表达减少，过表达MiR-320不能减轻心肌梗死的症状反而促进心肌凋亡和坏死的发生。产生这一现象的原因是MiR-320抑制的靶基因对心脏具有保护作用[2]。MiRNA在心肌梗死中的作用复杂，本研究以microarray芯片结果为切入点，寻找心肌梗死中新的LncRNA-miRNA-mRNA三元关系。

1 材料与方法

1.1 大鼠心肌梗死模型的制备8-10周龄雄性清洁级SD大鼠(275±25) g,6只，由南京大学医学院附属金陵医院提供，实验动物使用许可证号为SCXK(军)2012-0014。随机分为2组(n=6):心肌梗死(MI)组和假手术(sham)组。以50 mg·kg-1的剂量腹腔注射戊巴比妥钠(产品批号：11751，Sigma Aldrich，USA)麻醉大鼠并将其固定在无菌手术台上。心电图机(CardiMax FX-7202，锡尔摩电气有限公司，中国)电极连接在四肢皮下，记录标准导联心电图。气管插管并连接呼吸机(HX-100E小动物呼吸机，成都市泰盟软件有限公司，中国)，控制潮气量在3 mL·kg-1，呼吸比2 ∶1，呼吸频率60-70次/min。纵行切开心脏最强波动处的皮肤约4 cm，钝性分离，于第3或4肋间开胸，将心脏暴露。用6-0线在左心耳根部下方3-5 mm处连同少量心肌组织缝扎左冠状动脉前降支。当缺血区心肌变苍白，运动减弱，心电图I、avL导联ST段弓背向上抬高0.1 mV持续30 min表明结扎成功。逐层关胸，置于37 ℃恒温下使大鼠苏醒。假手术组进行相同的操作，但不结扎冠状动脉。术后4 h再次麻醉大鼠，快速开胸取心脏，用眼科剪(南京大学附属金陵医院普外科)小心沿左室剪下梗死心肌组织，立即置于液氮中冷冻，然后置于干冰中送样。

1.2 miRNA芯片由上海康成生物技术有限公司完成大鼠梗死心肌组织样品分析和miRNA芯片杂交。简要地说，用TRIzol试剂(Invitrogen)提取缺血心肌细胞中的RNA，纯化后用nanodropND-1000分光光度计(美国，赛默)，测量RNA浓度，并通过凝胶电泳确定RNA完整性。接下来使用miRCURYTMHy3TM/Hy5TM功率标记试剂盒标记样品，并在miRCURYTMLNA阵列(v.18.0)上杂交。在洗涤步骤之后，使用Agilent Microarray Scanner 和Agilent Feature Extraction软件进行miRNA阵列扫描，并进行标准化处理。

1.3 lncRNA芯片由上海康成生物技术有限公司完成大鼠梗死心肌组织样品分析和长链非编码RNA(long non-coding RNA，LncRNA)芯片杂交。简言之，在去除核糖体RNA(Ribosomal RNA，rRNA)后，从总RNA中纯化mRNA(mRNA-ONLYTM真核mRNA分离试剂盒，Epicentre)。然后，使用Quick Amp Labeling Kit将每个样品扩增并沿着转录物的全长转录为荧光cRNA。使用RNeasy Mini Kit(Qiagen)纯化标记的互补RNA(complementary RNA，cRNA)。标记的cRNA的浓度和比活性通过NanoDrop ND-1000测量。加入5 μL 10×封闭剂和1 μL 25×Fragmentation Buffer裂解1 μg每个标记的cRNA，然后将混合物在60 ℃加热30 min，最后加入25 μL 2×GE杂交缓冲液稀释标记的cRNA。将50 μL的杂交溶液分配到垫片载玻片中，并组装到LncRNA表达微阵列载玻片上。将载玻片在安捷伦杂交炉中于65 ℃孵育17 h。洗涤固定杂交阵列，使用Agilent Microarray Scanner(Agilent p/n G2565BA)进行扫描，用Agilent Feature Extraction软件进行数据处理。

1.4 mRNA芯片由上海康成生物技术有限公司完成大鼠梗死心肌组织样品分析和mRNA芯片杂交。根据制造商的说明，使用TRIzol(Invitrogen)和RNeasy试剂盒(Qiagen)，从大鼠心肌组织中提取总RNA，包括DNase消化步骤，并完成RNA质量评估。使用安Quick Amp Labeling Kit对样品进行扩增和标记。然后将标记的探针与全大鼠基因组寡核苷酸微阵列试剂盒(包括41 000个大鼠基因和转录本；安捷伦)在65 ℃杂交缓冲液中杂交。过17 h，在室温下，先用0.5×SSC/0.01% SDS洗涤5 min，再用0.06×SSC洗涤2 min。用安捷伦DNA微阵列扫描仪扫描标记物强度，用Agilent Feature Extraction Software(版本10.5.1.1)提取数据结果，用安捷伦GeneSpring GX软件(版本10.0)进行进一步分析。使用Agilent FE单色方案在GeneSpring GX中对微阵列数据集进行归一化处理。中对微阵列数据集进行归一化处理。

1.5 差异表达的miRNA分析使用MEV软件(v4.9)进行差异表达的miRNA的系统聚类，以显示可区分的基因。根据基因之间的距离将基因分组在一起。系统聚类的默认距离度量是皮尔逊相关性，选择的链接方法是平均链接聚类。采用在线网络工具Ehbio(http://www.ehbio.com/ImageGP/)绘制差异表达的miRNA的火山图。

1.6 lncRNA-miRNA相互作用预测对筛选出来的差异表达的lncRNA和mRNA进行Pearson相关分析，并计算Pearson相关系数(PCC)。LncTar是一款通过自由能(ndG)最小化预测lncRNA-RNA相互作用的在线网络工具(http://www.cuilab.cn/lnctar)。ndG被认为是确定配对的RNA是否相互作用的临界值，一般设为-0.1。可以用于预测各种类型的RNA分子(mRNA，lncRNA，miRNA)和其他类型的非编码RNA之间的相互作用。

1.7 miRNA靶基因预测及功能分析miRNA靶基因预测使用miRWalk3.0(http://mirwalk.umm.uni-heidelberg.de/)进行，该软件集成了mirTarBaseV7.0，TargetScanV7.1和miRDBV5.0的预测结果。Bindingp≥0.95被认为是作为miRWalk中预测分析的临界标准。将预测的上调的miRNA靶基因与mRNA芯片中下调的mRNA取交集，同时将预测的下调的miRNA的靶基因与mRNA芯片中上调的mRNA取交集，用KOBAS在线工具(http://kobas.cbi.pku.edu.cn)对这两组交集基因进行GO和KEGG Pathway分析。GO分析的结果显示了基因的细胞定位，生物学过程和分子功能。KEGG分析结果用于确定基因参与的重要途径。P<0.05的标准用作定义目标基因及其相关功能和途径之间存在明显相关性的阈值。

1.8 lncRNA-miRNA-mRNA网络构建综合lncRNA-miRNA相互作用和miRNA靶基因预测及功能分析结果，采用Cytoscape软件(3.72版)进行lncRNA-miRNA-mRNA网络整合。

1.9 统计学处理使用2-ΔCt方法将lncRNA，mRNA，miRNA的相对表达水平表示为倍数变化。Fisher的精确检验用于GO分析和pathway分析。

2 结果

2.1 miRNA芯片结果分析按Foldchange≥1.5，P≤0.05的标准进行筛选，508个miRNAs结果中，与对照组相比，大鼠心肌梗死组织miRNA芯片结果中，13个miRNAs明显上调(miR-132-3p、miR-146b-5p、miR-18a-5p、miR-21-3p、miR-21-5p、miR-222-5p、miR-223-3p、miR-328b-3p、miR-365-5p、miR-449a-5p、miR-490-5p、miR-6333、miR-877)，6个miRNAs明显下调(miR-122-3p、miR-146a-3p、miR-181b-1-3p、miR-30c-1-3p、miR-3591、miR-3596c)。如Fig 1A所示，明显上调表达的13个miRNAs位于火山图中右侧，而明显下调表达的6个miRNAs则位于左侧。系统聚类图(Fig 1B)中可见差异表达的miRNAs的变化情况。

Fig 1 Volcano plot ( A) and hierarchical clustering ( B) of differentially expressed miRNAs

2.2 lncRNA芯片结果分析按照Foldchange≥2，P≤0.05的标准进行筛选，在心肌梗死组织和正常组织之间有75个lncRNAs明显差异表达，其中30个lncRNAs上调和45个lncRNAs下调。数据库中已经存在的lncRNAs有14个，其中10个lncRNAs (NR_133672、ENSRNOT00000057772、ENSRNOT00000050986、ENSRNOT00000058957、ENSRNOT00000002846、ENSRNOT00000027501、ENSRNOT00000071997、ENSRNOT00000075631、NR_126574)下调，4个lncRNAs(ENSRNOT00000076118、NR_132625、ENSRNOT00000076620、ENSRNOT000000 48374)上调。

2.3 mRNA芯片结果分析按照Fold change≥2，P≤0.05的标准进行筛选，心肌梗死组与对照组相比，差异表达的mRNAs共有331个，其中189个mRNAs明显上调，142个mRNAs明显下调。

2.4 lncRNA-miRNA相互作用结果分析采用LncTar对19个差异表达的miRNAs(Tab 1)和筛选的14个数据库中已经存在的差异表达的lncRNAs(Tab 2)进行相互作用预测，发现它们之间都有相互结合的作用位点，这一结果印证了miRNA在体内具有广泛的作用靶点。对差异表达的lncRNA和差异表达的miRNA进行Pearson相关系数计算，Pearson相关系数绝对值≥0.95的lncRNA和miRNA之间可能存在着正向或负向调控的作用。为进一步寻找潜在的lncRNA-miRNA-mRNA调控网络，这里仅汇总Pearson相关系数≤-0.95的lncRNA-miRNA的结果，见Tab 3。

Tab 1 Differentially expressed miRNAs

Tab 2 Differentially expressed lncRNAs

Tab 3 Analysis of interaction between LncRNAs and mRNAs

2.5 miRNA靶基因预测及功能富集分析按照P≤0.05，Foldchange≥1.5的标准筛选mRNA芯片中的差异mRNA(附Tab 1),与差异表达的miRNA预测的靶基因取交集的结果如韦恩图(Fig 2)所示，交集基因见Tab 4。GO和KEGG Pathway分析结果(Fig 3)表明miRNA的靶基因主要定位在细胞质和细胞器，富集在Hippo 信号通路，HIF-1信号通路，FoxO信号通路，小分子代谢过程，辅酶生物合成过程，辅因子的生物合成，脱甲基酶活性等条目中。从众多GO和KEGG pathway条目中挑选出的凋亡相关基因有Rnf7、I118、Sort1、Igf1、Tgfb2、Serpine1、Tyrobp、Lrp1、Sqstm1、Sema5a、Atp7a、Deptor、Ace、Tsc22d1、Stat3；自噬相关基因有Stat3、Deptor、Sqstm1；线粒体相关基因有Coq5、Qrsl1、Rnf7、Acsl1、Igf1、Stat3、Acsl4、Sqstm1、Ywhaz、Tmem65、Atp7a。这些基因与miRNA的对应关系如Tab 5所示。

Tab 4 Genes in intersection of Venn diagram

Tab 5 Correspondence between apoptotic,autophagy, mitochondrial related genes and miRNAs

Fig 2 Intersection of down-regulated mRNAs with predicted target genes of up-regulated miRNAs(A) and intersection of up-regulated mRNAs and down-regulated miRNAs predicted target genes(B)

2.6 lncRNA-miRNA-mRNA网络构建结果综合miRNA与靶基因的对应关系，miRNA与lncRNA相互作用结果，以及凋亡，自噬或线粒体相关基因，找出多条lncRNA-miRNA-mRNA三元关系。调控网络如Fig 4所示。LncRNA ENSRNOT00000076620的下调可能引起miR-146b-5p的上调，从而抑制Stat3，Sh3pxd2b、Tmem167a、Nop16、Sele、Rnf7、Qrsl1相关基因的表达。lncRNAENSRNOT00000057772的上调可能引起miR-181b-1-3p的下调，从而使对Tnpo1、Sdc4、Pcsk5基因的抑制作用减弱。LncRNA ENSRNOT00000071991的上调可能引起miR-122-3p的下调，从而使对Cdkn1a、Rfx3、RGD1559622、Atp7a、Ckap4、Deptor、Cdc42ep4、St3gal1基因的抑制作用减弱。LncRNA ENSRNOT00000075631的上调可能引起miR-146a-3p、miR-3591的下调，从而使对Rfx3、St3gal1、Slc6a6、Sdc4、Mboat1、Rcn1、Agfg1、Sh3pxd2b、Ace、Adamts1、Add2、Adamts12、Ddx21、Nolc1、Uchl1、Ankar、Mboat1基因的抑制作用减弱。

Fig 3 KEGG analysis of intersection genes of down-regulated mRNAs and predicted target genes of up-regulated miRNAs(A); KEGG analysis of intersection genes of up-regulated mRNAs and predicted target genes of down-regulated miRNAs(B); GO analysis of intersection genes of down-regulated mRNAs and predicted target genes of up-regulated miRNAs(C); GO analysis of intersection genes of up-regulated mRNAs and predicted target genes of down-regulated miRNAs(D)

Fig 4 LncRNA-miRNA-mRNA ternary relationshipRed ovals represent down-regulated miRNAs,Yellow oval represents up-regulated miRNAs.Blue triangles represent target genes.Purple triangles represent apoptotic and autophagy as well as mitochondrial-related genes.Green arrows represent up-regulated lncRNAs.Orange arrows represent down-regulated lncRNAs.

3 讨论

心肌梗死是临床常见病和多发病，严重威胁人类健康。研究发现，miRNA 通过调节心肌细胞凋亡、新生血管形成及纤维化等表型发挥功能[2]。LncRNA-Xist 5'发挥海绵作用，通过调节miR-29b-1-5p对靶标Bcl2l2的作用，调节心肌细胞凋亡[4]。LncRNA- APF靶向miR-188-3p，后者导致ATG7上调，从而促进心肌缺血/再灌注损伤中的心肌细胞自噬[5]。

本研究成功制备了心肌梗死模型，共获得19个明显差异表达的miRNAs。其中一些已通过大鼠心梗模型实验得到验证。它们是miR-21、miR-132、miR-222、miR-223-3p、miR-146a/b、miR-181b、miR-449a-5p、miR-122。临床试验表明miR-328、miR-18a是治疗心肌梗死的潜在候选分子。目前没有相关报道表明miR-365-5p、miR-490-5p、miR-6333、miR-30c-1-3p、miR-3591、miR-3596c、miR-877与心肌梗死相关。有文献表明miR-3596c、miR-877经qPCR验证在内毒素血症大鼠心肌组织中差异表达[6]。它们可能成为新的心肌梗死分子标志物，且其功能有待进一步研究。

MiR-132能够抑制心肌细胞凋亡，改善心肌重塑，减少大鼠心肌梗死面积。miR-132上调，导致白细胞介素1β下调，Bax和Cleaved caspase-3的表达水平降低，左心室射血分数、最大心率增加，左心室压力降低[7]。LncRNA-Ang362与邻近基因(miR-221和miR-222)共同发挥调控功能。这两个miRNA在响应血管紧张素II时过表达，并与LncRNA-Ang362共转录以调节大鼠血管平滑肌细胞的增殖[8]。miR-146b与miR-146a在3′末端有两个不同的核苷酸，具有相似的生物学功能。miR-146a/b通过直接结合IRAK-1、TRAF6和MyD88，抑制核因子κB途径的激活，减少缺氧诱导的大鼠心肌细胞凋亡。另外，miR-146a/b还靶向内皮生长因子受体，这可能会增加缺氧引起的大鼠心肌细胞损伤。核糖核酸酶L是miR-146b的直接靶标，抑制miR-146b可促进激活核因子κB和具有心脏保护作用的信号转导与转录激活因子Stat3，miR-146b过表达保护心肌细胞免受缺氧诱导的细胞凋亡[9]。过表达的miR-223直接结合RASA1,促进MEK1/2、ERK1/2和AKT的磷酸化，防止大鼠心肌梗死后心脏功能恶化和心肌纤维化[10]。同时过表达miR-223能沉默PARP-1，激活Akt/mTOR途径，减轻大鼠心肌细胞因缺氧诱导的细胞凋亡和过度自噬[11]。PDCD4是一种促凋亡蛋白，miR-21靶向PDCD4基因调控心肌细胞凋亡，过表达心肌miR-21基因能够改善心脏功能[12]。miR-181b-5P通过与AKT3和PIK3R3的3'UTRs直接结合，调节大鼠缺血心肌的细胞凋亡和组织损伤的发生发展[13]。lncRNA-MIAT作为竞争性内源RNA通过使miR-181b海绵化而增加STAT3的表达促进增殖并阻止凋亡[14]。lncRNA-XIST通过结合靶向Notch1基因的miR-449调控大鼠心肌细胞凋亡[15]。miR-122-5p的下调可以造成Bax、caspase-3和Bcl-2表达水平的上调，诱导大鼠心肌细胞抗凋亡。目前miR-122-3p的作用还不清楚[16]。临床研究表明miR-328和miR-18a在心肌梗死发生时，对心肌细胞具有保护作用[17-18]。

为进一步研究 miRNAs 的上下游作用靶点以及它们在心肌缺血缺氧状态下的作用机制，我们结合筛选出的差异lncRNA和miRNA靶基因预测及功能分析结果，找到了一些可能存在的lncRNA-miRNA-mRNA三元调控网络。在心肌梗死发生时，有可能存在如下作用机制:(1) LncRNA ENSRNOT00000076620/NR_132625的下调导致miR-146b-5p的上调，上调的miR-146b-5p通过下调STAT3，Rnf7激活心肌细胞凋亡，或通过下调Qrsl1引起心肌细胞线粒体功能障碍，加剧心肌梗死相关症状。(2) LncRNA ENSRNOT00000071991的上调引起miR-122-3p、miR-122-5P及其前体miR-122的下调，下调的miR-122通过上调Deptor(mTOR抑制剂)，抑制大鼠心肌细胞自噬的发生，加剧心肌细胞的凋亡过程。(3) LncRNA NR_132625的下调导致miR-21-3P，miR-18a-5p的上调，线粒体相关基因Coq5、Acsl1、Tmem65表达因此降低，心肌细胞线粒体受损，心肌梗死症状加重。