基于基因分型技術的燕麥SNP標記研究

董艷輝，劉龍龍，溫 鑫，于宇鳳，楊 方，劉根科，崔 林，曹秋芬，秦永軍,

(1.山西省農業科學院 生物技術研究中心，山西 太原 030031；2.山西省農業科學院 農作物品種資源研究所，農業部黃土高原種質資源實驗室，山西 太原 030031；3.山西省農業科學院 農業科技信息研究所，山西 太原 030031；4.山西省農業科學院 右玉農業試驗站，山西 右玉 037200)

燕麥(AvenasativaL.)是世界各地廣泛栽培種植的一種重要的糧食兼飼草、飼料作物，其富含β-葡聚糖，被證明具有降血糖和血脂的功效[1]，是美國FDA和英國JHCI認定的功能性谷物[2-3]，隨著燕麥及其加工產品的世界消費量的增長，傳統燕麥育種技術已經不能有效滿足燕麥生產的需要，并且成為制約燕麥育種進程的重要因素，而分子輔助育種是加速燕麥育種的重要手段。

目前，燕麥的育種多依賴于常規育種和表型鑒定研究，已成功育成了多種適應不同需求的燕麥品種并用于大田生產[4]。早期的表型特征量化研究證明燕麥具有豐富的遺傳多樣性[5]，之后燕麥種質資源遺傳連鎖圖[6]、皮裸性基因相關聯SSR[7]和AFLP標記[8]也被相繼報道。單核苷酸多態性(SNP)是目前最具優勢的分子標記，被廣泛應用于作物遺傳圖譜構建、QTL定位、全基因組關聯分析(GWAS)、群體進化等[9-12]。近年來，基于下一代測序技術(Next generation sequencing,NGS)的基因分型技術(Genotyping by sequencing,GBS)是通過酶切處理去除基因組重復序列，顯著降低了測序量和測序成本，且能夠覆蓋整個基因組，不受參考基因組限制，因而被廣泛應用到開發高密度、高精度的農作物SNP分子標記研究中。目前，燕麥全基因組序列還沒有發表，基因分型技術對農作物遺傳多樣性研究、種質資源鑒定和分子輔助育種標記開發具有重要的意義[13-16]。基因分型技術從單酶切技術RAD[17]、2b-RAD[18]、GBS[19]到現在的雙酶切技術dd-GBS[20]、SLAF[21]均已被應用到玉米、小麥等[22-30]主要農作物分子標記開發研究。Winkler 等[31]采用GBS技術，證明了燕麥穎片顏色和穗型等表型特征和SNP位點相關。在另一項研究中，通過GBS技術構建的高密度SNP則能精準地區分起源于世界不同地區燕麥種質資源[32]，證實GBS-SNP技術是種質資源鑒定的有力工具。

本研究基于NGS的dd-GBS技術，以經過抗旱性鑒定的燕麥種質為試驗材料，選用適合燕麥基因組酶切的酶，運用適合燕麥簡化基因組參考序列的分析軟件，在構建燕麥簡化基因組參考序列的基礎上，研究燕麥相關SNP標記，旨在豐富當前燕麥的基因組數據庫，提供分子標記(SNP)基礎數據，為加速燕麥育種進程提供理論支撐。

1 材料和方法

1.1 試驗材料

供試材料由中國農業科學院品種資源研究所麥類室于1981-1983年進行抗旱鑒定的燕麥種質42份，包括11份高抗(High tolerance)、13份中抗(Moderate tolerance)、18份低抗(Low tolerance)和山西省農業科學院農作物品種資源研究所燕麥種質資源圑提供[33]的抗旱性未鑒定(Uncharacterized)燕麥育成品種(品系)12個(表1)，2016年種植于山西省農業科學院智能溫室。

1.2 試驗方法

GBS文庫構建與數據質控：燕麥種子發芽14 d，取鮮嫩葉片提取基因組DNA(Mag-MK Plant Genomic DNA Extraction Kit試劑盒，上海生工)，1.0%瓊脂糖凝膠電泳和紫外分光光度計分別檢測DNA質量。用限制性內切酶PstⅠ(CTGCAG)和MspⅠ(CCGG)對提取的DNA 進行酶切后，回收大小在220～450 bp的酶切片段，之后按照dd-GBS方法進行建庫[20]。采用llumina Hiseq測序平臺進行雙末端(Paired-end，PE)測序，獲得的下機數據以雙端FASTQ格式保存。

分別對每個樣品的下機數據運用 FastQC軟件進行質量控制，主要包括堿基質量分布(Per base sequence quality)、質量值(Per sequence quality scores)、GC分布(Per sequence GC content)，對反映下機數據質量的多項指標進行統計。進一步過濾數據包含一些帶接頭，雙末端reads1 5′端6 bp非酶切位點序列CTGCAG或reads2的5′端4 bp非酶切位點序列CCGG。采用Adapter Removal去除3′端的接頭污染[34]。采用滑動窗口進行質量過濾，窗口大小設置為5 bp，步長設置為1 bp，每一次往前移動1個堿基，取5個堿基計算窗口的平均Q值，若最后一個堿基的Q值≤2，則僅保留該位置之前的堿基；若窗口的平均Q值≤20，則僅保留該窗口倒數第2個堿基及之前的堿基。濾除小于50 bp 的reads。

表1 不同抗旱性的燕麥試驗材料Tab.1 Sequenced oats germplasm of different drought tolerance

1.3 數據分析

1.3.1 燕麥簡化基因組參考序列組裝與群體SNP檢測 應用Stacks軟件包[35]中的ustacks對每個樣品的reads進行聚類，同一個stack代表1個酶切位點(loci)，聚類參數-m設置為4，對每一個樣品的loci及loci的測序深度進行統計。下一步用cstacks將所有樣品的loci合并，不同樣品loci之間最多允許2個錯配，獲得每個loci的catalog consensus序列。

采用sstacks將每個樣品的loci序列與catalog consensus序列比對后populations過濾獲得群體SNP。主要參數包括：1個位點最少要在1個群體中出現；1個群體中檢測到同一位點的個體最低百分數50%(當群體中該位點的缺失率超過 50%，則去除該位點)；1個位點的最小等位基因頻率0.05。

1.3.2 燕麥SNP分析 采用 GCTA軟件利用群體SNP數據進行主成分聚類分析；對高、低抗材料的SNP數據進行Fisher Test統計分析獲取相關SNP標記；采用結合校正標準錯誤發現率(False Discovery Rate，FDR)小于0.05和0.001確定差異顯著性標記和極顯著性標記。

1.3.3 燕麥SNP功能預測 進一步對差異極顯著(FDR<0.001)的SNP所在源序列和小麥(TriticumaestivumL.)基因組序列進行Blast比對，篩查SNP相關基因。

2 結果與分析

2.1 燕麥GBS測序結果

本研究共測序54個樣本，注釋reads數在4 111 218～21 782 382；質控后保留的高質量數據(reads數4 111 218～19 019 296)用于后續的數據分析(圖1)。

圖1 燕麥GBS-SNP注釋序列Fig.1 Oat GBS-SNP annotation sequence

2.2 燕麥參考序列與群體SNP位點

采用Stacks數據包對燕麥reads進行聚類，結果顯示：每個loci平均測序深度為20.61×，合并所有樣本同一個loci，共獲得753 325個參考基因組序列(catalog consensus)。每個樣品的loci與其對應的catalog比對，去除群體位點缺失(./.)率超過50%和最小等位基因頻率小于0.05的SNP位點，共獲得74 657個群體SNP位點。由圖2可知，所有測序樣本中純合狀態的SNP位點(0/0或1/1)占少數，大部分呈現不同程度的雜合(0/1)或缺失狀態(./.)。

0/0 和1/1.純合體；0/1.雜合體；./..缺失。0/0 和1/1. Homozygote; 0/1.Heterozygote;./.. Deletion.

2.3 燕麥SNP研究

根據群體SNP數據，對54份燕麥種質的PCA結果顯示，SNP基因型明顯聚類為2簇，全部11份高抗、4份抵抗、8份中抗種質和6個未鑒定品種聚為一簇；而另外一簇則包含14份低抗和5份中抗種質；其余6個未鑒定的燕麥品種則分布在2簇之外(圖3)。

基于PCA基因型聚類結果，進一步對高抗(n=11)和低抗(n=18)種質2組進行SNP基因型差異顯著分析，共獲得2 937個燕麥抗旱性相關SNP(Fisher Test，FDR<0.05)。其中，差異極顯著(FDR<0.001)的55個SNP位點坐落于41個燕麥簡化基因組參考序列(同一參考序列包含2或3個SNP位點)。這55個相關SNP可分為3類，23個SNP 位點在高抗種質中缺失； 28個SNP位點在低抗種質中缺失；其余4個SNP位點在2種種質中都存在(表2)。

表2 高抗與低抗燕麥的差異顯著SNP位點Tab.2 Oat significant SNPs of high tolerance and low tolerance

表2(續)

注：▲.SNP位點缺失。

Note: ▲. Represents site deletion of SNP.

2.4 燕麥SNP的Blast功能預測

進一步對差異極顯著的55個相關SNP的源序列與小麥基因組序列進行了Blast比對，結果顯示，10個源序列(包括14個SNP位點)和多個基因轉錄序列聯配(表3)。其中，參與線粒體活性氧信號轉導、植物激素信號轉導的調控蛋白、自噬基因表達與植物抗逆、抗旱水平有關。

表3 燕麥SNP源序列聯配基因(比對小麥基因組)Tab.3 Oat SNPs alignment genes (against wheat genome)

3 結論與討論

3.1 GBS技術開發高密度SNP分子標記

高通量測序技術的面世是現代分子生物學技術的革命性飛躍，NGS技術促使生命科學研究進入基因組學時代，NGS在農作物上的重要應用之一就是GBS技術。GBS技術已經成為加速農作物種質資源鑒定篩選以及分子輔助育種的重要技術之一[16,19,46]，燕麥基因組相對復雜(2n=6x=42)，具有豐富的遺傳多樣性。本研究采用新一代ddGBS技術，成功組裝了基于54個樣本的燕麥簡化基因組參考序列。質控顯示，獲得的高質量reads通量和SNP雜合狀態與抗旱性沒有明顯關聯。GBS通過酶切過濾重復序列，降低基因組復雜程度，雖不代表全基因組序列結構，但非常適合于農作物分子輔助育種SNP標記開發。本研究之所以選擇經過多年抗旱性鑒定、遺傳穩定、抗旱梯度差異明顯(高抗、中抗、低抗)的燕麥種質資源進行測序，目的是通過抗旱性表型與SNP基因型關聯分析，開發可能與抗旱性有關聯的SNP。因本研究樣本數量相對偏小(n=54)，在構建簡化基因組參考序列中沒有采用報道的針對大樣本TASSEL分析流程[33]，而應用更適合于小群體的Stacks分析流程[35]。雖技術路線不同，但每項研究貢獻的簡化基因組參考序列無疑豐富了目前的燕麥基因組數據庫。

3.2 抗旱性和SNP基因型關聯分析

本研究的主要目的是挖掘可能與燕麥抗旱性有關聯的SNP，為燕麥種質資源鑒定和分子輔助育種服務。這些試驗材料根據多年的抗旱性鑒定劃分為高抗、中抗和低抗3組，但根據SNP基因型聚類分析則大致分為2簇，高抗材料SNP基因型和大多數低抗材料基因型各聚類為一簇，明確顯示抗旱性和SNP基因型是關聯的，也就是說SNP基因型很大程度上能反映抗旱性。中抗材料SNP基因型沒有聚類成簇，而是呈分散狀態，也從側面說明表型是基因型與環境共同作用的結果，同時作物抗旱性又是一個多指標的綜合性狀表現，但基因型的選擇是改良作物性狀的有效途徑[47-48]。燕麥抗旱表型是抗旱基因決定的性狀在環境作用下的具體表現，其抗旱性鑒定結果可能在不同的環境條件下會有所偏差。

因高、低抗材料的抗旱性和SNP基因型相關聯，本研究針對這2組材料，開發可能與燕麥抗旱性相關的SNP，大部分SNP屬于雜合或缺失狀態，只有極少數呈現純合狀態，可能是因為燕麥基因組高度雜合(六倍體)以及測序技術的原因。下一步仍需要應用精度更高的數字化 PCR 技術在更大群體驗證本研究發現的SNP。

利用GBS技術開發可能與燕麥抗旱性相關SNP標記在其他作物育種研究中已經報道[49-53]。對于燕麥來講，本研究采用GBS技術，第1個報道可能與燕麥抗旱性有關聯的SNP標記，以期為燕麥種質資源精準篩查和分子輔助育種提供數據基礎。

3.3 相關SNP生物學功能

基因組技術已被廣泛應用到抗旱基因挖掘[54-56]。GRLs[57]、OsMYB55[58]和GIF[59]被相繼報道與不同作物抗旱機制有關。本研究經Blast比對獲得的相關SNP源序列和不同生物通路基因轉錄子聯配。其中cat 186303的第83個堿基在高抗材料中位點缺失，而在低抗材料中則呈現GT雜合狀態。cat 186303源序列Blast 比對小麥基因組數據庫顯示，其序列與小麥硫氧還蛋白Trx序列相似。Trx蛋白參與細胞生長凋亡和基因轉錄[60]，多項研究證明，Trx基因和作物抗旱、抗逆功能有關[16, 61-65]。后續功能驗證包括研究 cat 186303對燕麥Trx基因表達調控，是解釋功能SNP的重要步驟。此外，其他抗旱性相關SNP源序列聯配的基因包括植物激素信號轉導的DUF221、雙特異蛋白磷酸酶和核孔蛋白，與植物抗逆有關。

3.4 本研究的局限性

燕麥是六倍體作物，本研究所開發的燕麥抗旱相關聯SNP不能定位，但證明不同抗性材料的燕麥在SNP存在表達差異，是作為燕麥材料分子篩選的潛在標記。因為抗旱性狀是一個綜合指數，單一參數不能全面代表抗旱能力。本試驗得出的特異性SNP 需要今后在大量的群體材料中進行驗證。此外，燕麥抗旱指標有待量化也是今后挖掘抗旱SNP 的重要基礎。本研究結果支持表型-基因型互作關系的同時也指出，抗旱性是基因型和環境因子綜合體現，初期開發的可能與燕麥抗旱性有關聯的SNP必須在更大種質群體中驗證以加強其SNP的廣譜性和實用性。

常規育種主要依賴于表型數據，而分子輔助育種則側重于基因型數據，構建燕麥抗旱性表型-基因型關聯則是精準育種的必經之路。應用GBS技術，本研究成功組裝了燕麥簡化基因組參考序列并初步研究了可能與燕麥抗旱性有關聯的SNP標記，研究結果不僅豐富了燕麥基因組數據庫，同時為燕麥種質資源早期篩選和燕麥基因數據庫提供數據基礎。

致謝：感謝荷蘭萊頓大學統計系Goeman J J教授的統計學方法指導；感謝美國加州戴維斯分校植物系羅明誠教授對本試驗給予的建議。