基于EST-SSR標(biāo)記的油松種子園遺傳多樣性分析

代劍峰，黨 磊，金雨晴，周先清，賈海寬，張鴻景

（1.河北省林業(yè)科學(xué)研究院，河北 石家莊 050061；2.河北省林木良種技術(shù)創(chuàng)新中心，河北 石家莊 050061；3.北京林業(yè)大學(xué)，北京 100083；4.平泉市國有七溝林場，河北 承德 067500；5.呼倫貝爾市紅花爾基林業(yè)局國家樟子松良種基地, 內(nèi)蒙古 呼倫貝爾 021112）

油松（Pinus tabulaeformis）是我國特有的針葉樹種，其天然林和人工林廣泛分布于國內(nèi)的 14個?。▍^(qū)）［1］，在木材生產(chǎn)和生態(tài)保護(hù)上占據(jù)重要地位。自20世紀(jì)70年代起，我國開始建立油松良種初級種子園，目前包括河北在內(nèi)的很多地區(qū)都已開展有關(guān)高世代油松種子園的研究工作［2-5］。

種子園的遺傳多樣性在一定程度上反映了育種群體配置的合理性以及遺傳改良的能力，隨著種子園向更高世代發(fā)展，如何建立有效的育種群體逐漸成為種子園建設(shè)的首要問題，對種子園遺傳多樣性的研究也越來越迫切［6-12］。利用分子標(biāo)記對種子園的遺傳多樣性進(jìn)行研究是目前較為常用的方法，其中簡單序列重復(fù)（Simple sequence repeat, SSR)標(biāo)記具有可靠性高、操作簡單且在同屬間具備較高通用性等優(yōu)點，廣泛應(yīng)用于各種林木遺傳分析中［13-15］。由于針葉樹基因組非常龐大且復(fù)雜,導(dǎo)致 SSR引物的開發(fā)成本較高［16-17］，利用SSR引物對油松種子園遺傳多樣進(jìn)行分析的研究較為少見，目前只在內(nèi)蒙和山西兩處種子園有過相關(guān)的報道［18-19］。近年來，隨著高通量測序技術(shù)的發(fā)展，利用植物物種的轉(zhuǎn)錄組序列來設(shè)計的微衛(wèi)星分子標(biāo)記引物，即EST-SSR引物，使引物開發(fā)的成本大幅降低［20-22］，能夠批量設(shè)計出大量有效引物用于遺傳多樣性研究，極大地提高了針葉樹遺傳分析的效率［23-26］。

河北省平泉七溝油松良種基地油松初級種子園始建于1984年，面積56 hm2，共使用217個優(yōu)樹無性系。優(yōu)樹接穗主要來源于遼河源大窩鋪、唐山遵化市清東陵、承德平泉市寬城縣和七溝鎮(zhèn)。本試驗以基地內(nèi)來源清晰的153個無性系為研究對象，利用EST-SSR分子標(biāo)記分析它們的遺傳多樣性，探討無性系之間的親緣關(guān)系和分化程度，對油松高世代種子園建設(shè)、遺傳改良策略的制定以及優(yōu)良家系的篩選提供參考。

1 材料方法

1.1 植物材料和基因組DNA提取

油松材料由河北省平泉國有林場良種基地提供，選取基地內(nèi)油松初級種子園中種源較為清晰的153個無性系作為試驗驗樣本，分別來自于4處天然油松群落：大窩鋪（Dawopu）（N41°19′，E118°32′）6個、東陵（Dongling）（N40°11′，E117°37′）97個、寬 城（Kuancheng）（N40°36′，E118°29′）11個 和七 溝（Qigou）（N41°02′，E118°32′）39個。采 集各無性系幼嫩針葉，裝于自封袋中冷藏帶回，樣品儲存于-20 ℃直至提取基因組DNA。

利用CTAB法提取各樣本針葉基因組DNA，將DNA溶于ddH2O后，利用Nano Drop 2000（Nano Drop 2000/2000c, Thermo Scientific, USA）檢測DNA濃度，最終將DNA樣本稀釋至20 ng/μL備用。

1.2 引物篩選及SSR-PCR反應(yīng)

本文利用油松及其近緣物種針葉組織轉(zhuǎn)錄組信息［21,27］，設(shè)計了70對EST-SSR引物用于多態(tài)性引物的篩選，對所有正向引物的5’端進(jìn)行M13修飾，以便與熒光標(biāo)簽（FAM，HEX，TAMRA，ROX）結(jié)合，利用瓊脂糖電泳法對引物進(jìn)行篩選。

PCR采 用20 μL體 系，包 括10 μL 2×TaqPCR master mix (Biomed Technologies, Beijing, China)，4 μL(4 pmol/L)正向/反向混合引物，4 μL (4 pmol/L) 熒光 標(biāo) 簽（Fluorescent-dye-labeled M13 primer），2 μL（20 μg）基因組DNA。PCR擴增條件如下： 94 ℃10 min，94 ℃ 30 s，退火60 ℃ 30 s，每個循環(huán)減0.5 ℃，延伸72 ℃ 45 s，進(jìn)行20個循環(huán)。95 ℃ 30 s，退火50 ℃ 30 s，延伸72 ℃ 45 s，進(jìn)行20個循環(huán)。最終延伸72 ℃ 5 min。

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

試驗采用毛細(xì)管電泳法，將擴增產(chǎn)物用ABI3730測序儀（3730XL, Applied Biosystems,USA）進(jìn)行檢測，檢測結(jié)果利用GeneMarker version 2.20（SoftGenetics, State College, Pennsylvania,USA）讀取，并記錄SSR位點信息（引物gp01對應(yīng)位點記為GP01，以此類推）。每個位點的等位基因數(shù)（Number of observed alleles,N）、近交系數(shù)（Coefficient of inbreeding,Fis） 等利用GenALEx version 6.5［28］計算獲得；每個位點無效等位基因頻率（Invalid allele frequency,Pn）通過FreeNA［29］計算得出；利用SPAGeDi［30］計算出等位基因總數(shù)（Total number of alleles,Nt）、每個位點的平均等位基因數(shù)（Average number of alleles,Na）、等位基因豐富度（Allelic richness,AR）、有效等位基因數(shù)（Effective number of alleles,Nea）、群體獨有等位基因數(shù)（Np）、觀測（Observed heterozygosity,Ho）和期望（Expected heterozygosity,He）雜 合 度；哈 迪—溫 伯 格 平衡（Hardy-Weinberg equilibrium, HWE）通過ARLEQUIN［31］進(jìn)行檢測；群體差異系數(shù)（Population differentiation coefficient,Fst）由GenALEx軟 件 利用多位點比較法（999 permutations）計算得出，群體間差異度（PairwiseFst）也同樣由該軟件進(jìn)行計算；包含無效等位基因以及剔除無效等位基因的群體差異系數(shù)Fst（ENA）通過FreeNA進(jìn)行計算，該參數(shù)用于評估群體間的遺傳分歧以及無效等位基因?qū)z傳差異的潛在影響；NJ Tree圖由SplitsTree軟件構(gòu)建并繪制。

2 結(jié)果與分析

2.1 引物篩選

利用隨機選取8個試驗樣本的混合DNA對70對引物進(jìn)行PCR初篩，共有32對引物擴增出了預(yù)期大小的片段，其余擴增沒有PCR產(chǎn)物。再用這32對引物分別對8個試驗樣本逐一進(jìn)行PCR擴增，其中20對引物擴增出了單態(tài)位點（或多態(tài)性低），另外12對為多態(tài)位點，圖1為引物gp08的擴增結(jié)果。由此，選用這12對引物進(jìn)行后續(xù)試驗，引物信息如表1。

圖1 引物gp08擴增產(chǎn)物檢測Fig.1 Amplified product detection of primer gp08

表1 引物信息Table 1 Information of primers

2.2 無性系及種源群體的遺傳多樣性

利用GeneMarker version 2.20讀取毛細(xì)管電泳結(jié)果（圖2），整理后發(fā)現(xiàn)153個無性系在12個SSR位點上，一共檢測到了62個等位基因，每個位點的等位基因數(shù)為2～10個，平均5.17個，其中在位點GP04檢測到了10個等位基因，位點GP02只檢測到了2個等位基因。各位點無效等位基因頻率（Pn）的檢測結(jié)果顯示，有5個位點（GP2、GP6、GP10、GP11、GP12）幾乎沒有無效等位基因，GP13的Pn值最高達(dá)到了0.162 4，其它位點的Pn值均未超過0.08。期望雜合度（He）的最大值出現(xiàn)在了位點GP01（0.757），所有位點的期望雜合度（He）和觀測雜合度（Ho）的平均值分別為0.414和0.403。有多個位點（GP02、GP04、GP06、GP10、GP11和GP12）出現(xiàn)了Ho大于He的情況，表明在這些位點上存在雜合子過量的情況。各位點的近交系數(shù)（Fis）較小，平均Fis為0.018 6（表2）。

圖2 引物gp01在部分樣本中的擴增條帶Fig. 2 Amplification bands of primer gp01 in partial samples

如表3所示，根據(jù)樣本來源，153個無性系被分為4個種源群體：大窩鋪、東陵、寬城和七溝。在3個群體的9個位點上發(fā)現(xiàn)了17個特殊等位基因（Np），出現(xiàn)頻率從0.005到0.045。4個群體檢測到的等位基因數(shù)（Nt）分別為28（大窩鋪）、55（東陵）、35（寬城）和44（七溝），平均等位基因數(shù)依次為2.33、4.58、2.92和3.67。大窩鋪、東陵和七溝3個群體的平均有效等位基因數(shù)（Nea）、等位基因豐富度（AR）、觀測雜合度（Ho）和期望雜合度（He）基本一致，都略高于寬城群體，這可能與群體樣本數(shù)量以及樣本選擇有關(guān)。

2.3 無性系種源差異及群體分化

153個無性系在12個位點上的群體差異系數(shù)（Fst）都很小，均值僅為0.015 6。去除無效等位基因后，F(xiàn)st）（ENA）也沒有發(fā)生顯著的變化，均值為0.017 2。

同時，對4處種源群體進(jìn)行兩兩比較，PairwiseFst最大值出現(xiàn)在大窩鋪群體和寬城群體之間（0.071 7），最小值在東陵群體和七溝群體之間（0.009 9），均未達(dá)到顯著水平，表明各群體間沒有明顯的遺傳分化（表4）。

表4 4處種源群體的遺傳分化Table 4 Genetic differentiation of groups between 4 provenances

圖3為依據(jù)153個樣本的遺傳距離繪制的NJ 樹圖，從圖中可以看到聚類的結(jié)果并沒有表現(xiàn)出明顯的種源相關(guān)性，大部分情況下都是多個種源的無性系相互混雜。表明雖然種源地不同，但這些無性系之間共享了大部分的遺傳信息。

依據(jù)個體基因型進(jìn)行主成分分析（PCA），前3個成分貢獻(xiàn)了79.54%的變異（分別為73.54%，3.09%和2.91%）。圖4為4個群體的無性系分別在PC1與PC2的二維坐標(biāo)圖（4A）和PC1與PC3的二維坐標(biāo)圖（4B）上的投影，從圖中可以看出4個群體是相互重疊的，很難進(jìn)行區(qū)分。這進(jìn)一步證明所有無性系的親緣關(guān)系較近，沒有明顯的群體界限。

3 討論與結(jié)論

3.1 EST-SSR引物篩選

隨著各物種cDNA文庫的不斷豐富，與傳統(tǒng)SSR引物相比利用轉(zhuǎn)錄組序列設(shè)計EST-SSR引物速度更快、成本更低且具備良好的屬間通用性，因此近年來被廣泛的應(yīng)用于種質(zhì)資源的評價和利用、遺傳多樣性檢測與分析以及種子園交配系統(tǒng)評價等領(lǐng)域［18,25,32］。然而，由于EST-SSR的有效擴增率低、多態(tài)性差，在使用前首先要進(jìn)行大量的引物篩選工作［33］。本試驗從油松及其近緣物種的70對ESTSSR引物中篩選出32對能夠有效擴增的引物，其中12對為多態(tài)性引物，引物多態(tài)率為17.14%。造成EST-SSR引物多態(tài)率較低的原因可能由于松科植物是比較原始的樹種，序列相對保守、變異程度相對較低，此外植物樣本之間的基因差異較小也可能造成EST-SSR引物多態(tài)率較低。本研究用這12對引物對油松種子園的育種群體進(jìn)行標(biāo)記，共檢測出62對等位基因，多樣性信息較為豐富，從基因水平上為油松種子園的多樣性分析提供了參考信息。此次篩選出的引物也可以為其他同類種子園多樣性的研究提供可靠的引物資源。

圖4 主成分分析圖Fig.4 Principal component analysis(PCA)

3.2 油松種子園的遺傳多樣性

種子園的遺傳多樣性在一定程度上反映了育種群體配置的合理性以及遺傳改良的能力，豐富的遺傳多樣性是種子園實現(xiàn)育種目的的基礎(chǔ)。本研究對油松種子園的153個無性系單株遺傳多樣性的分析結(jié)果表明，由這153個無性系單株構(gòu)成的群體在12個多態(tài)性位點上平均等位基因數(shù)為5.17，期望雜合度（He）和觀測雜合度（Ho）分別為0.414和0.403。與程翔等［18］利用EST-SSR對黑里河油松種子園遺傳多樣性進(jìn)行分析時得到的結(jié)果（He=0.366，Ho=0.446）基本一致，高于李悅等［34］利用同工酶標(biāo)記得到結(jié)果（He=0.329，Ho=0.285），又低于袁虎威等［5］利用SSR標(biāo)記得到的結(jié)果（He=0.570 3，Ho=0.595 5）。造成這種差異的原因可能與EST-SSR標(biāo)記的位點多態(tài)性高于同工酶標(biāo)記而低于SSR標(biāo)記有關(guān)。本文結(jié)果表明該處種子園有著較為豐富的遺傳多樣性，具備一定的育種潛力。

3.3 油松種子園無性系的遺傳分化

本試驗選用的153個無性系按來源可以分為4個種源群體：大窩鋪、東陵、寬城和七溝。通過分析發(fā)現(xiàn)它們的群體間平均差異系數(shù)（Fst）為0.015 6，該結(jié)果與馬尾松（0.018 1）［35］和側(cè)柏（0.011）［36］類似，遠(yuǎn)低于樟子松（0.219 4）和挪威云杉（0.749）等［37-38］針葉樹種，說明該種子園的油松群體遺傳分化水平較低，主要差異來源于群體內(nèi)部。同時，NJ樹圖聚類和主成分分析的結(jié)果進(jìn)一步印證了以上觀點。導(dǎo)致這種分化不明顯的原因可能是：該處種子園的種源地都分布于河北省東北部和北部山區(qū)，氣候條件類似，加之人類的頻繁引種活動，最終使得該區(qū)域內(nèi)的油松基因交流頻繁。為了避免近交的發(fā)生，油松種子園在繁育高世代良種時，一方面應(yīng)當(dāng)結(jié)合分子標(biāo)記對親本群體進(jìn)行篩選，另一方面應(yīng)當(dāng)從更廣泛的區(qū)域引進(jìn)油松良種擴充親本群體。

由于病蟲害問題導(dǎo)致種子園連年欠收，一直沒能獲得足夠的子代群體，因此結(jié)果尚需要子代數(shù)據(jù)的進(jìn)一步佐證，后續(xù)的研究中對種子園的交配系統(tǒng)進(jìn)行全面的分析，完善結(jié)論。另外，由于種子園的優(yōu)樹絕大部分都源于遵化市的清東陵，而源于其他種源地的資源較少，導(dǎo)致本文中各種源群體大小出現(xiàn)不平衡，這是否會影響結(jié)果產(chǎn)生偏差有待進(jìn)一步驗證。