楊銀盆 楊慧蘭 王銀肖 譚慧敏 陳詠霞

(河北大學生命科學與綠色發(fā)展研究院, 保定 071002)

生物多樣性及物種的分布格局都與冰期及期間形成的避難所密切相關(guān)[1]。對當今物種分布影響較大的是始于3—2 Ma前, 結(jié)束于0.02—0.01 Ma前的第四紀冰期[2]。第四紀冰期的最大特點是冰期和間冰期往復震蕩, 氣候出現(xiàn)冷暖周期性變化, 強烈地影響了種群的歷史動態(tài)和群落結(jié)構(gòu)[3]。冰期來臨時, 棲息于北部的溫熱生物向赤道方向遷移, 或沿著冰緣帶遷移, 或躲藏于避難所; 此時, 海平面下降, 大陸架裸露, 為不同區(qū)系物種交流提供了通道[4]。間冰期時, 隨著氣候變得溫暖和濕潤, 適應(yīng)于冷水的生物則向北遷移, 并擴大其分布范圍[5]。第四紀冰期的這種往復變動引起種群分布范圍及種群數(shù)量的擴張和收縮, 種群遺傳結(jié)構(gòu)受其歷史動態(tài)變動而發(fā)生改變。對現(xiàn)生生物種群遺傳結(jié)構(gòu)及其地理分布的研究, 了解第四紀冰期對現(xiàn)生生物種群結(jié)構(gòu)的影響以及生物避難所的確定是當前生物學領(lǐng)域研究的熱點問題, 可以為生物多樣性保護提供依據(jù)。

北鰍(Lefua costata)隸屬于鯉形目(Cypriniformes)、條鰍科(Nemacheilinae), 為古北區(qū)一種冷水性魚類, 在中國、韓國、日本、蒙古及俄羅斯東部等均有分布[6], 其中, 我國主要分布于淮河以北的東部地區(qū)。北鰍屬現(xiàn)知物種有7種, 我國僅北鰍1種。

北鰍為小型底棲魚類, 遷徙能力差, 對環(huán)境依存度高, 多見于溪流的上游, 喜低溫, 受人類活動影響小, 其分布格局可以很好地反映環(huán)境與地質(zhì)演變過程, 是探討冰期氣候變遷歷史的理想物種。已有研究探討了北鰍屬的起源、分化及屬內(nèi)物種的親緣關(guān)系?；诰€粒體全基因組序列, 北鰍屬大約在26 Ma前與高原鰍屬等條鰍科魚類發(fā)生分化[7]。Mihara等[1]認為北鰍屬魚類的祖先起源于亞歐大陸北部, 向南經(jīng)北海道擴散至朝鮮半島, 基于線粒體D-loop區(qū)序列分析認為北鰍屬的分歧時間在(1.91±0.64) Ma前; 其中, 北鰍和短體北鰍(Lefua nikkonis)的分歧時間在(1.50±0.58) Ma前, 且二者的親緣關(guān)系最近。而同樣基于D-loop區(qū)序列分析, Shedko等[8]則認為北鰍屬的分歧時間在11—12 Ma前, 北鰍的分歧時間在3 Ma前。目前, 對于北鰍遺傳結(jié)構(gòu)及生物地理學格局研究比較匱乏?；诰€粒體D-loop區(qū)序列, 本文分析我國北鰍的遺傳多樣性, 并對冰期避難所及北鰍間冰期可能的擴散路徑進行相關(guān)推測。

1 材料與方法

1.1 樣本采集

樣本于2016—2019年采自9個水系18樣點(圖1),共211尾, 所分析的樣本涵蓋了我國北鰍的主要分布范圍, 按地理分布劃分為7個組: 淮黃河(黃河水系和淮河水系, HH種群)2個樣點: 山東省西營鎮(zhèn)錦繡川水庫、河南省汝州市蟒川鎮(zhèn)汝河, 28尾; 海灤河(灤河水系和海河水系, HL種群)7個樣點: 河北省閃電河、河北省府河、河北省拒馬河、天津張莊水庫、河南省獅豹頭水庫、河南省小南海水庫、河南省塔崗水庫, 共89尾; 遼河水系(LH種群)3個樣點: 內(nèi)蒙古赤峰市西遼河、遼寧省彰武縣柳河、遼寧省朝陽市大凌河, 共29尾; 遼東半島水系(LD種群)2個樣點: 遼寧省瓦房店復州河、遼寧省普蘭店沙河, 共33尾; 大洋河水系(DY種群)1個樣點: 遼寧省圓臺大沙鎮(zhèn)大洋河, 共19尾; 鴨綠江水系(YL種群)2個樣點: 遼寧省鳳城和吉林省臨江, 共5尾; 松花江水系(SH種群)1個樣點: 黑龍江省樺南縣八虎力河, 共8尾, 后文中不同種群均簡稱表示。標本保存于河北大學博物館, 儲存于-20℃冰箱保存。

圖1 北鰍樣本采集地點Fig.1 Map of the sample collection sites of Lefua costata

1.2 DNA提取、PCR擴增和測序

樣本取自魚的胸鰭、背鰭或肌肉組織, 剪碎后使用基因組提取試劑盒DP304(北京天根生化科技有限公司)進行提取。D-loop區(qū)擴增的引物序列為DL1 (3′-ACCCCTGGCTCCGAAAGC-5′)和DH2 (3′-ATCTTAGCATCTTCAGTG-5′)[9]。PCR反應(yīng)條件為: 95℃預變性3min; 95℃變性30s, 退火30s溫度為51℃, 72℃延伸40s, 循環(huán)30次; 最后72℃延伸10min,4℃停止反應(yīng)。PCR產(chǎn)物用1.00%—2.00%瓊脂糖凝膠電泳對其條帶進行檢測, 送往通用生物系統(tǒng)(安徽)有限公司完成序列測定。

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

利用BioEdit 7.0.9[10]軟件進行DNA序列比對和拼接; 利用MEGA 5.05[11]軟件對序列變異情況、堿基組成進行分析, 并計算轉(zhuǎn)換/顛換比率等相關(guān)系數(shù), 基于Kimura-2參數(shù)模型計算各地理群體間的遺傳距離。利用DNAsp 5.1[12]軟件進行遺傳多樣性相關(guān)分析。利用Arlequin 3.5[13]軟件計算群體間分化系數(shù)(FST), 以及中性檢驗及錯配分析等種群歷史動態(tài)分析。

以細尾高原鰍(Tripiophysa stenura)、董氏須鰍(Barbatula toni)、短體北鰍(L.nikkonis)、斑北鰍(L.echigonia)和L.sp.作為外來群構(gòu)建北鰍的系統(tǒng)發(fā)育樹, 下載序列GenBank登錄號見表1。利用jModel-Test 0.1.1[14]軟件篩選構(gòu)建系統(tǒng)發(fā)育樹的最佳進化模型, 為HKY+G模型。利用PhyloSuite[16]中的MrBayes程序[17]構(gòu)建單倍型貝葉斯(BI)系統(tǒng)發(fā)育樹, 4條馬爾科夫鏈同時運行1×106代, 通過后驗概率評估每個節(jié)點, 每100代對系統(tǒng)樹進行抽樣。利用IQtree[15]構(gòu)建單倍型最大似然(ML)系統(tǒng)發(fā)育樹, 對每個節(jié)點的置信度對評估通過有1000個重復的bootstrap進行。利用PopART[18]軟件采用TCS法構(gòu)建單倍型網(wǎng)絡(luò)關(guān)系圖。利用BEAST 2.0[21]軟件推算分歧時間, 采用花鰍屬(Cobitis)魚類化石(約23 Ma前)[19]和高原鰍屬(Tripiophysa)魚類化石(約13 Ma前)[20]進行時間校準, 下載序列GenBank登錄號見表1。種群擴張時間的評估依據(jù)線粒體D-loop區(qū)每百萬年1.00%的平均進化速率計算[8], 公式為τ=2ut[22],u為每條序列每個世代的變異速率, 根據(jù)u=2μk計算,μ為平均進化速率,k為序列長度。北鰍屬魚類為一次產(chǎn)卵類型,因此按1世代/年計算[23]。

表1 GenBank下載序列Tab.1 List of the sequences from GenBank

2 結(jié)果

2.1 序列特征與遺傳多樣性

9條水系211條北鰍D-loop區(qū)序列平均長度為925 bp: 保守位點823個, 約占序列總長度的88.97%;變異位點102個, 約占序列總長度的11.02%; 簡約信息位點61個, 約占序列總長度的6.70%。序列轉(zhuǎn)換顛換比(Ti/Tv)為1.88。堿基含量T為32.50%, C為18.61%, A為35.90%, G為12.99%, 表現(xiàn)為G偏移。A+T含量(68.40%)明顯高于C+G含量(31.60%)。

在211個樣本中檢測到55個單倍型, 地理種群間共享單倍型有5個, 單倍型多樣性(h=0.9304)和核苷酸多樣性(π=0.0087)均較高(表2)。7個種群中LD種群單倍型多樣性和核苷酸多樣性均最高(h=0.8920,π=0.0074), SH種群(h=0.2500,π=0.0006)最低。7個種群的群體間遺傳分化指數(shù)(FST)為-0.0026—0.6870(表3), 除LH種群與SH種群之間FST不顯著外, 其余種群間均有一定程度的分化,FST均為顯著(P≤0.01), 其中LD種群與其他地理種群之間FST均較高, 為0.4853—0.6679。7個種群間遺傳距離(K2P)為0.0019—0.0158, 其中LD種群與其他地理種群間的遺傳距離均大于其他地理種群間的遺傳距離, 為0.0139—0.0158, 其次是YL種群與其他地理種群間遺傳距離較大, 為0.0079—0.0139,DY種群與其他地理種群間遺傳距離最小, 為0.0030—0.0158; 種內(nèi)遺傳距離最大的是LD種群 (0.0074),其次是YL種群(0.0053), SH種群最小(0.0006)。7個種群內(nèi)的遺傳差異(66.61%)大于種群間(33.39%);根據(jù)LD種群與其他6個種群間的分化程度, 將其與其他6個種群劃分為2個組分析顯示, 組間遺傳差異(30.12%)大于組內(nèi)遺傳差異(29.49%), 而種群內(nèi)遺傳差異最大, 為70.51%(表4)。

表2 北鰍遺傳多樣性、中性檢驗、錯配分析與擴張時間Tab.2 Genetic diversity, neutrality tests, mismatch distribution and population expansion time of Lefua costata

表3 北鰍群體間的遺傳分化系數(shù)FST(對角線下方)、群體間的遺傳距離(對角線上方)及群體內(nèi)遺傳距離(對角線上)Tab.3 Pairwise FST (below diagonal), genetic distances among populations (above diagonal) and genetic distances within populations (on diagonal) of Lefua costata

2.2 系統(tǒng)發(fā)育與單倍型網(wǎng)絡(luò)圖分析

以細尾高原鰍和董氏須鰍為外來群, 構(gòu)建的ML和BI單倍型系統(tǒng)發(fā)育樹, 其拓撲結(jié)構(gòu)基本一致(圖2)。我國北鰍9條水系的55個單倍型以高置信度(ML=96, BI=1)聚為一支, 并與短體北鰍形成姐妹群關(guān)系。9條水系的單倍型未能各自聚類, 而是嵌套分布, 形成3大分支(分支A、B和C)。分支A主要以LD種群的15個單倍型為主, 嵌套分布有HH種群單倍型Hap1和LH種群單倍型Hap3, 位于北鰍系統(tǒng)樹的根部; 其余單倍型以略高的置信度(ML=91,BI=0.75)聚為一支, 形成分支B和分支C兩個姐妹群,其中分支C主要以HH種群和HL種群為主, 嵌套分布于有LD種群單倍型Hap41和YL種群單倍型Hap5; 分支B則包含所有種群的單倍型。

表4 群體內(nèi)和群體間變異的分子方差分析(AMOVA)Tab.4 Analysis of molecular variance (AMOVA) among populations and within the population of Lefua costata

單倍型網(wǎng)絡(luò)圖(圖3)顯示55個單倍型整體呈復雜的網(wǎng)狀結(jié)構(gòu), LD種群Hap39、Hap36位于網(wǎng)狀中心, 表明7個種群以LD種群為中心經(jīng)歷過復雜的基因交流, 結(jié)合系統(tǒng)發(fā)育樹的結(jié)果, 將其分成3大分支。分支A以LD種群為主, 呈星狀或輻射狀, 單倍型Hap31位于中心位置, 表明此單倍型的祖先經(jīng)歷過種群擴張; 分支B包含所有種群的單倍型和種群間共享單倍型Hap3、Hap6、Hap8和Hap9, 形成復雜的網(wǎng)狀結(jié)構(gòu), 共享單倍型均處于網(wǎng)狀結(jié)構(gòu)中心,表明7個種群經(jīng)歷過復雜的基因交流; 分支C以HH種群和HL種群為主, 呈星狀或輻射, 單倍型Hap11位于中心位置。基于化石校準點推算出分支A與分支(B+C)的分歧時間約為2.6241 Ma (95%HPD:0.9264—4.6394 Ma)前, 分支C和分支B的分歧時間約為1.5808 Ma (95%HPD: 0.4686—3.0695 Ma)前(圖4)。

圖2 北鰍單倍型ML系統(tǒng)發(fā)育樹(節(jié)點上方為ML支持率, 下方為BI后驗概率)Fig.2 ML phylogenetic tree constructed by Lefua costata haplotype (bootstrap support for ML above the branch and posterior probability for Bayesian inferences below branch)

圖3 北鰍單倍型網(wǎng)絡(luò)進化圖Fig.3 Haplotype evolution network diagram of Lefua costata

2.3 種群動態(tài)分析

中性檢驗分析顯示, Tajima’sD值僅LH種群呈顯著負值(P＜0.05, 表2); Fu’sFS值僅LD種群和DY種群呈顯著負值, 表明LH種群、LD種群和DY種群可能經(jīng)歷過種群擴張事件。錯配分析中,Hri檢驗結(jié)果均為不顯著(P＞0.05),SSD檢驗結(jié)果僅LH種群顯著偏離了種群擴張模型, 表示7個種群均經(jīng)歷過種群擴張事件。錯配分布圖顯示(圖5), DY種群呈單峰分布曲線; HH種群呈波浪多峰曲線; HL種群、LH種群、LD種群、YL種群和SH種群呈鋸齒狀下降曲線, 表明DY種群經(jīng)歷過種群擴張事件,其他地理種群未出現(xiàn)種群擴張現(xiàn)象。參考北鰍屬魚類D-loop區(qū)序列1.00%/百萬年的平均進化速率,推算出北鰍不同種群擴散時間為2.3784—0.0477 Ma前(表2)。

圖4 基于松散分子鐘構(gòu)建的北鰍分歧時間Fig.4 Divergence times among Lefua costata derived from the Bayesian relaxed-molecular

3 討論

3.1 北鰍種群遺傳多樣性和種群分化

現(xiàn)生生物種群的遺傳結(jié)構(gòu)和分布區(qū)系與第四紀冰期和間冰期的往復更迭造成的冷暖氣候變化密切關(guān)聯(lián), 冰期和間冰期引起的海進和海退為不同區(qū)系物種遷徙以及基因交流提供通道。通常在冰期降溫時, 溫熱性魚類往往會縮小其分布范圍, 而冷水性魚類則會擴大其分布范圍, 而間冰期氣溫回暖時, 溫熱性魚類的分布范圍會擴大, 冷水性魚類的分布范圍則會縮小[24]。種群遺傳結(jié)構(gòu)分析能反映生物對自然選擇及環(huán)境變遷的適應(yīng)能力和演化潛力[25], 種群遺傳多樣性越高, 適應(yīng)環(huán)境能力越強,承受的生存壓力越小。群體遺傳多樣性衡量指標則有賴于核苷酸多樣性指數(shù)和單倍型多樣性指數(shù)[26]。種群經(jīng)歷過小種群事件, 如奠基者效應(yīng)或瓶頸效應(yīng),種群遺傳結(jié)構(gòu)通常表現(xiàn)為低的遺傳多樣性[27]; 而種群經(jīng)歷過長期積累遺傳分化, 如種群處于物種起源中心或冰期避難所, 種群遺傳結(jié)構(gòu)則表現(xiàn)為高的遺傳多樣性[7,28]。在本研究中, 核苷酸多樣性和單倍型多樣性最高的是遼東半島種群(LD), 距離遼東半島越遠的種群, 其種群的單倍型多樣性與核苷酸多樣性也普遍越低, 如距離遼東半島較近的大洋河種群(DY)和鴨綠江種群(YL)遺傳多樣性較高, 僅次于遼東半島種群(LD), 距離遼東半島較遠的松花江種群(SH)遺傳多樣性則最低。種群分歧時間推算顯示我國北鰍7個地理種群中, 遼東半島種群(LD)與其他地理種群分化時間最早, 這表明我國北鰍種群起源地之一為遼東半島或遼東半島種群(LD)受第四紀冰期影響較小, 該種群經(jīng)歷過長期遺傳分化的積累, 其他地理種群則是以遼東半島種群為中心在第四紀冰期和間冰期時向南或向北擴散類群。

從種群間的遺傳分化程度來看, 遼東半島種群(LD)與其他地理種群的遺傳分化程度最為顯著, 且距離遼東半島近的種群與其他地理種群的分化程度高, 而距離遼東半島遠的種群與其他地理種群的分化程度低, 這表明我國北鰍以遼東半島為中心在末次冰期向南擴散的種群和間冰期向北擴散的種群積累的遺傳變異最少。北鰍為底棲小型魚類, 遷移能力弱, 且為冷水性魚類, 多生活在山澗溪流中,對環(huán)境依存度高, 在末次冰期后不同地理群體之間無基因交流, 表現(xiàn)為群體內(nèi)遺傳變異大于群體間遺傳變異, 符合距離隔離模型(Isolation by distance,IBD)理論[29], 即遠距離分開的不同個體群由于有限距離的基因遷移而產(chǎn)生的遺傳分化[30]。

3.2 北鰍種群歷史動態(tài)

北鰍屬魚類的祖先種最早起源于亞歐大陸北部, 始于漸新世約26 Ma前[7], 經(jīng)東亞北部穿過對馬海峽和韃靼海峽遷移至日本相關(guān)島嶼, 向南經(jīng)北海道擴散至朝鮮半島[1]。在我國, 與朝鮮半島相鄰的吉林省延邊圖們江有較北鰍物種原始的新記錄種——史氏北鰍(Lefua pleskei) (待發(fā)表), 及我國北鰍種群中與朝鮮半島距離最近的遼東半島種群(LD)和鴨綠江種群(YL)具有較高的遺傳多樣性, 而越往北和往南北鰍種群遺傳多樣性降低, 表明我國北鰍屬魚類的起源地是由朝鮮半島遷徙而來。

圖5 核苷酸錯配分布圖Fig.5 The nucleotide of mismatch distribution

北鰍分化時間約在3 Ma前[9], 處于晚第三紀上新世晚期, 我國北鰍魚類主要分布于朝鮮半島鄰近地區(qū), 如鴨綠江和遼東半島地區(qū)等地。在第四紀冰期時, 由于遼東半島東側(cè)地勢由東北向西南漸低的長白山余脈千山山脈減輕了冰期氣候?qū)ι锏挠绊? 為生物提供了避難場所, 因此分布于遼東半島的北鰍魚類受冰期影響小, 能夠長期積累遺傳分化,具有較高的單倍型多樣性和核苷酸多樣性; 其他地理種群隨冰緣帶借助黃渤海平原[31]的河口向南及周圍擴散, 南遷至海河、黃河等, 西遷至大凌河、灤河等。文獻報道, 渤海是第四紀才下陷的, 在渤海下陷之前, 古黃河、古遼河及兩河中間水系曾為同一水系[32], 渤、黃海經(jīng)歷幾次海進海退, 為北鰍魚類在冰期和間冰期遷徙及基因交流提供了通道。本文研究與地質(zhì)歷史事件相吻合, 以遼東半島種群(LD)為主的分支A分化時間(4.6394—0.9264 Ma前)發(fā)生在第一次冰期期間[33], 為我國北鰍較原始的類群, 分支A所含有的淮黃河種群(HH)單倍型Hap4、遼河種群(LH)單倍型Hap27及海灤河種群(HL)共享單倍型Hap19, 表明遼東半島種群(LD)往南往西擴散至黃河、灤河和遼河干流等。

冰期避難所存在多個[34], 北京、河北等地在第四紀冰川期經(jīng)歷過四次持續(xù)性且分布范圍相似的冰期[32], 燕山山脈東側(cè)由灤河堆積形成的霧靈山為末次冰期的避難所之一[33,35]。本文研究支持霧靈山為冰期避難所之一。以海灤河種群(HL)為主的分支C分化時間(3.0695—0.4686 Ma前)約在第一次間冰期和第二次冰期, 結(jié)合海灤河種群(HL)與遼東半島種群(LD)的共享單倍型Hap19(位于分支A), 與大洋河種群(DY)、淮黃河種群(HH)的共享單倍型Hap3和Hap9(位于分支B), 及海灤河種群(HL)較高的核苷酸多樣性和較低的單倍型多樣性, 表明海灤河種群(HL)是在第一次冰期或間冰期由遼東半島種群(LD)擴散而來, 在以后的幾次冰期和間冰期中經(jīng)歷了瓶頸效應(yīng)后種群迅速擴張; 同時與大洋河種群(DY)、淮黃河種群(HH))的基因交流, 表明海灤河種群曾發(fā)生過不同地理種群的二次接觸。

末次冰期后, 北鰍往北擴大其分布范圍才遷徙至西遼河、東遼河以及松花江水系。因為在早更新世和中更新世期間, 遼河古水系(東遼河和西遼河)及松花江古水系流向松遼平原北部的松遼古大湖[36,37], 晚更新世時, 由于松遼分水嶺的上升和形成, 松遼古大湖消失[38—41], 末次冰期黃、渤海海平面下降, 大陸架裸露, 遼河干流襲奪了東遼河、西遼河, 形成遼河水系[31,42]; 另外松遼分水嶺的形成使得原本作為西遼河上源的嫩江與第一松花江合并后在三江平原與黑龍江干流合并[31], 最終形成現(xiàn)今松花江水系格局。在末次冰期后, 隨著遼河水系的形成和松花江的改道為北鰍往北遷徙提供了通道, 最終形成現(xiàn)今的遼河種群(LH)和松花江種群(SH)?；春铀档谋宾q起源于黃河水系, 是在黃河奪淮[43]事件后, 淮河水系才有了北鰍的分布, 由于隔絕時間短, 黃河水系和淮河水系北鰍仍存在基因交流。