中華蜜蜂幼蟲腸道的高表達基因分析

解彥玲，王鴻權，張 璐，侯志賢，刀 晨，江亮亮，熊翠玲，鄭燕珍，徐細建，黃枳腱，郭 睿，陳大福

(福建農林大學 蜂學學院，福建 福州 350002)

中華蜜蜂幼蟲腸道的高表達基因分析

解彥玲，王鴻權，張 璐，侯志賢，刀 晨，江亮亮，熊翠玲，鄭燕珍，徐細建，黃枳腱，郭 睿*，陳大福

(福建農林大學 蜂學學院，福建 福州 350002)

中華蜜蜂(中蜂)是中國特有的本土蜜蜂資源，也是養蜂生產中的常用蜂種之一，具有重要的生態和經濟價值。本研究利用RNA-seq技術對中蜂5和6日齡幼蟲腸道進行深度測序，得到的有效讀段(clean reads)數分別為29290630與26636038，Q20分別為98.45%和98.38%。GO富集分析結果顯示，Ac5中的高表達基因(HEGs)富集在44個GO term，其中基因富集數最多的為代謝進程(753 unigenes)，其次是催化活性(679 unigenes)與結合(657 unigenes)；Ac6中的HEGs富集在44個GO term，其中基因富集數最多的為代謝進程(792 unigenes)，其次是結合(760 unigenes)和催化活性(744 unigenes)。KEGG代謝通路(pathway)富集分析結果顯示，Ac5的HEGs富集在129個pathway，其中基因富集數最多的是核糖體(103 unigenes)、內質網蛋白質加工(68 unigenes)和碳代謝(67 unigenes)；Ac6的HEGs富集在129個代謝通路，其中基因富集數最多的是核糖體(103 unigenes)、氧化磷酸化(88 unigenes)及內質網蛋白質加工(68 unigenes)。進一步對HEGs進行Venn分析，結果顯示，二者共有的HEGs為2 015個，特有的HEGs分別為219和491個。研究結果不僅可提供中蜂幼蟲腸道發育過程中的基因表達譜數據，也為分子水平上深入研究中蜂幼蟲腸道發育提供重要信息。

RNA-seq；中華蜜蜂；幼蟲腸道；高表達基因

東方蜜蜂(Apiscerana)主要集中分布于亞洲的中國、日本、泰國、菲律賓和伊朗等國家，以及俄羅斯遠東地區[1]。中華蜜蜂(Apisceranacerana，簡稱中蜂)特指分布在我國的東方蜜蜂，其距今已有3 000多年的馴養歷史[2]。中蜂是我國重要的本土蜜蜂資源，也是養蜂生產中的常用蜂種之一，在長期的進化過程中已高度適應我國的氣候、蜜源等生態環境，具有諸多優點，例如善于利用零星蜜源、飼料消耗量少、對蜂螨和球囊菌等具有較強的抗性等[3]。

昆蟲腸道具有消化、吸收和免疫防御等重要功能。有關蜜蜂腸道的研究主要集中在腸道菌群，如Gilliam[4]曾通過傳統培養方法對西方蜜蜂(Apismellifera)腸道內的微生物進行分離鑒定，發現其腸道內有復雜而多樣的共生微生物群；Jeyaprakash等[5]利用基于16S rRNA序列分析的分子生物學技術研究了非洲兩種西方蜜蜂的腸道菌，發現其腸道內只有少數幾類特定的細菌；Li等[6]在東方蜜蜂中檢測到5種類型的細菌，其中Pasteurellaceae(Gamma-1)、Lactobacillus (Firm)、Neisseriaceae(Beta)和Bifidobacterium(Bifido)共占東方蜜蜂細菌總量的90%。然而，蜜蜂及其幼蟲腸道發育的相關研究非常滯后，腸道過程中的基因表達規律尚不明確。前期研究中，我們已對中蜂幼蟲腸道參考轉錄組進行了組裝和注釋(未發表數據)，可為其轉錄組學研究提供可靠的參考信息。本研究利用RNA-seq技術對中蜂5日齡和6日齡幼蟲腸道進行深度測序，通過GO(gene ontology)和KEGG(Kyoto encyclopedia of genes and genomes)代謝通路富集分析對幼蟲腸道高表達基因(HEGs)進行轉錄組學研究，研究結果可在轉錄組水平提供中蜂幼蟲腸道發育過程中的基因表達譜數據，也為在分子水平深入研究腸道發育提供重要信息。

1 材料與方法

1.1 生物材料

本研究中使用的中蜂幼蟲取自福建農林大學蜂學學院教學蜂場。

1.2 主要試劑及儀器

DNaseⅠ和Oligotex mRNA Kits Midi試劑盒購自德國Qiagen公司，Dynal M280磁珠購自Invitrogen公司，DNA ligase購自美國Thermo公司，RNA Reagent抽提試劑盒、ExTaqpolymerase及Superscript Ⅱ reverse transcriptase均購自日本TaKaRa公司。純化cDNA的Ampure beads為美國Agencourt產品，cDNA文庫構建試劑盒TruSeqTMDNA Sample Prep Kit-Set A為美國Illumina公司產品。其他試劑為國產分析純。

恒溫恒濕氣候箱購自中國寧波江南儀器廠，高速冷凍離心機購自德國Eppendorf公司，倒置顯微鏡為上海光學儀器五廠產品，超凈工作臺為蘇州安泰空氣技術有限公司產品，PCR儀為美國Bio Rad公司產品，凝膠成像系統為上海培清科技有限公司產品，超低溫冰箱為中科美菱低溫科技股份有限公司產品。

1.3 中蜂幼蟲的人工飼養

參照王倩等[7]的方法配制中蜂幼蟲飼料，其中將D-果糖和D-葡萄糖換為新鮮蜂蜜。預實驗中中蜂幼蟲7日齡成活率達到70%以上。選擇20群群勢較強的蜂群(外觀無白堊病癥狀)，對上述蜂群進行PCR檢測，結果陰性的蜂群作為本研究的實驗蜂群。將2日齡幼蟲移至無菌的24孔細胞培養板(每孔對應1只幼蟲，孔內加有35 ℃預溫的幼蟲飼料)，在35 ℃，相對濕度70%條件下培養至6日齡。每隔24 h更換飼料。

1.4 測序樣品準備

分別剖取中蜂5日齡和6日齡的幼蟲腸道，為盡量減少腸道RNA的降解，將從解剖取樣到樣品放入液氮速凍的時間控制在30 s以內，每剖取一組樣品，液氮速凍并迅速放入-80 ℃超低溫冰箱保存備用。

1.5 cDNA合成與RNA-seq

利用RNAiso Reagent試劑盒抽提上述2個幼蟲腸道的總RNA，然后用RNase-free DNaseⅠ去除基因組DNA殘留。RNA的質量通過瓊脂糖凝膠電泳和NanoDrop ND-1000(NanoDrop， Wilmington， DE， USA)進行檢測。cDNA文庫構建參照張曌楠等[8]的建庫方法。

上述幼蟲腸道樣品委托廣州基迪奧生物科技有限公司進行雙端測序，測序平臺為Illumina Hiseq 2500。測序數據已上傳NCBI SRA數據庫，SRA號為：SRA456721。

1.6 基因表達豐度計算與高表達基因分析

Illumina PE125測序下機的原始讀段(raw reads)，經去除含接頭(adapter)的、含N比例大于10%的和低質量的reads得到有效讀段(clean reads)，進而利用Tophat 軟件將有效讀段(clean reads)映射(mapping)至中蜂幼蟲腸道參考轉錄組(未發表數據)，表達量的計算使用FPKM(fragments per kilobase of transcript per million mapped reads)法，其計算公式為：FPKM=(106C×103)/NL。

選取各腸道樣品FPKM值大于15的高HEGs，利用在線分析工具OmicsShare(http://www.omicshare.com/tools/index.php/Home/Index/index.html)進行GO分類及KEGG pathway富集分析。

2 結果與分析

2.1 RNA-seq測序數據概述

對2個中蜂幼蟲腸道樣品進行深度測序，Ac5和Ac6分別得到的reads數為29893270與27182602，去除低質量reads后的clean reads數分別為29290630與26636038，Q20分別為98.45%和98.38%(表1)，說明本研究中的RNA-seq數據質量良好，可用于進一步分析。Ac5與Ac6的clean reads mapping到中蜂幼蟲腸道參考轉錄組的比例分別為88.98%和89.87%。Ac5和Ac6樣品表達基因數分別為25 987和28 004，分別占參考基因總數(43 556 unigenes)的59.66%和64.29%。

2.2 Ac5與Ac6中高表達基因的GO分類

一般認為FPKM值在0.10～3.75的是低表達基因，在3.75～15.00為中度表達基因，大于15.00的為高度表達基因(HEGs)[9]，根據FPKM值選取所有HEGs進行GO富集分析，結果顯示，Ac5中的HEGs富集在44個GO term，其中基因富集數最多的為代謝進程(753 unigenes)，其次是催化活性(679 unigenes)、結合(657 unigenes)、細胞進程(591 unigenes)和單組織進程(528 unigenes)；Ac6中的HEGs富集在44個GO term，其中基因富集最多的為代謝進程(792 unigenes)，其次是結合(760 unigenes)、催化活性(744 unigenes)、細胞進程(656 unigenes)和單組織進程(582 unigenes)(圖1)。

2.3 Ac5與Ac6中高表達基因的KEGG pathway富集分析

KEGG pathway 富集分析結果顯示，Ac5的HEGs富集在129個代謝通路(pathway)，其中基因富集數最多的是核糖體(103 unigenes)，其次是內質網蛋白質加工(68 unigenes)和碳代謝(67 unigenes)；Ac6的HEGs富集在129個代謝通路，其中基因富集數最多的是核糖體(103 unigenes)，其次是氧化磷酸化(88 unigenes)和內質網蛋白質加工(68 unigenes)(圖2)。上述結果說明，蛋白等基礎物質代謝和能量代謝在中蜂幼蟲腸道發育過程中非?；钴S，腸道發育對物質和能量的需求旺盛。

表1RNA-seq數據信息統計

Table1Overview of RNA-seq data in this study

樣本Sample原始數據Rawreads有效讀段數Cleanreads(%)接頭Adapter(%)低質量讀段Lowquality(%)99%堿基正確率(Q20)(%)Ac52989327029290630(97.98%)271218(0.91%)328722(1.1%)3604436630(98.45%)Ac62718260226636038(97.99%)224524(0.83%)319664(1.18%)3275561178(98.38%)

A， Ac5樣品HEGs的GO分類；B， Ac6樣品HEGs的GO分類；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， 轉運子活性A， GO enrichment analysis of HEGs in Ac5; B， GO enrichment analysis of HEGs in Ac6; 1， Behavior; 2， Biological adhesion; 3， Biological regulation; 4， Cellular component organization or biogenesis; 5， Cellular process; 6， Detoxification; 7， Developmental process; 8， Growth; 9， Immune system process; 10， Localization; 11， Locomotion; 12， Metabolic process; 13， Multi-organism process; 14， Multicellular organismal process; 15， Reproduction; 16， Reproductive process; 17， Response to stimulus; 18， Signaling; 19， Single-organism process; 20， Cell; 21， Cell junction; 22， Cell part; 23， Extracellular matrix; 24， Extracellular region; 25， Extracellular region part; 26， Macromolecular complex; 27， Membrane; 28， Membrane part; 29， Membrane-enclosed lumen; 30， Organelle; 31， Organelle part; 32， Virion; 33， Virion part; 34， Antioxidant activity; 35， Binding; 36， Catalytic activity; 37， Electron carrier activity; 38， Molecular function regulator; 39， Molecular transducer activity; 40， Nucleic acid binding transcription factor activity; 41， Signal transducer activity; 42， Structural molecule activity; 43， Transcription factor activity， protein blinding; 44， Transporter activity圖1 高表達基因的GO富集分析Fig.1 GO enrichment analysis of HEGs

2.4 Ac5與Ac6中高表達基因的Venn分析

對Ac5和Ac6中蜂幼蟲腸道樣品的Venn分析結果顯示(圖3)，二者的共有HEGs為2 015個，特有HEGs分別為219和491個，說明絕大多數HEGs在中蜂5和6日齡幼蟲腸道均有高表達，推測這些共有HEGs對中蜂幼蟲腸道的正常發育十分重要，而特有HEGs在幼蟲腸道的不同發育時期發揮特殊作用。

圖2 高表達基因的KEGG pathway富集分析Fig.2 KEGG pathway enrichment analysis of HEGs

圖3 Ac5和Ac6中高表達基因的Venn分析Fig.3 Venn analysis of HEGs in Ac5 and Ac6

3 討論

中蜂是我國特有的本土蜂種，也是我國生態系統不可或缺的組成部分，具有十分重要的生態價值[10]。此外，中蜂作為養蜂生產中的常用蜂種，具有較高的經濟價值。目前，有關中蜂及其幼蟲腸道發育的分子研究鮮有報道。近年來，高通量測序技術發展迅猛，因其通量大和準確度高而廣泛應用于動物、植物和微生物的研究領域，在蜜蜂的相關研究中也取得一些重要進展[11-13]。本研究利用RNA-seq技術對中蜂5日齡和6日齡腸道進行深度測序，發現二者的表達基因數分別占參考基因總數的59.66%和64.29%，而前期研究顯示4日齡中蜂幼蟲腸道表達基因數為36 602，占84.03%(未發表數據)，說明中蜂幼蟲腸道在不同發育時期所表達的基因數量差別較大，發育早期表達的基因數明顯多于發育后期。本研究還發現，Ac5和Ac6中均有大量HEGs富集在核糖體和內質網蛋白加工，說明中蜂幼蟲腸道的發育過程對蛋白質的需求大增，因而蛋白質代謝旺盛。還發現大量HEGs富集在氧化磷酸化，說明中蜂幼蟲腸道發育過程中伴隨物質代謝的增強，能量代謝也相應提高。

角質層和圍食膜是昆蟲免疫防御的第一道防線[1]，當病原微生物突破這道防線后，將遭遇昆蟲體內細胞免疫與體液免疫的抵抗，包括內吞作用、黑化作用、吞噬作用、蛋白的酶促水解及抗菌肽等[14]。值得注意的是，對于Ac5和Ac6，相當多的HEGs富集在細胞免疫通路，如內吞作用、溶酶體和吞噬體，且隨著發育時間的延長，富集在上述細胞免疫通路的HEGs數量呈增加趨勢，表明在中蜂幼蟲腸道發育后期，細胞免疫通路較為活躍，可能賦予中蜂幼蟲較強的免疫防御能力。

西方蜜蜂的基因組早在2016年就已公布[15]，為其分子研究奠定了重要基礎。前期研究中，我們已對健康及球囊菌脅迫的意大利蜜蜂(Apismelliferaligustica，簡稱意蜂)幼蟲腸道進行了HEGs分析(待發表)。較之意蜂，中蜂的分子研究進展緩慢，其中參考基因組的缺失是最大的制約因素。直到2015年，Park等[16]通過二代測序技術對雄蜂進行測序，組裝得到中蜂參考基因組，然而作者當時并未同時公布基因位置及注釋信息，導致該參考基因組無法被有效利用。腸道對于蜜蜂及其幼蟲而言具有十分重要的作用，它是重要的消化吸收場所和免疫防御器官。我們在前期研究中denovo組裝并注釋了中蜂幼蟲腸道的參考轉錄組，可為深入研究中蜂幼蟲腸道的發育機理、幼蟲響應病原(如球囊菌、東方蜜蜂微孢子蟲和囊狀幼蟲病毒)脅迫的應答機制提供重要的參考信息。

此前，中蜂幼蟲腸道發育的組學研究未見報道，有關中蜂幼蟲腸道發育機理的研究相對滯后。本研究在前期研究基礎上，對中蜂5日齡和6日齡幼蟲腸道進行HEGs分析，研究結果不僅可提供中蜂幼蟲腸道發育過程中的基因表達譜數據，也為在分子水平深入研究中蜂幼蟲腸道發育提供重要信息。本研究是根據基因表達量FPKM值篩選得到HEGs，進而預測其功能，若要在分子水平更加全面、深入地研究中蜂幼蟲腸道的發育，則需要對不同日齡的幼蟲腸道進行差異表達基因分析，進一步建立腸道發育過程中的分子調控網絡、挖掘關鍵調控基因。這將是我們下一步的工作重點。

[1] 龔一飛. 蜜蜂分類與進化[M]. 福建: 福建科學技術出版社， 2000.

[2] 楊冠煌. 中華蜜蜂[M]. 北京: 中國農業科技出版社， 2001.

[3] 張復興. 現代養蜂生產[M]. 北京: 中國農業大學出版社， 1998.

[4] GILLIAM M. Identification and roles of non-pathogenic microflora associated with honey bees[J].FEMSMicrobiologyLetters， 1997， 155(1): 1-10.

[5] JEYAPRAKASH A， HOY M A， ALLSOPP MH. Bacterial diversity in worker adults ofApismelliferacapensis andApismelliferascutellata (Insecta: Hymenoptera) assessed using 16S rRNA sequences[J].JournalofInvertebratePathology， 2003， 84(2): 96-103.

[6] LI JL， QIN HR， WU J， et al.The prevalence of parasites and pathogens in Asian honeybeesApisceranain China[J].PLoSOne， 2012， 7(11): e47955.

[7] 王倩， 孫亮先， 肖培新， 等. 室內人工培育中華蜜蜂幼蟲技術研究[J]. 山東農業科學， 2009(11): 113-116. WANG Q， SUN L X， XIAO P X， et al. Study on technology for indoor artificial feeding ofApisceranaceranalarvae[J].ShandongAgriculturalSciences， 2009(11): 113-116. (in Chinese with English abstract)

[8] 張曌楠， 熊翠玲， 徐細建， 等. 蜜蜂球囊菌的參考轉錄組denovo組裝及SSR分子標記開發[J]. 昆蟲學報， 2017，60(1): 34-44. ZHANG Z N， XIONG C L， XU X J， et al. De novo assembly of a reference transcriptome and development of SSR markers forAscosphaeraapis[J].ActaEntomologicaSinica， 2017， 60(1): 34-44. (in Chinese with English abstract)

[9] MORTAZAVI A， WILLIAMS BA， MCCUE K， et al. Mapping and quantifying mammalian transcriptomes by RNA-Seq[J].NatureMethods， 2008， 5(7): 621-628.

[10] 丁桂玲. 中華蜜蜂(Apiscerana)群體多樣性的研究[D]. 北京: 中國農業科學院， 2006. DING G L. Study on the population diversity ofApiscerana[D]. Beijing: Chinese Academy of Agricultural Sciences， 2006. (in Chinese with English abstract)

[11] CORNMAN R S， LOPEZ D， EVANS J D. Transcriptional response of honey bee larvae infected with the bacterial pathogenPaenibacilluslarvae[J].PLoSOne， 2013， 8(6): e65424.

[12] JULIE A， BARBARA MISME-AUCOUTURIER， BERNARD V， et al. Transcriptome analyses of the honeybee response toNosemaceranaeand insecticides[J].PloSOne， 2014， 9(3): e91686.

[13] CRISTINO A S， BARCHUUK A R， FREITAS F C， et al. Neuroligin-associated microRNA-932 targets actin and regulates memory in the honeybee[J].NatureCommunications， 2014， 5: 5529.

[15] HONEY BEE GENOME SEQUENCING CONSORTIUM. Insights into social insects form the genome of the honeybeeApismellifera[J].Nature， 2006， 443(7114): 931-949.

[16] PARK D， JUNE J W， CHOI B S， et al. Uncovering the novel characteristics of Asian honey bee，Apiscerana， by whole genome sequencing[J].BMCGenomics， 2015， 16: 1.

(責任編輯張 韻)

AnalysisofhighlyexpressedgenesoflarvalgutofApisceranacerana

XIE Yanling， WANG Hongquan， ZHANG Lu， HOU Zhixian， DAO Chen， JIANG Liangliang， XIONG Cuiling， ZHENG Yanzhen， XU Xijian， HUANG Zhijian， GUO Rui*， CHEN Dafu

(CollegeofBeeScience，FujianAgricultureandForestryUniversity，Fuzhou350002，China)

Apisceranaceranais not only a specific honeybee species resource in China， but also a frequently-used honeybee species in apiculture. In the present study， the 5-(Ac5) and 6-day-old (Ac6) larval guts ofA.c.ceranawere sequenced using RNA-seq. After filtration， 29290630 and 26636038 clean reads with a Q20 of 98.45% and 98.38% were obtained in Ac5 and Ac6. GO enrichment analysis results showed that the highly expressed genes (HEGs) in Ac5 were enriched in 44 GO terms， among them metabolic process (753 unigenes)， catalytic activity (679 unigenes) and binding (657 unigenes) were the largest groups; the HEGs in Ac6 were enriched in 44 GO terms， and the mostly enriched one was metabolic process (792 unigenes)， followed by binding (760 unigenes) and catalytic activity (744 unigenes). KEGG enrichment analysis results suggested that the HEGs of Ac5 were enriched in 129 pathways， among them ribosome (103 unigenes) was the largest group， followed by protein processing in endoplasmic reticulum (68 unigenes) and carbon metabolism (67 unigenes); the HEGs of Ac6 were enriched in 129 pathways， the mostly enriched ones were ribosome (103 unigenes)， oxidative phosphorylation (88 unigenes) as well as protein processing in endoplasmic reticulum (68 unigenes). Our findings can not only offer the gene expression profiles during the developmental process of the larval gut ofA.c.ceranaat transcriptome level， but also provide key information for further investigation of the larval gut’s development.

RNA-seq;Apisceranacerana; larval gut; highly expressed genes

S89

：A

：1004-1524(2017)09-1575-06

解彥玲，王鴻權，張璐，等. 中華蜜蜂幼蟲腸道的高表達基因分析[J].浙江農業學報，2017，29(9)： 1575-1580.

10.3969/j.issn.1004-1524.2017.09.22

2017-03-20

現代農業產業技術體系建設專項基金(CARS-45-KXJ7)；福建農林大學科技發展基金(KF2015123)；福建省大學生創新創業項目(201610389053)

解彥玲(1996—)，女，浙江湖州人，主要從事蜂學研究。E-mail: 1103430452@qq.com

*通信作者，郭睿，E-mail: fafu_ruiguo@126.com