基于高通量測序的真核微藻多樣性質(zhì)控分析

宋 倫，吳 景，李 楠，孫 明，杜 靜，王 鵬

( 1.哈爾濱工業(yè)大學 環(huán)境學院，城市水資源與水環(huán)境國家重點實驗室，黑龍江 哈爾濱 150090； 2.遼寧省海洋水產(chǎn)科學研究院，遼寧省海洋生物資源與生態(tài)學重點實驗室，遼寧 大連 116023 )

海洋微藻是海洋生態(tài)系統(tǒng)中物質(zhì)循環(huán)和能量流動的基礎，是生命和非生命系統(tǒng)聯(lián)系的關鍵環(huán)節(jié)，同時也是赤潮、褐潮暴發(fā)的致災種[1-2]。由于微藻采樣和顯微鏡精度限制，傳統(tǒng)的形態(tài)學鑒定對微藻多樣性的挖掘存在較大局限性。褐潮的暴發(fā)引起了諸多學者對微微型藻類多樣性的關注，早期的分子鑒定技術雖可增加目標種類的鑒定準確度，但對物種多樣性的種類挖掘尚顯遜色。高通量測序技術具有簡單快速、測序通量高、錯誤率低和成本低等特點，為微藻的多樣性檢測提供了新手段。18S rDNA是真核生物中具有信息和功能兩種作用的核糖體編碼基因[3]，由保守區(qū)和9個可變區(qū)交替排列組成[4]，可作為物種鑒定的分子標記。國內(nèi)外以18S rDNA作為分子標記使用高通量測序技術分析藻類多樣性的研究已有報道。

2000年，國內(nèi)學者針對微型和微微型浮游藻類建立了18S rDNA克隆文庫，分析了浮游藻類的多樣性[5-6]。2011年，于杰等[7]建立了微型浮游生物的18S rDNA可變區(qū)V9的克隆文庫并進行測序，發(fā)現(xiàn)褐潮是由抑食金球藻(Aureococcusanophagefferens)和多形微眼藻(Minutocelluspolymorphus)共同形成。2015年，也有學者對渤海海域[6]及我國東海靠近濟州島海域[8]真核浮游植物群落組成和多樣性進行了研究。Stoeck等[9]采用高通量測序方法分析了海水中真核生物群落，并比較了18S rDNA的V9區(qū)和V4區(qū)對群落多樣性研究的差異，結果顯示，V4區(qū)在生物多樣性研究的準確性中有一定優(yōu)勢。

以往454測序平臺以讀長優(yōu)勢廣泛應用于18S rDNA測序領域，隨著Hiseq平臺雙端PE250測序策略的發(fā)展，讀長不斷增長，使得基于HiSeq測序研究物種多樣性的準確性大幅度提高[5]。利用HiSeq平臺對18S rDNA的一個或多個高變區(qū)測序，具有測序深度高、利于鑒定低豐度群落物種以及費用低的特點，已成為研究微生物群落多樣性的首選之策[10-11]。目前有關微藻的高通量測序鑒定技術已有較多應用，但在鑒定序列聚類水平等方面尚存在爭議[12]。微生物多樣性研究中多采用97%以上的聚類水平[11,13]，真核微藻多樣性的研究中97%的聚類水平應用較多[8-9]。同時研究者對不同分類水平下微藻多樣性的差異進行了比較，并指出高的聚類水平可能會引入更多的誤差，所以應作適當篩選[9]。筆者結合遼東灣真核微藻的高通量測序結果，比較不同序列聚類水平對真核微藻多樣性研究的影響，并展示97%聚類水平下遼東灣海域浮游植物多樣性，以期進一步優(yōu)化真核微藻的分子鑒定技術。

1 材料與方法

1.1 樣品采集及制備

2014年5月在遼東灣西、中、東側海域分別網(wǎng)格化設置：靠近西岸LS01～05組(LS組)、遠岸OS01～05組(OS組)、靠近東岸RS01～07組(RS組)，3組17個調(diào)查站位(圖1)。每站采集表層海水1 L經(jīng)0.22 μm微孔濾膜收集浮游生物樣本，最后將濾膜轉(zhuǎn)移至1.5 mL無菌離心管中，置于-20 ℃或-80 ℃冷凍保存、運輸。同時采集0.5 L水樣5%甲醛固定用于微藻種類及數(shù)量鏡檢。

圖1 遼東灣采樣站位示意

1.2 分析方法

1.2.1 基因組DNA的提取

采用CTAB法提取真核微藻宏基因組，DNA含量及純度檢測合格后于-20 ℃冰箱保存?zhèn)溆?，具體檢測步驟參考文獻[14]。

1.2.2 18S rDNA可變區(qū)V4的PCR擴增

利用自行開發(fā)的真核微藻18S rDNA的V4區(qū)基因擴增引物進行PCR反應，產(chǎn)物檢測合格后構建文庫，再經(jīng)過Qubit定量和文庫檢測，合格后使用Hiseq2500 PE250進行上機測序，具體步驟參考文獻[14]。

1.2.3 數(shù)據(jù)分析

測序得到的原始數(shù)據(jù)使用FLASH軟件進行拼接，參照Qiime軟件質(zhì)量控制流程將拼接后的序列經(jīng)過截取、過濾得到有效數(shù)據(jù)。利用Uparse對有效數(shù)據(jù)進行可操作分類單元(OTUs)聚類和物種分類，采用RDP Classifier方法與SILVA數(shù)據(jù)庫對可操作分類單元代表序列進行物種注釋[15]。

公式計算、數(shù)據(jù)分析、方差檢驗、圖件繪制均通過Microsoft Excel 2007和SPSS 19.0、PRIMER 5.0軟件完成。

2 結 果

2.1 測序數(shù)據(jù)質(zhì)量分析

通過Illumina HiSeq 2500測序平臺進行PE250模式測序得到的原始數(shù)據(jù)，存在一定比例的干擾數(shù)據(jù)，為使信息分析的結果更加準確、可靠，首先使用FLASH軟件對原始數(shù)據(jù)進行拼接，參照Qiime軟件質(zhì)量控制流程，將拼接后的序列經(jīng)過截取、過濾得到有效數(shù)據(jù)。數(shù)據(jù)處理過程中各步驟得到的統(tǒng)計結果見表1。每個樣品獲得的原始序列經(jīng)過拼接和質(zhì)量過濾，高質(zhì)量序列數(shù)均達到87%以上，最終可用于后續(xù)分析的序列數(shù)中堿基質(zhì)量值大于20(測序錯誤率小于1%)和30(測序錯誤率小于0.1%)的堿基占比分別達到了98%和97%以上，表明測得的數(shù)據(jù)準確可靠。

2.2 序列聚類水平優(yōu)化

隨后基于有效數(shù)據(jù)進行可操作分類單元聚類和物種分類分析。為研究各樣本的物種組成，對所有樣本的有效序列進行可操作分類單元聚類注釋，分別以97%、98%、99%的一致性進行可操作分類單元聚類，研究不同聚類水平對真核微藻多樣性統(tǒng)計差異的影響。結果表明，隨著聚類水平的上升，獲得注釋的物種數(shù)量也隨之增加，97%聚類水平共獲得注釋的可操作分類單元777個，98%聚類水平共獲得注釋的可操作分類單元996個，99%聚類水平共獲得注釋的可操作分類單元1705個，99%聚類水平獲得注釋的可操作分類單元顯著高于97%和98%(P>0.01)，而97%和98%聚類水平獲得注釋的可操作分類單元差異不顯著(P<0.05)(圖2)。

為評估各聚類水平對目標物種豐度的鑒定準確度，將鏡檢豐度較高、各站位均出現(xiàn)的夜光藻作為評估對象，夜光藻在不同聚類水平下各站位豐度占比及與鏡檢值擬合度分析見圖3。結果表明，以各站位鏡檢夜光藻豐度占比為基準值，97%聚類水平獲得的夜光藻注釋豐度占比擬合度最高(r2=0.9239)，98%、99%聚類水平反而依次下降。由于PCR、測序中不可避免的錯誤，或無法追蹤基因組內(nèi)多態(tài)性，減少高聚類水平注釋物種多樣性的誤差，筆者應用了寬泛的聚類水平進行物種注釋分析(97%水平)。

表1 測序質(zhì)控數(shù)據(jù)

注：*Q20和Q30：序列中堿基質(zhì)量值大于20(測序錯誤率小于1%)和30(測序錯誤率小于0.1%)的堿基所占的百分比；**%：過濾低質(zhì)量和短長度后的序列數(shù)的數(shù)目與原始下機的序列數(shù)數(shù)目的百分比.

圖2 不同聚類水平下各站位注釋的可操作分類單元數(shù)量變化

圖3 夜光藻在不同聚類水平和鏡檢下豐度占比

2.3 測序數(shù)據(jù)量合理性

為評估物種測序數(shù)據(jù)量的合理性，從樣本中隨機抽取一定測序量的數(shù)據(jù)，統(tǒng)計它們所代表物種的數(shù)目(即可操作分類單元數(shù)目)，以抽取的測序數(shù)據(jù)量與對應的物種數(shù)來構建稀釋曲線。同時為考察所檢測物種豐度數(shù)據(jù)量的合理性，將樣本中的可操作分類單元按相對豐度由大到小排序得到對應的排序編號，再以可操作分類單元的排序編號為橫坐標，可操作分類單元中的相對豐度為縱坐標，將這些點用折線連接即繪制得到等級聚類曲線。在水平方向上，物種的豐富度由曲線的寬度來反映，物種的豐富度越高，曲線在橫軸上的跨度越大；在垂直方向上，曲線的平滑程度，反映了樣本中物種的均勻程度，曲線越平緩，物種分布越均勻[16]。結果表明，可操作分類單元相對豐度在10-4以下時，曲線趨于平緩，說明測序數(shù)據(jù)量已趨于合理，再多的數(shù)據(jù)量對群落整體豐度影響不大。按LS、OS、RS 3組進行稀釋曲線和等級聚類曲線分析見圖4。結果顯示，3組隨機抽取的測序條數(shù)達到15 000條時，曲線趨向平坦，說明測序數(shù)據(jù)量漸進合理，其中遼東灣西、東兩側的LS、RS站位獲得的可操作分類單元數(shù)量趨勢一致，均較中間海域的OS站位獲得的可操作分類單元數(shù)量高。同樣，3組站位可操作分類單元相對豐度在10-4以下時，曲線趨于平緩，說明測序數(shù)據(jù)量已趨于合理，其中RS站位曲線相對平緩，物種分布比較均勻(圖5)。

利用稀釋曲線和等級聚類曲線評估測序數(shù)據(jù)量的合理性。稀釋曲線主要關注物種α多樣性的挖掘程度，當曲線趨向平坦，說明測序數(shù)據(jù)量漸進合理，增加測序數(shù)據(jù)量浪費時間和成本。等級聚類曲線主要關注物種β多樣性，在垂直方向上，曲線趨于平緩，說明測序數(shù)據(jù)量已趨于合理，無需增加測序時間和成本。本研究對遼東灣的測序量無論站位評估還是組內(nèi)評估均達到合理水平。

2.4 遼東浮游植物名錄、分布及特征

根據(jù)高通量測序結果，筆者整理了種水平遼東灣浮游植物名錄、分布及特性(表 2)。名錄包含90種藻類，通過比較發(fā)現(xiàn)，遼東灣西岸(LS組)與東岸的浮游植物種類數(shù)更相近，比中部海域(OS組)站位種類豐富，包括古多甲藻(Archaeperidiniumsaanichi)、凱聰杜波斯克藻(Euduboscquellacachoni)、伊姆裸甲藻(Gymnodiniumimpudicum)、五刺多甲藻(Peridiniumquinquecorne)、渦輪共甲藻(Syndiniumturbo)、球半隱藻屬藻類(Hemiselmisbrunnescens)、恩諾巴夫藻(Pavlovaennorea)、硅鞭藻綱藻類(Ciliophrysinfusionum)、中華擬根管藻(Rhizosoleniasimiloides)、近囊胞藻(Paraphysomonasvestita)。遼東灣西岸與東岸的藻類分布也有差異，其中包括古多甲藻、多紋膝溝藻(Gonyaulaxpolygramma)、平行原多甲藻(Protoperidiniumpellucidum)、南極衣藻(Chlamydomonasraudensis)等10種藻類。本次調(diào)查中發(fā)現(xiàn)7種導致魚類死亡的藻類，25種引發(fā)過赤潮或者褐潮的藻類，21種含有毒素的藻類(表2)。

圖4 稀釋曲線

圖5 等級聚類曲線

分類地位學名站位LSOSRS毒害性所含毒素甲藻門甲藻綱紅哈卡藻(Akashiwosanguinea)+++B, FHEM相關亞歷山大藻(Alexandriumaffine)+++BPSP偽亞歷山大藻(A.fraterculus)+++BPSP李氏亞歷山大藻(A.leei)+++PSP假性亞歷山大藻(A. pseudogoniaulax)+++BPSP古多甲藻+凱聰杜波斯克藻++刺毛杜波斯克藻(E.crenulata)+++多紋膝溝藻+B, FHEM具刺膝溝藻(G.spinifera)+++BYTX底刺膝溝藻(G.verior)+++伊姆裸甲藻++BPSP無紋裸甲藻(G.instriatum)+++PSP微型裸甲藻(G.microreticulatum)+++PSP微小卡羅藻(Karlodiniummicrum)+++B,FHEM帕氏帆甲藻(Kofoidinium cf. pavillardii)+++多邊舌甲藻(Lingulodiniumpolyedrum)+++YTX夜光藻(Noctilucascintillans)+++B, F五刺多甲藻++B, F圓禿頂藻(Phalacromarotundatum)+++科夫多溝藻(Polykrikoskofoidii)+++斯氏多溝藻(P.schwartzii)+++B海洋原甲藻(Prorocentrummicans)+++微小原甲藻(P.minimum)+++BDSP網(wǎng)狀原角藻(P.reticulatum)+++YTX雙角多甲藻(Protoperidiniumbipes)+++單卷原多甲藻(P.monovelum)+++平行原多甲藻++B三角原多甲藻(P.tricingulatum)+++斯氏扁甲藻(Pyrophacussteinii)+++B渦輪共甲藻++綠藻門綠藻綱南極衣藻++單針藻(Monoraphidiumcontortum)+++綠藻類麥可屬(Mychonasteshomosphaera)+++綠枝藻綱塔形藻(Cymbomonastetramitiformis)+++微微型藻類(Dolichomastixtenuilepis)+++青綠藻(Mamiellagilva)+++細小微胞藻(Micromonaspusilla)+++B雙鞭綠藻(Nephroselmispyriformis)+++

續(xù)表2

分類地位學名站位LSOSRS毒害性所含毒素綠枝藻綱單細胞球狀綠藻(Ostreococcustauri)+++B具翅冠突藻(Pterospermacristatum)+++密球藻(Pycnococcusprovasolii)+++四嘴塔胞藻(Pyramimonastetrarhynchus)++共球藻綱尖粒藻(Amphikrikos cf. nanus)+++小球藻(Chlorellavulgaris)+++喜糖橢圓球藻(Chloroidiumsaccharophilum)+++石莼綱藻類(Dilabifilumarthropyreniae)+++球半隱藻屬藻類(H.brunnescens)+球半隱藻屬藻類(H.cryptochromatica)+++海洋白隱藻(Leucocryptosmarina)+++雙尖全溝藻(Teleaulaxamphioxeia)+++尖尾全溝藻(Teleaulaxgracilis)+++隱藻門隱藻綱恩諾巴夫藻++普林藻綱草莓定鞭金藻(Haptolinafragaria)+++球等鞭金藻(Isochrysisgalbana)++魚腥棕囊藻(Phaeocystiscordata)+++FHEM球形棕囊藻(Phaeocystisglobosa)+++B, FHEM浮生藻綱抑食金球藻+++BEPS硅鞭藻綱弧形無柄鐘藻(Apedinellaradians)+++硅鞭藻綱藻類(Ciliophrysinfusionum)++小鞭毛蟲(Florenciellaparvula)+++硅鞭藻綱藻類(Pseudochattonellafarcimen)+++彈性假柄鐘藻(Pseudopedinellaelastica)+++硅藻門中心硅藻綱冰河擬星桿藻(Asterionellopsisglacialis)+++B長角四棘藻(Attheyalongicornis)++大洋角管藻(Cerataulinapelagica)+++鈣質(zhì)角毛藻(Chaetoceroscalcitrans)+++并基角毛藻(C.decipiens)+++牟氏角毛藻(C.muellerii)+++布氏雙尾藻(Ditylumbrightwellii)+++B柔弱幾內(nèi)亞藻(Guinardiadelicatula)+++B泰晤士旋鞘藻(Helicothecatamesis)+++中華擬根管藻++曼氏骨條藻(Skeletonemamenzellii)+++塔形冠蓋藻(Stephanopyxisturris)+++海鏈藻屬藻類(Thalassiosiraconcaviuscula)+++大洋海鏈藻(T.oceanica)+++羽紋藻綱新月柱鞘藻(Cylindrothecaclosterium)+++B

續(xù)表2

分類地位學名站位LSOSRS毒害性所含毒素羽紋藻綱翼繭形藻(Amphiprora alata)+++奇異楔形藻(Licmophoraparadoxa)+++中型斜紋藻(Pleurosigmaintermedium)+++斜紋藻屬藻類(P.planktonicum)+++中華齒狀藻(Odontellasinensis)+++金藻門金藻綱黃金色鞭毛藻(Poterioochromonasmalhamensis)+++HEM近囊胞藻++針胞藻綱海洋卡盾藻(Chattonellamarina)++BHEM日本纖囊藻(Fibrocapsajaponica)+++BNSP真眼點藻綱微擬球藻(Nannochloropsisgaditana)++顆粒微擬球藻(N.granulata)+++黃藻門黃藻綱赤潮異彎藻(Heterosigmaakashiwo)+++BHEM

注：+.此站位出現(xiàn)； B.藻華暴發(fā)； F.致魚死； HEM.溶血素； PSP.麻痹性貝類中毒； DSP.腹瀉性貝毒； NSP.神經(jīng)性貝毒； YTX.蝦夷扇貝毒素.

3 討 論

3.1 測序數(shù)據(jù)質(zhì)量對真核微藻多樣性質(zhì)控的影響

自DNA樣本到最終數(shù)據(jù)獲得的過程中，樣本檢測、PCR、純化、建庫、測序每一個環(huán)節(jié)都會對數(shù)據(jù)質(zhì)量和數(shù)量產(chǎn)生影響，而數(shù)據(jù)質(zhì)量又會直接影響后續(xù)信息分析的結果。為了從源頭上保證測序數(shù)據(jù)的準確性、可靠性，應對樣本檢測、建庫、測序每一個試驗步驟嚴格把控，從根本上確保高質(zhì)量數(shù)據(jù)的產(chǎn)出。本研究在遼東灣17個站位每個樣品獲得的高質(zhì)量序列數(shù)均達到87%以上，用于注釋的序列數(shù)達97%以上，表明測得的數(shù)據(jù)準確可靠。

3.2 序列聚類水平對真核微藻多樣性質(zhì)控的影響

筆者研究了97%、98%、99%不同聚類水平對真核微藻多樣性統(tǒng)計差異的影響發(fā)現(xiàn)，雖然隨著聚類水平的上升，獲得注釋的物種數(shù)量也隨之增加，但增加的部分豐度很低，均無法獲得注釋(未獲培養(yǎng)的種類)。通過顯微鏡檢測方法檢測夜光藻的豐度，并以此為基準值，與高通量測序方法獲得的夜光藻豐度進行比較，發(fā)現(xiàn)97%聚類水平獲得的夜光藻注釋豐度占比擬合度最高。由于測序過程中存在一些不可避免的錯誤，如核苷酸誤讀，可能由于Taq聚合酶錯誤或測序技術錯誤導致[17]。在380 bp序列中無法追蹤到的基因組內(nèi)多態(tài)性和PCR重組的發(fā)生率也不可避免，這可能會導致對新物種和多樣性的高估[15]。另外，目前許多微藻分子測序鑒定文獻均使用97%一致性聚類水平，因此，為避免物種多樣性的高估，增加與相關研究文獻的可比性[18-19]，建議使用相對寬泛的聚類水平(97%)進行序列聚類分析。97%的聚類水平下，利用稀釋曲線和等級聚類曲線評估測序數(shù)據(jù)量的合理性。稀釋曲線主要關注物種α多樣性的挖掘程度，當曲線趨向平坦，說明測序數(shù)據(jù)量漸進合理，增加測序數(shù)據(jù)量浪費時間和成本。等級聚類曲線主要關注物種β多樣性，在垂直方向上，曲線趨于平緩，說明測序數(shù)據(jù)量已趨于合理，無需增加測序時間和成本。本研究對遼東灣的測序量無論站位評估還是組內(nèi)評估均達到合理水平。

3.3 遼東灣浮游植物名錄、分布及特征

根據(jù)3種聚類水平的比較結果，整理了97%水平下遼東灣浮游植物種類，名錄包含90種藻類。通過比較發(fā)現(xiàn)，遼東灣西岸(LS組)與東岸(RS組)的浮游植物種類更相近，比遠離岸邊的站位(OS組)種類要多，遼東灣西岸與東岸的藻類分布也有差異。浮游藻類分布的差異，可能是藻類本身特性影響，對營養(yǎng)物質(zhì)不同造成的，也可能是調(diào)查中的差異操作造成的，仍需要多次調(diào)查進行證實。該次調(diào)查中，發(fā)現(xiàn)7種導致魚類死亡的藻類，25種引發(fā)過赤潮或褐潮的藻類，21種含有毒素的藻類(表2)。本次在遼東灣海域發(fā)現(xiàn)的紅哈卡藻，是一種赤潮藻，曾在佛羅里達和南加州海岸附近引發(fā)赤潮[20]。抑食金球藻可產(chǎn)生一種抑制雙殼類動物攝食的胞外多糖，使硬殼蛤(Mercenariamercenaria)、海灣扇貝(Argopectenirradians)、紫貽貝(Mytilusgalloprovincialis)等纖毛活動受到抑制，進而停止攝食，最終大量死亡[21]。球形棕囊藻在福建、廣東沿海大量暴發(fā)，導致約6×104t魚的死亡[22]。本次調(diào)查發(fā)現(xiàn)，多種含有毒素的藻類，例如李氏亞歷山大藻、無紋裸甲藻、微型裸甲藻等含有麻痹性貝毒，具刺膝溝藻、微小卡羅藻、魚腥棕囊藻等含有溶血素，日本纖囊藻含有神經(jīng)性貝毒。這些含有毒素的藻類在沒有大規(guī)模暴發(fā)時，暫時不會對水環(huán)境和生物造成危害，但仍存在潛在的風險，需要連續(xù)的調(diào)查和監(jiān)測。

高通量測序檢測的藻類DNA信息大部分未能注釋到種，主要原因是由于很多種類尚未錄入數(shù)據(jù)庫，導致注釋到種的成功率較低。雖然擴增區(qū)域、引物、PCR、擴增效率、讀長解析和數(shù)據(jù)庫等因素會影響微藻分子鑒定準確度[23]，但相比傳統(tǒng)形態(tài)學鏡檢方式，分子鑒定對物種多樣性的發(fā)掘已有很大提升，尤其對微型、微微型藻類的發(fā)現(xiàn)提供了科學手段，同時大大降低了形態(tài)學鑒定的人為主觀誤差。隨著高通量測序技術的發(fā)展，測序通量不斷增大，測序片段長度也可逐漸增長，數(shù)據(jù)庫會逐漸完善，物種注釋的準確度也會逐漸提升。