貴州武陵山地區19份柚資源的遺傳多樣性

王小柯， 鄭乾明， 羅 懌， 李文云， 柏自琴， 李金強

(貴州省農業科學院 果樹科學研究所，貴州 貴陽 550006)

柚是一類重要的柑橘類果樹，廣泛分布于中國南方及東南亞地區。云貴地區作為內陸性系統柚的多樣化中心之一，擁有極其豐富的地方柚類資源[1]。貴州省內柚資源主要分布于黔東北武陵山區，目前調查發現有酸柚、紅肉柚、白柚等野生資源，經過長期的自然雜交和人工選育，這些野生資源形成了一系列具有優良性狀的地方品種。當前對黔東北武陵山區柚資源尚無系統的研究報道，對其遺傳多樣性缺乏了解。

SNP位點是目前廣泛應用的分子標記，具有分布廣泛、分辨率高、共顯性、易于高通量檢測和分型等優點，可應用于遺傳資源評價和遺傳圖譜構建等[2]。GBS技術是基于第2代測序技術大量快速挖掘SNP的簡化基因組測序方法，主要利用酶切降低基因組復雜度，多樣本高通量平行檢測酶切位點附近的SNP，從而達到基因分型的目的。GBS分析具有同時檢測多個樣本，通量大、費用低的優點，廣泛用于品種識別、遺傳多樣性評價和遺傳圖譜構建等[3-6]。筆者前期采集了貴州省武陵山地區6個區縣共19份柚資源，分析其形態學特征發現，在葉形、花色、果肉色澤上均存在差異。因此結合GBS測序技術，從基因組水平上挖掘SNP位點，進行核苷酸差異分析及群體遺傳多樣性分析。擬通過簡化基因組測序的方法，開發分子標記用于貴州地方柚品種及野生資源品種鑒定以及遺傳多樣性分析，以期為保護柚野生資源和地方品種資源、品種鑒定、挖掘優良性狀奠定基礎。

1 材料與方法

1.1 樣品采集

通過實地考察和走訪群眾，在貴州省銅仁市6個區縣共采集19份柚樣品(表1)。采集時選擇干凈、無病蟲害、幼嫩的春梢葉片。采集后用無菌水清洗2次并擦干，然后液氮速凍。

表1 武陵山地區不同區縣19份柚樣品的信息

1.2 基因組 DNA提取

19份柚樣品葉片總DNA的提取使用植物基因組DNA提取試劑盒[天根生化科技(北京)有限公司，DP305]，具體步驟參考試劑盒說明書。提取的DNA利用瓊脂糖凝膠電泳檢測其完整性，用Nanodrop 2000核酸測定儀檢測濃度。檢測合格的DNA樣品送北京諾禾致源生物信息科技有限公司進行后續檢測。

1.3 文庫構建與測序

對柚參考基因組(http://citrus.hzau.edu.cn/orange/download/index.php)進行電子酶切預測，最終選擇的限制性內切酶組合為MseⅠ+HaeⅢ+HaeⅡ。針對每個樣品，首先用限制性內切酶組合對基因組進行酶切，酶切后加上帶有條碼的接頭；接著對每個樣品進行PCR擴增，將擴增后的所有樣品混合；電泳后選擇大小為240～290 bp的片段，利用 Illumina HiSeq 4000測序平臺，進行雙末端各150 bp測序，用于文庫構建。

1.4 群體SNP檢測

測序獲得的原始數據，需要篩選以獲得高質量序列，篩選條件：去除原始測序序列中的接頭序列，去除單端測序中含N數量超過該序列長度10%的序列，去除單端測序中低質量(質量值Q≤5)堿基數超過該序列長度50%的序列。產生的高質量序列利用BWA(Burrows-Wheeler Aligner)程序與柚參考基因組比對(http://citrus.hzau.edu.cn/orange/download/index.php)。利用Samtools軟件進行SNP檢測，利用貝葉斯模型檢測群體中的多態性位點。然后過濾掉測序錯誤率>1%的SNPs，過濾掉兩者距離在5 bp范圍內的SNPs，選擇覆蓋深度在平均深度的[1/3,5]倍的SNP，最終獲得高質量的SNP。使用ANNOVAR軟件對過濾后的SNP檢測結果進行注釋。

1.5 系統進化樹構建及核苷酸多樣性(π)分析

基于每個SNP位點在19份樣品中的分型結果，采用鄰接法(neighbor-joining methods)構建系統進化樹。兩兩個體和之間的p-距離(Dij)通過如下公式計算：

式中，L為高質量SNPs區域長度，dij表示2個個體為不同基因型值的賦值表。

運用Treebest-1.9.2軟件計算距離矩陣，以此為基礎，通過鄰接法(neighbor-joining method)構建系統進化樹，使用FigTree軟件生成系統進化樹。引導值(bootstrap values)經過1 000次計算獲得。

核苷酸多樣性(π)是兩兩序列中具有不同核苷酸的比例。計算公式：

式中，pi是序列i的頻率，其中pj是序列j的頻率，πij是序列i與序列j的差異比例。使用軟件ngs Tools軟件(Fumagalli et al.,2014)進行計算。以200 kb長度為滑動窗口，100 kb為步長，進行核苷酸多樣性(π)分析。

2 結果與分析

2.1 樣本測序統計與質量評估

由表2看出，19個柚樣本的原始測序總堿基為5.93 Gbp，去除低質量的序列后，產生高質量序列堿基為5.61 Gbp。19個樣品測序堿基為255.11～388.38 Mbp，平均每個樣本堿基數為311.76 Mbp。19個樣品的Q20值為94.87%～96.28%，平均為95.76%；Q30值為87.44%～90.72%，平均為89.47%；GC含量為36.91%～37.62%，平均為37.3%，略高于柚參考基因組GC含量(34.99%)。

由表3可知，19個柚樣品的比對匹配率為97.56%～98.78%，平均為98.47%。測序深度范圍為9.39～12.54，平均為10.67。序列覆蓋度范圍為12.27%～13.84%，至少有1個堿基的平均覆蓋度為13.03%，至少有4個堿基的平均覆蓋度為7.24%。

表2 樣本GBS測序的數據統計

表3 與參考基因組比對的匹配程度

2.2 群體SNP的篩選及分布

利用 Samtools軟件檢測后共獲得608 466個SNP位點，過濾后共獲得79 825個高質量的SNP位點。由圖1可見，所有SNP位點在柚基因組9條染色體上基本均勻分布，分布密度為223.48～273.34個SNPs/Mbp。約5.52%(4 406個)的SNP位點未能有效定位到9條染色體上。對所有SNP位點在基因組上的位置統計表明(表4)，52.19%(41 674個)位于基因間區域，14.19%位于基因的上下游，20.19%位于內含子區域，13.38%位于外顯子區域，其中包括同義突變類型5.93%，非同義突變類型7.28%。

圖1 SNP位點在柚基因組染色體上的分布

Fig.1 Distribution of SNP sites in chromosomes ofC.maximagenome

表4柚基因組的SNP檢測及注釋結果

Table 4 SNP detection and annotation result ofC.maximagenome

注(Note)：外顯子(Exonic)：變異位于外顯子區域；內含子(Intronic)：變異位于內含子區域；剪接位點(Splicing)：變異位于剪接位點(內含子中靠近外顯子/內含子邊界的2 bp)；上游(Upstream)：基因上游1 kbp區域；下游(Downstream)：基因下游1 kbp區域；上游/下游(Upstream/Downstream)：基因上游1 kbp區域，同時也在另一基因的下游1 kbp區域；基因間區(Intergenic)：變異位于基因間區。

2.3 群體的系統進化樹

由圖2可見，19份樣品被聚類為3個主要的類群。A類群中大部分是采集于德江縣的柚子品種，另外還有少部分石阡縣、思南縣和印江縣的品種。這一類群又可分為3個亞類群，A1類群包括樣品1(石阡紅香柚)、樣品2(石阡紅橙)、樣品8(德江白橙4)、樣品9(德江白橙5)；A2類群包括樣品5(德江紅橙2)、樣品6(德江白橙2)、樣品7(德江白橙3)、樣品11(德江白橙6)、樣品16(思南紅橙)；A3類群包括樣品10(德江紅橙3)、樣品14(印江白橙)、樣品15(印江紅橙)。B類群包括樣品3(德江紅橙1)、樣品4(德江白橙1)、樣品13(印江淡紅橙)。C類群包括樣品12(印江紅香柚)、樣品17(利王白柚)、樣品18(利王紅柚)、樣品19(虎渡口紅橙)。

從進化樹可以看出，每個類群中的樣品基本分布在同一地理區域，一些變異材料與其親本也緊密聚類在一起，如樣品5(德江紅橙2)與6(德江白橙2)、3(德江紅橙1)與4(德江白橙1)。樣品6是在樣品5上發現的果肉呈現粉紅色的芽變枝條，樣品3是在樣品4上發現的果肉呈現粉紅色的芽變枝條，進一步驗證了結果的可靠性。

2.4 群體的核苷酸多樣性

從圖3可知，19份柚樣品的群體核苷酸多樣性值主要分布在4.84×10-6～4.0×10-4，最大值為1.5×10-3，平均值為1.09×10-4。

圖2 19份柚資源的系統進化樹

Fig.2 Phylogenetic tree of 19C.maximaresources

注：數值以200 kb為滑動窗口，100 kb為步長進行計算。

Note: These values were calculated in 200 kb sliding windows with 100 kb steps.

圖319份柚資源基因組的核苷酸多樣性

Fig.3 Nucleotide diversity of genome of 19C.maximaresources

3 結論與討論

柚多為單胚，且長期以實生繁殖，因而極易產生實生變異和天然雜種，目前對柚類種質資源的演化關系及品種間系統關系尚缺乏了解。陳鵬等[7]采用花粉形態學分析方法研究湖南省地方柚類分類和遺傳多樣性；劉勇等[8]則是結合SSR 引物和AFLP分子標記對110份柚類基因型、12份野生近緣種進行遺傳多樣性研究；陳巍等[9]利用SRAP標記對13份浙南柚類地方資源和琯溪蜜柚及芽變材料進行遺傳多樣性分析和鑒定。近年來，隨著高通量測序技術的發展，利用分子標記的方法對柚類進行遺傳多樣性研究變得越來越普遍，其中SSR和SNP 是廣泛使用的分子標記。而SNP 標記將成為未來基因型鑒定的主要標記類型[10]。基于第二代測序技術的GBS技術因其成本低、高效且不受參考基因組限制等優點，逐漸成為一種快捷有效的SNP標記開發手段[11]。

與常規SSR標記比較，GBS分析在獲得的標記數量上遠遠超過SSR標記[12-14]。并且，GBS分析獲得的標記在染色體上均勻分布，具有更高的分析精密度，這是SSR標記不可比擬的。同時，這種優勢也體現在試驗所需周期和成本上。因此，結合高通量測序和SNP標記優點的GBS技術，在未來必將廣泛用于資源多樣性評價相關研究。筆者利用GBS測序分型獲得了79 825個高質量SNP位點，這些位點廣泛分布于基因間區、內含子和外顯子中，在柚基因組9條染色體上均勻分布。群體核苷酸多樣性值主要分布在4.84×10-6～4.0×10-4，平均值為1.09×10-4，遺傳多樣性較低，可能與群體的樣本構成和選用的測序手段有關。利用這些SNP位點的分型結果構建系統進化樹，19份柚樣品被聚為3個主要類群，每個類群中的樣品基本分布在同一地理區域，一些變異材料與其親本也緊密聚類在一起，說明基于GBS測序分型的結果是準確可靠的。

研究中的柚群體核苷酸多樣性值主要分布于4.84×10-6～4.0×10-4，遠低于前人研究的柚群體19個樣品的核苷酸多樣性值(主要分布在1.0×10-3～2.0×10-3)[15]，說明19份柚樣品的遺傳多樣性較低。推測有以下2個方面的原因：首先是群體的樣品構成，試驗所用的19份樣品采自貴州武陵山地區6個區縣，盡管在植物學特征上存在差異，但樣品的采集范圍較窄；其次是檢測的基因組范圍不同，這與試驗使用的技術策略和目的不一致有關。筆者采用的GBS基因分型技術，屬于簡化基因組測序，目的在于獲得一些具有多態性的SNP位點，因此僅對酶切位點附近的SNP進行了檢測，最終僅覆蓋約13.03%的柚基因組，且試驗成本較低；WANG等[15]利用了全基因組重測序技術，目的在于評價整個基因組的核苷酸多樣性，其測序深度較高，基本覆蓋了整個基因組，檢測的位點更全面豐富，試驗成本也相對較高。

德江紅橙1(樣品3)和德江紅橙2(樣品5)均是實地調查中發現的果肉呈現粉紅或大紅色的變異材料，分別采自其親本德江白橙1(樣品4)和德江白橙2(樣品6)樹上出現的枝條變異株。在系統進化樹中，樣品3與樣品4、樣品5與樣品6分別聚在一起，其遺傳差異極小，這也證實了2份變異材料確實分別來自其親本上的枝條變異。這2對材料遺傳背景極為相似，僅在果肉色澤上表現差異，是研究果肉色澤變異機理和色素代謝調控的良好材料。對這2對材料的差異SNP進行了統計和注釋(結果未展示)，僅約6.8%的差異SNP位點位于編碼區。從注釋結果中也未發現這些位于編碼區的差異SNP與類胡蘿卜素代謝途徑基因相關，推測可能與采用的GBS技術所覆蓋的基因組范圍(僅覆蓋13.03%的基因組)和遺傳變異類型(僅檢測SNP)有關。因此，下一步將利用這2對材料，從全基因組水平探討色澤變異的遺傳機理，為解析果肉色澤變異和色素代謝調控提供研究基礎。