張 莉，張 勇，張淵海，張逸妍，解 濤，彭選明，楊 震*

（1. 湖南省核農(nóng)學(xué)與航天育種研究所，長沙 410125；2. 湖南省農(nóng)業(yè)生物輻照工程技術(shù)研究中心，長沙 410125）

水稻作為重要的糧食作物，它的豐產(chǎn)、穩(wěn)產(chǎn)對國民經(jīng)濟發(fā)展具有重要的戰(zhàn)略意義。誘變育種技術(shù)在農(nóng)作物新材料創(chuàng)制和優(yōu)良新品種的培育中已經(jīng)發(fā)揮了顯著作用［1-3］。常規(guī)水稻品種的誘變策略通常是改進一個或者一些容易鑒別的性狀，并利用突變體植株掃描突變位點和定位相關(guān)基因來研究突變的分子機制［4-5］。植物輻射誘變通過改變?nèi)旧w數(shù)目、染色體結(jié)構(gòu)或DNA的核苷酸構(gòu)成從而引起表型變異［6］。在重離子輻射誘變DRF缺失的煙草突變體中，NtDRF2整個基因缺失，NtDRF1發(fā)生1個堿基缺失使肽鏈編碼發(fā)生移框突變，最終形成白花煙草［7］。Wang等［8］通過12C6+輻射月牙藻獲得了5個葉綠素a缺陷的突變體，研究發(fā)現(xiàn)突變體中捕光色素復(fù)合體相關(guān)的Lhcb5、Lhcbm5和Lhcbm1基因的表達發(fā)生顯著變化。目前基因組水平上誘導(dǎo)突變的分子機制仍不清楚，隨著生物信息技術(shù)的高速發(fā)展，通過全基因組測序（whole-genome sequencing，WGS）技術(shù)可以大規(guī)模檢測許多植物的個體突變［9］，WGS的高效率和低成本使其成為發(fā)掘基因組突變的一種好方法。

遺傳變異包括序列變異和結(jié)構(gòu)變異。序列變異包括單核苷酸多態(tài)性（single nucleotide polymorphisms，SNPs）、插入缺失（insertions and deletions，In-Dels）、微衛(wèi)星或簡單序列重復(fù)（microsatellites or simple sequence repeats，SSRs）和轉(zhuǎn)座因子。這些序列多態(tài)性由于其低成本、穩(wěn)定性和高通量的應(yīng)用已被廣泛應(yīng)用于基因組選擇、數(shù)量性狀位點（quantitative trait locus，QTL）定位、單倍型和家系分析等領(lǐng)域［10］。如Zhang等［11］通過鑒定早熟三葉柑桔與其野生型的全基因組遺傳變異，開發(fā)出新的遺傳標(biāo)記進行柑桔重要性狀研究。結(jié)構(gòu)變異（structural variations，SVs）包括大片斷插入（large-scale insertions，INSs）、缺失（deletions，DELs）和倒置（inversions，INVs）以及基因組水平的染色體內(nèi)部和染色體間易位，這些類型的結(jié)構(gòu)變異統(tǒng)稱為拷貝數(shù)變異（copy number variants，CNVs），這些變異在形成基因組多樣性方面起著重要作用。Mu?oz-Amatriaín等［12］采用8個大麥栽培品種和6個野生大麥品種進行基因組雜交比較。與栽培大麥相比，野生大麥中存在更高水平的拷貝數(shù)變異多樣性?；蚩截悢?shù)占基因陣列編碼序列的9.5%，被拷貝數(shù)影響的基因標(biāo)記為抗病蛋白和蛋白激酶。栽培大麥Barke和Morex品種的CNVs序列比較表明，單鏈退火和合成依賴鏈退火的雙鏈斷裂DNA修復(fù)機制在大麥CNVs的發(fā)生過程中起著重要作用。Wallace等［13］對玉米5 000個近交系中的41種不同表型進行全基因組關(guān)聯(lián)分析，檢測到28 900 000個SNPs和800 000個CNVs，發(fā)現(xiàn)基因區(qū)和基因間區(qū)具有相反的富集模式、較小的等位基因頻率和效應(yīng)大小，同時全基因組關(guān)聯(lián)分析（genomewide association study，GWAS）標(biāo)記的基因具有豐富的調(diào)控功能，這表明基因調(diào)控和基因復(fù)制是表型變異的強大驅(qū)動因素。

目前對作物序列多態(tài)性的研究越來越多，不同作物中已經(jīng)開發(fā)了幾種模型的全基因組SNPs和InDels數(shù)據(jù)庫［14-15］，并廣泛應(yīng)用于遺傳多樣性研究［16］、功能和進化研究［17-18］以及馴化和基因組進化研究［9，19-20］。本研究利用60Co-γ輻照水稻品種黃華占獲得遺傳穩(wěn)定突變體湘輻1821，該突變體較黃華占株高增高，葉寬增寬，單產(chǎn)顯著提高。同時，利用高通量深度測序來研究其基因組變異，以期為輻射育種技術(shù)提供理論支撐。

1 材料與方法

1.1 試驗材料

水稻突變體湘輻1821（Xiangfu 1821，MT_HHZ，Oryza satiuaL.）系350 Gy60Co-γ射線輻照親本黃華占（Huang Hua Zhan，WT_HHZ，Oryza satiuaL.）干種子，經(jīng)多代大田選育而成。水稻品種黃華占購買于湖南農(nóng)豐種業(yè)有限公司。

1.2 誘變處理

2016年在湖南省核農(nóng)學(xué)與航天育種研究所輻照中心利用60Co-γ射線對黃華占進行誘變處理，輻照劑量為350 Gy，劑量率為5.25 Gy/min。

1.3 種植方法

M1代編號B16314當(dāng)季大田多本粗插種植，混收，只收獲其主穗種子，以期獲得更大的變異率。2017年隨機選取B16314種子，單本種植，群體量20 000株左右，在其田間選擇出綜合性狀良好的突變體，單株收獲，編號B1716。2017年冬在海南三亞繁殖B1716種子，并觀察其遺傳穩(wěn)定性，收獲其種子，編號D71157。2018年繼續(xù)種植D71157觀察其遺傳穩(wěn)定性，并進行聯(lián)合品比試驗，收獲其種子，編號B81247。2019年進行聯(lián)合品比試驗，并進行米質(zhì)分析，種子編號B1821，命名為湘輻1821。

1.4 測序幼苗樣品取樣

將黃華占和湘輻1821種子常規(guī)浸種3 d，發(fā)芽后的種子播種至帶土塑料桶于人工氣候箱中培養(yǎng)。至3葉期對地上部分取樣，將樣品送至北京諾禾致源生物信息科技有限公司進行后續(xù)工作。

1.5 DNA文庫構(gòu)建、測序和數(shù)據(jù)處理及分析

十六烷基三甲基溴化銨（cetyltriethylammnonium，CTAB）法提取DNA，質(zhì)檢合格的DNA樣品通過Covaris破碎機隨機打斷成長度為350 bp的片段。采用TruSeq Library Construction Kit 試劑盒進行建庫。DNA片段經(jīng)末端修復(fù)、加ployA尾、加測序接頭、純化、PCR擴增等步驟完成整個文庫制備。構(gòu)建好的文庫通過illumina進行測序。首先，對測序得到的原始序列進行過濾，去除帶接頭的序列對。其次，檢查測序錯誤率分布，每個堿基位置的測序錯誤率都應(yīng)低于1%。最終得到的過濾序列為有效測序數(shù)據(jù)。通過BWA［21］軟件比對到參考基因組，比對結(jié)果經(jīng) SAMTOOLS［22］軟件去除重復(fù) ；采用 SAMTOOLS軟件進行個體 SNPs 的檢測［SNPs的 reads 支持數(shù)不低于4；SNPs 的質(zhì)量值（mapping quality，MQ）不低于20］和InDels的檢測；通過CNVnator［23］檢測CNVs；利用 BreakDancer［24］軟件檢測SVs。基因本體（gene ontology，GO）分析首先把所有候選基因向Gene Ontology數(shù)據(jù)庫（http://www.geneontology.org/）的各個分類條目映射，計算映射到每個term的靶基因數(shù)目，然后應(yīng)用超幾何檢驗，找出與整個參考基因背景相比，在候選靶基因中顯著富集的GO條目，通過GO功能顯著性富集分析能確定候選靶基因行使的主要生物學(xué)功能。

2 結(jié)果與分析

2.1 測序數(shù)據(jù)質(zhì)量情況匯總

通過測序數(shù)據(jù)分析，本次測序共產(chǎn)生原始數(shù)據(jù)82.89 Gb，過濾后的有效數(shù)據(jù) 82.73 Gb，各樣品有效數(shù)據(jù)大小在39 884.828～43 004.680 Mb之間。GC含量在43.25%～43.63%之間，所有樣品的數(shù)據(jù)量足夠，測序質(zhì)量合格，GC分布正常，符合要求，可以進行后續(xù)分析。本文選用秈稻蜀恢498（R498）作為參考基因組，其大小為390 983 850 bp，所有樣本的比對率在96.25%～97.18%之間， 對參考基因組的平均覆蓋深度在87.79倍（X）～88.28倍（X）之間， 1倍（1X）覆蓋度（至少有1個堿基的覆蓋）在95.46% 以上。比對結(jié)果正常，可用于后續(xù)的變異檢測及相關(guān)分析。表1是測序數(shù)據(jù)和比對結(jié)果小結(jié)。

表1 黃華占（WT_HHZ）和湘輻1821（MT_HHZ）測序數(shù)據(jù)小結(jié)Tab. 1 Summary of illunima data in Huanghuazhan（WT_HHZ）and Xiangfu 1821（MT_HHZ）

2.2 SNPs變異分析

采用SAMTOOLS軟件進行個體SNPs的檢測。黃華占（WT_HHZ）和湘輻1821（MT_HHZ）分別得到758 215和799 434個SNPs，基因間隔區(qū)分布的SNPs最多，多達481 765和508 409個。其次為內(nèi)含子區(qū)域有69 806和72 355個，基因上游1 kb區(qū)域有64 274和68 140個。外顯子則只有47 603和50 561個，以非同義突變?yōu)橹饕愋?。湘?821的SNPs總數(shù)較黃華占多了41 219個，且不同分布區(qū)域的SNPs數(shù)值都高于黃華占。湘輻1821的SNPs的雜合率為0.236%，比黃華占高了0.021 %。根據(jù)堿基替換不同，SNPs分為轉(zhuǎn)換（transitions，Ts）和顛換（transversions，Tv）。湘輻1821轉(zhuǎn)換與顛換的比率為2.546，黃華占的轉(zhuǎn)換與顛換的比率為2.542，二者基本持平（表2）。這2種基因型的SNPs變異基本以轉(zhuǎn)換為主要類型。轉(zhuǎn)換即是G/A 和 C/T的2種變異，顛換則有A/C、C/G、G/T和T/A 4種變異類型。在黃華占中，轉(zhuǎn)換變異達到544 204個，顛換只有214 011個，湘輻1821中類似，轉(zhuǎn)換達到574 019個，顛換則只有225 415個，二者都以轉(zhuǎn)換為主要類型（圖1）。

表2 SNPs分布區(qū)間及數(shù)目Tab. 2 Summary of SNPs discovery via illunima sequencing

圖1 SNPs突變頻譜圖Fig. 1 The mutation spectrum of SNPs

2.3 InDels變異分析

利用SAMTOOLs軟件檢測長度小于50 bp的小片段的插入缺失（InDels），黃華占（WT_HHZ）/湘輻1821（MT_HHZ）一共得到142 313/147 686個InDels。基因間隔區(qū)分布的InDels最多，多達77 888/80 604個，其次為內(nèi)含子區(qū)域有19 174/19 977個，基因上游1 kb區(qū)域有15 632/16 319個。湘輻1821外顯子中的變異以移碼變異為主，共有2 468個，非移碼的有2 357個（表3）?；蚪MInDels長度從1 bp到21 bp不等，2種基因型都是1 bp的InDels為主要變異類型，占50%（圖2）。編碼區(qū)的InDels也是1～21 bp，1 bp的InDels在2種基因型中約占27%，2 bp的InDels約占11%，3 bp的InDels約占24%（圖3）。在黃華占和湘輻1821中基本以短序列InDels為主，說明點突變?yōu)橹饕蛔兎绞健?/p>

表3 InDels分布區(qū)間及數(shù)目統(tǒng)計Tab. 3 Summary of InDels discovery via illunima sequencing

2.4 SVs分析

SVs指基因組水平上INSs、DELs、INVs、染色體內(nèi)部遷移（intra-chromosomal translocations，ITXs）、染色體間的遷移（inter-chromosomal translocations，CTXs）。黃華占中一共得到16 775個SVs，包括29個 INSs、5 751個 DELs、1 275個 INVs、2 354個ITXs 和7 366個CTXs。大部分 SVs位于外顯子區(qū)域（2 542）和基因間隔區(qū)（2 917），基因上游1 kb區(qū)域有 527個， 下游1 kb區(qū)域有447個，內(nèi)含子區(qū)域有334個， 8個 SVs影響剪接位點。在湘輻1821中一共得到18 382個SVs，比黃華占多了1 607個，包括29個 INSs、6 251個 DELs、1 318個 INVs、2 596 個 ITXs 和 8 188 個 CTXs。多數(shù) SVs 位于外顯子區(qū)域（2 648）和基因間隔區(qū)（3 211），基因上游1 kb區(qū)域有598個， 下游1 kb區(qū)域有463個，內(nèi)含子區(qū)域有345個， 14個 SVs影響剪接位點。SVs的長度一般以100 bp為單位，長度大小從0～100 bp至 ＞1 200 bp。黃華占中58%的SVs和湘輻1821中57% 的SVs都是大于1 200 bp的序列（圖4）。

圖2 基因組InDels長度分布Fig. 2 The length distribution of the genome InDels

圖3 CDS區(qū)InDels長度分布Fig. 3 The length distribution of the CDS InDels

2.5 CNVs拷貝數(shù)變異分析

CNVs 可以在轉(zhuǎn)錄和翻譯水平上影響基因表達。利用CNVnator 軟件分別檢測到黃華占和湘輻1821中一共有16 614 和17 658個 CNVs。這些CNVs幾乎都分布在基因間隔區(qū)，黃華占11 430個，湘輻1821 12 203個。外顯子區(qū)域中黃華占2 090個，湘輻1821 2 238個?；蛏嫌? kb區(qū)域中黃華占1 092個，湘輻1821 1 118個。下游1 kb區(qū)域中黃華占826個，湘輻1821 866個。內(nèi)含子區(qū)域中黃華占有552個，湘輻1821有 614個（圖5）。黃華占中拷貝數(shù)增加個數(shù)是1 758個，拷貝數(shù)減少個數(shù)是14 856個，拷貝數(shù)總增加長度達到8 275 400 bp，拷貝數(shù)減少長度達到44 603 000 bp。湘輻1821中拷貝數(shù)增加個數(shù)是1 827個，拷貝數(shù)減少個數(shù)是15 831個，拷貝數(shù)總增加長度達到8 469 400 bp，拷貝數(shù)減少長度達到445 569 500 bp。2種基因型均以拷貝數(shù)減少為主要變異方式。

圖4 SVs長度分布圖Fig. 4 The length distribution of the SVs

圖5 CNVs位置分布圖Fig. 5 The CNVs distribution

2.6 變異的差異比較分析

SNPs和InDels變異可以直接地對應(yīng)到基因的變異，利于突變體性狀分析。通過對黃華占和湘輻1821中的SNPs和InDels進行差異分類統(tǒng)計，共找出差異SNPs 86 163個。其在12條染色體均有分布，其中第1染色體最多，達到11 389個，其次是第11染色體和12染色體。從分布區(qū)間來看，基因間隔區(qū)最多，有65 553個，外顯子區(qū)域只有3 375個（表4、表5）。對這些差異SNPs類型進行分類統(tǒng)計，湘輻1821中轉(zhuǎn)換類型的SNPs為38 121個，占總差異SNPs的44%；顛換類型的SNPs為13 100個，占總差異SNPs的15%。黃華占中轉(zhuǎn)換類型的SNPs為35 666個，占總差異SNPs的41%；顛換類型的SNPs為12 470個，占總差異SNPs的14%。差異InDels共88 777個，其在12條染色體上也均有分布，其中第1染色體最多，達到13 261個，其次是第5染色體和第2染色體。從分布區(qū)間來看，基因間隔區(qū)最多，有48 045個，外顯子只有2 843個（表4、表5）。

SNPs突變位點為單一堿基，對所涉及到的基因進行注釋，篩選候選基因可靠方便。通過對這些差異SNPs所涉及基因進行篩選，篩選得到的候選基因有3 092個，并對其進行GO分析。3 092個候選基因分為生物學(xué)過程（biological process）、分子功能（molecular function）和細胞成分（cellular component）3大類。如圖6所示：在生物學(xué)過程的17個小類中，大分子生物合成（macromolecule biosynthetic process）、細胞大分子生物合成（cellular macromolecule biosynthetic process）、大分子代謝調(diào)控（regulation of macromolecule metabolic process）、基因表達（gene expression）和DNA復(fù)制（DNA replication）類型基因較多，約占整個GO分析的70%；在分子功能的8個小類中，轉(zhuǎn)移酶活性（transferase activity）、鐵離子結(jié)合（iron ion binding）和DNA解旋酶活性（DNA helicase activity）的基因占主要類型；細胞成分包含5個小類，轉(zhuǎn)錄因子和RNA聚合酶相關(guān)基因較多。對這些差異候選基因的GO分析對湘輻1821表型分析可能具有重要參考意義。

表4 差異SNPs和InDels在染色體上的分布Tab. 4 The distribution of differential SNPs and InDels on chromosome

表5 差異SNPs和InDels在基因組上不同區(qū)間的分布Tab. 5 The distribution of differential SNPs and InDels on different genomic regions

3 討論

圖6 差異SNPs所在基因的GO分類Fig. 6 Gene ontology categories of differential SNPs associated gene

二代測序技術(shù)的高速發(fā)展推動了水稻基因組學(xué)的研究，高通量、高精度、低成本的超高優(yōu)勢使其迅速應(yīng)用于結(jié)構(gòu)基因組學(xué)、功能基因組學(xué)以及轉(zhuǎn)錄組學(xué)研究［25］。水稻基因組重測序是指在已知水稻基因組序列的前提下，對某個水稻品種基因組進行測序，對照參考基因組序列分析重測序品種變異情況，可挖掘大量的SNPs、InDels和SVs等，從而闡明該水稻品種的遺傳特征。結(jié)合二代測序技術(shù)，對輻射誘變機理的探索也廣泛開展起來，空間誘變獲得的突變體主要發(fā)生DNA水平上的突變［26］。重離子輻射改變植物體內(nèi)基因序列及表達方式，使突變體DNA甲基化、轉(zhuǎn)座子/反轉(zhuǎn)座子活性發(fā)生改變［27］。Cheng等［28］對從9311輻照來的突變體Red-1進行重測序，發(fā)現(xiàn)Red-1中9.19%的基因組序列發(fā)生改變，有381 403個SNPs，50 116個1～5 bp的InDels，在Red-1基因組中點突變變異方式是主要變異，這與本研究中的湘輻1821的變異方式一致，表明SNPs和InDels變異是γ射線輻射誘變的主要特征。

水稻育種已經(jīng)從高產(chǎn)育種目標(biāo)轉(zhuǎn)變?yōu)榧骖檭?yōu)質(zhì)、多抗、高氮利用率等綠色性狀。核誘發(fā)突變技術(shù)能夠誘發(fā)產(chǎn)生自然界稀有的新基因，大量的實踐證明利用人工誘發(fā)遺傳變異是豐富水稻種質(zhì)資源和選育新品種的重要手段之一［29］。輻射誘變可以改變一個或兩個主要性狀，并保持其他性狀不變。如經(jīng)γ輻照后，水稻和馬鈴薯中的淀粉含量發(fā)生改變［30-31］。本研究利用γ射線輻射處理黃華占干種子，獲得了綜合性狀良好、遺傳穩(wěn)定的突變新品系，較野生型其產(chǎn)量顯著提高，劍葉增寬，株高增加，且米質(zhì)經(jīng)農(nóng)業(yè)農(nóng)村部稻米檢測中心檢測達到經(jīng)農(nóng)業(yè)農(nóng)村部稻米檢測中心頒布的一等食用稻標(biāo)準(zhǔn)。通過二代高通量測序，在該突變體中發(fā)現(xiàn)了大量分子變異，湘輻1821中變異類型的個數(shù)均高于野生型黃華占。湘輻1821的SNPs轉(zhuǎn)換與顛換的比率為2.546，黃華占的為2.542，這2種基因型的SNPs變異基本以轉(zhuǎn)換為主要類型。在黃華占和湘輻1821的InDels中基本以短序列InDels為主，說明點突變?yōu)橹饕蛔兎绞健？截悢?shù)減少是黃華占和湘輻1821的主要拷貝數(shù)變異方式。通過比較分析，2種基因型的差異SNPs 共有86 163個，差異InDels共88 777個。差異SNPs所涉及基因篩選到的候選基因有3 092個。這些研究結(jié)果將為湘輻1821表型分析提供重要參考，同時為輻射誘變機理提供理論支撐。