基于溫帶和熱帶玉米群體全基因組選擇和雜種優(yōu)勢候選位點的鑒定

周玲　熊威　胡俏強　戴惠學(xué)　趙涵

摘要：基于187份種質(zhì)資源材料的全基因組重測序數(shù)據(jù)開發(fā)了120583個高質(zhì)量SNP變異位點，通過這些位點可以將具有不同遺傳背景的187份玉米種質(zhì)劃分為兩大類群，分別為包含100份材料的溫帶亞群和包含87份材料的熱帶亞群。通過對溫帶和熱帶玉米群體選擇信號的遺傳分化分析，檢測到3664個受到選擇的位點。選取187份材料中已報道能在溫/熱亞群形成雜種優(yōu)勢的135份代表性自交系，其中75份來自溫帶玉米自交系，60份為熱帶玉米系，基于兩大雜種優(yōu)勢群進行雜種優(yōu)勢性狀的全基因組關(guān)聯(lián)分析（GWAS），[JP2]結(jié)果鑒定出2407個雜種優(yōu)勢候選位點不均勻分布在玉米10條染色體上。整合選擇信號檢測和GWAS分析結(jié)果，共識別出1153個受到選擇的雜種優(yōu)勢相關(guān)位點，其中，619個位點與26個已報道的雜種優(yōu)勢相關(guān)QTLs一致。功能注釋發(fā)現(xiàn)與候選位點緊密連鎖的324個候選基因大部分都具有功能，其中包含61個重要的轉(zhuǎn)錄因子。根據(jù)GO富集分析發(fā)現(xiàn)這些候選基因主要參與了很多對雜種優(yōu)勢形成有貢獻的關(guān)鍵生化代謝途徑，包括氮化合物代謝、葉酸代謝、糖酵解、發(fā)育過程的負調(diào)控及轉(zhuǎn)錄調(diào)控等重要生物學(xué)途徑。

關(guān)鍵詞：溫熱帶玉米;全基因組選擇;雜種優(yōu)勢;關(guān)聯(lián)分析;候選位點

中圖分類號：S513.032文獻標志碼：A

文章編號：1002-1302（2021）04-0019-07

作者簡介：周玲（1987—），女，安徽合肥人，博士，副研究員，從事玉米遺傳育種學(xué)研究。E-mail：zlingxiaoyao@163.com。

通信作者：趙涵，博士，研究員，從事玉米遺傳育種學(xué)研究。E-mail：zhaohan@jaas.ac.cn。

玉米在我國是重要的糧食、飼料、經(jīng)濟作物，盡管目前培育和審定的玉米品種很多，但是雜交種親本趨同化嚴重，所用親本自交系主要集中在少數(shù)幾個骨干自交系或其衍生系，狹窄的玉米育種遺傳資源已成為制約當今玉米育種的瓶頸因素，且生產(chǎn)上又存在病蟲害大面積暴發(fā)的風險，擴增、改良和創(chuàng)新玉米種質(zhì)資源已經(jīng)迫在眉睫[1]。

雜種優(yōu)勢利用是玉米育種和生產(chǎn)的根本途徑，因此在種質(zhì)創(chuàng)新中應(yīng)優(yōu)先考慮創(chuàng)新種質(zhì)的雜種優(yōu)勢利用問題。根據(jù)生長環(huán)境的不同，玉米可劃分為溫帶和熱帶兩大類群。在溫帶最主要的雜優(yōu)利用模式是瑞德×蘭卡斯特，在熱帶、亞熱帶地區(qū)則是Tuxpeno×ETO。溫帶種質(zhì)主要表現(xiàn)出莖稈強韌、出籽率高、產(chǎn)量潛力大、農(nóng)藝性狀優(yōu)良等特點。熱帶、亞熱帶玉米種質(zhì)則表現(xiàn)出抗逆性強、根系發(fā)達、葉片濃綠、持綠期長等優(yōu)點[2]。由于熱帶和亞熱帶玉米種質(zhì)資源在病蟲、干旱等生物和非生物脅迫耐受性方面具有廣泛的遺傳變異，[JP2]對于玉米群體改良具有重要應(yīng)用價值，特別是對于溫帶玉米。因此利用熱帶、亞熱帶×溫帶雜優(yōu)模式，將有利于拓寬我國溫帶玉米種質(zhì)基礎(chǔ)，在我國玉米育種中發(fā)揮重要作用。

由于熱帶玉米在經(jīng)歷人工馴化和改良過程中逐漸適應(yīng)溫帶環(huán)境，因此，有必要利用選擇信號檢測法揭示在改良過程中發(fā)生選擇的位點，挖掘候選馴化基因，進一步從基因?qū)用嫔咸骄坑衩椎娜后w改良。例如Liu等利用溫帶、熱帶和亞熱帶玉米中代表性的260個玉米自交系，結(jié)合固定系數(shù)（FST）的計算方法，發(fā)現(xiàn)熱帶玉米相比于溫帶玉米具有更高的遺傳多樣性和更多的等位基因位點[3]。楊宇昕等通過對溫熱帶玉米群體的高測序深度的SNP進行分析，在熱帶玉米群體和溫帶玉米群體中鑒定到204752個SNP標記，通過FST選擇信號識別到557候選馴化基因，并鑒定到多個候選基因與玉米的開花調(diào)控密切相關(guān)[4]。[JP2]因此，利用選擇信號法探究溫熱帶玉米群體的基因組變化，可以鑒定在熱帶玉米適應(yīng)性改良過程中受到選擇的基因區(qū)段，并且深入挖掘更多與雜種優(yōu)勢關(guān)聯(lián)的基因，對于解析雜種優(yōu)勢群的遺傳機制和指導(dǎo)雜種優(yōu)勢群利用具有重要意義。

本研究利用來源不同的187份玉米自交系的全基因組重測序信息開發(fā)高質(zhì)量的SNP位點，通過這些遺傳變異位點研究自交系群體的遺傳結(jié)構(gòu)，劃分雜種優(yōu)勢群，通過選擇信號分析和關(guān)聯(lián)分析策略主要研究以下內(nèi)容：（1）鑒定溫帶和熱帶玉米在適應(yīng)環(huán)境和遺傳改良過程中，基因組上由于選擇遺留的選擇信號。（2）篩選出溫熱帶兩大雜種優(yōu)勢群之間的遺傳變異分化關(guān)聯(lián)的候選位點。（3）確定受到選擇的雜種優(yōu)勢關(guān)聯(lián)候選位點及其緊密連鎖的候選基因，并分析潛在選擇候選基因顯著富集途徑，為挖掘和利用控制雜種優(yōu)勢性狀的位點和候選基因提供信息參考。

1材料與方法

1.1研究材料與基因組測序數(shù)據(jù)

供試材料為187份玉米自交系，種質(zhì)基因組重測序數(shù)據(jù)下載于NCBI網(wǎng)站（www.ncbi.nlm.nih.gov）的SRA數(shù)據(jù)庫（PRJNA389800），對其進行過濾后利用BWA軟件比對到玉米B73（v3版本）基因組上并使用SAMtool對結(jié)果進行整理，最后利用GATK軟件進行SNP變異檢測，獲得變異位點后利用Plink軟件剔除缺失率大于10%和最小等位基因頻率低于0.05的SNP位點以及覆蓋率較低的樣本，采用Beagle軟件對缺失基因型進行填補后，最終獲得基因組的變異信息，SNP標記密度分布利用R軟件包CMplot繪制。

1.2雜種優(yōu)勢類群劃分

利用TASSELV5.0（www.maizegenetics.net/tassel）的鄰接算法（neighbor-joiningmethod，簡稱NJ）[5]計算187份玉米自交系之間的遺傳距離構(gòu)建遺傳樹，并進行雜種優(yōu)勢類群劃分。

1.3FST選擇信號檢測

基于不同群體的遺傳多樣性數(shù)據(jù)，可以利用FST進行溫熱帶玉米間選擇信號的檢測，其計算原理是依據(jù)染色體等位基因頻率變化。在中性進化條件下，F(xiàn)統(tǒng)計量的大小主要取決于遺傳漂變和遷移等因素的影響，若群體中一個等位基因因為對于特定環(huán)境的適合度較高而經(jīng)歷適應(yīng)性選擇，那么其頻率的升高會增大種群分化水平，反映在F統(tǒng)計量上就是有較高的FST值（0≤FST≤1，F(xiàn)ST<0.05表示種群間沒有遺傳分化，F(xiàn)ST>0.25則表示種群間分化程度非常高）[6]。為了降低因單位點SNP掃描受到遺傳漂變等因素影響而造成的假陽性，筆者利用VCFtools軟件計算滑動窗口10kb內(nèi)群體間的FST值來增加選擇信號的靈敏度[7]，并運用R包CMplot繪制全基因組水平上的FST結(jié)果圖。為了鑒定FST值的受選擇位點，選擇FST值高于0.25作為顯著閾值線，高于閾值線的SNP位點定義為“受選擇位點”。

1.4全基因組關(guān)聯(lián)分析（GWAS）

運用R軟件GAPIT軟件包中的壓縮混合線性模型（compressedmixedlinearmodel，簡稱CMLM）程序[8]，計算前5個主成分控制群體結(jié)構(gòu)，并根據(jù)VanRaden描述的方法[9]計算親緣關(guān)系K矩陣。在考慮研究材料群體結(jié)構(gòu)和親緣關(guān)系的情況下，對溫熱帶兩大類群進行全基因組關(guān)聯(lián)分析，溫帶的雜種優(yōu)勢表型值記作1，熱帶記作0。采用Bonferroni方法矯正P值，鑒定溫熱帶群體之間的特異性位點。

1.5功能注釋和GO富集分析

對基于FST和GWAS共識別的特異位點進行全基因組掃描，利用SnpEffv4.3軟件對變異信息進行功能注釋[10]，獲得與其緊密連鎖的候選基因，經(jīng)MaizeGDB（http：//www.maizegdb.org/blast.php）[JP2]數(shù)據(jù)庫比對尋找有同源序列的功能基因。為了更深刻解析選擇清除區(qū)域候選基因的生物學(xué)功能，利用在線平臺AgriGO的單一富集分析（SEA）功能（http：//systemsbiology.cau.edu.cn/agriGOv2/index.php）[11]對候選基因進行生物學(xué)過程（biologicalprocess）富集分析，篩選P值<005且錯誤發(fā)現(xiàn)率（FDR）<005的GO富集術(shù)語作為選擇候選基因。

2結(jié)果與分析

2.1基因型鑒定和進化樹構(gòu)建

本研究利用187份種質(zhì)資源材料的高通量重測序數(shù)據(jù)開發(fā)了120583個SNP位點，這些變異位點基本均勻分布于10條染色體上（圖1）。采用TASSEL5.0軟件的NJ法構(gòu)建了系統(tǒng)發(fā)育樹，結(jié)果（圖2）顯示，187份玉米材料被劃分為2個不同的類群上：100份材料被聚類為溫帶玉米亞群，另外87份材料被聚類為熱帶玉米亞群。

2.2基于全基因組選擇信號的檢測和關(guān)聯(lián)分析

基于全基因組檢測到的SNP標記利用VCFtools軟件計算溫熱帶亞群體之間的FST值，結(jié)果表明，全基因組水平上FST值高于閾值線0.25的受選擇位[CM（21][KG*8]點共有[KG*8]3[KG*3]664[KG*8]個（圖3-a），占總變異位點數(shù)的3%。其中染色體4受到選擇的顯著性位點最多（913個），染色體6（188個）和染色體7（189個）最少。第4染色體46380000～46390000區(qū)間內(nèi)含有最高的FST值，其值為0.66。

選取187份材料中具有廣泛代表性的135份玉米自交系（表1），其中75份是選自育種程序的溫帶玉米自交系，代表經(jīng)歷顯著改良并適應(yīng)溫帶長日照環(huán)境的溫帶玉米（如B73、Mo17等），而60份主要來源于國際玉米小麥改良中心（CIMMYT）的熱帶、亞熱帶玉米自交系（如CML103、NC350等），代表具有豐富的遺傳多樣性、改良程度較低、更適應(yīng)玉米馴化原始環(huán)境的熱帶玉米。在控制群體結(jié)構(gòu)和材料間親緣關(guān)系條件下，基于溫熱帶兩大雜種優(yōu)勢群體利用GAPIT軟件CMLM方法對其雜種優(yōu)勢性狀的GWAS，采用Bonferroni方法矯正P值，顯著水平閾值為α=0.05/120583=4.15×10-7，在-lgP>683的水平下，共檢測到2407個特異位點與雜種優(yōu)勢群劃分顯著關(guān)聯(lián)（圖3），其中顯著的P值的范圍是6.84～20.67，解釋了4.31%～11.97%的表型變異。在10條染色體的分布依次為214、180、568、576、41、51、97、274、122、284個，其中在染色體3（568個）和染色體4（576個）上面存在顯著富集。

整合全基因組選擇信號的檢測和GWAS分析的結(jié)果，發(fā)現(xiàn)受到選擇的雜種優(yōu)勢相關(guān)SNPs共有1153個，其中在染色體3（233個）、染色體4（433個）和染色體10（201個）上有顯著富集，在染色體1、2、5、6、7、8、9染色體上分別包含34、37、21、17、44、114、19個顯著關(guān)聯(lián)SNPs。與已報道的雜種優(yōu)勢相關(guān)QTL進行比較分析[12-15]，發(fā)現(xiàn)619個位點存在于24個已報道的QTL區(qū)段（表2），其中染色體3中4個QTL區(qū)段覆蓋了140個候選位點，染色體4中3個QTL區(qū)段覆蓋了146個位點，[JP2]染色體10中3個QTL區(qū)段覆蓋了178個位點，表明本研究識別的候選位點與已報道的雜種優(yōu)勢QTL結(jié)果高度一致。

2.3候選基因的富集分析和功能注釋

利用Perl腳本和SnpEff4.3軟件對受到選擇的雜種優(yōu)勢關(guān)聯(lián)的顯著SNP標記與玉米B73參考基因組進行比對，并進行候選基因預(yù)測及功能注釋，最終筆者獲得了324個候選基因，其中包含61個轉(zhuǎn)錄調(diào)控功能的基因（表3）。利用AgriGO在線軟件SEA對候選基因進行富集分析發(fā)現(xiàn)候選基因主要富集在55個GO富集術(shù)語生物過程（圖4）。這些GO富集術(shù)語主要涉及到氮化合物代謝（GO：0006807）、葉酸代謝（GO：0046655）、糖酵解（GO：0006096）、發(fā)育過程的負調(diào)控（GO：0051093）和轉(zhuǎn)錄調(diào)控（GO：0006355）等重要生物學(xué)途徑（圖4）。

3結(jié)論與討論

對玉米種質(zhì)遺傳多樣性深入認識是合理利用種質(zhì)的前提，隨著分子生物學(xué)的發(fā)展，分子標記已成為作物遺傳育種研究的重要手段，同時為評價玉米的遺傳變異提供了方便、快捷的研究方法。利用分子標記的方法對玉米種質(zhì)進行類群劃分，方法易行、結(jié)果可靠，已經(jīng)大規(guī)模用于玉米育種研究中。本研究以B73的參考基因組和187份玉米自交系的基因組信息為背景開發(fā)了120583個高質(zhì)量的SNP分子標記均勻分布于10條染色體上，基于這些標記將187份自交系劃分為2個雜種優(yōu)勢類群，分別為溫帶和熱帶兩大類群。合理的雜種優(yōu)勢群，是把遺傳關(guān)系較近的自交系劃分為同一雜種優(yōu)勢群，在雜交種選育過程中避免遺傳關(guān)系較近的自交系之間試組配工作，可以有效減少自交系之間的組配工作量，大大提高育種效率。

解析不同玉米類群遺傳差異，鑒定與雜種優(yōu)勢相關(guān)的基因有助于深入了解玉米分化過程，同時為玉米種質(zhì)的擴增和改良提供重要遺傳資源。在本研究中，筆者結(jié)合選擇信號分析和GWAS方法共鑒定出1153個受到選擇的雜種優(yōu)勢相關(guān)候選位點，前人已報道的26個雜種優(yōu)勢相關(guān)QTL區(qū)段覆蓋了本研究獲得的619個候選位點，表明本研究的結(jié)果與已報道的雜種優(yōu)勢QTL結(jié)果高度一致，這為溫熱帶種質(zhì)資源的雜交和改良奠定了基礎(chǔ)。功能注釋發(fā)現(xiàn)與候選位點緊密連鎖的324個候選基因中大部分都是具有功能的基因，利用GO富集分析發(fā)現(xiàn)，這些候選基因主要參與了很多對雜種優(yōu)勢形成有貢獻的關(guān)鍵生化代謝途徑，包括氮化合物代謝、葉酸代謝、糖酵解、發(fā)育過程的負調(diào)控及轉(zhuǎn)錄調(diào)控等重要生物學(xué)途徑。

此外，這些候選基因中還包括了61個轉(zhuǎn)錄因子，在馴化過程中，植物對各種環(huán)境、組織和發(fā)育信號做出反應(yīng)，需要各種功能基因的表達，從而進行精確的調(diào)控。而植物感受外界環(huán)境如干旱、病害、鹽堿等體內(nèi)細胞發(fā)育等信號時，需要一系列信號傳遞激發(fā)轉(zhuǎn)錄因子，從而啟動功能基因的轉(zhuǎn)錄表達，最后通過基因產(chǎn)物的作用對外界信號在生理生化等方面的變化做出適當?shù)恼{(diào)節(jié)反應(yīng)[16]?？梢?，植物的轉(zhuǎn)錄因子通過其功能域與DNA及其他蛋白質(zhì)間的相互作用，可以調(diào)節(jié)植物的生長發(fā)育及環(huán)境脅迫應(yīng)答等過程[17-18]。因此，這61個轉(zhuǎn)錄因子將成為研究玉米馴化的候選基因。

參考文獻：

[1]戴景瑞，鄂立柱.我國玉米育種科技創(chuàng)新問題的幾點思考[J].玉米科學(xué)，2010，18（1）：1-5.

[2]TallurySP，GoodmanMM.Experimentalevaluationofthepotentialoftropicalgermplasmfortemperatemaizeimprovement[J].TheoreticalandAppliedGenetics，1999，98（1）：54-61.

[3]LiuK，GoodmanM，MuseS，etal.GeneticstructureanddiversityamongmaizeinbredlinesasinferredfromDNAmicrosatellites[J].Genetics，2003，165（4）：2117-2128.

[4]楊宇昕，鄒棖.基于溫帶和熱帶玉米群體全基因組FST和XP-EHH的選擇信號檢測[J].中國農(nóng)業(yè)科學(xué)，2019，52（4）：579-590.

[5]SaitouN，NeiM.Theneighbor-joiningmethod：anewmethodforreconstructingphylogenetictrees[J].MolecularBiologyandEvolution，1987，4（4）：406-425.

[6]WrightS.Evolutionandgeneticsofpopulations[M].Chicago：UniversityofChicagoPress，1978.

[7]MaY，DingX，QanbariS，etal.Propertiesofdifferentselectionsignaturestatisticsandanewstrategyforcombiningthem[J].Heredity，2015，115（5）：426-436.

[8]ZhangZ，ErsozE，LaiCQ，etal.Mixedlinearmodelapproachadaptedforgenome-wideassociationstudies[J].NatureGenetics，2010，42（4）：355-360.

[9]VanRadenPM.Efficientmethodstocomputegenomicpredictions[J].JournalofDairyScience，2008，91（11）：4414-4423.

[10]CingolaniP，PlattsA，WangLL，etal.Aprogramforannotatingand[HT][HT][LM]

[KH*4D]

[HT8.][KG2]predictingtheeffectsofsinglenucleotidepolymorphisms，SnpEff：SNPsinthegenomeofDrosophilamelanogasterstrainw1118;iso-2;iso-3[J].Fly，2012，6（2）：80-92.

[11]TianT，LiuY，YanH，etal.AgriGOv2.0：aGOanalysistoolkitfortheagriculturalcommunity，2017update[J].NucleicAcidsResearch，2017，45（1）：122-129.

[12]StuberCW，LincolnSE，WolffDW，etal.Identificationofgeneticfactorscontributingtoheterosisinahybridfromtwoelitemaizeinbredlinesusingmolecularmarkers[J].Genetics，1992，132（3）：823-839.

[13]LuH，Romero-SeversonJ，BernardoR.Geneticbasisofheterosisexploredbysimplesequencerepeatmarkersinarandom-matedmaizepopulation[J].TheoreticalandAppliedGenetics，2003，107（3）：494-502.

[14]FrascaroliE，CanèMA，LandiP，etal.Classicalgeneticandquantitativetraitlocianalysesofheterosisinamaizehybridbetweentwoeliteinbredlines[J].Genetics，2007，176（1）：625-644.

[15][JP2]ThiemannA，F(xiàn)uJ，SeifertF，etal.Genome-widemeta-analysisofmaizeheterosisrevealsthepotentialroleofadditivegeneexpressionatpericentromericloci[J].BMCPlantBiology，2014，14：88.

[16]LedentV，VervoortM.Thebasichelix-loop-helixproteinfamily：comparativegenomicsandphylogeneticanalysis[J].GenomeResearch，2001，11（5）：754-770.

[17]VaragonaMJ，SchmidtRJ，RaikhelNV.Nuclearlocalizationsignal（s）requiredfornucleartargetingofthemaizeregulatoryproteinOpaque-2[J].ThePlantCell，1992，4（10）：1213-1227.

[18]劉欣，李云.轉(zhuǎn)錄因子與植物抗逆性研究進展[J].中國農(nóng)學(xué)通報，2006，22（4）：61-65.