【目的】通过对意大利蜜蜂（Apis mellifera ligustica，简称意蜂）工蜂中肠响应东方蜜蜂微孢子虫（Nosema ceranae）胁迫的差异表达基因（differentially expressed gene，DEG）及免疫通路进行深入细致的分析，揭示宿主的细胞和体液免疫应答，为关键免疫应答基因的筛选及功能研究打下基础。【方法】基于前期获得的正常及东方蜜蜂微孢子虫胁迫的意蜂7日龄和10日龄工蜂中肠（Am7CK、Am7T、Am10CK、Am10T）转录组数据，按照FDR≤1，P≤0.05和|log2 fold change|≥1的标准筛选出各组的DEG，利用相关生物信息学软件对DEG进行Pearson相关性分析、Venn分析、GO分类和KEGG代谢通路富集分析，进而对免疫通路富集基因进行统计和分析，利用实时荧光定量PCR（real-time quantitative PCR，RT-qPCR）验证转录组数据的可靠性。【结果】 差异分析结果显示，Am7CK vs Am7T比较组包含472个上调基因和385个下调基因；Am10CK vs Am10T比较组包含611个上调基因和360个下调基因。Venn分析结果显示，两个比较组特有的DEG分别为739和853个，共有的DEG为118个。GO分类结果显示，Am7CK vs Am7T比较组中上调和下调基因分别涉及23和29个功能条目，其中富集上调基因数最多的前5位分别为结合、催化活性、代谢进程、细胞进程和单组织进程；富集下调基因数最多的前5位分别为代谢进程、单组织进程、催化活性、细胞进程和结合。Am10CK vs Am10T比较组中上调和下调基因分别涉及36和26个功能条目，其中富集上调基因数最多的前5位分别为单组织进程、结合、细胞进程、催化活性和代谢活性；富集下调基因数最多的前5位分别为结合、细胞进程、催化活性、代谢进程和单组织进程。KEGG代谢通路富集分析结果显示，Am7CK vs Am7T比较组中上调和下调基因分别富集在38和33条代谢通路，富集上调基因数最多的前5条分别为胆汁分泌、内质网蛋白加工、泛素介导的蛋白水解、PI3K-Akt信号通路和神经营养因子信号通路；富集下调基因数最多的前5条分别为胞质DNA传感通路、嘌呤代谢、嘧啶代谢、RNA聚合酶和核糖体；涉及泛素介导的蛋白水解等3条细胞免疫通路，以及PI3K-Akt信号通路等7条体液免疫通路。Am10CK vs Am10T比较组中上调和下调基因分别富集在54和43条代谢通路，富集上调基因数最多的前5条分别为Hippo信号通路、药物代谢-细胞色素P450、P450对外源物质的代谢、泛素介导的蛋白水解和鞘脂类代谢；富集下调基因数最多的前5条分别为mRNA监测、鞘脂类信号通路、果糖和甘露糖代谢、半乳糖代谢和鞘脂类代谢；涉及泛素介导的蛋白水解等7条细胞免疫通路，以及NF-κB信号通路等2条体液免疫通路。RT-qPCR验证结果显示6个随机挑选的DEG的表达量变化趋势与测序结果一致，证实了本研究中测序数据的可靠性。进一步分析发现意蜂7日龄和10日龄工蜂中肠的NF-κB信号通路均被东方蜜蜂微孢子虫激活，随即启动了3种抗菌肽apidaecin、defensin-1和hymenoptaecin的合成，体现了它们在宿主抵御东方蜜蜂微孢子虫入侵中的重要性。【结论】在转录组水平解析了意蜂工蜂中肠对东方蜜蜂微孢子虫的免疫应答，揭示宿主在胁迫前期同时作出细胞和体液免疫应答，前者可能在抵御病原入侵方面发挥主要作用；宿主的细胞免疫在胁迫后期持续增强，但体液免疫较大程度地减弱；泛素介导的蛋白水解通路及富集DEG，NF-κB信号通路及富集DEG，以及apidaecin、defensin-1和hymenoptaecin编码基因可能在意蜂工蜂对东方蜜蜂微孢子虫的免疫应答中起到关键作用。

Honeybees are pivotal pollinators of crops and wild flora, and of great importance in supporting critical ecosystem balance. Nosema ceranae, a unicellular fungal parasite that infects midgut epithelial cells of honeybees, can dramatically reduce honeybee population and productivity. Here, midguts of Apis mellifera ligustica workers at 7 d and 10 d post inoculation (dpi) with sucrose solution (Ac7CK and Ac10CK) and midguts at 7 dpi and 10 dpi with sucrose solution containing N. ceranae spores (Ac7T and Ac10T) were sequenced using strand-specific cDNA library construction and next-generation sequencing. A total of 1956129858 raw reads were gained in this article, and after quality control, 1946489304 high-quality clean reads with a mean Q30 of 93.82% were obtained. The rRNA-removed clean reads were then aligned to the reference genome of Apis mellifera with TopHat2. For more insight please see “Genome-wide identification of long non-coding RNAs and their regulatory networks involved in Apis mellifera ligustica response to Nosema ceranae infection” [1]. Raw data were deposited in NCBI Sequence Read Archive (SRA) database under the BioProject number PRJNA406998. These data can be used for comparative analysis to identify differentially expressed coding RNAs and non-coding RNAs involved in host responses to N. ceranae stress, and for investigation of molecular mechanisms regulating N. ceranae-response.

Here, the expression profiles and differentially expressed miRNAs (DEmiRNAs) in the midguts of Apis cerana cerana workers at 7 d and 10 d post-inoculation (dpi) with N. ceranae were investigated via small RNA sequencing and bioinformatics. Five hundred and twenty nine (529) known miRNAs and 25 novel miRNAs were identified in this study, and the expression of 16 predicted miRNAs was confirmed by Stem-loop RT-PCR. A total of 14 DEmiRNAs were detected in the midgut at 7 dpi, including eight up-regulated and six down-regulated miRNAs, while 12 DEmiRNAs were observed in the midgut at 10 dpi, including nine up-regulated and three downregulated ones. Additionally, five DEmiRNAs were shared, while nine and seven DEmiRNAs were specifically expressed in midguts at 7 dpi and 10 dpi. Gene ontology analysis suggested some DEmiRNAs and corresponding target mRNAs were involved in various functions including immune system processes and response to stimulus. KEGG pathway analysis shed light on the potential functions of some DEmiRNAs in regulating target mRNAs engaged in material and energy metabolisms, cellular immunity and the humoral immune system. Further investigation demonstrated a complex regulation network between DEmiRNAs and their target mRNAs, with miR-598-y, miR-252-y, miR-92-x and miR-3654-y at the center. Our results can facilitate future exploration of the regulatory roles of miRNAs in host responses to N. ceranae, and provide potential candidates for further investigation of the molecular mechanisms underlying eastern honeybeemicrosporidian interactions.

Long non-coding RNAs (lncRNAs) have been proven to play pivotal roles in a wide range of biological processes in animals and plants. Here, we performed whole transcriptome strand-specific RNA sequencing of normal midguts of Apis mellifera ligustica workers (Am7CK, Am10CK) and Nosema ceranae-inoculated midguts (Am7T, Am10T), and 4749 conserved lncRNAs and 1604 novel lncRNAs were identified. These lncRNAs had minimal sequence similarities with other known lncRNAs in other species; however, their structural features were similar to counterparts in mammals and plants, including shorter exon and intron length, lower exon number, and lower expression level, compared with protein-coding transcripts. Further, 111 and 146 N. ceranae-responsive lncRNAs were identified from midguts at 7-days post-inoculation (dpi) and 10 dpi compared with control midguts. Twelve differentially expressed lncRNAs (DElncRNAs) were shared by Am7CK vs. Am7T and Am10CK vs. Am10T comparison groups, while the numbers of unique DElncRNAs were 99 and 134, respectively. Functional annotation and pathway analysis showed that the DElncRNAs may regulate the expression of neighboring genes by acting in cis and trans fashion. Moreover, we discovered 27 lncRNAs harboring eight known miRNA precursors and 513 lncRNAs harboring 2257 novel miRNA precursors. Additionally, hundreds of DElncRNAs and their target miRNAs were found to form complex competitive endogenous RNA (ceRNA) networks, suggesting that these DElncRNAs may act as miRNA sponges. Furthermore, DElncRNA-miRNA-mRNA networks were constructed and investigated, the results demonstrated that a portion of the DElncRNAs were likely to participate in regulating the host material and energy metabolism as well as cellular and humoral immune host responses to N. ceranae invasion. Our findings offer a rich genetic resource for further investigation of the functional roles of lncRNAs involved in the A. m. ligustica response to N. ceranae infection.

为探究中华蜜蜂（Apis cerana cerana）与意大利蜜蜂（Apis mellifera ligustica）幼虫的球囊菌抗性差异产生的原因，利用Venn分析、GO分类分析、KEGG代谢通路分析以及Ka/ks分析对二者的转录组进行比较研究。Venn分析结果显示，中蜂与意蜂的同源基因有17 656 个，归属于8 111 个基因家族，二者的特有基因分别有1 158和241 个，分别归属于458和86 个基因家族。GO分类结果显示，中蜂和意蜂的特有基因分别富集在38和28 个GO条目。KEGG代谢通路富集分析结果显示，中蜂和意蜂的特有基因分别富集于96和21 个代谢通路；进一步分析发现中蜂的特有基因富集在内吞作用、溶酶体、泛素介导的蛋白水解和黑化作用等4 个细胞免疫通路，以及MAPK信号通路和Toll-like受体信号通路等2 个体液免疫通路，而意蜂仅有1 个特有基因富集在MAPK信号通路。Ka/Ks分析结果显示受到强烈正向选择、弱正向选择、中性选择和负向选择的单拷贝同源基因分别有37、281、2 630和3 269 个。对受到强烈正向选择的基因进行GO分类和KEGG代谢通路富集分析，GO结果显示中蜂和意蜂的单拷贝同源基因富集的GO 条目类型相同；KEGG结果显示中蜂的单拷贝基因仅富集在氧化磷酸化。上述结果表明免疫相关基因数量的差异是二者的球囊菌抗性差异产生的重要原因之一，中蜂在与球囊菌的协同进化过程中可能通过提高能量利用从而限制球囊菌的增殖。本研究结果可为阐明中蜂及意蜂幼虫的球囊菌抗性差异产生的分子机制提供参考。