家蠶感染二分濃核病毒（鎮(zhèn)江株）的數(shù)字基因表達(dá)譜分析

高 坤，尚夢珂，錢荷英，覃光星，郭錫杰

（江蘇科技大學(xué)生物技術(shù)學(xué)院/中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所，江蘇鎮(zhèn)江 212018）

家蠶感染二分濃核病毒（鎮(zhèn)江株）的數(shù)字基因表達(dá)譜分析

高 坤，尚夢珂，錢荷英，覃光星，郭錫杰

（江蘇科技大學(xué)生物技術(shù)學(xué)院/中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所，江蘇鎮(zhèn)江 212018）

【目的】篩選家蠶與二分濃核病毒（Bombyx mori bidensovirus Zhenjiang strain，BmBDV-ZJ）感染相關(guān)的差異表達(dá)基因，鑒定家蠶與病毒感染有關(guān)的調(diào)控基因，為進(jìn)一步闡明家蠶抗BmBDV-ZJ的分子機(jī)制提供理論依據(jù)?！痉椒ā坎捎肐llumina高通量測序技術(shù)，構(gòu)建家蠶品種JS口服感染BmBDV-ZJ的數(shù)字基因表達(dá)譜，為排除個(gè)體間差異，以10頭蠶作為一個(gè)樣本用于DGE檢測。樣本中基因的差異表達(dá)檢測通過嚴(yán)格的運(yùn)算法則進(jìn)行，對差異檢驗(yàn)的P值（P value）作多重假設(shè)檢驗(yàn)校正，通過控制FDR（false discovery rate）來決定P值的域值。本研究中，差異表達(dá)基因定義為FDR≤0.001且差異倍數(shù)在2倍及以上（|log2ratio|≥1）的基因。采用基因本體論（GO）分類體系確定所有差異表達(dá)基因可能的功能。用GO計(jì)算P值和bonferoni校正。選用校正P值≤0.05作為基因組顯著富集的閾值。WEGO軟件用來視化、比較和繪制GO注釋結(jié)果。利用KEGG數(shù)據(jù)庫進(jìn)行通路富集分析,進(jìn)一步確定顯著富集代謝途徑或信號傳導(dǎo)途徑，Q值≤0.05的通路指定為DGEs中的顯著富集通路。通過qRT-PCR方法對部分差異表達(dá)基因進(jìn)行驗(yàn)證?！窘Y(jié)果】感染組和對照組分別得到4 850 663和4 875 307個(gè)原始標(biāo)簽,去除低質(zhì)量標(biāo)簽后，分別得到4 757 934和4 788 406個(gè)清潔標(biāo)簽，對應(yīng)的標(biāo)簽種類數(shù)量分別為62 436和63 680種。兩個(gè)文庫間的清潔標(biāo)簽和清潔標(biāo)簽種類的數(shù)量在不同拷貝區(qū)間分布類似，感染組和對照組樣本的測序量分別為3.5 M和3.7 M，測序深度符合試驗(yàn)的要求，兩樣本的DGE數(shù)據(jù)是可信的。將這兩個(gè)DGE數(shù)據(jù)庫的所有清潔標(biāo)簽與家蠶參考基因庫進(jìn)行比對，在對照組與感染組中，分別有36.39%和45.30%的清潔標(biāo)簽可以比對到基因。另有50.02%和43.34%的清潔標(biāo)簽可以比對到家蠶參考基因組，剩余的未知標(biāo)簽分別占清潔標(biāo)簽總數(shù)的13.59%和12.35%。共發(fā)現(xiàn)了447個(gè)差異表達(dá)基因，其中306個(gè)上調(diào)表達(dá)，141個(gè)下調(diào)表達(dá)。分別有218、147和179個(gè)差異表達(dá)基因涉及GO 3個(gè)本體中的分子功能、細(xì)胞組分和生物過程。利用KEGG公共數(shù)據(jù)庫進(jìn)行Pathway顯著性富集分析，注釋到的基因總數(shù)為8 473個(gè)。447個(gè)差異表達(dá)基因經(jīng)鑒定后，其中的330個(gè)基因被歸類到151個(gè)KEGG路徑中。差異表達(dá)基因顯著富集的Pathway（Q值≤0.05）有19個(gè)，其中最顯著富集的是細(xì)胞質(zhì)中DNA識別通路。挑選了24個(gè)差異表達(dá)基因進(jìn)行qRT-PCR驗(yàn)證，其中20個(gè)基因的差異表達(dá)趨勢與DGE的結(jié)果一致。其中在DNA識別通路中共檢測到9個(gè)差異表達(dá)基因，BGIBMGA009408-TA、BGIBMGA004913-TA、BGIBMGA011753-TA均為編碼RNA聚合酶III的基因，表達(dá)量均上調(diào)，是對照組的4.3、2.3、3.4倍。【結(jié)論】構(gòu)建了3齡家蠶JS感染BmBDV-ZJ后28 h感染組及對照組幼蟲的數(shù)字基因表達(dá)譜，Pathway顯著性富集分析和qRT-PCR驗(yàn)證顯示，家蠶感染BmBDV-ZJ后可能通過啟動胞質(zhì)內(nèi)DNA識別通路來感應(yīng)入侵病毒的異源DNA成分并迅速啟動天然免疫抵御BmBDV-ZJ病毒感染,為研究BmBDV-ZJ侵染家蠶和家蠶抵御BmBDV感染的分子機(jī)制打下了基礎(chǔ)。

家蠶；數(shù)字基因表達(dá)譜；二分濃核病毒；qRT-PCR

0 引言

【研究意義】家蠶二分濃核病毒鎮(zhèn)江株（Bombyx mori bidensovirus Zhenjiang strain，BmBDV-ZJ）是感染家蠶的一種重要病毒，主要通過食下感染，引起家蠶的濃核病?；疾⌒Q呈現(xiàn)空頭，下痢，吐液等癥狀，是嚴(yán)重影響蠶桑生產(chǎn)的病毒病之一，每年給養(yǎng)蠶業(yè)造成較大的經(jīng)濟(jì)損失。目前，對于該病毒感染家蠶的致病機(jī)理和家蠶抵御該病毒感染的應(yīng)答機(jī)制均不清楚。篩選鑒定家蠶與該病毒感染相關(guān)的差異表達(dá)基因，進(jìn)一步研究病毒對家蠶致病的分子機(jī)理，對家蠶抗病毒育種和家蠶濃核病的有效防治具有重要意義?！厩叭搜芯窟M(jìn)展】中國早在1959年就證實(shí)了生產(chǎn)上所發(fā)生的家蠶空頭性軟化病是由病毒引起，隨后各國科學(xué)家分離得到了該病毒的不同株系，分別命名為伊那株（DNV-1）[1]、佐久株（DNV-2）[2-3]、中國（鎮(zhèn)江）株（DNV-3）[4]。之前認(rèn)為該病毒均屬于細(xì)小病毒科（Parvoriridae）濃核病毒屬（Densovirus，DNV）的蠶濃核病毒（B. mori densonucleosis virus）；后來發(fā)現(xiàn)DNV-2和DNV-3基因組包含兩個(gè)節(jié)段的單鏈DNA且不采用滾環(huán)復(fù)制（細(xì)小病毒科常用的復(fù)制方式）[5]，因此將DNV-2和DNV-3重新劃分為一個(gè)新設(shè)定的二分DNA病毒科二分濃核病毒屬（Bidensovirus，BDV）[6]。本研究所用的病毒株是1981年由中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所（鎮(zhèn)江）錢元駿等分離得到的家蠶二分濃核病毒（鎮(zhèn)江）株（BmBDV-ZJ），與伊那株（DNV-1）在蠶品種感受性、病毒的物理化學(xué)性狀、血清學(xué)特性等方面明顯不同[7]。BDV-ZJ主要在家蠶幼蟲中腸柱狀上皮細(xì)胞的細(xì)胞核中復(fù)制，引起幼蟲的濃核病，不形成多角體[8]。其病毒基因組包括兩段不同的DNA單鏈（6 543和6 022 bp），一個(gè)病毒粒子只包含其中一條單鏈DNA，因此是兩種不同病毒粒子的混合物[9]。成熟的病毒粒子都會隨感染細(xì)胞的破裂而釋放，進(jìn)而感染鄰近細(xì)胞。大部分家蠶品種對BmBDV-ZJ是易感的，錢元駿等[10]通過對380多個(gè)家蠶品種的抗BmBDV-ZJ性能比較研究和不同抗性品種進(jìn)行雜交試驗(yàn)證實(shí)家蠶對BmBDV-ZJ抗性為隱性遺傳，同時(shí)也受微效多基因的影響。有的品種即使接種高濃度的病毒也完全不發(fā)病。日本科學(xué)家[11]研究發(fā)現(xiàn)了一個(gè)與DNV-2易感性有關(guān)的家蠶基因nsd-2，抗性品種中該基因開放閱讀框中大約6 000個(gè)堿基缺失，導(dǎo)致其編碼的具有12個(gè)跨膜域的氨基酸轉(zhuǎn)運(yùn)膜蛋白的缺失，且該蛋白僅在家蠶的中腸中表達(dá)；而通過轉(zhuǎn)基因技術(shù)修復(fù)該基因的缺失后可以使抗性品種對BmBDV-2病毒易感，該研究結(jié)果表明病毒識別的位點(diǎn)可能位于nsd-2基因缺失的膜蛋白部分，然而病毒是如何與這種膜蛋白相互作用仍不清楚。裘智勇等[12-13]研究表明，家蠶對BmBDV-ZJ的抵抗性可能與家蠶中腸細(xì)胞表面的某些特殊蛋白因子有關(guān)。BAO等[14]通過抑制消減雜交技術(shù)，研究了兩個(gè)不同抗性家蠶品種JS（對BmBDV-ZJ感受性）和NIL（對BmBDV-ZJ非感受性）感染BmBDV-ZJ后基因的相對表達(dá)變化，發(fā)現(xiàn)了在抗性品種中有11個(gè)基因明顯上調(diào)，推測這些基因與NIL的抗BmBDV-ZJ感染復(fù)制有關(guān)?！颈狙芯壳腥朦c(diǎn)】病毒生活史的每一步都受病毒和宿主間分子相互作用的調(diào)控，通過家蠶與BmBDV-ZJ病毒感染相關(guān)的基因差異表達(dá)分析，鑒定可能與病毒復(fù)制有關(guān)的宿主細(xì)胞內(nèi)調(diào)控因子，有助于闡明家蠶抗BmBDV-ZJ的抗病毒機(jī)制?！緮M解決的關(guān)鍵問題】通過對易感家蠶品種JS的差異基因表達(dá)譜和差異基因可能參與的信號通路分析，篩選和尋找更多可能與病毒感染相關(guān)的差異表達(dá)基因。

1 材料與方法

試驗(yàn)于2013年9月至2015年9月在中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所病理研究室完成。

1.1 供試材料與試劑

1.1.1 家蠶及病毒株 家蠶品種JS由中國農(nóng)業(yè)科學(xué)院國家蠶種質(zhì)資源保存中心提供；家蠶二分濃核病毒（BmBDV-ZJ）由中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所家蠶病理研究室繁殖保存；BmBDV-ZJ抗血清由中國農(nóng)業(yè)科學(xué)院蠶業(yè)研究所家蠶病理研究室制備保存。

1.1.2 主要試劑 PrimeScriptTMRT Reagent Kit、

SYBR Premix Ex TaqTMKit購自TaKaRa公司；總RNA提取試劑盒、Trizol、M-MLV反轉(zhuǎn)錄酶購自Invitrogen公司；DNA marker DL 2000購自TaKaRa公司；其他均為國產(chǎn)或進(jìn)口分析純試劑。

1.2 試驗(yàn)方法

1.2.1 家蠶幼蟲添毒感染 取二分濃核病毒感染蠶的中腸干粉0.4 g，加入8 mL 2% Na2CO3溶液研磨勻漿，靜止4—5 min，加5 mL蒸餾水至50 mL，3 000 r/min離心20 min，取上清加3倍體積0.2 mol·L-1醋酸，處理30 min，用Na2CO3調(diào)pH值至6—7，用蒸餾水稀釋至1/1 000備用，過濾除菌備用。

試驗(yàn)用家蠶幼蟲JS在標(biāo)準(zhǔn)溫度和光照條件下飼養(yǎng)至3齡起蠶。計(jì)數(shù)40頭蠶為感染組（添食BmBDV-ZJ），另40頭為對照組（添食滅菌水）。將制備病毒懸浮液均勻涂到桑葉上喂食感染組家蠶，對照組喂食同樣體積的滅菌水涂過的桑葉，12 h后全部喂飼正常桑葉?？紤]到個(gè)體間差異可能帶來的影響，每10條蠶收集一管提取RNA，共設(shè)置3組生物學(xué)重復(fù)。最后剩余10頭蠶用來檢測病毒感染發(fā)病情況。

1.2.2 病毒感染的確認(rèn) 在一塊7 cm×10 cm的瓊脂板（1%瓊脂糖）上打7個(gè)孔，中間孔加入抗血清，周圍6個(gè)孔加待測樣品，將瓊脂板放入補(bǔ)濕飯盒中，在20—30℃條件下孵育1 d觀察有無沉淀帶出現(xiàn)。

1.2.3 數(shù)字基因表達(dá)譜（digital gene expression, DGE）測序及分析 由華大基因公司提供技術(shù)支持，對28 h時(shí)間點(diǎn)的感染組及對照組JS家蠶幼蟲各取10 μg總RNA，經(jīng)過純化、反轉(zhuǎn)錄、酶切和PCR線性擴(kuò)增后，使用Illumina HiSeqTM2000進(jìn)行測序。差異表達(dá)基因的篩選及其Gene Ontology（GO）、Pathway顯著性富集分析參照文獻(xiàn)[15-16]進(jìn)行。

1.2.4 qRT-PCR分析 通過實(shí)時(shí)熒光定量PCR（qRTPCR）對DGE篩選得到的部分差異表達(dá)基因進(jìn)行驗(yàn)證。感染組和對照組各3管樣品，以每管樣品提取的總RNA作為模板，用PrimeScriptTMRT Reagent Kit（TaKaRa）進(jìn)行反轉(zhuǎn)錄合成第一鏈cDNA。以家蠶Actin 3作為內(nèi)參基因，使用SYBR? Premix Ex TaqTM（TaKaRa）試劑盒在ABI PRISM? 7300 檢測系統(tǒng)上進(jìn)行實(shí)時(shí)定量PCR，每個(gè)模板做3次重復(fù)。反應(yīng)條件：94℃ 10 min變性，94℃ 15 s，60℃ 31 s，40個(gè)循環(huán)。數(shù)據(jù)分析參考文獻(xiàn)[17]進(jìn)行，引物序列見表1。

2 結(jié)果

2.1 病毒感染的確定

通過雙向免疫擴(kuò)散法診斷病毒的感染。在BmBDV-ZJ感染28 h后，所有接種感染的家蠶幼蟲其中腸研磨液與血清孔之間均出現(xiàn)了白色沉淀帶，而與健康蠶（對照組）中腸研磨液及BmCPV、BmNPV感染蠶研磨液均不起反應(yīng)，表明感染組家蠶幼蟲全部被BmBDV-ZJ感染。

2.2 DGE數(shù)據(jù)庫分析

采用Illumina高通量測序技術(shù)，感染組和對照組分別得到4 850 663和4 875 307個(gè)原始標(biāo)簽，去除低質(zhì)量標(biāo)簽后，分別得到4 757 934和4 788 406個(gè)清潔標(biāo)簽，對應(yīng)的標(biāo)簽種類數(shù)量分別為62 436和63 680種（表2）。兩個(gè)文庫間的清潔標(biāo)簽和清潔標(biāo)簽種類的數(shù)量在不同拷貝區(qū)間分布類似（圖1），感染組和對照組樣本的測序量分別為3.5 M和3.7 M（圖2），測序深度符合試驗(yàn)的要求，兩樣本的DGE數(shù)據(jù)是可信的。

圖1 對照組和感染組的清潔標(biāo)簽和清潔標(biāo)簽種類在不同拷貝區(qū)的分布Fig. 1 Distribution of total clean tags and distinct clean tags in each library

表1 定量PCR引物Table 1 Primers for qRT-PCR

2.3 標(biāo)簽比對及標(biāo)準(zhǔn)化處理

將所有清潔標(biāo)簽與家蠶參考基因庫進(jìn)行比對，在對照組與感染組中，分別有36.39%和45.30%的清潔標(biāo)簽可以比對到基因，占總清潔標(biāo)簽種類數(shù)的30.56%和37.41%。另有50.02%和43.34%的清潔標(biāo)簽可以比對到家蠶參考基因組，剩余的未知標(biāo)簽分別占清潔標(biāo)簽總數(shù)的13.59%和12.35%（表2）。

表2 兩個(gè)樣本中標(biāo)簽分布表Table 2 Distributions of tags between the two libraries

圖2 數(shù)字基因表達(dá)譜的測序飽和趨勢圖Fig. 2 Trends of saturation of DGEs

2.4 差異基因分析及其驗(yàn)證

篩選了發(fā)現(xiàn)錯(cuò)誤率（false discovery rate）FDR≤0.001且差異倍數(shù)在2倍及以上（|log2ratio|≥1）的基因作為差異表達(dá)基因，共發(fā)現(xiàn)差異表達(dá)基因447個(gè)。在感染組中，上調(diào)表達(dá)基因306個(gè)，下調(diào)表達(dá)基因141個(gè)。為了驗(yàn)證DGE數(shù)據(jù)的準(zhǔn)確性，挑選了24個(gè)差異表達(dá)基因進(jìn)行RT-PCR驗(yàn)證，其中20個(gè)基因表達(dá)情況與DGE的結(jié)果一致（圖3），說明DGE數(shù)據(jù)基本可以準(zhǔn)確反映樣本的基因表達(dá)情況。

2.5 GO和KEGG分析

GO共有3個(gè)本體，分別描述基因的分子功能、所處的細(xì)胞組分、參與的生物過程。GO結(jié)果顯示分別有147、179、218個(gè)差異基因分配在這3個(gè)本體中，其中在核糖核蛋白復(fù)合體、核糖體、核糖體亞基、大分子復(fù)合物、細(xì)胞內(nèi)組分、蛋白酶體調(diào)控組分、代謝過程、基因表達(dá)和結(jié)構(gòu)分子活性中基因富集程度較高，校正后的P值≤0.05（表3）。

為了進(jìn)一步了解基因的生物學(xué)功能，利用KEGG公共數(shù)據(jù)庫進(jìn)行Pathway顯著性富集分析，注釋到的基因總數(shù)為8 473個(gè)，447個(gè)差異表達(dá)基因經(jīng)鑒定后，其中的330個(gè)基因被歸類到151個(gè)KEGG路徑中。定義為在差異表達(dá)基因中顯著富集的Pathway（Q值≤0.05）有19個(gè)（表4），其中最顯著富集的是細(xì)胞質(zhì)中DNA識別通路，其次是RNA聚合酶、核糖體、嘧啶代謝、蛋白酶體、抗原加工與呈遞、代謝途徑等。

圖3 候選基因在對照和感染BmBDV-ZJ的家蠶幼蟲中的差異表達(dá)Fig. 3 Differential expression levels of candidate genes in BmBDV-ZJ-infected and control B. mori larvae

3 討論

為了研究BmBDV-ZJ感染家蠶的可能分子機(jī)制和宿主的應(yīng)答反應(yīng)，基于本研究構(gòu)建JS家蠶感染組和對照組的兩個(gè)DGE文庫分析，篩選得到447個(gè)差異表達(dá)基因，其中感染組有306個(gè)上調(diào)表達(dá)基因和141個(gè)下調(diào)表達(dá)基因，該差異基因的上調(diào)和下調(diào)對比結(jié)果與之前家蠶感染質(zhì)型多角體病毒的DGE結(jié)果類似，都是上調(diào)表達(dá)的基因數(shù)目明顯高于下調(diào)表達(dá)的基因數(shù)目[15-16,18]。相關(guān)研究表明，家蠶在感染大部分病原微生物時(shí)的差異表達(dá)基因都有上調(diào)基因明顯多于下調(diào)基因的趨勢，LIU等[19]通過基因芯片研究發(fā)現(xiàn)家蠶感染革蘭氏陽性菌（Serratia marcescens）后有172個(gè)基因上調(diào)表達(dá)，61個(gè)基因下調(diào)表達(dá)；感染革蘭氏陰性菌（Staphylococcus aureus）后125個(gè)基因上調(diào)表達(dá)，104個(gè)基因下調(diào)表達(dá)；感染真菌（Beauveria bassiana）后有133個(gè)基因上調(diào)表達(dá)，24個(gè)基因下調(diào)表達(dá)；而被PBS處理后有201個(gè)基因上調(diào)表達(dá)和40個(gè)基因下調(diào)表達(dá)?？梢姛o論是病毒感染、細(xì)菌感染還是真菌感染都會引起家蠶更多基因的上調(diào)表達(dá)來應(yīng)對病原的入侵。

表3 基因本體分類Table 3 Gene Ontology classification

表4 顯著富集的PathwayTable 4 Significantly enriched pathways

利用KEGG公共數(shù)據(jù)庫進(jìn)行Pathway顯著性富集分析，其中最顯著富集的是細(xì)胞質(zhì)中DNA識別通路（cytosolic DNA-sensing pathway），該通路可以有效識別入侵病毒的異源DNA成分并迅速啟動天然免疫，以及隨后的特異性免疫應(yīng)答來對病毒進(jìn)行清除，是機(jī)體抵抗病毒感染的重要機(jī)制[20-21]。DNA感受器（DNA sensor）是宿主感受外源入侵DNA和免疫防御的橋梁，可以特異地識別病毒等外源DNA進(jìn)而激活下游的免疫信號途徑，并通過誘導(dǎo)表達(dá)抗病毒蛋白來抑制病毒的復(fù)制并向周圍的細(xì)胞示警[22-23]。目前已經(jīng)有超過10種DNA感受器被發(fā)現(xiàn)。例如，Toll樣受體9（TLR9）家族[24]、視黃酸誘導(dǎo)基因蛋白RIG-I相關(guān)受體（RLRs）家族[25]、DNA依賴的RNA聚合酶III[26]、干擾素誘導(dǎo)蛋白16（IFI16）[27]；DExD-H框解旋酶超家族[28-30]等。RNA聚合酶III可以識別侵入細(xì)胞內(nèi)的DNA病毒，然后以病毒DNA為模板，將病毒信息轉(zhuǎn)錄合成一種5′三磷酸化的特殊的雙鏈RNA，進(jìn)而被RIG-I分子所識別。RIG-I可以激活I(lǐng)型干擾素產(chǎn)生，也可以激活NF-κB等信號通路并誘導(dǎo)抗病毒基因的表達(dá)，從而抑制病毒等病原體復(fù)制[31-32]。家蠶感染BmBDV-ZJ后在DNA識別通路中共檢測到9個(gè)差異表達(dá)基因，其中BGIBMGA009408-TA、BGIBMGA004913-TA、BGIBMGA011753-TA均為編碼RNA聚合酶III的基因，表達(dá)量均上調(diào)，是對照組的4.3、2.3、3.4倍（圖2）。其功能可能作為家蠶細(xì)胞質(zhì)中的DNA感受器，通過識別病毒DNA激活家蠶抗BmBDV-ZJ感染的免疫應(yīng)答反應(yīng)。

干擾素-γ誘導(dǎo)的溶酶體巰基還原酶（IFN-γ inducible lysosomal thiol reductase，GILT）在脊椎動物適應(yīng)性免疫中的MHCII類抗原加工和呈遞過程中起著關(guān)鍵作用，可以催化未折疊的天然抗原蛋白二硫鍵的斷裂，進(jìn)而對其進(jìn)行酶解加工[33]。雖然適應(yīng)性免疫和干擾素的產(chǎn)生已明確只在脊椎動物中存在，但是GILT基因家族在脊椎動物和無脊椎動物中普遍存在，如在蝦[34]、蟹[35]、鮑[36]、果蠅[37]和線蟲[38]中均發(fā)現(xiàn)了該基因。無脊椎動物中的GILT基因在感染細(xì)菌和病毒后也都出現(xiàn)上調(diào)表達(dá)趨勢，推測其功能可能不同于脊椎動物中的抗原加工和呈遞，而是參與了無脊椎動物中的某種先天性免疫信號途徑。家蠶感染BmBDV-ZJ病毒后基因BGIBMGA003741-TA和BGIBMGA008360-TA均上調(diào)表達(dá)，且這2個(gè)基因都含有GILT保守結(jié)構(gòu)域，推測為家蠶中的GILT基因。其功能是否類似于接頭蛋白干擾素刺激基因，其上調(diào)表達(dá)是否可以進(jìn)一步誘導(dǎo)相關(guān)抗病毒蛋白的產(chǎn)生還需要進(jìn)一步試驗(yàn)驗(yàn)證。

BGIBMGA013054-TA作為缺血再灌注誘導(dǎo)蛋白的同源基因，其主要功能是調(diào)節(jié)細(xì)胞內(nèi)物質(zhì)裝配和相關(guān)轉(zhuǎn)運(yùn)蛋白的活性[39]。病毒感染家蠶后，會利用宿主細(xì)胞的蛋白合成系統(tǒng)來進(jìn)行病毒蛋白的復(fù)制，這些異源蛋白的產(chǎn)生，加速了細(xì)胞的物質(zhì)運(yùn)輸過程，以增加宿主細(xì)胞對病毒蛋白的多重耐受性。家蠶感染BmBDV-ZJ后該基因的上調(diào)表達(dá)可能與細(xì)胞內(nèi)一些蛋白的定位有關(guān)。

BmBDV-ZJ通常引起的是慢性病，感染的中腸上皮細(xì)胞不像感染BmBDV-1的中腸上皮細(xì)胞那樣容易脫落，而是通過增加細(xì)胞數(shù)目，引起上皮組織折疊，在10—20 d間死亡，少數(shù)幼蟲可以化蛹。因此本研究分析發(fā)現(xiàn)家蠶感染BmBDV-ZJ后，很多與細(xì)胞增殖相關(guān)的轉(zhuǎn)錄因子和蛋白都發(fā)生了明顯上調(diào)，如叉頭框轉(zhuǎn)錄因子G1（BGIBMGA010297-TA，forkfead transcription factor G1）、富含亮氨酸重復(fù)序列蛋白47（BGIBMGA013072-TA，leucine-rich repeat-containing protein 47-like protein）、賴氨酸去甲基化和組氨酸脫氫的雙功能酶（BGIBMGA014563-TA，bifunctional lysine-specific demethylase and histidyl-hydroxylase NO66）等基因。初步推測家蠶可以通過對這些基因的上調(diào)表達(dá)促進(jìn)細(xì)胞增殖來增加細(xì)胞的數(shù)目，進(jìn)而替代因感染BmBDV-ZJ而失去正常生理功能或凋亡的中腸柱狀細(xì)胞。

BGIBMGA007315-TA是家蠶的泛素特異蛋白酶7（Ubiquitin specific protease 7，USP7），定量結(jié)果顯示為下調(diào)，與DGE結(jié)果相反。該基因編碼的一種泛素化酶，可以水解Mdm2（p53的E3泛素連接酶）對p53蛋白進(jìn)行去泛素化，保護(hù)p53不被S26蛋白酶體降解，進(jìn)而調(diào)節(jié)p53與Mdm2的穩(wěn)定性[40]。相關(guān)研究表明，EB病毒的核抗原1蛋白可以通過與USP7結(jié)合，破壞p53的穩(wěn)定性，有利于病毒的潛伏感染[41]。BmBDV-ZJ感染家蠶后病毒是否會與家蠶的USP7相互作用，USP7下調(diào)表達(dá)是否也會影響到家蠶p53蛋白的穩(wěn)定性還需要進(jìn)一步驗(yàn)證。

綜上所述，差異表達(dá)基因最顯著富集的是細(xì)胞質(zhì)DNA識別通路，該通路涉及的基因BGIBMGA009408-TA、BGIBMGA004913-TA、BGIBMGA011753-TA、BGIBMGA003741-TA、BGIBMGA008360-TA和BGIBMGA013054-TA將進(jìn)一步進(jìn)行功能驗(yàn)證，明確DNA識別通路在家蠶抗BmBDV-ZJ病毒感染中的先天免疫機(jī)制。其他大部分差異基因功能未知，因此后期主要工作是對這些未知功能的差異表達(dá)基因做進(jìn)一步鑒定及其參與的抗病毒機(jī)制研究。

4 結(jié)論

構(gòu)建了3齡家蠶JS感染BmBDV-ZJ后28 h感染組及對照組幼蟲的數(shù)字基因表達(dá)譜，Pathway顯著性富集分析和qRT-PCR驗(yàn)證顯示，家蠶感染BmBDV-ZJ后可能通過啟動胞質(zhì)內(nèi)DNA識別通路來感應(yīng)入侵病毒的異源DNA成分并迅速啟動天然免疫抵御BmBDV-ZJ病毒感染。

[1] BANDO H, CHOI H, ITO Y, KAWASE S. Terminal structure of a densovirus implies a hairpin transfer replication which is similar to the model for AAV. Virology, 1990, 179(1): 57-63.

[2] BANDO H, CHOI H, ITO Y, NAKAGAKI M, KAWASE S. Structural analysis on the single-stranded genomic DNAs of the virus newly isolated from silkworm: the DNA molecules share a common terminal sequence. Archives of Virology, 1992, 124(1/2): 187-193.

[3] BANDO H, HAYAKAWA T, ASANO S, SAHARA K, NAKAGAKI M, IIZUKA T. Analysis of the genetic information of a DNA segment of a new virus from silkworm. Archives of Virology, 1995, 140(6): 1147-1155.

[4] IWASHITA Y, CHUN C Y. The development of a densonucleosis virus isolated from silkworm larvae, Bombyx mori, of China//AKAI H, KING R C, MOROHOSHI S. The Ultrastructure and Functioning of Insect Cell. Tokyo: Society for Insect Cells Japan, 1982: 161-164.

[5] HAYAKAWA T, ASANO S, SAHARA K, IIZUKA T, BANDO H. Detection of replicative intermediate with closed terminus of Bombyx densonucleosis virus. Archives of Virology,1997, 142: 393-399.

[6] ADAMS M J, CARSTENS E B. Ratification vote on taxonomic proposals to the International Committee on Taxonomy of Viruses. Archives of Virology, 2012, 157: 1411-1422.

[7] 錢元駿, 郭錫杰, 胡雪芳, 黃可威, 渡部仁. 我國和日本家蠶DNV的血清學(xué)關(guān)系. 蠶業(yè)科學(xué), 1985, 11(4): 241-242. QIAN Y J, GUO X J, HU X F, HUANG K W, WATANABE H. The serological relationship between China isolate densonucleosis virus and Japan isolate densonucleosis virus. Acta Sericologica Sinica, 1985, 11(4): 241-242. (in Chinese)

[8] 郭錫杰, 錢元駿, 胡雪芳, 王紅林. 我國家蠶濃核病毒(DNV)寄生部位研究. 蠶業(yè)科學(xué), 1985, 11(2): 93-98. GUO X J, QIAN Y J, HU X F, WANG H L. Studies on locations of Bombyx mori densonucleosis virus (China isolate) invasion. Acta Sericologica Sinica, 1985, 11(2): 93-98. (in Chinese)

[9] WANG Y J, YAO Q, CHEN K P, WANG Y, LU J, HAN X. Characterization of the genome structure of Bombyx mori densovirus (China isolate). Virus Genes, 2007, 35: 103-108.

[10] 錢元駿, 胡雪芳, 孫玉昆, 戴仁鳴. 家蠶濃核病毒的研究. 蠶業(yè)科學(xué), 1986, 12(2): 89-94. QIAN Y J, HU X F, SUN Y K, DAI R M. Studies on Bombyx mori densonucleosis virus. Acta Sericologica Sinica, 1986, 12(2): 89-94. (in Chinese)

[11] ITO K, KIDOKORO K, SEZUTSU H, NOHATA J, YAMAMOTO K, UCHINO K, KALYEBI A, EGUCHI R, HARA W, TAMURA T, KATSUMA S, MITA K, KADONO-OKUDA K. Deletion of a gene encoding an amino acid transporter in the midgut membrane causes resistance to a Bombyx parvo-like virus. Proceedings of the National Academy of Sciences of the United Dtates of America, 2008, 105(21): 7523-7527.

[12] 裘智勇, 李木旺, 沈興家, 郭錫杰. 家蠶對濃核病毒(鎮(zhèn)江株)抵抗性和感受性品種的中腸組織蛋白比較分析. 蠶業(yè)科學(xué), 2008, 34(2): 244-249. QIU Z Y, LI M W, SHEN X J, GUO X J. Comparative analysis of proteins extracted from midgut of silkworm strains susceptible and non-susceptible to Bomby mori densovirus (Zhenjiang strain). Acta Sericologica Sinica, 2008, 34(2): 244-249. (in Chinese)

[13] 裘智勇, 李木旺, 覃光星, 劉挺, 沈興家, 郭錫杰. 家蠶對濃核病毒中國鎮(zhèn)江株抵抗性機(jī)制的初步研究. 蠶業(yè)科學(xué), 2007, 33(4): 596-601. QIU Z Y, LI M W, QIN G X, LIU T, SHEN X J, GUO X J. Primary studies on mechanism of silkworm (Bombyx mori) resistance to densovirus China (Zhenjiang) strain. Acta Sericologica Sinica, 2007, 33(4): 596-601. (in Chinese)

[14] BAO Y Y, LI M W, ZHAO Y P, GE J Q, WANG C S, HUANG Y P, ZHANG C X. Differentially expressed genes in resistant and susceptible Bombyx mori strains infected with a densonucleosis virus. Insect Biochemistry and Molecular Biology, 2008, 38(9): 853-861.

[15] GAO K, DENG X Y, QIAN H Y, QIN G X, HOU C X, GUO X J. Cytoplasmic polyhedrosis virus-induced differential gene expression in two silkworm strains of different susceptibility. Gene, 2014, 539: 230-237.

[16] GAO K, DENG X Y, QIAN H Y, QIN G X, GUO X J. Digital gene expression analysis in the midgut of 4008 silkworm strain infected with cytoplasmic polyhedrosis virus. Journal of Invertebrate Pathology, 2014, 115(1): 8-13.

[17] 高坤, 鄧祥元, 裘智勇, 覃光星, 郭錫杰. 家蠶感染質(zhì)型多角體病毒 (BmCPV)后中腸組織差異蛋白質(zhì)分析. 中國農(nóng)業(yè)科學(xué), 2013, 46(13): 2796-2807. GAO K, DENG X Y, QIU Z Y, QIN G X, GUO X J. Comparative analysis of differential proteins from midgut of silkworm induced by cytoplasmic polyhedrosis virus infection. Scientia Agricultura Sinica, 2013, 46(13): 2796-2807. (in Chinese)

[18] GUO R, WANG S M, XUE R Y, CAO G L, HU X L, HUANG M L, ZHANG Y Q, LU Y H, ZHU L Y, CHEN F, LIANG Z, KUANG S L, GONG C L. The gene expression profile of resistant and susceptible Bombyx mori strains reveals cypovirus-associated variations in host gene transcript levels. Applied Microbiology and Biotechnology, 2015, 99: 5175-5187.

[19] LIU F, LING E, WU S. Gene expression profiling during early response to injury and microbial challenges in the silkworm, Bombyx mori. Archives of Insect Biochemistry and Physiology, 2009, 72(1): 16-33.

[20] YANAI H, SAVITSKY D, TAMURA T, TANIGUCHI T. Regulation of the cytosolic DNA-sensing system in innate immunity: a current view. Current Opinin in Immunology, 2009, 21(1): 17-22.

[21] MANSUR D S, SMITH G L, FERGUSON B J. Intracellular sensing of viral DNA by the innate immune system. Microbes and Infection, 2014, 16(12): 1002-1012.

[22] RATHINAM V A, FITZGERALD K A. Innate immune sensing of DNA viruses. Virology, 2011, 411(2): 153-162.

[23] 邢雅玲, 鄭洋, 王凱, 陳曉娟, 陳忠斌. 病原DNA識別及其誘導(dǎo)天然免疫調(diào)節(jié)機(jī)制研究進(jìn)展. 生物化學(xué)與生物物理進(jìn)展, 2011, 38(12): 1099-1105. XING Y L, ZHENG Y, WANG K, CHEN X J, CHEN Z B. The cellular recognition of pathogenic DNA and the related regulation of innate immunity. Progress in Biochemistry and Biophysics, 2011, 38(12): 1099-1105. (in Chinese)

[24] EBIHARA N, CHEN L, TOKURA T, USHIO H, IWATSU M, MURAKAMI A. Distinct functions between toll-like receptors 3 and 9 in retinal pigment epithelial cells. Ophthalmic Research, 2007, 39(3): 155-163.

[25] CHOI M K, WANG Z C, BAN T, YANAI H, LU Y, KOSHIBA R, NAKAIMA Y, HANGAI S, SAVITSKY D, NAKASATO M, NEGISHI H, TAKEUCHI O, HONDA K, AKIRA S, TAMURA T, TANIGUCHI T. A selective contribution of the RIG-I-like receptor pathway to type I interferon responses activated by cytosolic DNA. Proceedings of the National Academy of Sciences of the United States of America, 2009, 106(42): 17870-17875.

[26] CHIU Y H, MACMILLAN J B, CHEN Z J. RNA polymerase III detects cytosolic DNA and induces type I interferons through the RIG-I pathway. Cell, 2009, 138(3): 576-591.

[27] VEERANKI S, CHOUBEY D. Interferon-inducible p200-family protein IFI16, an innate immune sensor for cytosolic and nuclear double-stranded DNA: regulation of subcellular localization. Molecular Immunology, 2012, 49(4): 567-571.

[28] ZHANG Z Q, YUAN B, BAO M S, LU N, KIM T, LIU Y J. The helicase DDX41 senses intracellular DNA mediated by the adaptor STING in dendritic cells. Nature Immunology, 2011, 12(10): 959-965.

[29] KIM T, PAZHOOR S, BAO M S, ZHANG Z Q, HANABUCHI S, FACCHINETTI V, BOVER L, PLUMAS J, CHAPEROT L, QIN J, LIU Y J. Aspartate-glutamate-alanine-histidine box motif (DEAH)/ RNA helicase A helicases sense microbial DNA in human plasmacytoid dendritic cells. Proceedings of the National Academy of Sciences of the United States of America, 2010, 107(34): 15181-15186.

[30] ZHANG X, BRANN T W, ZHOU M, YANG J, OGUARIRI R M, LIDIE K B, IMAMICHI H, HUANG D W, LEMPICKI R A, BASELER M W, VEENSTRA T D, YOUNG H A, LANE H C, IMAMICHI T. Ku70 is a novel cytosolic DNA sensor that induces type III rather than type I IFN. Journal of Immunology, 2011, 186(8): 4541-4545.

[31] ABLASSER A, BAUERNFEIND F, HARTMANN G, LATZ E, FITZGERALD K A, HORNUNG V. RIG-I dependent sensing of poly (dA-dT) through the induction of an RNA polymerase III transcribed RNA intermediate. Nature Immunology, 2009, 10(10): 1065-1072.

[32] MELCHJORSEN J, RINTAHAKA J, S?BY S, HORAN K A, POLTAJAINEN A, ?STERGAARD L, PALUDAN S R, MATIKAINEN S. Early innate recognition of herpes simplex virus in human primary macrophages is mediated via the MDA5/MAVS- dependent and MDA5/MAVS/RNA polymerase III-independent pathways. Journal of Virology, 2010, 84(21): 11350-11358.

[33] ARUNACHALAM B, PHAN U T, GEUZE H J, CRESSWELL P. Enzymatic reduction of disulfide bonds in lysosomes: Characterization of a gamma interferon inducible lysosomal thiol reductase (GILT). Proceedings of the National Academy of Sciences of the United States of America, 2000, 97(2): 745-750.

[34] DE ZOYSA M, LEE J. Molecular cloning and expression analysis of interferon-γ inducible lysosomal thiol reductase (GILT)-like cDNA from disk abalone (Haliotis discus discus). Journal of Invertebrate Pathology, 2007, 96(3): 221-229.

[35] HUANG W S, DUAN L P, HUANG B, ZHOU L H, LIANG Y, TU CL, ZHANG F F, NIE P, WANG T. Identification of three IFN-gamma inducible lysosomal thiol reductase (GILT)-like genes in mud crab Scylla paramamosain with distinct gene organizations and patterns of expression. Gene, 2015, 570: 78-88.

[36] KONGTON K, PHONGDARA A, SRITHAWORN M T, WANNA W. Molecular cloning and expression analysis of the interferon-γinducible lysosomal thiol reductase gene from the shrimp Penaeus monodon. Molecular Biology Reports, 2011, 38: 3463-3470.

[37] KONGTON K, MCCALL K, PHONGDARA A. Identification of gamma-interferon-inducible lysosomal thiol reductase (GILT) homologues in the fruit fly Drosophila melanogaster. Developmental and Comparative Immunology, 2014, 44: 389-396.

[38] HASTINGS K T, CRESSWELL P. Disulfide reduction in the endocytic pathway: immunological functions of gamma-interferoninducible lysosomal thiol reductase. Antioxidants and Redox Signaling, 2011, 15(3): 657-668.

[39] PROKOPENKO O, MIROCHNITCHENKO O. Ischemia-reperfusioninducible protein modulates cell sensitivity to anticancer drugs by regulating activity of efflux transporter. American Journal of Physiology-Cell Physiology, 2009, 296: C1086-C1097.

[40] SHENG Y, SARIDAKIS V, SARKARI F, DUAN S, WU T, ARROWSMITH C H, FRAPPIER L. Molecular recognition of p53 and MDM2 by USP7/HAUSP. Nature Structral and Molecular Biology, 2006, 13(3): 285-291.

[41] HOLOWATY M N, FRAPPIER L. HAUSP/USP7 as an Epstein-Barr virus target. Biochemical Society Transactions, 2004, 32(5): 731-732.

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

Digital Gene Expression Analysis of Silkworm Infected by Bombyx mori Bidensovirus Zhenjiang Strain

GAO Kun, SHANG Meng-ke, QIAN He-ying, QIN Guang-xing, GUO Xi-jie

(College of Biotechnology, Jiangsu University of Science and Technology/Sericultural Research Institute, Chinese Academy of Agricultural Sciences, Zhenjiang 212018, Jiangsu)

【Objective】 The objective of this study is to screen differentially expressed genes in the Bombyx mori larvae infected with BmBDV-ZJ (B. mori bidensovirus Zhenjiang strain) and identify regulatory genes related to the virus infection and the host response so as to provide important clues for better understanding of the mechanism of B. mori resistance against BmBDV-ZJ infection. 【Method】 The differential gene expression profiles in JS B. mori larvae after oral infection with BmBDV-ZJ were constructed using Illumina Genome Analyzer platform. In order to exclude the effects of individual differences, 10 larvae were dissected and pooled as one sample for digital gene expression (DGE) analysis, respectively. The differential expression detection of genes across samples was performed using a rigorous algorithm method. False discovery rate (FDR) was used to determine the P value threshold in multiple tests and analyses. The significance of the gene expression difference was obtained through a FDR≤0.001 and the absolute value of log2ratio≥1. The gene ontology (GO) classification system was used to determine the possible functions of all differentially expressed genes. P value was calculated by GO (http://www.geneontology. org/) and corrected by Bonferoni. A corrected P value≤ 0.05 was selected as a threshold for significant enrichment of the gene sets. WEGO (web gene ontology annotation plot) software was used for visualizing, comparing and plotting GO annotation results. Pathway enrichment analysis was conducted to further identify the significantly enriched metabolic pathways or signal transduction pathways by using the KEGG database. Pathways with a Q value≤0.05 were designated as significantly enriched pathways in DGEs. Then some of the differentially expressed genes were verified by quantitative real-time PCR (qRT-PCR). 【Result】 Totally, 4 850 663 and 4 875 307 raw tags were generated in the control and BmBDV-ZJ infected DGE (digital gene expression) libraries, respectively. There were 4 757 934 and 4 788 406 clean tags corresponding to 62 436 and 63 680 distinct clean tags were filtered from the raw tags. The distribution of the total and distinct tags over the different tag abundance categories showed highly similar patterns in each DGE library. The sequencing depths reached approximately 3.5 and 3.7 million in the two DGE libraries, respectively, which satisfied the requirement for the experiment. So the two DGE libraries were reliable. The tag sequences of the two DGE libraries were mapped to the reference database of B. mori. In the control and BmBDV-ZJ-infected DGE library, 36.39% and 45.30% of the clean tags were mapped to a gene in the reference database, 50.02% and 43.34% of the clean tags could be mapped to genome of B. mori, while 13.59% and 12.35% of the clean tags were unknown tags. A total of 447 differentially expressed genes were detected, of which 306 were upregulated and 141 were downregulated. There were 218, 147, 179 differentially expressed genes have GO categories according to molecular function, cellular component and biological process, respectively. KEGG (http://www.genome.jp/kegg) ontology assignments were used to classify the functional annotations of the identified genes. Among the differentially expressed genes, 330 were mapped to 151 pathways in the KEGG database. Nineteen terms was significantly enriched (Q value≤0.05) and the cytosolic DNA-sensing pathway was significantly enriched. Moreover, 24 differentially expressed genes were verified using qRT-PCR, showing that 20 genes were concordant in the expression with DGE. Among the 9 differentially expressed genes related to cytosolic DNA-sensing pathway, BGIBMGA009408-TA, BGIBMGA004913-TA, BGIBMGA011753-TA, which were the DNA-directed RNA polymerase III genes in B. mori, were all up-regulated in the BmBDV-ZJ infected B. mori with 4.3, 2.3, 3.4-fold change, respectively. 【Conclusion】 The results of this study may serve as a basis for future research not only on the molecular mechanism of BmBDV-ZJ invasion but also on the mechanism of B. mori resistance against BmBDV-ZJ infection.

Bombyx mori; digital gene expression; Bombyx mori bidensovirus; quantitative real-time PCR

2016-04-08；接受日期：2016-06-02

國家自然科學(xué)基金（31402141）、江蘇省自然科學(xué)基金（BK20140508）

聯(lián)系方式：高坤，E-mail：gkunjn2002@126.com。通信作者郭錫杰，Tel：0511-84401328；E-mail：guoxijie@126.com