基于DNA條形碼基因和核基因的美麗小條鰍隱存多樣性研究

高 尚，向登高，李躍飛，李 捷，陳蔚濤

(1.中國水產(chǎn)科學研究院珠江水產(chǎn)研究所，廣州 510380；2.上海海洋大學海洋生態(tài)與環(huán)境學院，上海 201306；3.貴州大學動物科學學院，貴陽 550025)

美麗小條鰍(Micronemacheiluspulcher)隸屬于鯉形目(Cypriniformes)條鰍科(Nemacheilidae)小條鰍屬(Micronemacheilus)，是廣泛分布在我國華南地區(qū)的小型底棲魚類[1]。由于受到水利開發(fā)、過度捕撈、環(huán)境污染等諸多因素的影響,美麗小條鰍的資源量衰退嚴重,亟需給予關注和保護。運用遺傳進化分析手段了解美麗小條鰍的遺傳譜系及其空間分布格局對該物種遺傳種質(zhì)資源的保護和管理具有重要實際價值。已有研究表明，小型底棲魚類擴散能力弱，容易在小范圍內(nèi)形成遺傳分化，甚至形成高度分化的遺傳譜系或者新種[2,3]。我國華南地區(qū)河流網(wǎng)絡復雜，擁有珠江水系(包括東江、西江、北江三大干流)、海南島諸河流及其周邊小型陸封型河流，如南流江、漠陽江等。這些河流雖然相隔的地理位置不大，但卻相互被海水或者山脈隔離，極大地限制了魚類群體之間的遷移擴散。另外，華南地區(qū)經(jīng)歷的系列地質(zhì)事件，如海平面波動、山脈隆升等[4,5]也在一定程度上深刻影響著魚類的進化歷史[6]。鑒于上述背景資料，推測美麗小條鰍可能形成了相互獨立的遺傳譜系甚至新種。丘城鋒等[7]和慶寧等[8]分別利用線粒體Cytb基因和控制區(qū)片段發(fā)現(xiàn)華南西部及海南島的美麗小條鰍已經(jīng)形成了多個遺傳譜系。然而，由于早期研究涉及研究區(qū)域比較有限，并且僅使用了線粒體分子標記，因此很有必要利用高覆蓋度的樣本和不同類型基因來重新評估美麗小條鰍在遺傳層面上的多樣性水平以及是否已經(jīng)進化成多個物種。

DNA條形碼作為分類學中輔助物種鑒定的技術，已被廣泛應用于物種鑒定[9-11]或隱存種發(fā)現(xiàn)的相關研究中[12-14]。針對脊椎動物類群，DNA條形碼技術主要通過利用線粒體細胞色素c氧化酶Ⅰ中5′端約650 bp的堿基序列來實現(xiàn)物種的鑒定或者隱存種的挖掘。另外，隨著分析方法的改進和完善，越來越多的研究傾向于聯(lián)合線粒體基因和核基因等多種分子標記進行隱存種的驗證[2,15,16]。因此，本研究通過對華南地區(qū)多個水系19個站位的美麗小條鰍樣本的采集，利用DNA條形碼基因(COI基因)和糖基轉移酶基因(Glyt基因)評估美麗小條鰍的進化譜系組成，并驗證美麗小條鰍是否存在隱存種，旨在為美麗小條鰍后期的管理和保護提供科學支撐。

1 材料與方法

1.1 樣品采集

分別于2017年7月-2020年9月進行采集，在珠江水系(東江、西江、北江)、漠陽江、南流江、海南島的凌水河和南渡江的19個站位(圖1)共采集美麗小條鰍樣品107尾(表1)。剪取每尾樣本的少量胸鰭鰭條并保存在95%的酒精中，帶回實驗室用以DNA提取。

圖1 采樣示意圖

表1 樣本采集信息

1.2 基因組DNA提取、擴增與測序

采用高鹽法提取基因組總DNA[9]。利用通用引物FishF1(5′-TCAACCAACCACAAAGACATTGGCAC-3′)和FishR1(5′-TAGACTTCTGGGTGGCCAAAGAATCA-3′)[10]擴增線粒體COI基因5′端約650 bp的序列。PCR反應擴增體系為25 μL:2×PCR mix(0.1 U/μL)12.5 μL，正反向引物(10 mmol/L)1 μL,DNA模板1 μL，最后用ddH2O補齊至25 μL。PCR反應條件為：95 ℃預變性5 min，94 ℃變性30 s，54 ℃退火30 s，72 ℃延伸1 min，30個循環(huán)，最后再72 ℃延伸10 min。篩選部分樣本進行糖基轉移酶基因(Glyt)的擴增和測序。使用引物Glyt_F559(5′-GGACTGTCMAAGATGACCACMT-3′)和Glyt_R1562(5′-CCCAAGAGGTTCTTGTTRAAGAT-3′)[11]對Glyt基因進行擴增。PCR反應擴增提取與COI基因相同，反應條件為：95 ℃預變性5 min，95 ℃變性30 s，55 ℃退火1 min，72 ℃延伸1.5 min，35個循環(huán)，最后再72 ℃延伸10 min。擴增產(chǎn)物經(jīng)1%瓊脂糖凝膠電泳檢測后送往測序公司進行雙向測序。

1.3 數(shù)據(jù)處理與分析

測序獲得的DNA序列采用Lasergene v7.1軟件包中(DNASTAR,Inc.,USA)的Seqman軟件進行人工核查、校正及組裝拼接。拼接完成的序列用MEGA 6.06軟件[17]進行多重比對，剪切掉冗余片段獲得一致序列備用。

利用DnaSP 5.10[18]統(tǒng)計單倍型數(shù)，使用MEGA 6.06軟件分析序列特征。選取橫紋南鰍為外類群，采用鄰接法(Neighbor-Joining,NJ)并基于Kimura雙參數(shù)(K2P)模型[19]構建所有美麗小條鰍COI基因序列鄰接關系進化樹。另外，使用RaxML-VI-HPC[20]和MrBayes 3.1.2[21]構建基于美麗小條鰍COI基因單倍型數(shù)據(jù)和Glyt基因等位基因型數(shù)據(jù)的系統(tǒng)發(fā)育樹。利用MRMODELTEST 2.3基于赤池信息量原則(AIC)選擇最優(yōu)的核苷酸替代模型(COI:GTR+G;Glyt:K80)。最大似然樹自展重復數(shù)設置為1 000次；貝葉斯系統(tǒng)發(fā)育樹運算一共運行了2×108代，每1 000代取樣1次，舍棄前25%的運行代數(shù)。同時基于不同基因類型，分別利用MEGA估算了不同進化譜系內(nèi)部和不同進化譜系之間的遺傳距離。考慮到核基因進化速率慢的特點，本研究利用POPART 1.7.2繪制了不同譜系之間的等位基因網(wǎng)狀圖[22]，用以查看不同線粒體譜系之間在核基因水平上的進化關系。核基因數(shù)據(jù)的分型在DnaSP 5.10軟件中運行。

利用三種方法對美麗小條鰍物種是否存在隱存物種進行界定。首先，采用Automatic Barcode Gap Discovery(ABGD)[23]基于K2P和JC69距離模型查看了美麗小條鰍的可操作分類單元(Operational taxonomic unit)數(shù)目。此外，基于單倍型數(shù)據(jù)構建的最大似然樹，使用Poisson Tree Process(PTP)[24]算法評估美麗小條鰍有效可操作分類單元的數(shù)量。最后，利用the General Mixed Yule-coalescent(GMYC)[25]對物種有效性進行評估。因為GMYC分析需要超度量基因樹(Ultrametric gene tree)，因此利用BEAST version 1.8[26]基于單倍型數(shù)據(jù)構建了超度量基因樹。馬爾科夫鏈長10 000萬代，采用GTR+G模型、Speciation：Yule Process 假設，每1 000代建立一棵系統(tǒng)發(fā)育樹；舍棄前25%后將所有分枝平均后驗概率最高者選出，在FigTree中顯示最后所得樹形圖。PTP(http://mptp.h-its.org)和GMYC(http://species.h-its.org/gmyc/)分析分別在網(wǎng)頁版服務器上進行運算，參數(shù)采取默認設置參數(shù)。

為了進一步檢驗美麗小條鰍不同譜系是否達到物種級別，本研究聯(lián)合COI基因和Glyt基因，基于貝葉斯系統(tǒng)發(fā)育和系統(tǒng)地理的方法(Bayesian Phylogenetics and Phylogeography,BPP)綜合評估隱存物種有效性[27]。這個方法考慮了非單系性和不完全譜系分選等因素[27]。分別設置了三組祖先種群規(guī)模(?)和根齡(τ)，即?=(1,10)，τ=(1,10)；?=(2,2 000)，τ=(2,2 000)和?=(1,10)，τ=(2,2 000)。分析在BPP version 3.0軟件中運行，共運行100萬代，每5代取樣一次，舍棄前5萬代。選取基于COI基因單倍型數(shù)據(jù)獲得的系統(tǒng)發(fā)育拓撲結構作為引導樹。

2 結果

2.1 COI和Glyt序列特征分析

本研究成功獲得了107和47尾美麗小條鰍的COI和Glyt基因序列。在序列兩端截齊后，共獲得610 bp的COI基因序列，其中包含了60個變異位點和55個簡約信息位點；909 bp的Glyt基因序列，其中包含了17個變異位點和12個簡約信息位點。COI基因和Glyt基因分別界定了24個單倍型和27個等位基因型。

2.2 系統(tǒng)樹和遺傳距離

基于107尾美麗小條鰍COI基因序列構建的鄰接樹發(fā)現(xiàn)華南地區(qū)美麗小條鰍形成了2個進化譜系(Ⅰ和Ⅱ)，其中譜系Ⅰ由珠江水系、南流江和漠陽江群體組成(即大陸群體)，譜系Ⅱ由海南島群體組成(圖1；圖2)。另外，譜系Ⅰ擁有三個子譜系(Ⅰ-1，Ⅰ-2和Ⅰ-3)，其中Ⅰ-1主要包括珠江水系的西江和南流江群體，Ⅰ-2和Ⅰ-3分別由珠江水系的東江群體和漠陽江群體組成；譜系Ⅱ含有兩個子譜系，分別由陵水河群體和南渡江群體組成?；贑OI基因單倍型的最大似然樹和貝葉斯樹獲得了相近的系統(tǒng)發(fā)育結構，支持美麗小條鰍存在兩個主要譜系(Ⅰ和Ⅱ；圖3a)。然后，基于Glyt基因單倍型的最大似然樹和貝葉斯樹并未獲得譜系分化的拓撲結構(圖3b)，但是均強烈支持大陸群體、海南島的陵水河群體和南渡江群體均是獨立進化的。等位基因型網(wǎng)狀圖表明基于COI基因獲得的譜系Ⅰ和Ⅱ、譜系Ⅱ的兩個子譜系Ⅱ-1和Ⅱ-2均不存在共享等位基因型的現(xiàn)象(圖4)。

圖3 基于美麗小條鰍COI基因單倍型數(shù)據(jù)(a)和Glyt等位基因型數(shù)據(jù)(b)的最大似然樹和貝葉斯樹(僅展示最大似然樹)

圖4 Glyt等位基因型網(wǎng)狀圖

基于COI基因的美麗小條鰍種內(nèi)平均遺傳距離是1.75%，譜系Ⅰ和Ⅱ內(nèi)部的平均遺傳距離分別是0.61%和1.53%，譜系之間的平均遺傳距離是5.83%，其中子譜系Ⅱ-1和Ⅱ-2之間的平均遺傳也達到2.92%?；贕lyt基因的美麗小條鰍種內(nèi)平均遺傳距離是0.13%，譜系Ⅰ和Ⅱ內(nèi)部的平均遺傳距離分別是0.06%和0.38%，譜系之間的平均遺傳距離是0.37%。

2.3 基于DNA條形碼的物種界定

ABGD分析表明當種內(nèi)差異先驗值P在0.001 7～0.035 9之間時均支持3個可操作分類單元，分別是海南島的陵水河群體和南渡江群體以及珠江-漠陽江-南流江群體(圖5)。bPTP分析獲得結果與ABGD的結果一致，均支持3個可操作分類單元(圖5)。此外，單閾值和多閾值GMYC分析均顯示美麗小條鰍存在5個可操作分類單元(圖5)。

圖5 物種界定分析結果。系統(tǒng)發(fā)育樹拓撲結構由BEAST構建

2.4 貝葉斯物種界定

貝葉斯系統(tǒng)發(fā)育和系統(tǒng)地理(BPP)分析發(fā)現(xiàn)基于不同的祖先種群規(guī)模(?)和根齡(τ)，DNA條形碼基因界定的譜系Ⅰ、譜系Ⅱ的兩個子譜系(Ⅱ-1和Ⅱ-2)均獲得100%的后驗概率，表明這譜系Ⅰ與譜系Ⅱ的兩個子譜系可以上升為三個不同物種(圖6)。

圖6 貝葉斯物種界定結果

3 討論

3.1 美麗小條鰍多樣性的地理分布

已有研究表明，DNA條形碼不僅可以實現(xiàn)物種鑒定，還能夠劃分不同地理種群的分布界限[28-30]。本研究基于不同方法構建的系統(tǒng)發(fā)育樹均一致顯示美麗小條鰍存在兩個主要譜系(Ⅰ和Ⅱ)，其中譜系Ⅰ和譜系Ⅱ分別包含三個和兩個子譜系。丘城鋒等[8]和慶寧等[7]關于華南西部和海南島美麗小條鰍的遺傳變異研究也得到了類似的結果。另外，研究發(fā)現(xiàn)子譜系Ⅰ-1主要由西江、北江、南流江群體組成，Ⅰ-2由東江群體組成，Ⅰ-3由漠陽江群體組成。譜系Ⅱ包含兩個子譜系，其中子譜系Ⅱ-1和Ⅱ-2分別由陵水河群體和南渡江群體組成，呈現(xiàn)嚴格的地理分布格局(圖1和圖2)。上述研究結果表明美麗小條鰍擁有豐富的遺傳多樣性，同時也證實DNA條形碼能夠用來劃分美麗小條鰍不同的地理群體。美麗小條鰍的譜系地理分布格局主要受到華南地區(qū)山脈和海峽阻隔導致(圖1)。例如，譜系Ⅰ和Ⅱ的分布區(qū)被瓊州海峽阻隔，子譜系Ⅰ-3與其他兩個子譜系的分布區(qū)被云開山脈隔離，子譜系Ⅰ-1和Ⅰ-2的分布區(qū)被海水阻隔，子譜系Ⅱ-1和Ⅱ-2的分布區(qū)被五指山和鸚哥嶺阻隔。早期研究也證實華南地區(qū)的山脈和海峽塑造了該分布區(qū)其他魚類的遺傳多樣性地理分布格局，如三角魴(Megalobramaterminalis)[6]、海南鲌(Culterrecurviceps)[31]、福建紋胸(Glyptothoraxfokiensis)[3]、唐魚(Tanichthysalbonubes)[15]等。另外，鰍類較弱的擴散能力也是美麗小條鰍形成譜系分化的一個重要因素。研究表明，鰍類營底棲生活，游泳能力弱，容易在局域形成顯著的遺傳分化甚至進化成多個譜系[32]?？紤]到核基因進化速率慢的特點，基于Glyt基因的系統(tǒng)發(fā)育樹雖然未獲得高支持率的拓撲結構，但是等位基因網(wǎng)狀圖表明大陸群體與海南島群體、海南島的南渡江群體(DZ)和陵水河群體(SX)不存在共享等位基因型的現(xiàn)象，表明大陸群體與海南島群體以及海南島的陵水河群體與南渡江群體在線粒體和核基因水平上都是獨立進化的。

3.2 隱存種

分子測序技術的介入促使越來越多的隱存物種被發(fā)現(xiàn)，為生物多樣性的保護提供了重要數(shù)據(jù)[33]。DNA條形碼技術的一個重要作用是發(fā)現(xiàn)新種或隱存種，并且在近二十年得到了廣泛應用[34,35]。本研究發(fā)現(xiàn)譜系Ⅰ和譜系Ⅱ之間的平均遺傳距離達到了5.64%，子譜系Ⅱ-1和Ⅱ-2之間的平均遺傳距離也達到了2.92%，大于DNA條形碼研究中物種界定閾值(2%)[36]?；贒NA條形碼基因的物種界定分析結果也認為美麗小條鰍至少應該分為三個物種。另外，考慮了非單系性和不完全譜系分選等因素的BPP分析能夠綜合多個基因信息對不同譜系進行物種有效性的驗證，結果也支持美麗小條鰍應該上升為三個獨立的物種這個結論。加上嚴格的地理分布格局，認為華南地區(qū)的大陸群體、海南島的陵水河群體和南渡江群體應該看做三個獨立的物種。瓊州海峽、五指山與鸚哥嶺的地理隔離作用可能是導致美麗小條鰍隱性成種的重要因素。同樣的現(xiàn)象在華南地區(qū)的福建紋胸[3]和唐魚[15]也有所體現(xiàn)。例如，CHEN等[3]認為福建紋胸亞種(Glyptothoraxfokiensis)和海南紋胸亞種(Glyptothoraxfokiensishainanensis)應當上升為兩個種[8]，LI等[15]發(fā)現(xiàn)華南地區(qū)至少分布有8個唐魚物種。為尋找更多關于美麗小條鰍隱性成種的證據(jù)，后期研究應當開展傳統(tǒng)的形態(tài)度量學與基于更多分子標記的相關分析。

3.3 保護管理建議

通過分子進化手段揭示物種的遺傳多樣性的地理分布格局對于物種的精準管理與保護具有重要指導意義。本研究發(fā)現(xiàn)美麗小條鰍至少存在三個隱存種，應該當做三個物種有區(qū)別開展管理保護工作。另外，譜系Ⅰ擁有三個符合地理分布模式的子譜系，表明它們可能是三個獨立進化的種群，可以認為是三個不同的演化顯著單位(evolutionary significant unit，ESU)[37]。在未來的保護工作中，這三個演化顯著單元應該得到有針對性的管理和保護。由于美麗小條鰍廣泛分布于我國華南地區(qū)諸多水系江段，考慮到本研究調(diào)查站位僅覆蓋部分區(qū)域，因此該物種的遺傳多樣性可能被嚴重低估。為了更為科學地提出保護管理建議，更多樣本量、更高覆蓋度的相關研究亟待開展。

4 結論

本研究利用DNA條形碼COI基因和Glyt基因揭示了華南地區(qū)美麗小條鰍的遺傳多樣性分布以及隱存種問題，發(fā)現(xiàn)美麗小條鰍至少存在兩個主要進化譜系(Ⅰ和Ⅱ)，譜系Ⅰ存在三個具有嚴格地理分布的子譜系(Ⅰ-1，Ⅰ-2和Ⅰ-3)，譜系Ⅱ擁有兩個嚴格地理分布的子譜系(Ⅱ-1和Ⅱ-2)。同時研究證明譜系Ⅰ以及譜系Ⅱ的兩個譜系應當看做三個獨立物種。本研究認為美麗小條鰍的多樣性可能被嚴重低估，亟需使用形態(tài)度量學、高覆蓋度的樣品和多種類型的分子標記來全面系統(tǒng)揭示其隱藏的遺傳多樣性，為美麗小條鰍全面精準的管理保護提供數(shù)據(jù)支持。