中华蜜蜂6日龄幼虫肠道响应球囊菌胁迫的差异表达基因分析

郭 睿，张 璐，徐细建，史秀丽，熊翠玲，郑燕珍，付中民，黄枳腱，王鸿权，侯志贤，陈大福*

中华蜜蜂6日龄幼虫肠道响应球囊菌胁迫的差异表达基因分析

郭 睿1*，张 璐1*，徐细建1，史秀丽2，熊翠玲1，郑燕珍1，付中民1，黄枳腱1，王鸿权1，侯志贤1，陈大福1**

(1.福建农林大学蜂学学院，福州350002; 2.新疆维吾尔自治区蜂业技术管理总站, 乌鲁木齐 930001)

蜜蜂球囊菌特异性侵染蜜蜂幼虫而导致白垩病，严重危害养蜂生产。本研究利用RNA-seq技术对健康及球囊菌胁迫的中蜂幼虫肠道进行深度测序，进而对宿主的差异表达基因进行深入分析。本研究中，幼虫肠道样品的RNA-seq共得到191167730条原始读段(raw reads)，经过滤得到186284296条有效读段(clean reads)，差异表达基因(DEG)分析结果显示上调与下调基因的数量分别为4513和385个。Gene ontology(GO)富集分析结果显示，上调基因富集在45个GO条目(term)，富集基因数最多的是细胞进程、代谢进程及催化活性，下调基因富集在32个 GO term，富集基因数最多的是代谢进程、单组织进程及催化活性， KEGG代谢通路(pathway)富集分析结果显示上调基因富集在193个pathway，其中富集基因数最多的是核糖体、氨基酸的合成、碳代谢。下调基因富集在59个pathway，其中富集基因数最多的是甘氨酸、碳代谢以及二羧酸代谢。深入分析发现宿主的细胞免疫被显著激活，体液免疫中的Toll-like与Jak-STAT信号通路也被球囊菌所激活。研究结果为揭示中蜂幼虫在球囊菌入侵后期的胁迫应答机制提供了重要的信息，也为解析中蜂幼虫的球囊菌抗性机制奠定了基础。

中华蜜蜂；幼虫肠道；球囊菌；差异表达基因；免疫防御

蜜蜂是一种重要的社会学模式昆虫，因其在发育学、神经生物学、行为学和病原-宿主互作研究中的重要性而备受关注(Galiziaetal., 2012; Begnaetal., 2012; Zayed, 2012; Foretetal., 2012; Kurzeetal., 2016)。蜜蜂作为最重要的授粉昆虫，在农业生产和生态维持中也发挥着不可替代作用(Committee on the Status of Pollinators in North Acerica, 2007)。据报道，蜜蜂为全球70%的作物和野生植物授粉(Kleinetal., 2007; Elke, 2010)。蜜蜂因其群遭受细菌、真菌及病毒等病原的侵袭。蜜蜂白垩病是一种最具代表性的致死性真菌病，1913年Massen在德国首次报道发现白恶病(Aronstein and Murray, 2010)，中国大陆1990年发生白垩病(Liang and Chen, 2008)。近年来，随着养蜂活动及蜂产品全球贸易的快速发展，白垩病发病率逐年上升(Aizenetal., 2009)。白垩病是由蜜蜂球囊菌Ascosphaeraapis(简称球囊菌)特异性侵染蜜蜂幼虫而导致，可造成蜜蜂群势的大幅下降，从而严重影响蜂蜜等产品的产量(Bailey, 1963)，据报道，白垩病可造成蜂蜜产量下降5-37%(Zaghlouletal., 2005)。近二十年来，国内外学者在病原分类鉴定、形态学、病理学、流行病学、侵染过程、蜜蜂防御以及疾病防治等方面对白垩病开展了一系列研究。本课题组也在球囊菌的生化、检测及侵染过程等方面开展了较为系统的研究，如梁勤等从碳源、氮源、维生素、矿质元素等方面研究了营养生态条件对球囊菌生长及产孢的影响，结果表明营养生态条件的变化对球囊菌的影响极大(Liangetal., 2001)；郑志阳等对健康和患病蜜蜂幼虫血淋巴进行 SDS-PAGE电泳和蛋白酶、酯酶的活性染色，发现健康蜜蜂幼虫血淋巴中的蛋白含量丰富，主要由4种高分子质量的蛋白组成，而患病幼虫血淋巴中的蛋白含量很少，主要蛋白组分被降解，多种蛋白酶和酯酶的活性在患病幼虫血淋巴中检测到，但在健康幼虫中检测不到(Zhengetal., 2011)。

我国养蜂生产的主要蜂种是意大利蜜蜂和中华蜜蜂。中华蜜蜂Apisceranacerana(简称中蜂)基因组的公布(Parketal., 2015)，为中蜂的分子生物学研究提供了重要参考信息。Aronstein等利用cDNA-AFLP 对健康及球囊菌感染的西方蜜蜂Apismellifera幼虫进行了比较，结果表明差异表达基因(DEGs)参与了宿主的能量代谢和蛋白转运，其中的类几丁质编码基因很可能参与了蜜蜂幼虫对球囊菌的抵抗(Aronsteinetal., 2010)。Cornman等对来自培养基的球囊菌菌丝和来自蜜蜂幼虫感染组织的球囊菌菌丝进行了转录组测序，功能分析表明球囊菌的DEG参与了交配类型、细胞内信号转导和应激反应。

目前，白垩病的相关研究主要集中在意蜂，有关球囊菌侵染中蜂的研究报道极少。前期我们发现中蜂蜂群偶尔可见白垩病患病幼虫，从患病中蜂幼虫上分离培养真菌病原，经形态学、分子生物学以及交叉感染实验证实确为A.apis(未发表数据)。

本研究利用Illumina测序技术对对健康及球囊菌胁迫的中蜂6日龄幼虫肠道进行深度测序，得到宿主的差异表达基因(DEGs)，并通过Gene ontology(GO)和KEGG代谢通路(pathway)富集分析深入研究DEGs。研究结果可为揭示中蜂幼虫肠道响应球囊菌后期胁迫的应答机制提供重要信息，也能为关键应答基因的筛选及验证奠定基础。

1 材料与方法

1.1 生物材料

本研究中使用的中蜂幼虫取自福建农林大学蜂学学院教学蜂场，球囊菌菌株由福建农林大学蜂学学院蜜蜂保护实验室保存并活化。

1.2 主要实验试剂及仪器

DNaseI和Oligotex mRNA Kits Midi试剂盒购自德国Qiagen公司，Dynal M280磁珠购自Invitrogen公司，DNA ligase购自美国Thermo公司，RNA Reagent抽提试剂盒、Ex Taq polymerase及Superscript II reverse transcriptase均购自日本TaKaRa公司，纯化cDNA的Acpure beads为美国Agencourt产品，cDNA文库构建试剂盒TruSeqTMDNA SAcple Prep Kit-Set A为美国Illumina公司产品。其它试剂均为国产分析纯。

倒置显微镜为中国上海光学仪器五厂产品，超净工作台为中国苏州安泰空气技术有限公司产品，恒温恒湿气候箱购自中国宁波江南仪器厂，凝胶成像系统为中国上海培清科技有限公司产品，PCR仪为美国Bio Rad公司产品，超低温冰箱为中国中科美菱公司产品。

1.3 球囊菌活化、孢子纯化及计数

按照本实验室已建立的方法对球囊菌进行活化(Zhangetal., 2017)。按照Jensen等(2013)的方法纯化球囊菌孢子，将高浓度孢子溶液梯度稀释后用血球计数板对孢子进行计数。

1.4 中蜂幼虫的人工饲养及肠道样品准备

中蜂幼虫的人工饲养参照王倩等(2009)的方法。从蜂学学院蜂场群势较强的中蜂蜂群用移虫针挑取2日龄幼虫，放入无菌的24孔细胞培养板(每孔对应1只幼虫，孔内加有35℃预温的幼虫饲料)，将24孔板放入恒温恒湿培养箱，35℃，70%相对湿度(RH)条件下饲养。每隔24 h更换饲料。预先配制添球囊菌孢子的人工饲料，混匀后调整孢子终浓度至为1×107孢子/mL，饲喂处理组3日龄幼虫，对照组饲喂正常人工饲料。本次实验进行3次生物学重复。

本研究的分析重点是中锋幼虫肠道在球囊菌胁迫后期的DEGs。鉴于后续将对宿主在球囊菌入侵全过程的胁迫应答进行研究，为降低测序成本，拟将健康中蜂4日龄幼虫肠道作为唯一对照。因此，本研究中，分别剖取对照组4日龄幼虫肠道(AcCK)和处理组6日龄幼虫肠道(AcT)，AcCK与AcT的三个生物学重复分别为AcCK-1、AcCK-2、AcCK-3和AcT-1、AcT-2、AcT-3。每剖取一只幼虫肠道，迅速将肠道移至RNA Free的EP管，液氮速冻，待一组肠道样品(7只幼虫肠道)集齐后，转移保存于-80℃。

1.5 cDNA文库构建及Illumina测序

利用RNAiso Reagent试剂盒抽提处理组和对照组上述6头幼虫肠道的总RNA，然后用RNase-free DNaseI去除基因组DNA残留。RNA的质量通过琼脂糖凝胶电泳和NanoDrop ND-1000(NanoDrop, Wilmington, DE, USA)进行检测。cDNA文库构建参照张曌楠等的建库方法(Zhangetal., 2017)。委托广州基迪奥生物科技有限公司对上述幼虫肠道样品进行双端测序，测序平台为Illumina HiSeq2500。

1.6 数据分析

对于下机数据，利用Perl脚本去除含有adaptor、未知核苷酸比例大于5%和低质量reads，获得有效读段(clean reads)。利用R软件(Version 2.16.2)进行测序饱和度分析。使用短 reads 比对工具 bowtie(Langmeadetal., 2009)将clean reads映射(mapping)到核糖体数据库(最多允许5个错配)，去除比对上核糖体的 reads，将保留下来的数据用于转录组的组装及分析，进而利用SOAP aligner/soap2软件(Hurgobin, 2016)将未比对上核糖体的 reads mapping到中锋幼虫肠道参考转录组(前期已组装并注释，原始数据已上传NCBI SRA数据库，SRA号: SRA456721)。

利用FPKM(Fragments Per Kilobase of transcript per Million mapped reads)法计算基因表达量。利用R软件(version 2.16.2)计算各样品之间的相关性系数。利用edgeR软件(Robinsonetal., 2010)进行DEGs分析。DEGs的筛选标准为FDR ≤ 0.05且|log2Fold change | ≥ 1。将DEGs向GO数据库(http://www.geneontology.org/)的各条目(term)mapping，并计算每个term的基因数，从而得到具有某个GO功能的基因列表及基因数目统计，然后应用超几何检验，找出与整个基因组背景相比，在差异表达基因中显著富集的GO条目。KEGG(pathway)显著性富集分析以KEGG pathway为单位，应用超几何检验，找出与整个基因组背景相比，在DEGs中显著性富集的pathway。

2 结果与分析

2.1 RNA-seq数据质控与评估

中蜂幼虫肠道样品的转录组测序共测得191167730条raw reads，经过滤得到186284296条clean reads，各样品clean reads数均在26509638 (97.96%)以上(表1)。两端平均Q20为98.38%，两端平均Q30为为95.94%。随着测序量的增多，检测到的基因数也随之上升、增长速度趋于平缓，说明本研究的测序深度检测到的基因数趋于饱和(附图1)。AcCK与AcT的组内各生物学重复之间的相关性均在0.96以上，说明样本的重复性高(图1)。上述结果说明本研究的转录组数据质量良好，可用于进一步分析。

表1 RNA-seq数据统计Table 1 Overview of RNA-seq data

图1 各幼虫肠道样品不同生物学重复间的相关性.Fig.1 Pearson correlation between every two biological repeats within each Apis cerana cerana larval gut sample 注: A: AcCK-1与AcCK-2间的相关性; B: AcCK-2与AcCK-3间的相关性; C: AcCK-1与AcCK-3间的相关性; D: AcT-1与AcT-2间的相关性; E: AcT-2与AcT-3间的相关性; F: AcT-1与AcT-3间的相关性.Note: A: parson correlations between AcCK-1 and AcCK-2; B: pearson correlations between AcCK-2 and AcCK-3; C: pearson correlations between AcCK-1 and AcCK-3; D: pearson correlations between AcT-1 and AcT-2; E: pearson correlations between AcT-2 and AcT-3; F: pearson correlations between AcT-1 and AcT-3.

2.2 DEGs分析

DEGs分析结果显示，在AcCK vs AcT中共有4898个基因差异表达，其中，上调基因和下调基因的数量分别为4513和385个(表2)，上调基因的数量远远多余下调基因，说明在球囊菌胁迫后期，中蜂幼虫肠道的绝大多数基因被病原激活表达。

表2 差异表达基因统计Table2 Summary of DEGs

2.3 DEGs的GO分类

DEGs的GO分类结果显示，这些DEGs分为三类：生物学进程(biological process)、细胞组分(cellular component)和分子功能(molecular function)，上调基因分布于45个 GO term上，富集基因数最多的是细胞进程(cellular process)、代谢进程(metabolic process)、催化活性(catalytic activity)(图2)。下调基因分布于32个 GO term上，富集基因数最多的是代谢进程(metabolic process)、单组织进程(single-organism process)、催化活性(catalytic activity)(图2)。

2.4 DEGs的KEGG pathway富集分析

DEGs的KEGG pathway富集分析结果显示，上调基因富集在193个pathway，其中富集基因数最多的是核糖体(ribosome)、氨基酸生物合成(biosynthesis of Acino acids)以及碳代谢(carbon metabolism)。下调基因富集在59个pathway，其中富集基因数最多的是甘氨酸(biosynthesis of acino acids)，碳代谢(carbon metabolism)，二羧酸代谢(glyoxylate and dicarboxylate metabolism)(图3)。上述结果表明球囊菌侵染对中蜂幼虫肠道的物质代谢产生较大影响。进一步分析结果显示，有2、17、20、29和52个上调基因富集在凋亡(apoptosis)、溶酶体(lysosome)、泛素介导的蛋白水解(ubiquitin mediated proteolysis)、吞噬体(phagosome)和内吞作用(endocytosis)，说明宿主的细胞免疫在胁迫后期被显著激活；分别有1和2个上调基因富集在中的Toll-like受体信号通路与Jak-STAT信号通路，说明宿主的此二条体液免疫通路在球囊菌胁迫后期被激活(图4)。

3 结论与讨论

白垩病是养蜂生产的一大顽疾，每年给养蜂业造成巨大损失。前人研究主要集中在意蜂幼虫白垩病的诸多方面，而球囊菌侵染中蜂幼虫的研究进展几无报道。中蜂白垩病虫尸在自然蜂群中仅偶尔可见，本课题组前期已从自然蜂群中蜂白垩状虫尸上分离得到球囊菌，通过形态学和分子生物学手段证明该病原即为A.apis(未发表数据)。前期研究中，我们组装并注释了中蜂幼虫肠道的参考转录组并开发出15个SSR分子标记(Xiongetal., 2017)，在此基础上，本研究利用RNA-seq技术对健康及球囊菌胁迫的中蜂幼虫肠道进行转录组测序，进而对宿主响应球囊菌胁迫的DEGs进行深入分析。肠道是昆虫的重要免疫器官，在抵御病原微生物入侵过程中发挥重要作用。本研究选择中蜂幼虫肠道作为测序对象，其转录组变化能更为精确地反映宿主响应球囊菌胁迫的应答表现。

我国养蜂生产中的常用蜂种是意蜂和中蜂，中蜂作为我国养蜂生产的主要蜂种之一，较意蜂具有更强的球囊菌抗性。本研究发现球囊菌胁迫的中蜂6日龄幼虫肠道上调基因的数量(4513 unigenes)远远多于下调基因(385 unigenes)，宿主的绝大多数基因被球囊菌激活表达，说明宿主响应胁迫的应答活跃，这或许与中蜂幼虫的球囊菌抗性密切相关。球囊菌被 意幼虫摄入中肠后，在整个幼虫期因中肠缺氧而不萌发，至幼虫期结束进入蛹期后，此时蜜蜂的中后肠接通，球囊菌孢子伴随蛹便进入后肠，在氧气的刺激下在此迅速萌发生长，1-2 d内菌丝即突破体表而导致宿主死亡(Lietal., 2012)。本研究中，下调基因共富集在59个pathway，其中的49个pathway与新陈代谢相关，包括物质代谢(如氨基酸生物合成和半乳糖代谢)和能量代谢(如氮代谢)，说明球囊菌在胁迫后期通过病原-宿主互作对中蜂幼虫肠道的新陈代谢系统产生较强抑制。

图2 差异表达基因的GO分析Fig.2 GO analysis of DEGs between AcCK and AcT注: A: 上调基因; 1: 行为; 2: 生物附着; 3: 生物调控; 4: 细胞成分组织或生物合成; 5: 细胞进程; 6: 发展进程; 7: 生长; 8: 免疫系统进程; 9: 定位; 10: 运行; 11: 代谢进程; 12: 多组织进程; 13: 多细胞组织进程; 14: 生殖; 15: 生殖进程; 16: 应激; 17: 信号; 18: 单一有机体进程; 19: 细胞; 20: 细胞接合; 21: 细胞零件; 22: 细胞外基质; 23: 细胞外区域; 24: 细胞外区域组分; 25: 大分子复合物; 26: 细胞膜; 27: 细胞膜组分; 28: 细胞膜膜蛋白; 29: 细胞器; 30: 细胞器组分; 31: 突触; 32: 突触组分; 33: 病毒; 34: 病毒组分; 35: 抗氧化剂活性; 36: 结合; 37: 催化活性; 38: 电子载体活性; 39: 酶的调节; 40: 脒基核苷酸交换因子活性; 41: 分子功能调节; 42: 分子转导活性; 43: 核酸结合转化因素; 44: 分子结构活性; 45: 运输活性.B: 下调基因; 1: 生物附着; 2: 生物调控; 3 细胞成分组织或生物合成; 4: 细胞进程; 5: 发展进程; 6: 定位; 7: 运行; 8: 代谢进程; 9: 多细胞组织进程; 10: 应激; 11: 信号; 12: 单一有机体进程; 13: 细胞; 14: 细胞零件; 15: 细胞外基质; 16: 细胞外基质组成; 17: 细胞外区域; 18: 细胞外区域组分; 19: 大分子复合物; 20: 细胞膜; 21: 细胞膜组分; 22: 细胞器; 23: 细胞器组分; 24: 病毒; 25: 病毒组分; 26: 结合; 27: 催化活性; 28: 酶的调节; 29: 分子功能调节; 30: 分子转导活性; 31: 分子结构活性; 32: 运输活性。 Note: A: up-regulated genes; 1: behavior; 2: biological adhesion; 3: biological regulation; 4: cellular component organization or biogenesis; 5: cellular process; 6: developmental process; 7: growth; 8: immune system process; 9: localization; 10: locomotion; 11: metabolic process; 12: multi-organismal process; 13: multicellular organismal process; 14: reproduction; 15: reproductive process; 16: response to stimulus; 17: signaling; 18: single-organism process; 19: cell; 20: cell junction; 21: cell part; 22: extracellular matrix; 23: extracellular region; 24: extracellular region part; 25: macromolecular complex; 26: membrane; 27: membrane part; 28: membrane-enclosed lumen; 29: organelle; 30: organelle part; 31: synapse; 32: synapse part; 33: virion; 34: virion part; 35: antioxidant activity; 36: binding; 37: catalytic activity; 38: electron carrier activity; 39: enzyme regulator activity; 40: guanyl-nucleotide exchange factor activity; 41: molecular function regulator; 42: molecular transducer regulator; 43: nucleic acid binding transcription factor activity; 44: structural molecule activity; 45: transporter activity.B: down-regulated genes; 1: biological adhesion; 2: biological regulation; 3: cellular component organization or biogenesis; 4: cellular process; 5: developmental process; 6: localization; 7: locomotion; 8: metabolic process; 9: multicellular organismal process; 10: response to stimulus; 11: signaling; 12: single-organism process; 13: cell; 14: cell part; 15: extracellular matrix; 16: extracellular matrix part; 17: extracellular region; 18: extracellular region part; 19: macromolecular complex; 20: membrane; 21: membrane part; 22: organelle; 23: organelle part; 24: virion; 25: virion part; 26: binding; 27: catalytic activity; 28: enzyme regulator activity; 29: molecular function regulator; 30: molecular transducer regulator; 31: structural molecule activity; 32: transporter activity.

图3 差异表达基因的KEGG pathway富集分析Fig.3 KEGG enrichment analysis of DEGs between AcCK and AcT注: A: 上调基因; 1: 酪氨酸代谢; 2: 糖酵解; 3: 苯丙氨酸代谢; 4: 调节细胞骨架; 5: 半胱氨酸和蛋氨酸代谢; 6: 核糖体; 7: 卵母细胞减数分裂; 8: 苯丙氨酸、酪氨酸、色氨酸代谢; 9: 组氨酸代谢; 10: 脂肪酸退化; 11: 色氨酸代谢; 12: 氨基糖和核苷酸代谢; 13: 酮体的减数分裂; 14: 丁酮代谢; 15: 颉氨酸、亮氨酸、异亮氨酸的减数分裂; 16: 丙氨酸、天冬氨酸、谷氨酸代谢; 17: 柠檬酸循环; 18: 二氧代羧酸代谢; 19: 氨基酸生物合成; 20: 碳代谢.B: 下调基因; 1: 过氧物酶体; 2: 甘油酯代谢; 3: 淀粉与蔗糖的代谢; 4: 氰基氨基酸代谢 5: 磷酸戊糖途径; 6: 鞘糖脂生物合成7: 丙氨酸、天冬氨酸和谷氨酸代谢; 8: 其他葡聚糖降解; 9: 半乳糖代谢、果糖与甘露糖代谢; 10: 叶酸碳库; 11: 戊糖，葡萄糖醛酸转换; 12: 溶酶体; 13: 缬氨酸、亮氨酸和异亮氨酸降解; 14: 粘多糖的降解; 15: 糖酵解和糖异生; 16: 氨基酸生物合成; 17: 甘氨酸; 18: 丝氨酸和苏氨酸代谢; 19: 碳代谢作用; 20: 二羧酸代谢。 Note: A: up-regualted genes; 1: tyrosine metabolism; 2: glycolysis; 3: phenylalanine metabolism; 4: regulation of actin cytoskeleton; 5: cysteine and methionine metabolism; 6: ribosome; 7: oocyte meiosis; 8: phenylalanine, tyrosine and tryptophan biosynthesis; 9: histidine metabolism; 10: fatty acid degradation; 11: tryptophan metabolism; 12: amino sugar and nucleotide sugar metabolism; 13: synthesis and degradation of ketone bodies; 14: butanoate metabolism; 15: valine, leucine and isoleucine degradation; 16: alanine, aspartate and glutamate metabolism; 17: citrate cycle; 18: 2-oxocarboxylic acid metabolism; 19: biosynthesis of amino acids; 20: carbon metabolism.B: down-regulated genes; 1: peroxisome; 2: glycerolipid metabolism; 3: starch and sucrose metabolism; 4: cyanoamino acid metabolism; 5: pentose phosphate pathway; 6: glycosphingolipid biosynthesis-globo series; 7: alanine,aspartate and glutamate metabolism; 8: other glycan degradation; 9: galactose metabolism; 10: fructose and mammose metabolism; 11: one carbon pool by folate; 12: pentose and glucoronate interconversions; 13: lysosome; 14: valine,leucine and isoleucine degradation; 15: glycosaminoglycan degradation; 16: glycolysis; 17: biosynthesis of amino acids; 18: glycine,serine and threonine metabolism; 19: carbon metabolism; 20: glyoxylate and dicarboxylate metabolism.

图4 富集在Toll-like和Jak-STAT信号通路的上调基因的表达量聚类Fig.4 Expression cluster of the up-regulated unigenes enriched in Toll-like and Jak-STAT signaling pathway

前期研究发现球囊菌接种感染后，意蜂幼虫与中蜂幼虫的预蛹(7 d)死亡率分别为70.83%和16.67%，说明后者具有较强的球囊菌抗性，可通过某种机制抵御球囊菌入侵。昆虫中肠内侧有一层由几丁质和蛋白质构成的围食膜，它作为第一道物理屏障能够阻挡经口摄入的病原微生物的入侵(Vuocoloetal., 2001; Pengetal., 1999; Wang and Granados, 2001)。角质层蛋白是构成围食膜的主要成分之一。本研究中，有5个角质层蛋白编码基因表现为下调，说明球囊菌可通过抑制宿主角质层蛋白编码基因的表达来促进侵染。当物理屏障被病原突破后，昆虫随即启动细胞免疫和体液免疫，例如细胞内吞、黑化作用、吞噬作用、酶促降解及分泌抗菌肽等(Gliński and Jarosz, 2001; Glinski and Buczek, 2003)。本研究发现分别有2、17、20、29和52个上调基因富集在凋亡、溶酶体、泛素介导的蛋白水解、吞噬体和内吞作用，仅有1和8个下调基因分别富集在吞噬体和溶酶体，说明中蜂幼虫肠道的细胞免疫在球囊菌胁迫后期被显著激活。同时，我们还发现有1和2个上调基因分别富集在Toll-like及Jak-STAT信号通路上，推测中蜂幼虫的此二条免疫信号通路在抵御球囊菌入侵的过程中发挥关键作用。推测中蜂幼虫的细胞和体液免疫在一定程度上赋予其较强的球囊菌抗性。

目前尚无一种杀真菌剂被批准应用于养蜂生产(Galiziaetal., 2012)，一般通过选育抗病品系、改善养蜂管理和保持清洁卫生来防治白垩病(Gilliac Metal., 1988)，但是效果并不理想。蜜蜂幼虫在球囊菌胁迫过程中免疫应答机制及分子调控机制的缺失严重阻碍白垩病的有效治疗。本研究利用二代测序技术对健康及球囊菌胁迫的中蜂幼虫肠道进行测序，通过DEGs分析并对宿主的胁迫应答进行深入分析，研究结果不仅在转录组水平揭示了中蜂幼虫在球囊菌入侵后期的胁迫应答，也为解析中蜂幼虫的球囊菌抗性机制奠定了基础。

References)

Aizen MA, Garibaldi LA, Cunningh SA,etal.Lessons from long-term trends in crop production [J].AnnalsofBotany, 2009, 103(9): 1579-1588.

Aronstein KA, Murray KD.Chalkbrood disease in honey bees [J].JournalofinvertebratePathology, 2010, 103: S20-S29.

Aronstein KA, Murray KD, Saldivar E.Transcriptional responses in honey bee larvaeinfected with chalkbrood fungus [J].BMCGenomics, 2010, 11: 391.

Bailey L.Infectious Diseases of the Honeybee,1963; Wood M.Microbes help bees battle chalkbrood [M].Agricultural Research Washington D.c, 1998, 46 (8): 16-17.

Begna DHB, Feng M, Li J.Differential expression of nuclear proteomes between honeybee (ApismelliferaL.) queen and worker larvae: a deep insight into caste pathway decisions [J].ProteomeResearch, 2012, 11(2): 1317-1329.

Committee on the Status of Pollinators in North Acerica.Status of pollinators in north Acerica.National Academies Press, 2007.

Cornman RS, Lopez D, Evans JD.Transcriptional response of honey bee larvae infected with the bacterial pathogenPaenibacilluslarvae[J].PLoSOne, 2013,8(6): e65424.

Elke G.Honey bee pathology: current threats to honey bees and beekeeping [J].AppliedMicrobiologyandBiotechnology, 2010, 87(1): 87-97.

Evans JD, Aronstein KA, Chen YP,etal.Immune pathways and defense mechanisms in honey bees,Apismellifera[J].JournalofinvertebratePathology, 2006, 15: 645-656.

Evans JD, Spivak M.Socialized medicine: individual and communaldisease barriers in honey bees [J].JournalofInvertebratePathology, 2010, 103(1): S62-S72.

Foret S, Kucharski R, Pellegrini M ,etal.DNA methylation dynAcics, metabolic fluxes, gene splicing, and alternative phenotypes in honey bees [J].ProceedingsoftheNationalAcademyofSciencesoftheUnitedStatesofAmerica, 2012, 109(13): 4968-4973.

Galizia CG, Eisenhardt D, Giurfa M,etal.Honeybee neurobiology and behavior : a tribute to Randolf Menzel.Dordrecht Netherlands [M].New York: Springer, 2012, 2.

Gliński Z, Jarosz J.Infection and immunity in the honey beeApismellifera[J].Apiacta, 2001, 36: 12-24.

Glinski Z, Buczek K.Response of the Apoidea to fungal infections [J].Apiacta, 2003, 38: 183-189.

Gilliac M, Lii ST, Lorenz BJ,etal.Factors affecting development of chalkbrood disease in colonies of honey bee,Apismellifera, fed pollen contAcinated withAscosphaeraapis[J].JournalofinvertebratePathology, 1988, 52(2): 314-325.

Kurze C, Routtu J, Moritze.Parasite resistance and tolerance in honeybees at the individual and social level [J].Zoology,2016,119(4): 290-7.

KleinAC, Vaissiere BE, Cane JH,etal.Importance of pollinators in changing landscapes for world crops [J].ProceedingsoftheRoyalSocietyB-BiologicalSciences, 2007, 274(1608): 303-313.

Liang Q, Chen DF.Honeybee protection [M].China Agriculture Press, 2008.[梁勤, 陈大福 [M].蜜蜂保护学中国农业出版社, 2008.]

Liang Q, Chen DF, Wang JD.Effects on the mycelia growth and spore-forming ofAscosphaeraapisunder ecological condition of nutrients [J].ChineseJournalofEco-Agriculture, 2001, (9): 31-34.[梁勤, 陈大福, 王建鼎.营养生态条件对蜜蜂球囊菌生长及产孢的影响 [J].中国生态农业学报, 2001, (9): 31-34.]

Park D, Jung JW, Choi BS,etal.Uncovering the novel characteristics of Asian honey bee,Apiscerana, by whole genome sequencing [J].BMCGenomics, 2015, 16(1): 1-16.

Peng JX, Zhong J, Granados RR.A baculovirus enhancin alters the permeability of a mucosal midgut peritrophicmatrix from lepidopteran larvae [J].JournalofInsectPhysiology, 1999, 45(2): 159-166.

Vuocolo T, Eisemann CH, Pearson RD,etal.Identification and molecular characterisation of a peritrophin gene,Peritrophin-48, from the myiasis flyChrysomyabezziana[J].InsectBiochemistryandMolecularBiology, 2001, 31: 919-932.

Wang P, Granados RR.Molecular structure of the peritrophic membrane (PM): Identification of potential PM target sites for insect control [J].ArchivesofInsectBiochemistry&Physiology, 2001, 47(2): 110-118.

Wang Q, Sun LX, Xiao PXetal.Study on Technology for Indoor Artificial Feeding ofApisceranaceranaLarvae [J].ShandongAgricuituralSciences, 2009, 11:113-116.[王倩, 孙亮先, 肖培新, 等.室内人工培育中华蜜蜂幼虫技术研究 [J].山东农业科学, 2009, 11:113-116.]

Xiong CL, Zhang L, Fu ZM,etal.Large-scale development of SSR primers forApisceranaceranalarvae based on its RNA-seq datasets [J].JournalofEnvironmentalEntomology, 2017, 39(1)：68-74.[熊翠玲, 张璐, 付中民, 等.基于RNA-seq数据大规模开发中华蜜蜂幼虫的 SSR 分子标记 [J].环境昆虫学报, 2017, 39(1)：68-74]

Zaghloul OA,Mourad AK, Kady MBEetal.Assessment of losses in honey yield due to the chalkbrood disease, with reference to the determination of its economic injury levels in Egypt [J].CommunicationsinAgricultural&AppliedBiologicalSciences, 2005, 70 (4): 703-714.

Zayed ARG.Understanding the relationship between brain gene expression and social behavior: lessons from the honey bee [J].AnnualReviewofGenetics, 2012, 46(6): 591-615.

Zhang ZN, Xiong CL, Xu XJ,etal.Denovoassembly of a reference transcriptome and development of SSR markers forAscosphaeraapis[J].ActaEntomologicaSinica, 2017, 60 (1):34-44.[张曌楠, 熊翠玲, 徐细建, 等.蜜蜂球囊菌的参考转录组denovo组装及SSR分子标记开发 [J].昆虫学报, 2017, 60(1):34-44.]

Zheng ZY, Li JH, Liang Q.Ascosphaeraapissecretes multiple extracellular enzymes to infect honeybee larvae [J].JournalofFujianAgricultureandForestryUniversity(NaturalScienceEdition), 2011, (40): 280-284.[郑志阳, 李江红, 梁勤.蜜蜂球囊菌分泌多种胞外酶侵染蜜蜂幼虫 [J].福建农林大学学报(自然科学版), 2011, 40(3): 280-284.]

Analysis of the differentially expressed genes in the 6-day-old larval gut ofApisceranaceranaunder the stress ofAscosphaeraapis

GUO Rui1*, ZHANG Lu1*, XU Xi-Jian1, SHI Xiu-Li2, XIONG Cui-Ling1, ZHENG Yan-Zhen1, FU Zhong- Min1, HUANG Zhi-Jian1, WANG Hong-Quan1, HOU Zhi-Xian1, CHEN Da-Fu1**

(1.College of Bee Science, Fujian Agriculture and Forestry University, Fuzhou 350002, China; 2.Apicultural technology management station in the Xinjiang Uygur Autonomous Region, Wulumuqi 930001, China)

Ascosphaeraapisspecially infect honeybee larvae and leads to chalkbrood, which causes a huge loss for apiculture.Here, untreated andA.apis-treatedApisceranaceranalarval guts were sequenced using RNA-seq technology, followed by analysis of differentially expressed genes (DEGs).In this study, Illumina sequencing of larval guts generated a total of 191,1677,30 raw reads, and after filtration, 186,284,296 clean reads were obtained.DEGs analysis showed that 4513 genes were up-regulated and 385 were down-regulated.GO enrichment analysis displayed that the up-regulated genes were enriched in 45 terms, among them the mostly enriched ones were cell process, metabolic process and catalyitic activity; the down-regulated genes were enriched in 32 terms, and the mostly enriched ones were metabolic process, single organism process as well as catalyitic activity.Furthermore, KEGG enrichment analysis suggested that the up-regualted genes were engaged in 193 pathways, the mostly enriched one was ribosome, followed by biosynthesis of acino acids and carbon metabolism; the down-regualted genes were involved in 59 pathways, and the mostly enriched ones were glycine, carbon metabolism as well as glyoxylate and dicarboxylate metabolism.Further analysis showed that the host’s cellular immune was activated and the humoral immunity such as Toll-like and Jak-STAT signaling pathways were significantly induced to activation.Findings in the present study can not only provide key information for uncovering the mechanism regulating theA.c.ceranalarvae’s responses toA.apisduring the late stage of stress, but also lay a foundation for clarifying theA.apis-resistance mechanism ofA.c.ceranalarvae.

Apisceranacerana; larval gut;Ascosphaeraapis; differentially expressed gene; immune defense

郭睿，张璐，徐细建,等.中华蜜蜂6日龄幼虫肠道响应球囊菌胁迫的差异表达基因分析[J].环境昆虫学报，2017,39(3):539-547.

现代农业产业技术体系建设专项资金(CARS-45-KXJ7)；福建农林大学科技发展资金(KF2015123)

*共同第一作者简介：郭睿，男，1987年生，安徽六安人，讲师，研究方向为蜜蜂分子生物学，E-mail: fafu_ruiguo@126.com；张璐，女，1996年生，河南濮阳人，本科生，研究方向为蜂学，E-mail: m17805949180@163.com

**通讯作者Authors for correspondence, E-mail: dfchen826@163.com

Received: 2017-03-04; 接受日期Accepted: 2017-05-02

S891

A

1674-0858(2017)03-0539-09