參與調(diào)控意大利蜜蜂工蜂中腸基因表達的東方蜜蜂微孢子蟲miRNA的組學(xué)解析及其調(diào)控網(wǎng)絡(luò)

范小雪, 杜 宇, 張文德, 王 杰, 蔣海賓, 范元嬋, 馮睿蓉,萬潔琦, 周紫彧, 熊翠玲,2, 鄭燕珍,2, 陳大福,2,3, 郭 睿,2,3,*

(1. 福建農(nóng)林大學(xué)動物科學(xué)學(xué)院(蜂學(xué)學(xué)院), 福州 350002; 2. 福建農(nóng)林大學(xué)蜂療研究所, 福州 350002;3. 福建農(nóng)林大學(xué)蜂產(chǎn)品加工與應(yīng)用教育部工程研究中心, 福州 350002)

東方蜜蜂微孢子蟲Nosemaceranae是一種專性侵染成年蜜蜂中腸上皮細胞的單細胞真菌病原，可對蜜蜂宿主造成腸道破壞、能量脅迫、免疫抑制、消化紊亂及壽命縮短等不良影響(Higesetal., 2007; Mayack and Naug, 2009; Alauxetal., 2010; Botíasetal., 2013; Chenetal., 2013)，但到目前為止二者互作的分子機制尚未闡明。

微小RNA(microRNA, miRNA)是一類長度約為18～25 nt且高度保守的內(nèi)源性單鏈非編碼RNA(non-coding RNA, ncRNA)(Bartel, 2004)，可通過靶向結(jié)合mRNA使其降解或抑制其翻譯，從而在轉(zhuǎn)錄后水平調(diào)控基因表達，進而發(fā)揮廣泛而重要的生物學(xué)功能(Pillai, 2005)。此外，miRNA已被證實能夠作為媒介介導(dǎo)真菌病原及其宿主的互作(Linetal., 2016; Crostonetal., 2018)。2012年，張辰宇團隊研究發(fā)現(xiàn)植物來源的miR-168a可在動物體內(nèi)穩(wěn)定積累表達，并通過靶向抑制小鼠低密度脂蛋白LDLRAP1的表達水平，降低小鼠血漿LDL的去除率，首次證實miRNA的跨界調(diào)控作用(Zhangetal., 2012)。較多的研究證實結(jié)果表明小RNA(small RNA, sRNA)可在真菌病原及其宿主之間雙向傳播，真菌病原可將sRNA作為效應(yīng)因子轉(zhuǎn)移到宿主細胞內(nèi)，參與宿主的基因表達調(diào)控，從而影響自身毒力、宿主免疫力以及疾病進程等。例如，灰霉菌Botrytiscinerea和球孢白僵菌Beauveriabassiana來源的sRNA在細胞囊泡的包裹和保護作用下轉(zhuǎn)移到宿主體內(nèi)，通過與宿主Argonaute 1蛋白結(jié)合并劫持宿主RNAi機制，進而選擇性沉默宿主免疫基因，以達到有效感染的目的(Weibergetal., 2013; Cuietal., 2019)。對于miRNA介導(dǎo)的西方蜜蜂Apismellifera及東方蜜蜂微孢子蟲互作，前人開展了少量研究(Evans and Huang, 2018; Huangetal., 2019)，但總體進展緩慢。前期研究中，筆者所在團隊聯(lián)用鏈特異性建庫的RNA-seq和small RNA-seq(sRNA-seq)技術(shù)對東方蜜蜂微孢子蟲感染7 d和10 d的意大利蜜蜂Apismelliferaligustica工蜂中腸及未感染中腸進行測序，基于高質(zhì)量的全轉(zhuǎn)錄組數(shù)據(jù)系統(tǒng)解析了lncRNA, miRNA及ceRNA調(diào)控網(wǎng)絡(luò)介導(dǎo)的宿主免疫應(yīng)答(Chenetal., 2019; 付中民等, 2019)；此外，還在lncRNA, miRNA和mRNA組學(xué)水平深入分析和探討了東方蜜蜂微孢子蟲侵染意大利蜜蜂工蜂的分子機制(耿四海等, 2020a, 2020b)。近期，筆者所在團隊基于已獲得的東方蜜蜂微孢子蟲感染及未感染的意大利蜜蜂工蜂中腸轉(zhuǎn)錄組數(shù)據(jù)(付中民等, 2019)以及純凈的病原孢子全轉(zhuǎn)錄組數(shù)據(jù)(耿四海等, 2020a)篩選出宿主的差異表達miRNA(differentially expressed miRNA, DEmiRNA)和病原的差異表達mRNA(differentially expressed mRNA, DEmRNA)，通過生物信息學(xué)手段預(yù)測宿主DEmiRNA和病原DEmRNA的靶向關(guān)系，并對二者間的調(diào)控網(wǎng)絡(luò)進行了構(gòu)建和分析，進一步結(jié)合本團隊的前期研究結(jié)果和前人的相關(guān)研究結(jié)果探討了意大利蜜蜂工蜂的DEmiRNA對東方蜜蜂微孢子蟲基因和通路的調(diào)控作用(未發(fā)表數(shù)據(jù))。在長期的協(xié)同進化中，東方蜜蜂微孢子蟲的基因組大大簡化，多數(shù)物質(zhì)和能量代謝的途徑丟失，其增殖高度依賴宿主細胞提供物質(zhì)和能量(Burrietal., 2006)。同時東方蜜蜂微孢子蟲也進化出一些特殊的能力，能夠?qū)λ拗鞯纳砩δ苓M行操縱，如東方蜜蜂微孢子蟲侵染可對西方蜜蜂工蜂造成較強的能量脅迫，致使宿主的能量代謝速率加快且糖水的攝入量提高(Antúnezetal., 2009; Chenetal., 2013)。然而，對于東方蜜蜂微孢子蟲通過何種方式和機制對西方蜜蜂的基因表達和生理活動進行調(diào)控，相關(guān)研究仍然匱乏，迄今僅有一例報道(Evans and Huang, 2018)。此前，Evans和Huang(2018)通過高通量測序和生物信息學(xué)分析方法在東方蜜蜂微孢子蟲感染的西方蜜蜂工蜂中腸中預(yù)測出病原的6條miRNA，進一步分析發(fā)現(xiàn)其中5條miRNA不僅能靶向病原本身的1 545條mRNA，還可能被分泌到宿主細胞的胞質(zhì)內(nèi)，通過靶向918條宿主mRNA調(diào)控宿主的新陳代謝、細胞凋亡和免疫防御等過程。

本研究擬通過比較分析篩選出侵染意大利蜜蜂工蜂的東方蜜蜂微孢子蟲的DEmiRNA以及意大利蜜蜂工蜂響應(yīng)東方蜜蜂微孢子蟲侵染的意大利蜜蜂DEmRNA，通過生物信息學(xué)方法預(yù)測病原DEmiRNA靶向的宿主mRNA和DEmRNA，并通過注釋數(shù)據(jù)庫獲得相應(yīng)的功能和通路信息，進一步結(jié)合東方蜜蜂微孢子蟲侵染西方蜜蜂的生物學(xué)背景，筆者所在課題組的前期研究結(jié)果以及前人的相關(guān)研究結(jié)果深入探討東方蜜蜂微孢子蟲的DEmiRNA對意大利蜜蜂工蜂中腸基因表達的調(diào)控作用。研究結(jié)果可為探明東方蜜蜂微孢子蟲的侵染機制、東方蜜蜂微孢子蟲-意大利蜜蜂工蜂互作機制提供重要的參考信息和基礎(chǔ)。

1 材料與方法

1.1 東方蜜蜂微孢子蟲的miRNA組學(xué)數(shù)據(jù)來源

前期研究中，筆者所在課題組分別提取了3個東方蜜蜂微孢子蟲的純凈孢子樣品NcCK(NcCK1, NcCK2, NcCK3)的總RNA，并按照已建立的方法(耿四海等, 2020a)構(gòu)建了相應(yīng)的cDNA文庫，單端測序由廣州基迪奧生物科技有限公司完成，測序平臺為Illumina MiSeq。通過連續(xù)比對西方蜜蜂參考基因組(Assembly Amel_4.5)、GenBank數(shù)據(jù)庫、Rfam數(shù)據(jù)庫及東方蜜蜂微孢子蟲參考基因組(Assembly ASM98816v1)篩濾得到比對上東方蜜蜂微孢子蟲參考基因組的東方蜜蜂微孢子蟲純凈孢子miRNA組學(xué)數(shù)據(jù)(NcCK)，侵染意大利蜜蜂工蜂7 d(NcT1)和10 d(NcT2)的東方蜜蜂微孢子蟲miRNA組學(xué)數(shù)據(jù)(耿四海等, 2020a)。測序數(shù)據(jù)已上傳至NCBI SRA數(shù)據(jù)庫，NcCK的Bioproject登記號為PRJNA395264。通過對下機的原始數(shù)據(jù)進行質(zhì)控得到了最終的有效標簽序列(clean tags)(耿四海等, 2020a)。

1.2 東方蜜蜂微孢子蟲的DEmiRNA篩選

按照耿四海等(2020a)的方法進行東方蜜蜂微孢子蟲miRNA的鑒定和表達量計算以及NcCKvsNcT1和NcCKvsNcT2比較組顯著性DEmiRNA的篩選，篩選出的顯著性DEmiRNA(|log2fold change|≥1且P≤0.05)用于本研究中靶向意大利蜜蜂工蜂中腸mRNA和DEmRNA的預(yù)測和分析。

1.3 意大利蜜蜂工蜂中腸的mRNA組學(xué)數(shù)據(jù)來源

筆者所在團隊前期已對接種50%(w/v)蔗糖溶液7 d(AmCK1)和10 d(AmCK2)的意大利蜜蜂工蜂中腸樣品以及接種含東方蜜蜂微孢子蟲孢子的50%(w/v)蔗糖溶液7 d (AmT1)和10 d(AmT2)的意大利蜜蜂工蜂中腸樣品進行RNA抽提、cDNA文庫構(gòu)建、轉(zhuǎn)錄組測序以及測序數(shù)據(jù)質(zhì)控(付中民等, 2019)。測序數(shù)據(jù)已上傳至NCBI SRA數(shù)據(jù)庫，Bioproject號為PRJNA395264。質(zhì)控后的意大利蜜蜂工蜂中腸mRNA組學(xué)數(shù)據(jù)可用于本研究中東方蜜蜂微孢子蟲的顯著性DEmiRNA的靶標預(yù)測及分析。

1.4 東方蜜蜂微孢子蟲顯著性DEmiRNA靶向意大利蜜蜂工蜂中腸mRNA的預(yù)測

為探究1.2節(jié)篩選出的東方蜜蜂微孢子蟲DEmiRNA調(diào)控意大利蜜蜂工蜂中腸基因表達的廣泛性，利用TargetFinder軟件(Allenetal., 2005)預(yù)測NcCKvsNcT1和NcCKvsNcT2比較組中顯著性DEmiRNA靶向的1.3節(jié)的意大利蜜蜂工蜂中腸的全部mRNA，采用默認參數(shù)。

1.5 意大利蜜蜂工蜂中腸響應(yīng)東方蜜蜂微孢子蟲侵染的DEmRNA的篩選及分析

miRNA一般對基因進行負調(diào)控(Pillai, 2005)，故差異變化趨勢相反的miRNA-mRNA關(guān)系有重要意義。因此，本研究對東方蜜蜂微孢子蟲上調(diào)miRNA-意大利蜜蜂工蜂下調(diào)mRNA、東方蜜蜂微孢子蟲下調(diào)miRNA-意大利蜜蜂工蜂上調(diào)mRNA關(guān)系進行重點分析和探討。采用FPKM(fragments per kilobase of transcript per million fragments mapped)法計算和歸一化基因表達量。利用edgeR軟件(Robinsonetal., 2010)篩選AmCK1vsAmT1和AmCK2vsAmT2比較組的顯著性DEmRNA，篩選標準為：|log2fold change|≥1且P≤0.05。

利用OmicShare云平臺(www.omicshare.com)的GO分析軟件和KEGG通路分析軟件對1.4節(jié)的靶mRNA和上述靶DEmRNA進行GO數(shù)據(jù)庫(George, 2008)和KEGG(Kanehisa and Goto, 2006)通路注釋，采用默認參數(shù)。

1.6 東方蜜蜂微孢子蟲顯著性DEmiRNA與意大利蜜蜂工蜂中腸DEmRNA的調(diào)控網(wǎng)絡(luò)構(gòu)建及分析

根據(jù)1.4和1.5節(jié)預(yù)測出的東方蜜蜂微孢子蟲顯著性DEmiRNA與意大利蜜蜂工蜂中腸顯著性DEmRNA的靶向關(guān)系，構(gòu)建NcCKvsNcT1比較組中DEmiRNA與AmCK1vsAmT1比較組中DEmRNA以及NcCKvsNcT2比較組中DEmiRNA與AmCK2vsAmT2比較組中DEmRNA的調(diào)控網(wǎng)絡(luò)，并利用Cytoscape軟件(Smootetal., 2011)進行可視化，采用默認參數(shù)。

2 結(jié)果

2.1 東方蜜蜂微孢子蟲DEmiRNA靶向意大利蜜蜂工蜂中腸mRNA和DEmRNA的預(yù)測及分析

NcCKvsNcT1比較組中東方蜜蜂微孢子蟲的165條顯著性DEmiRNA整體上可靶向意大利蜜蜂工蜂中腸的11 711條mRNA，其中東方蜜蜂微孢子蟲的91條顯著上調(diào)和55條顯著下調(diào)miRNA分別靶向意大利蜜蜂工蜂中腸AmCK1vsAmT1的9 920和8 391條mRNA，其余19條DEmiRNA與意大利蜜蜂工蜂中腸mRNA之間無靶向關(guān)系。進一步分析發(fā)現(xiàn)東方蜜蜂微孢子蟲的77條顯著上調(diào)miRNA和52條顯著下調(diào)miRNA可分別靶向意大利蜜蜂工蜂中腸的118條顯著下調(diào)mRNA和135條顯著上調(diào)mRNA(圖1)。

NcCKvsNcT2比較組中東方蜜蜂微孢子蟲的124條顯著性DEmiRNA整體上可靶向意大利蜜蜂工蜂中腸的11 295條mRNA，其中東方蜜蜂微孢子蟲的56條顯著上調(diào)和49條顯著下調(diào)miRNA可分別靶向AmCK2vsAmT2比較組中意大利蜜蜂工蜂中腸的8 138和9 055條mRNA，其余19條東方蜜蜂微孢子蟲DEmiRNA與意大利蜜蜂工蜂中腸mRNA之間無靶向關(guān)系。進一步分析結(jié)果顯示，東方蜜蜂微孢子蟲的52條顯著上調(diào)miRNA可靶向AmCK2vsAmT2中意大利蜜蜂工蜂中腸的97條顯著下調(diào)mRNA，東方蜜蜂微孢子蟲的49條顯著下調(diào)miRNA可靶向意大利蜜蜂工蜂中腸的210條顯著上調(diào)mRNA(圖2)。

圖1 NcCK vs NcT1比較組中東方蜜蜂微孢子蟲的顯著性DEmiRNA與AmCK1 vs AmT1比較組中意大利蜜蜂工蜂中腸的顯著性DEmRNA之間的調(diào)控網(wǎng)絡(luò)

圖2 NcCK vs NcT2比較組中東方蜜蜂微孢子蟲的顯著性DEmiRNA與AmCK2 vs AmT2比較組中意大利蜜蜂工蜂中腸的顯著性DEmRNA之間的調(diào)控網(wǎng)絡(luò)

進一步分析發(fā)現(xiàn)，NcCKvsNcT1和NcCKvsNcT2中的11條共有的顯著上調(diào)miRNA(miR-3968-y, miR-5119-y和miR-317-y等)可靶向AmCK1vsAmT1和AmCK2vsAmT2中的6條皆顯著下調(diào)的mRNA (GenBank登錄號: XM_006564187.2, XM_006567501.2和XM_006558063.2等)；19條共有的顯著下調(diào)miRNA(novel-m0016-3p, let-7-x和miR-122-x等)靶向AmCK1vsAmT1和AmCK2vsAmT2中14條皆顯著上調(diào)的mRNA (GenBank登錄號: XM_016912123.1, XM_006558673.2和XM_016916036.1等)。

NcCKvsNcT1中共63條特有的DEmiRNA，其中36條顯著上調(diào)miRNA和18條顯著下調(diào)miRNA可分別靶向AmCK1vsAmT1中的97條顯著下調(diào)mRNA和91條顯著上調(diào)mRNA，其余9條東方蜜蜂微孢子蟲DEmiRNA與意大利蜜蜂工蜂中腸mRNA之間不存在靶向關(guān)系。NcCKvsNcT2中共22條特有的DEmiRNA，其中8條顯著上調(diào)miRNA和14條顯著下調(diào)miRNA可分別靶向AmCK2vsAmT2中51條顯著下調(diào)mRNA和139條顯著上調(diào)mRNA。

2.2 東方蜜蜂微孢子蟲顯著性DEmiRNA靶向意大利蜜蜂工蜂中腸顯著性DEmRNA的GO數(shù)據(jù)庫注釋

NcCKvsNcT1中的顯著上調(diào)miRNA靶向AmCK1vsAmT1中意大利蜜蜂工蜂中腸的顯著下調(diào)mRNA涉及細胞(25條下調(diào)mRNA)、細胞組分(25條下調(diào)mRNA)和細胞膜(24條下調(diào)mRNA)等10個細胞組分相關(guān)條目，細胞進程(75條下調(diào)mRNA)、代謝進程(60條下調(diào)mRNA)和單一組織進程(49條下調(diào)mRNA)等15個生物學(xué)進程相關(guān)條目，以及結(jié)合(59條下調(diào)mRNA)、催化活性(55條下調(diào)mRNA)和轉(zhuǎn)運活性(11條下調(diào)mRNA)等6個分子功能相關(guān)條目(圖3: A)。東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向的意大利蜜蜂工蜂中腸顯著上調(diào)mRNA涉及細胞膜(15條上調(diào)mRNA)、細胞膜組分(15條上調(diào)mRNA)和細胞(12條上調(diào)mRNA)等7個細胞組分相關(guān)條目，細胞進程(34條上調(diào)mRNA)、代謝進程(30條上調(diào)mRNA)和單一組織進程(21條上調(diào)mRNA)等10個生物學(xué)進程相關(guān)條目，以及結(jié)合(34條上調(diào)mRNA)、催化活性(24條上調(diào)mRNA)和核酸結(jié)合轉(zhuǎn)錄因子活性(5條上調(diào)mRNA)等8個分子功能相關(guān)條目(圖3: B)。

圖3 NcCK vs NcT1比較組中東方蜜蜂微孢子蟲顯著性DEmiRNA靶向AmCK1 vs AmT1比較組中意大利蜜蜂工蜂中腸顯著性DEmRNA的GO數(shù)據(jù)庫注釋

NcCKvsNcT2中東方蜜蜂微孢子蟲的顯著上調(diào)miRNA靶向AmCK2vsAmT2中意大利蜜蜂工蜂中腸的顯著下調(diào)mRNA涉及細胞膜組分(20條下調(diào)mRNA)、細胞膜(20條下調(diào)mRNA)和細胞(14條下調(diào)mRNA)等11個細胞組分相關(guān)條目，細胞進程(52條下調(diào)mRNA)、代謝進程(42條下調(diào)mRNA)和單一組織進程(35條下調(diào)mRNA)等10個生物學(xué)進程相關(guān)條目，以及結(jié)合(49條下調(diào)mRNA)、催化活性(41條下調(diào)mRNA)和轉(zhuǎn)運活性(9條下調(diào)mRNA)等6個分子功能相關(guān)條目(圖4: A)。東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向的意大利蜜蜂工蜂中腸顯著上調(diào)mRNA涉及細胞(24條上調(diào)mRNA)、細胞組分(24條上調(diào)mRNA)和細胞膜組分(22條上調(diào)mRNA)等9個細胞組分相關(guān)條目，細胞進程(59條上調(diào)mRNA)、代謝進程(47條上調(diào)mRNA)和單一組織進程(42條上調(diào)mRNA)等12個生物學(xué)進程相關(guān)條目，以及結(jié)合(57條上調(diào)mRNA)、催化活性(35條上調(diào)mRNA)和轉(zhuǎn)運活性(14條上調(diào)mRNA)等9個分子功能相關(guān)條目(圖4: B)。

此外，NcCKvsNcT1和NcCKvsNcT2中共同的顯著上調(diào)miRNA靶向AmCK1vsAmT1和AmCK2vsAmT2中共同的顯著下調(diào)mRNA可注釋到單一組織進程(2條下調(diào)mRNA)、細胞進程(2條下調(diào)mRNA)和結(jié)合(2條下調(diào)mRNA)等7個條目；NcCKvsNcT1和NcCKvsNcT2中共同的顯著下調(diào)miRNA靶向AmCK1vsAmT1和AmCK2vsAmT2中共同的顯著上調(diào)mRNA可注釋到結(jié)合(3條上調(diào)mRNA)、代謝進程(2條上調(diào)mRNA)和細胞(2條上調(diào)mRNA)等10個條目。

NcCKvsNcT1中特有的顯著上調(diào)miRNA靶向AmCK1vsAmT1中顯著下調(diào)mRNA涉及26個條目，包括細胞膜組分(9條下調(diào)mRNA)、細胞膜(9條下調(diào)mRNA)和細胞(8條下調(diào)mRNA)等7個細胞組分相關(guān)條目，細胞進程(27條下調(diào)mRNA)、單一組織進程(21條下調(diào)mRNA)和代謝進程(21條下調(diào)mRNA)等13個生物學(xué)進程相關(guān)條目，結(jié)合(24條下調(diào)mRNA)、催化活性(19條下調(diào)mRNA)、分子轉(zhuǎn)導(dǎo)活性(5條下調(diào)mRNA)等6個分子功能相關(guān)的條目。NcCKvsNcT1特有的東方蜜蜂微孢子蟲顯著下調(diào)miRNA靶向AmCK1vsAmT1中顯著上調(diào)mRNA涉及22個功能條目，包括9個生物學(xué)進程相關(guān)條目，7個細胞組分相關(guān)條目，以及6個分子功能相關(guān)條目。NcCKvsNcT2中特有的顯著上調(diào)miRNA靶向AmCK2vsAmT2中顯著下調(diào)mRNA可注釋到22個功能條目，包括細胞(10條下調(diào)mRNA)、細胞組分(10條下調(diào)mRNA)和細胞器(9條下調(diào)mRNA)等7個細胞組分相關(guān)條目，細胞進程(19條下調(diào)mRNA)、代謝進程(16條下調(diào)mRNA)和單一組織進程(10條下調(diào)mRNA)等8個生物學(xué)進程相關(guān)條目，結(jié)合(24條下調(diào)mRNA)、催化活性(13條下調(diào)mRNA)和核酸結(jié)合轉(zhuǎn)錄因子活性(4條下調(diào)mRNA)等7個分子功能相關(guān)條目；NcCKvsNcT2中特有的東方蜜蜂微孢子蟲顯著下調(diào)miRNA靶向AmCK2vsAmT2中意大利蜜蜂工蜂顯著上調(diào)mRNA涉及28個功能條目，包括9個細胞組分相關(guān)條目，11個生物學(xué)進程相關(guān)條目，以及8個分子功能相關(guān)條目。

2.3 東方蜜蜂微孢子蟲顯著性DEmiRNA靶向意大利蜜蜂工蜂中腸顯著性DEmRNA的KEGG通路注釋

NcCKvsNcT1中東方蜜蜂微孢子蟲的顯著上調(diào)miRNA靶向AmCK1vsAmT1中意大利蜜蜂工蜂中腸的顯著下調(diào)mRNA可注釋到剪接體(4條下調(diào)mRNA)、光傳導(dǎo)-果蠅(3條下調(diào)mRNA)和磷酸肌醇代謝(3條下調(diào)mRNA)等113條通路，其中包括胞吞作用(2條下調(diào)mRNA)等細胞免疫通路，以及泛素介導(dǎo)的蛋白水解(2條下調(diào)mRNA)、黑色素生成(3條下調(diào)mRNA)等體液免疫通路(圖5: A)。NcCKvsNcT1中東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向意大利蜜蜂工蜂中腸的顯著上調(diào)mRNA可注釋到Hedgehog信號通路-果蠅(3條上調(diào)mRNA)、粘著斑(3條上調(diào)mRNA)和生理節(jié)律-果蠅(2條上調(diào)mRNA)等107條通路，其中包括胞吞作用(4條上調(diào)mRNA)、泛素介導(dǎo)的蛋白水解(1條上調(diào)mRNA)等細胞免疫通路，黑色素生成(1條上調(diào)mRNA)等體液免疫通路，以及色氨酸代謝(1條上調(diào)mRNA)等3條氨基酸代謝通路，檸檬酸循環(huán)(TCA循環(huán))(2條上調(diào)mRNA)等4條碳水化合物代謝通路和氧化磷酸化(1條上調(diào)mRNA)等2條能量代謝通路(圖5: B)。注釋數(shù)量前20位的通路信息詳見圖5。

圖5 NcCK vs NcT1比較組中東方蜜蜂微孢子蟲顯著性DEmiRNA靶向AmCK1 vs AmT1比較組中意大利蜜蜂工蜂中腸顯著性DEmRNA注釋的前20位KEGG通路

NcCKvsNcT2中東方蜜蜂微孢子蟲的顯著上調(diào)miRNA靶向AmCK2vsAmT2中意大利蜜蜂工蜂中腸的顯著下調(diào)mRNA涉及代謝通路(11條下調(diào)mRNA)、鞘脂信號通路(5條下調(diào)mRNA)及硫代謝(1條下調(diào)mRNA)等97條通路，其中包括溶酶體(1條下調(diào)mRNA)、黑色素生成(1條下調(diào)mRNA)、胞吞作用(1條下調(diào)mRNA)等細胞免疫通路。NcCKvsNcT2中東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向意大利蜜蜂工蜂中腸的顯著上調(diào)mRNA涉及代謝通路(9條上調(diào)mRNA)、Hippo信號通路(5條上調(diào)mRNA)和其他聚糖降解(4條上調(diào)mRNA)等127條通路，其中包括溶酶體(4條上調(diào)mRNA)、胞吞作用(2條上調(diào)mRNA)等細胞免疫通路，以及Toll/Imd信號通路(1條上調(diào)mRNA)、Toll樣受體信號通路(1條上調(diào)mRNA)等體液免疫通路，嘧啶代謝(1條上調(diào)mRNA)等2條核苷酸代謝通路，角質(zhì)、亞黃素和蠟的生物合成(1條上調(diào)mRNA)等3條脂質(zhì)代謝通路，谷胱甘肽代謝(1條上調(diào)mRNA)等2條其他氨基酸的代謝通路，1條萜類和聚酮類化合物的代謝通路。注釋數(shù)量前20位的通路信息詳見圖6。

圖6 NcCK vs NcT2比較組中東方蜜蜂微孢子蟲顯著性DEmiRNA靶向AmCK2 vs AmT2比較組中意大利蜜蜂工蜂中腸顯著性DEmRNA注釋的前20位KEGG通路

此外，NcCKvsNcT1和NcCKvsNcT2中共同的顯著上調(diào)miRNA靶向AmCK1vsAmT1和AmCK2vsAmT2中共同的顯著下調(diào)mRNA未注釋到任何通路；共同的東方蜜蜂微孢子蟲顯著下調(diào)miRNA靶向共同的意大利蜜蜂工蜂中腸顯著上調(diào)mRNA可注釋到9條通路，包括生理節(jié)律(1條上調(diào)mRNA)、生理節(jié)律-果蠅(1條上調(diào)mRNA)、Hippo信號通路-多物種(1條上調(diào)mRNA)、Hedgehog信號通路(1條上調(diào)mRNA)、Hedgehog信號通路-果蠅(1條上調(diào)mRNA)、FoxO信號通路(1條上調(diào)mRNA)、Hippo信號通路(1條上調(diào)mRNA)、Hippo信號通路-果蠅(1條上調(diào)mRNA)和Wnt信號通路(1條上調(diào)mRNA)等。

NcCKvsNcT1中特有的顯著上調(diào)miRNA靶向AmCK1vsAmT1中特有的顯著下調(diào)mRNA涉及代謝途徑(5條下調(diào)mRNA)、剪接體(4條下調(diào)mRNA)和光傳導(dǎo)-果蠅(3條下調(diào)mRNA)等106條通路。NcCKvsNcT1東方蜜蜂微孢子蟲特有的顯著下調(diào)miRNA靶向意大利蜜蜂工蜂中腸特有的顯著上調(diào)mRNA涉及代謝途徑(4條上調(diào)mRNA)、Hippo信號通路(3條上調(diào)mRNA)和生理節(jié)律-果蠅(2條上調(diào)mRNA)等92條通路。

NcCKvsNcT2中特有的顯著上調(diào)miRNA靶向AmCK2vsAmT2中顯著下調(diào)mRNA涉及代謝通路(5條下調(diào)mRNA)、糖胺聚糖降解(1條下調(diào)mRNA)及Notch信號通路(1條下調(diào)mRNA)等36條通路。NcCKvsNcT2中東方蜜蜂微孢子蟲特有的顯著下調(diào)miRNA靶向意大利蜜蜂工蜂中腸特有的顯著上調(diào)mRNA涉及各種類型的N-聚糖的生物合成(4條上調(diào)mRNA)、溶酶體(4條上調(diào)mRNA)及cAMP信號通路(4條上調(diào)mRNA)等112條通路。

2.4 東方蜜蜂微孢子蟲顯著性DEmiRNA及其靶向意大利蜜蜂工蜂中腸的免疫防御相關(guān)DEmRNA的調(diào)控網(wǎng)絡(luò)

根據(jù)2.3節(jié)通路注釋的結(jié)果，東方蜜蜂微孢子蟲顯著性DEmiRNA靶向意大利蜜蜂工蜂中腸DEmRNA涉及的免疫防御相關(guān)通路包括：胞吞作用、黑色素生成、溶酶體、自噬、Toll樣受體信號通路、細胞凋亡、Ras信號通路、泛素介導(dǎo)的蛋白水解和MAPK信號通路。

圖7 東方蜜蜂微孢子蟲的顯著性DEmiRNA靶向意大利蜜蜂工蜂中腸免疫防御和能量代謝相關(guān)通路注釋的顯著性DEmRNA的調(diào)控網(wǎng)絡(luò)

NcCKvsNcT1中顯著上調(diào)的miR-14-y [log2(fold change)=16.41,P=3.76E-05]靶向意大利蜜蜂工蜂中腸的mRNA (GenBank登錄號: XM_395448.6) [log2(fold change)=-10.02,P=0.0076]； miR-283-x [log2(fold change)=5.50,P=0.049]和miR-216-x [log2(fold change)=16.33,P=1.43E-07]共同靶向意大利蜜蜂工蜂中腸的E3泛素蛋白連接酶Nedd-4亞型X8編碼基因的mRNA (GenBank登錄號: XM_016917135.1) [log2(fold change)=-10.57,P=0.0020]。NcCKvsNcT1中顯著下調(diào)的miR-181-x[log2(fold change)=-2.22,P=0.0051]和miR-128-y[log2(fold change)=-1.87,P=5.55E-04]共同靶向的意大利蜜蜂工蜂中腸mRNA (GenBank登錄號: XM_016915020.1) [log2(fold change)=1.04,P=0.019], miR-340-x[log2(fold change)=-1.12,P=0.024]靶向的mRNA(GenBank登錄號: XM_006570195.2)[log2(fold change)=9.52,P=0.043], miR-484-z[log2(fold change)=-1.64,P=0.029]靶向的mRNA (GenBank登錄號: XM_006571076.2)[log2(fold change)=9.40,P=0.033]和miR-320-y[log2(fold change)=-4.73,P=1.18E-04]靶向的mRNA (GenBank登錄號: XM_394836.6) [log2(fold change)=4.01,P=0.0086]均可注釋到胞吞作用(圖7: A)。顯著上調(diào)的miR-3202-x[log2(fold change)=12.47,P=1.75E-07], miR-317-y[log2(fold change)=19.06,P=1.72E-06]和miR-252-x[log2(fold change)=14.99,P=9.90E-05]等6條東方蜜蜂微孢子蟲DEmiRNA共同靶向的意大利蜜蜂工蜂中腸mRNA (GenBank登錄號: XM_016911563.1) [log2(fold change)=-8.28,P=0.0093]， miR-283-x, miR-8-y [log2(fold change)=7.03,P=0.032]和miR-216-x等7條DEmiRNA共同靶向的1-磷脂酰肌醇-4,5-二磷酸磷酸二酯酶I類和II類亞型X1編碼基因mRNA (GenBank登錄號: XM_003250942.3) [log2(fold change)=-3.41,P=0.0084]，miR-339-y [log2(fold change)=12.47,P=1.75E-07], miR-9226-y [log2(fold change)=21.32,P=1.98E-10]和miR-317-y 4條DEmiRNA共同靶向mRNA (GenBank登錄號: XM_006559228.1) [log2(fold change)=-2.45,P=0.045]以及顯著下調(diào)的novel-m0016-3p [log2(fold change)=-7.07,P=1.82E-05]和novel-m0028-5p [log2(fold change)=-14.19,P=3.86E-07]共同靶向mRNA (GenBank登錄號: XM_006567354.2) [log2(fold change)=2.08,P=9.21E-04]均可注釋到黑色素生成(圖7: A)。顯著上調(diào)的miR-14-y靶向的意大利蜜蜂工蜂中腸mRNA (GenBank登錄號: XM_395448.6) [log2(fold change)=-10.02,P=0.0076], miR-283-x和miR-216-x共同靶向的mRNA (GenBank登錄號: XM_016917135.1)以及顯著下調(diào)的miR-106-x [log2(fold change)=-2.30,P=0.026], miR-17-x [log2(fold change)=-1.26,P=0.036], let-7-y [log2(fold change)=-3.76,P=0.0070], miR-20-x [log2(fold change)=-3.70,P=9.29E-04)和miR-146-x [log2(fold change)=-2.86,P=0.037)共同靶向的mRNA (GenBank登錄號: XM_006564220.2) [log2(fold change)=10.41,P=0.010]均可注釋到泛素介導(dǎo)的蛋白水解(圖7: A)。顯著上調(diào)的miR-283-x, miR-216-x和miR-4796-y [log2(fold change)=15.51,P=1.21E-06]共同靶向的mRNA (GenBank登錄號: XM_016914371.1)，顯著下調(diào)的miR-126-y [log2(fold change)=-3.77,P=0.0098]和miR-128-y共同靶向的mRNA (GenBank登錄號: XM_394432.6) [log2(fold change)=7.97,P=0.035]均可注釋到MAPK信號通路(圖7: A)。

在NcCKvsNcT2中，東方蜜蜂微孢子蟲顯著下調(diào)的novel-m0016-3p [log2(fold change)=-4.37,P=3.24E-04]靶向意大利蜜蜂工蜂中腸的mRNA (GenBank登錄號: XM_016917619.1) [log2(fold change)=12.12,P=4.14E-05]和miR-320-y [log2(fold change)=-4.87,P=0.0016]靶向的mRNA (GenBank登錄號: XM_006572368.2) [log2(fold change)=8.59,P=0.0041]以及顯著上調(diào)的miR-5106-y [log2(fold change)=19.69,P=9.27E-09]靶向的mRNA (GenBank登錄號: XM_006561133.2)[log2(fold change)=-1.13,P=0.045]皆可注釋到胞吞作用(圖7: B)。顯著下調(diào)的novel-m0021-5p [log2(fold change)=-1.52,P=0.036]靶向意大利蜜蜂工蜂中腸的4條mRNAs (GenBank登錄號: XM_006567428.2) [log2(fold change)=1.52,P=0.0042]，(GenBank登錄號: XM_006567431.2) [log2(fold change)=1.34,P=0.039]，(GenBank登錄號: XM_016914010.1) [log2(fold change)=1.22,P=1.47E-04]和(GenBank登錄號: XM_016914011.1) [log2(fold change)=3.09,P=0.042]以及顯著上調(diào)的novel-m0022-5p [log2(fold change)=14.13,P=7.97E-06]和miR-216-x [log2(fold change)=16.33,P=1.43E-07]共同靶向的類α-乙酰氨基葡萄糖苷酶編碼基因mRNA (GenBank登錄號: XM_016913157.1) [log2(fold change)=-9.33,P=0.010]皆可注釋到溶酶體(圖7: B)。顯著下調(diào)的novel-m0002-3p [log2(fold change)=-5.85,P=6.42E-05], novel-m0021-5p和miR-374-x[log2(fold change)=-15.48,P=1.23E-07]共同靶向的mRNA (GenBank登錄號: XM_016911751.1) [log2(fold change)=1.65,P=0.0040]以及miR-26-x [log2(fold change)=-1.81,P=0.012]和miR-486-x [log2(fold change)=-8.86,P=0.028]共同靶向的mRNA (GenBank登錄號: XM_006557238.2) [log2(fold change)=3.73,P=0.0047]皆可注釋到自噬-動物(圖7: B)。顯著上調(diào)的bantam-y [log2(fold change)=16.57,P=1.84E-06], miR-8-y [log2(fold change)=5.53,P=0.044], miR-216-x和miR-5119-y [log2(fold change)=9.19,P=0.022]共同靶向的1-磷脂酰肌醇-4,5-二磷酸磷酸二酯酶I類和II類亞型X3編碼基因mRNA (GenBank登錄號: XM_006567225.2) [log2(fold change)=-8.75,P=0.016]可注釋到黑色素生成(圖7: B)。miR-424-x [log2(fold change)=13.86,P=9.14E-06], miR-7465-x [log2(fold change)=15.36,P=4.75E-07]和miR-547-x [log2(fold change)=13.98,P=5.96E-08]共同靶向的mRNA (GenBank登錄號: XM_006564132.2) [log2(fold change)=-8.92,P=0.038]可注釋到MAPK信號通路-果蠅。顯著下調(diào)的miR-26-x和miR-486-x共同靶向的mRNA (GenBank登錄號: XM_006557238.2)可同時注釋到細胞凋亡、MAPK信號通路和Toll樣受體信號通路(圖7: B)。顯著上調(diào)的miR-5106-y靶向的mRNA (GenBank登錄號: XM_006561133.2),顯著下調(diào)的miR-26-x和miR-486-x共同靶向的mRNA (GenBank登錄號: XM_006557238.2), novel-m0021-5p, miR-221-z [log2(fold change)=-2.89,P=0.0039], miR-374-x, miR-423-y [log2(fold change)=-2.04,P=0.012], miR-17-x [log2(fold change)=-1.92,P=0.039]和miR-16-y [log2(fold change)=-16.85,P=5.38E-08]共同靶向的mRNA (GenBank登錄號: XM_006569471.2) [log2(fold change)=5.64,P=0.027]以及顯著上調(diào)的miR-8-y, bantam-y和novel-m0022-5p [log2(fold change)=14.13,P=7.97E-06]等10條miRNA共同靶向的1-磷脂酰肌醇-4,5-二磷酸磷酸二酯酶ε-1編碼基因mRNA (GenBank登錄號: XM_016911813.1) [log2(fold change)=-7.89，P=0.023]皆可注釋到Ras信號通路(圖7: B)。

此外，根據(jù)2.3節(jié)通路注釋的結(jié)果，東方蜜蜂微孢子蟲顯著性DEmiRNA靶向的意大利蜜蜂工蜂中腸DEmRNA可注釋到氧化磷酸化和硫代謝等能量代謝相關(guān)通路。NcCKvsNcT1中顯著下調(diào)的miR-584-x [log2(fold change)=-1.32,P=0.042]和miR-151-y [log2(fold change)=-2.21,P=6.92E-04]等2條miRNA共同靶向的細胞色素c氧化酶亞基4亞型X1，線粒體亞型X1編碼基因mRNA (GenBank登錄號: XM_006572254.2) [log2(fold change)=-2.21,P=6.92E-04],顯著上調(diào)的miR-459-x[log2(fold change)=16.24,P=1.76E-07], miR-315-x[log2(fold change)=12.80,P=4.70E-06]和miR-5119-y分別靶向的mRNAs (GenBank登錄號: XM_001122063.4)[log2(fold change)=-1.07,P=0.022]、(GenBank登錄號: XM_006567193.2)[log2(fold change)=-1.02,P=8.31E-04]和(GenBank登錄號: XM_623229.5) [log2(fold change)=-2.86,P=0.0095]皆可注釋到氧化磷酸化(圖7: C)。NcCKvsNcT2中顯著上調(diào)的miR-301-y [log2(fold change)=15.58,P=6.07E-09]靶向的mRNA (GenBank登錄號: XM_396499.6) [log2(fold change)=-2.09,P=0.0048]可注釋到硫代謝(圖7: D)。

NcCKvsNcT1中特有的miR-186-x [log2(fold change)=1.02,P=0.047767213]靶向意大利蜜蜂工蜂中腸的mRNAs (GenBank登錄號: XM_006559228.1)和(GenBank登錄號: NM_001011611.2) [log2(fold change)=-4.44,P=0.0048]可分別注釋到黑色素生成通路和免疫系統(tǒng)進程條目，NcCKvsNcT2中特有的novel-m0022-5p靶向意大利蜜蜂工蜂中腸的mRNA (GenBank登錄號: XM_006559395.1) [log2(fold change)=-1.23,P=0.0088], mRNA(GenBank登錄號: XM_006560586.2) [log2(fold change)=-1.13,P=1.03E-07]和mRNA(GenBank登錄號: XM_016911813.1)可注釋到代謝途徑。

NcCKvsNcT1中特有的miR-22-x [log2(fold change)=-13.99,P=6.33E-06]靶向意大利蜜蜂工蜂中腸的mRNAs (GenBank登錄號: XM_006563625.2) [log2(fold change)=8.42,P=0.027]和(GenBank登錄號: XM_006559419.2) [log2(fold change)=1.30,P=0.0078]可分別注釋到生物附著和應(yīng)激反應(yīng)條目。NcCKvsNcT2中特有的miR-26-x靶向意大利蜜蜂工蜂中腸的mRNA (GenBank登錄號: XM_006557238.2)可同時注釋到Toll樣受體信號通路、Hippo信號通路、MAPK信號通路及生物進程調(diào)節(jié)功能條目。

3 討論

3.1 東方蜜蜂微孢子蟲DEmiRNA對意大利蜜蜂工蜂中腸的基因表達具有潛在的廣泛影響

本研究中，NcCKvsNcT1中東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向的意大利蜜蜂工蜂中腸顯著上調(diào)mRNA涉及7個細胞組分相關(guān)條目，10個生物學(xué)進程相關(guān)條目，8個分子功能相關(guān)條目(圖3: B)；以及3條氨基酸代謝通路，4條碳水化合物代謝通路，2條能量代謝通路(圖5: B)。NcCKvsNcT2中東方蜜蜂微孢子蟲的顯著下調(diào)miRNA靶向意大利蜜蜂工蜂中腸的顯著上調(diào)mRNA涉及9個細胞組分相關(guān)條目，12個生物學(xué)進程相關(guān)條目，9個分子功能相關(guān)條目(圖4: B)；以及2條核苷酸代謝通路，3條脂質(zhì)代謝通路，2條其他氨基酸的代謝通路，1條萜類和聚酮類化合物的代謝通路(圖6: B)。這些結(jié)果暗示東方蜜蜂微孢子蟲在侵染意大利蜜蜂工蜂中腸的不同時間點通過差異表達部分miRNA對相應(yīng)的意大利蜜蜂工蜂基因表達進行調(diào)控，進而影響上述與意大利蜜蜂工蜂物質(zhì)和能量代謝相關(guān)的多個方面。

此外，NcCKvsNcT1和NcCKvsNcT2比較組共同的顯著下調(diào)miRNA靶向AmCK1vsAmT1和AmCK2vsAmT2比較組共同的顯著上調(diào)mRNA可注釋到10個GO條目，以及9條KEGG通路。上述結(jié)果表明東方蜜蜂微孢子蟲在侵染意大利蜜蜂工蜂中腸的過程中可能通過持續(xù)差異表達一些miRNA對宿主的細胞活動、生化反應(yīng)、生理活動和信號轉(zhuǎn)導(dǎo)產(chǎn)生廣泛的影響。

3.2 東方蜜蜂微孢子蟲可能通過下調(diào)表達部分miRNA對意大利蜜蜂工蜂中腸的能量代謝進行調(diào)控

在長期的協(xié)同進化中，東方蜜蜂微孢子蟲的基因組大大簡化，多數(shù)物質(zhì)和能量代謝的途徑丟失，其增殖高度依賴宿主細胞提供物質(zhì)和能量(Burrietal., 2006)。Chen等(2013)研究發(fā)現(xiàn)東方蜜蜂微孢子蟲侵染可對西方蜜蜂工蜂造成較強的能量脅迫，致使宿主的能量代謝速率加快且糖水的攝入量提高。ATP/ADP轉(zhuǎn)移酶和ABC轉(zhuǎn)運蛋白參與微孢子蟲對宿主物質(zhì)和能量的竊取(Paldietal., 2010)。前期研究中，筆者所在課題組發(fā)現(xiàn)東方蜜蜂微孢子蟲的7個ABC轉(zhuǎn)運蛋白編碼基因在NcCKvsNcT1和NcCKvsNcT2中差異表達(耿四海等, 2020b)。本研究發(fā)現(xiàn)，NcCKvsNcT1中miR-151-y和miR-584-x顯著下調(diào)表達，且能靶向意大利蜜蜂工蜂中腸的mRNA (GenBank登錄號: XM_006572254.2)，該基因可注釋到能量代謝相關(guān)的氧化磷酸化信號通路。上述結(jié)果暗示東方蜜蜂微孢子蟲的DEmiRNA參與意大利蜜蜂工蜂中腸的關(guān)鍵能量代謝過程，控制宿主能量代謝增強以產(chǎn)生更多的ATP，從而掠奪更多的ATP利于自身增殖。

3.3 東方蜜蜂微孢子蟲可能通過上調(diào)表達部分miRNA調(diào)控意大利蜜蜂工蜂中腸的免疫防御

昆蟲的天然免疫系統(tǒng)包括細胞免疫和體液免疫，前者主要是指東方蜜蜂微孢子蟲體被血淋巴細胞等識別、包被和吞噬，而后者主要包括各種識別受體、胞內(nèi)信號因子、抗微生物肽、蛋白酶以及酶級聯(lián)反應(yīng)導(dǎo)致血淋巴凝結(jié)、黑化等(Leulieretal., 2000)。較多的研究證明，東方蜜蜂微孢子蟲在侵染宿主的過程中可通過多種方式跨界調(diào)控宿主的基因表達，并抑制宿主免疫防御(Wangetal., 2017; Shahidetal., 2018; Cuietal., 2019; Huangetal., 2019)。菟絲子將Bc-siR3.2裝載到AGO1蛋白，靶向擬南芥Arabidopsisthaliana的絲裂原活化蛋白激酶2基因(MPK2)和MPK1基因的轉(zhuǎn)錄本，從而抑制宿主的免疫反應(yīng)(Shahidetal., 2018)。球孢白僵菌將bba-milR-1載入囊泡并轉(zhuǎn)運到斯氏按蚊Anophelesstephensi的細胞內(nèi)，進而跨界調(diào)控宿主基因CLIPB9和Spz4的表達，抑制Toll樣受體信號通路和逃避黑化反應(yīng)(Cuietal., 2019)。本研究發(fā)現(xiàn)，NcCKvsNcT1和NcCKvsNcT2比較組中miR-216-x, miR-8-y, miR-4796-y, bantam-y和miR-5119-y等多個共同顯著上調(diào)的miRNA可靶向注釋到黑色素生成、MAPK信號通路、Ras信號通路、泛素介導(dǎo)的蛋白水解和溶酶體等通路的多條顯著下調(diào)mRNA。其中，NcCKvsNcT1中的miR-216-x同時靶向AmCK1vsAmT1中可注釋到宿主免疫相關(guān)的黑色素生成、MAPK信號通路及泛素介導(dǎo)的蛋白水解通路的mRNAs (GenBank登錄號分別為XM_003250942.3, XM_016914371.1和XM_016917135.1)(圖7: A)；NcCKvsNcT2中的miR-216-x還可同時靶向AmCK2vsAmT2中可注釋到意大利蜜蜂工蜂免疫相關(guān)的黑色素生成、溶酶體及Ras信號通路通路的mRNAs (GenBank登錄號分別為XM_006567225.2, XM_016913157.1和XM_016911813.1)(圖7: B)。NcCKvsNcT1中miR-8-y靶向AmCK1vsAmT1中的mRNA可注釋到黑色素生成通路(圖7: A)；NcCKvsNcT2中miR-8-y靶向AmCK2vsAmT2中的2個mRNA可分別注釋到黑色素生成和Ras信號通路(圖7: B)。NcCKvsNcT1中miR-4796-y靶向AmCK1vsAmT1中的mRNA可注釋到MAPK信號通路(圖7: A)；NcCKvsNcT2中miR-4796-y靶向AmCK2vsAmT2中的mRNA可注釋到Ras信號通路(圖7: B)。NcCKvsNcT1中bantam-y靶向AmCK1vsAmT1中的mRNA可注釋到黑色素生成(圖7: A)；NcCKvsNcT2中bantam-y同時靶向AmCK2vsAmT2中2條mRNA可分別注釋到Ras信號通路和黑色素生成(圖7: B)。NcCKvsNcT1中miR-5119-y靶向AmCK1vsAmT1中的mRNA可注釋到黑色素生成(圖7: A)；NcCKvsNcT2中miR-5119-y還同時靶向AmCK2vsAmT2中的2條mRNA可分別注釋到黑色素生成和MAPK信號通路(圖7: B)。

值得注意的是，筆者所在課題組前期研究發(fā)現(xiàn)，NcCKvsNcT2中bantam-y可靶向東方蜜蜂微孢子蟲的己糖激酶編碼基因mRNA (GenBank登錄號: XM_002995838.1)，該基因可注釋到東方蜜蜂微孢子蟲能量代謝相關(guān)的糖酵解/糖異生途徑；并且bantam-y和miR-5119-y可分別靶向東方蜜蜂微孢子蟲自身的ABC轉(zhuǎn)運蛋白編碼基因mRNA (GenBank登錄號分別為XM_002995835.1和XM_002996675.1)(耿四海等, 2020a)。據(jù)此猜測東方蜜蜂微孢子蟲在侵染意大利蜜蜂工蜂中腸的過程中通過不同程度地差異表達miR-216-x和miR-8-y等miRNA對意大利蜜蜂工蜂中腸的黑色素生成、溶酶體和Ras信號通路進行調(diào)控，此外東方蜜蜂微孢子蟲還能通過差異表達bantam-y和miR-5119-y等miRNA一方面調(diào)控意大利蜜蜂工蜂中腸的黑色素生成、Ras信號通路和MAPK信號通路，另一方面調(diào)節(jié)自身的物質(zhì)和能量運輸，從而促進增殖，體現(xiàn)了miRNA作為具有多面性的關(guān)鍵調(diào)控因子在東方蜜蜂微孢子蟲和意大利蜜蜂工蜂中腸的互作中扮演著特殊角色。

3.4 東方蜜蜂微孢子蟲在不同侵染時間點通過特異性差異表達部分miRNA對意大利蜜蜂工蜂中腸基因進行表達調(diào)控

NcCKvsNcT1中特有的miR-186-x靶向意大利蜜蜂工蜂中腸的mRNAs (GenBank登錄號分別為XM_006559228.1和NM_001011611.2)可分別注釋到黑色素生成通路和免疫系統(tǒng)進程條目，NcCKvsNcT2中特有的novel-m0022-5p靶向意大利蜜蜂工蜂可注釋到代謝途徑的mRNAs (GenBank登錄號分別為XM_006559395.1, XM_006560586.2和XM_016911813.1)。我們推測東方蜜蜂微孢子蟲在感染意大利蜜蜂工蜂7 d通過上調(diào)miR-186-x抑制意大利蜜蜂工蜂中腸的免疫防御相關(guān)基因的表達，在感染10 d時通過上調(diào)novel-m0022-5p對意大利蜜蜂工蜂中腸的能量代謝相關(guān)基因進行抑制，從而加強對宿主細胞內(nèi)能量的掠奪，供其增殖所需。此外，NcCKvsNcT1中特有的miR-22-x靶向意大利蜜蜂工蜂可分別注釋到生物附著和應(yīng)激反應(yīng)條目的mRNAs (GenBank登錄號分別為XM_006563625.2和XM_006559419.2)，NcCKvsNcT2中特有的miR-26-x靶向意大利蜜蜂工蜂的mRNA (GenBank登錄號: XM_006557238.2)可同時注釋到Toll樣受體信號通路、MAPK信號通路及生物進程調(diào)節(jié)功能條目。這些結(jié)果表明東方蜜蜂微孢子蟲在感染意大利蜜蜂工蜂7 d和10 d可能通過分別下調(diào)miR-22-x和miR-26-x對上述重要的功能條目和通路相關(guān)基因表達進行抑制，從而影響宿主的正常生命活動。

綜上所述，本研究結(jié)合前期獲得的高質(zhì)量miRNA和mRNA組學(xué)數(shù)據(jù)對東方蜜蜂微孢子蟲DEmiRNA靶向的意大利蜜蜂工蜂中腸mRNA和DEmRNA進行預(yù)測、分析和探討，揭示病原的DEmiRNA對宿主的基因表達調(diào)控具有廣泛性，病原可能通過上調(diào)部分miRNA調(diào)控宿主的免疫防御以促進增殖，通過下調(diào)部分miRNA跨界調(diào)控宿主的能量代謝以加強能量竊取并促進增殖。研究結(jié)果為闡明東方蜜蜂微孢子蟲跨界調(diào)控意大利蜜蜂工蜂中腸的分子機制提供了理論依據(jù)和數(shù)據(jù)基礎(chǔ)，也為深入理解二者的互作機制提供了有價值的參考信息。