馮睿蓉，付中民,2，杜宇，張文德，范小雪，王海朋，萬潔琦，周紫彧，康育欣，陳大福,2，郭睿,2，史培穎

中華蜜蜂幼蟲腸道中微小RNA的鑒定及分析

馮睿蓉1，付中民1,2，杜宇1，張文德1，范小雪1，王海朋1，萬潔琦1，周紫彧1，康育欣1，陳大福1,2，郭睿1,2，史培穎1

1福建農(nóng)林大學動物科學學院（蜂學學院），福州 350002；2福建農(nóng)林大學蜂療研究所，福州 350002

【】利用small RNA-seq（sRNA-seq）技術和生物信息學方法對中華蜜蜂（，簡稱中蜂）幼蟲腸道的微小RNA（microRNA，miRNA）進行全轉錄組鑒定和分析，旨在豐富中蜂的miRNA信息，并為深入研究miRNA調(diào)控中蜂幼蟲腸道發(fā)育的分子機理提供依據(jù)。利用sRNA-seq技術對中蜂4、5和6日齡幼蟲腸道樣品（Ac1、Ac2和Ac3）進行測序，通過數(shù)據(jù)質控獲得有效標簽序列（clean tags）。采用Blast工具將clean tags連續(xù)比對東方蜜蜂（）基因組和miRBase數(shù)據(jù)庫，以鑒定保守miRNA和新miRNA。采用TPM法對miRNA的表達量進行歸一化處理。通過GraphPad Prism 7軟件統(tǒng)計各組腸道樣品中sRNA占比、miRNA長度分布及首位堿基偏向性。利用相關軟件預測上述miRNA靶向的mRNA并進行GO和KEGG數(shù)據(jù)庫注釋。進一步根據(jù)靶向結合關系構建和分析注釋到發(fā)育和免疫相關通路的基因及其靶向miRNA的調(diào)控網(wǎng)絡，并利用Cytoscape軟件進行可視化。利用莖環(huán)反轉錄PCR（Stem-loop RT-PCR）、分子克隆和Sanger測序驗證miRNA的表達和序列的真實性。共鑒定到中蜂的371個保守miRNA和64個新miRNA；這些miRNA的長度介于18—25 nt且首位堿基主要偏向于U；上述miRNA共靶向14 750條mRNA，涉及離子結合、金屬離子結合、細胞膜、細胞膜組件和單一有機體進程等2 270個GO條目，以及內(nèi)吞作用、細胞凋亡、mTOR信號通路、RNA轉運和昆蟲激素的生物合成等332條KEGG通路。進一步分析結果顯示156個miRNA與注釋到Wnt、Hippo、Notch和mTOR等生長發(fā)育相關通路的67個靶基因存在調(diào)控關系，145個miRNA與注釋到Toll、Imd/JNK、Jak-STAT和抗菌效應因子等免疫相關途徑的21個靶基因存在調(diào)控關系。Stem-loop RT-PCR 結果顯示miR-8-y、miR-9-z、miR-14-y、miR-281-y、miR-283-x和miR-306-x均能擴增出預期的特異性片段；Sanger測序結果顯示上述6個miRNA的序列與深度測序結果一致。提供了中蜂miRNA的數(shù)量、結構特征和表達譜；揭示中蜂幼蟲腸道m(xù)iRNA潛在調(diào)控諸多生命進程與細胞活動；中蜂幼蟲腸道的部分miRNA可通過靶向結合相應的mRNA參與調(diào)節(jié)發(fā)育和免疫相關途徑。

東方蜜蜂；中華蜜蜂；幼蟲；腸道；微小RNA；調(diào)控網(wǎng)絡

0 引言

【研究意義】中華蜜蜂（，簡稱中蜂）是高度適應我國自然環(huán)境的東方蜜蜂（）亞種，具有特殊的生態(tài)價值和經(jīng)濟價值。微小RNA（microRNA，miRNA）作為一類重要的基因表達調(diào)控因子，通過在轉錄后水平調(diào)控基因表達而發(fā)揮諸多重要的生物學功能[1]。目前，中蜂的miRNA研究較為滯后。對中蜂幼蟲腸道的miRNA進行全轉錄組鑒定和分析，可豐富中蜂的miRNA信息，也為深入開展中蜂miRNA的功能研究打下基礎?！厩叭搜芯窟M展】研究表明，miRNA參與調(diào)控昆蟲的生長、發(fā)育、生殖、晝夜節(jié)律、學習記憶、免疫應答及抗藥性等重要生物學過程[2-3]。目前，miRBase數(shù)據(jù)庫（http://www.mirbase.org/）已收錄黑腹果蠅（）、家蠶（）、赤擬谷盜（）等26種昆蟲的3 000余個miRNA，包括西方蜜蜂（）的262個miRNA，但尚未收錄東方蜜蜂miRNA。前人對西方蜜蜂的miRNA進行了較多研究[4-5]。例如，CRISTINO等[6]研究發(fā)現(xiàn)，miR-932可靶向調(diào)控西方蜜蜂體內(nèi)肌動蛋白基因的表達，從而影響蜜蜂大腦的記憶產(chǎn)生過程。筆者團隊前期對意大利蜜蜂（，簡稱意蜂）4、5和6日齡幼蟲腸道進行深度測序和分析，共預測到560個miRNA，并發(fā)現(xiàn)ame-bantam和ame-miR-927a等潛在調(diào)控Wnt、Hippo和Notch等信號通路進而參與意蜂幼蟲腸道的發(fā)育[7]；差異表達譜分析揭示ame-miR-6001-3p和miR-342-y等差異表達miRNA（differentially expressed miRNA，DEmiRNA）通過調(diào)控生長發(fā)育及物質能量代謝相關的多條通路潛在影響意蜂工蜂中腸的發(fā)育過程[8]。相比于西方蜜蜂，東方蜜蜂的miRNA研究較為滯后，相關信息有限[9-11]。筆者團隊前期發(fā)現(xiàn)中蜂工蜂中腸內(nèi)miR-8、miR-14和bantam等miRNA維持較高表達水平，而miR-1、miR-92和miR-194等miRNA在宿主響應東方蜜蜂微孢子蟲（）脅迫的過程中不同程度的差異表達，并潛在調(diào)節(jié)應激反應、能量代謝和細胞凋亡等途徑進而參與宿主免疫應答[11]?！颈狙芯壳腥朦c】目前，東方蜜蜂的miRNA信息總體有限，中蜂的miRNA信息更為匱乏?！緮M解決的關鍵問題】利用small RNA-seq（sRNA-seq）技術對中蜂4、5和6日齡幼蟲腸道樣品進行深度測序，通過生物信息學方法對中蜂miRNA進行全轉錄組鑒定，并構建和分析腸道發(fā)育和免疫相關調(diào)控網(wǎng)絡，進而對部分miRNA進行表達和序列驗證，為東方蜜蜂和中蜂的miRNA信息提供有益補充，并為進一步研究miRNA調(diào)控中蜂幼蟲腸道發(fā)育的分子機理打下基礎。

1 材料與方法

試驗于2018年11月至2020年1月在福建農(nóng)林大學動物科學學院（蜂學學院）蜜蜂保護實驗室完成。

1.1 供試蜜蜂幼蟲

中蜂幼蟲取自福建農(nóng)林大學動物科學學院（蜂學學院）教學蜂場。

1.2 中蜂幼蟲的人工飼養(yǎng)與腸道樣品制備

參照筆者所在實驗室前期建立的技術流程[12-13]，從蜂群群勢較強且整體健康的中蜂蜂群中提取已限王產(chǎn)卵的巢脾，將2日齡幼蟲用移蟲針移至已預置50 μL飼料的48孔細胞培養(yǎng)板中，并將培養(yǎng)板放置于溫度為35℃、相對濕度（RH）為90%的培養(yǎng)箱內(nèi)。幼蟲飼料按照蜂王漿63%、無菌水30%、蜂蜜6%和酵母提取物1%的質量比例配制；每24 h更換一次飼料，并及時移除死亡幼蟲。分別剖取4、5和6日齡幼蟲腸道（分別命名為Ac1、Ac2和Ac3），置于RNA-Free的EP管，放入液氮速凍，然后-80℃超低溫冰箱保存?zhèn)溆?。Ac1-1、Ac1-2和Ac1-3為Ac1組的3個生物學重復；Ac2-1、Ac2-2和Ac2-3為Ac2組的3個生物學重復；Ac3-1、Ac3-2和Ac3-3為Ac3組的3個生物學重復。

1.3 sRNA-seq與數(shù)據(jù)質控

按照陳華枝等[14]的方法，利用Trizol法提取上述9個幼蟲腸道樣品的總RNA；通過瓊脂糖凝膠電泳、切膠回收并連接3′接頭序列；15%的變性聚丙烯酰胺凝膠電泳進行分離，切膠回收連接5′接頭序列；反轉錄的產(chǎn)物經(jīng)3.5%的瓊脂糖凝膠電泳分離，進而切膠回收140—160 bp區(qū)域條帶，產(chǎn)物即為終文庫。委托廣州基迪奧生物科技有限公司進行單端測序，測序結果已上傳至NCBI SRA數(shù)據(jù)庫，BioProject號：PRJNA395108。參照杜宇等[10]的方法對下機的原始讀段（raw reads）進行過濾，進一步對得到的有效讀段（clean reads）進行嚴格質控，將得到有效標簽序列（clean tags）用于后續(xù)分析。

1.4 miRNA的生物信息學預測與分析

參照熊翠玲等[15]的方法進行miRNA的預測和分析，先利用Blast工具將各樣品的clean tags序列比對到GenBank（https://www.ncbi.nlm.nih.gov/genbank/）及Rfam（http://rfam.xfam.org/）數(shù)據(jù)庫，去除可注釋上的clean tags，再將未注釋上的clean tags比對到東方蜜蜂的參考基因組（Assembly ACSNU-2.0），以獲得參考基因組上的位置信息。進一步將比對上參考基因組的clean tags與miRBase數(shù)據(jù)庫中已知miRNA的前體序列進行比對，以鑒定保守miRNA（known miRNA），保守miRNA后綴“-x”和“-y”分別代表從miRNA前體的5′臂或3′臂加工而來，“-z”代表來源位置尚不確定；根據(jù)miRNA是否具有經(jīng)典的莖環(huán)結構鑒定新miRNA（novel miRNA）。根據(jù)公式TPM=T×106/N（T表示miRNA的tags，N表示總miRNA的tags）對表達量進行歸一化處理，從而得到各組樣品中全部miRNA的表達譜。通過GraphPad Prism 7軟件統(tǒng)計各組腸道樣品中sRNA的占比、miRNA的長度分布和首位堿基偏向性。利用OmicShare在線工具集合（www.omicshare.com）對各組腸道樣品中的miRNA進行Venn分析。

1.5 miRNA的靶向預測及調(diào)控網(wǎng)絡構建與分析

參照耿四海等[16]和陳華枝等[17]的方法，聯(lián)用TargetFinder/RNAhybrid（v2.1.2）+svm_light（v6.01）、Miranda（v3.3a）、TargetScan（v7.0）軟件預測上述miRNA靶向的mRNA，均采用默認參數(shù)，將3種預測結果的交集作為可靠的靶標集合。利用Blast工具將預測得到的靶mRNA序列與GO和KEGG數(shù)據(jù)庫比對，獲得靶mRNA的功能注釋信息。

鑒于中蜂的全部miRNA與靶向mRNA之間的調(diào)控網(wǎng)絡十分復雜，本研究根據(jù)KEGG數(shù)據(jù)庫及Nr數(shù)據(jù)庫注釋結果篩選出注釋到發(fā)育和免疫相關通路的mRNA，將上述mRNA對應的基因命名為靶基因，根據(jù)miRNA與靶基因之間的結合關系構建相關調(diào)控網(wǎng)絡，進而通過Cytoscape軟件（v3.6.1）進行可視化。

1.6 miRNA的Stem-loop RT-PCR驗證

參照杜宇等[10]的方法，隨機挑選6個miRNA（miR-8-y、miR-9-z、miR-14-y、miR-281-y、miR-283-x和miR-306-x）進行表達驗證。根據(jù)相應的核酸序列，設計特異性Stem-loop引物和上游引物以及通用下游引物，委托上海生工生物工程股份有限公司合成。利用RNA抽提試劑盒（Promega，美國）分別提取中蜂4、5和6日齡幼蟲腸道的總RNA，然后各0.5 μg進行混合作為反轉錄模板，再利用Stem-loop引物進行反轉錄得到相應的cDNA，作為模板進行PCR擴增（翊圣，中國）。PCR產(chǎn)物經(jīng)1.5%瓊脂糖凝膠電泳檢測。膠回收100 bp附近的片段，連接pMD-19T載體（Takara，中國），轉化DH5大腸桿菌（天根，中國）后挑斑搖菌，取菌液PCR呈陽性的菌液送上海生工進行單端測序。

2 結果

2.1 數(shù)據(jù)質控與評估

Ac1、Ac2和Ac3組幼蟲腸道樣品測得的clean reads數(shù)分別平均為11 273 306、11 349 964和11 122 092條，經(jīng)質控得到的clean tags數(shù)分別平均為9 791 926、9 402 531和9 394 648條，占比均≥82.77%（表1）。上述結果說明本研究中sRNA-seq數(shù)據(jù)質量良好，可滿足后續(xù)分析需要。

2.2 中蜂幼蟲腸道m(xù)iRNA結構特征與表達譜

共預測出435個miRNA，包含371個保守miRNA和64個新miRNA，長度分布介于18—25 nt，分布數(shù)量最多的長度為18 nt和22 nt（圖1-A）；此外上述miRNA的首位堿基主要偏向于U，其次偏向于A（圖1-B）。

表1 sRNA-seq數(shù)據(jù)信息概覽

A：miRNA的長度分布Length distribution of miRNAs；B：miRNA的首位堿基偏向性First nucleotide bias of miRNAs

表達譜分析結果顯示，上述435個中蜂miRNA的表達量介于0.08—192 644，表達量最高的前10位miRNA分別為miR-8-y（TPM=192 644）、miR-750-y（TPM=104 019）、miR-184-y（TPM=70 814）、bantam-y（TPM=47 727）、miR-275-y（TPM=44 336）、miR-14-y（TPM=43 636）、miR-9-z（TPM=39 580）、miR-281-y（TPM=36 366）、miR-306-x（TPM=33 676）和miR-283-x（TPM=31 896）。

2.3 中蜂幼蟲腸道m(xù)iRNA的Stem-loop RT-PCR驗證

Stem-loop RT-PCR結果顯示，隨機選擇的6個miRNA（miR-8-y、miR-9-z、miR-14-y、miR-281-y、miR- 283-x和miR-306-x）均真實表達。進一步的分子克隆和Sanger測序結果與其測序結果中的序列完全一致（圖2）。上述結果證明了測序數(shù)據(jù)和預測結果的可靠性。

2.4 中蜂幼蟲腸道m(xù)iRNA的功能注釋

靶向預測結果顯示，中蜂幼蟲腸道的435個miRNA共靶向14 750條mRNA，涉及分子功能、細胞組分和生物學進程相關的2 270個GO條目；其中，有941條靶mRNA可注釋到離子結合（GO:0043167）和565條靶mRNA可注釋到金屬離子結合（GO:0046872）等570個分子功能相關條目，有874條靶mRNA可注釋到細胞膜（GO:0016020）和826條靶mRNA可注釋到細胞膜組件（GO:0044425）等241個細胞組分相關條目，有1 509條靶mRNA可注釋到單一有機體進程（GO: 0044699）和281條靶mRNA可注釋到單一有機體定位（GO:1902578）等1 459個生物學進程相關條目（圖3）。

A：miRNA的Stem-loop RT-PCR Stem-loop RT-PCR of miRNAs；B：miRNA的Stem-loop RT-PCR產(chǎn)物的Sanger測序結果 Sanger sequencing result of Stem-loop RT-PCR product from miRNAs；黃色區(qū)域表示相應miRNA的成熟序列 Yellow regions indicate mature sequences of corresponding miRNAs

圖3 中蜂幼蟲腸道m(xù)iRNA靶向mRNA的GO條目富集圈圖

上述靶mRNA還可注釋到332條通路，包括內(nèi)吞作用（ko04144）和細胞凋亡（ko04210）等細胞進程相關通路，F(xiàn)oxO信號通路（ko04068）和mTOR信號通路（ko04150）等環(huán)境信息加工相關通路，泛素介導的蛋白水解（ko04120）和RNA轉運（ko03013）等遺傳信息進程相關通路，Toll和Imd信號通路（ko04624）等有機系統(tǒng)相關通路，2-氧代羧酸代謝（ko01210）和昆蟲激素的生物合成（ko00981）等與代謝相關通路（圖4）。

2.5 中蜂幼蟲腸道m(xù)iRNA及其靶基因的調(diào)控網(wǎng)絡

共有156個miRNA靶向注釋在Wnt、Hippo、Notch和mTOR等發(fā)育相關通路的67個靶基因；在miRNA-靶基因調(diào)控網(wǎng)絡中，有95個miRNA靶向的33個靶基因涉及Wnt信號通路，76個miRNA靶向的22個靶基因涉及Hippo信號通路，22個miRNA靶向的7個靶基因涉及Notch信號通路，31個miRNA靶向的7個靶基因涉及mTOR信號通路（圖5）。進一步分析發(fā)現(xiàn)，miR-466-x和miR-1277-x等5個miRNA的靶基因同時注釋到Wnt、Hippo、Notch和mTOR信號通路，miR-9394-y和miR-466-y等6個miRNA的靶基因同時注釋到Wnt、Hippo和mTOR信號通路，miR-3759-y、miR-6873-y和miR-6528-x的靶基因同時注釋到Wnt、Hippo和Notch信號通路，miR-6001-y和miR-8109-x等17個miRNA的靶基因同時注釋到Wnt和Hippo信號通路，bantam-y和miR-195-x等9個miRNA的靶基因同時注釋到Wnt和Notch信號通路，miR-4787-x和novel-m0035-3p等5個miRNA的靶基因同時注釋到Wnt和mTOR信號通路，miR-228-x、miR-263-x和miR-183-x的靶基因同時注釋到Hippo和mTOR信號通路，miR-6538-x的靶基因同時注釋到Notch和mTOR信號通路（圖5）。

此外，共有145個miRNA靶向Toll、Imd/JNK、Jak-STAT和抗菌效應因子等免疫通路相關的21個靶基因，88個miRNA靶向Toll信號通路相關的7個基因，60個miRNA靶向Imd/JNK信號通路相關的7個基因，11個miRNA靶向Jak-STAT信號通路相關的2個基因，14個miRNA靶向抗菌效應因子信號通路相關的5個基因（圖6）。進一步分析發(fā)現(xiàn)，miR-466-y和miR-1277-x等4個miRNA的靶基因同時注釋到Toll、Imd/JNK和抗菌效應因子信號通路，miR-235-y和miR-92-y的靶基因同時注釋到Toll、Imd/JNK和Jak-STAT信號通路，miR-283-x和miR-3759-y等16個miRNA的靶基因同時注釋到Toll和Imd/JNK信號通路，miR-4206-y和miR-6516-y等6個miRNA的靶基因同時注釋到Toll和Jak-STAT信號通路，miR-2480-x和novel-m0032-5p等4個miRNA的靶基因同時注釋到Toll和抗菌效應因子信號通路，miR-25-y和miR-653-x的靶基因同時注釋到Imd/JNK和Jak-STAT信號通路，novel-m0033-5p的靶基因同時注釋到Imd/JNK和抗菌效應因子信號通路（圖6）。

圖4 中蜂幼蟲腸道m(xù)iRNA靶向mRNA的KEGG通路富集圈圖

圖5 中蜂幼蟲腸道m(xù)iRNA及其靶向發(fā)育相關基因的調(diào)控網(wǎng)絡

圖6 中蜂幼蟲腸道m(xù)iRNA及其靶向免疫相關通路基因的調(diào)控網(wǎng)絡

3 討論

3.1 鑒定到的中蜂miRNA豐富了東方蜜蜂的miRNA信息

東方蜜蜂與西方蜜蜂同屬于蜜蜂科蜜蜂屬，是親緣關系較近的兩個蜂種。西方蜜蜂的miRNA研究起步較早且研究較多，但有關東方蜜蜂的miRNA研究比較滯后且信息匱乏。筆者團隊前期在意蜂4、5和6日齡幼蟲腸道中共鑒定到560個miRNA，包括515個保守miRNA和45個新miRNA；其中有331個miRNA為3個日齡幼蟲腸道共有且穩(wěn)定表達[7]。本研究利用sRNA-seq技術對中蜂4、5和6日齡幼蟲腸道進行測序，通過生物信息學分析方法鑒定到371個保守miRNA和64個新miRNA，略少于意蜂幼蟲腸道內(nèi)鑒定得到的保守miRNA[7]，其中257個miRNA為中蜂4、5和6日齡幼蟲腸道共有且穩(wěn)定表達。進一步對本研究鑒定到的中蜂幼蟲腸道m(xù)iRNA與前期鑒定到的意蜂幼蟲腸道m(xù)iRNA進行比較分析，發(fā)現(xiàn)有多達174個中蜂幼蟲腸道m(xù)iRNA（如miR-8-y、bantam-y、miR-275-y、miR-9-z和miR-12-z等）的核酸序列與意蜂幼蟲腸道內(nèi)鑒定到的保守miRNA核酸序列[7]完全一致，表明這些miRNA在親緣關系較近的兩個蜂種中具有較強的保守性。這暗示上述保守miRNA在中蜂幼蟲腸道和意蜂幼蟲腸道中可能發(fā)揮相近的調(diào)控功能。此外，筆者發(fā)現(xiàn)中蜂miRNA長度介于18—25 nt，并主要富集在18 nt和22 nt，且首位堿基主要偏向于U，與意蜂及其他昆蟲的miRNA結構特征相似[14,18]。前期研究發(fā)現(xiàn)意蜂幼蟲腸道m(xù)iRNA潛在調(diào)控Wnt、Hippo、FoxO和Hedgehog信號通路相關靶mRNA[7]。本研究中，中蜂miRNA與涉及內(nèi)吞作用、細胞凋亡、泛素介導的蛋白水解、Toll和Imd、FoxO和mTOR信號通路等相關的靶mRNA具有潛在調(diào)控關系，可形成較復雜的調(diào)控網(wǎng)絡。這說明在意蜂幼蟲腸道和中蜂幼蟲腸道中miRNA均發(fā)揮廣泛的調(diào)控作用。本研究鑒定到的miRNA較好地補充了中蜂的miRNA信息，也為東方蜜蜂其他亞種的miRNA研究提供了寶貴的數(shù)據(jù)資源。

3.2 miR-306在中蜂和其他昆蟲中具有很高的保守性

較多的研究表明miRNA在不同物種中具有較高的保守性[19]。通過與miRBase數(shù)據(jù)庫中收錄的miRNA相比較，筆者發(fā)現(xiàn)中蜂miR-306-x的5′端第2—8位種子區(qū)堿基序列（CAGGTAC）與黑腹果蠅dme-miR-306-5p、家蠶bmo-miR-306a-5p、西方蜜蜂ame-miR-306-5p和豌豆長管蚜（）api-miR-306的種子區(qū)堿基序列完全一致，表明中蜂和上述幾種昆蟲的miR-306具有很高的保守性。研究表明miR-306負向調(diào)控，并影響黑腹果蠅的翅脈形成過程[20]；miR-306在意蜂處女蜂王和產(chǎn)卵蜂王卵巢中差異表達，并可能參與卵巢激活過程[21]。miR-306-x在中蜂體內(nèi)的生物學功能尚不清楚，可作為候選靶標用于下一步的功能研究。

3.3 中蜂幼蟲腸道m(xù)iRNA潛在參與調(diào)節(jié)生長發(fā)育相關生命進程

昆蟲等真核生物可通過復雜而精密的信號轉導途徑調(diào)控各類生命活動[22]。Wnt、Hippo、mTOR及Notch等關鍵信號通路上的多個基因已被證實參與昆蟲的生長發(fā)育及器官形成等過程[23-25]。飛蝗（）miR-2/13/71簇通過抑制靶基因的表達，降低卵黃原蛋白（Vg）mRNA的水平，進而影響卵巢的發(fā)育和生殖[24]。植物花粉來源的miR-162a可靶向抑制西方蜜蜂的和黑腹果蠅的，進而延緩蟲體的生長發(fā)育過程[25]。本研究發(fā)現(xiàn)共有156個miRNA靶向Wnt、Hippo、mTOR及Notch信號通路相關的67個基因，此外許多miRNA潛在調(diào)控兩條及以上發(fā)育相關通路。在前期關于意蜂幼蟲腸道和工蜂中腸的miRNA相關研究中，發(fā)現(xiàn)ame-bantam、ame-miR-14和ame-miR-8等大量miRNA參與對意蜂幼蟲腸道內(nèi)Wnt信號通路相關基因的調(diào)控[7]。ame-miR-6001-3p在意蜂工蜂腸道發(fā)育過程中差異表達，并與Hippo信號通路相關基因存在靶向調(diào)控關系[8]。這表明miRNA作為關鍵調(diào)控因子，可能通過調(diào)控靶基因表達參與關鍵信號通路的調(diào)節(jié)，進而影響中蜂幼蟲腸道發(fā)育過程。

3.4 中蜂幼蟲腸道m(xù)iRNA潛在參與調(diào)控免疫途徑應答病原侵染

在與病原的長期協(xié)同進化與相互作用中，蜜蜂進化出多條細胞和體液免疫途徑以應答和抵御病原侵染，其中Toll、Imd、JNK和Jak-STAT等體液免疫途徑可通過激活抗菌肽基因的表達增強defensin、hymenoptaecin、apidaecin、abaecin和apisimin的合成與釋放，并與體內(nèi)多酚氧化酶（PPO）和溶菌酶（Lys）共同參與對多種病原的免疫過程[26-27]。LOUREN?O等[28]研究發(fā)現(xiàn)，ame-miR-2、ame-miR-184和ame-bantam等miRNA在被黏質沙雷氏菌（）和藤黃微球菌（）接種感染的西方蜜蜂工蜂體內(nèi)顯著差異表達，可能通過共同調(diào)控Imd、JNK和Toll等體液免疫途徑，促進抗菌肽的合成以應答細菌的感染。本研究中，中蜂幼蟲腸道的145個miRNA與Toll、Imd/JNK、Jak-STAT和抗菌效應因子等免疫途徑相關的21個基因存在靶向結合關系，例如miR-9-z靶向Imd/JNK信號通路上的（），miR-184-y靶向調(diào)控Toll信號通路上的（），miR-283-x同時靶向通路Toll信號通路上的（）和（）；此外，許多miRNA潛在調(diào)控兩條及以上免疫途徑。上述結果表明相關miRNA與中蜂幼蟲腸道免疫防御的潛在調(diào)控關系。鑒于蜜蜂幼蟲腸道不僅是重要的免疫器官，也是其與蜜蜂球囊菌（）等病原互作的主要部位[29]。上述潛在調(diào)控中蜂幼蟲腸道免疫防御的miRNA值得進一步深入研究。

4 結論

在中蜂幼蟲腸道中鑒定到371個保守miRNA和64個新miRNA，這些miRNA具有與其他昆蟲miRNA類似的結構特征；通過分子生物學手段驗證了6個miRNA的表達和序列真實性；上述miRNA潛在調(diào)控Wnt、Hippo、mTOR和Notch等發(fā)育相關通路，以及Toll、Imd/JNK、Jak-STAT和抗菌效應因子等免疫途徑。

[1] HUSSAIN M, ASGARI S. MicroRNAs as mediators of insect host-pathogen interactions and immunity. Journal of Insect Physiology, 2014, 70: 151-158.

[2] ASGARI S. MicroRNA functions in insects. Insect Biochemistry and Molecular Biology, 2013, 43(4): 388-397.

[3] 楊婕, 謝苗, 徐雪嬌, 白建林, 尤民生. 昆蟲miRNA研究進展. 昆蟲學報, 2021, 64(2): 259-280.

YANG J, XIE M, XU X J, BAI J L, YOU M S. Research progress of insect miRNAs. Acta Entomologica sinica, 2021, 64(2): 259-280. (in Chinese)

[4] 施騰飛, 余林生, 劉芳, 宗超, 汪天澍. 蜜蜂microRNA的研究進展. 昆蟲學報, 2014, 57(5): 601-606.

SHI T F, YU L S, LIU F, ZONG C, WANG T S. Progress in microRNAs in honey bees. Acta Entomologica sinica, 2014, 57(5): 601-606. (in Chinese)

[5] LIU F, SHI T, YIN W, SU X, QI L, HUANG Z Y, ZHANG S, YU L. The microRNA ame-miR-279a regulates sucrose responsiveness of forager honey bees (). Insect Biochemistry and Molecular Biology, 2017, 90: 34-42.

[6] CRISTINO A S, BARCHUK A R, FREITAS F C, NARAYANAN R K, BIERGANS S D, ZHAO Z, SIMOES Z L, REINHARD J, CLAUDIANOS C. Neuroligin-associated microRNA-932 targets actin and regulates memory in the honeybee. Nature communications, 2014, 5: 5529.

[7] 熊翠玲, 杜宇, 陳大福, 鄭燕珍, 付中民, 王海朋, 耿四海, 陳華枝, 周丁丁, 吳素珍, 石彩云, 郭睿. 意大利蜜蜂幼蟲腸道的miRNAs的生物信息學預測及分析. 應用昆蟲學報, 2018, 55(6): 1023-1033.

XIONG C L, DU Y, CHEN D F, ZHENG Y Z, FU Z M, WANG H P, GENG S H, CHEN H Z, ZHOU D D, WU S Z, SHI C Y, GUO R. Bioinformatic prediction and analysis of miRNAs in thelarval gut. Chinese Journal of Applied Entomology, 2018, 55(6): 1023-1033. (in Chinese)

[8] 杜宇, 范小雪, 蔣海賓, 王杰, 范元嬋, 祝智威, 周丁丁, 萬潔琦, 盧家軒, 熊翠玲, 鄭燕珍, 陳大福, 郭睿. 微小RNA及其介導的競爭性內(nèi)源RNA調(diào)控網(wǎng)絡在意大利蜜蜂工蜂中腸發(fā)育過程中的潛在作用. 中國農(nóng)業(yè)科學, 2020, 53(12): 2512-2526.

DU Y, FAN X X, JIANG H B, WANG J, FAN Y C, ZHU Z W, ZHOU D D, WAN J Q, LU J X, XIONG C L, ZHENG Y Z, CHEN D F, GUO R. The potential role of microRNAs and microRNA-mediated competing endogenous networks during the developmental process ofworker’s midgut. Scientia Agricultura Sinica, 2020, 53(12): 2512-2526. (in Chinese)

[9] 吳小波, 王子龍, 石元元, 張飛, 曾志將. 婚飛對中華蜜蜂性成熟處女蜂王sRNAs表達的影響. 中國農(nóng)業(yè)科學, 2013, 46(17): 3721-3728.

WU X B, WANG Z L, SHI Y Y, ZHANG F, ZENG Z J. Effects of mating flight on sRNAs expression in sexual matured virgin queens (). Scientia Agricultura Sinica, 2013, 46(17): 3721-3728. (in Chinese)

[10] 杜宇, 童新宇, 周丁丁, 陳大福, 熊翠玲, 鄭燕珍, 徐國鈞, 王海朋, 陳華枝, 郭意龍, 隆琦, 郭睿. 中華蜜蜂幼蟲腸道響應球囊菌脅迫的microRNA應答分析. 微生物學報, 2019, 59(9): 1747-1764.

DU Y, TONG X Y, ZHOU D D, CHEN D F, XIONG C L, ZHENG Y Z, XU G J, WANG H P, CHEN H Z, GUO Y L, LONG Q, GUO R. MicroRNA responses in the larval gut oftostress. Acta Microbiologica Sinica, 2019, 59(9):1747-1764. (in Chinese)

[11] CHEN D F, DU Y, CHEN H Z, FAN Y C, FAN X X, ZHU Z W, WANG J, XIONG C L, ZHENG Y Z, HOU C S, DIAO Q Y, GUO R. Comparative identification of microRNAs inworkers’ midguts in response toinvasion. Insects, 2019, 10(9): 258.

[12] 郭睿, 杜宇, 熊翠玲, 鄭燕珍, 付中民, 徐國鈞, 王海朋, 陳華枝, 耿四海, 周丁丁, 石彩云, 趙紅霞, 陳大福. 意大利蜜蜂幼蟲腸道發(fā)育過程中的差異表達microRNA及其調(diào)控網(wǎng)絡. 中國農(nóng)業(yè)科學, 2018, 51(21): 4197-4209.

GUO R, DU Y, XIONG C L, ZHENG Y Z, FU Z M, XU G J, WANG H P, CHEN H Z, GENG S H, ZHOU D D, SHI C Y, ZHAO H X, CHEN D F. Differentially expressed microRNA and their regulation networks during the developmental process oflarval gut. Scientia Agricultura Sinica, 2018, 51(21): 4197-4209. (in Chinese)

[13] 郭睿, 杜宇, 童新宇, 熊翠玲, 鄭燕珍, 徐國鈞, 王海朋, 耿四海, 周丁丁, 郭意龍, 吳素珍, 陳大福. 意大利蜜蜂幼蟲腸道在球囊菌侵染前期的差異表達microRNA及其調(diào)控網(wǎng)絡. 中國農(nóng)業(yè)科學, 2019, 52(1): 166-180.

GUO R, DU Y, TONG X Y, XIONG C L, ZHENG Y Z, XU G J, WANG H P, GENG S H, ZHOU D D, GUO Y L, WU S Z, CHEN D F. Differentially expressed microRNAs and their regulation networks inlarval gut during the early stage ofinfection. Scientia Agricultura Sinica, 2019, 52(1): 166-180. (in Chinese)

[14] 陳華枝, 熊翠玲, 祝智威, 王杰, 范小雪, 蔣海賓, 范元嬋, 萬潔琦, 盧家軒, 鄭燕珍, 付中民, 徐國鈞, 陳大福, 郭睿. 基于small RNA組學分析揭示意大利蜜蜂響應東方蜜蜂微孢子蟲脅迫的免疫應答機制. 微生物學報, 2020, 60(7): 1458-1478.

CHEN H Z, XIONG C L, ZHU Z W, WANG J, FAN X X, JIANG H B, FAN Y C, WAN J Q, LU J X, ZHENG Y Z, FU Z M, XU G J, CHEN D F, GUO R. Unraveling the mechanism underlying the immune responses oftostress based on small RNA omics analyses. Acta Microbiologica Sinica, 2020, 60(7): 1458-1478. (in Chinese)

[15] 熊翠玲, 杜宇, 馮睿蓉, 蔣海賓, 史小玉, 王海朋, 范小雪, 王杰, 祝智威, 范元嬋, 陳華枝, 周丁丁, 鄭燕珍, 陳大福, 郭睿. 侵染中華蜜蜂6日齡幼蟲的蜜蜂球囊菌的微小RNA差異表達譜及調(diào)控網(wǎng)絡. 微生物學報, 2020, 60(5): 992-1009.

XIONG C L, DU Y, FENG R R, JIANG H B, SHI X Y, WANG H P, FAN X X, WANG J, ZHU Z W, FAN Y C, CHEN H Z, ZHOU D D, ZHENG Y Z, CHEN D F, GUO R. Differential expression pattern and regulation network of microRNAs ininvading6-day-old larvae. Acta Microbiologica Sinica, 2020, 60(5): 992-1009. (in Chinese)

[16] 耿四海, 石彩云, 范小雪, 王杰, 祝智威, 蔣海賓, 范元嬋, 陳華枝, 杜宇, 王心蕊, 熊翠玲, 鄭燕珍, 付中民, 陳大福, 郭睿. 微小RNA介導東方蜜蜂微孢子蟲侵染意大利蜜蜂工蜂的分子機制. 中國農(nóng)業(yè)科學, 2020, 53(15): 3187-3204.

GENG S H, SHI C Y, FAN X X, WANG J, ZHU Z W, JIANG H B, FAN Y C, CHEN H Z, DU Y, WANG X R, XIONG C L, ZHENG Y Z, FU Z M, CHEN D F, GUO R. The mechanism underlying microRNAs-mediatedinfection toworker. Scientia Agricultura Sinica, 2020, 53(15): 3187-3204. (in Chinese)

[17] 陳華枝, 祝智威, 蔣海賓, 王杰, 范元嬋, 范小雪, 萬潔琦, 盧家軒, 熊翠玲, 鄭燕珍, 付中民, 陳大福, 郭睿. 蜜蜂球囊菌菌絲和孢子中微小RNA及其靶mRNA的比較分析. 中國農(nóng)業(yè)科學, 2020, 53(17): 3606-3619.

CHEN H Z, ZHU Z W, JIANG H B, WANG J, FAN Y C, FAN X X, WAN J Q, LU J X, XIONG C L, ZHENG Y Z, FU Z M, CHEN D F, GUO R. Comparative analysis of microRNAs and corresponding target mRNAs inmycelium and spore. Scientia Agricultura Sinica, 2020, 53(17): 3606-3619. (in Chinese)

[18] XIAO S, WANG B, LI K, XIONG S, YE X, WANG J, ZHANG J, YAN Z, WANG F, SONG Q, STANLEY D W, YE G, FANG Q. Identification and characterization of miRNAs in an endoparasitoid wasp,. Archives of insect biochemistry and physiology, 2020, 103(2): e21633.

[19] ROUSH S, SLACK F J. The let-7 family of microRNAs. Trends in cell biology, 2008, 18(10): 505-516.

[20] SIMOES DA SILVA C J, SOSPEDRA I, APARICIO R, BUSTURIA A. The microRNA-306/abrupt regulatory axis controls wing and haltere growth in. Mechanisms of Development, 2019, 158: 103555.

[21] MACEDO L M, NUNES F M, FREITAS F C, PIRES C V, TANAKA E D, MARTINS J R, PIULACHS M D, CRISTINO A S, PINHEIRO D G, SIM?ES Z L. MicroRNA signatures characterizing caste- independent ovarian activity in queen and worker honeybees (L.). Insect Molecular Biology, 2016, 25(3): 216-226.

[22] KLEINO A, SILVERMAN N. Regulation of theImd pathway by signaling amyloids. Insect biochemistry and molecular biology, 2019, 108: 16-23.

[23] GOU J, LIN L, OTHMER H G. A model for the Hippo pathway in thewing disc. Biophysical journal, 2018, 115(4): 737-747.

[24] SONG J, LI W, ZHAO H, ZHOU S. Clustered miR-2, miR-13a, miR-13b and miR-71 coordinately targetgene to regulate oogenesis of the migratory locust. Insect biochemistry and molecular biology, 2019, 106: 39-46.

[25] ZHU K, LIU M, FU Z, ZHOU Z, KONG Y, LIANG H, LIN Z, LUO J, ZHENG H, WAN P,. Plant microRNAs in larval food regulate honeybee caste development. PLoS genetics, 2017, 13(8): e1006946.

[26] Hu Y T, Wu T C, Yang E C, Wu P C, Lin P T, Wu Y L. Regulation of genes related to immune signaling and detoxification inby an inhibitor of histone deacetylation. Scientific reports, 2017, 7: 41255.

[27] BRUTSCHER L M, DAUGHENBAUGH K F, FLENNIKEN M L. Antiviral defense mechanisms in honey bees. Current opinion in insect science, 2015, 10: 71-82.

[28] LOUREN?O A P, GUIDUGLI-LAZZARINI K R, FREITAS F C, BITONDI M M, SIM?ES Z L. Bacterial infection activates the immune system response and dysregulates microRNA expression in honey bees. Insect Biochemistry and Molecular Biology, 2013, 43(5): 474-482.

[29] 李江紅, 鄭志陽, 陳大福, 梁勤. 影響蜜蜂球囊菌侵染蜜蜂幼蟲的因素及侵染過程觀察. 昆蟲學報, 2012, 55(7): 790-797.

LI J H, ZHENG Z Y, CHEN D F, LIANG Q. Factors influencinginfection on honeybee larvae and observation on the infection process. Acta Entomologica sinica, 2012, 55(7): 790-797. (in Chinese)

Identification and analysis of microRNAs in the larval gut of

FENG RuiRong1, FU ZhongMin1,2, DU Yu1, ZHANG WenDe1, FAN XiaoXue1, WANG HaiPeng1, WAN JieQi1, ZHOU ZiYu1, KANG YuXin1, CHEN DaFu1,2, GUO Rui1,2, SHI PeiYing1

1College of Animal Sciences (College of Bee Science), Fujian Agriculture and Forestry University, Fuzhou 350002;2Apitherapy Research Institute, Fujian Agriculture and Forestry University, Fuzhou 350002

【】In this study, transcriptome-wide identification and analysis of miRNAs in the larval guts ofwas conducted using a combination of small RNA-seq (sRNA-seq) technology and bioinformatic method, aiming to enrich the information ofmiRNAs and offer a basis for further investigation of miRNA-regulated molecular mechanism underlyinglarval gut development.【】Gut samples of4-, 5-, and 6-day-old larvae (Ac1, Ac2, and Ac3 ) were sequenced using sRNA-seq technology, and clean tags were obtained after quality control. By using Blast tool, clean tags were continuously mapped togenome and miRBase database to identify known miRNAs and novel miRNAs. TPM method was used to perform normalization of miRNAs’ expression. The ratio of sRNAs, length distribution of miRNAs and first base bias were calculated with GraphPad Prism 7 software. Using related software, target mRNAs of above-mentioned miRNAs were predicted followed by GO and KEGG database annotation. Further, regulatory networks between genes associated with development and immune-related pathways and corresponding target miRNAs were constructed and analyzed, followed by visualization of regulatory networks with Cytoscape software. The authenticity of miRNA expression and sequence was verified by using Stem-loop RT-PCR, molecular cloning and Sanger sequencing.【】In total, 371 known miRNAs and 64 novel ones ofwere identified; their length was distributed among 18-25 nt, and the first base had an U bias. The aforementioned miRNAs could target 14 750 mRNAs, involving 2 270 GO terms such as ion binding, metal ion binding, membrane, membrane part and single-organism process, as well as 332 KEGG pathways including endocytosis, apoptosis, mTOR signaling pathway, RNA transport and insect hormone biosynthesis. Further investigation suggested that 156 miRNAs could target 67 genes relative to development-associated pathways such as Wnt, Hippo, Notch and mTOR signaling pathways, while 145 miRNAs could target 21 genes relevant to immune-associated pathways such as Toll, Imd/JNK, Jak-STAT signaling pathways and antimicrobial effectors. Stem-loop RT-PCR result indicated that specific fragments with expected sizes could be amplified from miR-8-y, miR-9-z, miR-14-y, miR-281-y, miR-283-x and miR-306-x; Sanger sequencing result demonstrated that sequences of above-mentioned six miRNAs were in accordance with those in deep sequencing result.【】Our findings provide number, structural characteristics and expression profile ofmiRNAs, and unraveled that miRNAs inlarval gut potentially regulate a lot of life processes and cellular activities, part of miRNAs can participate in regulation of development-related and immune-related pathways by targeting corresponding mRNAs.

;; larva; gut;microRNA; regulatory network

10.3864/j.issn.0578-1752.2022.01.017

2021-04-28；

2021-06-07

國家現(xiàn)代農(nóng)業(yè)產(chǎn)業(yè)技術體系建設專項（CARS-44-KXJ7）、福建農(nóng)林大學杰出青年科研人才計劃（xjq201814）、福建農(nóng)林大學碩士生導師團隊項目（郭睿）、福建省大學生創(chuàng)新創(chuàng)業(yè)訓練計劃（202010389012，202010389016）

馮睿蓉，E-mail：frr18705911634@163.com。付中民，E-mail：369699776@qq.com。馮睿蓉和付中民為同等貢獻作者。通信作者史培穎，E-mail：peiyshi@126.com。通信作者郭睿，E-mail：ruiguo@fafu.edu.cn

（責任編輯 岳梅）