武鵬雨，李冬，吳翠云，包建平，金強(qiáng)，虎海防 ，張銳*

（1塔里木盆地生物資源保護(hù)利用兵團(tuán)重點(diǎn)實(shí)驗(yàn)室，新疆 阿拉爾 843300）

（2南疆特色果樹高效優(yōu)質(zhì)栽培與深加工技術(shù)國(guó)家地方聯(lián)合工程實(shí)驗(yàn)室，新疆 阿拉爾 843300）

(3塔里木大學(xué)園藝與林學(xué)學(xué)院，新疆 阿拉爾 843300)

（4新疆佳木果樹學(xué)國(guó)家長(zhǎng)期科研基地，新疆 溫宿 843100）

核桃（Juglans regia L.）隸屬于胡桃科（Juglandaceae）、核桃屬（Juglans），是世界4大堅(jiān)果（核桃、扁桃、板栗、腰果）之一，是不可或缺的堅(jiān)果和木本油料樹種［1］。我國(guó)核桃栽培歷史悠久，是核桃屬（Juglans）植物的起源和分布中心之一［2］，新疆作為我國(guó)核桃重要產(chǎn)區(qū)之一，核桃栽培面積達(dá)3.5×105hm2，其中塔里木盆地周圍綠洲產(chǎn)量最高，位居全國(guó)前列［3］。新疆栽培核桃資源豐富，擁有野核桃和紙皮核桃等種質(zhì)資源，也有薄皮核桃如‘扎343’‘新新2號(hào)’‘新萃豐’等品種，殼薄如紙，透光可見殼內(nèi)果仁形狀，薄皮核桃硬殼的平均厚度一般小于1.5 mm［4］，在薄皮不露仁栽培種中殼厚最厚的是‘新萃豐’，厚度為1.45 mm［5］；而‘新露’是天然內(nèi)果皮發(fā)育不全的露仁核桃品種，其內(nèi)果皮的殼面在生長(zhǎng)發(fā)育過程中發(fā)生自動(dòng)破裂使核桃仁露出，果仁裸露雖方便食用，但不耐運(yùn)輸儲(chǔ)藏。對(duì)擁有早實(shí)、薄殼、抗逆性強(qiáng)的‘新萃豐’和天然果殼裸露的‘新露’進(jìn)行研究，可以使新疆核桃內(nèi)果皮發(fā)育的遺傳背景清晰化，為新疆核桃的保護(hù)與科學(xué)利用提供理論依據(jù)。

隨著測(cè)序成本的降低和已知基因組序列物種的增多，基因組測(cè)序已成為研究植物分子育種、群體進(jìn)化中最為迅速有效的方法之一［6］?；蚪M測(cè)序技術(shù)為不同物種提供了高質(zhì)量的參考基因組序列，基因組重測(cè)序技術(shù)根據(jù)供參考的高質(zhì)量序列進(jìn)行整個(gè)基因組、區(qū)域或基因的重新測(cè)序，每個(gè)個(gè)體的測(cè)序數(shù)據(jù)與高質(zhì)量的參考數(shù)據(jù)進(jìn)行比對(duì)，不同樣本之間的遺傳變異可以被檢測(cè)為高度可信的序列差異，并在此基礎(chǔ)上對(duì)個(gè)體或群體進(jìn)行比較基因組學(xué)分析，可以從基因組水平對(duì)重要性狀的候選功能基因或QTL進(jìn)行定位與分析［7］，通過與參考基因組比較獲得多態(tài)性插入與缺失（InDel）標(biāo)記，使用這些分析標(biāo)記用于基因鑒定［8］，構(gòu)建雜交群體使用基因組重測(cè)序單核苷酸多態(tài)性位點(diǎn)（SNP）數(shù)據(jù)庫(kù)篩選與目的形狀相關(guān)基因［9］。本研究通過Illumina NovaSeq 6000測(cè)序平臺(tái)對(duì)新疆栽培核桃中的‘新露’和‘新萃豐’進(jìn)行全基因組重測(cè)序分析，對(duì)其SNP和InDel等變異位點(diǎn)進(jìn)行深度挖掘，從基因組層面分析新疆核桃內(nèi)果皮發(fā)育過程中出現(xiàn)的果仁裸露這一性狀的機(jī)理，通過檢測(cè)變異位點(diǎn)所在區(qū)域，為后續(xù)內(nèi)果皮發(fā)育等性狀的相關(guān)基因的挖掘、品種選育等研究提供一定的基礎(chǔ)。

1 材料與方法

1.1 試驗(yàn)材料

試驗(yàn)所選‘新露’與‘新萃豐’品種材料種植于新疆佳木果樹學(xué)國(guó)家長(zhǎng)期科研基地，2019年選取生長(zhǎng)情況一致的4年生核桃植株各4株，在5月份新出嫩葉時(shí)采摘葉片，放入液氮中保存帶回實(shí)驗(yàn)室。

1.2 試驗(yàn)方法

‘新露’與‘新萃豐’的DNA提取方法按照DNA提取試劑盒的使用說(shuō)明操作（TianGen）。DNA樣品的質(zhì)量控制通過0.8%瓊脂糖凝膠電泳檢測(cè)，同時(shí)采用紫外分光光度計(jì)對(duì)DNA進(jìn)行定量分析。

1.2.1 測(cè)序

對(duì)提取的DNA進(jìn)行超聲隨機(jī)打斷，隨后進(jìn)行末端損傷修復(fù)及連接接頭，在與測(cè)序接頭連接后進(jìn)行PCR擴(kuò)增，對(duì)測(cè)序文庫(kù)模板進(jìn)行富集并純化文庫(kù)產(chǎn)物，連接產(chǎn)物通過2%瓊脂糖凝膠電泳進(jìn)行最終片段選擇與純化，選擇插入片段大約為400 bp回收后根據(jù)制造商建庫(kù)流程構(gòu)建測(cè)序文庫(kù)。在Illumina Nova-Seq 6000測(cè)序平臺(tái)上對(duì)構(gòu)建的文庫(kù)進(jìn)行雙端（PE 150 bp）高通量測(cè)序。

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

獲得原始數(shù)據(jù)后采用FastQC將數(shù)據(jù)進(jìn)行質(zhì)量控制，進(jìn)而獲得高質(zhì)量數(shù)據(jù)。使用Burrows-Wheeler-Alignment Tool軟件將高質(zhì)量數(shù)據(jù)比對(duì)到已測(cè)序完成的核桃參考基因組，參考基因組數(shù)據(jù)為：JrSerr v1.0(GCA_004785585.1)［10］，利用 picard 軟件將比對(duì)得到的sam結(jié)果文件轉(zhuǎn)化為bam文件，使用GATK軟件對(duì)SNP和InDel進(jìn)行檢測(cè)，并用Unified Genotyper對(duì)變異位點(diǎn)位置進(jìn)行獲取，使用ANNOVAR軟件對(duì)獲取的SNP和InDel進(jìn)行功能注釋。獲取‘新露’與‘新萃豐’相同位置的SNP及InDel后，再根據(jù)非同義SNP（nonsynonymous，nsSNP）/InDel獲取相關(guān)基因。候選基因的富集分析遞交于AgriGO軟件用于富集Gene ontology terms。

2 結(jié)果與分析

2.1 測(cè)序數(shù)據(jù)預(yù)處理

以Illumina Nova Seq測(cè)序平臺(tái)提供的初始測(cè)序數(shù)據(jù)為原始數(shù)據(jù)（Raw Data），即本次測(cè)序中獲得的各樣本得到的短序列數(shù)（Reads）和堿基總數(shù)（Bases），對(duì)2個(gè)核桃樣本產(chǎn)出的數(shù)據(jù)進(jìn)行匯總（見表1）。2個(gè)品種共得到232 715 892個(gè)短序列，其中‘新露’核桃短序列數(shù)目127 572 756個(gè)，Q30比為91.30%；‘新萃豐’核桃中Reads數(shù)目105 143 136個(gè)，Q30比為91.40%。為剔除Illumina平臺(tái)錯(cuò)誤率對(duì)結(jié)果的影響，需對(duì)原始數(shù)據(jù)進(jìn)行質(zhì)量控制，包括去除低質(zhì)量序列，去除接頭，質(zhì)量過濾和長(zhǎng)度過濾等獲得純化數(shù)據(jù)（Clean Data），以進(jìn)行后續(xù)工作。質(zhì)量控制后‘新露’與‘新萃豐’分別獲得121 633 052個(gè)和100 378 460個(gè)短序列，測(cè)序總有效短序列與初始下機(jī)測(cè)序的比例分別為95.34%，95.46%。

表1 測(cè)序數(shù)據(jù)匯總

2.2 短序列匹配統(tǒng)計(jì)

將過濾后獲得的高質(zhì)量數(shù)據(jù)比對(duì)到核桃參考基因組，統(tǒng)計(jì)結(jié)果表明：‘新露’中總匹配的短序列數(shù)為121 332 322個(gè)，占所有短序列數(shù)的98.1%，重復(fù)序列占比15.15%，覆蓋全基因的深度為26.76 x，當(dāng)覆蓋深度≥20時(shí)，覆蓋全基因組的百分比為75.32%；‘新萃豐’中總匹配的短序列數(shù)為99 895 457個(gè)，占所有短序列數(shù)的97.96%，重復(fù)序列占比14.02%，覆蓋全基因的深度為22.45 x，覆蓋深度≥20時(shí)，覆蓋全基因組的百分比為59.92%。

2.3 SNP統(tǒng)計(jì)分析

‘新露’與‘新萃豐’核桃樣品檢測(cè)后獲得SNP位點(diǎn)統(tǒng)計(jì)信息（見表2），得到變異位點(diǎn)所在區(qū)域及變異影響。‘新露’樣本的總SNP數(shù)量為4 700 220個(gè)，雜合SNPs數(shù)量為2 645 454個(gè)，純合SNPs數(shù)量為2 054 766個(gè)，編碼區(qū)非同義突變數(shù)量（nsSNP）為81 580個(gè)，密度占比1.74%，SNP位于內(nèi)含子區(qū)域內(nèi)（Intronic）總數(shù)為535 307個(gè)，在基因間區(qū)域內(nèi)總數(shù)（Intergenic）為3 469 438個(gè)。

表2 ‘新露’與‘新萃豐’SNPs數(shù)目與類型檢測(cè)結(jié)果

‘新萃豐’樣本經(jīng)質(zhì)量過濾后獲得4 675 346個(gè)高質(zhì)量的SNP位點(diǎn)，編碼區(qū)內(nèi)非同義突變總數(shù)為83 811個(gè)，SNP位于內(nèi)含子區(qū)域內(nèi)總數(shù)為541 875個(gè)，位于基因間區(qū)域的SNP占比73.2%。

2.4 InDel統(tǒng)計(jì)分析

在基因組中檢測(cè)長(zhǎng)度小于50 bp的小片段插入或缺失（InDel）并進(jìn)行變異位置及類型的注釋，對(duì)‘新露’與‘新萃豐’的InDels變異位點(diǎn)信息進(jìn)行統(tǒng)計(jì)（見表3）?！侣丁目侷nDel數(shù)量為648 603個(gè)，編碼區(qū)內(nèi)移碼插入總數(shù)（Frame shift insertion）和移碼缺失（Frameshift deletion）數(shù)量分別為3 560個(gè)、2090個(gè)，InDel在基因5’UTR內(nèi)總數(shù)（5’UTR），3’UTR內(nèi)總數(shù)（3’UTR）和InDel位于不同基因的5’UTR和3’UTR內(nèi)總數(shù)（5’UTR/3’UTR）為3 095個(gè)、5 912個(gè)、2個(gè)，InDel位于內(nèi)含子區(qū)域內(nèi)總數(shù)為106 672個(gè)，在基因間區(qū)域內(nèi)總數(shù)為420 118個(gè)。

表3 ‘新露’和‘新萃豐’InDels數(shù)目與類型檢測(cè)結(jié)果

‘新萃豐’鑒定出606 188個(gè)InDel位點(diǎn)，其中插入297 690個(gè)，缺失311 822個(gè)，編碼區(qū)內(nèi)移碼插入總數(shù)和移碼缺失數(shù)量為3 449個(gè)、2 015個(gè)。InDels大部分位于基因間區(qū)域，有389 525個(gè)（64.2%），次之為內(nèi)含子 202 381個(gè)（16.89%），外顯子 8 048個(gè)（1.33%）。

2.5 DNA水平變異的基因分析

發(fā)生在外顯子區(qū)的變異可能會(huì)引起基因功能的異常改變，通過與參考基因組比對(duì)，檢測(cè)到‘新露’和‘新萃豐’基因組間發(fā)生非同義突變、同義突變和移碼突變等信息，分別獲得16 686個(gè)、19 829個(gè)變異SNPs以及1 687個(gè)、2 143個(gè)InDels，將這些變異位點(diǎn)比對(duì)到Gene ontology數(shù)據(jù)庫(kù)中。

2.5.1 ‘新露’核桃變異基因的GO富集分析

‘新露’核桃的SNP和InDel變異位點(diǎn)共注釋到4 865個(gè)基因，富集分析注釋結(jié)果分為生物過程、細(xì)胞組分和分子功能三類（如圖1）。注釋結(jié)果中，草酸代謝過程、囊泡運(yùn)輸與錳離子結(jié)合分別在生物過程、細(xì)胞組分和分子功能類別中被顯著富集。在這些被顯著富集的三類反應(yīng)中草酸代謝過程在植物體內(nèi)廣泛存在，作用于植物抗逆過程中；高爾基體的囊泡運(yùn)輸參與到植物細(xì)胞壁的形成，說(shuō)明露仁性狀可能和次生細(xì)胞壁形成有關(guān)；而錳離子則是植物體細(xì)胞中許多酶的催化劑。

圖1 ‘新露’變異基因的GO注釋分類圖

2.5.2 ‘新萃豐’核桃變異基因的GO富集分析

‘新萃豐’核桃共注釋到6 131個(gè)基因，GO terms富集分析如圖2所示，在生物過程和細(xì)胞組分中微管過程和囊泡運(yùn)輸被顯著富集，用于維持細(xì)胞骨架和信號(hào)傳導(dǎo)。在分子功能中富集到了天冬氨酸和苯丙氨酸的2種轉(zhuǎn)移酶，天冬氨酸轉(zhuǎn)氨酶，苯丙氨酸轉(zhuǎn)氨酶在植物體內(nèi)分別作用于苯丙烷途徑前后，最終目的為合成木質(zhì)素的單體。為進(jìn)一步對(duì)變異基因的生物學(xué)功能進(jìn)行了解，對(duì)富集到的天冬氨酸轉(zhuǎn)氨酶和苯丙氨酸氨基轉(zhuǎn)移后的裂解反應(yīng)的基因進(jìn)行挖掘，篩選出7個(gè)相關(guān)基因，結(jié)果見表4。

表4 木質(zhì)素合成相關(guān)變異基因

圖2 ‘新萃豐’變異基因的GO注釋分類圖

3 討論

隨著測(cè)序技術(shù)的快速發(fā)展，蘋果［11］、柑橘［12］、梨［13］等多種果樹的全基因組序列陸續(xù)公布，不同果樹物種的基因組信息逐漸被完整破譯。通過全基因組重測(cè)序可以找到大量的SNP、拷貝數(shù)變異（CNV）、InDel、結(jié)構(gòu)變異（SV）等變異信息，基于檢測(cè)到的變異位點(diǎn)能進(jìn)一步研究該物種的特性［14］、群體進(jìn)化［15］和定位目標(biāo)性狀基因位點(diǎn)［16］。

SNP在基因組中廣泛存在，并含有豐富的遺傳信息，隨著測(cè)序技術(shù)的發(fā)展以及生物信息學(xué)分析能力的提升，SNP作為第三代分析標(biāo)記應(yīng)用于相關(guān)性狀的基因定位研究中［17］。SNP的最早研究用于人類基因組，最終繪制成142萬(wàn)個(gè)SNP的人類基因組圖譜，這種高密度的SNP圖譜為整個(gè)基因組單倍型變異提供了依據(jù)［18］。將2個(gè)橄欖品種雜交鑒定獲得10 941個(gè)SNP，并使用其中多態(tài)性的SNP位點(diǎn)構(gòu)建覆蓋基因組3 049 cM，25個(gè)連鎖群體的高密度遺傳連鎖圖譜，遺傳連鎖圖譜的構(gòu)建可更有效地定位數(shù)量性狀和其他重要經(jīng)濟(jì)性狀［19］。使用基因測(cè)序方法進(jìn)行基因分型，對(duì)180種棗進(jìn)行全基因組關(guān)聯(lián)分析（GWAS），鑒定出4 651個(gè)高質(zhì)量SNP和45個(gè)與SNP性狀關(guān)聯(lián)基因，為標(biāo)記輔助育種提供了數(shù)據(jù)基礎(chǔ)［20］。將SNP作為新型分子標(biāo)記也應(yīng)用于核桃分子育種中，將14 761個(gè)SNP進(jìn)行全基因組關(guān)聯(lián)，確定了葉芽萌發(fā)和花芽分化等性狀的標(biāo)記位點(diǎn)，為后續(xù)QTL鑒定奠定了基礎(chǔ)［21］。在自然適應(yīng)和選擇性育種過程中，最常見且研究最多的SNP外還有InDel，InDel變異是不同個(gè)體間基因組同一位點(diǎn)的序列發(fā)生不等大小片段的插入與缺失，從而使同源序列比對(duì)產(chǎn)生空位（gap）的現(xiàn)象［22］。有研究者將SNP、InDel、SV以及簡(jiǎn)單序列重復(fù)(SSR)組裝成核桃基因組變異圖譜，利用遺傳信息豐富的標(biāo)記位點(diǎn)對(duì)中國(guó)五個(gè)核桃品種進(jìn)行系統(tǒng)發(fā)育研究，將5個(gè)核桃種以100%的bootstrap(BS)值劃分為兩個(gè)已知區(qū)段核桃/胡桃組(Juglans/Dioscaryon 和核桃楸組 Cardiocaryon)［23］。在核桃上使用454焦磷酸測(cè)序技術(shù)開發(fā)鑒定獲得48 165個(gè)SNPs和1 037個(gè)InDels，轉(zhuǎn)換/顛換（ts/tv）比例為2.79:1，略高于本次新疆核桃測(cè)序的比例［24］。本研究選用具有代表性的新疆露仁品種與薄皮中的厚殼品種的‘新露’與‘新萃豐’，利用基因組重測(cè)序技術(shù)檢測(cè)核桃內(nèi)果皮果殼發(fā)育不全，果仁裸露這一性狀的各種遺傳變異?！侣丁c‘新萃豐’的總SNP數(shù)目為4 700 220個(gè)，4 675 346個(gè)，總InDel數(shù)目為648 603個(gè)，606 188個(gè)，變異位點(diǎn)在基因編碼區(qū)，基因間區(qū)域和內(nèi)含子區(qū)域都有分布，在全基因組水平上有著顯著的遺傳變異，2者的SNP數(shù)量顯著多于InDel的變異數(shù)量，說(shuō)明外界環(huán)境對(duì)基因組的變化主要以單核苷酸多態(tài)性變異為主，這與梨［25］、馬鈴薯［26］等測(cè)序結(jié)果一致。本研究將檢測(cè)到的‘新露’與‘新萃豐’核桃基因組間所發(fā)生的非同義突變、同義突變等基因比對(duì)到GO數(shù)據(jù)庫(kù)中，發(fā)現(xiàn)在薄皮不露仁系列中的厚殼性狀品種‘新萃豐’富集到與天冬氨酸轉(zhuǎn)氨酶和苯丙氨酸轉(zhuǎn)氨酶的相關(guān)過程，通過與基因組序列比對(duì)分析，篩選到7個(gè)與木質(zhì)素生物合成有關(guān)的基因。本試驗(yàn)通過全基因組重測(cè)序技術(shù)對(duì)新疆2個(gè)核桃栽培種進(jìn)行變異檢測(cè)分析，通過比對(duì)變異位點(diǎn)，挖掘變異的關(guān)鍵基因，為基因定位和克隆功能驗(yàn)證以及為核桃分子育種的利用奠定了基礎(chǔ)。