陸地棉抗黃萎病性狀與SSR標記的關(guān)聯(lián)分析

王龍，李黎貝，馬啟峰，耿洪偉，陳全家，曲延英，范術(shù)麗

(1.新疆農(nóng)業(yè)大學(xué)農(nóng)學(xué)院，烏魯木齊 830052；2.棉花生物學(xué)國家重點實驗室/中國農(nóng)業(yè)科學(xué)院棉花研究所，河南安陽 455000)

0 引 言

【研究意義】棉花黃萎病是以大麗輪枝菌浸染為主的維管束病害，對棉花的產(chǎn)量和纖維品質(zhì)等農(nóng)藝性狀有很大影響。1914年在美國發(fā)現(xiàn)棉花黃萎病，隨后在澳大利亞等主要產(chǎn)棉國也陸續(xù)發(fā)現(xiàn)了棉花黃萎病[1]。1935年，棉花黃萎病由美國傳入中國，20世紀90年代，棉花黃萎病在我國主要產(chǎn)棉區(qū)逐年蔓延，其中西北內(nèi)陸植棉區(qū)、黃河流域植棉區(qū)、長江流域植棉區(qū)等地的范圍不斷擴大。棉花黃萎病成為阻礙棉花實現(xiàn)優(yōu)質(zhì)、高產(chǎn)、穩(wěn)產(chǎn)等育種目標的主要因素，已成為提高棉花的產(chǎn)量、纖維品質(zhì)的主要障礙[2]。【前人研究進展】有研究認為，棉花黃萎病抗性由單基因控制的質(zhì)量性狀[3-5]。多數(shù)學(xué)者認為，棉花黃萎病抗性由多個微效基因共同控制的數(shù)量性狀[6-9]。產(chǎn)生兩種結(jié)論的原因，可能是研究人員所選不同供試材料而造成的，并且棉花抗黃萎病性狀受環(huán)境影響較大。棉花黃萎病抗性QTL(Quantitative trait locus)定位在國外的相關(guān)研究較少，國內(nèi)相關(guān)研究較多；SSR分子標記已成為研究黃萎病的主要技術(shù)，研究群體由先前的陸陸群體轉(zhuǎn)變?yōu)楹ｊ懭后w，因為后者的多態(tài)性較好。作圖群體也由初級的F2、BC1群體向更適合精細定位的RIL(Recombinant inbred lines)、NIL(Near isogenic lines)等群體過度。高玉千等[10]利用海陸組合的F2群體、通過RAPD(Random amplified polymorphic DNA)和SSR(Simple sequence repeat)標記找到2個主效QTLs和一個微效QTL；王紅梅等[11]利用F2/F2∶3群體定位到6個與黃萎病抗性相關(guān)的QTLs；昝偉等[12]研究結(jié)果表明：在F2群體中，通過SSR分子標記檢測到與棉花黃萎病抗性相關(guān)的2個主效QTLs，并且定位到第一連鎖群上；李磊[13]研究結(jié)果表明：在構(gòu)建的兩個F2群體上進行定位，找到6個與黃萎病抗性相關(guān)的QTLs，其中包括3個抗病QTLs，3個感病QTLs。目前，關(guān)聯(lián)分析已經(jīng)成為一項技術(shù)體系比較完善、相對比較成熟的研究功能基因定位的一種方法。關(guān)聯(lián)分析要比常規(guī)QTL定位技術(shù)優(yōu)勢更明顯，(1)以自然群體為材料，無需構(gòu)建專門的作圖群體，而QTL作圖至少需要花費2年以上的時間來完成群體的構(gòu)建；(2)可實現(xiàn)對其作圖群體(自然群體)一個基因座上所有等位基因的考察，而絕大部分QTL作圖所利用的群體是兩親本雜交、重組自交后代，每一基因座一般只能涉及兩個等位基因；(3)定位精度高，可以達到單基因水平。關(guān)聯(lián)分析的方法在玉米、大麥、大豆、小麥等多個作物的一些功能基因定位研究中得到廣泛應(yīng)用[14-24]。近年來，通過關(guān)聯(lián)分析定位與棉花農(nóng)藝性狀相關(guān)功能基因的報道已有一些。Abdurakhmonov[18-19]利用285份陸地棉材料構(gòu)建自然群體，并結(jié)合使用202個SSR標記，研究了該群體纖維品質(zhì)性狀的關(guān)聯(lián)分析，研究結(jié)果表明：與纖維長度顯著關(guān)聯(lián)的標記有23個，與纖維比強度顯著關(guān)聯(lián)的標記有19個。賀道華等[20]用92份栽培棉構(gòu)建自然群體，對其進行掃描，并結(jié)合均勻分布于26條染色體上的132個 SSR標記，獲得21個與纖維品質(zhì)性狀相關(guān)的SSR位點。梁冰等[21]構(gòu)建68份陸地棉材料為自然群體，以160對SSR引物對該群體進行基因型檢測，發(fā)現(xiàn)了與棉花纖維品質(zhì)等主要農(nóng)藝性狀相關(guān)的SSR位點。Kantartzi等[22]以56份亞洲棉品種構(gòu)建自然群體，通過98個SSR標記對該研究群體的纖維品質(zhì)性狀進行關(guān)聯(lián)分析，研究結(jié)果表明：與纖維長度顯著關(guān)聯(lián)的標記3個，與纖維比強度顯著關(guān)聯(lián)的標記有2個。Wang等[23]以55個海島棉品系構(gòu)建自然群體，使用258個SRAP標記、170個SSR標記與該研究群體的產(chǎn)量、纖維品質(zhì)等性狀進行關(guān)聯(lián)分析， 研究結(jié)果表明，與產(chǎn)量性狀相關(guān)聯(lián)的標記有26個，與纖維品質(zhì)性狀相關(guān)聯(lián)的標記有46個?！颈狙芯壳腥朦c】近年來，關(guān)于研究棉花黃萎病抗性方面的關(guān)聯(lián)分析報道較少。Zhao等[24]構(gòu)建自然群體以158個陸地棉材料為研究對象，利用與棉花抗黃萎病性狀相關(guān)的212個SSR標記進行關(guān)聯(lián)分析，研究結(jié)果表明，與黃萎病抗性相關(guān)聯(lián)的標記位點有42個，其中新發(fā)現(xiàn)的有32個標記位點。郭志軍等[25]構(gòu)建自然群體以178份材料為研究對象，利用74對SSR引物與該群體進行關(guān)聯(lián)分析，檢測到位于NAU1514- NAU2754區(qū)間內(nèi)與棉花抗黃萎病性狀相關(guān)的1個QTL?！緮M解決的關(guān)鍵問題】以186份陸地棉材料組成的自然群體為研究對象，對其抗黃萎病性狀進行4個播期的表型鑒定，利用篩選出的140個SSR標記進行基因分型，分析群體結(jié)構(gòu)和品種間親緣關(guān)系，利用TASSEL4.0軟件的廣義線性模型(General line model，GLM)(P<0.001)和混合線性模型(Mixed linear model，MLM)(P<0.01)軟件，對棉花農(nóng)藝性狀和標記進行關(guān)聯(lián)分析，發(fā)掘棉花農(nóng)藝性狀顯著關(guān)聯(lián)的SSR標記，為棉花抗黃萎病分子標記輔助選擇育種提供依據(jù)。

1 材料與方法

1.1 材 料

試驗以186份陸地棉材料構(gòu)建自然群體。表1

表1 186份陸地棉材料來源
Table 1 The origin regions of 186 upland cotton accessions

中國農(nóng)業(yè)科學(xué)院棉花所試驗基地位于安陽市白璧鎮(zhèn)，將供試材料種植在試驗基地，試驗設(shè)計為2年(2015和2016年)、2個不同播期(春播和夏播)一共4個環(huán)境的種植模式。186份供試材料播種期分別標記為 M(2015年4月26日)、 W(2015年5月23日) 、 Y(2016年5月1日)和 Z(2016年5月25日)，一共4個環(huán)境。采用隨機區(qū)組排列，設(shè)置3次重復(fù)，行長6 m，行距0.7 m，株距0.3 m，常規(guī)大田農(nóng)事管理。

在棉花黃萎病發(fā)病高峰期進行黃萎病表型病級鑒定，分別為2015年8月1日、2016年8月2日進行，試驗地棉花黃萎病的表型病級鑒定采用5級制的劃分方法[26]。表2

表2 棉花葉片黃萎病級劃分標準
Table 2 The classification standard of cotton leaf verticillium wilt

1.2 方 法

1.2.1 SSR標記

以棉花幼苗八葉期葉片為試驗材料，提取DNA參考宋國立等[27]改進的CTAB法。選取覆蓋陸地棉26條染色體的140個SSR標記位點，對186份供試材料進行基因型檢測，在該群體中檢測到355個等位變異。所用SSR引物由上海生工生物工程股份有限公司合成。

PCR反應(yīng)體系為10 μL： 包括模板DNA 30 ng，10×PCR Buffer 1.0 μL，10 mmol/L dNTP 0.2 μL，10 mmol/L SSR引物2 μL，5 U/μLTaqDNA聚合酶0.1μL，加ddH2O補足至10 μL。試驗所用的Taq聚合酶和dNTP均由北京全式金生物技術(shù)有限公司提供。PCR反應(yīng)程序：95℃預(yù)變性5 min；30個循環(huán)(94℃變性30 s，58℃退火溫度45 s，延伸45 s)；72℃延伸10 min；10℃保存。使用濃度為8的聚丙烯酰胺凝膠電泳對PCR擴增產(chǎn)物進行分離，待跑完電泳以后，銀染方法依據(jù)張軍等[28]的方法進行并照相保存。

1.3 數(shù)據(jù)處理

帶型統(tǒng)計采用0、1統(tǒng)計法，在同一遷移率位置上，有帶記為1，無帶記為0，缺失記為‘？’。根據(jù)DNA Marker結(jié)合Quantity One軟件的方法，以此來判斷多態(tài)性條帶的大小。

利用Power Maker V3.25軟件計算引物的多態(tài)性信息含量(Polymorphism information Content，PIC)[29]。

參照Henderson[30]提出的BLUP(Best Linear Unbiased Prediction)法，通過R語言中的“l(fā)me4”包進行四個環(huán)境下的育種值計算。

對186份陸地棉材料構(gòu)建的自然群體基于數(shù)學(xué)模型進行的類群劃分，通過使用 STRUCTURE2.3.4分析軟件[31]，進一步對該群體的群體結(jié)構(gòu)進行估算，并計算得出的186份研究群體中每個個體所相應(yīng)的Q值。使用該軟件程序設(shè)定K為1至10，將MCMC(Markov Chain Monte Carlo)開始時的不作數(shù)迭代(Length of burn- in period)設(shè)為10 000次，再將不作數(shù)迭代后的MCMC設(shè)為100 000次，參照文自翔等[32]的設(shè)置對程序進行處理，并進行10次重復(fù)計算[24]。通過SAS運算出最大似然值，根據(jù)運行結(jié)果確定合適的K值，若似然值LnPD隨K值的增大而持續(xù)增大，則可參照Evanno等[33]的方法根據(jù)ΔK確定合適的K值。

通過使用TASSEL4.0軟件模型，利用廣義線性模型(General linear model，Q模型)是以個體Q值作為協(xié)變量，而混合線性模型(Mixed liner model，Q+K模型)是以Q值親緣關(guān)系(Kinship)作為協(xié)變量，運用兩種模型對4個環(huán)境下抗黃萎病表型數(shù)據(jù)與SSR位點進行回歸分析。利用SPAGeDi-1.3[34]軟件對親緣關(guān)系進行計算。通過運用Bresghellohe和Sorrells等[35]提出的“無效等位變異”方法，分析了棉花抗黃萎病表型效應(yīng)值情況。

2 結(jié)果與分析

2.1 田間鑒定

供試材料種植在中國農(nóng)業(yè)科學(xué)院棉花所試驗基地，在棉花黃萎病發(fā)病高峰期進行黃萎病表型鑒定，分別為2015年8月1日、2016年8月2日進行，黃萎病的表型病級鑒定采用馬存[26]提出的5級制劃分方法。圖1

圖1 田間棉花葉片病級劃分標準
Fig.1 The standard phenotype of cotton leaf in classifying verticillium wilt resistance in the field

2.2 陸地棉材料的群體結(jié)構(gòu)

研究186份供材料的群體結(jié)構(gòu)，通過使用Structure version2.3.4軟件對其進行分析，結(jié)果顯示：LnP(D)值隨K值的增大而增大，發(fā)現(xiàn)沒有一個拐點，因此，無法確定該研究群體結(jié)構(gòu)的K值；隨后采用通過ΔK法來確定群體結(jié)構(gòu)的K值[33]。當K=2時，平均似然值的斜率改變，將186份研究材料分成2個亞群。一個亞群有96份材料，占總數(shù)的51.60% ；另一亞群有90份材料，占總數(shù)的48.40%。

利用Structure軟件，計算出當K=2時，可將186份研究材料劃分為2個亞群，得到對應(yīng)的Q值，為關(guān)聯(lián)分析協(xié)變量的分析提供依據(jù)。圖2

圖2 186份陸地棉材料的群體結(jié)構(gòu)
Fig.2 The population structure 186 cultivated upland cotton

2.3 單個環(huán)境下關(guān)聯(lián)

通過使用140對SSR引物做分子標記，在該群體中檢測到355個等位變異，平均每對引物的等位變異為2.54，其中平均值為0.76，引物多態(tài)性信息含量變化范圍：0.50～0.99；將186份供試材料兩年、兩個播期共四個環(huán)境的棉花黃萎病指數(shù)，結(jié)合140個SSR標記，通過GLM(P<0.001)模型和MLM(P<0.01)模型同時進行關(guān)聯(lián)分析，在檢測到的355個位點中，同時能在兩個模型下檢測到且與黃萎病抗性相關(guān)的位點有22個。有3個標記位點與M播期相關(guān)聯(lián)，有6個標記位點與W播期相關(guān)聯(lián)，有3個標記位點與Y播期相關(guān)聯(lián)，有11個標記位點與Z播期相關(guān)聯(lián)。其中，NAU998和CGR5258兩個位點不僅與Y播期、Z播期的黃萎病病指相關(guān)聯(lián)，還分別于M播期、W播期的黃萎病病指相關(guān)聯(lián)。NAU998和CGR5258兩個位點表型變異解釋率分別為5.53%和12.07%。表3

表3 陸地棉抗黃萎病性狀關(guān)聯(lián)標記的貢獻率及P值

2.4 多點環(huán)境下關(guān)聯(lián)

將兩年、兩個播期下共四個環(huán)境的表型數(shù)據(jù)構(gòu)建BLUP(Best linear unbiased prediction)模型，同時進行該模型與單個環(huán)境下的GLM模型以及MLM模型下的關(guān)聯(lián)分析。研究表明，在GLM模型下關(guān)聯(lián)到14個位點，在MLM模型下關(guān)聯(lián)到3個位點，最終發(fā)現(xiàn)兩個位點同時存在于2個環(huán)境以上，分別為NAU998和CGR5258，該結(jié)果與單個環(huán)境下關(guān)聯(lián)分析的結(jié)果相吻合，說明這兩個標記可能在棉花黃萎病抗性的分子鑒定方面有較好的效果，對棉花分子標記輔助育種有一定的應(yīng)用價值。圖3

圖3 關(guān)聯(lián)分析流程
Fig.3 The flow chart of association analysis

2.5 抗黃萎病性狀表型效應(yīng)值

與棉花抗黃萎病性狀相關(guān)的兩個位點，分別為NAU998和CGR5258。依據(jù)文自翔[32]提出的“無效等位變異”判斷方法，得出SSR位點對表型效應(yīng)的影響。位點NAU998在四個環(huán)境下的表型效應(yīng)分別為0.03/0.08/-0.005/0.033，增效效應(yīng)更為顯著，說明這個位點的存在有利于提高棉花黃萎病的抗病性；位點CGR5258在四個環(huán)境下的表型效應(yīng)都表現(xiàn)為-0.001/-0.001/-0.001/-0.001，均為減效效應(yīng)，這個位點的存在會使棉花黃萎病的感病性提高。在育種過程中，可選擇具有NAU998位點同時缺失CGR5258位點的中間材料，以提高選育材料的黃萎病抗性。表4

表4 抗黃萎病性狀的表型效應(yīng)
Table 4 Phenotypic effects of resistance of Verticillum wilt

3 討 論

目前，有關(guān)棉花關(guān)聯(lián)分析的研究大多數(shù)于集中于主要農(nóng)藝性狀等方面，然而關(guān)于棉花抗黃萎病性狀的報道仍然較少[21,25]。陸地棉黃萎病抗性的遺傳方式為數(shù)量性狀，由多個微效基因共同控制，且致病菌有變化快的特點，而且棉田的發(fā)病易受濕度和溫度等條件的影響。研究利用兩年、兩個播期的四環(huán)境數(shù)據(jù)，發(fā)掘與陸地棉抗黃萎病性狀相關(guān)聯(lián)的SSR位點，為培育棉花抗黃萎病新品種提供可靠的分子標記。在檢測到與抗黃萎病性狀相關(guān)聯(lián)的22個位點中，兩個環(huán)境中同時被關(guān)聯(lián)到的位點NAU998和CGR5258，并且在前人研究結(jié)果中尚未出現(xiàn)。因此，可通過構(gòu)建分離群體對這2個位點進行驗證，并應(yīng)用于改良棉花抗黃萎病育種實踐。

研究在進行關(guān)聯(lián)分析過程中運用了GLM和MLM兩種模型，GLM通過群體結(jié)構(gòu)Q值作協(xié)變量；而使用MLM模型在進行關(guān)聯(lián)分析時，不但使用了Q值作協(xié)變量，還引用了親緣關(guān)系K；因此，在一定程度上能降低關(guān)聯(lián)結(jié)果的假陽性，使結(jié)果更為準確，無論是單個環(huán)境下的關(guān)聯(lián)結(jié)果還是多個環(huán)境下的關(guān)聯(lián)結(jié)果，MLM模型所關(guān)聯(lián)到的位點始終少于GLM模型。研究在進行關(guān)聯(lián)分析時使用了Henderson提出了BLUP法，即最佳線性無偏預(yù)測法，這種方法最初運用于動物育種當中，而20世紀80年代中期以后，BLUP方法開始應(yīng)用于牧草、農(nóng)作物及園林植物的遺傳改良研究，在關(guān)聯(lián)分析研究中用來構(gòu)建模型處理多年多點的數(shù)據(jù)，可進一步消除多年多點對數(shù)據(jù)關(guān)聯(lián)結(jié)果的影響，且已在大豆，玉米等[36-38]的關(guān)聯(lián)分析研究中得到應(yīng)用。

結(jié)合前人對棉花黃萎病抗性QTL定位的研究結(jié)果，利用e-PCR程序[39]，與前人研究結(jié)果相比，NAU998和CGR5258不能比對到Zhang等[40]構(gòu)建的連鎖群中。但是在最近的研究中，分別與抗病性QTLs位點qDI-C3-3的近距離標記MUSS59相近，與QTLs位點qDI-C15-3的近距離標記NAU3337相近[41]，且分別位于NAU3775與BNL0663和NAU3684與BNL3971之間[42]。出現(xiàn)這種差異的原因可能是由于定位的材料不同而導(dǎo)致的，NAU998和CGR5258標記的發(fā)現(xiàn)，將對棉花抗黃萎病基因的克隆和功能驗證奠定堅實的基礎(chǔ)[43]。

4 結(jié) 論

將四個環(huán)境下的表型數(shù)據(jù)進行BLUP之后進行關(guān)聯(lián)分析，無論是GLM模型還是MLM模型所關(guān)聯(lián)到的位點都少于單個環(huán)境下所關(guān)聯(lián)到的位點，BLUP的精準性，研究得到的NAU998和CGR5258兩個位點具有較好的可靠性。

近年來，有關(guān)棉花關(guān)聯(lián)分析主要以SSR標記居多，傳統(tǒng)的研究提供基因組多態(tài)性信息相對較少，未來的發(fā)展更趨向于通過基于簡化基因組測序技術(shù)和重測序技術(shù)進行全基因組關(guān)聯(lián)分析[44-46]。研究通過使用140對SSR做分子標記，檢測到355個等位變異，進行關(guān)聯(lián)分析，標記數(shù)量需要進一步加大，研究結(jié)果只是粗略意義上的全基因組關(guān)聯(lián)分析， 因此，還需增加標記的密度，可以結(jié)合SNP標記，開展更加詳實的棉花抗黃萎病GWAS研究；第二，要對已發(fā)掘的SSR位點進一步功能驗證，通過構(gòu)建分離群體進行標記與性狀的相關(guān)性分析，可進行QTL分析，以此來提高相應(yīng)位點的準確性。通過對棉花抗黃萎病性狀的關(guān)聯(lián)分析，發(fā)掘了一些與黃萎病抗性相關(guān)聯(lián)的SSR位點，為陸地棉抗黃萎病材料進行全基因組掃描提供了更多的分子標記。