意大利蜜蜂工蜂中腸發(fā)育過程中的差異表達(dá)環(huán)狀RNA及其調(diào)控網(wǎng)絡(luò)分析

郭睿，陳華枝，熊翠玲，鄭燕珍，付中民，徐國鈞，杜宇，王海朋，耿四海，周丁丁，劉思亞，陳大福

?

意大利蜜蜂工蜂中腸發(fā)育過程中的差異表達(dá)環(huán)狀RNA及其調(diào)控網(wǎng)絡(luò)分析

郭睿，陳華枝，熊翠玲，鄭燕珍，付中民，徐國鈞，杜宇，王海朋，耿四海，周丁丁，劉思亞，陳大福

（福建農(nóng)林大學(xué)蜂學(xué)學(xué)院，福州 350002）

【目的】環(huán)狀RNA（circular RNA，circRNA）在可變剪接、轉(zhuǎn)錄調(diào)控和來源基因的表達(dá)調(diào)控等方面具有重要功能。本研究旨在探究意大利蜜蜂（，簡(jiǎn)稱意蜂）工蜂中腸發(fā)育過程中circRNA的表達(dá)譜及其發(fā)育過程中的差異表達(dá)circRNA（differentially expressed circRNA，DEcircRNA），進(jìn)而在轉(zhuǎn)錄組水平探究DEcircRNA在中腸發(fā)育中的作用?！痉椒ā炕谇捌讷@得的意蜂7和10日齡工蜂中腸樣品（Am7和Am10）的全轉(zhuǎn)錄組數(shù)據(jù)，利用find_circ軟件從質(zhì)控后的數(shù)據(jù)中預(yù)測(cè)circRNA。采用RPM算法歸一化處理得到circRNA的表達(dá)量。利用DEGseq軟件對(duì)circRNA進(jìn)行差異分析，按照差異倍數(shù)（fold change）≥2、＜0.05及錯(cuò)誤發(fā)現(xiàn)率（false discovery rate，F(xiàn)DR）＜0.05條件篩選DEcircRNA。通過BLAST比對(duì)GO和KEGG數(shù)據(jù)庫，對(duì)DEcircRNA的來源基因進(jìn)行功能和代謝通路注釋。利用TargetFinder軟件預(yù)測(cè)DEcircRNA-miRNA及DEcicRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)，通過Cytoscape v.3.2.1軟件對(duì)調(diào)控網(wǎng)絡(luò)進(jìn)行可視化。通過實(shí)時(shí)熒光定量PCR（RT-qPCR）驗(yàn)證測(cè)序數(shù)據(jù)的可靠性?！窘Y(jié)果】意蜂中腸各樣品比對(duì)上參考基因組的短序列讀段數(shù)平均為19 616 356條。Am7與Am10的組內(nèi)Pearson相關(guān)系數(shù)均≥0.950。共預(yù)測(cè)出256個(gè)DEcircRNA，包括105個(gè)上調(diào)circRNA和151個(gè)下調(diào)circRNA。Novel_circ_009675和novel_circ_013879分別在Am7和Am10中高量表達(dá)。DEcircRNA的來源基因可注釋到包括結(jié)合、單組織進(jìn)程及細(xì)胞進(jìn)程在內(nèi)的32個(gè)GO條目，其中分別有35、35和7個(gè)來源基因注釋到催化活性、代謝進(jìn)程和應(yīng)激反應(yīng)。上述來源基因還可注釋到35條KEGG代謝通路，其中分別有5、5和4個(gè)來源基因注釋到Hippo信號(hào)通路、內(nèi)吞作用和吞噬體；進(jìn)一步分析發(fā)現(xiàn)分別有1、2和2個(gè)來源基因注釋到磷酸肌醇代謝、淀粉和蔗糖代謝和半乳糖代謝等物質(zhì)代謝通路，5、4、3、1和1個(gè)來源基因注釋到內(nèi)吞作用、吞噬體、溶酶體、泛素介導(dǎo)的蛋白水解和MAPK信號(hào)通路等免疫通路。上述結(jié)果表明相應(yīng)的DEcircRNA廣泛參與意蜂工蜂中腸的生長(zhǎng)發(fā)育、新陳代謝和免疫防御。DEcircRNA-miRNA調(diào)控網(wǎng)絡(luò)分析結(jié)果顯示，141個(gè)DEcircRNA可靶向結(jié)合107個(gè)miRNA，其中多數(shù)僅能結(jié)合1—2個(gè)miRNA，但novel_circ_011577和novel_circ_010719結(jié)合的靶miRNA數(shù)可達(dá)32和28個(gè)；此外，mir-136-y、ame-miR-6001-3p及mir-136-y結(jié)合的circRNA數(shù)最多，分別為15、14和14個(gè)，表明相應(yīng)的DEcircRNA可作為競(jìng)爭(zhēng)性內(nèi)源RNA在意蜂中腸發(fā)育過程發(fā)揮作用。進(jìn)一步構(gòu)建DEcircRNAs-ame-miR-6001-3p-mRNA調(diào)控網(wǎng)絡(luò)，分析結(jié)果顯示14個(gè)DEcircRNA可共同靶向結(jié)合ame-miR-6001-3p，表明它們可能通過調(diào)控ame-miR-6001-3p對(duì)意蜂中腸干細(xì)胞的分裂和分化進(jìn)行間接調(diào)控。隨機(jī)選取6個(gè)DEcircRNA進(jìn)行RT-qPCR驗(yàn)證，結(jié)果顯示5個(gè)DEcircRNA的表達(dá)量變化趨勢(shì)與測(cè)序結(jié)果一致，證實(shí)了本研究測(cè)序結(jié)果的可靠性。【結(jié)論】通過對(duì)意蜂工蜂中腸發(fā)育過程中的DEcircRNA深入分析，提供了circRNA在意蜂工蜂中腸發(fā)育過程中的表達(dá)譜和差異表達(dá)信息，揭示了DEcircRNA在中腸發(fā)育過程中的作用，為中腸發(fā)育相關(guān)的關(guān)鍵circRNA的篩選和功能研究打下了基礎(chǔ)。

意大利蜜蜂；中腸；環(huán)狀RNA；調(diào)控網(wǎng)絡(luò)；發(fā)育

0 引言

【研究意義】蜜蜂是自然界最重要的授粉昆蟲，也是社會(huì)行為學(xué)模式昆蟲，具有非常重要的經(jīng)濟(jì)和生態(tài)價(jià)值[1-2]。意大利蜜蜂（，簡(jiǎn)稱意蜂）屬于西方蜜蜂（），具有優(yōu)越的采集能力、造脾能力和分泌蜂王漿能力，在世界各地的養(yǎng)蜂生產(chǎn)中廣泛使用[3]。目前，有關(guān)蜜蜂腸道的發(fā)育機(jī)理及調(diào)控機(jī)制的研究十分滯后。對(duì)于環(huán)狀RNA（circular RNA，circRNA）在蜜蜂腸道發(fā)育過程中的作用，相關(guān)信息仍然缺失。利用circRNA-seq技術(shù)對(duì)意蜂工蜂中腸進(jìn)行測(cè)序，并對(duì)中腸發(fā)育過程中的差異表達(dá)circRNA（DEcircRNA）及其調(diào)控網(wǎng)絡(luò)進(jìn)行深入分析，可揭示DEcircRNA在中腸發(fā)育過程中的作用，為關(guān)鍵circRNA的篩選和功能研究打下基礎(chǔ)。【前人研究進(jìn)展】CircRNA是新近發(fā)現(xiàn)的一類非編碼RNA（non-coding RNA，ncRNA），通過外顯子或/和內(nèi)含子的反向剪切形成環(huán)化RNA分子[4]。CircRNA大量存在于真核細(xì)胞中，在不同物種中具有保守性、豐富性、穩(wěn)定性、組織表達(dá)特異性和時(shí)序表達(dá)特異性等特點(diǎn)[5]。CircRNA已被證明具有多種生物學(xué)調(diào)控功能，包括作為微小RNA（microRNA，miRNA）“海綿”，形成RNA-蛋白質(zhì)復(fù)合物（RNA-binding protein，RBP），以及調(diào)控靶基因的轉(zhuǎn)錄和可變剪接等[6]。最新的研究結(jié)果表明含有核糖體進(jìn)入位點(diǎn)[7]的circRNA能夠翻譯蛋白[8]。與線性RNA相比，circRNA缺少5′帽子和3′尾巴結(jié)構(gòu)，穩(wěn)定性更高且能抵抗核糖核酸外切酶RNase R的消化，因而成為理想的生物標(biāo)志物[9]。利用最新的高通量測(cè)序技術(shù)和生物信息學(xué)分析方法，已在人類[10-12]、動(dòng)物[13]、植物[14]及微生物[8,15-16]中預(yù)測(cè)和鑒定出circRNA。SALZMAN等[12]對(duì)人類15種細(xì)胞類型中的circRNA進(jìn)行預(yù)測(cè)分析，發(fā)現(xiàn)circRNA具有組織特異性、保守性及潛在的調(diào)控功能；SHEN等[13]對(duì)斑馬魚肌肉、卵巢及眼睛等不同組織進(jìn)行高通量測(cè)序，利用3種算法對(duì)circRNA進(jìn)行預(yù)測(cè)，共預(yù)測(cè)到3 868個(gè)circRNA，進(jìn)而對(duì)具有較高可信度的176個(gè)circRNA進(jìn)行實(shí)驗(yàn)驗(yàn)證，發(fā)現(xiàn)其中84%真實(shí)存在；LU等[14]對(duì)水稻進(jìn)行深度測(cè)序及生物信息學(xué)分析，共預(yù)測(cè)出2 354個(gè)circRNA，進(jìn)一步分析發(fā)現(xiàn)circRNA具有相當(dāng)數(shù)量的亞型；GUO等[16]對(duì)蜜蜂球囊菌（）的菌絲及孢子混合樣品進(jìn)行高通量測(cè)序，通過生物信息學(xué)分析發(fā)掘出551個(gè)長(zhǎng)度介于200—600 nt的circRNA，并發(fā)現(xiàn)它們與miRNA存在復(fù)雜的調(diào)控關(guān)系。較之人類和哺乳動(dòng)物，昆蟲的circRNA研究還處于初級(jí)階段，相關(guān)信息極為有限[17-19]。WESTHOLM等[17]在全基因組水平對(duì)果蠅（）的circRNA進(jìn)行分析，發(fā)現(xiàn)其circRNA具有組織表達(dá)特異性，并且可作為果蠅衰老的潛在生物標(biāo)志物；GAN等[18]通過對(duì)家蠶（）中腸和后腸進(jìn)行高通量測(cè)序，共預(yù)測(cè)出來源于1 727個(gè)基因的3 155個(gè)circRNA，同樣具有組織表達(dá)特異性，進(jìn)一步分析發(fā)現(xiàn)中腸和后腸中的來源基因的功能和代謝通路注釋信息類似；CHEN等[19]對(duì)西方蜜蜂蜂王的卵巢組織中的circRNA進(jìn)行預(yù)測(cè)和分析，發(fā)掘出12 211個(gè)circRNA，并通過對(duì)剛產(chǎn)卵蜂王、處女王和限制產(chǎn)卵蜂王進(jìn)行比較分析推測(cè)DEcircRNA可通過競(jìng)爭(zhēng)性結(jié)合miRNA在蜂王卵巢組織活化及產(chǎn)卵過程中發(fā)揮作用。較多的研究表明circRNA可作為競(jìng)爭(zhēng)性內(nèi)源RNA（competing endogenous RNA，ceRNA）與mRNA競(jìng)爭(zhēng)性結(jié)合miRNA，間接調(diào)控基因表達(dá)[20-22]。CHENG等[21]對(duì)椎間盤退行性病變的病人髓核細(xì)胞和組織通過微距陣芯片分析，發(fā)現(xiàn)circVMA21能夠通過吸附miR-200c抑制細(xì)胞凋亡相關(guān)的X蛋白基因，將circVMA21注射到患有椎間盤退行性病變的大鼠模型可使疾病得到緩解；ZHENG等[22]通過對(duì)6個(gè)正常組織和7個(gè)癌變組織進(jìn)行深度測(cè)序，發(fā)現(xiàn)與癌變相關(guān)的CircHIPK3包含結(jié)合9個(gè)miRNA的18個(gè)潛在結(jié)合位點(diǎn)，可靶向結(jié)合miR-124抑制腫瘤形成?！颈狙芯壳腥朦c(diǎn)】蜜蜂中腸既是消化食物、吸收營(yíng)養(yǎng)的重要場(chǎng)所[23]，也是多種病原微生物寄生的主要部位。筆者所在課題組前期已在mRNA組學(xué)水平探究了意蜂4、5和6日齡幼蟲腸道發(fā)育過程中的基因表達(dá)譜和差異表達(dá)規(guī)律[24]，并在長(zhǎng)鏈非編碼RNA（long non-coding RNA，lncRNA）組學(xué)水平探究了意蜂7和10日齡工蜂中腸的差異表達(dá)lncRNA及其調(diào)控網(wǎng)絡(luò)[25]。目前，蜜蜂circRNA的信息極為有限，circRNA在蜜蜂中腸發(fā)育中的作用尚不明確?！緮M解決的關(guān)鍵問題】通過生物信息學(xué)方法對(duì)意蜂7和10日齡工蜂中腸發(fā)育過程中的差異表達(dá)circRNA（DEcircRNA）及其調(diào)控網(wǎng)絡(luò)進(jìn)行深入分析，提供意蜂工蜂中腸的circRNA表達(dá)譜及差異表達(dá)信息，并在組學(xué)水平探究DEcircRNA在中腸發(fā)育過程中的作用。

1 材料與方法

試驗(yàn)于2017年在福建農(nóng)林大學(xué)蜂學(xué)學(xué)院蜜蜂保護(hù)實(shí)驗(yàn)室完成。

1.1 生物材料

意蜂工蜂取自福建農(nóng)林大學(xué)蜂學(xué)學(xué)院教學(xué)蜂場(chǎng)。

1.2 測(cè)序樣品及全轉(zhuǎn)錄組數(shù)據(jù)來源

筆者所在課題組前期已利用二代測(cè)序技術(shù)對(duì)意蜂工蜂中腸進(jìn)行深度測(cè)序，獲得了高質(zhì)量的全轉(zhuǎn)錄組數(shù)據(jù)[25]。其中，測(cè)序樣品的制備過程簡(jiǎn)述如下：從外觀健康且群勢(shì)較強(qiáng)的意蜂蜂群中選擇蜂子狀況良好的封蓋子脾，迅速提至實(shí)驗(yàn)室，放入（34±0.5）℃培養(yǎng)箱，將剛出房的工蜂（記為0 d）放入四周打孔的干凈塑料盒（10只/盒），每個(gè)塑料盒上方插入一支裝有50%（w/v）無菌糖水的飼喂器。（34±0.5）℃飼養(yǎng)工蜂至10 d。每日檢查工蜂存活情況，及時(shí)清理死亡工蜂。進(jìn)行3次生物學(xué)重復(fù)。待意蜂工蜂7日齡和10日齡時(shí)，在超凈臺(tái)用干凈鑷子拉取工蜂中腸，放入RNA-Free的EP管，經(jīng)液氮速凍后迅速轉(zhuǎn)移至-80℃超低溫冰箱保存?zhèn)溆?。意?日齡工蜂中腸樣品的3個(gè)重復(fù)為Am7-1、Am7-2和Am7-3；意蜂10日齡工蜂中腸樣品的3個(gè)重復(fù)為Am10-1、Am10-2和Am10-3。

上述6個(gè)中腸樣品的全轉(zhuǎn)錄組測(cè)序由廣州基迪奧生物科技有限公司完成，測(cè)序平臺(tái)為Illumina HiSeq 4000，測(cè)序策略為雙端（paired-end）測(cè)序，測(cè)序數(shù)據(jù)已上傳NCBI SRA數(shù)據(jù)庫，BioProject號(hào)：PRJNA406998。cDNA建庫過程簡(jiǎn)述如下：利用RNA抽提試劑盒AxyPrepTMMultisource Total RNA Miniprep Kit（TaKaRa公司，中國）抽提蜜蜂中腸樣品的總RNA，為最大限度地保留所有非編碼RNA（ncRNA），去除核糖體RNA后的mRNA和ncRNA用裂解緩沖液隨機(jī)打斷為小片段，作為模板用六堿基隨機(jī)引物、緩沖液、dNTPs、RNase H和DNA polymerase I合成cDNA第2鏈；經(jīng)過QiaQuick PCR試劑盒（Qiagen公司，德國）純化并加EB緩沖液洗脫經(jīng)末端修復(fù)、加堿基A，加測(cè)序接頭，然后通過尿嘧啶-N-糖基化酶（UNG）降解cDNA第2鏈。

1.3 測(cè)序數(shù)據(jù)質(zhì)控及circRNA預(yù)測(cè)

對(duì)于測(cè)序得到的原始讀段（raw reads），去除N的比例大于10%的、質(zhì)量值Q≤10的和堿基數(shù)占整條讀段數(shù)的50%以上的讀段，過濾得到的有效讀段（clean reads）用于后續(xù)的數(shù)據(jù)分析。利用TopHat軟件[26]將各中腸樣品的有效讀段與西方蜜蜂參考基因組（Amel_4.5）[27]進(jìn)行比對(duì)，從比對(duì)結(jié)果中提取未比對(duì)上讀段（unmapped reads），然后截取每一條未比對(duì)上讀段的兩端（默認(rèn)20 nt），得到短序列讀段（anchors reads），繼而再次比對(duì)參考基因組，將得到的比對(duì)結(jié)果提交給find_circ軟件[11]，從而對(duì)circRNA進(jìn)行篩選和鑒定，篩選條件包括：（1）breakpoints=1，即只保留有且只有1個(gè)清晰breakpoint的環(huán)狀RNA；（2）anchor overlap≤2，即每條reads的兩個(gè)anchors reads比對(duì)到基因組上的位置overlap不能超過2 bp；（3）edit≤2，即只允許2 bp錯(cuò)配；（4）n uniq＞2，即uniq reads大于2條；（5）best qual A＞35或best qual B＞35，即每條reads的其中一個(gè)anchors reads比對(duì)到基因組上最好的mapping結(jié)果要比其排第二的結(jié)果的分值高35分以上；（6）n uniq＞int (samples/2)，即支持該circRNA的uniq reads要大于總樣品數(shù)的1/2；（7）circRNA的長(zhǎng)度小于100 kb。

1.4 CircRNA的表達(dá)量、差異表達(dá)及來源基因分析

采用RPM（mapped back-splicing junction reads per million mapped reads）法計(jì)算circRNA的表達(dá)量。先計(jì)算Am7和Am10中所有circRNA的lg (RPM+1)值，再利用基迪奧在線工具平臺(tái)Omicshare（www. omicshare.com/tools）繪制Am7和Am10的circRNA的箱線圖，從而比較二者的circRNA的表達(dá)量，參數(shù)為默認(rèn)參數(shù)。利用DEGseq軟件[28]對(duì)circRNA進(jìn)行差異表達(dá)分析，DEcircRNA的篩選條件包括：（1）差異倍數(shù)（fold change，F(xiàn)C）≥2；（2）＜0.05；（3）錯(cuò)誤發(fā)現(xiàn)率（false discovery rate，F(xiàn)DR）＜0.05。來源于內(nèi)含子或由內(nèi)含子和外顯子環(huán)化的circRNA能夠調(diào)節(jié)親本基因的表達(dá)[22,29]，使用短reads比對(duì)工具Bowtie[30]分別將這些circRNA兩端的anchors reads比對(duì)西方蜜蜂參考基因組[27]，兩端都比對(duì)上同一個(gè)基因即為該circRNA的來源基因，再將來源基因通過BLAST注釋到GO和KEGG數(shù)據(jù)庫，按照值≤0.5篩選顯著富集GO條目或KEGG代謝通路的DEcircRNA的來源基因注釋。

1.5 DEcircRNA的靶miRNA預(yù)測(cè)及調(diào)控網(wǎng)絡(luò)構(gòu)建

利用TargetFinder軟件[31]預(yù)測(cè)DEcircRNA靶向結(jié)合的miRNA，得到DEcircRNA-miRNA的靶向調(diào)控關(guān)系。按照≤0.05，自由能≤35的標(biāo)準(zhǔn)，從結(jié)果中提取DEcircRNA潛在的靶向結(jié)合miRNA的位點(diǎn)信息。進(jìn)一步預(yù)測(cè)DEcircRNA的靶miRNA靶向結(jié)合的mRNA，根據(jù)DEcircRNA與靶miRNA、miRNA與靶mRNA的結(jié)合關(guān)系，得到DEcircRNA-miRNA-mRNA的靶向調(diào)控關(guān)系，最后通過Cytoscape v.3.2.1軟件[32]對(duì)各調(diào)控網(wǎng)絡(luò)進(jìn)行可視化，參數(shù)采用默認(rèn)參數(shù)。

1.6 DEcircRNA的實(shí)時(shí)熒光定量PCR（RT-qPCR）驗(yàn)證

為了驗(yàn)證circRNA-seq數(shù)據(jù)的可靠性，隨機(jī)選取6個(gè)DEcircRNA（novel_circ_000663、novel_circ_005547、novel_circ_014049、novel_circ_002507、novel_circ_ 012440及novel_circ_001915）進(jìn)行RT-qPCR驗(yàn)證。參照GUO等[16]的方法，利用DNAMAN軟件根據(jù)所選circRNA的堿基序列設(shè)計(jì)相應(yīng)的特異性反向引物，委托上海生工生物工程有限公司進(jìn)行合成，相關(guān)引物信息詳見表1。利用RNA抽提試劑盒（TaKaRa公司，中國）分別提取Am7和Am10樣品的總RNA，分為兩份，其中一份總RNA用3 U/mg RNase R（吉賽生物公司，中國）進(jìn)行消化以去除線性RNA，37℃處理15 min后以隨機(jī)引物進(jìn)行反轉(zhuǎn)錄，得到circRNA的cDNA；另一份總RNA直接作為模板，以O(shè)ligo (dT)23作為引物進(jìn)行反轉(zhuǎn)錄，得到所有含polyA的RNA對(duì)應(yīng)的cDNA，作為內(nèi)參基因的qPCR模板[33]。反應(yīng)體系（20 μL）包含正、反向引物（10.0 μmol·L-1）各1 μL，cDNA模板DNA 1 μL，SYBR Green Dye 10 μL，DEPC水7 μL。反應(yīng)在QuantStudio熒光定量PCR儀（ThermoFisher公司，美國）上進(jìn)行，按照SYBR Green Dye試劑盒（Vazyme公司，中國）操作說明書上的方法，每個(gè)反應(yīng)進(jìn)行3次重復(fù)。反應(yīng)條件：95℃預(yù)變性1 min，95℃變性15 s，60℃延伸30 s，共45個(gè)循環(huán)，最后72℃延伸45 s。利用2?ΔΔCt法對(duì)上述DEcircRNAs的相對(duì)表達(dá)量進(jìn)行計(jì)算。通過Graph Prism軟件進(jìn)行相關(guān)數(shù)據(jù)分析及繪圖。

表1 RT-qPCR驗(yàn)證的引物信息

2 結(jié)果

2.1 高通量測(cè)序數(shù)據(jù)的質(zhì)控及評(píng)估

前期研究中的意蜂工蜂中腸樣品Am7及Am10測(cè)序分別得到134 802 058條raw reads，經(jīng)過濾得到clean reads數(shù)平均為147 051 470條，各樣品的平均Q20和Q30分別為97.34%和93.86%[25]。此外，Am7與Am10的組內(nèi)Pearson相關(guān)系數(shù)都≥0.950（圖1）。上述結(jié)果表明本研究中的測(cè)序數(shù)據(jù)質(zhì)量良好，可用于下一步分析。本研究中，各樣品的anchors reads平均為181 583 825條，比對(duì)上參考基因組的anchors reads數(shù)平均為19 616 356條（表2）。

表2 未比對(duì)上核糖體數(shù)據(jù)庫的有效讀段比對(duì)參考基因組的信息統(tǒng)計(jì)

2.2 意蜂工蜂中腸circRNA的表達(dá)譜分析

對(duì)Am7和Am10中所有circRNA進(jìn)行表達(dá)量分析，結(jié)果顯示前者中的circRNA的整體表達(dá)水平較高（圖2-A）。Am7中表達(dá)量最高的circRNA為novel_ circ_003183、novel_circ_010717、novel_circ_012530、novel_circ_000476及novel_circ_011750（表3）；Am10中表達(dá)量最高的circRNA為novel_circ_003183、novel_ circ_010717、novel_circ_012530、novel_circ_000896、novel_circ_000476（表4）。進(jìn)一步分析發(fā)現(xiàn)，novel_ circ_003183、novel_circ_010717、novel_circ_012530及novel_circ_000476在Am7和Am10中均高量表達(dá)，推測(cè)它們?cè)谝夥涔し渲心c的發(fā)育過程中發(fā)揮重要作用。Novel_circ_009675和novel_circ_013879分別在Am7和Am10中高量表達(dá)，推測(cè)二者分別在意蜂工蜂中腸發(fā)育的不同時(shí)期發(fā)揮特殊功能。差異分析結(jié)果顯示，Am7 vs Am10中共有256個(gè)DEcircRNA，包含105個(gè)上調(diào)和151個(gè)下調(diào)circRNA（圖2-B）。

2.3 意蜂工蜂中腸DEcircRNA來源基因的GO數(shù)據(jù)庫注釋

GO分類結(jié)果顯示，Am7 vs Am10的DEcircRNA來源基因可注釋到32個(gè)GO條目（term），涉及細(xì)胞組分、生物學(xué)進(jìn)程和分子功能3大功能分類（圖3）。細(xì)胞組分中，注釋基因數(shù)最多的前5位分別是細(xì)胞（21個(gè)）、細(xì)胞組件（21個(gè)）、細(xì)胞膜組件（15個(gè)）、細(xì)胞膜（15個(gè)）及細(xì)胞器（15個(gè)）；生物學(xué)進(jìn)程中，注釋基因數(shù)最多的前5位分別是單組織進(jìn)程（38個(gè)）、細(xì)胞進(jìn)程（38個(gè)）、代謝進(jìn)程（35個(gè)）、定位（13個(gè)）及生物調(diào)控（10個(gè)），此外應(yīng)激反應(yīng)和發(fā)育進(jìn)程注釋基因數(shù)分別為7和2個(gè)；分子功能中，注釋基因數(shù)最多的前5位分別是結(jié)合（50個(gè)）、催化活性（35個(gè)）、轉(zhuǎn)運(yùn)器活性（8個(gè)）、核苷酸結(jié)合轉(zhuǎn)錄因子活性（4個(gè)）及分子傳感器活性（3個(gè)）。上述結(jié)果說明DEcircRNA在意蜂工蜂中腸的生長(zhǎng)、發(fā)育、代謝、免疫和細(xì)胞生命活動(dòng)中發(fā)揮功能。

圖1 各意蜂工蜂中腸樣品的不同生物學(xué)重復(fù)間的Pearson相關(guān)性

A：Am7和Am10中circRNA表達(dá)量的箱線圖Box plot of the circRNAs’ expression quantity in Am7 and Am10；B：Am7 vs Am10中的DEcircRNA DEcircRNAs in Am7 vs Am10

表3 Am7樣品中表達(dá)量最高的前10位circRNA

圖3 DEcircRNA來源基因的GO數(shù)據(jù)庫注釋

2.4 意蜂工蜂中腸DEcircRNA來源基因的KEGG數(shù)據(jù)庫注釋

KEGG數(shù)據(jù)庫注釋結(jié)果顯示（圖4），意蜂工蜂中腸DEcircRNA的來源基因可注釋到35條代謝通路（pathway），其中注釋基因數(shù)最多的前10位分別是Hippo信號(hào)通路（5個(gè)，圖5-A）、內(nèi)吞作用（5個(gè)，圖5-B）、吞噬體（4個(gè)）、FoxO信號(hào)通路（3個(gè)）、mRNA監(jiān)測(cè)通路（3個(gè)）、溶酶體（3個(gè)）、RNA轉(zhuǎn)導(dǎo)（3個(gè)）、背腹軸形成（2個(gè)）、淀粉和蔗糖代謝（2個(gè)）及半乳糖代謝（2個(gè)），說明DEcircRNA通過參與多條信號(hào)通路、免疫通路對(duì)意蜂工蜂中腸的生長(zhǎng)發(fā)育、免疫防御進(jìn)行調(diào)節(jié)。內(nèi)吞作用、Hippo信號(hào)通路的概貌詳見圖5-A和圖5-B。

2.5 意蜂工蜂中腸的DEcircRNA-miRNA、DEcircRNA- ame-miR-6001-3p-mRNA調(diào)控網(wǎng)絡(luò)構(gòu)建及分析

通過TargetFinder和Cytoscape軟件構(gòu)建及可視化意蜂工蜂中腸的DEcircRNA-miRNA調(diào)控網(wǎng)絡(luò)，88個(gè)上調(diào)circRNA靶向結(jié)合71個(gè)miRNA，其中novel_ circ_011577、novel_circ_010719及novel_circ_ 009951結(jié)合的miRNA數(shù)最多，分別為32、28和18個(gè)；53個(gè)下調(diào)circRNA可靶向結(jié)合36個(gè)miRNA，其中novel_ circ_013731、novel_circ_002319及novel_circ_006352結(jié)合的miRNA數(shù)最多，分別為23、6和5個(gè)。多數(shù)DEcircRNA（70.21%）僅能結(jié)合1—2個(gè)miRNA。此外，mir-136-y、ame-miR-6001-3p及mir-136-y結(jié)合的circRNA數(shù)最多，分別為15、14和14個(gè)circRNA（圖6）。影響miRNA與靶基因的結(jié)合，從而對(duì)下游的基因表達(dá)進(jìn)行間接調(diào)控。進(jìn)一步構(gòu)建DEcircRNA-ame- miR-6001-3p-mRNA的調(diào)控網(wǎng)絡(luò)（圖7），分析結(jié)果顯示共有14個(gè)DEcircRNA（novel_circ_011733、novel_circ_001055、novel_circ_000661、novel_circ_ 002196、novel_circ_014839、novel_circ_000663、novel_ circ_005547、novel_circ_010719、novel_circ_006612、novel_circ_004676、novel_circ_006035、novel_circ_ 002507、novel_circ_014049及novel_circ_010227）靶向結(jié)合同一個(gè)miRNA（ame-miR-6001-3p），說明這些DEcircRNA可通過競(jìng)爭(zhēng)性結(jié)合ame-miR-6001-3p而抑制其與靶基因的結(jié)合，從而影響靶基因的表達(dá)水平。上述結(jié)果說明意蜂工蜂中腸發(fā)育過程中DEcircRNA通過作為ceRNA充當(dāng)miRNA的海綿，間接調(diào)控基因表達(dá)。

表4 Am10中表達(dá)量最高的前10位circRNA

圓圈大小代表富集在某一通路的基因數(shù)多少，越大代表基因數(shù)越多；圓圈顏色代表某一通路的富集基因的顯著性高低，越紅代表顯著性越高

紅色方框代表注釋在該通路的DEcircRNA的來源基因Red boxes indicate DEcircRNAs’ source genes annotated in certain pathway

圖6 意蜂工蜂中腸的DEcircRNA-miRNA調(diào)控網(wǎng)絡(luò)

圖7 意蜂工蜂中腸的DEcircRNA-ame-miR-6001-3p-mRNA調(diào)控網(wǎng)絡(luò)

2.6 DEcircRNA的RT-qPCR驗(yàn)證

隨機(jī)挑取6個(gè)DEcircRNA進(jìn)行RT-qPCR驗(yàn)證，結(jié)果顯示其中5個(gè)DEcircRNA的表達(dá)量變化趨勢(shì)和轉(zhuǎn)錄組數(shù)據(jù)的結(jié)果一致（圖8），說明本研究中的測(cè)序數(shù)據(jù)真實(shí)可靠。

3 討論

CircRNA在可變剪接、轉(zhuǎn)錄調(diào)控和來源基因的表達(dá)調(diào)控等方面具有重要功能[6]。隨著高通量測(cè)序技術(shù)和生物信息學(xué)的發(fā)展，人們已在人類[10-12]、牛[34]、豬[35]、雞[36]、老鼠[37]、斑馬魚[13]、草魚[38]、水稻[14]、大麥[39]、擬南芥[40]、線蟲[41]、蜜蜂球囊菌[16]、古細(xì)菌[15]和丁型肝炎病毒[8]等動(dòng)植物及微生物中預(yù)測(cè)和鑒定出circRNA。2018年，CHEN等[19]通過高通量測(cè)序技術(shù)和生物信息學(xué)方法從西方蜜蜂蜂王的卵巢組織中預(yù)測(cè)出12 211個(gè)circRNA，這是迄今唯一有關(guān)蜜蜂circRNA的研究報(bào)道。前期研究中，筆者通過對(duì)意蜂7和10日齡工蜂中腸進(jìn)行全轉(zhuǎn)錄組測(cè)序和全面分析，提供了中腸發(fā)育過程中l(wèi)ncRNA的表達(dá)譜及差異表達(dá)信息，構(gòu)建并分析了差異表達(dá)lncRNA（DElncRNA）的調(diào)控網(wǎng)絡(luò)，在組學(xué)水平揭示了DElncRNA在意蜂工蜂中腸發(fā)育過程中的作用[25]。本研究在此基礎(chǔ)上，進(jìn)一步通過生物信息學(xué)方法對(duì)意蜂工蜂中腸的circRNA進(jìn)行挖掘，從意蜂7和10日齡工蜂中腸中分別預(yù)測(cè)出7 341和9 092個(gè)circRNA，為蜜蜂的ncRNA信息提供了重要補(bǔ)充。

本研究發(fā)現(xiàn)，novel_circ_003183、novel_circ_010717、novel_circ_012530、novel_circ_000896、novel_circ_ 000476、novel_circ_011750、novel_circ_011749、novel_ circ_006484、novel_circ_012115及novel_circ_009675在Am7和Am10中皆高量表達(dá)，說明它們?cè)谝夥渲心c發(fā)育過程中具有重要功能，但它們是否在其他日齡的意蜂工蜂中腸中高量表達(dá)，有待于進(jìn)一步研究。Novel_circ_009675和novel_circ_013879分別在Am7和Am10中高量表達(dá)，二者在意蜂中腸發(fā)育的不同時(shí)期發(fā)揮特殊功能。Am7 vs Am10中包含105個(gè)上調(diào)和151個(gè)下調(diào)circRNA，推測(cè)這些DEcircRNA與意蜂中腸發(fā)育關(guān)系密切。

RT-qPCR組中，*表示p＜0.05，**表示p＜0.01 In RT-qPCR group, * indicates p＜0.05, ** indicates p＜0.01

來源于外顯子和內(nèi)含子組成的circRNA（exon- intron circRNA，EIciRNA）可通過與RNA聚合酶II、U1小核核糖核蛋白及基因啟動(dòng)子相互作用對(duì)其來源基因的轉(zhuǎn)錄進(jìn)行調(diào)控[29]。蜜蜂的天然食物主要包括含葡萄糖、果糖和淀粉等成分的蜂蜜及含蛋白質(zhì)的花粉[42]。蜜蜂中腸作為主要的消化和吸收器官，含有大量的消化酶如淀粉酶、蔗糖酶和蛋白酶等[43]。本研究中，DEcircRNA的來源基因中有35個(gè)可注釋到催化活性；分別有1、2和2個(gè)來源基因可注釋到磷酸肌醇代謝、淀粉和蔗糖代謝和半乳糖代謝等消化吸收相關(guān)通路。上述結(jié)果表明相應(yīng)的DEcircRNA參與意蜂工蜂中腸的食物消化、營(yíng)養(yǎng)吸收等過程。此外，分別有2和35個(gè)來源基因注釋到發(fā)育進(jìn)程和代謝進(jìn)程，暗示DEcircRNA在此過程扮演重要角色。Hippo信號(hào)通路可通過抑制細(xì)胞增殖和促進(jìn)細(xì)胞凋亡對(duì)器官大小進(jìn)行調(diào)節(jié)[44]，也能和其他信號(hào)通路相互作用共同調(diào)節(jié)腸道組織的穩(wěn)態(tài)[45]，還在腸道結(jié)構(gòu)的維持特別是上皮細(xì)胞的分化過程起重要作用[46]。Wnt信號(hào)通路與哺乳動(dòng)物胚胎形成、卵巢發(fā)育、平面細(xì)胞極化等生理學(xué)進(jìn)程密切相關(guān)[47]，還能夠與Hippo、Notch和TGF-beta信號(hào)通路共同影響昆蟲的體節(jié)形成、色素沉淀、附肢發(fā)育以及翅等器官的發(fā)育[48]。本研究發(fā)現(xiàn)，分別有5和2個(gè)來源基因注釋到Hippo和Wnt信號(hào)通路，表明相應(yīng)的DEcircRNA可能通過調(diào)節(jié)Hippo和Wnt信號(hào)通路對(duì)意蜂工蜂中腸的細(xì)胞生長(zhǎng)、分化及凋亡以及結(jié)構(gòu)和穩(wěn)態(tài)維持進(jìn)行調(diào)控。還發(fā)現(xiàn)有7個(gè)來源基因注釋到應(yīng)激反應(yīng)，表明相應(yīng)的DEcircRNA參與調(diào)控意蜂工蜂中腸對(duì)外界環(huán)境的適應(yīng)過程。

蜜蜂的免疫系統(tǒng)包括群體免疫與個(gè)體免疫，后者又分為細(xì)胞免疫和體液免疫。其中，蜜蜂的細(xì)胞免疫包括內(nèi)吞作用、吞噬作用、黑化作用等，體液免疫主要為抗菌肽的合成與釋放[49]。MAPK信號(hào)通路能夠通過產(chǎn)生免疫效應(yīng)因子刺激淋巴細(xì)胞活性，對(duì)外界入侵的病原微生物進(jìn)行免疫防御[50]。本研究中，對(duì)于意蜂工蜂中腸發(fā)育過程中DEcircRNA的來源基因，分別有5、4、3和1個(gè)注釋到內(nèi)吞作用、吞噬體、溶酶體和泛素介導(dǎo)的蛋白水解等細(xì)胞免疫通路，僅有1個(gè)注釋在MAPK信號(hào)通路。上述結(jié)果表明相應(yīng)的DEcircRNA參與意蜂工蜂中腸的細(xì)胞免疫、體液免疫的調(diào)控過程，并且可能在細(xì)胞免疫調(diào)控方面發(fā)揮更重要的作用。

2011年，Salmena等[51]提出了“競(jìng)爭(zhēng)性內(nèi)源RNA”假說，即任何包含miRNA反應(yīng)元件（miRNA response element）的RNA，例如mRNA、假基因、lncRNA和circRNA，可以競(jìng)爭(zhēng)性地結(jié)合miRNA。此后，越來越多的研究結(jié)果都證實(shí)了ceRNA假 說[20-22,52-54]。Deng等[52]通過對(duì)患急性心肌梗塞病人和正常人的DEcircRNA進(jìn)行分析，發(fā)現(xiàn)circRNA_081881能夠靶向結(jié)合miR-548，使急性心肌梗塞相關(guān)的表達(dá)量上調(diào)；Xu等[53]分別比較分析了患惡性腫瘤的組織和正常組織的DEcircRNA，發(fā)現(xiàn)hsa_circ_000984通過靶向結(jié)合miR-106b，有效上調(diào)細(xì)胞生長(zhǎng)調(diào)節(jié)相關(guān)基因的表達(dá)水平。本研究中，調(diào)控網(wǎng)絡(luò)分析結(jié)果顯示141個(gè)DEcircRNA可靶向結(jié)合107個(gè)miRNA，其中大部分DEcircRNA（70.21%）僅能結(jié)合1—2個(gè)miRNA，少數(shù)DEcircRNA可結(jié)合多個(gè)miRNA，例如novel_circ_011577和novel_circ_010719結(jié)合的miRNA數(shù)可達(dá)32和28個(gè)，表明二者作為ceRNA的重要性。還發(fā)現(xiàn)多個(gè)DEcircRNA可共同靶向結(jié)合同一個(gè)miRNA，例如分別有15和14個(gè)DEcircRNA共同結(jié)合mir-136-y和ame-miR-6001-3p。上述結(jié)果表明DEcircRNA可作為ceRNA吸附結(jié)合miRNA，減少其對(duì)mRNA的抑制和降解，從而影響意蜂工蜂腸道的生長(zhǎng)發(fā)育。目前，circRNA的功能研究尚處起步階段，絕大多數(shù)circRNA的功能還不明確。Liu等[54]通過對(duì)骨關(guān)節(jié)炎相關(guān)circRNA的深入研究，發(fā)現(xiàn)CircRNA-CER可通過吸附mir-136對(duì)的表達(dá)進(jìn)行調(diào)控，從而參與軟骨細(xì)胞外基質(zhì)（EMC）的降解。本研究中，novel_circ_000976、novel_circ_005547和novel_circ_006344等23個(gè)DEcircRNA靶向結(jié)合的mir-136-y和mir-136-x與mir-136具有高度同源性，推測(cè)相應(yīng)的DEcircRNA可能通過對(duì)mir-136-y、mir-136-x表達(dá)水平的調(diào)節(jié)，參與意蜂工蜂中腸細(xì)胞外基質(zhì)的降解。蛻皮激素主要調(diào)控昆蟲幼蟲老舊組織的程序化細(xì)胞死亡與成蟲干細(xì)胞的分裂分化[55]，Mello等[56]通過敲除變態(tài)發(fā)育過程中的蜜蜂個(gè)體蛻皮激素受體編碼基因，發(fā)現(xiàn)ame-miR-6001-3p在敲除個(gè)體中上調(diào)表達(dá)，表明ame-miR-6001-3p能夠抑制蛻皮激素的分泌，從而影響中腸干細(xì)胞的分裂和分化。本研究中，DEcircRNA-ame-miR-6001-3p-mRNA調(diào)控網(wǎng)絡(luò)分析結(jié)果顯示，共有novel_circ_000661、novel_circ_000663和novel_circ_001055等14個(gè)DEcircRNA能夠共同靶向結(jié)合ame-miR-6001-3p，推測(cè)上述DEcircRNA可通過競(jìng)爭(zhēng)性結(jié)合ame-miR-6001-3p調(diào)控蛻皮激素受體編碼基因的表達(dá)，從而影響中腸干細(xì)胞的分裂及分化。進(jìn)一步通過設(shè)計(jì)跨反向剪切位點(diǎn)的引物對(duì)隨機(jī)挑選DEcircRNA進(jìn)行RT-qPCR驗(yàn)證，結(jié)果顯示83.33%（5/6）的DEcircRNA相對(duì)表達(dá)量與測(cè)序結(jié)果一致，證實(shí)了本研究中circRNA差異表達(dá)信息的可靠性。

本研究?jī)H對(duì)意蜂7和10日齡工蜂中腸的circRNA進(jìn)行相關(guān)分析，預(yù)測(cè)出的circRNA在其他日齡工蜂中腸、同一日齡的不同組織是否表達(dá)以及表達(dá)水平的高低仍需要進(jìn)一步研究。此外，若要全面解析意蜂工蜂中腸的發(fā)育機(jī)理及調(diào)控機(jī)制，則需對(duì)更多日齡的工蜂中腸進(jìn)行測(cè)序，進(jìn)而在全局水平進(jìn)行更加深入的分析，此為下一步的工作重點(diǎn)。

4 結(jié)論

通過對(duì)意蜂7和10日齡工蜂中腸的深度測(cè)序和生物信息學(xué)分析，提供了中腸發(fā)育過程的circRNA表達(dá)譜及差異表達(dá)信息。DEcircRNA可能通過調(diào)控來源基因的表達(dá)水平和作為ceRNA在意蜂工蜂中腸發(fā)育過程中發(fā)揮重要功能。靶向結(jié)合ame-miR-6001-3p的14個(gè)DEcircRNA可能參與意蜂中腸干細(xì)胞的分裂及分化。

[1] PARK D, JUNG J W, CHOI B S, JAYAKODI M, LEE J, LIM J, Yu Y, Choi Y S, Lee M L, Park Y, Choi I Y, Yang T J, Edwards O R, Nah G. Uncovering the novel characteristics of Asian honey bee,, by whole genome sequencing., 2015, 16: 1.

[2] National Research Council.. the National Academies Press, 2006.

[3] 周冰峰. 蜜蜂飼養(yǎng)管理學(xué). 廈門: 廈門大學(xué)出版社, 2002.

ZHOU B F.. Xiamen: Xiamen University Publishing Company, 2002. (in Chinese)

[4] MENG S J, ZHOU H C, FENG Z Y, XU Z H, TANG Y, LI P Y, WU M H. CircRNA: functions and properties of a novel potential biomarker for cancer., 2017, 16: 94.

[5] HE J, XIE Q C, XU H L, LI J T, LI Y S. Circular RNAs and cancer., 2017, 396: 138-144.

[6] QU S B, YANG X S, LI X L, WANG J L, GAO Y, SHANG R Z, SUN W, DOU K F, LI H M. Circular RNA: a new star of noncoding RNAs., 2015, 365(2): 141-148.

[7] ASHWAL-FLUSS R, MEYER M, PAMUDURTI N R, IVANOV A, BARTOK O, HANAN M, EVANTAL N, MEMCZAK S, RAJEWSKY N, KADENER S. CircRNA biogenesis competes with pre-mRNA splicing., 2014, 56(1): 55-66.

[8] KOS A, DIJKEMA R, ARNBERG A C, VAN DER MEIDE P H, SCHELLEKENS H. The hepatitis delta (delta) virus possesses a circular RNA., 1986,323(6088): 558-560.

[9] PERKEL J M. Assume nothing: the tale of circular RNA., 2013, 55(2): 55-57.

[10] JECK W R, SORRENTINO J A, WANG K, SLEVIN M K, BURD C E, LIU J, MARZLUFF W F, SHARPLESS N E. Circular RNAs are abundant, conserved, and associated with ALU repeats., 2013, 19(2): 141-157.

[11] MEMCZAK S, JENS M, ELEFSINIOTI A, Torti F, Krueger J, Rybak A, Maier L, Mackowiak S D, Gregersen L H, Munschauer M, Loewer A, Ziebold U, Landthaler M, Kocks C, le Noble F, Rajewsky N. Circular RNAs are a large class of animal RNAs with regulatory potency., 2014, 495(7441): 333-338.

[12] SALZMAN J, CHEN R E, OLSEN M N, WANG P L, BROWN P O. Cell-type specific features of circular RNA expression., 2013, 9(9): e1003777.

[13] SHEN Y D, GUO X W, WANG W M. Identification and characterization of circular RNAs in zebrafish., 2017, 591(1): 213-220.

[14] LU T T, CUI L L, ZHOU Y, ZHU C R, FAN D L, GONG H, ZHAO Q, ZHOU C C, ZHAO Y, LU D F, LUO J H, WANG Y C, TIAN Q L, FENG Q, HUANG T, HAN B. Transcriptome-wide investigation of circular RNAs in rice., 2015, 21(12): 2076-2087.

[15] DANAN M, SCHWARTZ S, EDELHEIT S, SOREK R. Transcriptome- wide discovery of circular RNAs in archaea., 2012, 40(7): 3131-3142.

[16] GUO R, CHEN D F, CHEN H Z, FU Z M, XIONG C L, HOU C S, ZHENG Y Z, GUO Y L, WANG H P, DU Y, DIAO Q Y. Systematic investigation of circular RNAs in, a fungal pathogen of honeybee larvae., 2018, 678: DOI: 10.1016/j.gene. 2018.07.076.

[17] WESTHOLM J O, MIURA P, OLSON S, SHENKER S, JOSEPH B, SANFILIPPO P, CELNIKER S E, GRAVELEY B R, LAI E C. Genome-wide analysis ofcircular RNAs reveals their structural and sequence properties and age-dependent neural accumulation., 2014, 9(5): 1966-1980.

[18] GAN H, FENG T, WU Y,Liu C, Xia Q, Cheng T. Identification of circular RNA in thesilk gland., 2017, 89: 97-106.

[19] CHEN X, SHI W, CHEN C. Differential circular RNAs expression in ovary during oviposition in honeybees., 2018: doi. org/10.1016/j.ygeno.2018.03.015.

[20] DU W W, YANG W, LIU E, YANG Z, DHALIWAL P, YANG B B. Foxo3 circular RNA retards cell cycle progression via forming ternary complexes with p21 and CDK2., 2016, 44(6): 2846-2858.

[21] CHENG X, ZHANG L, ZHANG K, ZHANG G, HU Y, SUN X, ZHAO C, LI H, LI Y M, ZHAO J. Circular RNA VMA21 protects against intervertebral disc degeneration through targeting miR-200c and X linked inhibitor-of-apoptosis protein., 2018, 77(5): 770-779.

[22] ZHENG Q P, BAO C Y, GUO W J, LI S Y, CHEN J, CHEN B, LUO Y T, LYU D B, LI Y, SHI G H, LIANG L H, GU J R, HE X H, HUANG S L. Circular RNA profiling reveals an abundant circHIPK3 that regulates cell growth by sponging multiple miRNAs., 2016, 7: 11215.

[23] 李兆英. 意大利蜜蜂胚后發(fā)育過程中中腸上皮組織細(xì)胞的更替. 昆蟲學(xué)報(bào), 2011, 54(10): 1127-1132.

LI Z Y. Replacement of midgut epithelium in(Hymenoptera: Apidae) during postembryonic development., 2011, 54(10): 1127-1132. (in Chinese)

[24] 郭睿, 解彥玲, 熊翠玲, 尹偉軒, 鄭燕珍, 付中明, 陳大福. 意大利蜜蜂4、5和6日齡幼蟲腸道發(fā)育過程中差異表達(dá)基因的趨勢(shì)分析. 上海交通大學(xué)學(xué)報(bào), 2018, 36(4): 14-21, 29.

GUO R, XIE Y L, XIONG C L, YI W X, ZHENG Y Z, FU Z M, CHEN D F. Trend analysis for differentially expressed genes in developmental process of 4-, 5- and 6-day-old larval guts of., 2018, 36(4): 14-21, 29. (in Chinese)

[25] 郭睿, 耿四海, 熊翠玲, 鄭燕珍, 付中民, 王海朋, 杜宇, 童新宇, 趙紅霞, 陳大福. 意大利蜜蜂工蜂中腸發(fā)育過程中長(zhǎng)鏈非編碼RNA的差異表達(dá)分析. 中國農(nóng)業(yè)科學(xué), 2018, 51(18): 3600-3613.

GUO R, GENG S h, XIONG C l, ZHENG Y z, FU Z m, WANG H p, DU Y, TONG X y, ZHAO H x, CHEN D f. Differential expression analysis of long non-coding RNAs during the developmental process ofworker’s midgut., 2018, 51(18): 3600-3613. (in Chinese)

[26] KIM D, PERTEA G, TRAPNELL C, PIMENTEL H, KELLEY R, SALZBERG S L. Tophat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions., 2013, 14(4): R36.

[27] The Honeybee Genome Sequencing Consortium. Insights into social insects from the genome of the honeybee., 2006, 443(7114): 931-949.

[28] WANG L K, FENG Z X, WANG X W, WANG X, ZHANG X G. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data., 2010, 26(1): 136-138.

[29] LIZ Y, HUANGC, BAO C, CHEN L, LIN M, WANG X L, ZHONG G L, YU B, HU W C, DAI L M, ZHU P F, CHANG Z X, WU Q F, ZHAO Y, JIA Y, XU P, LIU H J, SHAN G. Exon-intron circular RNAs regulate transcription in the nucleus., 2015, 22(3): 256-264.

[30] LANGMEAD B, TRAPNELL C, POP M, SALZBERG S L. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome., 2009, 10(3): R25.

[31] ALLEN E, XIE Z, GUSTAFSON A M, CARRINGTON J C. MicroRNA-directed phasing during-acting siRNA biogenesis in plants., 2005, 121(2): 207-221.

[32] SMOOT M E, ONO K, RUSCHEINSKI J, WANG P L, IDEKER T. Cytoscape 2.8: new features for data integration and network visualization., 2011, 27(3): 431-432.

[33] HU X L, ZHU M, ZHANG X, LIU B, LIANG Z, HUANG L X, XU J, YU L, LI K, ZAR M S, XUE R Y, CAO G L, GONG C L. Identification and characterization of circular RNAs in the silkworm midgut followingcytoplasmic polyhedrosis virus infection., 2018, 15(2): 292-301.

[34] ZHANG C L, WU H, WANG Y H, ZHU S Q, LIU J Q, FANG X T, CHEN H.Circular RNA of cattle casein genes are highly expressed in bovine mammary gland., 2016, 99(6): 4750-4760.

[35] HUANG M J, SHEN Y F, MAO H G, CHEN L X, CHEN J C, GUO X L, XU N G. Circular RNA expression profiles in the porcine liver of two distinct phenotype pig breeds., 2018, 31(6): 812-819.

[36] ZHANG X H, YAN Y M, LEI X Y, LI A J, ZHANG H M, DAI Z K, LI X J, CHEN W G, LIN W C, CHEN F, MA J Y, XIE Q M. Circular RNA alterations are involved in resistance to avian leukosis virus subgroup-J-induced tumor formation in chickens., 2017, 8(21): 34961-34970.

[37] FAN X Y, ZHANG X N, WU X L, GUO H S, HU Y Q, TANG F C, HUANG Y Y. Single-cell RNA-seq transcriptome analysis of linear and circular RNAs in mouse preimplantation embryos., 2015, 16: 148.

[38] HE L B, ZHANG A D, XIONG L, LI Y M, HUANG R, LIAO L J, ZHU Z Y, WANG Y P. Deep circular RNA sequencing provides insights into the mechanism underlying grass carp reovirus infection., 2017, 18(9): 1977.

[39] DARBANI B, NOEPARVAR S, BORG S. Identification of circular RNAs from the parental genes involved in multiple aspects of cellular metabolism in barley., 2016, 7: 776.

[40] SUN X Y, WANG L, DING J C, WANG Y R, WANG J S, ZHANG X Y, CHE Y L, LIU Z W, ZHANG X R, YE J Z, WANG J, SABLOK G, DENG Z P, ZHAO H W. Integrative analysis oftranscriptomics reveals intuitive splicing mechanism for circular RNA., 2016, 590(20): 3510-3516.

[41] 劉駿武, 陳玲玲. 線蟲環(huán)狀RNA分析. 計(jì)算生物學(xué), 2015, 5(2): 17-28.

LIU J W, CHEN L L. Analysis of circular RNA in., 2015, 5(2): 17-28. (in Chinese)

[42] 劉春蕾, 胥保華, 劉振國, 王穎, 王紅芳. 不同越冬飼料對(duì)蜜蜂中腸消化酶活性、組織發(fā)育狀態(tài)以及抗氧化酶基因表達(dá)的影響. 動(dòng)物營(yíng)養(yǎng)學(xué)報(bào), 2017, 29(4): 1183-1190.

LIU C L, XU B H, LIU Z G, WANG Y, WANG H F. Effects of different overwintering feeds on midgut digestive enzyme activities, tissue development status and antioxidant enzyme gene expression of honeybees., 2017, 29(4): 1183-1190. (in Chinese)

[43] 劉彩珍. 中華蜜蜂(Fabricius)中腸消化酶活性的探討[D]. 福州: 福建農(nóng)林大學(xué), 2001.

LIU C Z. Inquisition of the digestive enzyme activity in the midgut of the honeybee (Fab.)[D]. Fuzhou: Fujian Agriculture and Forestry University, 2001. (in Chinese)

[44] HALDER G, JOHNSON R L. Hippo signaling: growth control and beyond., 2011, 138(1): 9-22.

[45] Camargo F D, Gokhale S, JOHNNIDIS J B, FU D, BELL G W, Jaenisch R, Brummelkamp T R.increases organ size and expands undifferentiated progenitor cells., 2007, 17(23): 2054-2060.

[46] FEVR T, ROBINE S, LOUVARD D, HUELSKEN J. Wnt/beta- catenin is essential for intestinal homeostasis and maintenance of intestinal stem cells., 2007, 27(21): 7551-7559.

[47] 孫曉陽, 王雁玲. Wnt 信號(hào)通路與哺乳動(dòng)物生殖. 生物化學(xué)與生物物理進(jìn)展, 2003, 30(2): 180-184.

SUN X Y, WANG Y L. Wnt signaling pathways in mammalian reproduction., 2003, 30(2): 180-184. (in Chinese)

[48] BARRY E R, CAMARGO F D. The Hippo superhighway: signaling crossroads converging on the Hippo/Yap pathway in stem cells and development., 2013, 25(2): 247-253.

[49] ARONSTEIN K A, MURRAY K D. Chalkbrood disease in honey bees.2010, 103(Suppl. 1): S20-S29.

[50] SARAAV I, SINGH S, SHARMA S. Outcome of mycobacterium tuberculosis and Toll-like receptor interaction: immune response or immune evasion?, 2014, 92(9): 741-746.

[51] SALMENA L, POLISENO L, TAY Y, KATS L, PANDOLFI P P. Ahypothesis: the rosetta stone of a hidden RNA language?, 2011, 146(3): 353-358.

[52] DENG Y Y, ZHANG W P, SHE J Q, ZHANG L S, CHEN T, ZHOU J, YUAN Z Y. Circular RNA related tofunction as ceRNA of microRNA in human acute myocardial infarction., 2016, 68(16): C51-C52.

[53] XU X W, ZHENG B A, HU Z M, Qian Z Y, Huang C J, Liu X Q, Wu W D. Circular RNA hsa_circ_000984 promotes colon cancer growth and metastasis by sponging miR-106b., 2017, 8(53): 91674-91683.

[54] LIU Q, ZHANG X, HU X Q, Dai L H, Fu X, Zhang J Y, Ao Y F. Circular RNA related to the chondrocyte ECM regulatesexpression by functioning as a mir-136 ‘Sponge’ in human cartilage degradation., 2016, 6: 22572.

[55] WU Y, PARTHASARATHY R, BAI H, PALLI S R. Mechanisms of midgut remodeling: juvenile hormone analog methoprene blocks midgut metamorphosis by modulating ecdysone action., 2006, 123(7): 530-547.

[56] MELLO T R, ALEIXO A C, PINHEIRO D G, NUNES F M, BITONDI M M, HARTFELDER K, BARCHUK A R, SIMOES Z L. Developmental regulation of ecdysone receptor (EcR) and EcR-controlled gene expression during pharate-adult development of honeybees ()., 2014, 5: 445.

Analysis of differentially expressed circular RNAs and their regulation networks during the developmental process ofworker’s midgut

GUO Rui, CHEN HuaZhi, XIONG CuiLing, ZHENG YanZhen, FU ZhongMin, XU GuoJun, DU Yu, WANG HaiPeng, GENG SiHai, ZHOU DingDing, LIU SiYa, CHEN DaFu

(College of Bee Science, Fujian Agriculture and Forestry University, Fuzhou 350002)

【Objective】Circular RNA (circRNA) plays a primary role in alternative splicing, transcription regulation and expression regulation of parental gene. The objective of this study is to investigate the profile expression of circRNAs and differentially expressed circRNAs (DEcircRNAs) during the developmental process of the midguts ofworkers, and to explore the role of DEcircRNAs in the development of midgut at the transcriptional level. 【Method】 On basis of the whole transcriptome data from 7- and 10-day-old worker’s midguts of(Am7 and Am10), find_circ software was used to predict circRNAs based on the filtered sequencing data.The circRNA expression level was normalized by RPM algorithm. Differential expression analysis for circRNAs was conducted via DEGseq software following standards fold change≥2.0,＜0.05 and false discovery rate (FDR)＜0.05. Source genes of DEcircRNAs were annotated to GO and KEGG databases to gain function and pathway annotations by using BLAST. DEcircRNAs-miRNAs and DEcircRNAs-miRNAs-mRNAs networks were predicted with TargetFinder software and visualized using Cytoscape v.3.2.1 software. RT-qPCR was conducted to verify the reliability of sequencing data.【Result】 On average, 19 616 356 anchors reads were obtained from eachworker’s midgut sample.Pearson correlations between different biological repeats within Am7 and Am10 groups were ≥0.950. In total, 256 DEcircRNAs including 105 up-regulated circRNAs and 151 down-regulated circRNAs were predicted. Novel_circ_009675 and novel_circ_013879 were highly expressed in Am7 and Am10, respectively. Source genes of DEcircRNAs could be annotated to 32 GO terms including binding, single-organism process and cellular process, among them 35, 35 and 7 source genes were involved in catalytic activity, metabolic process and stress response. Additionally, these source genes could be annotated to 35 KEGG pathways, in which 5, 5 and 4 source genes were associated with Hippo signaling pathway, endocytosis and phagosome, respectively; further investigation showed that 1, 2 and 2 source genes could be annotated to material metabolisms such as phosphoinositol metabolism, starch and sucrose metabolism and galactose metabolism; 5, 4, 3, 1 and 1 source genes could be annotated to immune signaling pathways including endocytosis, phagosome, lysosome, ubiquitin-mediated proteolysis and MAPK signaling pathway, respectively.These results suggested that the corresponding DEcircRNA was involved in the development, metabolism and immune defense of the midgut of.DEcircRNA-miRNA regulation network analysis showed that 141 DEcircRNAs could link to 107 miRNAs, most of these DEcircRNAs could only bind to 1-2 miRNAs, but novel_circ_011577 and novel_circ_010719 could respectively bind to 32 and 28 miRNAs.In addition, the number of DEcircRNAs combined with mir-136-y, ame-miR-6001-3p and mir-136-y was the highest (15, 14 and 14, respectively), which indicated that the corresponding DEcircRNA could play roles during the developmental process ofworker’s midgut as competing endogenous RNAs. Furthermore, DEcircRNAs-ame-miR-6001-3p-mRNA network was constructed and analyzed, and the result indicated that 14 DEcircRNAs could jointly link to ame-miR-6001-3p, implying they were likely to indirectly regulate division and differentiation of stem cells inworker’s midgut via regulation of ame-miR-6001-3pSix DEcircRNAs were randomly selected for RT-qPCR assay, the result showed the alteration trend of expression levels of 5 DEcircRNAs was consistent with that of the sequencing data, which proved the reliability of trancriptome data.【Conclusion】Through the deep investigation of DEcircRNAs during the developmental process ofworker’s midgut, the expression profile and differential expression information of circRNAs in the development of midgut of worker bee were provided, and the role of DEcircRNAs in the development of midgut was revealed. It provides a basis for the screening and functional study of key circRNAs associated with the development of the midgut.

; midgut; circRNAs; regulation network; development

10.3864/j.issn.0578-1752.2018.23.015

2018-07-16；

2018-09-10

國家現(xiàn)代農(nóng)業(yè)產(chǎn)業(yè)技術(shù)體系建設(shè)專項(xiàng)資金（CARS-44-KXJ7）、福建省科技計(jì)劃項(xiàng)目（2018J05042）、福建省教育廳中青年教師教育科研項(xiàng)目（JAT170158）、福建農(nóng)林大學(xué)科技創(chuàng)新專項(xiàng)基金（CXZX2017342，CXZX2017343）、福建省大學(xué)生創(chuàng)新創(chuàng)業(yè)訓(xùn)練計(jì)劃（3165602032，3155006018）

郭睿，E-mail：ruiguo@fafu.edu.cn。陳華枝，E-mail：18965015689@163.com。郭睿和陳華枝為同等貢獻(xiàn)作者。

陳大福，E-mail：dfchen826@fafu.edu.cn

（責(zé)任編輯 岳梅）