張志毅 安小平 莊 璐 馬麥卷 米志強(qiáng) 黃 勇 范 航 劉 瑋 童貽剛

(軍事醫(yī)學(xué)科學(xué)院微生物流行病研究所，病原微生物生物安全國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100071)

近年來(lái)，由于氣候變化和環(huán)境破壞等因素，蟲媒傳染病疫情不時(shí)發(fā)生，其控制變得日益復(fù)雜(鄭學(xué)禮，2011)。除了已熟知的蟲媒傳染病，20世紀(jì)70年代以來(lái)，多數(shù)年份都有一種或一種以上的新傳染病被發(fā)現(xiàn)，新發(fā)傳染病至今達(dá)40余種，而且其病原微生物種類復(fù)雜，有病毒、細(xì)菌、(包括立克次體、衣原體、螺旋體)及寄生蟲等(寶福凱等，2009)。蟲媒傳染病病原體的本底調(diào)查，對(duì)于蟲媒傳染病防控有重要意義。

傳統(tǒng)的病原體檢測(cè)方法費(fèi)時(shí)費(fèi)力，尤其對(duì)于新發(fā)傳染病，早期發(fā)現(xiàn)及診斷更為困難。測(cè)序方法是確定微生物種類的重要方法之一(秦楠等，2011)。傳統(tǒng)的Sanger序列測(cè)定方法耗時(shí)長(zhǎng)、成本高、通量低，無(wú)法滿足對(duì)于未知病原的高通量分析。由于未知病原體的核酸序列未知，也不能夠直接采用PCR 技術(shù)進(jìn)行擴(kuò)增和測(cè)序(李建彬等，2011)。近幾年來(lái)，高通量測(cè)序技術(shù)，又稱新一代測(cè)序技術(shù)(next-generation sequencing)的出現(xiàn)，使得測(cè)序費(fèi)用相對(duì)于Sanger測(cè)序方法大幅下降，已廣泛應(yīng)用于動(dòng)植物全基因組測(cè)序、基因組重測(cè)序、轉(zhuǎn)錄組測(cè)序、小RNAs測(cè)序和表觀基因組測(cè)序等領(lǐng)域(岳桂東等，2012)，對(duì)于整個(gè)生命科學(xué)研究都產(chǎn)生巨大的影響，同時(shí)也為蟲媒傳染病病原體的快速檢測(cè)提供了一條新的有效途徑。

宏基因組(Metagenomics)高通量測(cè)序方法篩查病原體是一種新型的有效方法(Nakamuraetal.，2011)，但對(duì)于病毒病原而言，由于其基因組相對(duì)于原核生物和高等生物的基因組而言非常小，相對(duì)含量很低，其序列在測(cè)序數(shù)據(jù)分析時(shí)很容易被宿主序列所掩蓋。小分子RNA(Small RNA，sRNA)是長(zhǎng)度一般為20～30個(gè)核苷酸的非編碼RNA 分子，無(wú)脊椎動(dòng)物如蠶、蚊蟲、線蟲等應(yīng)對(duì)RNA病毒入侵時(shí)可產(chǎn)生與其對(duì)應(yīng)的小干擾RNA(siRNA)，昆蟲siRNA免疫機(jī)制對(duì)入侵的RNA病毒的特定序列進(jìn)行加工和放大，并利用這些小RNA分子抑制病毒基因的表達(dá)，從而達(dá)到控制病毒復(fù)制和致病的目的(Wuetal.，2010)。這種免疫機(jī)制使昆蟲體內(nèi)小RNA中病毒序列相對(duì)的比例更加突出。研究發(fā)現(xiàn)sRNA在真核細(xì)胞和原核細(xì)胞中對(duì)基因的轉(zhuǎn)錄后調(diào)控發(fā)揮著重要作用，因此動(dòng)植物體內(nèi)都有數(shù)量巨大的小分子RNA，如何從巨量的sRNA中確定病原體RNA或病原體相關(guān)sRNA，是一個(gè)值得探討的科學(xué)問(wèn)題，然而傳統(tǒng)測(cè)序法操作復(fù)雜，花費(fèi)大，測(cè)序深度有限，效率較低，難以解決這一問(wèn)題，新一代高通量測(cè)序技術(shù)具有速度快、成本低、覆蓋度深、產(chǎn)出巨大等優(yōu)點(diǎn)，非常適合小分子RNA 測(cè)序(衛(wèi)波等，2009)。已有研究報(bào)道小分子RNA高通量測(cè)序可以用來(lái)發(fā)現(xiàn)媒介昆蟲攜帶病毒的報(bào)道(Janetal.，2009; Shietal.，2009；陳斌等，2011)。對(duì)于其他病原體，還未見(jiàn)相關(guān)報(bào)道。目前，具有傳播媒介能力的節(jié)肢動(dòng)物中以蚊類最多。已登記的535種蟲媒病毒中，從蚊類分離到的病毒占近50%，其次是從蜱分離到的病原體(116種)，占蟲媒病毒總數(shù)的21.68%(李文剛等，2011)。因此本文選取蚊和蜱為例，嘗試了基于小分子RNA高通量測(cè)序篩查媒介昆蟲攜帶病原體的方法。結(jié)果發(fā)現(xiàn)，小RNA病原體篩查范圍其實(shí)不僅限于RNA病毒和DNA病毒，還可以發(fā)現(xiàn)原核以及真核生物病原體。現(xiàn)報(bào)道如下。

1 材料與方法

1.1 實(shí)驗(yàn)材料

蚊蟲從云南西雙版納地區(qū)野外捕獲，分別為3個(gè)種類：中華按蚊Anophelessinensis、三帶喙庫(kù)蚊Culextritaeniorhynchus和致倦庫(kù)蚊Cx.quinquefasciatus，各100只左右。蜱采自黑龍江牡丹江市林區(qū)，包括全溝硬蜱Ixodespersulcatus18只、邊緣革蜱Dermacentormarginatus19只和長(zhǎng)角血蜱Haemaphysalislongicornis21只。

1.2 總RNA 的提取

液氮冷凍，同種的所有蟲體取出研磨后各以500 μL 的Trizol (購(gòu)自Invitrogen) 裂解，加入蛋白酶以65℃處理30 min后提取總RNA。具體步驟：在室溫15～30℃放置5 min，加入氯仿，用力震蕩約15 s，室溫下放置2～3 min， 12 000 r/min(4℃) 離心15 min；取上層水相加入異丙醇，在室溫下放置10 min， 12 000 r/min(2℃～8℃) 離心10 min；棄上清，加75%乙醇洗滌，渦旋混合，7 500 r/min(2～8℃) 離心5 min，棄上清; 讓沉淀的RNA在室溫下自然干燥5～10 min，用DEPC 水溶解。RNA 沉淀，55～60℃孵育10 min。

1.3 樣品的制備與Solexa 測(cè)序

siRNA 分離與高通量測(cè)序由深圳華大基因公司完成，具體流程如下：首先用聚丙烯酰胺電泳分離siRNA，切膠回收18～30 nt 大小的片段，兩端分別連接RNA 接頭，然后進(jìn)行反轉(zhuǎn)錄和PCR 擴(kuò)增，最后將擴(kuò)增產(chǎn)物上樣進(jìn)行Solexa 法測(cè)序。

1.4 數(shù)據(jù)過(guò)濾

初始數(shù)據(jù)采用軟件Solexa pipeline進(jìn)行以下幾步過(guò)濾：1)去掉低質(zhì)量的reads；2)去掉5′接頭和3′接頭序列；3)過(guò)濾無(wú)插入接頭序列；4)過(guò)濾ployA reads。

1.5 判別方法

1.5.1序列比對(duì)。使用NCBI的Blast(2.22)軟件中的blastn作核酸比對(duì)，下載NCBI核酸數(shù)據(jù)庫(kù)(nt庫(kù)2013年4月)建立本地化核酸數(shù)據(jù)庫(kù)，以過(guò)濾后的測(cè)序reads作為查詢項(xiàng)，進(jìn)行全面比對(duì)搜索。只給出得分在前10位的結(jié)果。

1.5.2blast結(jié)果過(guò)濾。通過(guò)以下幾步對(duì)結(jié)果進(jìn)行過(guò)濾：每個(gè)結(jié)果只保留得分最高的，如果有并列最高的兩項(xiàng)屬于不同物種則認(rèn)為屬于不同物種間保守序列而予以忽略；去掉核糖體RNA(因?yàn)楹颂求wRNA有較高保守性且數(shù)量較多，因此特異性差)；去掉匹配長(zhǎng)度小于21 bp的結(jié)果(因?yàn)閞eads序列較短，匹配短的很可能是隨機(jī)匹配)。

1.5.3物種分析。因?yàn)镹CBI只給出了gi號(hào)與種屬分類編號(hào)的對(duì)應(yīng)文件，從nt庫(kù)的標(biāo)題可以得到gi號(hào)與具體數(shù)據(jù)庫(kù)編號(hào)(如ref、emb、dbj等)的對(duì)應(yīng)關(guān)系。而blast結(jié)果只給出具體數(shù)據(jù)庫(kù)編號(hào)，所以由具體數(shù)據(jù)庫(kù)編號(hào)得到對(duì)應(yīng)gi號(hào)，再得到種屬分類編號(hào)，進(jìn)行統(tǒng)計(jì)得到與各屬匹配的reads數(shù)和匹配總長(zhǎng)度，最后轉(zhuǎn)化為屬的名稱。全部過(guò)程用自行編寫的python程序完成。

1.5.4校正排序。因?yàn)樾NA序列較短，比對(duì)搜索的結(jié)果中必然有較多的非特異的隨機(jī)匹配?？紤]到nt庫(kù)中搜集的某一種屬的序列越多，被隨機(jī)匹配到的可能性越大，因此我們又將最后得到每個(gè)種屬的匹配總長(zhǎng)度除以該種屬在nt庫(kù)中的總長(zhǎng)度，得到檢出比值，按照這個(gè)比值從大到小排序。這個(gè)比值越大的，表明該種屬有更大的檢出可能性。

2 結(jié)果

2.1 高通量測(cè)序統(tǒng)計(jì)結(jié)果

6只樣品分別經(jīng)過(guò)RNA提取，小RNA分離、純化和測(cè)序，并經(jīng)過(guò)數(shù)據(jù)過(guò)濾，得到有效小RNA序列的數(shù)量和測(cè)序的平均覆蓋倍數(shù)估計(jì)值在表1中列出。

表1 各蟲種測(cè)序reads數(shù)和平均覆蓋倍數(shù)Tab.1 Number of reads and coverage of the 6 vector insects

2.2 Blast及結(jié)果過(guò)濾

將各樣品的有效reads作為查詢項(xiàng)，輸入到nt庫(kù)做Blast比對(duì)搜索。只保留得分在前10的結(jié)果。然后將每一條比對(duì)結(jié)果中得分最高的一項(xiàng)匯集起來(lái)進(jìn)行統(tǒng)計(jì)。因?yàn)槊恐昀ハx樣品中，核糖體RNA的數(shù)量都占有較高比例，所以首先將結(jié)果中的核糖體RNA的條目去掉。這一步可去除11%～18%的匹配項(xiàng)。

小RNA長(zhǎng)度較短，從表1可以看出，22 nt長(zhǎng)的reads占了大多數(shù)。而nt庫(kù)含有1 600多萬(wàn)條序列，總堿基長(zhǎng)度達(dá)到37 594 899 229 nt，有許多研究詳盡的物種如人、鼠、果蠅等物種的序列，基因組很大且收錄得有冗余，因此搜索時(shí)，當(dāng)長(zhǎng)度小于等于20時(shí)，得到隨機(jī)匹配的幾率非常大。我們?cè)鲞^(guò)測(cè)試，設(shè)計(jì)了10 000條22 bp的隨機(jī)序列，去和nt庫(kù)比對(duì)，發(fā)現(xiàn)匹配在18～20 bp的條目占了接近60%的比例。因此，選取Blast結(jié)果中匹配堿基數(shù)大于20的條目能夠進(jìn)一步去除非特異隨機(jī)匹配，使其不至于干擾真實(shí)結(jié)果。結(jié)果表明，可進(jìn)一步去掉約2/3的條目。

2.3 病原體篩查結(jié)果

對(duì)剩下的Blast結(jié)果進(jìn)行統(tǒng)計(jì)，首先根據(jù)每一條的數(shù)據(jù)庫(kù)編號(hào)找到對(duì)應(yīng)的GenBank序號(hào)(gi號(hào))，并記下其匹配堿基數(shù)，這一步完成后，將gi號(hào)相同的條目的匹配堿基數(shù)相加，得到每個(gè)gi號(hào)對(duì)應(yīng)匹配總長(zhǎng)度之和。而后將gi號(hào)變?yōu)閷?duì)應(yīng)的物種分類編號(hào)，就得到所包含的每個(gè)物種的匹配總長(zhǎng)度。結(jié)果中可以發(fā)現(xiàn)，由于人、小鼠等物種在nt庫(kù)中總長(zhǎng)度很大，所以其隨機(jī)匹配幾率很高，在結(jié)果中排名很靠前。為了排除這類非隨機(jī)匹配的干擾，用每個(gè)種屬的匹配長(zhǎng)度除以該種屬在nt庫(kù)中的總長(zhǎng)度的比值大小來(lái)排序，就可以使檢出比超出背景值的物種排在前面。

在結(jié)果中，排在最前面的是所測(cè)昆蟲的物種及其近緣物種，其數(shù)量占有絕對(duì)優(yōu)勢(shì)。此外一些環(huán)境微生物、植物、真菌也占有相當(dāng)比例。除去此類一般意義的物種外，一些重要的病原體是我們最關(guān)注的。

根據(jù)以上算法，我們得到的初步結(jié)果中(表2)，中華按蚊可能攜帶有乙型腦炎病毒。乙型腦炎又稱“日本腦炎”，主要通過(guò)蚊蟲叮咬傳播，庫(kù)蚊、伊蚊、按蚊中均有一些種可以傳播此病，主要流行區(qū)域?yàn)闁|南亞、西太平洋地區(qū)，我國(guó)除新疆、西藏、青海外均為乙腦疫區(qū)(郭楊等，2008)。三帶喙庫(kù)蚊則可能攜帶的病原體是庫(kù)蚊濃核病毒(Brevidensovirus)。蚊濃核病毒屬于細(xì)小病毒科，特異性感染蚊類，引發(fā)特異性的細(xì)胞核致密增生病變并可最終導(dǎo)致宿主發(fā)病或死亡(顧金保等，2008)。根據(jù)顧金保的統(tǒng)計(jì)，2007年以前發(fā)現(xiàn)的蚊濃核病毒，有兩株是在中國(guó)發(fā)現(xiàn)的，分別是淡色庫(kù)蚊濃核病毒(Culexpipienspallensdensovirus, CpDNV)和C6/36濃核病毒(C6/36 DNV)。我們樣品中的病毒應(yīng)該更接近于前者。庫(kù)蚊黃病毒為只感染蚊蟲的黃病毒(Kimetal.，2009)，有研究報(bào)道該病毒可與黃病毒屬的其他病毒共同感染三帶喙庫(kù)蚊(Kentetal.，2010)，但是否相互影響還有待研究。致倦庫(kù)蚊中發(fā)現(xiàn)的可能病原體則是果蠅X病毒，該病毒由兩條核酸鏈組成，A鏈序列于1996年測(cè)定(Chungetal.，1996)，B鏈于2002年測(cè)定(Shwedetal.，2002)，但該病毒在研究昆蟲的病毒免疫機(jī)制中應(yīng)用較多，對(duì)人的致病性還待深入研究。全溝硬蜱最可能攜帶的病原體有：peacockii立克次體(Rickettsiapeacockii，Rp)。邊緣革蜱最可能攜帶的病原體有：立氏立克次體(Rickettsiarickettsii，Rr)。長(zhǎng)角血蜱最可能攜帶的病原體是：伯氏柯克斯體Coxiellaburnetii。R.peacockii和R.rickettsii都屬于斑點(diǎn)熱群立克次體，目前NCBI物種分類數(shù)據(jù)庫(kù)收錄的斑點(diǎn)熱群立克次體不少于幾十種。2005的一項(xiàng)調(diào)查研究，總結(jié)了我國(guó)之前引起立克次體病的10種病原體(張麗娟等，2005)，其中只有伯氏柯克斯體出現(xiàn)血蜱中，R.peacockii和R.rickettsii在國(guó)內(nèi)較少有報(bào)道。伯氏柯克斯體習(xí)慣上稱為Q熱立克次體，它所引起的Q熱是一種人獸共患病。

表2 蜱與蚊蟲中可能含有的病原體Tab.2 Possible pathogens contained in ticks and mosquitos

3 討論

常規(guī)病原體高通量檢測(cè)需要對(duì)不同的病原體采用不同的純化方案和測(cè)序方案，如細(xì)菌和原蟲作為細(xì)胞生物，需要按照細(xì)胞生物的分離方法將它們與人體血液和組織成分分離，并提取DNA進(jìn)行測(cè)序；而DNA病毒及支原體、衣原體需要按照可濾過(guò)性病原體進(jìn)行分離和DNA樣品提取；RNA病毒則需要采用濾過(guò)性提取方法提取和進(jìn)行RNA測(cè)序。對(duì)于未知病原體，上述3種方案需要分別考慮，選擇使用，最可靠的辦法是同時(shí)采用上述3種方案進(jìn)行高通量測(cè)序分析，這些操作在時(shí)間、人力和經(jīng)濟(jì)上都會(huì)造成很大負(fù)擔(dān)。如果采用RNA高通量測(cè)序技術(shù)，則可以用一種測(cè)序方案檢測(cè)所有類型的病原體，這不僅能夠節(jié)省很多的人工和試劑費(fèi)用，更重要的是可以節(jié)省大量時(shí)間，尤其是在有暴發(fā)疫情的時(shí)候，病原體篩查的速度將顯得十分關(guān)鍵。小RNA高通量測(cè)序篩選媒介昆蟲攜帶病毒已經(jīng)被證明是一種有效的策略(Jan，2009; Wuetal.，2010)。本研究結(jié)果表明，用這一方法也不僅可以檢測(cè)媒介昆蟲攜帶的病毒，而且對(duì)細(xì)菌及原蟲等病原體的檢測(cè)也具有很大的潛力。因此用小RNA高通量測(cè)序篩查病原體，只需一次測(cè)序，一次分析就得到樣品可能含有的病毒、原核及真核的各類病原體，具有快速、簡(jiǎn)便、節(jié)省費(fèi)用的突出優(yōu)點(diǎn)。

小RNA由于其在調(diào)控、免疫等方面的重要作用，當(dāng)前受到格外關(guān)注，各類研究已積累了大量的小RNA高通量測(cè)序數(shù)據(jù)。因此用本文所述的方法對(duì)這些小RNA數(shù)據(jù)進(jìn)行分析，可獲得樣品中可能含有的各類物種的信息，這對(duì)于充分利用實(shí)驗(yàn)數(shù)據(jù)獲取有用信息具有一定意義，也可能對(duì)疾病的發(fā)生機(jī)制、發(fā)展進(jìn)程與預(yù)后提供一些新的線索。

通過(guò)上述的生物信息學(xué)分析，大大縮小了后期實(shí)驗(yàn)驗(yàn)證的范圍。但由于小分子RNA核苷酸序列長(zhǎng)度較短，因此在對(duì)大的核酸數(shù)據(jù)庫(kù)比對(duì)搜索時(shí)，非特異性的隨機(jī)匹配可能性較大，因此有必要進(jìn)一步優(yōu)化生物信息學(xué)分析參數(shù)來(lái)盡可能排除掉比對(duì)結(jié)果中的各種非特異匹配，并需進(jìn)一步做實(shí)驗(yàn)驗(yàn)證。