馬鈴薯資源晚疫病抗性的全基因組關聯分析

蔣 偉 潘哲超 包麗仙 周福仙 李燕山 隋啟君,* 李先平,*

1 云南省農業科學院經濟作物研究所, 云南昆明 650205; 2 云貴高原馬鈴薯與油菜科學觀測實驗站, 云南昆明 650200

晚疫病是馬鈴薯的主要病害, 每年可造成我國馬鈴薯減產約10%～15%, 經濟損失高達20億美元。西南地區低溫高濕的氣候環境, 晚疫病發生更嚴重,馬鈴薯減產可達 15%～40%[1-2]。提高品種晚疫病抗性是馬鈴薯育種重要目標之一, 而挖掘持久、高抗的晚疫病抗性基因對馬鈴薯抗病育種具有重要的意義。

20世紀50年代, 育種家將野生種Solanum demissum的抗性基因導入到栽培種中, 但抗性很快被新出現的晚疫病生理小種克服[3-4]。目前已經鑒定并克隆到一些對晚疫病表現出廣譜抗性的R基因。通過構建BAC文庫和長片段PCR擴增的方法, 從野生種S. bulbocastanum中克隆了第1個馬鈴薯晚疫病廣譜抗性基因RB或Rpi-blb1[5-6]。通過構建S. bulbocastanum作圖群體, 利用AFLP分子標記, 克隆到另1個廣譜抗性基因Rpi-blb2[4]。利用野生種S. venturii的 5個作圖群體, 結合分離群體分組分析法(bulked segregant analysis, BSA)找到與抗性基因連鎖的AFLP分子標記, 通過構建BAC文庫從而克隆了廣譜抗性基因Rpi-vnt1[7]。對含有廣譜抗性基因R8分子標記或QTL位點的群體構建BAC文庫, 通過酶切PCR克隆到了該基因[8-9]。上述這些基因的鑒定主要是通過構建遺傳群體、分子標記定位獲得, 而通過自然群體全基因組關聯分析定位晚疫病抗性基因的研究較少。

本研究選擇的自然群體引自國際馬鈴薯中心,該中心從19世紀80年代起以淘汰馬鈴薯晚疫病垂直抗性而保留其水平抗性為目的, 篩選了一批群體材料。因此, 該群體中很可能保留著晚疫病廣譜抗性基因。通過對四倍體馬鈴薯群體連續2年的田間晚疫病抗性評價, 結合群體簡化基因組測序獲得大量的SNP位點, 進行全基因組關聯分析挖掘晚疫病抗性相關的遺傳位點和候選基因, 旨在為晚疫病抗性品種選育和抗病機理研究提供一定的理論和材料基礎。

1 材料與方法

1.1 材料

引自國際馬鈴薯中心的288份四倍體馬鈴薯晚疫病抗性群體材料(表1), 感病對照為Desiree (D214)。群體A中晚疫病抗性來源為改良的S. demissum衍生抗性, 包括來自安第斯栽培種S. tuberosumgroupsandigena、phureja和stenotomum, 以及野生種S.acaule和S. bulbocastanum抗性的材料 8份; 群體B3來源于群體A, 但篩選保留表現出晚疫病水平抗性的資源材料 82份; 群體 B1來自S. tuberosumgroupsandigena的材料15份; 群體LTVR具有馬鈴薯病毒病(PVY、PVX和PLRV)抗性、生長期短和適應溫暖環境特性的材料111份; 群體B3-HT結合群體B3的晚疫病抗性、北美和歐洲品種的耐熱性以及群體LTVR特性的材料34份; 群體B3-LTVR包括群體B3和群體LTVR的雜交后代22份; 群體BW具有細菌性枯萎病抗性的材料 12份; 群體PREBRED具有從野生種導入到四倍體B3或LTVR中的晚疫病抗性材料1份; VARIETY表示馬鈴薯品種 3份[10]。

1.2 田間試驗和表型數據分析

2015—2016年, 群體材料種植于云南省曲靖市會澤縣試驗基地(26o06'N, 103o22'E), 海拔 2650 m。每個材料種植1行10株, 行距為0.7 m, 株距為0.3 m, 種植3次重復。采用Alpha試驗設計(CIP提供軟件), 整個生長季進行正常中耕、施肥管理等, 但不噴灑殺菌類農藥。

調查田間晚疫病發病比率。發現初始病株, 每7 d調查1次, 共調查8次, 計算AUDPC和sAUDPC值。利用Microsoft Excel對所觀測到的數據進行統計分析。

式中,n表示總的調查次數,xi表示第i次調查的晚疫病嚴重程度,ti表示第i次調查的時間。

1.3 DNA提取和測序方法

利用 CTAB法提取 DNA, 檢測合格后利用Illumina HiSeq2500測序儀(Illumina, Inc; San Diego,CA, USA)測序。以馬鈴薯S. chacoenseM6的基因組測序信息為參考[11], 利用北京百邁客生物有限公司自主研發的簡化基因組測序流程 SLAF-seq[12]進行電子酶切預測, 根據最佳酶切方案選擇RsaI+HaeIII進行酶切, 酶切長度在314～414的序列被定義為SLAF標簽, 獲得多態性標簽1,032,778個。以馬鈴薯S. chacoenseM6為參考基因組, 利用BWA軟件[13]將 reads比對到參考基因組上, 利用 GATK[14]和SAMtools[15]軟件篩選SNP。根據完整度大于0.8, 次要基因頻率大于0.05過濾后, 得到353,753個SNP位點進行后續的分析。

1.4 晚疫病抗性全基因組關聯分析

基于SNP位點, 利用TASSEL軟件[16]的混合線性模型(compressed mixed linear model, CMLM), 通過公式y= Xα+Qβ+Kμ+e得到關聯值。公式中, 通過Admixture軟件[17]計算樣品群體結構 Q矩陣,SPAGeDi軟件[18]計算樣品間親緣關系K矩陣, X為基因型矩陣, y為表型, 最終獲得每個SNP位點的關聯值, 并計算它們的P值。以 AUDPC_2015、AUDPC_2016、sAUDPC_2015、sAUDPC_2016, 以及這 2年的 AUDPC和 sAUDPC平均值即AUDPC_Mean和 AUDPC_Mean為表型數據, 共分析6個性狀。對這些性狀利用GLM、MLM、CMLM、FASTLMM和EMMAX五種模型進行關聯分析。利用 GGplot2軟件[19]繪制 Quantile-Quantile散點圖(Q-Qplot), 利用 QQman軟件繪制曼哈頓圖(Manhattan), 展示關聯分析顯著的SNP位點。

1.5 候選基因預測

對6個性狀在5種模型下獲得的顯著水平是0.1和0.01的關聯SNP位點進行注釋。選取這些位點前后100 kb范圍內的基因, 利用BLAST軟件分別與NR、SwissProt、GO、COG和KEGG數據庫進行比對, 獲得關聯區域內基因注釋信息, 并根據注釋信息尋找與晚疫病抗性相關的候選基因。

表1 自國際馬鈴薯中心引進的288份馬鈴薯群體材料Table 1 288 potato population resources from CIP

(續表 1)

(續表 1)

2 結果與分析

2.1 晚疫病抗性的表型分析

通過對288份馬鈴薯資源材料連續2年的晚疫病田間抗性評價發現, 晚疫病抗性評價值 AUDPC和 sAUDPC呈連續分布。2015年所有材料的AUDPC平均值為2244.1, 范圍是121.3～4740.8, 變異系數為 61.73%, sAUDPC平均值為 3.9, 范圍是0.2～8.3, 變異系數為 61.56%; 2016年所有材料的AUDPC平均值為1228.5, 范圍是72.7～3371.3, 變異系數為80.47%, sAUDPC平均值為3.6, 范圍是0.2～9.5,變異系數為80.21%。表明該群體表現出晚疫病抗性廣泛的遺傳差異, 且受多基因控制, 屬于數量性狀遺傳。相關性分析中, 2015—2016年AUDPC值相關系數為0.8988, 高度線性相關; sAUDPC值相關系數為0.8974,高度線性相關(圖1)。表明2年間群體材料田間晚疫病抗性表現一致, 抗性穩定。

2.2 馬鈴薯晚疫病抗性的全基因組關聯分析

5種模型下對6個性狀分析的QQ圖結果表明,在所有模型下關聯到顯著相關的SNP位點都是可靠的(圖2是2個性狀AUDPC_mean和sAUDPC_mean在 5種模型下的 QQ圖)。在 5種分析模型下, 共關聯到與6種晚疫病抗性性狀顯著相關(顯著水平為 0.1和 0.01)的 SNP位點 82個(表 2和圖 3)。其中, 6個SNP位點分布在3號染色體(圖3-A, B); 3個 SNP位點分布在 4號染色體(圖 3-A, B); 1個SNP位點分布在5號染色體(圖3-B, C); 分別有2個SNP位點分布在6號、7號染色體(圖3-A, B, D);66個 SNP位點分布在 9號染色體(圖 3-A～D); 各有1個SNP位點分布在11號和12號染色體(圖3-A,B, D)。在82個SNP位點中, 20個SNP位點在所有模型分析中, 2種顯著水平情況下, 與所有性狀都相關。它們是31,375,673 (第9個)、31,707,832(第 23 個)、32,382,888 (第 25 個)、32,973,114 (第35 個)、33,118,463 (第 36 個)、33,145,927 (第 38個)、33,148,740 (第 39 個)、33,156,262 (第 40 個)、33,211,574 (第 43 個)、33,244,698 (第 48 個)、33,283,731 (第 49 個)、33,283,949 (第 50 個)、33,451,096 (第 53 個)、33,480,937 (第 55 個)、34,959,324 (第 62 個)、35,060,347 (第 66 個)、35,065,836 (第 67 個)、35,123,523 (第 72 個)、35,212,534 (第 77 個)、35,256,410 (第 80 個)。同一性狀下, GLM模型分析定位到的SNP位點數目最多,而 MLM 模型分析定位到的 SNP位點數目最少。在GLM 模型下分析, 顯著水平為 0.1時定位到與AUDPC_Mean相關的SNP位點數目(74個)最多(表2)。

2.3 候選基因預測

本研究在82個關聯SNP位點上下游區域100 kb,共檢測到922個候選基因。其中, 460個候選基因在所有分析(60種)中都存在, 其余基因在2～56種分析中出現。根據基因組注釋, 篩選出54個已知或可能與晚疫病抗性相關的基因(表3)。其中, 23個基因為抗性基因包括晚疫病抗性基因R1同源基因、Sw-5同源基因(R8)和Rpi-vnt1以及編碼多效性耐藥蛋白基因; 5個基因編碼MAPK蛋白和WRKY轉錄因子;1個基因參與茉莉酸途徑; 3個基因與水楊酸途徑相關; 6個基因是病程相關的基因; 3個基因參與苯基丙酸類合成途徑; 其他與晚疫病抗性相關的基因如HMGR基因(2個)、細胞色素P450 (21個)。

表2 6個馬鈴薯晚疫病抗性性狀在不同分析模型和顯著水平下的顯著關聯標記Table 2 Significant correlation markers of six phenotype of potato resistance to late blight under different analytical models and significant levels

(續表 2)

(續表 2)

3 討論

3.1 不同模型分析對全基因關聯分析結果的影響

本研究采用了5種模型對基因組數據和表型數據進行關聯分析。GLM 模型(general linear model)采用Q矩陣[16], MLM模型(mixed linear model)采用Q+K矩陣[16], CMLM模型(compressed mixed linear model)算法介于GLM和MLM之間, 尋找最優算法壓縮計算時間[20], EMMAX模型(efficient mixed model association eXpedited)在MLM模型基礎上校正關聯群體的群體結構和遺傳親緣關系并加速這種運算過程[21-22], FASTLMM (factored spectrally transformed linear mixed model)同樣修正了MLM模型的運算公式, 提高運算速度[23]。盡管不同的模型存在差異, 但是在所有模型下分析的表型性狀得到的結果都比較可靠(圖2)。5種模型下共關聯到82個SNP位點, 而同一性狀、同樣的顯著水平下, GLM模型分析關聯到的SNP位點數量明顯多于其他模型(表2)。與晚疫病抗性相關的候選基因中, 關聯到48個 SNP位點, 其中 13個(第 6、14、16、31、32、56、60、65、68、69、70、75和81號) SNP位點只在GLM模型下關聯到, 16個SNP位點在所有模型到關聯到, 剩余的19個SNP位點在2～4個模型下獲得。因此, 利用多種模型分析能關聯到更多可能與晚疫病抗性相關的候選基因。

3.2 馬鈴薯晚疫病水平抗性群體中可能存在的抗性基因

栽培馬鈴薯晚疫病抗性主要來自安第斯栽培種S. tuberosumAndigenum group (包括andigena、phureja和stenotomun)[24-26], 并通過導入野生種如S.demissum[27]、S. venturii[7]、S. bulbocastanum[28]等的抗性基因不斷改良馬鈴薯栽培種的晚疫病抗性。盡管國際馬鈴薯中心通過不斷的輪回選擇創制了一批晚疫病水平抗性群體[29], 但是這批材料中可能仍存在垂直抗性基因(R1-R11)[30]。本研究所用的試驗材料即是來自國際馬鈴薯中心篩選的水平抗性群體,并且關聯到12個R1類似基因, 3個番茄斑萎病抗性基因, 4個Rpi-vnt1抗性基因以及4個多效性耐藥蛋白基因。

植物中抗性基因編碼的最大的一類蛋白具有核苷酸結合位點和亮氨酸重復結構(nucleotide-binding site plus leucine-rich repeat, NB-LRR)[31], 而目前已知的馬鈴薯晚疫病抗性基因編碼的都是NB-LRR蛋白[32]。R1基因定位在5號染色體上, 是植物抗性基因亮氨酸zipper/NBS/LRR中的一員[33]。根據前人的研究, 至少有4個抗性基因位于9號染色體, 分別是R8[8-9]、R9a[34]、Rpi-vnt1[7]和Rpi-moc1[35], 這 4 個基因均位于 9號染色體長臂末端, 而且遺傳距離比較近[34,36]。本研究所關聯到的12個R1類似基因位于9號染色體, 其長度為 35.34 Mb, 而關聯到的 SNP的位置在 34.8 Mb左右, 因此, 這些抗性基因肯定不是R1, 而可能是已克隆的R8或其他抗性基因。Sw-5是番茄中一類廣譜抗性的番茄斑萎病(ToMV)抗性基因[37], 位于 9號染色體上的抗性基因R8與Sw-5同源, 而Rpi-vnt1與番茄中另1個ToMV持久抗性基因 Tm-22同源[7,38]。本研究中另外關聯到與RPi-vnt1基因類似的基因, 其位置在32.5 Mb附近,而關聯到的ToMV抗性基因位于34.7 Mb左右, 根據R8與Rpi-vnt1的遺傳位置[36],Rpi-vnt1在R8的上方, 因此關聯到位于32.5 Mb的基因可能是Rpi-vnt1,而位于34.7 Mb的基因可能是R8。多效性耐藥蛋白是一類 ABC轉運蛋白, 可參與到擬南芥[39](PEN3)和小麥[40](Lr34)中非宿主抗性和對病原菌的持久抗性中。馬鈴薯中的多效性耐藥蛋白也可參與到晚疫病菌的侵染過程中, 但是具體怎么調節這個過程尚不清楚[41]。本研究關聯到的多效性耐藥蛋白基因位于 31.7 Mb左右, 其對馬鈴薯晚疫病抗性貢獻的大小, 有待進一步研究。

?

?

?

3.3 其他可能影響馬鈴薯晚疫病抗性的相關基因

本研究所用材料主要是表現出晚疫病水平抗性,即非R基因主導的持久抗性[29]。還存在其他途徑影響晚疫病抗性[42-43], 包括參與信號轉導和基因調控、茉莉酸途徑、水楊酸途徑、病程相關基因、苯基丙酸類合成途徑、HMGR基因和細胞色素P450。馬鈴薯中MAPK激酶和WRKY轉錄因子可參與植物對晚疫病菌的防御過程如StMPK1[44]、StWRKY8[45]、StWRKY1[46]。本研究在 7號染色體上關聯到 1個MAPK基因和WRKY基因, 在3號染色體上關聯到3個WRKY基因, 它們對晚疫病水平抗性的貢獻有待進一步研究。激素信號在植物病害防御中發揮著重要作用, 植物激素茉莉酸和水楊酸代謝途徑中的因子不同程度地參與到馬鈴薯晚疫病水平抗性中[47-49]。本研究中關聯到 1個與茉莉酸代謝途徑可能相關和3個與水楊酸代謝途徑可能相關的基因, 在其他馬鈴薯晚疫病抗性全基因組關聯分析和 QTL定位研究中也定位到與茉莉酸和水楊酸代謝途徑相關的基因[50-51]。有趣的是, 本研究中關聯到 3個參與苯基丙酸類合成途徑的查爾酮合成酶基因(CHS), 位于9號染色體。馬鈴薯中CHS基因至少有6個拷貝,研究較多的是位于5號染色體的CHS基因參與馬鈴薯塊莖中花青素合成[52-54]。但是, 查爾酮合成酶能夠參與植物病害防御反應[55-56], 表明位于 9號染色體的CHS基因可能是構成馬鈴薯晚疫病水平抗性因素之一。其他基因包括與病程相關的基因、HMGR基因和細胞色素 P450都有報道可以參與植物病害防御過程, 但是在馬鈴薯晚疫病水平抗性中所起的作用是不清楚的, 而本研究關聯到的基因可以為馬鈴薯晚疫病水平抗性的形成機理提供參考[42,57-59]。

4 結論

基于5種模型的GWAS結果, 共檢測到82個與晚疫病抗性顯著相關的SNP位點, 其中48個位點關聯到可能與晚疫病抗性相關的候選基因, 包括R1同源基因、Sw-5同源基因(R8)、Rpi-vnt1基因、與信號轉導相關的MAPK基因、WRKY基因、參與植物激素茉莉酸和水楊酸代謝相關的基因、參與苯基丙酸類合成途徑的查爾酮酶基因和其他病程相關的基因、HMGR基因以及細胞色素P450基因。