瑞士乳杆菌对泥鳅肠道转录组基因表达的影响

崔晓琳，丁志雯，米浩宇，孟文蓉，刘 姝，杨 光，侯晓月，董志国，房耀维

( 1.江苏海洋大学 海洋科学与水产学院，江苏 连云港 222005； 2.江苏海洋大学 食品科学与工程学院，江苏 连云港 222005 )

泥鳅(Misgurnusanguillicaudatus)分布于中国各地，主要分布于长江流域。由于其优越的营养价值和独特的风味，在东亚地区广受欢迎[1-2]。近年来，泥鳅在国内外的需求量年年攀升，是我国外贸出口的重要水产品之一[3]。基于此，泥鳅在我国一直是重要的淡水养殖品种。随着养殖技术的不断发展和进步，集约化和工厂化的养殖模式发展迅速，养殖规模不断扩大，但因养殖过程中出现的病害较多，抗生素使用广泛，耐药性和药物残留的问题也日益突出，严重影响农民的增产增收，一定程度上阻碍了泥鳅产业的健康发展[4-5]。

饲喂益生菌是一种持久的、环境友好型的且能够替代抗生素的鱼类疾病解决方案之一，使水产养殖业朝着绿色健康养殖的方向发展。益生菌作为水产养殖的安全添加剂，通过促进生长、提供营养、调节微生物定殖、提高免疫反应、提高饲料利用率、增加消化酶活性和消化率、改善水质和控制疾病等来提高宿主的健康水平[6-7]。目前报道的在水产养殖中应用的益生菌主要有酵母菌(Saccharomycetes)、芽孢杆菌(Bacillus)、乳杆菌(Lactobacillus)及光合细菌等[8]。Kim等[9]证明乳酸乳杆菌(L.lactis)BFE920对链球菌(Streptococcus)具有很强的抗菌活性，在褐牙鲆(Paralichthysolivaceus)的养殖中提高了鱼体的消化率和质量增加率。芽孢杆菌能够抑制潜在病原体在虾肠道内的定殖，从而提高凡纳滨对虾的免疫力和抗病性[10-12]。益生菌的作用机制除了直接的抑菌作用外，还有对消化道和免疫系统的影响。据报道，有益微生物可在养殖对象体内产生多种消化酶，从而提高饲料利用率[13-14]。乳杆菌以某种免疫调节因子的形式起作用，刺激肠道某种局部性免疫反应[15-16]。但是，与抗生素相比，益生菌的作用机理在理论研究上还不够深入，此外，乳酸菌在泥鳅的健康养殖中的研究较少，相关机制研究鲜有报道。

笔者对投喂瑞士乳杆菌(L.helveticus)HML037的泥鳅转录组基因的表达影响进行探索，以期解释瑞士乳杆菌HML037对泥鳅的生长和免疫反应的分子机制，为益生菌的促生长和免疫作用机制奠定理论基础，为乳酸菌在泥鳅的健康养殖提供理论支撑。

1 材料与方法

1.1 材料

泥鳅购于江苏连云港赣榆养殖技术服务中心，体质量(3.14±0.17) g。瑞士乳杆菌HML037自马奶中筛得，由内蒙古农业大学乌云达来教授所赠。

1.2 泥鳅的养殖

选取120尾泥鳅随机分到6个70 L聚氯乙烯水族缸(20尾/箱)中，试验分为2组：瑞士乳杆菌HML037组和对照组，每组3个平行。每日7：00和18：00对泥鳅进行投喂，投喂量为其体质量的3%，饲养8周。

1.3 组织收集和RNA制备

泥鳅喂养8周后，对照组和瑞士乳杆菌HML037组分别取泥鳅6尾，解剖取其肠道组织，并-80 ℃保存用于RNA制备。根据总RNA提取试剂盒(Trizol)说明书，对采集的泥鳅肠道组织样品进行RNA提取。用Qubit 2.0荧光计检测RNA浓度，琼脂糖凝胶检测RNA完整性以及基因组污染情况。最后将总RNA送至北京诺禾致源生物信息科技有限公司测序。

1.4 文库构建和RNA测序

通过Oligo(dT)磁珠(Illumina, 美国)富集带有polyA尾的mRNA，随后在NEB Fragmentation Buffer中用二价阳离子将得到的mRNA随机打断。以片段化的mRNA为模版，随机寡核苷酸为引物，在M-MuLV逆转录酶体系中合成cDNA第一条链，随后用RNaseH降解RNA链，并在DNA polymerase Ⅰ体系下,以dNTPs为原料合成cDNA第二条链。纯化后的双链cDNA经过末端修复、加A尾并连接测序接头，用AMPure XP beads筛选250～300 bp的cDNA，进行PCR扩增并再次使用AMPure XP beads纯化PCR产物，最终获得文库。cDNA文库在Illumina HiSeqTM2500上测序。

1.5 转录组De nove组装、质量评估和注释

转录组测序在Illumina HiSeq 2000平台上进行，该平台产生约100 bp的配对末端(PE)原始读数(北京诺禾致源生物信息科学技术有限公司)。去除衔接序列、不明确的“N”核苷酸(N的比率大于10%)和低质量序列(质量分数小于5%)后，使用Trinity软件组装剩余的干净读数，如无参基因组的从头转录组组装所述[17]。

将unigenes与公共数据库进行比较，并基于基因的相似性进行功能注释。使用NCBI公共数据库 (http://www.ncbi.nlm.nih.gov/)将unigenes与Nr、Nt、Swiss-Prot、GO、COG和KEGG数据库中的序列进行比较。如果不同数据库的结果有冲突，则遵循Nr、Nt、KEGG、Swiss-Prot、GO和COG数据库的比对优先顺序。

1.6 差异基因表达

使用DESeqR包(1.10.1)进行两个条件/组的差异表达分析。DESeq提供统计程序，用于使用基于负二项分布的模型确定数字基因表达数据中的差异表达。使用Benjamini和Hochberg的方法调整得到的P值以控制错误发现率。由DESeq发现的具有调整的P<0.05的基因被指定为差异表达。采用GOseq和KOBAS软件对差异基因集进行GO功能富集分析、KEGG通路富集分析等。富集分析基于超几何分布原理，其中差异基因集为差异显著分析所得差异基因并注释到GO或KEGG数据库的基因集，背景基因集为所有进行差异显著分析的基因并注释到GO或KEGG数据库的基因集。富集分析结果是对每个差异比较组合的所有差异基因集、上调差异基因集、下调差异基因集进行富集。

1.7 实时荧光定量PCR(qRT-PCR)

参照说明书，使用SYBR Premix Ex Taq试剂盒(Takara，日本)选择9个基因(TRP、TMSB10、APC11、COLEC12、TNFRSF10、IL-10、IL-1B、FGL和IL-5)，β-actin作为内参基因，通过实时荧光定量PCR验证RNA序列结果。根据转录组测序数据用Primer premier 5.0软件设计引物(表1)。用优化的比较Ct(2-ΔCt)值法估计基因表达水平。数据用3个平行试验的平均值表示[18]。

表1 实时荧光定量PCR引物

2 结果与分析

2.1 转录组测序和De nove组装

对对照组和瑞士乳杆菌HML037组的泥鳅进行转录组测序，分别获得21 626 795和23 112 691 个原始读长。两组数据Q20碱基比例超过95.6%，表明数据的质量适合后续分析。对原始数据过滤，即去除含有带接头的、低质量的序列，分别获得21 228 708和22 748 744个有效读长(表 2)。使用Trinity将有效数据De novo组装，获得82 944个转录本。将转录本去除冗余后获得40 877个unigenes。 unigenes中70.58%的长度大于500 bp，43.63% 的长度大于1000 bp，其长度为301～31 819 bp，平均长度为1218 bp(表3)。后续分析均以unigenes序列作为参考序列。

表2 转录组测序数据统计

表3 De nove组装结果统计

2.2 基因注释

基于BLAST算法将unigenes基因序列分别于Nr、Nt、KOG、PFAM、SwissProt、GO和KEGG数据库进行比对，取相似度>30%，且E<10-5的注释，合并基因得到的所有注释详细信息，并统计各数据库中unigenes基因注释数量及百分比(表4)。对于7个主要的数据库Nr、Nt、KOG、PFAM、SwissProt、GO和KEGG，其分别获得61.97%、69.63%、23.00%、47.24%、53.48%、49.16%和35.82%的注释，其中7691个unigenes在5个数据库中均有注释(图1)。

表4 基因注释统计

图1 蛋白数据库注释

GO是一套国际标准化的基因功能描述的分类系统。基于转录本与Swissprot、TrEMBL的注释结果，446 636个unigenes被分配到生物过程、分子功能和细胞组分中，其数量分别为297 873、52 533和96 230个(图2)。基于基因知悉同源关系，使用KOG数据库进一步对功能进行分类，9403个unigenes被分配到26个特定的KOG类别，其中信号传导机制最多(1582个unigenes)，其次是一般功能预测基因(1476个unigenes)，翻译后修饰、蛋白质转运及伴侣(993个unigenes)(图3)。KEGG是系统分析基因产物和化合物在细胞中的代谢途径以及这些基因产物功能的数据库。根据KO功能注释与Pathway的联系利用KAAS系统对13 580个unigenes进行KEGG代谢通路分类，获得34条通路途径注释。将34条通路分为5类(A～E)，其中unigenes的数量分别为3070、3147、1624、2977、和5160个(图4)。

图2 GO功能分类注释

图3 KOG功能分类注释

图4 KEGG分类注释

2.3 表达差异分析

引入TPM计算样品中基因表达丰度，根据表达丰度使用DEGseq软件进行差异表达分析，并通过设置Q<0.05且∣Fold Change∣>2来筛选显著差异基因。通过对对照组与瑞士乳杆菌HML037组比较，共鉴定出604个显著差异基因，其中，326个基因显著上调，278个基因显著下调。比较两组的基因表达水平(图5)，横坐标表示基因在不同样本中的表达倍数变化[log2(Fold Change)]，纵坐标表示表达差异的显著性水平[-log10(padj)]，上调基因用红色点表示，下调基因用绿色点表示。在这些差异基因中，注释到几种与生长和免疫相关的基因，如生长素(GH)、胰岛素生长因子-Ⅰ(IGH-Ⅰ)、Toll样受体、MyD88和主要组织相容性复合体-Ⅱ(MCH-Ⅱ)基因，具有调控细胞增殖、分化及蛋白合成，抗肿瘤和抗病菌感染等功能。

图5 差异基因表达水平

2.4 GO和KEGG数据库对差异基因富集分析

利用GO和KEGG富集技术进一步分析差异基因，从而确定其潜在功能和代谢途径。将差异基因映射到GO数据库进行功能富集分析，差异基因总共显著富集到436个，显著富集到5条通路中，其中生物过程包含0个，细胞组分包含1个，分子功能包含4个(图6)。细胞成分功能类型中富集在氧化还原过程。分子功能富集在氧化还原酶活性、G蛋白偶联受体结合、趋化因子活性和趋化因子受体结合。

图6 GO功能富集分析

将差异基因映射到KEGG数据库进行通路富集分析，当矫正后的P值Q<0.05时，认为该功能存在显著富集情况，且Q值越小越显著富集。结果显示，26 545个差异表达基因总共富集到215条通路上，其中差异表达基因主要集中在核糖体、吞噬体、丙酮酸代谢、过氧化物酶体增殖激活受体信号通路、糖酵解/糖质新生、脂肪酸延长和脂肪消化和吸收，其中核糖体、丙酮酸代谢、脂肪酸延长和脂肪消化吸收通路显著富集(图7)。

图7 KEGG通路富集分析

2.5 qRT-PCR结果

为了进一步验证高通量测序结果的可靠性，通过qRT-PCR验证泥鳅9个重要功能基因的表达，即TRP、TMSB10、APC11、COLEC12、TNFRSF10、IL-10、IL-1B、FGL和IL-5。将对照组与瑞士乳杆菌HML037组的差异基因与转录组测序结果相比较(图8)，结果显示，被检测的基因表达趋势与转录组表达分析结果基本一致，即证明转录组测序结果可靠。

图8 功能基因qRT-PCR与转录组之间的相对表达水平

3 讨 论

笔者通过转录组分析鉴定了泥鳅对瑞士乳杆菌HML037的应答基因，为鱼类的生长和免疫机制的研究提供了分子模型。从对照组和瑞士乳杆菌HML037组共获得 40 877个平均大小为1218 bp 的unigenes，且被注释到公共蛋白数据库中。RNA-seq结果的相关性分析表明RNA-seq数据的可靠性和理论分析的可行性。

3.1 瑞士乳杆菌 HML037对泥鳅GO注释的影响

GO注释显示，瑞士乳杆菌HML037组上调表达的基因明显多于对照组。瑞士乳杆菌 HML037对泥鳅的氧化还原酶活性、G蛋白偶联受体结合、趋化因子活性和趋化因子受体结合等影响显著。根据表达方式和免疫系统中的作用，可以把趋化因子分为内环境稳定性趋化因子和炎性趋化因子，前者有维持内环境稳态的功能，对淋巴细胞归巢、成熟起着重要作用；后者主要功能是募集效应细胞，在协调先天免疫和获得性免疫反应中起重要作用[19]。趋化因子通过与细胞表面的特异性受体结合发挥相应的生物学作用，趋化因子受体属G蛋白偶联受体[20]。趋化因子首先作为趋化性介质被发现，但随着对趋化因子研究的深入，他们不仅对白细胞具有趋化作用，在免疫细胞和器官的发育、免疫应答过程、炎症反应、病原体感染、创伤修复及肿瘤形成和转移等方面也发挥广泛的生理和病理作用[21]。在对照组和瑞士乳杆菌HML037组的差异基因比较中显示，Toll样受体、MyD88和MCHⅡ等免疫相关基因均呈显著上调。有研究表明，乳杆菌可以上调MHC、CD40或CD80基因的表达[22]。这与本研究结果相一致，表明瑞士乳杆菌HML037能够激活Toll/MyD88和MCH信号通路，非特异性免疫应答和免疫基因的表达很可能是因为免疫系统的激活增强，使其具有较高的免疫状态，更好保护其被病原体入侵，减少养殖中抗生素使用的风险。

3.2 瑞士乳杆菌HML037对泥鳅KEGG通路的影响

KEGG数据库中包含丰富的通路信息，这将有助于在系统水平去了解基因的生物学功能[23]。对差异基因进行KEGG富集分析，结果显示，差异表达基因主要集中在核糖体、吞噬体、丙酮酸代谢、过氧化物酶体增殖激活受体信号通路、糖酵解/糖质新生、脂肪酸延长和脂肪消化和吸收等过程。脂肪酸延长酶是高不饱和脂肪酸合成的关键酶之一，高不饱和脂肪酸是一类重要的生物学功能物质，在鱼类的生长发育过程中发挥重要的作用[24-26]。糖酵解是机体重要的糖代谢途径，丙酮酸激酶是糖酵解过程中重要的限速和调节酶[27]。过氧化物酶体增殖物激活受体信号通路可能是线粒体代谢和生长在基因转录水平的共同调控通路[28]。此外，笔者筛选出了生长素、胰岛素生长因子-Ⅰ、Δ6脂肪酸去饱和酶等重要的差异基因。胰岛素生长因子基因的表达增加，直接调控草鱼(Ctenopharyngodonidella)肌肉的细胞分化，从而使其体质量增加，生长加快[29-31]。相关基因表达的增加，可能与泥鳅机体的蛋白质合成、糖代谢和脂质代谢呈正相关。这些基因表达增加表明，瑞士乳杆菌HML037具有调控生理功能的作用，从而使泥鳅的生长加快和营养调节引起的内部应答加快。但是，导致这些变化的分子机制尚不清楚，需进一步的科学研究。

4 结 论

测序获得泥鳅的转录组数据，成为泥鳅转录组学研究的宝贵资源，也将有助于其他重要水产经济品种的相关研究。饲料中添加瑞士乳杆菌HML037对泥鳅的生长和免疫相关的基因和通路具有上调和激活作用，可能提高泥鳅的生长速度和抗病力。由此可知，瑞士乳杆菌HML037在泥鳅养殖过程中可以替代抗生素的部分作用，减少抗生素的使用，进而达到泥鳅健康养殖的目的。