侵染中华蜜蜂6日龄幼虫的蜜蜂球囊菌的微小RNA差异表达谱及调控网络

【目的】蜜蜂球囊菌(Ascosphaeraapis,简称球囊菌)是专性侵染蜜蜂幼虫的致死性真菌病原。MicroRNA(miRNA)作为一类重要的基因表达调控因子,能够广泛参与真菌及其宿主的相互作用过程。本研究通过比较分析球囊菌孢子(AaCK)和侵染中华蜜蜂(Apis cerana cerana,简称中蜂) 6日龄幼虫肠道内的球囊菌(AaT)的smallRNA(sRNA)组学数据对球囊菌的差异表达miRNA(differentiallyexpressed miRNA,DEmiRNA)、靶mRNA及二者间的调控网络进行全面解析,旨在揭示miRNA介导的球囊菌对中蜂幼虫的侵染机制。【方法】对于球囊菌侵染的中蜂6日龄幼虫肠道的small RNA-seq (sRNA-seq)数据,利用BLAST工具连续比对东方蜜蜂(Apiscerana)和球囊菌的参考基因组筛滤得到AaT的sRNA组学数据。分别将AaCK和AaT的sRNA组学数据比对miRBase数据库,对球囊菌侵染宿主前后miRNA的数量和结构特征进行分析。联用RNAhybrid+svm_light、Miranda和TargetScan软件预测AaCK vs AaT比较组中DEmiRNA的靶mRNA,进而利用相关生物信息学软件对上述靶mRNA进行GO分类和KEGG代谢通路富集分析。通过Cytoscape软件对DEmiRNA-mRNA调控网络进行可视化。利用Stem-loop RT-PCR、RT-qPCR和分子克隆验证测序结果的可靠性。【结果】在AaCK和AaT中分别鉴定到380和387个miRNA。结构特征分析结果显示,AaCK和AaT的mi RNA皆集中分布在18–25 nt,且首位碱基主要偏向于U。AaCKvsAaT比较组共有270个DEmiRNA,包含155个上调miRNA和115个下调miRNA,分别靶向结合6091和6145个mRNA。GO分类结果显示,上述靶mRNA主要涉及代谢进程、细胞进程、应激反应等15个生物学进程;细胞、细胞组分、细胞器等12个细胞组分;催化活性、结合、转运子活性等11个分子功能。KEGG代谢通路富集分析结果显示,上述靶mRNA富集在123条代谢通路,参与对氨基酸代谢、碳水化合物代谢以及核苷酸代谢等物质代谢,氧化磷酸化、硫代谢、氮代谢等能量代谢,以及MAPK和Hippo等信号通路的调控。球囊菌DEmiRNA与靶mRNA之间存在复杂的调控关系,其中miR-29-x、miR-250-x、miR-4968-y、miR-11200-x、novel-m0023-5p、novel-m0130-5p和novel-m0135-5p等DEmiRNA可靶向结合与球囊菌的半胱氨酸蛋白酶、DNA甲基化转移酶以及几丁质酶相关的mRNA;此外,miR-7-x、miR-9-z、miR-319-y和miR-5951-y等同时参与调控MAPK信号通路;进一步分析发现,miR-250-x同时参与对DNA甲基化转移酶、MAPK信号通路及其他酶类合成与代谢途径的调控,并可能参与球囊菌与中蜂6日龄幼虫之间的跨界调控。通过Stem-loopRT-PCR和RT-qPCR验证了4个DEmiRNA的差异表达,并利用分子克隆和Sanger测序证实miR-7-x的序列与测序结果一致。【结论】本研究解析了侵染中蜂6日龄幼虫的球囊菌的miRNA差异表达谱及DEmiRNA的调控网络,揭示了球囊菌DEmiRNA可能通过调控病原的物质和能量代谢、增殖、毒力、信号通路及相关mRNA参与对中蜂幼虫的侵染过程。miR-7-x、miR-250-x、novel-m0023-5p等关键DEmiRNA有望作为白垩病治疗的新型分子靶点。     

参与调控意大利蜜蜂工蜂中肠基因表达的东方蜜蜂微孢子虫miRNA的组学解析及其调控网络

【目的】本研究结合前期已获得的miRNA和mRNA组学数据对东方蜜蜂微孢子虫Nosema ceranae的差异表达miRNA(differentially expressed miRNA, DEmiRNA)靶向意大利蜜蜂Apis mellifera ligustica工蜂中肠的mRNA和差异表达mRNA(differentially expressed mRNA, DEmRNA)进行生物信息学预测、数据库注释和调控网络分析,旨在解析东方蜜蜂微孢子虫对意大利蜜蜂工蜂中肠的基因表达调控。【方法】比较东方蜜蜂微孢子虫感染7 d(AmT1)和10 d(AmT2)的意大利蜜蜂工蜂中肠和未感染中肠(分别为AmCK1和AmCK2)的mRNA数据,筛选意大利蜜蜂工蜂的显著性DEmRNA;通过比较侵染AmT1和AmT2的东方蜜蜂微孢子虫(分别为NcT1和NcT2)和东方蜜蜂微孢子虫纯净孢子(NcCK)的miRNA数据筛选出东方蜜蜂微孢子虫的DEmiRNA。利用TargetFinder软件预测东方蜜蜂微孢子虫DEmiRNA靶向的意大利蜜蜂工蜂中肠DEmRNA。利用相关生物信息学工具对上述意大利蜜蜂工蜂中肠靶DEmRNA进行GO和KEGG数据库注释。根据KEGG数据库注释信息筛选出意大利蜜蜂工蜂中肠免疫防御和能量代谢通路相关DEmRNA,并构建和分析上述DEmRNA与相应的东方蜜蜂微孢子虫DEmiRNA之间的调控网络。【结果】NcCK vs NcT1比较组中东方蜜蜂微孢子虫的77条显著上调miRNA和52条显著下调miRNA可分别靶向AmCK1 vs AmT1比较组中意大利蜜蜂工蜂中肠的118条显著下调mRNA和135条显著上调mRNA,这些mRNA可分别注释到31和25个GO条目,以及113和107条KEGG通路。NcCK vs NcT2比较组中东方蜜蜂微孢子虫的52条显著上调miRNA和49条显著下调miRNA可分别靶向AmCK2 vs AmT2比较组中意大利蜜蜂工蜂中肠的97条显著下调mRNA和210条显著上调mRNA,这些mRNA可分别注释到27和30个GO条目以及97和127条KEGG通路。NcCK vs NcT1和NcCK vs NcT2比较组中的11条共同显著上调miRNA和19条共同显著下调miRNA分别靶向AmCK1 vs AmT1和AmCK2 vs AmT2比较组的6条共同显著下调和14条共同显著上调mRNA,可分别注释到7和10个GO条目以及0和9条KEGG通路。NcCK vs NcT1和NcCK vs NcT2比较组中东方蜜蜂微孢子虫的DEmiRNA可靶向AmCK1 vs AmT1和AmCK2 vs AmT2比较组中意大利蜜蜂工蜂中肠的氧化磷酸化和硫代谢等能量代谢通路相关DEmRNA,以及胞吞作用、黑色素生成、溶酶体、自噬、Toll样受体信号通路、细胞凋亡、Ras信号通路、泛素介导的蛋白水解和MAPK信号通路等免疫防御通路相关DEmRNA。进一步分析发现,miR-216-x, miR-5119-y, bantam-y和miR-8-y在NcCK vs NcT1和NcCK vs NcT2比较组中皆显著上调表达,且靶向意大利蜜蜂工蜂中肠的溶酶体、黑色素生成、泛素介导的蛋白水解、MAPK信号通路及Ras信号通路等免疫防御通路相关的显著下调mRNA。【结论】在侵染过程中,东方蜜蜂微孢子虫的DEmiRNA对意大利蜜蜂工蜂的基因表达具有广泛的潜在影响,东方蜜蜂微孢子虫可能通过上调部分miRNA跨界调控意大利蜜蜂工蜂的免疫防御以促进侵染,通过下调部分miRNA跨界调控意大利蜜蜂工蜂的能量代谢以加强能量窃取并促进增殖。     

蜜蜂球囊菌菌丝和孢子中环状RNA的鉴定及比较分析

[目的] 本研究旨在明确蜜蜂球囊菌（Ascosphaera apis,简称球囊菌）菌丝和孢子中环状RNA（circular RNA,circRNA）的数量、种类和表达谱差异,并探讨共有circRNA、特有circRNA和差异表达circRNA（differentially expressed circRNA,DEcircRNA）在菌丝与孢子中的潜在作用。[方法] 基于前期获得的球囊菌菌丝（AaM）和孢子（AaS）的高质量RNA-seq数据,利用find_circ软件预测circRNA。通过Venn分析筛选AaM和AaS的共有circRNA和特有circRNA。根据P ≤ 0.05且|log₂ fold change|≥ 1的标准筛选AaMvs.AaS比较组的DEcircRNA。通过比对GO和KEGG数据库对circRNA的来源基因进行功能和通路注释。利用TargetFinder软件预测circRNA靶向结合的miRNA及miRNA靶向结合的mRNA。采用Cytoscape软件对竞争性内源RNA（competing endogenous RNA,ceRNA）调控网络进行可视化。通过RT-qPCR对DEcircRNA进行验证。[结果] AaM和AaS中分别含有13210156和19011000条短序列读段（anchors reads）,其中分别有6124922和11392886条能够比对上球囊菌参考基因组。在AaM和AaS中分别鉴定到1868个和2225个circRNA,二者共有的circRNA为1098个,AaM的特有circRNA为770个,AaS的特有circRNA为1127个。AaM和AaS的circRNA长度主要介于1000-2000 nt,基因间区circRNA为主要环化类型。AaMvs.AaS比较组包含456个上调circRNA和97个下调circRNA。共有circRNA的来源基因可注释到29个功能条目和14类通路;AaM的特有circRNA的来源基因可注释到31个功能和17类通路;AaS的特有circRNA的来源基因可注释到34个功能条目和16类通路;DEcircRNA的来源基因可注释到29个功能条目和40条通路。调控网络分析结果显示,36个共有circRNA靶向4个miRNA进而调控6个与内吞作用通路相关的靶mRNA;4（255）个AaM（AaS）的特有circRNA靶向2（2）个miRNA进而调控8（2）个次级代谢产物生物合成通路相关的靶mRNA;9个DEcircRNA靶向2个DEmiRNA进而调控3个MAPK信号通路相关的DEmRNA。RT-qPCR结果显示10个DEcircRNA的表达趋势与测序数据一致,证实了测序数据的可靠性。[结论] 球囊菌菌丝和孢子的共有circRNA、特有circRNA和DEcircRNA可能通过调控来源基因表达和充当ceRNA的方式调节球囊菌的物质和能量代谢、内吞作用、次级代谢产物生物合成和MAPK信号通路,进而影响球囊菌菌丝生长、孢子萌发和致病性。     

转录组分析揭示东方蜜蜂微孢子虫侵染意大利蜜蜂的分子机制

【目的】本研究旨在通过差异表达基因(differentially expressed gene, DEG)分析以及毒力因子和其他侵染相关因子分析,在转录组水平揭示东方蜜蜂微孢子虫Nosema ceranae侵染意大利蜜蜂Apis mellifera ligustica的分子机制。【方法】基于前期已获得高质量的东方蜜蜂微孢子虫纯化孢子(NcCK)及侵染意大利蜜蜂工蜂7和10 d的东方蜜蜂微孢子虫(分别为NcT1和NcT2)转录组数据,根据P≤0.05且|log_2(Fold change)|≥1的标准,通过比较分析筛选出NcCK vs NcT1, NcCK vs NcT2和NcT1 vs NcT2比较组的DEG。通过相关生物信息学软件对上述DEG进行Venn分析、GO分类和KEGG代谢通路富集分析。根据Nr和KEGG数据库注释信息和相关文献进行对东方蜜蜂微孢子虫的毒力因子和侵染相关因子的统计和分析。通过RT-qPCR验证转录组数据及DEG表达趋势。【结果】从NcCK vs NcT1, NcCK vs NcT2和NcT1 vs NcT2比较组分别鉴定出1 397, 1 497和52个DEG。Venn分析结果显示各比较组共有的上调和下调基因分别为10和1个。GO分类结果显示,NcCK vs NcT1和NcCK vs NcT2中DEG富集数最多的功能条目为代谢进程、细胞进程、单组织进程、细胞、细胞组件、细胞器、催化活性和结合,而NcT1 vs NcT2中DEG富集数最多的是代谢进程、细胞进程、单组织进程、催化活性和结合。KEGG代谢通路富集分析结果显示,NcCK vs NcT1和NcCK vs NcT2中DEG分别富集到80和79条通路;富集在糖酵解/糖异生和MAPK信号通路的上调基因数量多于下调基因。毒力因子分析结果显示,孢壁蛋白9基因和孢壁蛋白12基因在NcCK vs NcT1和NcCK vs NcT2中均下调表达,孢壁蛋白8基因仅在NcCK vs NcT1中表达量下调;此外孢壁蛋白前体基因、孢壁和锚定盘复合蛋白基因、几丁质合酶基因、极管蛋白基因、蓖麻毒素B凝集素基因的表达水平在NcCK vs NcT1和NcCK vs NcT2中表现为上调。侵染相关因子分析结果表明,糖酵解途径的3个关键酶基因在NcCK vs NcT1和NcCK vs NcT2中上调表达;3个涉及ATP/ADP移位酶的基因在NcCK vs NcT1和NcCK vs NcT2中上调表达,但有1个表达量下调;2个涉及ABC转运蛋白的基因在NcCK vs NcT1和NcCK vs NcT2中上调表达,另有4个下调表达。RT-qPCR结果证实了本研究中转录组数据及DEG表达趋势的真实可靠性。【结论】本研究通过比较分析解析东方蜜蜂微孢子虫侵染意大利蜜蜂工蜂过程的转录组动态,揭示了孢壁蛋白、孢壁和锚定盘复合蛋白、几丁质酶、极管蛋白和蓖麻毒素B凝集素等毒力因子编码基因,以及己糖激酶、丙酮酸激酶、6-磷酸果糖激酶、ATP/ADP移位酶和ABC转运蛋白等侵染相关因子编码基因在病原增殖中扮演重要角色,为阐明东方蜜蜂微孢子虫的侵染机制提供了基础。     

基于small RNA组学分析揭示意大利蜜蜂响应东方蜜蜂微孢子虫胁迫的免疫应答机制

【目的】通过深度测序和组学分析在small RNA组学层面揭示意大利蜜蜂(Apis mellifera ligustica,简称意蜂)响应东方蜜蜂微孢子虫(Nosema ceranae)的免疫应答机制。【方法】利用small RNA-seq技术对正常及N. ceranae胁迫7 d和10 d的意蜂工蜂中肠(Am7CK、Am7T、Am10CK和Am10T)进行深度测序。利用相关生物信息学软件对测序数据进行质控、已知micro RNA (mi RNA)鉴定、新mi RNA预测及mi RNA的结构特征分析。通过Stem-loop RT-PCR验证新mi RNA的表达。按照|log2(Fold change)|≥1和P≤0.05的标准筛选出Am7CKvs.Am7T和Am10CKvs.Am10T比较组的差异表达mi RNA(differentially expressed mi RNA,DEmi RNA)。利用软件预测DEmi RNA靶向结合的m RNA并进行GO和KEGG数据库注释,根据注释信息对细胞和体液免疫相关通路及富集靶m RNA进行统计和分析。根据靶向结合关系构建DEmi RNA和免疫通路相关差异表达m RNA (differentially expressed m RNA,DEm RNA)的调控网络。采用RT-q PCR对数据的可靠性和DEmi RNA的差异表达进行验证。【结果】共获得165895574条原始读段和132028990条有效序列标签,各组的组内Pearson相关性平均在87.92%及以上。共鉴定到928个已知mi RNA和56个新mi RNA。这些mi RNA的长度介于18–28nt,多数的长度为18nt和22nt且首位碱基主要偏向U。验证了12个新mi RNA的真实表达。Am7CKvs.Am7T比较组包含48个上调mi RNA和36个下调mi RNA;Am10CK vs. Am10T比较组包含56个上调mi RNA和51个下调mi RNA。两个比较组的DEmi RNA可分别靶向结合9827个和10720个m RNA。这些靶m RNA可分别注释到50和47条功能条目,以及138和135条KEGG通路。DEmi RNA与免疫通路相关靶m RNA的调控网络分析结果显示,Am7CK vs. Am7T中有26个DEmi RNA靶向与内吞作用等免疫通路相关的10个DEm RNA;Am10CK vs. Am10T中有15个DEmi RNA靶向与MAPK信号通路等免疫通路相关的10个DEm RNA。验证了测序数据和4个DEmi RNA差异表达的可靠性。【结论】研究结果揭示了宿主DEmi RNA可能通过调控物质和能量代谢、细胞和体液免疫对N. ceranae产生应答,但DEmi RNA不参与抗菌肽基因的表达调控;mi R-1-z可能参与宿主的细胞增殖、细胞凋亡和免疫进程;氧化磷酸化通路可能在宿主免疫应答及宿主-病原互作中发挥特殊作用。     

中华蜜蜂6日龄幼虫响应蜜蜂球囊菌胁迫的环状RNA应答

【目的】蜜蜂球囊菌(Ascosphaera apis,简称球囊菌)专性侵染蜜蜂幼虫而导致白垩病,危害蜜蜂健康和养蜂生产。本研究旨在探究中华蜜蜂(Apis cerana cerana,简称中蜂) 6日龄幼虫响应球囊菌胁迫的环状RNA(circular RNA,circRNA)差异表达谱及差异表达circRNA (differentially expressed circRNA,DEcircRNA)在宿主胁迫应答中的潜在功能。【方法】利用去除线性RNA的circRNA-seq技术对正常和球囊菌侵染的中蜂6日龄幼虫肠道(AcCK和AcT)进行测序。利用find_circ软件鉴定circRNA,统计circRNA的长度和环化类型。根据|log_2(Fold change)|≥1和P≤0.05的标准筛选DEcircRNA。将DEcircRNA的来源基因比对Gene ontology (GO)数据库和Kyoto Encyclopedia of Genes and Genomes (KEGG)数据库,从而获得功能及通路(pathway)注释。随机挑选3个DEcircRNA进行RT-qPCR验证。【结果】AcCK和AcT的circRNA-seq分别得到76342570和68269362条原始读段(raw reads),经严格质控得到74524108和66974392条有效读段(clean reads),Q30分别为92.75%和94%,GC含量分别为54.31%和54.90%。比对上东方蜜蜂(Apis cerana)参考基因组的短序列读段(anchor reads)共计23648400条。AcCK和AcT中分别鉴定到805和702个circRNA,长度均介于201–1000 nt,数量最多的环化类型均为已注释外显子circRNA,但分布在不同长度、不同环化类型的circRNA数量存在差异。AcCK vs AcT比较组共有494个DEcircRNA,包括257个上调circRNA和237个下调circRNA;上调和下调幅度最大的circRNA分别为novel_circ_000123和novel_circ_000726。上述DEcircRNA的来源基因可注释到11条生物学进程相关条目,9条分子功能相关条目,9条细胞组分相关条目,以及73条通路。进一步分析发现,部分DEcircRNA的来源基因注释到7条细胞免疫通路和3条体液免疫通路。【结论】中蜂6日龄幼虫响应球囊菌胁迫的过程中可能通过改变分布在不同长度和环化类型的circRNA数量,以及特异性表达一些circRNA和调节部分circRNA的表达量对病原产生应答;novel_circ_000027、novel_circ_000127、novel_circ_000312等DEcircRNA在宿主的胁迫应答过程中可能通过调控氧化磷酸化、细胞和体液免疫等通路发挥特殊作用。研究结果为深入理解中蜂幼虫对球囊菌的胁迫应答机制及二者的相互作用机制提供了新见解。     

东方蜜蜂微孢子虫基因的可变剪接及可变腺苷酸化解析

本研究旨在基于已获得的第三代纳米孔全长转录组数据对东方蜜蜂微孢子虫Nosema ceranae基因的可变剪接（alternative splicing,AS）和可变多聚腺苷酸化（alternative polyadenylation,APA）进行分析。通过Astalavista软件鉴定东方蜜蜂微孢子虫基因的AS事件类型。采用TAPIS pipeline对东方蜜蜂微孢子虫基因的APA位点进行鉴定。利用MEME软件分析所有全长转录本的poly (A)剪接位点上游50bp的序列特征并鉴定motif。共鉴定出发生于5个基因的5次AS事件,包括1次可变供体位点（alternative donor site,ADS）和4次内含子保留（intron retention,IR）。鉴定出233个基因发生APA,其中含有1个poly (A)剪接位点的基因数量最多（143,61.37%）。序列特征分析结果显示,东方蜜蜂微孢子虫全长转录本的3’ UTR的上下游表现出明显的碱基倾向性,U在3’ UTR的上游富集,而A在3’ UTR的下游富集。此外,在东方蜜蜂微孢子虫全长转录本的poly (A)剪接位点上游鉴定到3个motif（AAUAAA、UGAUGC和GCGACG）。研究结果揭示了东方蜜蜂微孢子虫转录组的复杂性,完善了现有的参考基因组和转录组注释,也为深入探究不同剪接异构体的分子功能提供了基础。此外,本研究为其他真菌的AS和APA研究提供了有益的思路和方法借鉴。     

基于第三代长读段测序数据解析蜜蜂球囊菌基因的可变剪切与可变腺苷酸化

【目的】蜜蜂球囊菌(Ascosphaera apis,简称球囊菌)是一种专性侵染蜜蜂幼虫的真菌病原,导致的白垩病是严重影响养蜂生产的顽疾,每年给养蜂业造成较大损失。本研究旨在基于已获得的第三代长读段测序数据对球囊菌菌丝(Aam)和孢子(Aas)中基因的可变剪切(alternative splicing,AS)和可变多聚腺苷酸化(alternative polyadenylation,APA)进行深入分析。【方法】利用Astalavista软件鉴定Aam和Aas中基因的AS事件类型。利用IGV浏览器对部分剪切异构体(isoform)的结构进行可视化。通过TAPIS pipeline对Aam和Aas中基因的APA位点进行鉴定。通过MEME软件对Aam和Aas中全长转录本的APA位点上游50bp的序列特征进行分析并对motif进行鉴定。【结果】在Aam中共鉴定到286次AS事件,包括162次RI(Retained intron),87次A3(Alternative 3'splice-site)、32次A5(Alternative 5'splice-site)和5次SE (Skipping exon);在Aas中共鉴定到559次AS事件,包括305次RI、155次A3、85次A5、13次SE和1次MEE (Mutually exclusive exon)。进一步分析发现,现有参考基因组上的多数注释基因结构并不完整;部分注释基因在菌丝和孢子中转录形成的isoform在数量和结构方面均存在差异;部分isoform在参考基因组中没有对应的注释基因。Aam中共鉴定到2748个基因含有1个及以上的APA位点,其中含有1个APA位点的基因数量最多(726,26.42%);Aas中共鉴定到2768个基因含有1个及以上的APA位点,其中含有5个以上APA位点的基因数量最多(1180,42.63%)。部分基因在菌丝和孢子中含有不同的APA位点数。序列特征分析结果显示,球囊菌全长转录本的3'UTR的上下游表现出明显的碱基倾向性,U和A分别富集在3'UTR的上游和下游。此外,在球囊菌全长转录本的APA位点上游鉴定到4个motif,分别是UCUCCU、UCUUCU、CCCACC和CCCCCU。【结论】本研究通过对球囊菌菌丝和孢子中基因的AS和APA进行深入分析,揭示了球囊菌转录组的复杂性,为完善现有的基因组和转录组注释提供了宝贵信息,也为探究AS和APA在球囊菌的基因表达调控中的作用提供了关键基础。