低氧胁迫下斑马鱼鳃microRNAs差异分析

林 枫，贾若南，王法祥，许强华,2

1. 上海海洋大学 海洋科学学院，上海 201306

2. 上海海洋大学/大洋渔业资源可持续开发教育部重点实验室/国家远洋渔业工程技术研究中心/远洋渔业协同创新中心，上海 201306

鱼类在水体中的生命活动会受到不同环境的影响，水中的溶解氧浓度变化对鱼类有至关重要的影响。水中溶解氧质量浓度在4.0 mg·L−1以上时，鱼类可以正常生长发育；低于1.0 mg·L−1时，大部分鱼类会出现浮头现象[1]。水中溶解氧浓度降低，不仅使鱼类呼吸和摄氧能力下降[2]，还会影响鱼类细胞的存活和信号传递[3]，并直接影响其产卵、交配、生长、发育等一系列生命活动[4]，严重时还会出现死亡，破坏种群内部的动态平衡。

斑马鱼 (Danio rerio) 具有易繁殖、发育快等优点[5]，一直作为模式生物用于生物医学研究[6-7]。作为模式生物，斑马鱼的应激调节功能、心血管系统功能等与人类高度相似[8-9]。鳃是鱼类重要的黏膜免疫器官和呼吸器官，在鱼类低氧适应研究中常作为主要的研究对象，如低氧胁迫下鳃组织中热休克蛋白的研究[10]，草鱼 (Ctenopharyngodon idellus)在低氧胁迫下鳃的差异蛋白质组学及热休克诱导[11]，低氧胁迫对卵形鲳鲹 (Trachinotus ovatus) 鱼体鳃器官的影响[12]等。

miRNAs可以通过调控其靶基因的表达水平从而参与细胞的各类进程，并参与许多关键的生理进程与病理过程[13]。已有多项研究表明，miRNAs是植物和动物面临环境胁迫作出响应的关键调节剂，一些miRNAs可以通过调控基因表达，恢复或重建新的表达程序，从而增强细胞对胁迫的耐受性[14]。低氧是生物常常面临的一种胁迫，目前已有一些应对低氧胁迫miRNAs的报道。如缺氧性神经胶质瘤来源的外泌体通过信号转导和转录活化因子 3 (STAT3) 和核因子 κB (NF-κB) 途径靶向端粒重复结合因子2 (terf2ip)，传递microRNA-1246诱导M2巨噬细胞极化[15]。miR-204可以作用于血管内皮生长因子 (vegf) 的3'-UTR区，是一个内源性的vegf表达调控因子。miR-204通过基因网络调控应答低氧胁迫[16]。miR-126-5p可以作为一种新型的miRNA，在缺氧条件下靶向白细胞介素17 (IL-17A)来调节大鼠心肌细胞 (H9c2) 的活力和凋亡[17]。

热休克蛋白 (Heat shock proteins, HSPs) 是生物在面对环境中的物理、化学、生物等刺激发生应激反应后大量产生的，常被称为应激蛋白[18]。热休克蛋白是生物体内最古老的分子之一。生物体为抵御环境变化所带来的刺激，会减少其他正常蛋白的合成，同时增加HSP的合成以应对环境的变化，帮助生物体恢复正常[10]。因此本文在对常氧和低氧胁迫下的斑马鱼鳃组织进行小RNA组测序的基础上，进一步筛选可能与低氧胁迫相关的micRNAs，并用其对热休克蛋白基因进行了靶基因预测，以期进一步挖掘斑马鱼的低氧适应机制。

1 材料与方法

1.1 材料

1.1.1 实验试剂

small RNA 分离试剂盒，DEPC 水 (生物生工有限公司)，氯仿 (24∶1)，Trizol 试剂，异丙醇，无水乙醇 (吉泰生物公司)。

1.1.2 实验器材

研磨棒、EP管、超净工作台、冰块、移液器、各类枪头、剪刀、镊子、锡纸均购买于上海生物生工有限公司，低氧驯化箱 (长沙华晓电子科技有限公司定制)，离心机，旋涡仪 (德国Eppendorf有限公司)。

1.2 方法

1.2.1 样本低氧处理

将本实验室培养的多代斑马鱼分别置于1.0 mg·L−1的低氧驯化箱中，保持该浓度进行低氧胁迫，同时保留6.7 mg·L−1的溶解氧质量浓度作为常氧对照组。在使用1.0 mg·L−1的溶解氧质量浓度进行低氧胁迫时，发现2周后斑马鱼不再出现浮头的现象，推测经2周低氧胁迫后其逐渐适应了低氧环境。故采用低氧驯化2周后的斑马鱼与常氧下的斑马鱼进行比对。低氧胁迫2周后，取出常氧/低氧条件下的斑马鱼，每个溶解氧质量浓度下取样本15尾 (体长3～4 cm)，解剖取出两组样本的鳃组织，备用。

1.2.2 高通量测序和生信分析

提取RNA样本，按照NucleoZOL试剂盒进行。将所提取的RNA进行电泳检测。检测后将EP管放入−80 ℃冰箱进行保存。将上述常氧和低氧下斑马鱼鳃组织提取后检测合格的RNA[19]，按照small RNA分离试剂盒的方法，分别进行small RNA的分离。将提取的small RNA进行纯化，运用荧光光度计进行定量，上机检测测序得到原始测序结果。为保证数据的质量，对原始数据进行处理。去除低质量的reads (N比例大于10%的reads、有5'接头污染的reads、无3'接头序列和插入片段的reads、3'接头序列以及polyA/T/G/C的reads)后得到的 clean small RNA reads数。其中大多数小RNA的长度为21～23 nt。对于该物种的ncRNA注释，若有该物种小RNA的注释信息，就用该物种ncRNA注释所测的small RNA。若没有该物种的信息，则选择Rfam数据库中rRNA、tRNA、snRNA 和 snoRNA 来注释测序所得的 small RNA。

1.2.3 靶基因预测

本实验室前期通过对低氧、常氧条件下鳃组织的转录组比较分析发现，常氧低氧条件下斑马鱼鳃中一共筛选获得28个显著差异表达的热休克蛋白基因，包括表达量显著下调基因12个，显著上调基因16个[12]。对miRNAs 测序和斑马鱼鳃转录组进行关联分析。针对低氧胁迫和常氧条件下斑马鱼鳃中显著差异表达的miRNAs，对实验室前期筛选获得的28个差异热休克蛋白基因进行靶基因的预测分析。miRNAs的靶基因预测使用TargetScan-Fish (http://www.targetscan.org/fish_62/)、miRanda(http://www.microrna.org/microrna/home.do) 2 个网站同时进行。

2 结果

2.1 低氧与常氧下斑马鱼鳃的小RNA测序

常氧和低氧斑马鱼鳃样本的小RNA组测序分别产生 6 995 009 和 6 662 504 bp 的原始数据。去除低质量的数据后分别得到 6 585 748 和 5 941 304 bp的 clean data。

常氧与低氧小RNA序列比对后获得了相应结果，获得饼状图 (图1)。根据饼图比例分析，可以看出，其中常氧的miRNA约占整个small RNA总数量的40%。而低氧的约占20%。提示与常氧相比，低氧胁迫下的miRNA有降低的趋势 (图1)。

图1 常氧与低氧中鳃组织miRNA占小RNA的比例Fig. 1 Proportion of miRNA in whole small RNAs sequence in normoxic and hypoxic gill tissues

2.2 低氧/常氧下鳃组织的差异 miRNAs 分析

根据测序结果，首先排除tRNA、rRNA等小分子的RNA。利用归一化法对比低氧胁迫斑马鱼与常氧斑马鱼中同一个miRNA的表达差异量。利用火山图呈现差异miRNAs的整体分布情况。结果显示，低氧/常氧条件下的斑马鱼鳃间一共筛选出差异表达的miRNAs共32个 (图2)。使用校正后的显著水平 (P) 和差异倍数 (Fold change) 2 个水平进行评估，设置显著差异表达miRNAs的筛选条件为P<0.01 和 |log2(fold change)|>1。

图2 常氧和低氧鳃组织差异miRNAs火山图Fig. 2 miRNAs volcano map of difference between normoxic and hypoxic gill tissues

针对32个差异miRNAs，同时使用 |log2FC|≥1，P<0.05，且表达量≥50作为临界值，从低氧和常氧的比较中，一共鉴定获得低氧与常氧条件下显著差异表达的15个miRNAs。其中，13个miRNAs在低氧胁迫斑马鱼鳃中的表达量显著上调、2个miRNAs (miR-455-3p、miR-125b-5p) 的表达量显著下调。图3为显著差异表达miRNAs的聚类。

图3 常氧和低氧鳃组织差异miRNAs聚类图Fig. 3 miRNAs clustering map of difference between normoxic and hypoxic gill tissues

2.3 差异 miRNAs 靶向热休克蛋白基因的预测结果

针对前期筛选获得的低氧胁迫与常氧条件下显著差异表达的28个热休克蛋白基因 (包括12个表达量显著下调的hsp基因、16个显著上调的hsp基因)[11]进行靶基因预测，结果显示，7个显著差异表达的miRNAs可以靶向9个热休克蛋白基因。图4显示出单个miRNA靶向热休克蛋白基因的数量。

图4 差异miRNAs靶基因数量统计分析Fig. 4 Statistical analysis of number of differential miRNAs target genes

表达量显著下调的miR-455-3p可以靶向2个显著上调的热休克蛋白基因 (表1)。表达量显著上调的 miRNAs (dre-miR-194a、dre-miR-155、dremiR-130c、dre-miR-9、dre-miR-29a、dre-miR-96-5p)可以靶向7个显著下调的热休克蛋白基因 (表2)。

表1 斑马鱼低氧与常氧鳃中下调的miRNAs靶基因预测Table 1 Prediction of down-regulated miRNAs target genes in hypoxic and normoxic gills of D. rerio

表2 斑马鱼低氧与常氧鳃中上调miRNAs靶基因预测Table 2 Prediction of up-regulated miRNAs target genes in hypoxic and normoxic gills of D. rerio

2.4 差异 miRNAs 靶向的热休克蛋白基因富集分析

针对差异miRNAs靶向的热休克蛋白基因进行富集分析发现，富集到的生物过程功能前6条通路主要与发育相关，包括dnajb6b和hspg2两个靶基因 (图5)。低氧胁迫的斑马鱼鳃中表达量显著上调的miR-455-3p靶向dnajb6b，而hspg2同时受到2 个上调的 miRNAs (miR-155 和 miR-29a) 调控。差异miRNAs靶向的热休克蛋白基因富集的前20条通路主要涉及的基因除了上面提及的2个miRNAs外，还包括miR-194a和miR-130c同时靶向的hspb7。

图5 差异miRNAs靶向的热休克蛋白基因GO富集Fig. 5 GO enrichment of differential miRNAs-targeted heat shock protein genes

3 讨论

本研究发现了在缺氧反应中差异表达的15个miRNAs，13个miRNAs在低氧胁迫的斑马鱼鳃中表达量上调，2个miRNAs的表达显著下调。也有研究表明，在下调的miRNAs中，miR-125b-5p的靶基因rps3a通过在翻译机制中发挥调节作用以抑制细胞凋亡[20]；而在上调的miRNAs中，mir-192可以通过靶向E盒结合锌指蛋白基因 (Zeb2) 来保护肝脏免受氧化应激诱导的损伤，增强肝脏的低氧耐受性[21-22]。在人类细胞中，miR-29b的上调可以靶向TNF受体相关因子5 (TRAF5) 保护心肌细胞免受缺氧诱导的细胞凋亡[23]。miR-29b在调节细胞凋亡中显示出重要作用[24]。另外，小鼠中miR-216b可以作用于自噬相关蛋白13基因 (Atg13)，且使缺氧条件下细胞的自噬减少，并减少细胞的凋亡[25]。

高海拔人群中的miR-210-3p水平与红细胞计数以及血红蛋白和血细胞比容显示出强正相关性，被认为是人类适应高海拔地区生活的重要miRNA[26]。miR-194过表达可以保护缺氧诱导的人肾皮质近曲小管上皮细胞 (HK-2) 损伤[27]。miR-155被发现在低氧条件下促进内皮细胞的血管生成[28]。人参皂甙 (GS-Rb1) 可以增加mir-29a的表达量，保护缺氧的心肌细胞[29]。大鼠中Hif-1α诱导的miR-9上调在缺氧期间有助于肺动脉平滑肌细胞的表型调节[30]，miR-96-5p已知有抑制细胞凋亡的功能[31]。miR-1可能在转录后水平直接或间接调控HSP90aa1和HSP90b1。过表达miR-1后缺氧复氧的HSP90蛋白及其亚型90aa1和90b1表达水平更低，结合之前的结果提示miR-1可能在心肌缺氧复氧中调控HSP90[32]。在高血压心肌肥厚的早期代偿阶段，心肌miR-378表达的下降使其对内源性HSF1转录后抑制作用减弱，进而对HSF1的代偿性升高发挥了重要的调控作用[33]。可见，miRNAs在低氧条件下在其他器官起着抑制细胞凋亡，增强器官的低氧耐受性，保护细胞免受缺氧带来的损伤等一系列功能。因此，本研究筛选出的低氧和常氧之间差异表达的miRNAs很可能在低氧适应机制中起重要作用。已知生物体在面对环境变化时会通过改变蛋白或者调节mRNAs的翻译以适应环境，本实验室已经发表过低氧胁迫下鳃组织相关热休克差异蛋白基因，验证了一些热休克蛋白基因应对低氧胁迫的作用[11]。因此，本研究利用筛选到的斑马鱼低氧与常氧状态下差异表达的15条miRNAs对28个热休克蛋白基因[11]进行靶基因预测，再结合miRNAs与预测的mRNAs的表达呈负相关这一特性，对预测的靶基因进行分析。

热休克蛋白被认为与正常和异常的胚胎发育密切相关。低氧胁迫下，低表达的miR-455-3p通过同时靶向提高hspa14和dnajb6b的表达，有可能增强了生物体的发育和机体保护，进而增强对低氧的适应。本研究发现，miR-194a同时靶向5个热休克蛋白基因 (hspa12a、dnajc5aa、hspb7、hsp70.3、dnajb2)；而miR-155可以同时靶向4个热休克蛋白基因 (hspa12a、hspg2、hspa13、dnajb2)。本研究还发现，热休克蛋白基因dnajb2同时受到5个miRNAs调控。基因dnajb2是HSP70的伴侣调节剂，主要在神经系统中表达[34]。等距遗传性运动神经病 (dHMN) 是一组罕见的遗传性神经肌肉疾病，其特征是在没有感觉症状的情况下会影响腓骨肌肉的萎缩，dnajb2是23个与dHMN有关的基因之一，主要从dnajb2起始[35]。可以推测，5个miRNAs靶向抑制dnajb2基因的表达，有可能减轻了对低氧环境下生物体神经系统的伤害。

Hspa14可能是肢体发育的相关基因[36]，在成年斑马鱼中进行无偏见的诱变筛选，确定了dnajb6b是心肌病的新型遗传修饰剂[37]。缺乏热休克蛋白基因hspg2可减轻在低氧中引起的动脉高压[38]。另外，在靶向的热休克蛋白基因GO功能富集前20条通路中，均有hspg2的参与。miR-155和miR-29a同时靶向hspg2，很可能通过抑制hspg2的表达来降低低氧胁迫下引起的动脉高压，增强对低氧环境的适应。

综上，本研究结合常氧与低氧下差异表达的miRNAs和热休克蛋白基因的关联分析，为探究鱼类低氧适应机制提供了新的研究思路。