大豆抗旱遺傳位點及候選基因發(fā)掘

田蕊， 張華， 黃玫紅， 邵振啟， 李喜煥， 張彩英

（河北農(nóng)業(yè)大學(xué)農(nóng)學(xué)院，教育部華北作物種質(zhì)資源研究與利用重點實驗室，河北 保定 071001）

大豆是重要的糧食和油料作物，是植物蛋白和食用油的主要來源[1]。干旱嚴(yán)重影響大豆生長發(fā)育，干旱條件下大豆植株生長緩慢，葉片萎蔫，生物量降低，以至減產(chǎn)甚至絕收。葉片是大豆植株重要光合器官，其葉綠素含量與光合速率及產(chǎn)量密切相關(guān)，易受干旱影響，是反映品種抗旱性的重要指標(biāo)[2]。已有研究對作物干旱條件下的葉綠素含量進(jìn)行定位分析，沈波等[3]利用水稻重組自交系（recombinant inbred lines， RIL）群體，在水分脅迫和正常條件下對葉綠素含量進(jìn)行連鎖分析，定位到13 個控制葉綠素含量的數(shù)量性狀基因座（quantitative trait locus，QTL），其中水分脅迫條件下6 個，解釋表型變異4.74%～19.86%；胡頌平等[4]以水稻F9重組自交系（珍汕97B×IRAT109）為材料，在正常與水分脅迫下對葉綠素含量進(jìn)行定位，檢測到13 個葉綠素含量QTLs，其中正常條件下檢測到7個，聯(lián)合貢獻(xiàn)率為56.19%；干旱處理下檢測到6 個，聯(lián)合貢獻(xiàn)率為47.39%。Kumar 等[5]利用C306×HUW206 構(gòu)建的小麥RIL 群體在干旱條件下定位到1 個控制葉綠素含量QTL（QChl.ksu-3B），該QTL 解釋表型變異的14.2%。Peleg 等[6]以硬粒小麥和野生二粒小麥構(gòu)建的RIL 群體為材料，在干旱及正常水分條件下，對葉綠素含量進(jìn)行QTL定位，獲得8個葉綠素含量QTLs，解釋表型變異的2%～10%。

崔世友等[7]以‘波高’ב南農(nóng)94-156’構(gòu)建RIL群體為材料，在4 個不同生育期對大豆葉綠素含量進(jìn)行連鎖分析，定位到與葉綠素含量相關(guān)的10 個位點，分別位于1、7、12、13、18 和19 號染色體，解釋表型變異的6.9%～23.4%。Fang 等[8]利用809 份大豆自然群體重測序，對葉綠素含量進(jìn)行關(guān)聯(lián)分析，獲得9 個大豆葉綠素含量相關(guān)位點，分別位于11、13、18、20 號染色體。Kato 等[9]利用‘Minsoy’בT225H’構(gòu)建的F2:3家系在14 及17 號染色體定位到2 個葉綠素相關(guān)位點。由此可見，相較其他農(nóng)作物，目前大豆多就正常水分條件下葉綠素含量進(jìn)行研究。因此，本研究以199 份大豆品種資源構(gòu)建的自然群體為材料，在干旱和正常澆水條件下，測定其葉綠素含量，并計算抗旱指數(shù)，結(jié)合20×重測序基因型數(shù)據(jù)進(jìn)行關(guān)聯(lián)分析，發(fā)掘基于葉綠素含量的抗旱遺傳位點和相關(guān)基因，為大豆品種抗旱分子遺傳改良提供依據(jù)。

1 材料與方法

1.1 供試材料

以本課題組前期構(gòu)建的大豆自然群體為材料，包括152個育成品種和47個地方品種[10]。

1.2 大豆種植方法

將蛭石和土以體積比1∶1 進(jìn)行混合，裝入花盆（直徑35 cm，高40 cm，花盆下鋪塑料布），每盆裝土20 kg，施復(fù)合肥5 g，澆水6 L；當(dāng)基質(zhì)含水量為20%時播種，每盆播種20 粒；待2 片真葉完全展開時定苗，每盆保留10株。試驗分別于2021年6月（E1）、7月（E2）和8月（E3）在河北農(nóng)業(yè)大學(xué)作物育種中心試驗基地旱棚內(nèi)進(jìn)行，設(shè)置對照組和干旱脅迫處理組，2 次重復(fù)，其中干旱脅迫處理組于大豆V3 期開始進(jìn)行干旱處理，對照組正常澆水，期間利用托普云農(nóng)土壤水分含量測定儀對土壤水分含量進(jìn)行測定，待干旱處理組的土壤含水量低于10%（植株受到干旱脅迫）時測定葉綠素含量[11-12]。

1.3 葉綠素含量測定

利用葉綠素儀SPAD-502Plus 分別測定干旱處理組和對照組的葉綠素含量。每盆隨機(jī)選擇3 個單株，分別測定其最上部完全展開三出復(fù)葉的葉綠素含量，其中每個葉片選擇非葉脈處測定3次。

1.4 葉綠素含量抗旱指數(shù)計算

利用干旱處理及正常澆水條件下的葉綠素含量計算抗旱指數(shù)（drought resistance index，DRI），公式如下。

式中，CCD為干旱條件下葉綠素含量；CCW為正常條件下葉綠素含量；CCDM為干旱條件下供試大豆自然群體葉綠素含量平均值。

1.5 葉綠素抗旱指數(shù)遺傳變異分析方法

利用SPSS 25.0 對基于葉綠素含量計算的抗旱指數(shù)進(jìn)行描述統(tǒng)計及方差分析，利用公式計算其廣義遺傳力（broad-sense heritability，h2）[10]。

1.6 葉綠素抗旱指數(shù)關(guān)聯(lián)分析方法

利用GEMMA（genome-wide efficient mixed model association algorithm）軟件，采用混合線性模型（mixed linear model, MLM）對葉綠素抗旱指數(shù)進(jìn)行全基因組關(guān)聯(lián)分析（genome-wide association study, GWAS），其中供試大豆自然群體的基因型數(shù)據(jù)包含566萬個單核苷酸多態(tài)性（single nucleotide polymorphism，SNP）標(biāo)記[10]；此外，參照已有報道根據(jù)LD（linkage disequilibrium）衰減距離，同一LD區(qū)域顯著關(guān)聯(lián)SNP視為1個位點[10,13]。

1.7 候選基因的篩選

在顯著關(guān)聯(lián)SNP 上、下游100 kb 篩選候選基因，結(jié)合大豆公共數(shù)據(jù)庫Soybase（https://www.soybase.org/）注釋基因功能，初步篩選候選基因；然后利用已有大豆干旱脅迫條件下的轉(zhuǎn)錄組數(shù)據(jù)[14]，經(jīng)標(biāo)準(zhǔn)化處理分析各基因表達(dá)量，根據(jù)基因差異表達(dá)結(jié)果確定候選基因。

2 結(jié)果與分析

2.1 大豆自然群體葉綠素含量抗旱指數(shù)遺傳變異分析

對199 份供試大豆品種資源進(jìn)行葉綠素含量抗旱指數(shù)方差分析（表1）發(fā)現(xiàn)，每個環(huán)境條件下的供試品種資源抗旱指數(shù)均具有極顯著差異，說明干旱脅迫處理有效，品種資源抗旱能力存在較大差異。同時發(fā)現(xiàn)，不同環(huán)境下的葉綠素含量抗旱指數(shù)也存在差別，其平均值由大到小依次為E3、E2、E1，變化范圍為1.35～1.48。進(jìn)一步分析不同品種資源葉綠素抗旱指數(shù)遺傳變異（表1）發(fā)現(xiàn)，變異系數(shù)為8.71%～11.64%，以E2 環(huán)境條件下的變異系數(shù)最大。分析葉綠素含量抗旱指數(shù)廣義遺傳力發(fā)現(xiàn)，其值為81.44%，表明主要由遺傳因素決定。另外，分析葉綠素含量抗旱指數(shù)的偏度和峰度系數(shù)發(fā)現(xiàn)，絕對值均趨近于0，說明葉綠素抗旱指數(shù)屬多基因控制數(shù)量性狀（圖1）。

表1 大豆自然群體抗旱指數(shù)遺傳變異分析Table 1 Genetic variation analysis of drought resistance index in soybean natural population

2.2 大豆葉綠素含量抗旱指數(shù)GWAS分析

對3 種環(huán)境條件下的葉綠素抗旱指數(shù)及其平均值、最佳線性無偏預(yù)測（best linear unbiased prediction，BLUP）值進(jìn)行關(guān)聯(lián)分析（圖1）發(fā)現(xiàn)，關(guān)聯(lián)到237 個SNPs，分別位于1、2、3、7、10、16、17 及18 號染色體，其中，11 個SNPs 位于基因外顯子并引起非同義突變，28 個位于基因上游或下游，2 個位于基因3’UTR，這些SNPs位點可能影響基因功能或調(diào)控基因表達(dá)。此外，依據(jù)LD 衰減距離（100 kb）可將237 個SNPs 劃分為18 個位點，其中位點1、5、6、11、12、13、14、15、17、18 在E1、E2、E3、平均值及BLUP 值中均被關(guān)聯(lián)到，屬多環(huán)境關(guān)聯(lián)穩(wěn)定位點（表2）。

進(jìn)一步分析發(fā)現(xiàn)，在16號染色體關(guān)聯(lián)到181個SNPs（分布于Locus 10～Locus 14），占總SNPs 的76.37%，其中4 個位點在上述3 種環(huán)境以及平均值和BLUP 值同時關(guān)聯(lián)到，說明這4 個位點不僅是多環(huán)境穩(wěn)定位點，而且是控制葉綠素抗旱指數(shù)熱點區(qū)域；并且，上述181 個SNPs 中，49 個位于基因內(nèi)，其中9 個在內(nèi)含子，17 個在外顯子（6 個屬于非同義突變，9 個屬于同義突變，2 個屬于提前終止），11 個在基因上游，10 個在下游，2 個在3’UTR。此外，在10 號染色體關(guān)聯(lián)到27 個SNPs（位于Locus 5～Locus 9），其中Locus 5 和Locus 6屬多環(huán)境穩(wěn)定關(guān)聯(lián)位點，且在上述27 個SNPs 中，2 個位于基因上游或下游，2 個在內(nèi)含子區(qū)域。

2.3 葉綠素抗旱指數(shù)顯著關(guān)聯(lián)SNP驗證

為驗證顯著關(guān)聯(lián)SNPs及其位點準(zhǔn)確性，選取18 個位點的最顯著P值關(guān)聯(lián)SNPs，對參試品種含有的葉綠素抗旱指數(shù)優(yōu)異等位SNPs 數(shù)目及其抗旱指數(shù)進(jìn)行回歸分析，結(jié)果（圖2）表明，隨優(yōu)異等位變異數(shù)目增加，葉綠素抗旱指數(shù)逐漸增加（R2=0.76）。依據(jù)參試品種葉綠素抗旱指數(shù)，分別選取高、低極端品種各5個，對每個品種含有的18個位點中的11 個純合SNPs 數(shù)目進(jìn)行分析，結(jié)果發(fā)現(xiàn)（表3），5個抗旱品種中，3個具有所有的11個優(yōu)異等位變異，另2 個分別含有9 和10 個優(yōu)異等位變異，而5個敏感品種則不含任何優(yōu)異等位變異。

圖2 大豆葉綠素抗旱指數(shù)與優(yōu)異等位變異數(shù)量回歸分析Fig. 2 Regression analysis of drought resistance index and elite allele numbers of soybean accessions

表3 大豆抗旱特異種質(zhì)顯著關(guān)聯(lián)SNP分析Table 3 Analysis of associated SNPs in different soybean accessions with different drought resistance index

2.4 大豆抗旱相關(guān)基因發(fā)掘

在顯著關(guān)聯(lián)SNPs上、下游100 kb范圍內(nèi)尋找到57 個候選基因（表4），其中10 個基因外顯子存在等位變異，10 個基因內(nèi)含子存在等位變異，15 個基因上游或下游存在等位變異，2 個基因3’UTR存在等位變異?；诖耍鶕?jù)Soybase網(wǎng)站（https://www.soybase.org）Wm82 基因組注釋信息，從57 個候選基因中初步篩選出16 個相關(guān)基因，分別參與植物激素信號轉(zhuǎn)導(dǎo)、蛋白磷酸化、抗旱反應(yīng)等；利用已有大豆轉(zhuǎn)錄組數(shù)據(jù)[14]分析候選基因表達(dá)量（圖3），最終篩選出3 個抗旱相關(guān)基因，即Glyma.16G063600、Glyma.10G007000和Glyma.17G143900。

圖3 大豆抗旱候選基因表達(dá)分析Fig. 3 Expressions of candidate genes for soybean drought tolerance

表4 大豆抗旱指數(shù)候選基因Table 4 Candidate genes of drought resistance index in soybean

Glyma.16G063600存在1 個非同義突變（C/A，天冬氨酸/谷氨酸）、3 個同義突變、2 個下游SNPs、17 個基因間SNPs（圖4）。該基因編碼F-box 蛋白，可與小熱激蛋白結(jié)合，進(jìn)而調(diào)控植物耐旱性[15]。依據(jù)已發(fā)表轉(zhuǎn)錄組數(shù)據(jù)[14]發(fā)現(xiàn)，該基因受干旱脅迫誘導(dǎo)顯著上調(diào)表達(dá)。利用該基因的非同義突變分析供試群體發(fā)現(xiàn)，含有優(yōu)異等位變異A與含有非優(yōu)異等位變異C的品種間的葉綠素抗旱指數(shù)差異極顯著。

圖4 候選基因Glyma.16G063600分析Fig. 4 Identification of candidate gene Glyma.16G063600

Glyma.10G007000編碼AP2-EREBP 轉(zhuǎn)錄因子。AP2/EREBP 轉(zhuǎn)錄因子含有1 個或多個AP2/ERF 結(jié)構(gòu)域，該類轉(zhuǎn)錄因子可響應(yīng)植物鹽、旱、冷等多種非生物脅迫[16]。轉(zhuǎn)錄組數(shù)據(jù)分析發(fā)現(xiàn)，干旱脅迫條件下，Glyma.10G007000表達(dá)顯著上調(diào)；利用其上游SNP 分析群體發(fā)現(xiàn)，含有優(yōu)異等位變異的品種其葉綠素抗旱指數(shù)顯著高于含有非優(yōu)異等位變異的品種（圖5）。

圖5 候選基因Glyma.10G007000分析Fig. 5 Identification of candidate gene Glyma.10G007000

Glyma.17G143900編碼AP2/ERF 結(jié)構(gòu)域蛋白。利用基因間最顯著SNP 分析群體發(fā)現(xiàn)，含有優(yōu)異等位變異T的品種其葉綠素抗旱指數(shù)顯著高于含有非優(yōu)異等位變異C的品種（圖6）。

圖6 候選基因Glyma.17G143900分析Fig. 6 Identification of candidate gene Glyma.17G143900

2.5 大豆抗旱優(yōu)異基因型篩選

利用Glyma.10G007000、Glyma.16G063600及Glyma.17G143900等位變異，對供試品種資源進(jìn)行分析，結(jié)果（表5）發(fā)現(xiàn)，可將其分為8種類型，其中基因型Ⅶ（AAT）含3 個優(yōu)異等位變異，其葉綠素抗旱指數(shù)最高；基因型Ⅳ（CAT）、基因型Ⅵ（ACT）及基因型Ⅷ（AAC）含有2 個優(yōu)異等位變異，葉綠素抗旱指數(shù)次之；而不含任何優(yōu)異等位變異的基因型Ⅰ（CCC）葉綠素抗旱指數(shù)最低，說明3 個候選基因與干旱脅迫下的抗旱指數(shù)相關(guān)，為大豆抗旱育種提供了基礎(chǔ)材料。

表5 供試大豆自然群體基于入選抗旱基因等位變異分析Table 5 Analysis of soybean natural population based on the allele variations of three selected causal genes

3 討論

葉綠素是光能轉(zhuǎn)化為化學(xué)能的關(guān)鍵色素[17]。對干旱條件下的葉綠素含量進(jìn)行研究發(fā)現(xiàn)，隨干旱脅迫增加，葉綠素含量呈現(xiàn)先增加后降低的趨勢，說明葉綠素含量受干旱脅迫的影響[18-19]?？购抵笖?shù)是指性狀抗旱系數(shù)與干旱條件下該性狀指數(shù)的乘積[20]。劉桂茹等[21]對小麥抗旱系數(shù)、干旱敏感系數(shù)、單株產(chǎn)量降低指數(shù)分析發(fā)現(xiàn)，抗旱指數(shù)能反應(yīng)品種的抗旱性，可作為小麥抗旱篩選的重要指標(biāo)。楊玉敏等[22]以100 份小麥資源為材料，對抗旱系數(shù)、抗旱指數(shù)、耐旱指數(shù)進(jìn)行分析發(fā)現(xiàn)，抗旱指數(shù)綜合考慮了各品種在干旱條件下的敏感度和干旱脅迫下目標(biāo)性狀在整個群體中的地位，比較適合作為作物抗旱性評價的標(biāo)準(zhǔn)。基于此，本研究利用干旱處理和正常澆水處理下的葉綠素含量計算抗旱指數(shù)，用于大豆品種資源抗旱種質(zhì)篩選以及遺傳位點與候選基因挖掘；并依據(jù)群體供試品種資源的葉綠素抗旱指數(shù)，篩選出優(yōu)異種質(zhì)5 份，分別是‘冀10B5’‘滄0627’‘邯15-685’‘冀11B9’和‘HN1030’，其抗旱指數(shù)為1.70～1.83；進(jìn)一步分析發(fā)現(xiàn)，其含有9～11 個優(yōu)異等位變異，且均含有3 個入選相關(guān)基因的優(yōu)異等位變異，為大豆抗旱育種奠定了重要物質(zhì)基礎(chǔ)。

目前，有關(guān)水稻、小麥等作物的干旱脅迫葉綠素含量遺傳位點發(fā)掘研究已有報道，但大豆中開展相關(guān)研究較少[23-26]。本研究利用SPAD-502Plus測定干旱及正常條件下供試大豆自然群體葉綠素含量，計算抗旱指數(shù)，并依據(jù)本團(tuán)隊前期重測序數(shù)據(jù)對群體的葉綠素含量抗旱指數(shù)進(jìn)行關(guān)聯(lián)分析，結(jié)果獲得18 個位點，通過與已有研究進(jìn)行比較發(fā)現(xiàn)，Locus3 與李廣軍等[27]定位到的葉綠素含量QTLLeaflet chlorophyll 1-4重合，其加性效應(yīng)值為0.97，解釋表型變異5%。同時發(fā)現(xiàn)，本研究定位到的Locus15與Du等[28]定位到的干旱敏感指數(shù)相關(guān)QTLDrought index 1-9距離700 kb，該QTL 加性效應(yīng)值為0.128 1，可解釋表型變異6.42%；且QTLDrought index 1-5與本研究中的Locus10 相距4.2 Mb，該QTL 解釋表型變異9.36%。另外，本研究獲得的Locus3 與Carpentieri 等[29]定位到的耐旱相關(guān)QTLDrought tolerance 6-1及Shi 等[30]定位到的葉綠素相關(guān)QTLLeaflet chlorophyll 3-1相距2.2 Mb。上述結(jié)果證實了本研究結(jié)果的可靠性，為大豆抗旱遺傳改良提供了選擇標(biāo)記。

為篩選大豆抗旱相關(guān)基因，依據(jù)關(guān)聯(lián)區(qū)域內(nèi)的基因SNP突變、基因注釋及其表達(dá)量篩選出3個抗旱相關(guān)基因，其中候選基因Glyma.16G063600編碼F-box 蛋白，該類蛋白既可調(diào)控植物生長發(fā)育，也可對植物多種生物及非生物脅迫作出響應(yīng)[31-32]。Qu 等[33]發(fā)現(xiàn)，與擬南芥F-box 基因突變體fof2相比，過表達(dá)FOF2株系表現(xiàn)為ABA 含量增加及耐旱性提高；Lim 等[34]發(fā)現(xiàn)，辣椒F-box 類基因CaDIF1沉默植株表現(xiàn)出對干旱脅迫的敏感性，而過表達(dá)CaDIF1植株則表現(xiàn)為耐旱；Xu 等[15]發(fā)現(xiàn)，在干旱條件下，過表達(dá)大豆F-box 類基因GmFBL144可提高擬南芥過氧化氫及丙二醛含量。本研究篩選的F-box蛋白基因Glyma.16G063600與上述已報道的GmFBL144屬同一亞族，說明其可能在大豆抵御干旱脅迫逆境中具有重要作用。

本研究還篩選出AP2/EREBP（APETALA2/ethylene-responsive element binding proteins）轉(zhuǎn)錄因子家族基因Glyma.10G007000，該類轉(zhuǎn)錄因子通常包含1 或2 個AP2/ERF 結(jié)構(gòu)域，依據(jù)結(jié)構(gòu)域數(shù)目，可將其分為2 個亞族，即EREBP 亞族（包含1 個AP2/ERF 結(jié)構(gòu)域）和AP2 亞族（包含2 個AP2/ERF 結(jié)構(gòu)域）；AP2 亞族可調(diào)控花、胚珠及種子發(fā)育過程，EREBP 亞族（包括DREB 和ERF 類）參與植物激素調(diào)節(jié)，并響應(yīng)非生物脅迫[16]。研究發(fā)現(xiàn)，過表達(dá)蘋果AP2/EREBP 基因MdSHINE2可增加擬南芥植株莖和葉表面蠟質(zhì)，并提高其抗旱性[35]。大豆DREB 同源基因GmDREB2在干旱、高鹽及低溫等非生物脅迫條件下被誘導(dǎo)表達(dá)，且過表達(dá)該基因會激活下游轉(zhuǎn)錄因子進(jìn)而提高抗旱性[36]。因此，推斷AP2/EREBP 轉(zhuǎn)錄因子基因Glyma.10G007000可能參與大豆抗旱反應(yīng)。

另外，本研究還獲得了諸如MYB 轉(zhuǎn)錄因子、Dof （DNA binding with one finger）結(jié)構(gòu)域蛋白、CCR-4（carbon catabolite repressor 4）、LRR（Leucinerich repeat）候選基因，已有研究發(fā)現(xiàn)這些蛋白均可響應(yīng)植物干旱脅迫，Cominelli 等[37]發(fā)現(xiàn)，AtMYB60通過調(diào)控葉片保衛(wèi)細(xì)胞的開關(guān)影響植物抗旱；Ma 等[38]在干旱脅迫下對白菜9 個Dof 亞族基因表達(dá)模式進(jìn)行分析發(fā)現(xiàn)，不同亞族Dof基因均可響應(yīng)干旱脅迫；Li 等[39]發(fā)現(xiàn)，大豆LRR-RLK基因GmSARK可通過調(diào)控葉綠素降解基因Gmlls1及葉綠素合成基因Gmgtr1調(diào)控葉片衰老，說明該類基因位于Gmlls1和Gmgtr1上游，可間接調(diào)控葉綠素合成或降解。因此，本研究篩選出的這些候選基因可作為今后開展大豆抗旱分子育種研究的重要目標(biāo)基因。